MathGroup Archive 2010

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

Search the Archive

Re: NDSolve and hybrid dynamics (Differential Algebraic Equation


> Thank you, Bob.  I guess I oversimplified the original problem.  The 
> variable y is state-dependent such as
> y[x] := Piecewise[{{1, x < 1/10}}, 2]
> I tried the following system of equations, but it failed:
> eqns = {x'[t] == y[t] + Sin[t],
>     y[t] == Piecewise[{{1, x[t] < 1/10}}, 2], x[-Pi] == 1/2};
> sol = {x[t], y[t]} /. NDSolve[eqns, {x[t], y[t]}, {t, -5, 5}]
> NDSolve::ndsz: At t == -3.33251, step size is effectively zero; 
> singularity or stiff system suspected.
> I presume this is what you meant by breaking the problem into two 
> regions, but that transition is dependent on x.
> With your suggestion of using Piecewise, I was able to run the same 
> system by creating a function y
> y//Clear
> y[xVar_?NumericQ] := Piecewise[{{1, xVar < 1/10}}, 2]
> eqns = {x'[t] == y[x[t]] + Sin[t], x[-Pi] == 1/2};
> xSoln = NDSolve[eqns, {x[t]}, {t, -5, 5}][[1, 1, 2]]
> Plot[xSoln, {t, -5, 5}, Frame -> True, Axes -> False,
>   PlotStyle -> Thick]
> Plot[y[xSoln], {t, -5, 5}, Frame -> True, Axes -> False,
>   PlotStyle -> Thick]
> Is there a way that I can get NDSolve to manage the y variable for me as 
> part of the state because my system will have multiple such non-smooth 
> variables?

As you have seen, in principle Mathematica can handle DAEs
(Differential-Algebraic Equation), but my experience is that if you can
formulate the problem as a pure DE (or set of those), then NDSolve will
in general perform better. My guess is that it then has more Methods to
choose from and most often makes a good choice. In your case the message
about the stiff system might be a hint: while the pure DE methods will
often handle stiff systems very well (presumably by using the stiffness
switching method) the DAE methods will not be able to do that. Whether
that is a problem in principle or just caused by the way the DAE solving
methods are implemented I can't say, but for the moment the situation is
like that: NDSolve does work better for pure DEs than for DAEs in many

Fortunately, for those cases where a reformulation is possible,
Mathematica can help you do the symbolic manipulations in an automated
way by using replacement rules or whatever fits for you. Here is a very
simple example for your case:

rules = {y[x_] :> Piecewise[{{1, x < 1/10}}, 2]};

eqns = {x'[t] == y[x[t]] + Sin[t], x[-Pi] == 1/2} /. rules;

xSoln = NDSolve[eqns, {x[t]}, {t, -5, 5}][[1, 1, 2]]



  • Prev by Date: How can I generate this graphic?
  • Next by Date: Re: How can I generate this graphic?
  • Previous by thread: Re: How can I generate this graphic?
  • Next by thread: Adding labels to a GraphicsGrid object (2, with code)