NAME

txdodeint - integrate a given ordinary differential equation system along a coordinate

SYNOPSIS

pipe | txdodeint [parameters] [<ode> [<initval> [initval ...]]] | pipe

DESCRIPTION

By default, the first column in the data is used as time coordinate for which the expression(s) of the ODE are defined. The derivative can be a system of ODEs. The time derivative values shall be stored in the array members D0 .. Dn, with the current corresponding variable values available in V0 .. Vn and the corresponding values of the (interpolated) auxilliary timeseries from the piped data in [1] .. [m]. There are also the auxillary array values A0 .. Ak (to be used at your leisure to store/load values) and the constants C0 .. Cl (initialised from program parameters).

For convenience, the basic setup of time column, ODE, and intial values can be given directly on the command line without mentioning parameter names. The integrated values are written out at the points in time given by the input data. For each variable (size of the array V), a column is appended.

Example for a simple system (constant acceleration):

D0=V1; D1=A0; D2=C0*[1]

Assuming the time in column 1 and the constant acceleration in C0, this computes the evolution of the covered distance V0 via the numerically accelerated speed V1 and also directly from the analytically accelerated speed as variable V2, to give a hint about the accuracy of the numerical integration.

The employed integration method is a standard Runge-Kutta scheme with up to 4 stages, which should be fine for any application where you consider a humble Perl script for your numerical integration. A comparison of the fully numerical with the fully analytical solution can be constructed via the pipeline

txdconstruct -n=11 '[1]=C0-1' \
| txdodeint --rksteps=1 'D0=V1; D1=3' 0 0 \
| txdcalc '[4]=3/2*[1]**2; [5]=[4] ? ([2]-[4])/[4] : 0' \
| txdfilter -N=%g 'integration test' 't/s' 's/m' 'v/m' 's_ref/m' 'error'

With rksteps>1, you will not see any difference in this example. In general, the choice of integration stages, time step and interpolation might have significant influence on your results. A simple test of the quality of the chosen integration employs trivial polynomials:

txdconstruct -n=11 '[1]=C0-1' \
| txdodeint --rksteps=2 --timediv=1 \
  'D0=1; D1=2*[1]; D2=3*[1]**2; D3=4*[1]**3' \
     0     0         0            0 \
| txdcalc '[2]-=[1]; [3]-=[1]**2; [4]-=[1]**3; [5]-=[1]**4' \
| txdfilter -N=g 'integration order test' x err0 err1 err2 err3

These polynomials can actually be solved exactly to machine precision (depending on rksteps value) and smaller time steps would introduce rounding errors here from the summation. Finally, another classic example, the Lorenz attractor:

txdconstruct -n=5001 '[1]=(C0-1)/100' | txdodeint --timediv=1 \
  'D0=10*(V1-V0); D1=28*V0-V1-V0*V2; D2=-8/3*V2+V0*V1' \
  20 -20 1

More practical applications actually have some more data columns besides the time in the input data (measurements) and involve derivative expressions that make use of this data-driven time-dependence.

PARAMETERS

These are the general rules for specifying parameters to this program:

txdodeint -s -xyz -s=value --long --long=value [--] [files/stuff]

You mention the options to change parameters in any order or even multiple times. They are processed in the oder given, later operations overriding/extending earlier settings. Using the separator "--" stops option parsing An only mentioned short/long name (no "=value") means setting to 1, which is true in the logical sense. Also, prepending + instead of the usual - negates this, setting the value to 0 (false). Specifying "-s" and "--long" is the same as "-s=1" and "--long=1", while "+s" and "++long" is the sames as "-s=0" and "--long=0".

There are also different operators than just "=" available, notably ".=", "+=", "-=", "*=" and "/=" for concatenation / appending array/hash elements and scalar arithmetic operations on the value. Arrays are appended to via "array.=element", hash elements are set via "hash.=name=value". You can also set more array/hash elements by specifying a separator after the long parameter line like this for comma separation:

--array/,/=1,2,3  --hash/,/=name=val,name2=val2

The available parameters are these, default values (in Perl-compatible syntax) at the time of generating this document following the long/short names:

black (scalar)

0

ignore whitespace at beginning and end of line (disables strict mode) (from Text::NumericData)

comchar (scalar)

undef

comment character (if not set, deduce from data or use #) (from Text::NumericData)

comregex (scalar)

'[#%]*[^\\S\\015\\012]*'

regex for matching comments (from Text::NumericData)

config, I (array)

[]

Which configfile(s) to use (overriding automatic search in likely paths); special: just -I or --config causes printing a current config file to STDOUT

const, n (array)

[]

array of constants to use in ODE

debug, d (scalar)

0

give some info that may help debugging, >1 increasing verbosity

empty (scalar)

0

treat empty lines as empty data sets, preserving them in output (from Text::NumericData)

fill (scalar)

undef

fill value for undefined data (from Text::NumericData)

help, h (scalar)

0

Show the help message. Value 1..9: help level, par: help for paramter par (long name) only.

Additional fun with negative values, optionally followed by comma-separated list of parameter names: -1: list par names, -2: list one line per name, -3: -2 without builtins, -10: dump values (Perl style), -11: dump values (lines), -100: print POD.

interpolate (scalar)

'spline'

type of interpolation to use for intermediate points: spline or linear

lineend (scalar)

undef

line ending to use: (DOS, MAC, UNIX or be explicit if you can, taken from data if undefined, finally resorting to UNIX) (from Text::NumericData)

numformat, N (array)

[]

printf formats to use (if there is no "%" present at all, one will be prepended) (from Text::NumericData)

numregex (scalar)

'[\\+\\-]?\\d*\\.?\\d*[eE]?\\+?\\-?\\d*'

regex for matching numbers (from Text::NumericData)

ode, e (scalar)

'D0 = 1'

The ordinary differential equation system. The return value of the generated function does not matter, only that you set the values of the D array.

outsep (scalar)

undef

use this separator for output (leave undefined to use input separator, fallback to TAB) (from Text::NumericData)

plainperl (scalar)

0

Use plain Perl syntax for formula for full force without confusing the intermediate parser.

quote (scalar)

undef

quote titles (from Text::NumericData)

quotechar (scalar)

undef

quote character to use (derived from input or ") (from Text::NumericData)

rksteps, k (scalar)

'4'

steps (stages) of the RK integration scheme, only 4 supported right now

separator (scalar)

undef

use this separator for input (otherwise deduce from data; TAB is another way to say "tabulator", fallback is ) (from Text::NumericData)

strict, S (scalar)

0

strictly split data lines at configured separator (otherwise more fuzzy logic is involved) (from Text::NumericData)

text, T (scalar)

1

allow text as data (not first column) (from Text::NumericData)

timecol, t (scalar)

1

Column for the (time) coordinate the ODE shall be advanced on. In the ODE, you can access it via [1] if the column is 1 (just like any other variable of the (interpolated) input data).

timediv (scalar)

10

Divide input time intervals by that to get the integration time step. If timestep is set to non-zero, still an integer division is used, but one that yields a step close to the desired one (subject to rounding).

timestep, s (scalar)

0

desired time step size (see timediv)

varinit, i (array)

[]

array of initial values; must match number of derivatives from ODE

vartitle (array)

[]

array of column titles for the integrated variables

version (scalar)

0

print out the program version

AUTHOR

Thomas Orgis <thomas@orgis.org>

LICENSE AND COPYRIGHT

Copyright (c) 2005-2023 Thomas Orgis, Free Software licensed under the same terms as Perl 5.10