1######################################################################## 2## 3## Copyright (C) 2019-2021 The Octave Project Developers 4## 5## See the file COPYRIGHT.md in the top-level directory of this 6## distribution or <https://octave.org/copyright/>. 7## 8## This file is part of Octave. 9## 10## Octave is free software: you can redistribute it and/or modify it 11## under the terms of the GNU General Public License as published by 12## the Free Software Foundation, either version 3 of the License, or 13## (at your option) any later version. 14## 15## Octave is distributed in the hope that it will be useful, but 16## WITHOUT ANY WARRANTY; without even the implied warranty of 17## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 18## GNU General Public License for more details. 19## 20## You should have received a copy of the GNU General Public License 21## along with Octave; see the file COPYING. If not, see 22## <https://www.gnu.org/licenses/>. 23## 24######################################################################## 25 26## -*- texinfo -*- 27## @deftypefn {} {} streamline (@var{x}, @var{y}, @var{z}, @var{u}, @var{v}, @var{w}, @var{sx}, @var{sy}, @var{sz}) 28## @deftypefnx {} {} streamline (@var{u}, @var{v}, @var{w}, @var{sx}, @var{sy}, @var{sz}) 29## @deftypefnx {} {} streamline (@dots{}, @var{options}) 30## @deftypefnx {} {} streamline (@var{hax}, @dots{}) 31## @deftypefnx {} {@var{h} =} streamline (@dots{}) 32## Plot streamlines of 2-D or 3-D vector fields. 33## 34## Plot streamlines of a 2-D or 3-D vector field given by 35## @code{[@var{u}, @var{v}]} or @code{[@var{u}, @var{v}, @var{w}]}. The vector 36## field is defined over a rectangular grid given by @code{[@var{x}, @var{y}]} 37## or @code{[@var{x}, @var{y}, @var{z}]}. The streamlines start at the seed 38## points @code{[@var{sx}, @var{sy}]} or @code{[@var{sx}, @var{sy}, @var{sz}]}. 39## 40## The input parameter @var{options} is a 2-D vector of the form 41## @code{[@var{stepsize}, @var{max_vertices}]}. The first parameter 42## specifies the step size used for trajectory integration (default 0.1). A 43## negative value is allowed which will reverse the direction of integration. 44## The second parameter specifies the maximum number of segments used to 45## create a streamline (default 10,000). 46## 47## If the first argument @var{hax} is an axes handle, then plot into this axes, 48## rather than the current axes returned by @code{gca}. 49## 50## The optional return value @var{h} is a graphics handle to the hggroup 51## comprising the field lines. 52## 53## Example: 54## 55## @example 56## @group 57## [x, y] = meshgrid (-1.5:0.2:2, -1:0.2:2); 58## u = - x / 4 - y; 59## v = x - y / 4; 60## streamline (x, y, u, v, 1.7, 1.5); 61## @end group 62## @end example 63## 64## @seealso{stream2, stream3, streamtube, ostreamtube} 65## @end deftypefn 66 67function h = streamline (varargin) 68 69 [hax, varargin, nargin] = __plt_get_axis_arg__ ("streamline", varargin{:}); 70 71 if (nargin == 0) 72 print_usage (); 73 endif 74 75 nd = ndims (varargin{1}); 76 if (nd > 3) 77 error ("streamline: input data must be 2-D or 3-D"); 78 endif 79 80 if (isempty (hax)) 81 hax = gca (); 82 else 83 hax = hax(1); 84 endif 85 86 h = []; 87 if (nd == 2) 88 xy = stream2 (varargin{:}); 89 for i = 1 : length (xy) 90 sl = xy{i}; 91 if (! isempty (sl)) 92 htmp = line (hax, "xdata", sl(:, 1), "ydata", sl(:, 2), "color", "b"); 93 h = [h; htmp]; 94 endif 95 endfor 96 else 97 xyz = stream3 (varargin{:}); 98 for i = 1 : length (xyz) 99 sl = xyz{i}; 100 if (! isempty (sl)) 101 htmp = line (hax, 102 "xdata", sl(:, 1), "ydata", sl(:, 2), "zdata", sl(:, 3), 103 "color", "b"); 104 h = [h; htmp]; 105 endif 106 endfor 107 endif 108 109endfunction 110 111 112%!demo 113%! clf; 114%! [x, y] = meshgrid (-2:0.5:2); 115%! u = - y - x / 2; 116%! v = x - y / 2; 117%! [sx, sy] = meshgrid (-2:2:2); 118%! h = streamline (x, y, u, v, sx, sy); 119%! set (h, "color", "r"); 120%! hold on; 121%! quiver (x, y, u, v); 122%! scatter (sx(:), sy(:), 20, "filled", "o", "markerfacecolor", "r"); 123%! title ("Spiral Sink"); 124%! grid on; 125%! axis equal; 126 127%!demo 128%! clf; 129%! [x, y, z] = meshgrid (-3:3); 130%! u = - x / 2 - y; 131%! v = x - y / 2; 132%! w = - z; 133%! [sx, sy, sz] = meshgrid (3, 0:1.5:1.5, 0:1.5:3); 134%! h = streamline (x, y, z, u, v, w, sx, sy, sz); 135%! set (h, "color", "r"); 136%! hold on; 137%! quiver3 (x, y, z, u, v, w); 138%! scatter3 (sx(:), sy(:), sz(:), 20, "filled", "o", "markerfacecolor", "r"); 139%! view (3); 140%! title ("Spiral Sink"); 141%! grid on; 142%! axis equal; 143 144%!demo 145%! clf; 146%! [x, y, z] = meshgrid (-1:0.4:1, -1:0.4:1, -3:0.3:0); 147%! a = 0.08; 148%! b = 0.04; 149%! u = - a * x - y; 150%! v = x - a * y; 151%! w = - b * ones (size (x)); 152%! hold on; 153%! sx = 1.0; 154%! sy = 0.0; 155%! sz = 0.0; 156%! plot3 (sx, sy, sz, ".r", "markersize", 15); 157%! t = linspace (0, 12 * 2 * pi(), 500); 158%! tx = exp (-a * t).*cos (t); 159%! ty = exp (-a * t).*sin (t); 160%! tz = - b * t; 161%! plot3 (tx, ty, tz, "-b"); 162%! h = streamline (x, y, z, u, v, w, sx, sy, sz); 163%! set (h, "color", "r"); 164%! view (3); 165%! title ("Heuns Scheme (red) vs. Analytical Solution (blue)"); 166%! grid on; 167%! axis equal tight; 168 169## Test input validation 170%!error streamline () 171%!error <Invalid call to streamline> 172%! hf = figure ("visible", "off"); 173%! unwind_protect 174%! hax = axes (); 175%! streamline (hax); 176%! unwind_protect_cleanup 177%! close (hf); 178%! end_unwind_protect 179%!error <input data must be 2-D or 3-D> streamline (ones (2,2,2,2)) 180