# DIFFEQ: Solve (Systems of) Differential Equations

### Numeric solution of ordinary (ODE) or parabolic (PDE) differential equations in a 1-line statement with an additional callback subroutine.

• convert any higher order differential equations to systems of first order equations ⇾ example, e.g.:
• y1" + a*y1' + k*y1 + F = 0
• convert to
• y1' = y2
• y2' = -a*y2 - k*y1 - F
• DIFFEQ(CallBack=subname, T=t, Y=z, DY=dz, T0=t0, Go=go, TOLerance=err, DT1=Step1, itrOK=ok)
• Update any f(y,t,...)-dependencies in callback subroutine sub:
•  required keywords type mini sample Keyword sequence is irrelevant CallBack SUB CB=sub subroutine sub called back from DIFFEQ T NUM T=x the independent variable Y NUM Y=y if only 1 ODE: the dependent variable arr Y=vec dimension(N) : N ODEs, N dependent variables DY arr DY=dy same dimension as Y, to be set in callback subroutine sub optional keywords Go num G=50 max iterations (default is -1: from T0 to T) T0 num T0=3 initial value of T (default = 0) DT1 num DT1=1e-8 initial stepsize (default = 0.001) TOLerance num TOL=0.01 error limit for each iteration (def = 1E-5) NonZero num NZ=1e-8 y > 0 error control threshold to avoid numeric problems caused by small values of the dependent variable(s) 0 NZ=0 y error control OFF (allow y < 0 results again) itrOK NUM OK=L L==1: after a successful iteration, else 0 use for output, or to update any f(y,t,..) work ERror LBL ER=99 on error jump to label 99
• solve systems of _or_
• ∂y/∂t = ∂ 2(E⋅y)/∂x 2 + ∂(v⋅y)/∂x + R(x, y, y2, ..., t)
• Eddy diffusion E, may be a function of t and y, physical dimension e.g. m 3/m/sec, or miles 2/hour
• linear velocity v, may be a function of t and y, physical dimension e.g. m 3/m 2/sec, or miles/hour
• reaction rate R mole/m 3/sec or cars/mile/hour.
R usually couples multiple reaction-diffusion equations
• called without DY ...)
• boundary conditions at
• x = xfeed: d(E⋅y)/dx = (v⋅y) - (Vfeed⋅Yfeed)
• x = xout: d(E⋅y)/dx = 0
 required keywords type mini sample Keyword sequence is irrelevant CallBack SUB CB=sub subroutine sub called back from DIFFEQ T NUM T=t the independent variable. Physical dimension e.g. sec Y arr Y=vec array of length N: N/M coupled PDEs (N == M if only 1 PDE). Physical dimension e.g. mole/m 3 or cars/mile X arr X=xvec array of length M: x pillar positions (variable space OK). Physical dimension e.g. m, evaluated at M positions YFeed arr YF=inp array of length N/M (input at x=0 or x=xmax, see Velocity) Velocity arr V=flow array of length N (space dependent), N/M (space constant) flow > 0 enters at x=0, flow < 0 at x=xmax Eddy arr E=dif array of length N (space dependent), N/M (space constant) Rateatx arr R=r array of length N: differential in/out along X optional keyword BCFdiff num bcf=0 differences: -1:backward, 0=central, +1=forward
• REAL :: y(2), dy(2)
• ! set t and t0 (and optional parameters)
• DIFFEQ(CallBack=subname, T=t, Y=y, DY=dy, T0=t0)
• ...
• SUBROUTINE subname
• dy(1) = y(2)
• dy(2) = -y2 - y1 - SIN(t)
• END Support HicEst   ⇾ Impressum