Numerical integration "backwards"

32 visualizaciones (últimos 30 días)
José Goulart
José Goulart el 2 de Abr. de 2011
Hello everyone,
I have a system of nonlinear equations on the form
x1' = f1(x1,x2,x3)
x2' = f2(x1,x2,x3)
x3' = f3(x1,x2,x3)
Let the flow F(x,t) be the solution to the system, that is, F(x0,t) as a function of t is an integral curve that goes through x0 at t=0.
I need the value of the inverse flow to some particular values of x0 and t0, that is, the value of F^-1(x0,t0) = F(x0,-t0).
Wen I use Matlab's odeXX functions to integrate F() numerically, I get well-behaved integral curves. Now, the problem is, when I try to use them to integrate F^-1() by reversing the time span parameter (or by exchanging the signs of the right hand sides of the differential equations), then apparently a numerical unstability occurs, and the integral curve grows extremely fast.
Does any one know how I can get the curve for F^-1(x0,t0), for particular t0 and x0?
Code used so far:
%code to compute the solution for F(x0,t0):
[t1,flow] = ode15s(@f1f2f3,[0 t0],x0);
%code to compute the solution for F^-1(x0,t0)=F(x0,-t0):
[t2,flowinv] = ode15s(@f1f2f3,[t0 0],x0);
  1 comentario
Teja Muppirala
Teja Muppirala el 4 de Abr. de 2011
It'd be interesting to see the actual equations though.

Iniciar sesión para comentar.

Respuesta aceptada

Andrew Newell
Andrew Newell el 3 de Abr. de 2011
You seem to be talking about equations for some physical phenomenon like fluid flow. In general, you can't integrate physical equations backwards. You can only do it if you don't have dissipative terms like friction. The numerical blowup is probably telling you that you have such a term. This is not really a MATLAB issue.
EDIT: You could get an intuitive feel for whether it is feasible or not by plotting a phase portrait of the system, i.e., choose a grid of reasonable starting values, calculate the flow, and then plot all the curves together. Do you have attractors?
  7 comentarios
Andrew Newell
Andrew Newell el 8 de Abr. de 2011
Yes. The phase portrait is just a picture of flow lines while the Poincaré map is a less intuitive plot that looks at intersections of trajectories in some traversal section in phase space. Sorry, I can't really explain it in a comment. The point is that if the orbits are periodic they'll keep intersecting at the same points, but if it's chaotic you'll get a much more complex pattern.
José Goulart
José Goulart el 8 de Abr. de 2011
Even when I start the simulation in a state x0 very far from the origin (picking state components with even rather unrealistic values for the quantities they represent), the computation of the time-reversed trajectory only works correctly for a very small interval of time (some miliseconds). After this, a numerical blowup still occurs.
Thinking about your previous comment, I'm starting to think it will be really difficult to get this time-reversed trajectory for my system, because the equations have very strong dissipative terms. Computing the conventional (forward) trajectory numerically, I see that even for initial states quite far from the origin, the "decaying time" (i.e., the time required for the trajectory to get very close to the origin) is quite short (the observed times were never longer than about 30 ms).
Moreover, when I compute the Jacobian of the equations around the initial state, its eigenvalues have a really large magnitude. So, I conclude that by inverting the signs of the equations, I have a very unstable set of equations that would need and infinitesimal step size to be numerically integrated.

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Programming en Help Center y File Exchange.

Productos

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by