MathGroup Archive 2011

[Date Index] [Thread Index] [Author Index]

Search the Archive

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 =&nbsp;<br>(-((4 E^(-4 t) =
Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^2)/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Lambda^4) - =
(12 g[t]^2)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(3/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2)) =
factorp[t]<br>Out[24]= -((E^(5 t)<br>&nbsp;&nbsp;Lambda^5 (-((E^(-4 t) =
Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^4) =
-&nbsp;<br>&nbsp;&nbsp;&nbsp;3/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(3/2)))/(24 Pi^2))<br><br>In[25]:= vacuumcuadraticrhsg =
=&nbsp;<br>((24 E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^4)/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Lambda^6 + (216 =
g[t]^4)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(5/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2)) =
factorg[t]&nbsp;<br>Out[25]= (2 E^(5 t) Lambda^5 g[<br>&nbsp;t]^4 =
((E^(-6 t) Sqrt[E^(2 t) Lambda^2] (-1 + n))/Lambda^6 =
+&nbsp;<br>&nbsp;&nbsp;9/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(5/2)))/Pi^2<br><br>In[26]:= vacuumcubicrhsh =
=&nbsp;<br>[(-((160 E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + n) =
g[t]^6)/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Lambda^8) - =
(4320 g[t]^6)/(E^(2 t) Lambda^2 + 16 g[t]^2 =
p0[t])^(7/<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2)) =
factorh[t]<br>Out[26]= -((40 E^(5 t)<br>&nbsp;&nbsp;Lambda^5 =
g[<br>&nbsp;&nbsp;t]^6 (-((E^(-8 t) Sqrt[E^(2 t) Lambda^2] (-1 + =
n))/Lambda^8) -&nbsp;<br>&nbsp;&nbsp;&nbsp;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 =&nbsp;-17090/10000;&nbsp;<br>msig = =
400;&nbsp;<br><br>mp = 0;&nbsp;<br>fp = 93;&nbsp;<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;&nbsp;<br><br><br>sol = NDSolve[{D[p0[t], t] ===
 vacuumlinearrhsp,&nbsp;<br>&nbsp;&nbsp;&nbsp;4 D[g[t]^2, t] + 3 h[t] =
D[p0[t], t] == =
vacuumcuadraticrhsg,&nbsp;<br>&nbsp;&nbsp;&nbsp;D[h[t], t] == =
vacuumcubicrhsh, p0[L] == pk0, g[L]^2 == =
gk0^2,&nbsp;<br>&nbsp;&nbsp;&nbsp;h[0] == 0}, {p0, g, h}, {t, L, =
0},&nbsp;<br>&nbsp;&nbsp;Method -&gt; =
{StiffnessSwitching,&nbsp;<br>&nbsp;&nbsp;&nbsp;&nbsp;Method -&gt; =
{ExplicitRungeKutta, Automatic}}];<br><br>with =
output:<br><br>Power::infy: Infinite expression 1/0. encountered. =
&gt;&gt;<br><br>Infinity::indet: Indeterminate expression (0. Sqrt[5] =
ComplexInfinity)/\[Pi]^2 encountered. &gt;&gt;<br><br>Power::infy: =
Infinite expression 1/0.^2 encountered. &gt;&gt;<br><br>Infinity::indet: =
Indeterminate expression (0. Sqrt[5] ComplexInfinity)/\[Pi]^2 =
encountered. &gt;&gt;<br><br>Power::infy: Infinite expression 1/0. =
encountered. &gt;&gt;<br><br>General::stop: Further output of =
Power::infy will be suppressed during this calculation. =
&gt;&gt;<br>Infinity::indet: Indeterminate expression 0. ComplexInfinity =
encountered. &gt;&gt;<br>General::stop: Further output of =
Infinity::indet will be suppressed during this calculation. =
&gt;&gt;<br>NDSolve::ndnum: Encountered non-numerical value for a =
derivative at t == -1.709. &gt;&gt;<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>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;24 Pi^2) ( (n - 1)/(Lambda =
Exp[t])^3 =
+&nbsp;<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;3/((Lambda =
Exp[t])^2 + 16 p0[t] =
g[t]^2)^(3/2)),&nbsp;<br>&nbsp;&nbsp;&nbsp;D[g[t]^2, t] == ((Lambda =
Exp[t])^5 * g[t]^4 )/(<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;2 Pi^2) ((n - =
1)/(Lambda Exp[t])^5 =
+&nbsp;<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;9/((Lambda =
Exp[t])^2 + 16 p0[t] g[t]^2)^(5/2)), p0[L] == =
pk0,&nbsp;<br>&nbsp;&nbsp;&nbsp;g[L]^2 == gk0^2}, {p0, g}, {t, L, =
0},&nbsp;<br>&nbsp;&nbsp;Method -&gt; =
{StiffnessSwitching,&nbsp;<br>&nbsp;&nbsp;&nbsp;&nbsp;Method -&gt; =
{ExplicitRungeKutta, Automatic}}];<br><br>For which I get nice =
solutions.&nbsp;<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--


  • Prev by Date: Re: How to plot y=f(x) or given y interval?
  • Next by Date: Re: Position
  • Previous by thread: Re: Again : Is there a BNF for Mathematica?
  • Next by thread: Re: NDSolve issues with initial and boundary conditions (corrected characters)