2
$\begingroup$

I'm trying to solve an ODE system with the NDSolve method. This my ODE system, with the BC and the functions V and dV defined:

ODE system, with BC in (4) and the function V(B) and dV(B) defined in (5)

    s = NDSolve[{F'[
  t] == -(2 /Sqrt[3]) (F[t]^2/2 + a[t]^2*  V[B[ t]])^(1/2)*F[t] - 
  a[t]^2*dV[B[t]], 
a'[t] == a[t]/Sqrt[3] (F[t]^2/2 + a[t]^2*  V[B[ t]])^(1/2), 
B'[t] == F[t], F[tmin] == Subscript[v, 0], a[tmin] == a0 , 
B[tmin] == Subscript[B, 0] }, {F, a, B}, {t, tmin, tmax}, Method->"ExplicitRungeKutta"]

Where the solution functions are:

F(t), a(t), B(t)

I need the solutions for a big interval, like:

t = [-10^(10), -10^(-10)]

So:

tmin = -10^(10)
tmax = -10^(-10)

But the system seems to be stiff. I'm trying to replace the "ExplicitRungeKutta" method for "StiffnessSwitching method", replacing the method with:

Method -> {"StiffnessSwitching", Method -> {"ExplicitRungeKutta", Automatic}}

But the system seems to remain stiff.

The complete code is:

Needs["DifferentialEquations`NDSolveProblems`"];
Needs["DifferentialEquations`NDSolveUtilities`"];
Needs["NumericalCalculus`"];

V[B_] := 3;
dV[B_] :=  V'[B];

B0 = 1;
F0 = 0;
a0 = 10^(-10);

tmin = -10^(10);
tmax = -10^(-10);

  s = NDSolve[{F'[
  t] == -(2 /Sqrt[3]) (F[t]^2/2 + a[t]^2*  V[B[ t]])^(1/2)*F[t] - 
  a[t]^2*dV[B[t]], 
a'[t] == a[t]/Sqrt[3] (F[t]^2/2 + a[t]^2*  V[B[ t]])^(1/2), 
B'[t] == F[t], F[tmin] == F0, a[tmin] == a0 , B[tmin] == B0 }, {F,
 a, B}, {t, tmin, tmax},  Method->"ExplicitRungeKutta"]

Question: Is there any way to obtain the solutions of these equations without stiffness? What method should I use?

I'd appreciate much all the answers.

$\endgroup$
6
  • $\begingroup$ What if you just use Method -> "StiffnessSwitching" without invoking RK at all? $\endgroup$ Jun 21, 2016 at 23:03
  • 2
    $\begingroup$ (1) What's v0? (2) What if it's a singularity and not stiffness? $\endgroup$
    – Michael E2
    Jun 22, 2016 at 1:02
  • $\begingroup$ Are tmin & tmax switched, or is tmin supposed to be greater than tmax? $\endgroup$
    – Michael E2
    Jun 22, 2016 at 3:09
  • $\begingroup$ V0 was F0, sorry for the mistake. And both variables tmin and tmax were negative (another mistake writting the code here). $\endgroup$ Jun 22, 2016 at 23:57
  • $\begingroup$ And the stiffnessSwitching method alone does not solve the (aparent) stiffness problem. $\endgroup$ Jun 22, 2016 at 23:59

1 Answer 1

3
$\begingroup$

Let's read the error message:

NDSolve`Iterate::ndsz: At t == -1.*10^10, step size is effectively zero; singularity or stiff system suspected.

The "step size is effectively zero" means that in floating-point arithmetic t + dt is equal to t for the computed time step dt. If this is the problem, increasing working precision might help.

de = {F'[t] == -(2/Sqrt[3]) (F[t]^2/2 + a[t]^2*V[B[t]])^(1/2)*F[t] - 
     a[t]^2*dV[B[t]], 
   a'[t] == a[t]/Sqrt[3] (F[t]^2/2 + a[t]^2*V[B[t]])^(1/2), 
   B'[t] == F[t]};
ics = {F[tmin] == F0, a[tmin] == a0, B[tmin] == B0};
s = 
 NDSolve[{de, ics}, {F, a, B}, {t, tmin, tmax},
  PrecisionGoal -> 8, AccuracyGoal -> 8, WorkingPrecision -> 20,
  Method -> "StiffnessSwitching"]

Mathematica graphics

Note: The first of the OP's differential equations has

F'[t] == <..stuff..> * F[t]

because dV[B] is zero and the second terms goes away. Since F[t] starts out at F[tmin] == 0, that means F[t] is a constant 0 (zero). Thus the third differential equation B'[t] == F[t] and its initial condition imply B[t] is a constant 1 (one).

$\endgroup$

Your Answer

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

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