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--