Re: [OctDev] formula2mat : fill a matrix with values of formula

[Paul K. <pki...@ja...>]

Etienne,

I would model this on fplot, which does this for functions of one variable. E.g.,

   fmesh('fn') - plots the bivariate function or expression.  It is a function
	if it is a single identifier otherwise it is an expression.

   fmesh(x,y,'fn') - use the given x and y grid spacing

   fmesh(...,[xmin xmax ymin ymax]) - use the given axis limits

   fmesh(...,n) - use the given number of steps in each dimension x, y

   z = fmesh(...) - return the mesh rather than plotting

   [x,y,z] = fmesh(...) - return the mesh and the mesh spacing

   fmesh(...,'loop') - don't vectorize the expression, but instead evaluate
	it once for each point in x cross y.

E.g.,

   sombrero = fmesh('sinc(x^2+y^2)',[-8 8 -8 8])


To distinguish between expression and function use:

   if all(isalnum(expr)), z=feval(expr,x,y); else eval(["z=",expr,";"]); end

To vectorize we will need a function which converts bare ^ and * to .^ and .*.
This would be something like:

   expr = strrep(expr,'*','.*');
   expr = strrep(expr,'^','.^');
   expr = strrep(expr,'..','.');

It may be faster to do things directly rather than using strrep:

   # insert a '.' before each * and ^
   idx = expr=="*" | expr== "^";
   n = length(expr);
   target = [1:n] + cumsum(idx);
   vec = repmat('.',1,target(n));
   vec(target) = expr;

   # remove duplicate '..'
   n = length(vec);
   idx = find(vec(1:n-1)=='.' & vec(2:n)=='.');
   vec(idx) = []; 

I would put it in main/plot.

What do you think?
   
- Paul

On Tue, Sep 10, 2002 at 08:37:40PM +0100, Etienne Grossmann wrote:
> 
>   Hello,
> 
> what about a function like this (in main/miscellaneous)? For people
> too lazy to call linspace() and meshgrid().
> 
>   Cheers,
> 
>   Etienne
> 
> ======================================================================
> ## z = formula2mat (f, ...)    - Fill matrix 'z' with values of formula 'f'
> ##
> ## INPUT ------
> ## f    :string: A formula written in octave code.
> ##
> ## OPTIONS ----
> ## bndx : 2x1 : Boundary  of x in which f is evaluated (overrides bnd)
> ## bndy : 2x1 : Boundary  of y in which f is evaluated (overrides bnd)
> ## bnd  : 2x2 : bnd = [bndx; bndy] Boundary in which f is evaluated
> ##                                                     Default = [-1,1;-1,1]
> ##                                                             
> ## sz   : 2x1 : Number of x and y values at which f is evaluated
> ##                                                     Default =   [11,11]
> ## 
> ## x    : vec : Values of x at which f is evaluated (overrides bndX and szX)
> ## y    : vec : Values of y at which f is evaluated (overrides bndY and szY)
> ##
> ## vec  : bool: If true, assume that z can be obtained by evaluating f when
> ##              x and y are matrices, rather than scalars.
> ## 
> ## OUTPUT -----
> ## z    : sz-sized matrix containing the values of f.
> ##
> ## EXAMPLES ---
> ##                                ## Fills a matrix z containing the values
> ##                                ## "x^2 - y^2" for x and y varying in
> ##                                ## [-1,1] by intervals of 0.2.
> ## z = formula2mat ("x^2 - y^2");
> ## 
> ##                                ## Same thing, quicker.
> ## z = formula2mat ("x.^2 + y.^2","vec",1);
> ##
> function zz = formula2mat (f, ...)
> 
> op1 = " x y bndx bndy bnd sz szx szy vec ";
> 
> default = struct ("x", nan,\
> 		  "y", nan,\
> 		  "bnd",[-1,1;-1,1],\
> 		  "bndx",nan,\
> 		  "bndy",nan,\
> 		  "sz",[11,11],\
> 		  "szx",nan,\
> 		  "szy",nan,\
> 		  "vec",0);
> 
> if nargin > 1,
>   opts = read_options (list (all_va_args), "default", default, \
> 		       "op1", op1);
> else opts = default;
> end
> 
> [x,y,bnd, bndx, bndy, sz, szx, szy, vec] = \
>     getfield (opts, "x","y","bnd" ,"bndx", "bndy", "sz", "szx", "szy", "vec");
> 
> ## szx overwrites sz (same for szy)
> if ! isnan (szx), sz(1) = szx; end
> if ! isnan (szy), sz(2) = szy; end
> 
> if length (sz) != 2
>   error ("length of sz should be 2, not %i\n", length (sz));
> end
> bnd
> ## bndx overwrites bnd
> if ! isnan (bndx(1)), bnd(1,1) = bndx(1); end
> if length (bndx)>1 && ! isnan (bndx(2)), bnd(1,2) = bndy(2); end
> if ! isnan (bndy(1)), bnd(2,1) = bndy(1); end
> if length (bndy)>1 && ! isnan (bndy(2)), bnd(2,2) = bndy(2); end
> bnd
> if any (size (bnd) != [2,2])
>   error ("size of bnd should be [2,2], not [%i,%i]\n", size (bnd));
> end
> 
> ## If x or y are not def'd, define them now
> if isnan (x), x = linspace (bnd(1,1),bnd(1,2),sz(1)); end
> if isnan (y), y = linspace (bnd(2,1),bnd(2,2),sz(2)); end
> 
> 
> [xx,yy] = meshgrid (x,y);
> 
> 
> ## x and y will become scalar variables below
> 
> 
> if ! vec
>   zz = zeros (sz);
>   eval (sprintf (["for ii = 1:sz(1),\n",\
> 		  "  for jj = 1:sz(2),\n",\
> 		  "    x = xx(ii,jj); y = yy(ii,jj);\n",\
> 		  "    zz(ii,jj) = %s;\n",\
> 		  "  end;\n",\
> 		  "end;\n"],\
> 		 f));
> else
>   x = xx;
>   y = yy;
>   eval (["zz = ",f,";"]);
> end
> 
> -- 
> Etienne Grossmann ------ http://www.isr.ist.utl.pt/~etienne
> 
> 
> -------------------------------------------------------
> This sf.net email is sponsored by: OSDN - Tired of that same old
> cell phone?  Get a new here for FREE!
> https://www.inphonic.com/r.asp?r=sourceforge1&refcode1=vs3390
> _______________________________________________
> Octave-dev mailing list
> Oct...@li...
> https://lists.sourceforge.net/lists/listinfo/octave-dev