home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / e / e032 / 3.ddi / FILES / GRAPHICS.PAK / PLOTFIE1.M < prev    next >
Encoding:
Text File  |  1992-07-29  |  7.1 KB  |  239 lines
  1.  
  2. (* :Title: 3D Vector Fields *)
  3.  
  4. (* :Context: Graphics`PlotField3D` *)
  5.  
  6. (* :Author: John M. Novak *)
  7.  
  8. (* :Summary:
  9.     Plots vector fields in 3D
  10. *)
  11.  
  12. (* :Mathematica Version: 2.0 *)
  13.  
  14. (* :Package Version: 1.0 *)
  15.  
  16. (* :History:
  17.     V 1.0 April 91 by John M. Novak - based extensively on
  18.         PlotField.m by Kevin McIsaac, Mike Chan, ECM, and John Novak
  19.         VectorField3D.m by Wolfram Research `90 and ECM
  20. *)
  21.  
  22. (* :Keywords:
  23.     vector fields, gradient field, 3D graphics
  24. *)
  25.  
  26. (* :Limitations: *)
  27.  
  28. BeginPackage["Graphics`PlotField3D`","Utilities`FilterOptions`"]
  29.  
  30. ListPlotVectorField3D::usage =
  31.     "ListPlotVectorField3D[{{pt,vec},...},(options)] plots a
  32.     list of vectors in three dimensions, each vector based at
  33.     a corresponding point pt."
  34.  
  35. PlotVectorField3D::usage =
  36.     "PlotVectorField3D[{xfunc,yfunc,zfunc},xrange,yrange,zrange]
  37.     plots a vector field designated by the given functions,
  38.     over the given ranges, where a range is described as
  39.     {variable,min,max,(increment)}.  Also accepts options like
  40.     those of ListPlotVectorField3D."
  41.  
  42. PlotGradientField3D::usage =
  43.     "PlotGradientField3D[function,xrange,yrange,zrange,(options)]
  44.     plots the gradient of the given scalar function, over the
  45.     designated ranges, where a range is given as {variable,
  46.     min,max,(increment)}."
  47.  
  48. ScaleFactor::usage=
  49.     "ScaleFactor is an option for the PlotField3D functions that
  50.     scales the vectors to a specified length.  Default is Automatic;
  51.     at this setting, those functions that use a coordinate grid
  52.     (PlotVectorField3D, etc.) have the vectors scaled to this grid.
  53.     This scaling is applied after ScaleFunction and MaxArrowLength.";
  54.  
  55. ScaleFunction::usage=
  56.     "ScaleFunction rescales each vector to a length determined by applying
  57.     a pure function to the current length of that vector.  It will ignore
  58.     vectors of 0 magnitude. Note that because this is applied before the
  59.     ScaleFactor, this is most useful for resizing the relative lengths of the
  60.     vectors.  This is also applied before MaxArrowLength."
  61.  
  62. MaxArrowLength::usage=
  63.     "MaxArrowLength is an option for the PlotField3D functions
  64.     that determines the largest vector to be drawn. The 
  65.     value is compared to the magnitudes of all 
  66.     the vectors and causes all longer vectors to not be
  67.     drawn. This is applied after the ScaleFunction but before the 
  68.     ScaleFactor. Default is None (no maximum.)";
  69.  
  70. ColorFunction::usage=
  71.     "ColorFunction is an option for the PlotField3D functions that
  72.     determines the color and style used to display the vectors. It
  73.     is a pure function that accepts a value of 0 to 1; 0 corresponds
  74.     to the shortest vector, 1 the longest.";
  75.  
  76. VectorHeads::usage =
  77.     "VectorHeads is an option for the PlotField3D functions that
  78.     determines whether the vectors will be displayed with heads.
  79.     Default is VectorHeads->False."
  80.  
  81. Begin["`Private`"]
  82.  
  83. cross3[{a1_, a2_, a3_}, {b1_, b2_, b3_}] := 
  84.     {-(a3 b2) + a2 b3, a3 b1 - a1 b3,  -(a2 b1) + a1 b2}
  85.  
  86. mag[a_] := Sqrt[Apply[Plus, a^2]]
  87.  
  88. automatic[x_, value_] :=
  89.     If[x === Automatic, value, x]
  90.  
  91. vector3D[point:{x_, y_, z_}, grad:{dx_, dy_, dz_},False] :=
  92.     Line[{point, point + grad}]
  93.  
  94. vector3D[point:{x_,y_,z_}, grad:{dx_,dy_,dz_},True] :=
  95.     Point[{x,y,z}]/;grad == {0,0,0}
  96.  
  97. vector3D[point:{x_, y_, z_}, grad:{dx_, dy_, dz_},True] :=
  98.     Module[{endpoint, perp, perpm, offsetPoint,
  99.            arrowA, arrowB, arrowC, arrowD},
  100.       endpoint = point + grad;
  101.  
  102.       perp = cross3[grad, {0,0,1}];
  103.       perpm = mag[perp];
  104.       If[perpm == 0,
  105.         perp = cross3[grad, {0,1,0}];
  106.         perpm = mag[perp]
  107.       ];
  108.       perp = perp mag[grad]/(7 perpm);
  109.       
  110.       offsetPoint = point + 4/5 grad;
  111.       arrowA = offsetPoint + perp;
  112.       
  113.       perp = cross3[grad, perp];
  114.       perp = perp mag[grad]/(7 mag[perp]);
  115.       arrowB = offsetPoint + perp;
  116.       
  117.       perp = cross3[grad, perp];
  118.       perp = perp mag[grad]/(7 mag[perp]);
  119.       arrowC = offsetPoint + perp;
  120.       
  121.       perp = cross3[grad, perp];
  122.       perp = perp mag[grad]/(7 mag[perp]);
  123.       arrowD = offsetPoint + perp;
  124.       
  125.       {Line[{point, endpoint}],             (* 3D arrow shaft *)
  126.        Line[{arrowA, endpoint, arrowC}],         (* point of arrow *)
  127.        Line[{arrowB, endpoint, arrowD}],         (* point of arrow *)
  128.        Line[{arrowA, arrowB, arrowC, arrowD, arrowA}] (* base of point *)
  129.       }
  130.     ]
  131.  
  132. Options[ListPlotVectorField3D] = 
  133.     {ScaleFactor->Automatic, 
  134.      ScaleFunction->None,
  135.      MaxArrowLength->None,
  136.      ColorFunction->None,
  137.      VectorHeads->False};
  138.  
  139.  
  140. ListPlotVectorField3D[vects:{{_?VectorQ,_?VectorQ}..},opts___] :=
  141.     Module[{maxsize,scale,scalefunct,colorfunct,heads,points,
  142.             vectors,mags,colors,scaledmag,allvecs,vecs=N[vects]},
  143.         {maxsize,scale,scalefunct,colorfunct,heads} =
  144.             {MaxArrowLength,ScaleFactor,ScaleFunction,
  145.             ColorFunction,VectorHeads}/.{opts}/.
  146.             Options[ListPlotVectorField3D];
  147.         
  148.         (* option checking *)
  149.         If[Not[NumberQ[maxsize]] && maxsize != None,
  150.             maxsize = None,
  151.             maxsize = N[maxsize]];
  152.         If[Not[NumberQ[scale]] && scale != Automatic,
  153.             scale = Automatic,
  154.             scale = N[scale]];
  155.         heads = TrueQ[heads];
  156.         
  157.         vecs = Cases[vecs,{_,_?(VectorQ[#,NumberQ]&)}];
  158.         {points, vectors} = Transpose[vecs];
  159.         mags = Map[mag,vectors];
  160.         If[colorfunct == None, colorfunct = {}&];
  161.         If[Max[mags - Min[mags]] == 0,
  162.             colors = Map[colorfunct,Table[0,{Length[mags]}]],
  163.             colors = Map[colorfunct,
  164.                 (mags - Min[mags])/Max[mags - Min[mags]]]
  165.         ];
  166.  
  167.         If[scalefunct =!= None,
  168.              scaledmag = (If[# == 0, 0, scalefunct[#]]&) /@ mags;
  169.              vectors = MapThread[If[#2 == 0, {0,0,0}, #1 #2/#3]&,
  170.                 {vectors,scaledmag,mags}];
  171.              mags = scaledmag
  172.            ];
  173.  
  174.         allvecs = Transpose[{colors,points,vectors,mags}];  
  175.  
  176.         If[maxsize =!= None,
  177.              allvecs = Select[allvecs, (#[[4]]<=maxsize)&]
  178.            ];
  179.         
  180.         If[Max[mags] != 0,
  181.             scale = automatic[scale,Max[mags]]/Max[mags];
  182.              allvecs = Map[{#[[1]],#[[2]],scale #[[3]]}&,
  183.                  allvecs]
  184.          ];
  185.  
  186.         (* alternate method of vector generation requires pr.
  187.         pr = PlotRange[ Graphics3D[
  188.                 Flatten[Apply[Line[{#2,#2+#3}]&,allvecs,{1}]]]];
  189.         *)
  190.         
  191.         Show[Graphics3D[
  192.                  Flatten[Apply[{#1,vector3D[#2,#3,heads]}&,
  193.                      allvecs,{1}]],
  194.                  FilterOptions[Graphics, opts]]]
  195.     ]/; Last[Dimensions[vects]] === 3
  196.  
  197. Options[PlotVectorField3D] =
  198.     Join[Options[ListPlotVectorField3D],{PlotPoints->7}]
  199.  
  200. SetAttributes[PlotVectorField3D, HoldFirst]
  201.  
  202. PlotVectorField3D[f_, {u_, u0_, u1_, du_:Automatic},
  203.              {v_, v0_, v1_, dv_:Automatic},
  204.              {w_,w0_,w1_,dw_:Automatic},opts___] :=
  205.     Module[{plotpoints,dua,dva,dwa,vecs},
  206.         {plotpoints} = {PlotPoints}/.{opts}/.
  207.             Options[PlotVectorField3D];
  208.         dua = automatic[du,(u1 - u0)/(plotpoints-1)];
  209.         dva = automatic[dv,(v1 - v0)/(plotpoints-1)];
  210.         dwa = automatic[dw,(w1 - w0)/(plotpoints-1)];
  211.         vecs = Flatten[Table[{N[{u,v,w}],N[f]},
  212.             Evaluate[{u,u0,u1,dua}],Evaluate[{v,v0,v1,dva}],
  213.             Evaluate[{w,w0,w1,dwa}]],2];
  214.         ListPlotVectorField3D[vecs,
  215.             FilterOptions[ListPlotVectorField3D,opts],
  216.             FilterOptions[Graphics,opts],
  217.             ScaleFactor->N[Min[dua,dva,dwa]]]
  218.     ]
  219.  
  220. PlotGradientField3D[function_, 
  221.         {u_, u0__}, 
  222.         {v_, v0__},
  223.         {w_, w0__},
  224.         options___] :=
  225.     PlotVectorField3D[Evaluate[{D[function, u],
  226.                     D[function, v],D[function,w]}],
  227.             {u, u0},
  228.             {v, v0},
  229.             {w, w0},
  230.             options]
  231.  
  232. End[]
  233.  
  234. EndPackage[]
  235.  
  236. (* :Tests: *)
  237.  
  238. (* :Examples: *)
  239.