Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Topic: Problems with NDSolve
Replies: 0

 Search Thread: Advanced Search

 Steven Wilkinson Posts: 5 Registered: 12/7/04
Problems with NDSolve
Posted: Jul 25, 1996 10:39 AM
 Plain Text Reply

Mathgroupers,

I have a system of ODE, dx/dt = f(t,x,y), dy/dt = g(t,x,y), x(0) = x0, y(0)
= y0, for which NDSolve is giving unexpected results. Moreover, I have run
it with no problems with a fixed stepsize Runge-Kutta package so I think
the math is correct.

The functions f and g are piecewise defined (they are both essentially
continuous), and it is around where the pieces fit together that NDSolve
has problems. If you have some expertise with NDSolve and/or variable
stepsize algorithms could you take a look at the system and see if you can
spot what is happening?

The system:

f(t,x,y) =
If[ Sqrt[9x^4+12x^3+4x^2+4y^2] >= 0.1,
-Sign[x]*2y/Sqrt[9x^4+12x^3+4x^2+4y^2],
1/Sqrt[3x+2]
]

g(t,x,y) =
If[ Sqrt[9x^4+12x^3+4x^2+4y^2] >= 0.1,
-Sign[x]*(3x^2+2x)/Sqrt[9x^4+12x^3+4x^2+4y^2],
-Sqrt[3x+1]/Sqrt[3x+2]
]

x0 = -1, y0 = 0

The trajectory of this lies on the algebraic curve y^2 = x^2 + x^3. This
curve has a double point at (0,0) which is where the second definitions in
f and g kick in. To keep this simple I have given you the system so that
smoothness of solution is guaranteed only for t >= 0.

The problem with NDSolve is that it "bounces back" from the origin at
(0,0), instead of passing through it as the fixed step Runge-Kutta method
does. What follows is the NDSolve version with all the variables printed
out so that you can see it "bouncing back" at the origin. (Sorry about the
length of this, but I have found that NDSolve is picky about how you input
a piecewise defined equation. I need the flexibility of using a Module in
order to obtain the desired solution for all t.)

f[t_,x_,y_,defn_]:=
Module[{a},
Print["Using piecewise definition ",defn,". Variable t = ",t];
Print["(x,y) = ",{x,y}];
If[ defn==1,
a=-Sign[x]*2y/Sqrt[9x^4+12x^3+4x^2+4y^2],
a=1/Sqrt[3x+2]
];
Print["dx/dt = ",a];
Return[a];
];

g[t_,x_,y_,defn_]:=
Module[{a},
If[ defn==1,
a=-Sign[x]*(3x^2+2x)/Sqrt[9x^4+12x^3+4x^2+4y^2],
a=-Sqrt[3x+1]/Sqrt[3x+2]
];
Print["dy/dt = ",a];Print[""];
Return[a];
];

dx[t_,x_,y_]:=
If[ Sqrt[9x^4+12x^3+4x^2+4y^2 ] < 0.1,
f[t,x,y,2],
f[t,x,y,1]
];

dy[t_,x_,y_]:=
If[ Sqrt[9x^4+12x^3+4x^2+4y^2 ] < 0.1,
g[t,x,y,2],
g[t,x,y,1]
];

NDSolve[
{x'[t] == dx[t,x[t],y[t]], y'[t] == dy[t,x[t],y[t]], x[0] == -1,
y[0]== 0},
{x,y},{t,0,2}]

Thanks for any help you can give me.

Steve Wilkinson
wilkinson@nku.edu

© The Math Forum at NCTM 1994-2018. All Rights Reserved.