mfpic 1.06 (21mar11)

git-svn-id: svn://tug.org/texlive/trunk@21800 c570f23f-e606-0410-a88d-b1316a301751

1 files changed, 75 insertions, 21 deletions

diff --git a/Master/texmf-dist/metafont/mfpic/grafbase.mf b/Master/texmf-dist/metafont/mfpic/grafbase.mf
index 10f23e9688f..87492e00849 100644
--- a/Master/texmf-dist/metafont/mfpic/grafbase.mf
+++ b/Master/texmf-dist/metafont/mfpic/grafbase.mf
@@ -8,7 +8,7 @@
 %% 
 %% -------------------------------------------------------------------
 %% 
-%% Copyright 2002--2010, Daniel H. Luecking
+%% Copyright 2002--2011, Daniel H. Luecking
 %% 
 %% Mfpic may be distributed and/or modified under the conditions of the
 %% LaTeX Project Public License, either version 1.3c of this license or (at
@@ -21,9 +21,6 @@
 %% is Daniel H. Luecking. There are several Base Interpreters associated
 %% with mfpic: plain TeX, LaTeX, plain Metafont and plain MetaPost.
 %% 
-%%% File: grafbase.dtx
-%%% A part of mfpic 1.05 2010/06/10
-%%%
 if (known grafbaseversion) or (known grafbase):
   message "Grafbase (" & jobname & "): You have loaded grafbase more "
     & "than once! Please make sure that it is loaded only once.";
@@ -32,7 +29,7 @@ fi
 boolean grafbase; grafbase := true;
 
 string fileversion, filedate;
-fileversion := "1.05";  filedate := "2010/06/10";
+fileversion := "1.06";  filedate := "2011/02/25";
 
 message " Loading grafbase macros, version " & fileversion & ", " &
   filedate & ".";
@@ -60,7 +57,7 @@ def checkversions (expr g)=
   fi
 enddef;
 
-checkversions (105);
+checkversions (106);
 
 if unknown base_name       :  input plain;
 elseif not string base_name:  input plain;
@@ -1642,7 +1639,7 @@ enddef;
 vardef gettime (suffix arr, ct) (expr lngth) =
   setnumeric (_gtl) emax (arr[ct], emin (arr[arr], lngth));
   setsplit (_s) segment_split;
-  forever:  exitif ( (arr[ct] <= _gtl) and (_gtl <= arr[ct+1]) );
+  forever:  exitif inrange (arr[ct], arr[ct+1]) (_gtl);
     next ct;
   endfor
   if arr[ct] = arr[ct+1]:  ct
@@ -1869,11 +1866,17 @@ vardef borderrect =
   rect((xlow,ylow),(xhigh,yhigh))
 enddef;
 
+vardef between (expr A, B, X) = (A <  X) and (X <  B)  enddef;
+vardef inrange (expr A, B, X) = (A <= X) and (X <= B)  enddef;
+
 vardef inbounds (expr Z) =
-  (xpart Z >= xlow ) and (ypart Z >= ylow ) and
-  (xpart Z <= xhigh) and (ypart Z <= yhigh)
+  inrange (xlow, xhigh) (xpart Z) and inrange (ylow, yhigh) (ypart Z)
 enddef;
 
+tertiarydef X isbetween P =  between (xpart P, ypart P, X)  enddef;
+tertiarydef X isinrange P =  inrange (xpart P, ypart P, X)  enddef;
+tertiarydef P contains X  =  between (xpart P, ypart P, X)  enddef;
+
 numeric inside, outside, centered, onleft, onright, ontop, onbottom;
 inside   := -2;
 outside  := -1;
@@ -1955,7 +1958,7 @@ def plrvectorfield (expr len, rsp, tsp) (text fcn) (text cond) =
   vardef _vf (expr r,t) =  ((0,0)--(fcn)) shifted (r*dir t)  enddef;
   vardef _is_OK (expr r,t) =
     save _X, _Y;  _X := r*cosd t;  _Y := r*sind t;
-    (cond) and (_A < _X) and (_X < _B) and (_C < _Y) and (_Y < _D)
+    (cond) and between (_A, _B) (_X) and between (_C, _D) (_Y)
   enddef;
   mkplrvectorfield (len, rsp, tsp) (_vf, _is_OK);
 enddef;
@@ -2053,7 +2056,7 @@ def getpolarbounds =
       if rmax < r[j]:  rmax := r[j];  fi
     endfor
     rmin := 0;
-    if (xneg < 0) and (xpos > 0) and (yneg < 0) and (ypos > 0):
+    if between (xneg, xpos) (0) and between (yneg, ypos) (0):
       tmin := 0;  tmax := 360;
     elseif (p0 = origin):  tmin :=    0;  tmax :=  90;
     elseif (p1 = origin):  tmin :=  -90;  tmax :=   0;
@@ -2066,9 +2069,9 @@ def getpolarbounds =
         if tmax < t:  tmax := t;  fi
         if tmin > t:  tmin := t;  fi
       endfor
-      if (xneg < 0) and (0 < xpos):
+      if between (xneg, xpos) (0):
         rmin := emin (abs(yneg), abs(ypos));
-      elseif (yneg < 0) and (0 < ypos):
+      elseif between (yneg, ypos) (0):
         rmin := emin (abs(xneg), abs(xpos));
       else:
         rmin := min (r0, r1, r2, r3);
@@ -2147,8 +2150,10 @@ vardef anglefromto (expr u, v) =
 enddef;
 
 vardef cornerangle (expr A, B, C) =
-  if (A = B) and (B = C) :  60
-  elseif (A=B) or (A=C) :  90
+  if (A = B) or (A = C) :
+    if (B = C) :  60
+    else:         90
+    fi
   else:  anglefromto (B - A, C - A)
   fi
 enddef;
@@ -2704,14 +2709,14 @@ vardef pshcircle (expr disk, ctr, rad) =
   if disk:
     if rad >= 1 :
       if rad > 1:
-        GBerrmsg ("Impossible pseudohyperbolic circle.")
+        GBerrmsg ("Impossible radius of pseudohyperbolic circle.")
           "The radius of a pseudohyperbolic circle can be at most 1.";
       fi
       circle ((0,0),1)
     elseif abs(ctr) >= 1 :
       if abs(ctr) > 1:
-        GBerrmsg ("Impossible pseudohyperbolic circle.")
-          "The center of a pseudohyperbolic circle must be in"
+        GBerrmsg ("Impossible center of pseudohyperbolic circle.")
+          "The center of a pseudohyperbolic circle must be in "
           & "the unit disk.";
       fi
       onepointpath (true,ctr)
@@ -2731,7 +2736,7 @@ vardef pshcircle (expr disk, ctr, rad) =
     elseif ypart ctr <= 0:
       if ypart ctr < 0:
         GBerrmsg ("Impossible pseudohyperbolic circle.")
-          "The center of a pseudohyperbolic circle must be in"
+          "The center of a pseudohyperbolic circle must be in "
           & "the upper half-plane.";
       fi
       onepointpath (true,ctr)
@@ -2829,8 +2834,7 @@ tolerancefactor := .02;
 vardef mklevelset (expr sm, tens, X, Y, t, a, b, c, d) =
   save _inside_;
   vardef _inside_ (expr U, V) =
-    inside_levelset (U, V) and (a < U) and (U < b)
-      and (c < V) and (V < d)
+    inside_levelset (U, V) and between (a, b) (U) and between (c, d) (V)
   enddef;
   if not _inside_ (X, Y):
     GBwarn "Invalid seed point for levelset.";
@@ -2881,6 +2885,56 @@ vardef tlevelset (expr smth, tens, seed, seg) (text cond) =
     (xneg, xpos, yneg, ypos)
 enddef;
 
+def RKIV (expr sm) =  tRKIV (sm, default_tension)
+enddef;
+vardef tRKIV (expr sm, tens, zstart, ds, N) (text _RHS_) =
+  save _trj, _ztr, _dz, _ztmp, _ctm;
+  pair _trj[],          % The trajectory
+       _ztr,            % current point
+       _dz[],           % array[4] of displacements
+       _ztmp;           % current point for calculating velocity
+  _trj := N+1;          % ultimate size of _trj array
+  _trj1 := _ztr := zstart;
+  save _tt,             % current time
+       _dt,             % current time step
+       _th;             % current time plus half a step
+  _tt := 0;
+  for _idx := 2 upto _trj:
+    _dt   := ds/emax(1,abs(_RHS_(_tt,_ztr)));
+    _th   := _tt + .5_dt;
+    _dz1  := _dt*_RHS_(_tt, _ztr);  % displacement based on current point
+    _ztmp := _ztr + .5_dz1;         % 1st midpoint
+  % use _th instead of twice calculating (_tt + .5_dt)
+    _dz2  := _dt*_RHS_(_th, _ztmp); % displacement based on 1st midpoint
+    _ztmp := _ztr + .5_dz2;         % 2nd midpoint
+    _dz3  := _dt*_RHS_(_th, _ztmp); % displacement based on 2nd midpoint
+    _ztmp := _ztr + _dz3;           % temporary end point
+  % get time for next loop now since we need it right away in next line:
+    _tt   := _tt + _dt;
+    _dz4  := _dt*_RHS_(_tt, _ztmp); % displacement based on end point
+  % get next point
+    _ztr  := _ztr + (_dz1 + 2_dz2 + 2_dz3 + _dz4)/6;
+    _trj[_idx] := _ztr;
+  endfor
+  mkpath (sm, tens, false, _trj)
+enddef;
+
+def xyRKIV (expr sm) = txyRKIV (sm, default_tension)
+enddef;
+vardef txyRKIV (expr sm, tens, zstart, ds, N) (text _RHS_) =
+  save _fgxy, __fgxy;
+  vardef __fgxy (expr t, x, y) =  _RHS_  enddef;
+  vardef _fgxy (expr t, Z) = __fgxy(t, xpart Z, ypart Z) enddef;
+  tRKIV (sm, tens, zstart, ds, N) (_fgxy)
+enddef;
+
+def odeRKIV (expr sm) = todeRKIV (sm, default_tension)
+enddef;
+vardef todeRKIV (expr sm, tens, xstart, ystart, ds, N)
+(text _fxy) =
+  txyRKIV (sm, tens, (xstart, ystart), ds, N) ((1, _fxy))
+enddef;
+
 vardef lclosed expr f =
   f
   if not cycle f: