[pstricks] Searching for a numeric solution for coupled differential equation systems of order 2

[Juergen Gilg]

Dear Alexander and Manfred,

>
> What follows is the (untested) PS code for the set of ODEs to be
> integrated with \odesolve from the animate package documentation. Anything
> marked as ?.*?  must be replaced by meaningful values. Note that angles
> are measured in degrees. ps2pdf must be invoked as ps2pdf -dNOSAFER,
> because the result is written into a file.
>
> %define m_1, m_2, l_1, l_2, g
> \pstVerb{
> 	/m1 ?? def
> 	/m2 ?? def
> 	/l1 ?? def
> 	/l2 ?? def
> 	/g  ?? def
> }
>
> %implement ODE set in RP notation as a PS function
> \pstVerb{
>    /setofodes {
>      %pop elements of current solution vector from operand stack, note the
>      %reverse order of popping
>      /phi22 exch def  %\dot\varphi_2
>      /phi21 exch def  %\varphi2
>      /phi12 exch def  %\dot\varphi_1
>      /phi11 exch def  %\varphi1
>
>      %auxiliary variables a_1, b_1, c_1, a_2, b_2, c_2
>      /a1 m2 m1 m2 add div l2 mul phi22 phi22 mul mul phi11 phi21 sub sin mul def
>      /b1 g phi11 sin mul def
>      /c1 m2 m1 m2 add div l2 mul phi11 phi21 sub cos mul def
>
>      /a2 l1 phi12 mul phi11 phi21 sub sin mul def
>      /b2 g phi21 sin mul def
>      /c2 l1 phi11 phi21 sub cos mul def
>
>      %compute right hand side of odes (4 values), leave them on the
> 	%operand stack
>      phi12 %\dot\varphi_1
>
>      %\ddot\varphi_1
>      a1 b1 sub l1 div c1 l1 div a2 b2 sub l2 div mul sub 1 c1 c2 mul l1 l2 mul
> 	div sub
>
>      phi22 %\dot\varphi_2
>
>      %\ddot\varphi_2
>      a2 b2 sub l2 div c2 l2 div a1 b1 sub l1 div mul sub 1 c1 c2 mul l1 l2 mul
> 	div sub
>    } def
> }
>
> \odesolve{result.dat}{% output file name
>    %output format: integr. param. `t'&  complete solution vector (4 elements)
>    %on one table row
>    (t) 0 1 2 3
> }{?ta?}{?tb?}{% start and end of integration parameter `t'
>    ?number of output points?
> }{
>    % initial values \varphi1, \dot\varphi_1, \varphi2, \dot\varphi_2 at t=0
>    ?? ?? ?? ??
> }{setofodes} % name of the PS function for right hand side of ODE set
>
many thanks for your help and your hints.

Regards,

Jürgen

-- 
Jürgen Gilg
Austr. 59
70376 Stuttgart

tel   0711 59 27 88
email gilg at acrotex.net