%% %% This is file `mfpic4ode.tex', %% generated with the docstrip utility. %% %% The original source files were: %% %% mfpic4ode.dtx (with options: `tex') %% %% This is a generated file. %% %% Copyright (C) 2007 by Robert Marik %% %% This file may be distributed and/or modified under the conditions of %% the LaTeX Project Public License, either version 1.2 of this license %% or (at your option) any later version. The latest version of this %% license is in: %% %% http://www.latex-project.org/lppl.txt %% %% and version 1.2 or later is part of all distributions of LaTeX version %% 1999/12/01 or later. %% \catcode`\@=11 \newif\ifcolorODEarrow %%%\colorODEarrowfalse \colorODEarrowtrue %%% The line from one point to another \def\ODEline#1#2{\lines{#1,#2}} %%% The variable ODErhs is used to store the function from the right %%% hand side of ODE in the form y'=f(x,y). We use command %%% ODEdefineequation to set up this variable. \def\ODEdefineequation#1{\fdef{ODErhs}{x,y}{#1}} %%% Integral curve using Euler method. The step of this method is %%% ODEstep and the number of steps is ODEstepcount. The points are %%% stored in metapost variables x1,y1. \def\trajectory#1#2{ \begin{connect} \mfsrc{x1:=#1;y1:=#2; for i=1 upto ODEstepcount: x2:=x1+ODEstep; y2:=y1+ODEstep*ODErhs(x1,y1);} \ODEline{z1}{z2} \mfsrc{ if ((y2>yneg) and (y2yneg) and (y2yneg) and (y2y1:ODEcolorarrow:=blue else: ODEcolorarrow:=red fi; } \ifcolorODEarrow \drawcolor{ODEcolorarrow} \headcolor{ODEcolorarrow} \fi \draw\arrow\lines{z1,z3} } \def\ODEarrows#1{\ODE@cycle@points#1;,;} \def\trajectories#1{\ODE@cycle@IC#1;,;} \def\ODE@last@point{} \def\ODE@cycle@points#1,#2;{\def\temp{#1}\ifx\temp\ODE@last@point\let\next\relax \else\ODEarrow{#1}{#2}\relax\let\next\ODE@cycle@points\fi\next} \def\ODE@cycle@IC#1,#2;{\def\temp{#1}\ifx\temp\ODE@last@point\let\next\relax \else \trajectoryRKF{#1}{#2}\relax\let\next\ODE@cycle@IC\fi\next} \mfsrc{path p,q;color ODEcolorarrow;} %%% One-dimensional autonomous systems y'=f(y) where '=d/dx \def\ODEharrow#1{ \mfsrc{x1:=#1; if ODErhs(0,x1)>0: x3:=x1+ODEarrowlength else: x3:=x1-ODEarrowlength fi; if ODErhs(0,x1)*ODErhs(0,x3)<0: x1:=-100;x3:=-100 fi; if x3>x1:ODEcolorarrow:=blue else: ODEcolorarrow:=red fi; } \ifcolorODEarrow \drawcolor{ODEcolorarrow} \headcolor{ODEcolorarrow} \fi \pen{1.5pt} \draw\arrow\lines{(x1,0),(x3,0)} } \def\ODEvarrow#1{ \mfsrc{x1:=#1; if ODErhs(0,#1)>0: x3:=x1+(ODEarrowlength/yscale) else: x3:=x1-(ODEarrowlength/yscale) fi; if ODErhs(0,x1)*ODErhs(0,x3)<0: x1:=-100;x3:=-100 fi; if x3>x1:ODEcolorarrow:=blue else: ODEcolorarrow:=red fi; } \ifcolorODEarrow \drawcolor{ODEcolorarrow} \headcolor{ODEcolorarrow} \fi \pen{1.5pt} \draw\arrow\lines{(0,x1),(0,x3)} } %%% Two-dimensional autonomous systems x'=f(x,y), y'=g(x,y) where '=d/dt \def\ASdefineequations#1#2{\fdef{ASf}{x,y}{#1}\fdef{ASg}{x,y}{#2}} \def\AStrajectory#1#2{ \begin{connect} \mfsrc{x1:=#1;y1:=#2; for i=1 upto ODEstepcount: x2:=x1+ODEstep*ASf(x1,y1); y2:=y1+ODEstep*ASg(x1,y1);} \ODEline{z1}{z2} \mfsrc{ if ((y2>yneg) and (y2y1:ODEcolorarrow:=blue else: ODEcolorarrow:=red fi; } \ifcolorODEarrow \drawcolor{ODEcolorarrow} \headcolor{ODEcolorarrow} \fi \draw\arrow\lines{z1,z3} } \def\ASarrows#1{\AS@cycle@points#1;,;} \def\AS@cycle@points#1,#2;{\def\temp{#1}\ifx\temp\ODE@last@point\let\next\relax \else\ASarrow{#1}{#2}\relax\let\next\AS@cycle@points\fi\next} \def\AStrajectories#1{\AS@cycle@IC#1;,;} \def\AS@cycle@IC#1,#2;{\def\temp{#1}\ifx\temp\ODE@last@point\let\next\relax \else \AStrajectoryRKF{#1}{#2}\relax\let\next\AS@cycle@IC\fi\next} \def\AStrajectoryRKF#1#2{ \begin{connect} \mfsrc{x1:=#1;y1:=#2; TIMEsteps:=abs(TIMEend/TIMEstep); TIME:=0; for i=1 upto TIMEsteps: k1:=ASf(x1,y1); l1:=ASg(x1,y1); k2:=ASf(x1+(TIMEstep*k1/2),y1+(TIMEstep*l1/2)); l2:=ASg(x1+(TIMEstep*k1/2),y1+(TIMEstep*l1/2)); k3:=ASf(x1+(TIMEstep*k2/2),y1+(TIMEstep*l2/2)); l3:=ASg(x1+(TIMEstep*k2/2),y1+(TIMEstep*l2/2)); k4:=ASf(x1+(TIMEstep*k3),y1+(TIMEstep*l3)); l4:=ASg(x1+(TIMEstep*k3),y1+(TIMEstep*l3)); k5:=((k1)/6)+((k2)/3)+((k3)/3)+((k4)/6); l5:=(l1/6)+(l2/3)+(l3/3)+(l4/6); x2:=x1+(TIMEstep*k5); y2:=y1+(TIMEstep*l5);} \ODEline{z1}{z2} \mfsrc{ if ((y2>yneg) and (y2xneg)): x1:=x2; y1:=y2 fi; endfor } \end{connect} } \catcode`\@12\relax \endinput %% %% End of file `mfpic4ode.tex'.