1
$\begingroup$

Problem Setting

I have defined the following functions.

atomNumber[R_, d_, n_] := 2*n*volume[R, d]

volume[R_, d_ ] := Pi R^2 d

\[Gamma]Fn[x_] := 2/x^2 (1 - Exp[-x^2] (BesselI[0, x^2] + BesselI[1, x^2]))

\[Eta]Eval[\[Lambda]_, d_,  R_, rc_, n_] := 
\[Lambda] *((atomNumber[R, d, n])^2/d^2) \[Gamma]Fn[R/(Sqrt[2] rc)] (1 - Exp[-d^2/(4*rc^2)])

\[Del]Fn[\[Omega]_, m_] := Sqrt[Quantity["ReducedPlanckConstant"]/(\[Omega]*m)]

The inputs to all these functions should eventually be physical units. I then combine these functions towards the following:

coherence[\[Omega]_, d_, R_, rc_, \[Lambda]_, n_, T_, m_] := 
 4*\[Eta]Eval[\[Lambda], d,  R, rc, n]*T*(\[Del]Fn[\[Omega], m])^2

If I provide units and evaluate the functions separately, everything works fine, so there is no issue with any definition here.

Goal

I two variables: $r_c$ (rc) and $\lambda$ (\[Lambda]), all other inputs are parameters. I want to find all tuples $(r_c, \lambda)$ for which the following inequality holds:

coherence[\[Omega]_, d_, R_, rc_, \[Lambda]_, n_, T_, m_] <= 1

In particular, I want to use Manipulate to find these tuples for different parameter values, i.e. I want to have slides where I can change $\omega, d, R$ given the right physical units and then get a plot for the tuples.

My solution

I use RegionPlot for the inequality and Manipulate for the parameter slide. That leads to

Manipulate[
 RegionPlot[
  coherence[\[Omega]Paper Quantity[40, "Terahertz"], 
    d  Quantity[0.25, "Millimeters"], R Quantity[3.6, "Micrometers"], 
    rc, \[Lambda], n Quantity[176.2 * 10^(27), "Meters"^(-3)], 
    T Quantity[350, "Femtoseconds"], 
    m Quantity[6, "AtomicMassUnit"]] < 1, {rc , 
   Quantity[10^-9, "Meters"], Quantity[10^-1, "Meters"]}, {\[Lambda], 
   Quantity[10^-10, ("Seconds")^-1], 
   Quantity[10^-1, ("Seconds")^-1]}, 
  PlotRange -> All], {\[Omega]Paper, 20, 60}, {d, 0.01, 1}, {R, 0.01, 
  10}, {n, 100, 200}, {T, 300, 400}, {m, 1, 10}]

Here, I specified my variables $r_c$ and $\lambda$ and provided slides for all other parameters. However, I get an empty plot not showing anything.

Can anybody help me with this? I suspect that the problem has something to do with the Quantities package and my misuse thereof.

Additional Input

This is the list of hard-coded parameters I use:

\[Omega] = Quantity[40, "Terahertz"]
R = Quantity[3.6, "Micrometers"]
d = Quantity[0.25, "Millimeters"]
n = Quantity[176.2 * 10^(27), "Meters"^(-3)]
rc = Quantity[10^(-7), "Meters"]
\[Lambda] = Quantity[10^(-17), "Seconds"^(-1)]
T = Quantity[350, "Femtoseconds"]
m = Quantity[6, "AtomicMassUnit"]

This shows that coherence works:

testCoherence = coherence[\[Omega], d, R, rc, \[Lambda], n, T, m]
Out: 2.28009*10^-15

This is an example of RegionPlot without Manipulate, which also shows an empty plot:

RegionPlot[
 Evaluate[coherence[\[Omega]Paper, d, R, rc, \[Lambda], n, T, m]] < 
  1, {rc , Quantity[10^-9, "Meters"], 
  Quantity[10^-1, "Meters"]}, {\[Lambda], 
  Quantity[10^-10, ("Seconds")^-1], Quantity[10^-1, ("Seconds")^-1]}, 
 PlotRange -> Automatic]
$\endgroup$
3
  • $\begingroup$ Can you present an example of a working RegionPlot for your data with hard coded parameters of your choice as a starting point? $\endgroup$ Commented Mar 26, 2022 at 14:58
  • $\begingroup$ You have not provided the definition of the function volume $\endgroup$ Commented Mar 26, 2022 at 15:13
  • $\begingroup$ @BobHanlon Fixed both issues $\endgroup$ Commented Mar 26, 2022 at 15:18

1 Answer 1

4
$\begingroup$

You should also add some precision control.

Clear["Global`*"]

volume[R_, d_] := Pi R^2 d

atomNumber[R_, d_, n_] := 2*n*volume[R, d]

γFn[x_] := 2/x^2 (1 - Exp[-x^2] (BesselI[0, x^2] + BesselI[1, x^2]))

ηEval[λ_, d_, R_, rc_, n_] := 
 UnitConvert[λ*((atomNumber[R, d, n])^2/d^2) γFn[
    R/(Sqrt[2] rc)] (1 - Exp[-d^2/(4*rc^2)])]

∇Fn[ω_, m_] := 
 UnitConvert[Sqrt[Quantity["ReducedPlanckConstant"]/(ω*m)]]

coherence[ω_, d_, R_, rc_, λ_, n_, T_, m_] := 
 UnitConvert[4*ηEval[λ, d, R, rc, n]*T*(∇Fn[ω, m])^2]

Plotting,

Manipulate[
 RegionPlot[
  Evaluate[coherence[
     ωPaper Quantity[40, "Terahertz"],
     d Quantity[1/4, "Millimeters"],
     R Quantity[18/5, "Micrometers"],
     rc, λ,
     n Quantity[1762*10^(26), "Meters"^(-3)],
     T Quantity[350, "Femtoseconds"],
     m Quantity[6, "AtomicMassUnit"]] < 1],
  {rc, Quantity[10^-9, "Meters"], Quantity[10^-1, "Meters"]},
  {λ, Quantity[10^-10, ("Seconds")^-1], 
   Quantity[10^-1, ("Seconds")^-1]},
  PlotPoints -> 75,
  MaxRecursion -> 3,
  PlotRange -> All],
 {{ωPaper, 30}, 20, 60, 0.25, Appearance -> "Labeled"},
 {{d, 0.25}, 0.01, 1, 0.01, Appearance -> "Labeled"},
 {{R, 2}, 0.01, 10, 0.05, Appearance -> "Labeled"},
 {{n, 125}, 100, 200, 1, Appearance -> "Labeled"},
 {{T, 325}, 300, 400, 1, Appearance -> "Labeled"},
 {{m, 3}, 1, 10, 1, Appearance -> "Labeled"},
 SynchronousUpdating -> False,
 TrackedSymbols :> All]
$\endgroup$
3
  • $\begingroup$ Thanks. I see that you specified {n, 125}. In what units is now n? Is it like in n Quantity[1762*10^(26), "Meters"^(-3)], i.e. we have $125 * 10^{26} m^{-3}$? $\endgroup$ Commented Mar 26, 2022 at 16:17
  • $\begingroup$ Like all of your control variables, n is a number (dimensionless) which is multiplied by a Quantity inside the arguments of coherence. The arguments to coherence have the units of the associated quantities. You originally defined n as a number in the interval {100, 200}. I gave it an initial value (125) within that interval and specified a step size of 1. $\endgroup$ Commented Mar 26, 2022 at 16:36
  • $\begingroup$ Got it. I now changed {{n, 125}, 100, 200, 1, Appearance -> "Labeled"} to {{n, 176.2 * 10^(27)}, 170 * 10^(27), 180 * 10^(27), 10^26, Appearance -> "Labeled"} and again see no plot. But I will look into it, thanks. $\endgroup$ Commented Mar 26, 2022 at 16:42

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.