Problem-Solving Example #4
Implicit Differential Equations
An implicit (first order) ordinary differential equation cannot be reduced to the normal form
y' = f( x, y)
but rather has the general form
g( x, y, y') = 0
or g( x, y, y', y", y"', ...) = 0 for general Implicit ODEs
Systems of such equations, which may involve higher order derivatives, are solved in Fotran Calculus by a combined calculus process. The integration process is executed exactly as if the model involved explicit equations. Within the model, however, the equations are actually calculated by a general algebraic equation solver, such as AJAX.
As far as the specification of the integration process is concerned, there is absolutely no difference between the solution of implicit and explicit equations. In a given model some equations may be implicit while others are explicit. The calculus model itself, of course, is quite different, reflecting the different character of the mathematics.
Any Fotran Calculus integration solver may be used to solve implicit equations and all of the preceding discussions of integration apply to such equations. The preferred solvers, however, are JANUS and MERCURY since the fact that they execute the model only twice per integration step provides a very significant computational efficiency.
The procedure involved in solving this combined process is best described in terms of the following representative example. For descriptions of the solvers employed, refer to appropriate sections of this manual.
Example:
The implicit differential equations:
dx/dt + [1 + (dx/dt + dy/dt) / (x+y)] ^ .5 = 0
dy/dt + [1 + (dx/dt + dy/dt) / (x + y)] ^ .5 = 0
may be solved by the following program.
Program IDE
common /imp/ dxdt, x, dydt, y
read *, x, y, dydt, dydt ! initial conditions & initial guesses
read *, t, dt, tprint, tfinal
initiate Janus; for Deq;
equations dxdt / x, dydt / y; of t; step dt; to tf
tf = tprint
do while (tf .lt. tfinal) ! print every 'tprint' units of t
integrate deq; by Janus
print *, t, dxdt, x, dydt, y
tf = tf + tprint
end do
End
Controller Set (Ajax)
Summary = 0
End
Model Deq
common /imp/ dxdt, x, dydt, y
Find dxdt, dydt; In IDE(gx, gy); By Ajax(Set); To Match gx, gy
End
Model IDE(gx, gy) ! Implicit Equations
common /imp/ dxdt, x, dydt, y
rad = sqrt(1 + (dxdt + dydt) / (x + y))
gx = dxdt + rad
gy = dydt + rad
End