Re: NDsolve problem
- To: mathgroup at smc.vnet.net
- Subject: [mg110721] Re: NDsolve problem
- From: Peter Pein <petsie at dordos.net>
- Date: Sat, 3 Jul 2010 08:18:39 -0400 (EDT)
- References: <i0kif0$183$1@smc.vnet.net>
Am Fri, 2 Jul 2010 11:27:28 +0000 (UTC)
schrieb "ntina.86 at libero.it" <ntina.86 at libero.it>:
> Hello,
> I haven't been using Mathematica for some years and now I need to ask
> you for help with NDSolve.
> I'm using version 6.0.
>
> I want to solve the radial Schrodinger equation with a Yukawa
> potential, and I started trying to solve this:
>
> a = 500
> b = 10;
> s = NDSolve[{y''[x] + (2/x) y'[x] + (1 + (2*a/x)*E^(-b*x)) (y[x]) ==
> 0, y [0.000001] == 1, y'[0.000001] == -a}, y, {x, 0.000001,
> 60},PrecisionGoal -> \ [Infinity], MaxSteps -> Infinity]
>
> This gives me results that might be correct, but the problem comes
> when I try to plot:
>
> Plot[Evaluate[(x*y[x])^2 + (((x - 0.5 Pi)*y[x - 0.5
> Pi]))^2 /.s[[1]]], {x, 10, 60}, PlotRange -> All]
>
> It should be a constant (i.e. the sum of the square of a cosine and
> of a sine, which is the solution of the Schrod.) at least for x>10,
> instead I get a badly oscillating function.
>
> Is it a problem of the method of integration or I'm doing something
> wrong? I tried to change it, using ExplicitRK or the
> SymplecticPartitionedRK, but the latter gives me these errors:
>
>
> ---------------------------------------------------------------------------------------------------------
> NDSolve::nlnum1: The function value {-21065.3 (1.23594*10^7+82.8327 \
> TemporaryVariable[2][1.1462*10^-6])} is not a list of numbers with \
> dimensions {1} when the arguments are {1.1462*10^-6,0.494377}.
>
> InterpolatingFunction::dmval: Input value {28.4292} lies outside the \
> range of data in the interpolating function. Extrapolation will be \
> used. >>
>
> InterpolatingFunction::dmval: Input value {30} lies outside the range
> \ of data in the interpolating function. Extrapolation will be used.
> >>
>
> NDSolve::nlnum1: The function value {-21193.5 (1.23496*10^7+82.3232 \
> TemporaryVariable[2][1.14631*10^-6])} is not a list of numbers with \
> dimensions {1} when the arguments are {1.14631*10^-6,0.493983}.
>
> InterpolatingFunction::dmval: Input value {28.4292} lies outside the \
> range of data in the interpolating function. Extrapolation will be \
> used. >>
>
> General::stop: Further output of InterpolatingFunction::dmval will be
> \ suppressed during this calculation. >>
>
> NDSolve::nlnum1: The function value {-21321.9 (1.23397*10^7+81.8195 \
> TemporaryVariable[2][1.14643*10^-6])} is not a list of numbers with \
> dimensions {1} when the arguments are {1.14643*10^-6,0.493589}.
>
> General::stop: Further output of NDSolve::nlnum1 will be suppressed \
> during this calculation. >>
>
>
> -------------------------------------------------------------------------------------------------------
>
>
> Thank you in advance.
> Regards,
> Valentina
>
>
>
Hi Valentina,
it's strange. With Mathematica 7 but without the Options for
NIntegrate, I get a plot-function which is (almost) constant a little
bit below 1/1000:
In[1]:= a=500;
b=10;
\[CurlyEpsilon]=10^-6;
s=NDSolve[{y''[x]+(2/x) y'[x]+(1+(2*a/x)*E^(-b*x))
(y[x])==0,y[\[CurlyEpsilon]]==1,y'[\[CurlyEpsilon]]==-a},y,{x,\[CurlyEpsilon],60}];
In[5]:= Plot[Evaluate[(x*y[x])^2+(((x-0.5 Pi)*y[x-0.5
Pi]))^2/.s[[1]]],{x,\[CurlyEpsilon]+Pi/2,60},PlotRange->All]
(* see: http://dl.dropbox.com/u/3030567/Mathematica/schroed.png *)
In[6]:=
Evaluate[(x*y[x])^2+(((x-0.5 Pi)*y[x-0.5 Pi]))^2/.s[[1]]]/.x->{10,60}
Out[6]= {0.0000979763,0.0000976444}
does that help?
Peter
