0
$\begingroup$

I have a function like this Exp[-(\[Psi]+\[Delta]* x)^2/2]. I need to numerically integrate this function and represent using ListContourPlot. Hence I Use the code like this

t2 = Flatten[ Table[{\[Psi], \[Delta], If[0 < NIntegrate[ Exp[-(\[Psi]+\[Delta]* x)^2/2], {x, 0, \[Infinity]}] < 3, 0, 1]}, {\[Psi], -0.01, -2, -01}, {\[Delta], 0.1, 5, 01}], 1] ListContourPlot[t2, PlotTheme -> "Detailed", Evaluate@plotset2, PlotRange -> {{-0.01, -1}, {0, 4}}, PlotLegends -> None, FrameLabel -> {"\[Psi]", "\[Delta]"}, Contours -> 1]

and I obtain the plot as:

enter image description here

Now, I need to extract the line that separates these two regions (blue and cream region). I know that the boundary indicates the region where the function has value not within the interval. So I think If there is a way to extract the values of psi and delta when the value of the integral changes will serve the purpose. I am not sure.

Please help me out. Thanks in advance.

Improve this question
$\endgroup$
1
  • $\begingroup$ the picture does not match the updated code. $\endgroup$ Commented Jan 28, 2021 at 6:34

2 Answers 2

Reset to default
3
$\begingroup$ 
Clear["Global`*"]

The integral can be done analytically

int[ψ_, δ_] = Assuming[Element[{ψ, δ}, Reals], Integrate[Exp[-(ψ*δ*x)^2/2], {x, 0, ∞}]] (* Sqrt[π/2]/(Abs[δ] Abs[ψ]) *) f[ψ_, δ_] = Piecewise[{{1, int[ψ, δ] <= 0 || int[ψ, δ] >= 3}}] // Simplify (* Piecewise[ {{1, Sqrt[2*Pi]/(Abs[δ]*Abs[ψ]) >= 6}}, 0] *)

The boundary is int[ψ, δ] == 3

Simplify[Solve[{int[ψ, δ] == 3, ψ < 0, δ > 0}, δ][[1]], ψ < 0] (* {δ -> -(Sqrt[(π/2)]/(3 ψ))} *) Show[ ContourPlot[f[ψ, δ], {ψ, -0.01, -1}, {δ, 0, 4}, PlotTheme -> "Detailed", PlotRange -> {{-0.01, -1}, {0, 4}}, PlotLegends -> None, FrameLabel -> (Style[#, 14, Bold] & /@ {"ψ", "δ"})], Plot[-Sqrt[Pi/2]/(3 ψ), {ψ, -0.01, -1}, PlotStyle -> {Thick, Red}]]

enter image description here

EDIT: For the revised function

t2 = Flatten[ Table[{ψ, δ, If[0 < NIntegrate[Exp[-(ψ + δ*x)^2/2], {x, 0, ∞}] < 3, 0, 1]}, {ψ, -0.01, -2, -0.025}, {δ, 0.01, 5, 0.025}], 1];

Select the {ψ, δ} pairs on the boundary

boundary = First[MaximalBy[#, Last]] & /@ GatherBy[Most /@ Select[t2, #[[3]] == 1 &], First];

Use FindFormula to find a function that fits the boundary values

b[ψ_] = FindFormula[boundary, ψ] (* 0.403272 - 0.369685 ψ - 0.0845795 ψ^2 *) Show[ ListContourPlot[t2, PlotTheme -> "Detailed", PlotRange -> {{-0.01, -2}, {0, 2}}, PlotLegends -> None, FrameLabel -> {"ψ", "δ"}], Plot[b[ψ], {ψ, -2, -0.01}, PlotStyle -> {Red, Thick}]]

enter image description here

Improve this answer
$\endgroup$
1
  • $\begingroup$ Thanks.. But what if my function is not solvable analytically?? For example, Please see the edited question... $\endgroup$ Commented Jan 28, 2021 at 5:50
1
$\begingroup$

You can get the desired contour using the analytical solution for Integrate[Exp[-(ψ + δ x)^2/2], {x, 0, ∞}] under the assumptions ψ and δ are real and δ is non-zero:

ClearAll[ψδcontour] ψδcontour[ψ_, δ_] := FullSimplify[Integrate[Exp[-(ψ + δ x)^2/2], {x, 0, ∞}], (ψ | δ) ∈ Reals && δ != 0] ψδcontour[ψ, δ] // TeXForm

$$\frac{\sqrt{\frac{\pi }{2}} \left(\left| \delta \right| -\delta \text{erf}\left(\frac{\psi }{\sqrt{2}}\right)\right)}{\delta ^2}$$

ContourPlot[Evaluate[ψδcontour[ψ, δ]], {ψ, -2, 0}, {δ, 0, 3}, Contours -> {3}, ContourLabels -> True, PlotRange -> All, ContourStyle -> Directive[Red, Thick], ContourShading -> {LightBlue, Yellow}, BaseStyle -> FontSize -> 16, FrameLabel -> {ψ, δ}]

enter image description here

Plotting multiple contours:

ContourPlot[Evaluate[ψδcontour[ψ, δ]], {ψ, -2, 0}, {δ, 0, 2}, PlotRange -> All, Contours -> Range[4], ContourLabels -> True, ContourShading -> None, ContourStyle -> (Directive[Thick, #] & /@ { Red, Green, Blue, Orange}), BaseStyle -> FontSize -> 16, FrameLabel -> {ψ, δ}]

enter image description here

"what if my function is not solvable analytically?"

You can extract the contour line coordinates from ListContourPlot output (using Cases) to construct a BSplineFunction and use it with ParametricPlot:

Using Bob Hanlon's setup to get t2:

t2 = Flatten[Table[{ψ, δ, If[0 < NIntegrate[Exp[-(ψ + δ*x)^2/2], {x, 0, ∞}] < 3, 0, 1]}, {ψ, -0.01, -2, -0.025}, {δ, 0.01, 5, 0.025}], 1]; lcp = ListContourPlot[t2, FrameLabel -> {"ψ", "δ"}, PlotRange -> All, BaseStyle -> FontSize -> 16, FrameLabel -> {ψ, δ}, ImageSize -> Medium]; bSF = Cases[Normal[lcp], Line[x_] :> BSplineFunction[x], All][[1]]; Row[{lcp, Show[lcp, ParametricPlot[bSF[t], {t, 0, 1}, PlotStyle -> Directive[Thick, Red]]]}, Spacer[10]]

enter image description here

Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.