2
$\begingroup$
   we = 30*10^-6; P0 = 3; \[Gamma] = 0.070; g = 10; \[Rho] = 1000; c = 
   3*10^8;
   \[Kappa][r_] := (z''[r] + z'[r]/r); \[CapitalPi][r_] := (4 P0)/(
    Pi*we^2*c) Exp[-((2*r^2)/we^2)];
   Solsol = NDSolve[{-\[Gamma]*\[Kappa][r] + \[Rho]*g*z[r] - \[CapitalPi][
   r] == 0, z[50*we] == 0, z'[0] == 0}, z[r], {r, 0, 10 we}]
   g = Plot[Evaluate[{z[r], z'[r], z''[r]} /. sol], {r, 0, 10 we}, 
    PlotStyle -> Automatic]
$\endgroup$
1
  • $\begingroup$ You have a differential equation, therefore you need NDSolve instead of NSolve. Further, Infinity is not a number and can not be used as function argument. $\endgroup$ Dec 25, 2022 at 9:47

1 Answer 1

3
$\begingroup$

You have singularity at r=0 due to (z''[r] + z'[r]/r) so can't start at zero. Try $MachineEpsilon

ClearAll["Global`*"]
we = 30*10^-6;
P0 = 3;
γ = 0.070;
g = 10;
ρ = 1000;
c = 3*10^8;
ode0 = (z''[r] + z'[r]/r)
h = (4 P0)/(Pi*we^2*c) Exp[-((2*r^2)/we^2)];
ode = -γ*ode0 + ρ*g*z[r] - h == 0
ic = {z[50*we] == 0, z'[$MachineEpsilon] == 0};
sol = NDSolve[{ode, ic}, z[r], {r, $MachineEpsilon, 10 we}]

Mathematica graphics

g = Plot[Evaluate[{z[r], z'[r], z''[r]} /. sol], {r, $MachineEpsilon, 10 we}, 
  PlotStyle -> Automatic]

Mathematica graphics

Update

To plot the derivatives, you can do this. You can actually have NDSolve solve for these.

sol = NDSolve[{ode, ic}, {z, z', z''}, {r, $MachineEpsilon, 10 we}]

Then do

g = Plot[
  Evaluate[{z[r], z'[r], z''[r]} /. sol], {r, $MachineEpsilon, 10 we}, 
  PlotStyle -> Automatic, PlotRange -> All]

Now 3 curves show up. But the scales are so different, hard to see them all on same plot. You can plot each on separate plots? Something like

opt = {PlotRange -> All, ImageSize -> 300, GridLines -> Automatic, 
   GridLinesStyle -> LightGray, PlotStyle -> Red};
p1 = Plot[Evaluate[z[r] /. sol], {r, $MachineEpsilon, 10 we}, 
   Evaluate@opt, PlotLabel -> "z(r)"];
p2 = Plot[Evaluate[z'[r] /. sol], {r, $MachineEpsilon, 10 we}, 
   Evaluate@opt, PlotLabel -> "z'(r)"];
p3 = Plot[Evaluate[z''[r] /. sol], {r, $MachineEpsilon, 10 we}, 
   Evaluate@opt, PlotLabel -> "z''(r)"];
Column[{p1, p2, p3}]

Mathematica graphics

I am trying to vary [Gamma] using sol = ParametricNDSolve[{ode, ic}, z, {r, $MachineEpsilon, 200 we}, {[Gamma]}] then I want to calculate A0 = -[Gammaa^2( (z''[[Gamma]][r] - z'[[Gamma]][r]/r for two differnt value of [Gamma] ) 0.030 and 0.070.

You can use ParametricNDSolve

ClearAll["Global`*"]
we = 30*10^-6;
P0 = 3;
g = 10;
ρ = 1000;
c = 3*10^8;
ode0 = (z''[r] + z'[r]/r)
h = (4 P0)/(Pi*we^2*c) Exp[-((2*r^2)/we^2)];
ode = -γ*ode0 + ρ*g*z[r] - h == 0
ic = {z[50*we] == 0, z'[$MachineEpsilon] == 0};
sol=ParametricNDSolve[{ode,ic},z,{r,$MachineEpsilon,10 we},γ]

z1 = z[0.030] /. sol;
z2 = z[0.070] /. sol;
Plot[{z1[r], z2[r]}, {r, $MachineEpsilon, 10 we}, 
 PlotLegends -> {"0.03", "0.07"}]

Mathematica graphics

γ = 0.030;
z1 = z[γ] /. sol
a = 1;
A0 = -γ*a^2*(z1''[r] - z1'[r]/r)

γ = 0.070;
z2 = z[γ] /. sol
A0 = -γ*a^2*(z2''[r] - z2'[r]/r)
$\endgroup$
8
  • $\begingroup$ Thanks, but it it showing only z(r) its derivatives z'(r) z''(r) not showing. $\endgroup$ Dec 25, 2022 at 11:22
  • $\begingroup$ thanks! derivative z''(r) near 0 is not so clear, Please also try to plot (z''(r)+z'(r)/r)/2 which is the mean curvature. $\endgroup$ Dec 25, 2022 at 12:13
  • $\begingroup$ Please suggest how to use ParametricNDSolve so that I can vary "we" for two different values. $\endgroup$ Jan 12, 2023 at 11:27
  • $\begingroup$ @GopalVerma I do not think you can use the location of initial conditions as a parameter. But why not simply make a function and pass it we value, and solve the same equation each time you want to change we by making a call to this function with different we? This is much simpler solution and you do not have to deal with complexity of ParametricNDSolve. $\endgroup$
    – Nasser
    Jan 12, 2023 at 23:00
  • $\begingroup$ Thanks, I got it. $\endgroup$ Jan 13, 2023 at 6:41

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.