From e0c6872cf40896c7be36b11dcc744620f10adf1d Mon Sep 17 00:00:00 2001 From: Norbert Preining Date: Mon, 2 Sep 2019 13:46:59 +0900 Subject: Initial commit --- macros/latex/contrib/mfpic4ode/README | 16 + macros/latex/contrib/mfpic4ode/demo/demo-plain.pdf | Bin 0 -> 213804 bytes macros/latex/contrib/mfpic4ode/demo/demo-plain.tex | 368 ++++++++++++++++++ macros/latex/contrib/mfpic4ode/demo/demo.pdf | Bin 0 -> 221236 bytes macros/latex/contrib/mfpic4ode/demo/demo.tex | 381 ++++++++++++++++++ macros/latex/contrib/mfpic4ode/mfpic4ode.dtx | 429 +++++++++++++++++++++ macros/latex/contrib/mfpic4ode/mfpic4ode.ins | 59 +++ macros/latex/contrib/mfpic4ode/mfpic4ode.pdf | Bin 0 -> 231163 bytes 8 files changed, 1253 insertions(+) create mode 100644 macros/latex/contrib/mfpic4ode/README create mode 100644 macros/latex/contrib/mfpic4ode/demo/demo-plain.pdf create mode 100644 macros/latex/contrib/mfpic4ode/demo/demo-plain.tex create mode 100644 macros/latex/contrib/mfpic4ode/demo/demo.pdf create mode 100644 macros/latex/contrib/mfpic4ode/demo/demo.tex create mode 100644 macros/latex/contrib/mfpic4ode/mfpic4ode.dtx create mode 100644 macros/latex/contrib/mfpic4ode/mfpic4ode.ins create mode 100644 macros/latex/contrib/mfpic4ode/mfpic4ode.pdf (limited to 'macros/latex/contrib/mfpic4ode') diff --git a/macros/latex/contrib/mfpic4ode/README b/macros/latex/contrib/mfpic4ode/README new file mode 100644 index 0000000000..12f93d2aa0 --- /dev/null +++ b/macros/latex/contrib/mfpic4ode/README @@ -0,0 +1,16 @@ +The mfpic4ode is a set of macros for drawing direction fields, phase +portraits and trajectories of differential equations and two +dimensional autonomous systems. The Euler, Runge-Kutta and 4th order +Runge-Kutta algorithms are available to solve the ODE's. The picture +is translated into mfpic macros and MetaPost is used to create the +final drawing. The package is intended for LaTeX, but it can be used +in plainTeX as well. + +Online demonstration of the mfpic4ode macros is available on +http://wood.mendelu.cz/math/mfpicpreviewer/ as Example 6. + +You can get the documentation but compiling mfpic4ode.dtx file and the +files for installation by compiling mfpic4ode.ins by LaTeX. + +The author is supported by the grant FRVS 99/2008 (Fund for +Developement of Czech Universities). \ No newline at end of file diff --git a/macros/latex/contrib/mfpic4ode/demo/demo-plain.pdf b/macros/latex/contrib/mfpic4ode/demo/demo-plain.pdf new file mode 100644 index 0000000000..0d96e43985 Binary files /dev/null and b/macros/latex/contrib/mfpic4ode/demo/demo-plain.pdf differ diff --git a/macros/latex/contrib/mfpic4ode/demo/demo-plain.tex b/macros/latex/contrib/mfpic4ode/demo/demo-plain.tex new file mode 100644 index 0000000000..f8e758f3d6 --- /dev/null +++ b/macros/latex/contrib/mfpic4ode/demo/demo-plain.tex @@ -0,0 +1,368 @@ +% compile with pdftex demo-plain; mpost portret; pdftex demo-plain +\input multido +\input mfpic + +\usemetapost +\opengraphsfile{portret} +\tlabelsep{3pt} +\usemplabels + + +\def\frac#1#2{{#1 \over #2}} +\mfpicunit=1cm + +\input mfpic4ode + + + + + +\centerline{\bf Test file for mfpic4ode package} +\centerline{Robert Ma\v r\'\i k} +\centerline{January 3, 2008} + +\bigskip + +See the source file {\tt demo-plain.tex} for comments in the \TeX{} code. + +\clipmfpic +\bigskip +\centerline{\bf Logistic equation} +Here we draw a simple picture which describes stability of stationary +points of teh equation and then draw phase portrait of the equation. + +$$ x'= + {{r}\cdot\left(1-{x\over K}\right)x} +$$ + +% We set parameters for logistic equation +\mfsrc{r:=1;K:=0.98;} + +% We set parameters for drawing and for the ODE solver +\mfsrc{ODEarrowlength:=0.07; ODEstep:=0.02; ODEstepcount:=500;} + +% We define the equation +\ODEdefineequation{r*y*(1-(y/K))} + +Stability and sign of the right--hand side. + +\mfpic[5][3]{-0.1}{1.5}{-0.1}{0.5} + \axes + \xmarks{K} + \tlabel[tc](K,0){$K$} + \tlabel[bc](0,ypos){$\smash{f(x)}$} + \tlabel[cl](xpos,0){$x$} + + % This code draws arrows on x axis, the arrow points to the right if + % $f(x)$ is positive and to the left if f(x) is negative. If the + % starting point of the arrow is x_0 and the final point x_1 and if + % the function changes sign between x_0 and x_1, the arrow from x_0 + % is not drawn + \multido{\r=0.01+0.2}{8}{\ODEharrow{\r}} + + \pen{1pt} + \draw[rgb(0,0.5,0)]\parafcn{0,1.5,0.1}{(t,ODErhs(1,t))} + \draw[red]\parafcn{0,1.5,0.1}{(t,ODErhs(1,t))} + \draw[blue]\parafcn{0,K,0.1}{(t,ODErhs(1,t))} +\endmfpic + + +\mfsrc{ODEarrowlength:=0.3;} + + +Phase portrait + +\mfpic[1][4]{-0.1}{10}{-0.1}{1.5} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \penwd{1pt} + \tlabel[tc](xpos,0){$t$} + \tlabel[bc](0,ypos){$x$} + \pen{0.7pt} + + % This code draws arrows in 19x15 points and draws three integral + % curves + + \multido{\r=0.0+0.1}{15} + {\multido{\R=0.0+0.5}{19} + {\ODEarrow{\R}{\r}}} + \trajectories{0,0.3;0,0.01;0,1.4} + +\endmfpic + +\bigskip +\centerline{\bf Logistic equation with harvesting} + +Similar to the previous picture, but both pictures are drawn together +to see the relations between them. + +$$ x'= + {{r}\cdot\left(1-\frac x{{K}}\right)x}-{p} +$$ + +% We set parameters for logistic equation with harvesting and define +% this equation. +\mfsrc{r:=1;K:=0.98;lov:=0.15;} +\ODEdefineequation{r*y*(1-(y/K))-lov} + +% If the equation possesses stationary points, we store them into +% variables meza and mezb. If not, we set these variables to negative +% values (and the are not drawn in view of mfpicclip option. +\mfsrc{if (r**2*(K**2)-4*r*lov*K)<0: meza:=-1;mezb:=-1.1 + else: meza:=(r*K-sqrt(r**2*(K**2)-4*r*lov*K))/(2*r); + mezb:=(r*K+sqrt(r**2*(K**2)-4*r*lov*K))/(2*r) + fi;} + + +\hbox to \hsize{\hss\mfpic[1][4]{-0.1}{10}{-0.1}{1.3} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \penwd{1pt} + \tlabel[tc](xpos,0){$t$} + \tlabel[bc](0,ypos){$x$} + + % here we draw lines - stationary solutions stored in meza and mezb + % variables + \draw[gray(0.5)]\dashed\lines{(0,meza),(xpos,meza)} + \draw[gray(0.5)]\dashed\lines{(0,mezb),(xpos,mezb)} + + % We draw direction field using metapost cycle. Another option is to + % use multido command as in the previous example. + \mfsrc{for j=0 step 0.07 until 1.2: + for i:=0 step 0.5 until 10:} + \ODEarrow{i}{j} + \mfsrc{endfor;endfor;} + + % We draw trajectories using black color. + \drawcolor{black} + \trajectories{0,0.1;0,0.2;0,0.4;0,0.6;0,0.8;0,1.1} +\endmfpic\qquad +% On the right we draw the right hand side of the equation +\mfpic[3][4]{-0.15}{0.6}{-0.1}{1.3} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \tlabel[br](xpos,0){$f(x)$} + \tlabel[bc](0,ypos){$x$} + \pen{1pt} + + % This code draws the graph of right-hand side of logistic equation + % without harvesting. + \drawcolor{gray(0.7)}\draw\parafcn{0,1.2,0.1}{(ODErhs(1,t)+lov,t)} + + % This code draws the graph of right-hand side. We use the blue + % color for positivce and red color for negative parts. We draw also + % arrows which are determined by the sigh of the right-hand side of + % the equation. + \drawcolor{red}\draw\parafcn{0,1.2,0.1}{(ODErhs(1,t),t)} + \drawcolor{blue}\parafcn{meza,mezb,0.05}{(ODErhs(1,t),t)} + \multido{\r=0.1+0.15}{7}{\ODEvarrow{\r}} +\endmfpic\hss} + +\bigskip +\centerline{\bf Three numerical methods for ODEs} + +Here we draw solution of ODE using all three available numerical +methods. We use big step to see the difference between Euler, +Runge--Kutta and fourth order Runge--Kutta method. +$$ y'=x+y^3\qquad y(0)=1 +$$ + \leavevmode + \mfpic[20][5]{0}{0.5}{0.9}{2.4} + % We set-up parameters + \mfsrc{ODEarrowlength:=0.5;} + \ODEdefineequation{x+(y**3)} + \pen{1pt}\tlabelsep{1pt} + + % We set up parameters for small step + \nomplabels + \drawcolor{green} + \mfsrc{ODEstep:=0.02; ODEstepcount:=30;} + \trajectoryRKF{0}{1} + \tlabel[tr](0.39,2.4){\bf{Exact solution}} + + % We use bigger step to see the difference between various + % methods. + \mfsrc{ODEstep:=0.2; ODEstepcount:=2;} + + % We draw trajectory by Euler method + \drawcolor{black} + \trajectory{0}{1} + \tlabel[bl](0.4,1.6){\bf{Euler}} + \tlabel[tl](0,1){\bf{$k_1$}} + \tlabel[tl](0.2,1.2){\bf{$k_1$ for second step}} + + % We draw trajectory by Runge-Kutta method + \drawcolor{rgb(0.5,0.5,0.5)} + \trajectoryRK{0}{1} + \tlabel[cl](0.4,2.05){\bf{RK}} + \tlabel[tl](0.1,1.1){\bf{$k_2$}} + + % We draw trajectory by fourth order Runge-Kutta method + \drawcolor{rgb(1,0,0)} + \trajectoryRKF{0}{1} + \tlabel[tl](0.4,2.15){\bf{RK4}} + + \tlabelsep{3pt} + + % We draw direction field using blue arrows and metapost cycle + \penwd{1pt} + \drawcolor{blue}\headcolor{blue} + \mfsrc{for j=0.9 step 0.1 until 2.3: + for i:=0 step 0.05 until 0.5:} + \ODEarrow{i}{j} + \mfsrc{endfor;endfor;} + + + \drawcolor{black} + \doaxes{lbrt} + \bmarks{0,0.2,0.4} + \tmarks{0,0.2,0.4} + \lmarks{1,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,2,2.2} + \pointcolor{red} + \pointfilltrue + \point[4pt]{(0,1)} + + % We draw some other and bigger arrows in the direction field. The + % slopes of hese arrows are important for the first step by Euler + % and Ruge-Kutta method and the second step by Euler method. + \pen{2pt} + \mfsrc{ODEarrowlength:=1;} + \colorODEarrowfalse + \drawcolor{red}\headcolor{red} + \ODEarrows{0,1;0.1,1.1} + \ODEarrow{0.2}{1.2} + + \axislabels{b}[tc]{{$0$}0,{$0.2$}0.2,{$0.4$}0.4} + \axislabels{l}[cr]{{$0.8$}0.8,{$1$}1,{$1.2$}1.2,{$1.4$}1.4,{$1.6$}1.6,{$1.8$}1.8,{$2$}2,{$2.2$}2.2,{$2.4$}2.4} + \endmfpic + + +\bigskip +\centerline{\bf Autonomous systems} + +% The color arrows have no sense in the phase portrait of autonomous +% system. +\colorODEarrowfalse + +We draw the phase portrait of autonomous system, nulclines, invariant +set between nulclines, trajectories. We draw arrows in regular grid +and add few more arrows on nulclines and outside the regular grid. + +\mfsrc{TIMEstep:=0.05; TIMEend:=30;} + +\mfpic + [0.5]{-2}{15}{-2}{15} + \nomplabels + \tlabel[cc](8,15.5){\bf Competing species} + \usemplabels + + % We set up parameters, define equations, define functions which + % describe nulclines and store stationary points into variables z10, + % z11, z12, z1. + + \mfsrc{a:=11;b:=1;c:=0.8;k:=10;l:=1.1;m:=1.2;} + \mfsrc{ODEarrowlength:=0.3;} + \ASdefineequations{x*(a-b*x-c*y)}{y*(k-l*x-m*y)} + \fdef{xnulklina}{x}{(a-b*x)/c} + \fdef{ynulklina}{x}{(k-l*x)/m} + \mfsrc{z10=(0,a/c);z11=(0,k/m);z12=(a/b,0);z13=(k/l,0);} + + % Here we draw a gray polygon - invariant set fot the system. + \pen{0.3pt} + \gfill[gray(0.7)]\lclosed\lines{z10,z11,z13,z12} + \axes + \tlabel[bc](0,ypos){$y$} + \tlabel[cl](xpos,0){$x$} + + \pointsize=3pt + + \pointfilltrue\pointcolor{red}\point{(a/b,0)} + \draw[red]\function{0,a/b,1}{xnulklina(x)} + \draw[red]\lines{(0,0),(0,ypos)} + \tlabel[cr](x10,y10){$a\over c$} + \tlabel[tc](x12,y12){$a\over b$} + + \pointfilltrue\pointcolor{green}\point{(0,k/m)} + \draw[green]\function{0,k/l,1}{ynulklina(x)} + \draw[green]\lines{(0,0),(xpos,0)} + \tlabel[cr](x11,y11){$\alpha\over \beta$} + \tlabel[tc](x13,y13){$\alpha\over \gamma$} + + + \penwd{1.5pt} + \drawcolor{gray(0.25)}\headcolor{gray(0.25)} + \ASarrows{0,6;0,11;7,0;9,0;4.5,0;13,0;2,0;0,14;0,2} + \ASarrows{3,3;8,7;13,2;3,14} + \ASarrows{4,ynulklina(4);5,ynulklina(5);6,ynulklina(6); + 1.7,ynulklina(1.7)} + \ASarrows{0.5,xnulklina(0.5);1.8,xnulklina(1.8);6,xnulklina(6);4,xnulklina(4); + 7.5,xnulklina(7.5)} + + \drawcolor{black} + \AStrajectories{12,12;8,12;4,12;0.2,3;1,3;6,12;1,1;1,0.12;0.3,14} + \mfsrc{TIMEstep:=-0.05; TIMEend:=5;} + \AStrajectories{12,12;8,12;4,12;0.2,3;1,3;6,12;1,1;1,0.12;0.3,14} + \penwd{0.5pt} + \drawcolor{gray(0.5)} + \headcolor{gray(0.5)} + \multido{\r=0.5+1}{15}{\multido{\R=0.5+1}{15}{\ASarrow{\R}{\r}}} + + +\endmfpic + + +\break +\bigskip +\centerline{\bf Predator prey system with HollingII response function} + + + + +\mfpic[2]{-0.1}{4}{-0.1}{3} + % we define functions and parameters, right hand sides of the system + % and a function which defines nulcline + \mfsrc{r:=1;K:=3;a:=1;k:=0.8;P:=1;Alfa:=0.42;} + \fdef{funkceV}{x}{a*x/(x+P)} + \ASdefineequations{r*x*(1-(x/K))-funkceV(x)*y}{(-Alfa+k*funkceV(x))*y} + \fdef{xnulklina}{x}{r*(1-(x/K))*(x+P)/a} + \mfsrc{ODEarrowlength:=0.2;} + + % Here we draw axes and nulclines + \axes + \tlabel[bc](0,ypos){$y$} + \tlabel[cl](xpos,0){$x$} + \draw[red]\function{0,K,0.1}{xnulklina(x)} + \draw[green]\lines{(P/((k*a/Alfa)-1),0),(P/((k*a/Alfa)-1),ypos)} + + % Here we draw some arrows on nulclines and then arrows in the plane + \penwd{0.5pt} + \drawcolor{gray(0.25)}\headcolor{gray(0.25)} + \ASarrows{P/((k*a/Alfa)-1),1;P/((k*a/Alfa)-1),2;P/((k*a/Alfa)-1),0.5;P/((k*a/Alfa)-1),1.5} + \ASarrows{0,xnulklina(0);0.5,xnulklina(0.5);1,xnulklina(1);1.5,xnulklina(1.5);2,xnulklina(2);2.25,xnulklina(2.25)} + \multido{\r=0.1+0.25}{20}{\multido{\R=0.1+0.25}{20}{\ASarrow{\R}{\r}}} + + %% We draw trajectory with IC x=2, y=2 + \drawcolor{black} + \mfsrc{TIMEstep:=0.05; TIMEend:=10;} + \AStrajectoryRKF{2}{2} + %% We continue the trajectory (spiral) from the last point + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + + %% We continue backwards + \mfsrc{TIMEstep:=-0.05;} + \AStrajectoryRKF{2}{2} + +\endmfpic + + +\closegraphsfile +\end diff --git a/macros/latex/contrib/mfpic4ode/demo/demo.pdf b/macros/latex/contrib/mfpic4ode/demo/demo.pdf new file mode 100644 index 0000000000..256f443f0b Binary files /dev/null and b/macros/latex/contrib/mfpic4ode/demo/demo.pdf differ diff --git a/macros/latex/contrib/mfpic4ode/demo/demo.tex b/macros/latex/contrib/mfpic4ode/demo/demo.tex new file mode 100644 index 0000000000..0a56345461 --- /dev/null +++ b/macros/latex/contrib/mfpic4ode/demo/demo.tex @@ -0,0 +1,381 @@ +% compile with pdflatex demo; mpost portret; pdflatex demo +\documentclass{article} + +\usepackage{multido} +\usepackage{amsmath} +%\usepackage{color} +\usepackage%[clip] +{mfpic} +\usemetapost +\opengraphsfile{portret} +\tlabelsep{3pt} +\usemplabels + + +\mfpicunit=1cm + + + +\usepackage{mfpic4ode} + + +\begin{document} + +\title{Test file for mfpic4ode package} +\author{Robert Ma\v r\'\i k} +\date\today +\maketitle + +See the source file \texttt{demo.tex} for comments in the \TeX{} code. + +\clipmfpic +\section{Logistic equation} +Here we draw a simple picture which describes stability of stationary +points of teh equation and then draw phase portrait of the equation. + +\begin{center} + \begin{equation} + x'= + {{r}\cdot\left(1-\frac x{{K}}\right)x} + \end{equation} +\end{center} + +% We set parameters for logistic equation +\mfsrc{r:=1;K:=0.98;} + +% We set parameters for drawing and for the ODE solver +\mfsrc{ODEarrowlength:=0.07; ODEstep:=0.02; ODEstepcount:=500;} + +% We define the equation +\ODEdefineequation{r*y*(1-(y/K))} + +Stability and sign of the right--hand side. + +\begin{mfpic}[5][3]{-0.1}{1.5}{-0.1}{0.5} + \axes + \xmarks{K} + \tlabel[tc](K,0){$K$} + \tlabel[bc](0,ypos){$\smash{f(x)}$} + \tlabel[cl](xpos,0){$x$} + + % This code draws arrows on x axis, the arrow points to the right if + % $f(x)$ is positive and to the left if f(x) is negative. If the + % starting point of the arrow is x_0 and the final point x_1 and if + % the function changes sign between x_0 and x_1, the arrow from x_0 + % is not drawn + \multido{\r=0.01+0.2}{8}{\ODEharrow{\r}} + + \pen{1pt} + \draw[rgb(0,0.5,0)]\parafcn{0,1.5,0.1}{(t,ODErhs(1,t))} + \draw[red]\parafcn{0,1.5,0.1}{(t,ODErhs(1,t))} + \draw[blue]\parafcn{0,K,0.1}{(t,ODErhs(1,t))} +\end{mfpic} + + +\mfsrc{ODEarrowlength:=0.3;} + + +Phase portrait + +\begin{mfpic}[1][4]{-0.1}{10}{-0.1}{1.5} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \penwd{1pt} + \tlabel[tc](xpos,0){$t$} + \tlabel[bc](0,ypos){$x$} + \pen{0.7pt} + + % This code draws arrows in 19x15 points and draws three integral + % curves + + \multido{\r=0.0+0.1}{15} + {\multido{\R=0.0+0.5}{19} + {\ODEarrow{\R}{\r}}} + \trajectories{0,0.3;0,0.01;0,1.4} + +\end{mfpic} + + +\section{Logistic equation with harvesting} + +Similar to the previous picture, but both pictures are drawn together +to see the relations between them. + +\begin{center} + \begin{equation} + x'= + {{r}\cdot\left(1-\frac x{{K}}\right)x}-{p} + \end{equation} +\end{center} + +% We set parameters for logistic equation with harvesting and define +% this equation. +\mfsrc{r:=1;K:=0.98;lov:=0.15;} +\ODEdefineequation{r*y*(1-(y/K))-lov} + +% If the equation possesses stationary points, we store them into +% variables meza and mezb. If not, we set these variables to negative +% values (and the are not drawn in view of mfpicclip option. +\mfsrc{if (r**2*(K**2)-4*r*lov*K)<0: meza:=-1;mezb:=-1.1 + else: meza:=(r*K-sqrt(r**2*(K**2)-4*r*lov*K))/(2*r); + mezb:=(r*K+sqrt(r**2*(K**2)-4*r*lov*K))/(2*r) + fi;} + + +\hbox to \hsize{\hss\begin{mfpic}[1][4]{-0.1}{10}{-0.1}{1.3} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \penwd{1pt} + \tlabel[tc](xpos,0){$t$} + \tlabel[bc](0,ypos){$x$} + + % here we draw lines - stationary solutions stored in meza and mezb + % variables + \draw[gray(0.5)]\dashed\lines{(0,meza),(xpos,meza)} + \draw[gray(0.5)]\dashed\lines{(0,mezb),(xpos,mezb)} + + % We draw direction field using metapost cycle. Another option is to + % use multido command as in the previous example. + \mfsrc{for j=0 step 0.07 until 1.2: + for i:=0 step 0.5 until 10:} + \ODEarrow{i}{j} + \mfsrc{endfor;endfor;} + + % We draw trajectories using black color. + \drawcolor{black} + \trajectories{0,0.1;0,0.2;0,0.4;0,0.6;0,0.8;0,1.1} +\end{mfpic}\qquad +% On the right we draw the right hand side of the equation +\begin{mfpic}[3][4]{-0.15}{0.6}{-0.1}{1.3} + \axes + \ymarks{K} + \tlabel[cr](0,K){$K$} + \tlabel[br](xpos,0){$f(x)$} + \tlabel[bc](0,ypos){$x$} + \pen{1pt} + + % This code draws the graph of right-hand side of logistic equation + % without harvesting. + \drawcolor{gray(0.7)}\draw\parafcn{0,1.2,0.1}{(ODErhs(1,t)+lov,t)} + + % This code draws the graph of right-hand side. We use the blue + % color for positivce and red color for negative parts. We draw also + % arrows which are determined by the sigh of the right-hand side of + % the equation. + \drawcolor{red}\draw\parafcn{0,1.2,0.1}{(ODErhs(1,t),t)} + \drawcolor{blue}\parafcn{meza,mezb,0.05}{(ODErhs(1,t),t)} + \multido{\r=0.1+0.15}{7}{\ODEvarrow{\r}} +\end{mfpic}\hss} + +\section{Three numerical methods for ODEs} + +Here we draw solution of ODE using all three available numerical +methods. We use big step to see the difference between Euler, +Runge--Kutta and fourth order Runge--Kutta method. +\begin{center} + \begin{equation*} + \boxed{\begin{aligned}y'&=x+y^3\\y(0)&=1 + \end{aligned}}\qquad + \boxed{\begin{aligned}x_{n+1}&=x_n+h\\ + y_{n+1}&=y_n+kh\\h&=0.2 + \end{aligned}} + \end{equation*} + \leavevmode + \begin{mfpic}[20][5]{0}{0.5}{0.9}{2.4} + % We set-up parameters + \mfsrc{ODEarrowlength:=0.5;} + \ODEdefineequation{x+(y**3)} + \pen{1pt}\tlabelsep{1pt} + + % We set up parameters for small step + \nomplabels + \drawcolor{green} + \mfsrc{ODEstep:=0.02; ODEstepcount:=30;} + \dashed\trajectoryRKF{0}{1} + \tlabel[tr](0.39,2.4){\bf{Exact solution}} + + % We use bigger step to see the difference between various + % methods. + \mfsrc{ODEstep:=0.2; ODEstepcount:=2;} + + % We draw trajectory by Euler method + \drawcolor{black} + \trajectory{0}{1} + \tlabel[bl](0.4,1.6){\bf{Euler}} + \tlabel[tl](0,1){\bf{$k_1$}} + \tlabel[tl](0.2,1.2){\bf{$k_1$ for second step}} + + % We draw trajectory by Runge-Kutta method + \drawcolor{rgb(0.5,0.5,0.5)} + \trajectoryRK{0}{1} + \tlabel[cl](0.4,2.05){\bf{RK}} + \tlabel[tl](0.1,1.1){\bf{$k_2$}} + + % We draw trajectory by fourth order Runge-Kutta method + \drawcolor{rgb(1,0,0)} + \trajectoryRKF{0}{1} + \tlabel[tl](0.4,2.15){\bf{RK4}} + + \tlabelsep{3pt} + + % We draw direction field using blue arrows and metapost cycle + \penwd{1pt} + \drawcolor{blue}\headcolor{blue} + \mfsrc{for j=0.9 step 0.1 until 2.3: + for i:=0 step 0.05 until 0.5:} + \ODEarrow{i}{j} + \mfsrc{endfor;endfor;} + + + \drawcolor{black} + \doaxes{lbrt} + \bmarks{0,0.2,0.4} + \tmarks{0,0.2,0.4} + \lmarks{1,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,2,2.2} + \pointcolor{red} + \pointfilltrue + \point[4pt]{(0,1)} + + % We draw some other and bigger arrows in the direction field. The + % slopes of hese arrows are important for the first step by Euler + % and Ruge-Kutta method and the second step by Euler method. + \pen{2pt} + \mfsrc{ODEarrowlength:=1;} + \colorODEarrowfalse + \drawcolor{red}\headcolor{red} + \ODEarrows{0,1;0.1,1.1} + \ODEarrow{0.2}{1.2} + + \axislabels{b}[tc]{{$0$}0,{$0.2$}0.2,{$0.4$}0.4} + \axislabels{l}[cr]{{$0.8$}0.8,{$1$}1,{$1.2$}1.2,{$1.4$}1.4,{$1.6$}1.6,{$1.8$}1.8,{$2$}2,{$2.2$}2.2,{$2.4$}2.4} + \end{mfpic} + +\end{center} + + +\section{Autonomous systems} + +% The color arrows have no sense in the phase portrait of autonomous +% system. +\colorODEarrowfalse + +We draw the phase portrait of autonomous system, nulclines, invariant +set between nulclines, trajectories. We draw arrows in regular grid +and add few more arrows on nulclines and outside the regular grid. + +\mfsrc{TIMEstep:=0.05; TIMEend:=30;} + +\begin{mfpic} + [0.5]{-2}{15}{-2}{15} + \nomplabels + \tlabel[cc](8,15.5){\bf Competing species} + \usemplabels + + % We set up parameters, define equations, define functions which + % describe nulclines and store stationary points into variables z10, + % z11, z12, z1. + + \mfsrc{a:=11;b:=1;c:=0.8;k:=10;l:=1.1;m:=1.2;} + \mfsrc{ODEarrowlength:=0.3;} + \ASdefineequations{x*(a-b*x-c*y)}{y*(k-l*x-m*y)} + \fdef{xnulklina}{x}{(a-b*x)/c} + \fdef{ynulklina}{x}{(k-l*x)/m} + \mfsrc{z10=(0,a/c);z11=(0,k/m);z12=(a/b,0);z13=(k/l,0);} + + % Here we draw a gray polygon - invariant set fot the system. + \pen{0.3pt} + \gfill[gray(0.7)]\lclosed\lines{z10,z11,z13,z12} + \axes + \tlabel[bc](0,ypos){$y$} + \tlabel[cl](xpos,0){$x$} + + \pointsize=3pt + + \pointfilltrue\pointcolor{red}\point{(a/b,0)} + \draw[red]\function{0,a/b,1}{xnulklina(x)} + \draw[red]\lines{(0,0),(0,ypos)} + \tlabel[cr](x10,y10){$a\over c$} + \tlabel[tc](x12,y12){$a\over b$} + + \pointfilltrue\pointcolor{green}\point{(0,k/m)} + \draw[green]\function{0,k/l,1}{ynulklina(x)} + \draw[green]\lines{(0,0),(xpos,0)} + \tlabel[cr](x11,y11){$\alpha\over \beta$} + \tlabel[tc](x13,y13){$\alpha\over \gamma$} + + + \penwd{1.5pt} + \drawcolor{gray(0.25)}\headcolor{gray(0.25)} + \ASarrows{0,6;0,11;7,0;9,0;4.5,0;13,0;2,0;0,14;0,2} + \ASarrows{3,3;8,7;13,2;3,14} + \ASarrows{4,ynulklina(4);5,ynulklina(5);6,ynulklina(6); + 1.7,ynulklina(1.7)} + \ASarrows{0.5,xnulklina(0.5);1.8,xnulklina(1.8);6,xnulklina(6);4,xnulklina(4); + 7.5,xnulklina(7.5)} + + \drawcolor{black} + \AStrajectories{12,12;8,12;4,12;0.2,3;1,3;6,12;1,1;1,0.12;0.3,14} + \mfsrc{TIMEstep:=-0.05; TIMEend:=5;} + \AStrajectories{12,12;8,12;4,12;0.2,3;1,3;6,12;1,1;1,0.12;0.3,14} + \penwd{0.5pt} + \drawcolor{gray(0.5)} + \headcolor{gray(0.5)} + \multido{\r=0.5+1}{15}{\multido{\R=0.5+1}{15}{\ASarrow{\R}{\r}}} + + +\end{mfpic} + + +\begin{center} + \textbf{Pedator prey system with HollingII response function} +\end{center} + + + +\begin{mfpic}[2]{-0.1}{4}{-0.1}{3} + % we define functions and parameters, right hand sides of the system + % and a function which defines nulcline + \mfsrc{r:=1;K:=3;a:=1;k:=0.8;P:=1;Alfa:=0.42;} + \fdef{funkceV}{x}{a*x/(x+P)} + \ASdefineequations{r*x*(1-(x/K))-funkceV(x)*y}{(-Alfa+k*funkceV(x))*y} + \fdef{xnulklina}{x}{r*(1-(x/K))*(x+P)/a} + \mfsrc{ODEarrowlength:=0.2;} + + % Here we draw axes and nulclines + \axes + \tlabel[bc](0,ypos){$y$} + \tlabel[cl](xpos,0){$x$} + \draw[red]\function{0,K,0.1}{xnulklina(x)} + \draw[green]\lines{(P/((k*a/Alfa)-1),0),(P/((k*a/Alfa)-1),ypos)} + + % Here we draw some arrows on nulclines and then arrows in the plane + \penwd{0.5pt} + \drawcolor{gray(0.25)}\headcolor{gray(0.25)} + \ASarrows{P/((k*a/Alfa)-1),1;P/((k*a/Alfa)-1),2;P/((k*a/Alfa)-1),0.5;P/((k*a/Alfa)-1),1.5} + \ASarrows{0,xnulklina(0);0.5,xnulklina(0.5);1,xnulklina(1);1.5,xnulklina(1.5);2,xnulklina(2);2.25,xnulklina(2.25)} + \multido{\r=0.1+0.25}{20}{\multido{\R=0.1+0.25}{20}{\ASarrow{\R}{\r}}} + + %% We draw trajectory with IC x=2, y=2 + \drawcolor{black} + \mfsrc{TIMEstep:=0.05; TIMEend:=10;} + \AStrajectoryRKF{2}{2} + %% We continue the trajectory (spiral) from the last point + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + \AStrajectoryRKF{x1}{y1} + + %% We continue backwards + \mfsrc{TIMEstep:=-0.05;} + \AStrajectoryRKF{2}{2} + +\end{mfpic} + + + +\end{document} diff --git a/macros/latex/contrib/mfpic4ode/mfpic4ode.dtx b/macros/latex/contrib/mfpic4ode/mfpic4ode.dtx new file mode 100644 index 0000000000..f8e337b3df --- /dev/null +++ b/macros/latex/contrib/mfpic4ode/mfpic4ode.dtx @@ -0,0 +1,429 @@ +% \iffalse meta-comment +% +% Copyright (C) 2006-2008 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. +% +% \fi +% +% \iffalse +%<*driver> +\ProvidesFile{mfpic4ode.dtx} +% +%\NeedsTeXFormat{LaTeX2e}[1999/12/01] +%\ProvidesPackage{mfpic4ode} +%<*sty> + [2009/04/15 v0.3 mfpic4ode.dtx file] +% +% +%<*driver> +\documentclass{ltxdoc} +\EnableCrossrefs +\CodelineIndex +\RecordChanges +\begin{document} + \DocInput{mfpic4ode.dtx} + \PrintChanges + \PrintIndex +\end{document} +% +% \fi +% +% \CheckSum{185} +% +% \CharacterTable +% {Upper-case \A\B\C\D\E\F\G\H\I\J\K\L\M\N\O\P\Q\R\S\T\U\V\W\X\Y\Z +% Lower-case \a\b\c\d\e\f\g\h\i\j\k\l\m\n\o\p\q\r\s\t\u\v\w\x\y\z +% Digits \0\1\2\3\4\5\6\7\8\9 +% Exclamation \! Double quote \" Hash (number) \# +% Dollar \$ Percent \% Ampersand \& +% Acute accent \' Left paren \( Right paren \) +% Asterisk \* Plus \+ Comma \, +% Minus \- Point \. Solidus \/ +% Colon \: Semicolon \; Less than \< +% Equals \= Greater than \> Question mark \? +% Commercial at \@ Left bracket \[ Backslash \\ +% Right bracket \] Circumflex \^ Underscore \_ +% Grave accent \` Left brace \{ Vertical bar \| +% Right brace \} Tilde \~} +% +% +% \changes{v0.2}{2008/02/03}{First public version} +% \changes{v0.3}{2009/04/15}{Added connect environment to paths} +% \changes{v0.4}{2010/04/07}{Updated documentations} +% +% \GetFileInfo{mfpic4ode.dtx} +% +% \DoNotIndex{\newcommand,\newenvironment} +% +% \def\fileversion{0.4} +% \def\filedate{2010/04/07} +% +% \title{The \textsf{mfpic4ode} package\thanks{This document +% corresponds to \textsf{mfpic4ode}~\fileversion, dated +% \filedate.}} \author{Robert Marik \\ \texttt{marik@mendelu.cz} +% \thanks{Supported by the grant FRV\v{S} 99/2008 (Fund for +% Developement of Czech Universities).}} +% \maketitle +% +% \section{Introduction} +% +% The package |mfpic4ode| is a set of macros for drawing phase +% portraits and integral curves of differential equations and +% autonomous systems using |mfpic| macros. These macros have been used +% by the author to prepare some pictures for classrooms and the +% results seem to be acceptable for this purpose, but always remember +% that due to the fixed points arithmetics in Metapost, the error in +% computations could be significant. Another excellent tool which can +% be used to draw trajectories is Sage\TeX{} which gives you full +% power of computer algebra system Sage in \LaTeX. +% +% \section{Usage} +% You can load the package in \LaTeX{} using standard +% |\usepackage{mfpic4ode}| command, or you can use the macros in +% plain\TeX{} and load by |\input mfpic4ode.tex| command. +% +% \subsection{First order differential equation} +% +% To draw phase portrait of first order ordinary differential equation +% $$ +% y'=f(x,y) +% $$ +% we define commands |\ODEarrow| for drawing element of direction +% field and |\trajectory|, |\trajectoryRK| and |\trajectoryRKF| for +% drawing integral curves using Euler, second order Runge-Kutta and +% fourth order Runge-Kutta methods, respectively. Some important +% parameters, such as the number of steps, the length of step or the +% function from the right-hand side of the equations are stored in +% MetaPost variables and to keep the package simple and short, these +% variables are accessible using |\mfsrc| command. +% +% \DescribeMacro{\ifcolorODEarrow} \DescribeMacro{\colorODEarrowtrue} +% \DescribeMacro{\colorODEarrowfalse} If the \TeX{} boolean variable +% |\ifcolorODEarrow| is true, then the arrows from direction field are +% blue if the solution is increasing and red if decreasing. If +% |\ifcolorODEarrow| is false, the mfpic color from |\drawcolor| and +% |\headcolor| macros is used. More precisely, +% \begin{itemize} +% \item if |\ifcolorODEarrow| is true and $f(x_0,y_0)>0$, then the +% arrow at the point $(x_0,y_0)$ is blue +% \item if |\ifcolorODEarrow| is true and $f(x_0,y_0)\leq 0$, then the +% arrow at the point $(x_0,y_0)$ is red +% \item if |\ifcolorODEarrow| is false, then color from |\drawcolor| +% is used to draw the body of an arrow and color |\headcolor| is +% used to draw the head. +% \end{itemize} +% Arrows are drawn using mfpic |\draw\arrow\lines{...}| command and +% hence the parameters for customizing shape and size of the head from +% mfpic are also available. The MetaPost variable |ODEarrowlength| is +% used to customize the length of each arrow. If the arrow is +% horizontal, then the length of the arrow in mfpic coordinates is +% equal to |ODEarrowlength/xscale|. (This fixes the case when +% different |xscale| and |yscale| are used. All arrows have the same +% length.) You can set this variable using |\mfsrc| command, you can +% write e.g. +% \begin{verbatim} +% \mfsrc{ODEarrowlength:=0.07;} +% \end{verbatim} +% in your document. +% +% To draw arrows in regular rectangular grid you should use the +% |\ODEarrow| macro in a double cycle such as +% \begin{verbatim} +% \mfsrc{for j=0 step 0.07 until 1.2: for i:=0 step 0.5 until 10:} +% \ODEarrow{i}{j} +% \mfsrc{endfor;endfor;} +% \end{verbatim} +% or using the |multido| package +% \begin{verbatim} +% \multido{\r=0.0+0.1}{15}{\multido{\R=0.0+0.5}{19}{\ODEarrow{\R}{\r}}} +% \end{verbatim} +% +% \DescribeMacro{\ODEdefineequation{f(x,y)}} The macro |\ODEdefineequation| is +% used to save the right hand side of the ODE, i. e. the function +% $f(x,y)$. You should write the expression in the MetaPost format, +% the independent variable is supposed to be $x$, the dependent +% variable is $y$. +% +% \DescribeMacro{\trajectory\{x0\}\{y0\}} +% \DescribeMacro{\trajectoryRK\{x0\}\{y0\}} +% \DescribeMacro{\trajectoryRKF\{x0\}\{y0\}} The macros |\trajectory|, +% |\trajectoryRK| and |\trajectoryRKF| are used for drawing integral +% curves with initial condition $y(x_0)=y_0$ using Euler, second order +% Runge-Kutta and fourth order Runge-Kutta methods, respectively. +% The length of each step is stored in MetaPost variable +% |ODEstepcount|, the length of each step is in the MetaPost variable +% |ODEstep|. You can set these variables using |\mfsrc| macro as +% follows +% \begin{verbatim} +% \mfsrc{ODEstep:=0.02; ODEstepcount:=500;} +% \end{verbatim} +% The integral curve is drawn from short linear parts using |\ODEline| +% command which expands to |\lines| command from |mfpic| package by +% default. These linear parts are connected in connect environment and +% this allows to use prefixes like |\dotted| or |\dashed| to the +% trajectories. A simple test is used to keep the arithmetics in +% reasonable bounds: if after the step the curve leaves the horizontal +% strip between |yneg| and |ypos| variables, then the evaluation is +% stopped (in fact, in this case we do not change the independent +% variable and we do the remaining steps with the same last point). +% Recall that |yneg| and |ypos| variables are set when you call +% |mfpicture| environment. If you call the environment as follows +% \begin{verbatim} +% \begin{mfpic}[5][3]{-0.1}{1.5}{-0.1}{0.5} +% ........... +% \end{mfpic} +% \end{verbatim} +% then no more than one short linear part of the integral curve is +% outside the horizontal strip between $y=-0.1$ and $y=0.5$. +% +% \DescribeMacro{\trajectories}\DescribeMacro{\ODEarrows} To draw more +% trajectories you can use |\trajectories| command. The command +% |\trajectories{x1,y1;x2,y2;x3,y3;....;xn,yn}| expands to $n$ +% |\trajectoryRKF| commands with initial conditions $y(x_i)=y_i$ for +% $i=1..n$. In a similar way, +% |\ODEarrows{x1,y1;x2,y2;x3,y3;....;xn,yn}| expands into $n$ +% |\ODEarrow| commands. +% +% \subsection{Two-dimensional autonomous systems} +% Trajectories for two-dimensional autonomous system +% \begin{eqnarray*} +% x'&=&f(x,y)\\y'&=&g(x,y) +% \end{eqnarray*} +% are drawn using a very simple method based on the direction field. +% This could be improved in the next release of the package, but till +% now the results obtained in this way are qualitatively correct and +% sufficiently accurate (remember that you cannot expect accurate +% approximation due to the limitation of arithmetics in MetaPost). +% Anyway, some users may prefer the fourth order Runge--Kutta method. +% +% The macros |\ASdefineequation| |\ASarrow|, |\AStrajectory|, +% |\AStrajectoryRKF|, |\ASarrows| and |\AStrajectories| are +% counterparts to their |\ODE....| versions. The last point of each +% trajectory is stored in the |x1| and |x2| MetaPost variables. Hence, +% you can say |\AStrajectory{2}{2}| to draw trajectory with initial +% conditions $x(0)=2$, $y(0)=2$ and then you can continue this +% trajectory using |\AStrajectory{x1}{y1}| command. The macro +% |\AStrajectory| uses |ODEstep| and |ODEstepcount| variables, the +% macro |\AStrajectoryRKF| uses |TIMEstep| and |TIMEend| variables do +% perform the steps in the numerical solution. The number of steps is +% in the latter case evaluated as absolute value of the quotient +% |TIMEend/TIMEstep|. You can use negative value for |TIMEstep| to +% continue the trajectory backwards. +% +% \section{Troubleshooting} +% \subsection{The catcode of @ is messed} We set the category of @ to +% 11 (letter) when we load the package and at the end of definitions +% for mfpic4ode we set the category to 12. This could be a source of +% rare problems, if you use different value in your document. +% +% \subsection{Metapost: Not implemented: (unknown numeric) \dots } +% You have to set |ODEstep|, |ODEstepcount|, |TIMEstep| and |TIMEend| +% other variables using |\mfsrc| command (depending on the type of the +% problem). +% +% \StopEventually{} +% +% \section{Implementation} +% \begin{macrocode} +%<*tex> +\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 +% +%\input mfpic4ode.tex\relax +% \end{macrocode} +% +% \Finale +\endinput + + + diff --git a/macros/latex/contrib/mfpic4ode/mfpic4ode.ins b/macros/latex/contrib/mfpic4ode/mfpic4ode.ins new file mode 100644 index 0000000000..a1fe49991f --- /dev/null +++ b/macros/latex/contrib/mfpic4ode/mfpic4ode.ins @@ -0,0 +1,59 @@ +%% +%% 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. +%% + +\input docstrip.tex +\keepsilent + +\usedir{tex/latex/mfpic4ode} + +\preamble + +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. + +\endpreamble + +\generate{\file{mfpic4ode.tex}{\from{mfpic4ode.dtx}{tex}} +} +\generate{\file{mfpic4ode.sty}{\from{mfpic4ode.dtx}{sty}} +} + +\obeyspaces +\Msg{*************************************************************} +\Msg{* *} +\Msg{* To finish the installation you have to move the following *} +\Msg{* files into a directory searched by TeX: *} +\Msg{* *} +\Msg{* mfpic4ode.tex *} +\Msg{* mfpic4ode.sty *} +\Msg{* *} +\Msg{* To produce the documentation run the file *} +\Msg{* mfpic4ode.dtx through LaTeX. *} +\Msg{* *} +\Msg{* Happy TeXing! *} +\Msg{* *} +\Msg{*************************************************************} + +\endbatchfile diff --git a/macros/latex/contrib/mfpic4ode/mfpic4ode.pdf b/macros/latex/contrib/mfpic4ode/mfpic4ode.pdf new file mode 100644 index 0000000000..ff9a173961 Binary files /dev/null and b/macros/latex/contrib/mfpic4ode/mfpic4ode.pdf differ -- cgit v1.2.3