NDSolve issues with initial and boundary conditions (corrected characters)
- To: mathgroup at smc.vnet.net
- Subject: [mg119220] NDSolve issues with initial and boundary conditions (corrected characters)
- From: Arturo Amador <arturo.amador at ntnu.no>
- Date: Wed, 25 May 2011 06:59:12 -0400 (EDT)
Hi,
Sorry for the previous message, it had some weird characters I have =
corrected it and resend it.
I am having some issues when trying to solve a system of three coupled
differential equations numerically using NDSolve. I am trying to specify
boundary conditions for two of the variables in the same point (point at
L) and a boundary condition for the last variable at zero. The
mathematica code is:
IIn[24]:= vacuumlinearrhsp =
(-((4 E^(-4 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) g[t]^2)/
Lambda^4) - (12 g[t]^2)/(E^(2 t) Lambda^2 + 16 g[t]^2 p0[t])^(3/
2)) factorp[t]
Out[24]= -((E^(5 t)
Lambda^5 (-((E^(-4 t) Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^4) -
3/(E^(2 t) Lambda^2 + 16 g[t]^2 p0[t])^(3/2)))/(24 Pi^2))
In[25]:= vacuumcuadraticrhsg =
((24 E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) g[t]^4)/
Lambda^6 + (216 g[t]^4)/(E^(2 t) Lambda^2 + 16 g[t]^2 p0[t])^(5/
2)) factorg[t]
Out[25]= (2 E^(5 t) Lambda^5 g[
t]^4 ((E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^6 +
9/(E^(2 t) Lambda^2 + 16 g[t]^2 p0[t])^(5/2)))/Pi^2
In[26]:= vacuumcubicrhsh =
[(-((160 E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) g[t]^6)/
Lambda^8) - (4320 g[t]^6)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(7/
2)) factorh[t]
Out[26]= -((40 E^(5 t)
Lambda^5 g[
t]^6 (-((E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^8) -
27/(E^(2 t) Lambda^2 + 16 g[t]^2 p0[t])^(7/2)))/(3 Pi^2))
(*Declarations*)
Lambda = Sqrt[5] msig;
L = -17090/10000;
msig = 400;
mp = 0;
fp = 93;
lambda = 2 (msig^2 - mp^2)/fp^2;
gk0 = (lambda/24)^(1/2);
pk0 = 1/2 fp^2;
n = 4;
sol = NDSolve[{D[p0[t], t] == vacuumlinearrhsp,
4 D[g[t]^2, t] + 3 h[t] D[p0[t], t] == vacuumcuadraticrhsg,
D[h[t], t] == vacuumcubicrhsh, p0[L] == pk0, g[L]^2 == =
gk0^2,
h[0] == 0}, {p0, g, h}, {t, L, 0},
Method -> {StiffnessSwitching,
Method -> {ExplicitRungeKutta, Automatic}}];
with output:
Power::infy: Infinite expression 1/0. encountered. >>
Infinity::indet: Indeterminate expression (0. Sqrt[5] =
ComplexInfinity)/\[Pi]^2 encountered. >>
Power::infy: Infinite expression 1/0.^2 encountered. >>
Infinity::indet: Indeterminate expression (0. Sqrt[5] =
ComplexInfinity)/\[Pi]^2 encountered. >>
Power::infy: Infinite expression 1/0. encountered. >>
General::stop: Further output of Power::infy will be suppressed during =
this calculation. >>
Infinity::indet: Indeterminate expression 0. ComplexInfinity =
encountered. >>
General::stop: Further output of Infinity::indet will be suppressed =
during this calculation. >>
NDSolve::ndnum: Encountered non-numerical value for a derivative at t ===
-1.709. >>
I am sure there is no singularity. I am getting this output no matter =
what value I am giving for h[0], as long as I specify the boundary =
condition in a point that is not L, it gives me this same error message. =
I have tried h[0.9L] and still the same. When h[t]==0 The system =
reduces to this:
sol = NDSolve[{D[p0[t], t] == (Lambda Exp[t])^5/(
24 Pi^2) ( (n - 1)/(Lambda Exp[t])^3 +
3/((Lambda Exp[t])^2 + 16 p0[t] g[t]^2)^(3/2)),
D[g[t]^2, t] == ((Lambda Exp[t])^5 * g[t]^4 )/(
2 Pi^2) ((n - 1)/(Lambda Exp[t])^5 +
9/((Lambda Exp[t])^2 + 16 p0[t] g[t]^2)^(5/2)), p0[L] == pk0,=
g[L]^2 == gk0^2}, {p0, g}, {t, L, 0},
Method -> {StiffnessSwitching,
Method -> {ExplicitRungeKutta, Automatic}}];
For which I get nice solutions.
Thanks in advance
--
Arturo Amador Cruz
Stipendiat Norges =
Teknisk-Naturvitenskapelige Universitet (NTNU)
Institutt for Fysikk
Fakultet for =
Naturvitenskap og Teknologi
H=F8gskoleringen 5, =
Realfagbygget, E5-108
Telefon (735) 93366 NO-7491 =
TRONDHEIM
=
NORWAY
arturo.amador at ntnu.no
--Apple-Mail-2--555294862
Content-Type: text/html; charset="iso-8859-1"
Content-Transfer-Encoding: quoted-printable
X-Sun-Content-Length: 7914
<html><head></head><body style="word-wrap: break-word; =
-webkit-nbsp-mode: space; -webkit-line-break: after-white-space; =
">Hi,<br>Sorry for the previous message, it had some weird characters I =
have corrected it and resend it.<div><br>I am having some issues when =
trying to solve a system of three coupled differential equations =
numerically using NDSolve. I am trying to specify boundary conditions =
for two of the variables in the same point (point at L) and a boundary =
condition for the last variable at zero. The mathematica code =
is:<br><br>IIn[24]:= vacuumlinearrhsp = <br>(-((4 E^(-4 t) =
Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^2)/<br> Lambda^4) - =
(12 g[t]^2)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(3/<br> 2)) =
factorp[t]<br>Out[24]= -((E^(5 t)<br> Lambda^5 (-((E^(-4 t) =
Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^4) =
- <br> 3/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(3/2)))/(24 Pi^2))<br><br>In[25]:= vacuumcuadraticrhsg =
= <br>((24 E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^4)/<br> Lambda^6 + (216 =
g[t]^4)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(5/<br> 2)) =
factorg[t] <br>Out[25]= (2 E^(5 t) Lambda^5 g[<br> t]^4 =
((E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^6 =
+ <br> 9/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(5/2)))/Pi^2<br><br>In[26]:= vacuumcubicrhsh =
= <br>[(-((160 E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^6)/<br> Lambda^8) - =
(4320 g[t]^6)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(7/<br> 2)) =
factorh[t]<br>Out[26]= -((40 E^(5 t)<br> Lambda^5 =
g[<br> t]^6 (-((E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + =
n))/Lambda^8) - <br> 27/(E^(2 t) Lambda^2 + 16 =
g[t]^2 p0[t])^(7/2)))/(3 Pi^2))<br><br><br>(*Declarations*)<br>Lambda = =
Sqrt[5] msig;<br><br>L = -17090/10000; <br>msig = =
400; <br><br>mp = 0; <br>fp = 93; <br><br>lambda = =
2 (msig^2 - mp^2)/fp^2;<br>gk0 = (lambda/24)^(1/2);<br>pk0 = 1/2 =
fp^2;<br><br>n = 4; <br><br><br>sol = NDSolve[{D[p0[t], t] ===
vacuumlinearrhsp, <br> 4 D[g[t]^2, t] + 3 h[t] =
D[p0[t], t] == =
vacuumcuadraticrhsg, <br> D[h[t], t] == =
vacuumcubicrhsh, p0[L] == pk0, g[L]^2 == =
gk0^2, <br> h[0] == 0}, {p0, g, h}, {t, L, =
0}, <br> Method -> =
{StiffnessSwitching, <br> Method -> =
{ExplicitRungeKutta, Automatic}}];<br><br>with =
output:<br><br>Power::infy: Infinite expression 1/0. encountered. =
>><br><br>Infinity::indet: Indeterminate expression (0. Sqrt[5] =
ComplexInfinity)/\[Pi]^2 encountered. >><br><br>Power::infy: =
Infinite expression 1/0.^2 encountered. >><br><br>Infinity::indet: =
Indeterminate expression (0. Sqrt[5] ComplexInfinity)/\[Pi]^2 =
encountered. >><br><br>Power::infy: Infinite expression 1/0. =
encountered. >><br><br>General::stop: Further output of =
Power::infy will be suppressed during this calculation. =
>><br>Infinity::indet: Indeterminate expression 0. ComplexInfinity =
encountered. >><br>General::stop: Further output of =
Infinity::indet will be suppressed during this calculation. =
>><br>NDSolve::ndnum: Encountered non-numerical value for a =
derivative at t == -1.709. >><br><br>I am sure there is no =
singularity. I am getting this output no matter what value I am giving =
for h[0], as long as I specify the boundary condition in a point that is =
not L, it gives me this same error message. I have tried h[0.9L] and =
still the same. When h[t]==0 The system reduces to this:<br><br>sol =
= NDSolve[{D[p0[t], t] == (Lambda =
Exp[t])^5/(<br> 24 Pi^2) ( (n - 1)/(Lambda =
Exp[t])^3 =
+ <br> 3/((Lambda =
Exp[t])^2 + 16 p0[t] =
g[t]^2)^(3/2)), <br> D[g[t]^2, t] == ((Lambda =
Exp[t])^5 * g[t]^4 )/(<br> 2 Pi^2) ((n - =
1)/(Lambda Exp[t])^5 =
+ <br> 9/((Lambda =
Exp[t])^2 + 16 p0[t] g[t]^2)^(5/2)), p0[L] == =
pk0, <br> g[L]^2 == gk0^2}, {p0, g}, {t, L, =
0}, <br> Method -> =
{StiffnessSwitching, <br> Method -> =
{ExplicitRungeKutta, Automatic}}];<br><br>For which I get nice =
solutions. <br><br><br><br>Thanks in =
advance<br><br><br>--<br>Arturo Amador Cruz<br>Stipendiat<span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span>Norges =
Teknisk-Naturvitenskapelige Universitet (NTNU)<br><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span>Institutt =
for Fysikk<br><span class="Apple-tab-span" style="white-space: pre; =
"> </span><span class="Apple-tab-span" style="white-space: pre; =
"> </span><span class="Apple-tab-span" style="white-space: pre; =
"> </span><span class="Apple-tab-span" style="white-space: pre; =
"> </span><span class="Apple-tab-span" style="white-space: pre; =
"> </span><span class="Apple-tab-span" style="white-space: pre; =
"> </span>Fakultet for Naturvitenskap og Teknologi<br><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> =
</span>H=F8gskoleringen 5, Realfagbygget, E5-108<br>Telefon<span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span><span =
class="Apple-tab-span" style="white-space: pre; "> </span>(735) =
93366<span class="Apple-tab-span" style="white-space: pre; "> =
</span><span class="Apple-tab-span" style="white-space: pre; "> =
</span>NO-7491 TRONDHEIM<br><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span><span class="Apple-tab-span" =
style="white-space: pre; "> </span>NORWAY<br><a =
href="mailto:arturo.amador at ntnu.no">arturo.amador at ntnu.no</a><br><br></d=
iv></body></html>=
--Apple-Mail-2--555294862--