changeset 30893:e1788b1a315f

maint: Use "fcn" as preferred abbreviation for "function" in m-files. * accumarray.m, accumdim.m, quadl.m, quadv.m, randi.m, structfun.m, __is_function__.m, uigetfile.m, uimenu.m, uiputfile.m, doc_cache_create.m, colorspace_conversion_input_check.m, imageIO.m, argnames.m, vectorize.m, vectorize.m, normest1.m, inputname.m, nthargout.m, display_info_file.m, decic.m, ode15i.m, ode15s.m, ode23.m, ode23s.m, ode45.m, odeset.m, check_default_input.m, integrate_adaptive.m, ode_event_handler.m, runge_kutta_23.m, runge_kutta_23s.m, runge_kutta_45_dorpri.m, runge_kutta_interpolate.m, starting_stepsize.m, __all_opts__.m, fminbnd.m, fminsearch.m, fminunc.m, fsolve.m, fzero.m, sqp.m, fplot.m, plotyy.m, __bar__.m, __ezplot__.m, flat_entry.html, profexport.m, movfun.m, bicg.m, bicgstab.m, cgs.m, eigs.m, gmres.m, pcg.m, __alltohandles__.m, __sprand__.m, qmr.m, tfqmr.m, dump_demos.m: Replace "func", "fun", "fn" in documentation and variable names with "fcn".
author Rik <rik@octave.org>
date Mon, 04 Apr 2022 18:14:56 -0700
parents 1a3cc2811090
children 9c686fa8132e
files scripts/general/accumarray.m scripts/general/accumdim.m scripts/general/quadl.m scripts/general/quadv.m scripts/general/randi.m scripts/general/structfun.m scripts/gui/private/__is_function__.m scripts/gui/uigetfile.m scripts/gui/uimenu.m scripts/gui/uiputfile.m scripts/help/doc_cache_create.m scripts/image/private/colorspace_conversion_input_check.m scripts/image/private/imageIO.m scripts/legacy/@inline/argnames.m scripts/legacy/@inline/vectorize.m scripts/legacy/vectorize.m scripts/linear-algebra/normest1.m scripts/miscellaneous/inputname.m scripts/miscellaneous/nthargout.m scripts/miscellaneous/private/display_info_file.m scripts/ode/decic.m scripts/ode/ode15i.m scripts/ode/ode15s.m scripts/ode/ode23.m scripts/ode/ode23s.m scripts/ode/ode45.m scripts/ode/odeset.m scripts/ode/private/check_default_input.m scripts/ode/private/integrate_adaptive.m scripts/ode/private/ode_event_handler.m scripts/ode/private/runge_kutta_23.m scripts/ode/private/runge_kutta_23s.m scripts/ode/private/runge_kutta_45_dorpri.m scripts/ode/private/runge_kutta_interpolate.m scripts/ode/private/starting_stepsize.m scripts/optimization/__all_opts__.m scripts/optimization/fminbnd.m scripts/optimization/fminsearch.m scripts/optimization/fminunc.m scripts/optimization/fsolve.m scripts/optimization/fzero.m scripts/optimization/sqp.m scripts/plot/draw/fplot.m scripts/plot/draw/plotyy.m scripts/plot/draw/private/__bar__.m scripts/plot/draw/private/__ezplot__.m scripts/profiler/html/flat_entry.html scripts/profiler/profexport.m scripts/signal/movfun.m scripts/sparse/bicg.m scripts/sparse/bicgstab.m scripts/sparse/cgs.m scripts/sparse/eigs.m scripts/sparse/gmres.m scripts/sparse/pcg.m scripts/sparse/private/__alltohandles__.m scripts/sparse/private/__sprand__.m scripts/sparse/qmr.m scripts/sparse/tfqmr.m scripts/testfun/private/dump_demos.m
diffstat 60 files changed, 973 insertions(+), 972 deletions(-) [+]
line wrap: on
line diff
--- a/scripts/general/accumarray.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/accumarray.m	Mon Apr 04 18:14:56 2022 -0700
@@ -26,9 +26,9 @@
 ## -*- texinfo -*-
 ## @deftypefn  {} {@var{A} =} accumarray (@var{subs}, @var{vals})
 ## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz})
-## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{func})
-## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{func}, @var{fillval})
-## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{func}, @var{fillval}, @var{issparse})
+## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{fcn})
+## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{fcn}, @var{fillval})
+## @deftypefnx {} {@var{A} =} accumarray (@var{subs}, @var{vals}, @var{sz}, @var{fcn}, @var{fillval}, @var{issparse})
 ##
 ## Create an array by accumulating the elements of a vector into the
 ## positions defined by their subscripts.
@@ -49,14 +49,14 @@
 ##
 ## The default action of @code{accumarray} is to sum the elements with
 ## the same subscripts.  This behavior can be modified by defining the
-## @var{func} function.  This should be a function or function handle
+## @var{fcn} function.  This should be a function or function handle
 ## that accepts a column vector and returns a scalar.  The result of the
 ## function should not depend on the order of the subscripts.
 ##
 ## The elements of the returned array that have no subscripts associated
 ## with them are set to zero.  Defining @var{fillval} to some other value
 ## allows these values to be defined.  This behavior changes, however,
-## for certain values of @var{func}.  If @var{func} is @code{@@min}
+## for certain values of @var{fcn}.  If @var{fcn} is @code{@@min}
 ## (respectively, @code{@@max}) then the result will be filled with the
 ## minimum (respectively, maximum) integer if @var{vals} is of integral
 ## type, logical false (respectively, logical true) if @var{vals} is of
@@ -126,7 +126,7 @@
 ## The complexity of accumarray in general for the non-sparse case is
 ## generally O(M+N), where N is the number of subscripts and M is the
 ## maximum subscript (linearized in multi-dimensional case).  If
-## @var{func} is one of @code{@@sum} (default), @code{@@max},
+## @var{fcn} is one of @code{@@sum} (default), @code{@@max},
 ## @code{@@min} or @code{@@(x) @{x@}}, an optimized code path is used.
 ## Note that for general reduction function the interpreter overhead can
 ## play a major part and it may be more efficient to do multiple
@@ -135,7 +135,7 @@
 ## @seealso{accumdim, unique, sparse}
 ## @end deftypefn
 
-function A = accumarray (subs, vals, sz = [], func = [], fillval = [], issparse = [])
+function A = accumarray (subs, vals, sz = [], fcn = [], fillval = [], issparse = [])
 
   if (nargin < 2)
     print_usage ();
@@ -163,10 +163,10 @@
     endif
   endif
 
-  if (isempty (func))
-    func = @sum;
-  elseif (! is_function_handle (func))
-    error ("accumarray: FUNC must be a function handle");
+  if (isempty (fcn))
+    fcn = @sum;
+  elseif (! is_function_handle (fcn))
+    error ("accumarray: FCN must be a function handle");
   endif
 
   if (isempty (fillval))
@@ -204,7 +204,7 @@
       error ("accumarray: in the sparse case, values must be numeric or logical");
     endif
 
-    if (func != @sum)
+    if (fcn != @sum)
 
       ## Reduce values.  This is not needed if we're about to sum them,
       ## because "sparse" can do that.
@@ -216,7 +216,7 @@
       jdx = find (any (diff (subs, 1, 1), 2));
       jdx = [jdx; n];
 
-      vals = cellfun (func, mat2cell (vals(:)(idx), diff ([0; jdx])));
+      vals = cellfun (fcn, mat2cell (vals(:)(idx), diff ([0; jdx])));
       subs = subs(jdx, :);
       mode = "unique";
     else
@@ -267,7 +267,7 @@
 
     ## Some built-in reductions handled efficiently.
 
-    if (func == @sum)
+    if (fcn == @sum)
       ## Fast summation.
       if (isempty (sz))
         A = __accumarray_sum__ (subs, vals);
@@ -283,7 +283,7 @@
         mask(subs) = false;
         A(mask) = fillval;
       endif
-    elseif (func == @max)
+    elseif (fcn == @max)
       ## Fast maximization.
 
       if (isinteger (vals))
@@ -310,7 +310,7 @@
         mask(subs) = false;
         A(mask) = fillval;
       endif
-    elseif (func == @min)
+    elseif (fcn == @min)
       ## Fast minimization.
 
       if (isinteger (vals))
@@ -358,8 +358,8 @@
       ## Optimize the case when function is @(x) {x}, i.e., we just want
       ## to collect the values to cells.
       persistent simple_cell_str = func2str (@(x) {x});
-      if (! strcmp (func2str (func), simple_cell_str))
-        vals = cellfun (func, vals);
+      if (! strcmp (func2str (fcn), simple_cell_str))
+        vals = cellfun (fcn, vals);
       endif
 
       subs = subs(jdx);
@@ -458,9 +458,9 @@
 %!   assert (accumarray (zeros (0, 1), [], [] , funcs{idx}), zeros (0, 1));
 %! endfor
 
-## Matlab returns an array of doubles even though FUNC returns cells.  In
+## Matlab returns an array of doubles even though FCN returns cells.  In
 ## Octave, we do not have that bug, at least for this case.
 %!assert (accumarray (zeros (0, 1), [], [0 1] , @(x) {x}), cell (0, 1))
 
-%!error <FUNC must be a function handle>
+%!error <FCN must be a function handle>
 %! accumarray ([1; 2; 3], [1; 2; 3], [3 1], '@(x) {x}')
--- a/scripts/general/accumdim.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/accumdim.m	Mon Apr 04 18:14:56 2022 -0700
@@ -27,8 +27,8 @@
 ## @deftypefn  {} {@var{A} =} accumdim (@var{subs}, @var{vals})
 ## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim})
 ## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim}, @var{n})
-## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim}, @var{n}, @var{func})
-## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim}, @var{n}, @var{func}, @var{fillval})
+## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim}, @var{n}, @var{fcn})
+## @deftypefnx {} {@var{A} =} accumdim (@var{subs}, @var{vals}, @var{dim}, @var{n}, @var{fcn}, @var{fillval})
 ## Create an array by accumulating the slices of an array into the
 ## positions defined by their subscripts along a specified dimension.
 ##
@@ -43,7 +43,7 @@
 ##
 ## The default action of @code{accumdim} is to sum the subarrays with the
 ## same subscripts.  This behavior can be modified by defining the
-## @var{func} function.  This should be a function or function handle
+## @var{fcn} function.  This should be a function or function handle
 ## that accepts an array and a dimension, and reduces the array along
 ## this dimension.  As a special exception, the built-in @code{min} and
 ## @code{max} functions can be used directly, and @code{accumdim}
@@ -69,7 +69,7 @@
 ## @seealso{accumarray}
 ## @end deftypefn
 
-function A = accumdim (subs, vals, dim, n = 0, func = [], fillval = 0)
+function A = accumdim (subs, vals, dim, n = 0, fcn = [], fillval = 0)
 
   if (nargin < 2)
     print_usage ();
@@ -108,7 +108,7 @@
     error ("accumdim: dimension mismatch");
   endif
 
-  if (isempty (func) || func == @sum)
+  if (isempty (fcn) || fcn == @sum)
     ## Fast summation case.
     A = __accumdim_sum__ (subs, vals, dim, n);
 
@@ -137,10 +137,10 @@
   subsc{dim} = idx;
   vals = mat2cell (vals(subsc{:}), szc{:});
   ## Apply reductions.  Special case min, max.
-  if (func == @min || func == @max)
-    vals = cellfun (func, vals, {[]}, {dim}, "uniformoutput", false);
+  if (fcn == @min || fcn == @max)
+    vals = cellfun (fcn, vals, {[]}, {dim}, "uniformoutput", false);
   else
-    vals = cellfun (func, vals, {dim}, "uniformoutput", false);
+    vals = cellfun (fcn, vals, {dim}, "uniformoutput", false);
   endif
   subs = subs(jdx);
 
--- a/scripts/general/quadl.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/quadl.m	Mon Apr 04 18:14:56 2022 -0700
@@ -28,7 +28,7 @@
 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol})
 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace})
 ## @deftypefnx {} {@var{q} =} quadl (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace}, @var{p1}, @var{p2}, @dots{})
-## @deftypefnx {} {[@var{q}, @var{nfun}] =} quadl (@dots{})
+## @deftypefnx {} {[@var{q}, @var{nfev}] =} quadl (@dots{})
 ##
 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b} using
 ## an adaptive @nospell{Lobatto} rule.
@@ -55,7 +55,7 @@
 ##
 ## The result of the integration is returned in @var{q}.
 ##
-## The optional output @var{nfun} indicates the total number of function
+## The optional output @var{nfev} indicates the total number of function
 ## evaluations performed.
 ##
 ## Reference: @nospell{W. Gander and W. Gautschi}, @cite{Adaptive Quadrature -
@@ -72,7 +72,7 @@
 ## 2003-08-05 Shai Ayal
 ##   * permission from author to release as GPL
 
-function [q, nfun] = quadl (f, a, b, tol = [], trace = false, varargin)
+function [q, nfev] = quadl (f, a, b, tol = [], trace = false, varargin)
 
   if (nargin < 3)
     print_usage ();
@@ -97,17 +97,17 @@
   endif
 
   y = feval (f, [a, b], varargin{:});
-  nfun = 1;
+  nfev = 1;
 
   fa = y(1);
   fb = y(2);
 
   h = b - a;
 
-  [q, nfun, hmin] = adaptlobstp (f, a, b, fa, fb, Inf, nfun, abs (h),
+  [q, nfev, hmin] = adaptlobstp (f, a, b, fa, fb, Inf, nfev, abs (h),
                                  tol, trace, varargin{:});
 
-  if (nfun > 10_000)
+  if (nfev > 10_000)
     warning ("quadl: maximum iteration count reached -- possible singular integral");
   elseif (any (! isfinite (q(:))))
     warning ("quadl: infinite or NaN function evaluations were returned");
@@ -117,13 +117,13 @@
 
 endfunction
 
-function [q, nfun, hmin] = adaptlobstp (f, a, b, fa, fb, q0, nfun, hmin,
+function [q, nfev, hmin] = adaptlobstp (f, a, b, fa, fb, q0, nfev, hmin,
                                         tol, trace, varargin)
 
   persistent alpha = sqrt (2/3);
   persistent beta = 1 / sqrt (5);
 
-  if (nfun > 10_000)
+  if (nfev > 10_000)
     q = q0;
     return;
   endif
@@ -136,7 +136,7 @@
   mrr = m + alpha*h;
   x = [mll, ml, m, mr, mrr];
   y = feval (f, x, varargin{:});
-  nfun += 1;
+  nfev += 1;
   fmll = y(1);
   fml  = y(2);
   fm   = y(3);
@@ -150,25 +150,25 @@
   endif
 
   if (trace)
-    disp ([nfun, a, b-a, i1]);
+    disp ([nfev, a, b-a, i1]);
   endif
 
-  ## Force at least one adaptive step (nfun > 2 test).
-  if ((abs (i1-i2) < tol || mll <= a || b <= mrr) && nfun > 2)
+  ## Force at least one adaptive step (nfev > 2 test).
+  if ((abs (i1-i2) < tol || mll <= a || b <= mrr) && nfev > 2)
     q = i1;
   else
     q = zeros (6, 1, class (x));
-    [q(1), nfun, hmin] = adaptlobstp (f, a  , mll, fa  , fmll, q0/6, nfun, hmin,
+    [q(1), nfev, hmin] = adaptlobstp (f, a  , mll, fa  , fmll, q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
-    [q(2), nfun, hmin] = adaptlobstp (f, mll, ml , fmll, fml , q0/6, nfun, hmin,
+    [q(2), nfev, hmin] = adaptlobstp (f, mll, ml , fmll, fml , q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
-    [q(3), nfun, hmin] = adaptlobstp (f, ml , m  , fml , fm  , q0/6, nfun, hmin,
+    [q(3), nfev, hmin] = adaptlobstp (f, ml , m  , fml , fm  , q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
-    [q(4), nfun, hmin] = adaptlobstp (f, m  , mr , fm  , fmr , q0/6, nfun, hmin,
+    [q(4), nfev, hmin] = adaptlobstp (f, m  , mr , fm  , fmr , q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
-    [q(5), nfun, hmin] = adaptlobstp (f, mr , mrr, fmr , fmrr, q0/6, nfun, hmin,
+    [q(5), nfev, hmin] = adaptlobstp (f, mr , mrr, fmr , fmrr, q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
-    [q(6), nfun, hmin] = adaptlobstp (f, mrr, b  , fmrr, fb  , q0/6, nfun, hmin,
+    [q(6), nfev, hmin] = adaptlobstp (f, mrr, b  , fmrr, fb  , q0/6, nfev, hmin,
                                       tol, trace, varargin{:});
     q = sum (q);
   endif
@@ -189,11 +189,11 @@
 
 ## test different tolerances.
 %!test
-%! [q, nfun1] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.5, []);
+%! [q, nfev1] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.5, []);
 %! assert (q, (60 + sin (4) - sin (64))/12, 0.5);
-%! [q, nfun2] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.1, []);
+%! [q, nfev2] = quadl (@(x) sin (2 + 3*x).^2, 0, 10, 0.1, []);
 %! assert (q, (60 + sin (4) - sin (64))/12, 0.1);
-%! assert (nfun2 > nfun1);
+%! assert (nfev2 > nfev1);
 
 %!test  # test single input/output
 %! assert (class (quadl (@sin, 0, 1)), "double");
--- a/scripts/general/quadv.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/quadv.m	Mon Apr 04 18:14:56 2022 -0700
@@ -28,7 +28,7 @@
 ## @deftypefnx {} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol})
 ## @deftypefnx {} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace})
 ## @deftypefnx {} {@var{q} =} quadv (@var{f}, @var{a}, @var{b}, @var{tol}, @var{trace}, @var{p1}, @var{p2}, @dots{})
-## @deftypefnx {} {[@var{q}, @var{nfun}] =} quadv (@dots{})
+## @deftypefnx {} {[@var{q}, @var{nfev}] =} quadv (@dots{})
 ##
 ## Numerically evaluate the integral of @var{f} from @var{a} to @var{b}
 ## using an adaptive Simpson's rule.
@@ -57,7 +57,7 @@
 ##
 ## The result of the integration is returned in @var{q}.
 ##
-## The optional output @var{nfun} indicates the total number of function
+## The optional output @var{nfev} indicates the total number of function
 ## evaluations performed.
 ##
 ## Note: @code{quadv} is written in Octave's scripting language and can be
@@ -70,7 +70,7 @@
 ## Algorithm: See https://en.wikipedia.org/wiki/Adaptive_Simpson%27s_method
 ## for basic explanation.  See NOTEs and FIXME for Octave modifications.
 
-function [q, nfun] = quadv (f, a, b, tol = [], trace = false, varargin)
+function [q, nfev] = quadv (f, a, b, tol = [], trace = false, varargin)
 
   if (nargin < 3)
     print_usage ();
@@ -89,7 +89,7 @@
 
   if (trace)
     ## Print column headers once above trace display.
-    printf ("  nfun          a            (b - a)         q_interval\n");
+    printf ("  nfev          a            (b - a)         q_interval\n");
   endif
 
   ## NOTE: Split the interval into 3 parts which avoids problems with periodic
@@ -110,46 +110,46 @@
   fb2 = feval (f, b2, varargin{:});
   fc3 = feval (f, c3, varargin{:});
   fb  = feval (f, b,  varargin{:});
-  nfun = 7;
+  nfev = 7;
 
   ## NOTE: If there are edge singularities, move edge point by eps*(b-a) as
   ## discussed in Shampine paper used to implement quadgk.
   if (any (isinf (fa(:))))
     fa = feval (f, a + eps * (b-a), varargin{:});
-    nfun++;
+    nfev++;
   endif
   if (any (isinf (fb(:))))
     fb = feval (f, b - eps * (b-a), varargin{:});
-    nfun++;
+    nfev++;
   endif
 
   ## Region 1
   h = (b1 - a);
   q1 = h / 6 * (fa + 4*fc1 + fb1);
 
-  [q1, nfun, hmin1] = simpsonstp (f, a, b1, c1, fa, fb1, fc1, q1, tol,
-                                  nfun, abs (h), trace, varargin{:});
+  [q1, nfev, hmin1] = simpsonstp (f, a, b1, c1, fa, fb1, fc1, q1, tol,
+                                  nfev, abs (h), trace, varargin{:});
 
   ## Region 2
   h = (b2 - b1);
   q2 = h / 6 * (fb1 + 4*fc2 + fb2);
 
-  [q2, nfun, hmin2] = simpsonstp (f, b1, b2, c2, fb1, fb2, fc2, q2, tol,
-                                  nfun, abs (h), trace, varargin{:});
+  [q2, nfev, hmin2] = simpsonstp (f, b1, b2, c2, fb1, fb2, fc2, q2, tol,
+                                  nfev, abs (h), trace, varargin{:});
 
   ## Region 3
   h = (b - b2);
   q3 = h / 6 * (fb2 + 4*fc3 + fb);
 
-  [q3, nfun, hmin3] = simpsonstp (f, b2, b, c3, fb2, fb, fc3, q3, tol,
-                                  nfun, abs (h), trace, varargin{:});
+  [q3, nfev, hmin3] = simpsonstp (f, b2, b, c3, fb2, fb, fc3, q3, tol,
+                                  nfev, abs (h), trace, varargin{:});
 
   ## Total integral over all 3 regions and verify results
   q = q1 + q2 + q3;
 
   hmin = min ([hmin1, hmin2, hmin3]);
 
-  if (nfun > 10_000)
+  if (nfev > 10_000)
     warning ("quadv: maximum iteration count reached -- possible singular integral");
   elseif (any (! isfinite (q(:))))
     warning ("quadv: infinite or NaN function evaluations were returned");
@@ -159,10 +159,10 @@
 
 endfunction
 
-function [q, nfun, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q0, tol,
-                                       nfun, hmin, trace, varargin)
+function [q, nfev, hmin] = simpsonstp (f, a, b, c, fa, fb, fc, q0, tol,
+                                       nfev, hmin, trace, varargin)
 
-  if (nfun > 10_000)   # stop endless recursion
+  if (nfev > 10_000)   # stop endless recursion
     q = q0;
     return;
   endif
@@ -171,7 +171,7 @@
   e = (c + b) / 2;
   fd = feval (f, d, varargin{:});
   fe = feval (f, e, varargin{:});
-  nfun += 2;
+  nfev += 2;
   q1 = (c - a) / 6 * (fa + 4*fd + fc);
   q2 = (b - c) / 6 * (fc + 4*fe + fb);
   q = q1 + q2;
@@ -184,7 +184,7 @@
 
   if (trace)
     printf ("%5d   %#14.11g   %16.10e   %-16.11g\n",
-            nfun,  a,         b-a,      q + delta/15);
+            nfev,  a,         b-a,      q + delta/15);
   endif
 
   ## NOTE: Not vectorizing q-q0 in the norm provides a more rigid criterion
@@ -194,10 +194,10 @@
     ## each bisection interval should use tol/2.  However, Matlab does not
     ## do this, and it would also profoundly increase the number of function
     ## evaluations required.
-    [q1, nfun, hmin] = simpsonstp (f, a, c, d, fa, fc, fd, q1, tol,
-                                   nfun, hmin, trace, varargin{:});
-    [q2, nfun, hmin] = simpsonstp (f, c, b, e, fc, fb, fe, q2, tol,
-                                   nfun, hmin, trace, varargin{:});
+    [q1, nfev, hmin] = simpsonstp (f, a, c, d, fa, fc, fd, q1, tol,
+                                   nfev, hmin, trace, varargin{:});
+    [q2, nfev, hmin] = simpsonstp (f, c, b, e, fc, fb, fe, q2, tol,
+                                   nfev, hmin, trace, varargin{:});
     q = q1 + q2;
   else
     q += delta / 15;   # NOTE: Richardson extrapolation correction
--- a/scripts/general/randi.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/randi.m	Mon Apr 04 18:14:56 2022 -0700
@@ -188,9 +188,9 @@
 %! max_single = double (flintmax ("single"));
 
 ## Test that no warning thrown if IMAX is exactly on the limits of the range
-%!function test_no_warning (func, varargin)
+%!function test_no_warning (fcn, varargin)
 %!  lastwarn ("");
-%!  func (varargin{:});
+%!  fcn (varargin{:});
 %!  assert (lastwarn (), "");
 %!endfunction
 %!test test_no_warning (@randi, max_int8, "int8");
--- a/scripts/general/structfun.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/general/structfun.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,16 +24,16 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{A} =} structfun (@var{func}, @var{S})
-## @deftypefnx {} {@var{A} =} structfun (@dots{}, "ErrorHandler", @var{errfunc})
+## @deftypefn  {} {@var{A} =} structfun (@var{fcn}, @var{S})
+## @deftypefnx {} {@var{A} =} structfun (@dots{}, "ErrorHandler", @var{errfcn})
 ## @deftypefnx {} {@var{A} =} structfun (@dots{}, "UniformOutput", @var{val})
 ## @deftypefnx {} {[@var{A}, @var{B}, @dots{}] =} structfun (@dots{})
 ##
 ## Evaluate the function named @var{name} on the fields of the structure
-## @var{S}.  The fields of @var{S} are passed to the function @var{func}
+## @var{S}.  The fields of @var{S} are passed to the function @var{fcn}
 ## individually.
 ##
-## @code{structfun} accepts an arbitrary function @var{func} in the form of an
+## @code{structfun} accepts an arbitrary function @var{fcn} in the form of an
 ## inline function, function handle, or the name of a function (in a character
 ## string).  In the case of a character string argument, the function must
 ## accept a single argument named @var{x}, and it must return a string value.
@@ -57,16 +57,16 @@
 ## @end group
 ## @end example
 ##
-## Given the parameter @qcode{"ErrorHandler"}, @var{errfunc} defines a function
-## to call in case @var{func} generates an error.  The form of the function is
+## Given the parameter @qcode{"ErrorHandler"}, @var{errfcn} defines a function
+## to call in case @var{fcn} generates an error.  The form of the function is
 ##
 ## @example
-## function [@dots{}] = errfunc (@var{se}, @dots{})
+## function [@dots{}] = errfcn (@var{se}, @dots{})
 ## @end example
 ##
 ## @noindent
-## where there is an additional input argument to @var{errfunc} relative to
-## @var{func}, given by @nospell{@var{se}}.  This is a structure with the
+## where there is an additional input argument to @var{errfcn} relative to
+## @var{fcn}, given by @nospell{@var{se}}.  This is a structure with the
 ## elements @qcode{"identifier"}, @qcode{"message"} and @qcode{"index"},
 ## giving respectively the error identifier, the error message, and the index
 ## into the input arguments of the element that caused the error.  For an
@@ -75,7 +75,7 @@
 ## @seealso{cellfun, arrayfun, spfun}
 ## @end deftypefn
 
-function varargout = structfun (func, S, varargin)
+function varargout = structfun (fcn, S, varargin)
 
   if (nargin < 2)
     print_usage ();
@@ -105,7 +105,7 @@
   endif
 
   varargout = cell (max ([nargout, 1]), 1);
-  [varargout{:}] = cellfun (func, struct2cell (S), varargin{:});
+  [varargout{:}] = cellfun (fcn, struct2cell (S), varargin{:});
 
   if (! uniform_output)
     varargout = cellfun ("cell2struct", varargout, {fieldnames(S)}, {1}, ...
--- a/scripts/gui/private/__is_function__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/gui/private/__is_function__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,13 +24,14 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{result} =} __is_function__ (@var{func})
-## Undocumented internal function.
+## @deftypefn {} {@var{tf} =} __is_function__ (@var{fcn})
+## Return true if @var{fcn} is an m-file function, @file{.oct} or @file{.mex}
+## function, or a built-in function.
 ## @end deftypefn
 
-function result = __is_function__ (func)
+function tf = __is_function__ (fcn)
 
-  existval = exist (func);
-  result = (existval == 2 || existval == 3 || existval == 5);
+  existval = exist (fcn);
+  tf = (existval == 2 || existval == 3 || existval == 5);
 
 endfunction
--- a/scripts/gui/uigetfile.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/gui/uigetfile.m	Mon Apr 04 18:14:56 2022 -0700
@@ -41,7 +41,7 @@
 ## If a filename is given then the file extension is extracted and used as
 ## filter.  In addition, the path is selected as current path in the dialog and
 ## the filename is selected as default file.
-## Example: @code{uigetfile ("myfun.m")}
+## Example: @code{uigetfile ("myfcn.m")}
 ##
 ## @item A single file extension @qcode{"*.ext"}
 ## Example: @code{uigetfile ("*.ext")}
--- a/scripts/gui/uimenu.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/gui/uimenu.m	Mon Apr 04 18:14:56 2022 -0700
@@ -40,9 +40,9 @@
 ##
 ## @item @qcode{"callback"}
 ## Is the function called when this menu entry is executed.  It can be either a
-## function string (e.g., @qcode{"myfun"}), a function handle (e.g., @@myfun)
+## function string (e.g., @qcode{"myfcn"}), a function handle (e.g., @@myfcn)
 ## or a cell array containing the function handle and arguments for the
-## callback function (e.g., @{@@myfun, arg1, arg2@}).
+## callback function (e.g., @{@@myfcn, arg1, arg2@}).
 ##
 ## @item @qcode{"checked"}
 ## Can be set @qcode{"on"} or @qcode{"off"}.  Sets a mark at this menu entry.
--- a/scripts/gui/uiputfile.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/gui/uiputfile.m	Mon Apr 04 18:14:56 2022 -0700
@@ -37,7 +37,7 @@
 ## @item @qcode{"/path/to/filename.ext"}
 ## If a filename is given the file extension is extracted and used as filter.
 ## In addition the path is selected as current path in the dialog and the
-## filename is selected as default file.  Example: @code{uiputfile ("myfun.m")}
+## filename is selected as default file.  Example: @code{uiputfile ("myfcn.m")}
 ##
 ## @item @qcode{"*.ext"}
 ## A single file extension.
--- a/scripts/help/doc_cache_create.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/help/doc_cache_create.m	Mon Apr 04 18:14:56 2022 -0700
@@ -147,8 +147,8 @@
   endif
 
   ## create cache
-  func = @(s_) create_cache (__list_functions__ (s_));
-  cache = cellfun (func, directory, "UniformOutput", false);
+  f_lsfcn = @(dir) create_cache (__list_functions__ (dir));
+  cache = cellfun (f_lsfcn, directory, "UniformOutput", false);
 
   ## concatenate results
   cache = [cache{:}];
--- a/scripts/image/private/colorspace_conversion_input_check.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/image/private/colorspace_conversion_input_check.m	Mon Apr 04 18:14:56 2022 -0700
@@ -29,7 +29,7 @@
 ## by the complementary private function colorspace_conversion_revert()
 
 function [in_arg, sz, is_im, is_nd] ...
-            = colorspace_conversion_input_check (func, arg_name, in_arg)
+            = colorspace_conversion_input_check (fcn, arg_name, in_arg)
 
   cls = class (in_arg);
   sz = size (in_arg);
@@ -38,11 +38,11 @@
   if (! iscolormap (in_arg))
     if (! any (strcmp (cls, {"uint8", "int8", "int16", "uint16", ...
                              "single", "double"})))
-      error ("%s: %s of invalid data type '%s'", func, arg_name, cls);
+      error ("%s: %s of invalid data type '%s'", fcn, arg_name, cls);
     elseif (size (in_arg, 3) != 3)
-      error ("%s: %s must be a colormap or %s image", func, arg_name, arg_name);
+      error ("%s: %s must be a colormap or %s image", fcn, arg_name, arg_name);
     elseif (! isreal (in_arg) || ! isnumeric (in_arg))
-      error ("%s: %s must be numeric and real", func, arg_name);
+      error ("%s: %s must be numeric and real", fcn, arg_name);
     endif
     is_im = true;
 
@@ -60,7 +60,7 @@
       is_nd = true;
       in_arg = permute (in_arg, [1 2 4 3]);
     elseif (nd > 4)
-      error ("%s: invalid %s with more than 4 dimensions", func, arg_name);
+      error ("%s: invalid %s with more than 4 dimensions", fcn, arg_name);
     endif
     in_arg = reshape (in_arg, [numel(in_arg)/3 3]);
   else
--- a/scripts/image/private/imageIO.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/image/private/imageIO.m	Mon Apr 04 18:14:56 2022 -0700
@@ -35,8 +35,8 @@
 ##
 ## Usage:
 ##
-## func      - Function name to use on error message.
-## core_func - Function handle for the default function to use if we can't
+## fcn       - Function name to use on error message.
+## core_fcn  - Function handle for the default function to use if we can't
 ##             find the format in imformats.
 ## fieldname - Name of the field in the struct returned by imformats that
 ##             has the function to use.
@@ -46,7 +46,7 @@
 ## varargin  - Followed by all the OTHER arguments passed to imread and
 ##             imfinfo.
 
-function varargout = imageIO (func, core_func, fieldname, filename, varargin)
+function varargout = imageIO (fcn, core_fcn, fieldname, filename, varargin)
 
   ## First thing: figure out the filename and possibly download it.
   ## The first attempt is to try the filename and see if it exists.  If it
@@ -77,7 +77,7 @@
   endif
 
   if (isempty (fn))
-    error ([func ": unable to find file '" filename "'"]);
+    error ([fcn ": unable to find file '" filename "'"]);
   endif
 
   ## unwind_protect block because we may have a file to remove in the end
@@ -89,7 +89,7 @@
     ## try to guess the format from the file extension.  Finally, if
     ## we still don't know the format, we use the Octave core functions
     ## which is the same for all formats.
-    foo = []; # the function we will use
+    hfcn = []; # the function we will use
 
     ## We check if the call to imformats (ext) worked using
     ## "numfields (fmt) > 0 because when it fails, the returned
@@ -99,13 +99,13 @@
     if (! isempty (varargin) && ischar (varargin{1}))
       fmt = imformats (varargin{1});
       if (numfields (fmt) > 0)
-        foo = fmt.(fieldname);
+        hfcn = fmt.(fieldname);
         varargin(1) = []; # remove format name from arguments
       endif
     endif
 
     ## try extension from filename
-    if (isempty (foo))
+    if (isempty (hfcn))
       [~, ~, ext] = fileparts (fn);
       if (! isempty (ext))
         ## remove dot from extension
@@ -113,16 +113,16 @@
       endif
       fmt = imformats (ext);
       if (numfields (fmt) > 0)
-        foo = fmt.(fieldname);
+        hfcn = fmt.(fieldname);
       endif
     endif
 
     ## use the core function
-    if (isempty (foo))
-      foo = core_func;
+    if (isempty (hfcn))
+      hfcn = core_fcn;
     endif
 
-    [varargout{1:nargout}] = foo (fn, varargin{:});
+    [varargout{1:nargout}] = hfcn (fn, varargin{:});
 
   unwind_protect_cleanup
     if (file_2_delete)
--- a/scripts/legacy/@inline/argnames.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/legacy/@inline/argnames.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,18 +24,18 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{args} =} argnames (@var{fun})
+## @deftypefn {} {@var{args} =} argnames (@var{fcn})
 ## Return a cell array of character strings containing the names of the
-## arguments of the inline function @var{fun}.
+## arguments of the inline function @var{fcn}.
 ## @seealso{inline, formula, vectorize}
 ## @end deftypefn
 
-function args = argnames (obj)
+function args = argnames (fcn)
 
-  if (nargin < 1)
+  if (nargin != 1)
     print_usage ();
   endif
 
-  args = obj.args;
+  args = fcn.args;
 
 endfunction
--- a/scripts/legacy/@inline/vectorize.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/legacy/@inline/vectorize.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {} vectorize (@var{fun})
-## Create a vectorized version of the inline function @var{fun} by
+## @deftypefn {} {@var{vfcn} =} vectorize (@var{fcn})
+## Create a vectorized version of the inline function @var{fcn} by
 ## replacing all occurrences of @code{*}, @code{/}, etc., with
 ## @code{.*}, @code{./}, etc.
 ##
@@ -47,15 +47,15 @@
 ## The following function was translated directly from the original C++
 ## version.  Yes, it will be slow, but the use of inline functions is
 ## strongly discouraged anyway, and most expressions will probably be
-## short.  It may also be buggy.  Well, don't use this object!  Use
+## short.  It may also be buggy.  Well, don't use this function!  Use
 ## function handles instead!
 
-function fcn = vectorize (obj)
+function vfcn = vectorize (fcn)
 
   if (nargin < 1)
     print_usage ();
   endif
 
-  fcn = inline (__vectorize__ (obj.expr));
+  vfcn = inline (__vectorize__ (fcn.expr));
 
 endfunction
--- a/scripts/legacy/vectorize.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/legacy/vectorize.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,9 +24,9 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{vfun} =} vectorize (@var{fun})
+## @deftypefn {} {@var{vfcn} =} vectorize (@var{fcn})
 ## Create a vectorized version of the anonymous function or expression
-## @var{fun} by replacing all occurrences of @code{*}, @code{/}, etc.,
+## @var{fcn} by replacing all occurrences of @code{*}, @code{/}, etc.,
 ## with @code{.*}, @code{./}, etc.
 ##
 ## Note that the transformation is extremely simplistic.  Use of this
@@ -39,7 +39,7 @@
 ## anyway, and most expressions will probably be short.  It may also be
 ## buggy.  Well, don't use this function!
 
-function vfun = vectorize (fun)
+function vfcn = vectorize (fcn)
 
   persistent warned = false;
   if (! warned)
@@ -52,21 +52,21 @@
     print_usage ();
   endif
 
-  if (isa (fun, "function_handle"))
-    finfo = functions (fun);
+  if (isa (fcn, "function_handle"))
+    finfo = functions (fcn);
     if (! strcmp (finfo.type, "anonymous"))
-      error ("vectorize: FUN must be a string or anonymous function handle");
+      error ("vectorize: FCN must be a string or anonymous function handle");
     endif
     expr = finfo.function;
     idx = index (expr, ")");
     args = expr(1:idx);
     expr = expr(idx+1:end);
     new_expr = __vectorize__ (expr);
-    vfun = str2func ([args, new_expr]);
-  elseif (ischar (fun))
-    vfun = __vectorize__ (fun);
+    vfcn = str2func ([args, new_expr]);
+  elseif (ischar (fcn))
+    vfcn = __vectorize__ (fcn);
   else
-    error ("vectorize: FUN must be a string or anonymous function handle");
+    error ("vectorize: FCN must be a string or anonymous function handle");
   endif
 
 endfunction
@@ -92,4 +92,4 @@
 
 ## Test input validation
 %!error <Invalid call> vectorize ()
-%!error <FUN must be a string or anonymous function handle> vectorize (1)
+%!error <FCN must be a string or anonymous function handle> vectorize (1)
--- a/scripts/linear-algebra/normest1.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/linear-algebra/normest1.m	Mon Apr 04 18:14:56 2022 -0700
@@ -27,7 +27,7 @@
 ## @deftypefn  {} {@var{nest} =} normest1 (@var{A})
 ## @deftypefnx {} {@var{nest} =} normest1 (@var{A}, @var{t})
 ## @deftypefnx {} {@var{nest} =} normest1 (@var{A}, @var{t}, @var{x0})
-## @deftypefnx {} {@var{nest} =} normest1 (@var{Afun}, @var{t}, @var{x0}, @var{p1}, @var{p2}, @dots{})
+## @deftypefnx {} {@var{nest} =} normest1 (@var{Afcn}, @var{t}, @var{x0}, @var{p1}, @var{p2}, @dots{})
 ## @deftypefnx {} {[@var{nest}, @var{v}] =} normest1 (@var{A}, @dots{})
 ## @deftypefnx {} {[@var{nest}, @var{v}, @var{w}] =} normest1 (@var{A}, @dots{})
 ## @deftypefnx {} {[@var{nest}, @var{v}, @var{w}, @var{iter}] =} normest1 (@var{A}, @dots{})
@@ -39,7 +39,7 @@
 ## estimate of the 1-norm of a linear operator @var{A} when matrix-vector
 ## products @code{@var{A} * @var{x}} and @code{@var{A}' * @var{x}} can be
 ## cheaply computed.  In this case, instead of the matrix @var{A}, a function
-## @code{@var{Afun} (@var{flag}, @var{x})} is used; it must return:
+## @code{@var{Afcn} (@var{flag}, @var{x})} is used; it must return:
 ##
 ## @itemize @bullet
 ## @item
@@ -69,7 +69,7 @@
 ## @end example
 ##
 ## The parameters @var{p1}, @var{p2}, @dots{} are arguments of
-## @code{@var{Afun} (@var{flag}, @var{x}, @var{p1}, @var{p2}, @dots{})}.
+## @code{@var{Afcn} (@var{flag}, @var{x}, @var{p1}, @var{p2}, @dots{})}.
 ##
 ## The default value for @var{t} is 2.  The algorithm requires matrix-matrix
 ## products with sizes @var{n} x @var{n} and @var{n} x @var{t}.
@@ -134,8 +134,8 @@
     Aismat = false;
     Aisreal = A ("real", [], varargin{:});
     n = A ("dim", [], varargin{:});
-    Afun = @(x) A ("notransp", x, varargin{:});
-    A1fun = @(x) A ("transp", x, varargin{:});
+    Afcn = @(x) A ("notransp", x, varargin{:});
+    A1fcn = @(x) A ("transp", x, varargin{:});
   else
     error ("normest1: A must be a square matrix or a function handle");
   endif
@@ -175,7 +175,7 @@
     if (Aismat)
       Y = A * X;
     else
-      Y = Afun (X);
+      Y = Afcn (X);
     endif
     iter(2)++;
     [nest, j] = max (sum (abs (Y), 1), [], 2);
@@ -222,7 +222,7 @@
       if (Aismat)
         Z = A' * S;
       else
-        Z = A1fun (S);  # (3) of Algorithm 2.4
+        Z = A1fcn (S);  # (3) of Algorithm 2.4
       endif
       iter(2)++;
       h = max (abs (Z), [], 2);
@@ -254,7 +254,7 @@
 endfunction
 
 
-%!function z = afun_A (flag, x, A, n)
+%!function z = afcn_A (flag, x, A, n)
 %!  switch (flag)
 %!  case {"dim"}
 %!    z = n;
@@ -266,7 +266,7 @@
 %!    z = A * x;
 %!  endswitch
 %!endfunction
-%!function z = afun_A_P (flag, x, A, m)
+%!function z = afcn_A_P (flag, x, A, m)
 %!  switch (flag)
 %!  case "dim"
 %!    z = length (A);
@@ -297,17 +297,17 @@
 %! X = [1,1,-1;1,1,1;1,1,-1;1,-1,-1]/3;
 %! assert (normest1 (A, 3, X), norm (A, 1), 2 * eps);
 
-## test Afun
+## test Afcn
 %!test
 %! A = rand (10);
 %! n = length (A);
-%! Afun = @(flag, x) afun_A (flag, x, A, n);
-%! assert (normest1 (Afun), norm (A, 1), 2 * eps);
+%! Afcn = @(flag, x) afcn_A (flag, x, A, n);
+%! assert (normest1 (Afcn), norm (A, 1), 2 * eps);
 
-## test Afun with parameters
+## test Afcn with parameters
 %!test
 %! A = rand (10);
-%! assert (normest1 (@afun_A_P, [], [], A, 3), norm (A ^ 3, 1), 1000 * eps);
+%! assert (normest1 (@afcn_A_P, [], [], A, 3), norm (A ^ 3, 1), 1000 * eps);
 
 ## test output
 %!test
--- a/scripts/miscellaneous/inputname.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/miscellaneous/inputname.m	Mon Apr 04 18:14:56 2022 -0700
@@ -59,9 +59,9 @@
 ##
 ##   c = {'a', 'b'}
 ##   y = 1; z = 2;
-##   func (c, y, z)
+##   fcn (c, y, z)
 ##   % inputname() would return 'c', 'y', 'z' for the inputs.
-##   func (c{1}, y, z)
+##   fcn (c{1}, y, z)
 ##   % inputname() would return '', '', '' for the inputs.
 ##
 ## 3) If inputname is not called from a function, Matlab walks up the stack
--- a/scripts/miscellaneous/nthargout.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/miscellaneous/nthargout.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,13 +24,13 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{arg} =} nthargout (@var{n}, @var{func}, @dots{})
-## @deftypefnx {} {@var{arg} =} nthargout (@var{n}, @var{ntot}, @var{func}, @dots{})
+## @deftypefn  {} {@var{arg} =} nthargout (@var{n}, @var{fcn}, @dots{})
+## @deftypefnx {} {@var{arg} =} nthargout (@var{n}, @var{ntot}, @var{fcn}, @dots{})
 ## Return the @var{n}th output argument of the function specified by the
-## function handle or string @var{func}.
+## function handle or string @var{fcn}.
 ##
-## Any additional arguments are passed directly to @var{func}.  The total
-## number of arguments to call @var{func} with can be passed in @var{ntot}; by
+## Any additional arguments are passed directly to @var{fcn}.  The total
+## number of arguments to call @var{fcn} with can be passed in @var{ntot}; by
 ## default @var{ntot} is @var{n}.  The input @var{n} can also be a vector of
 ## indices of the output, in which case the output will be a cell array of the
 ## requested output arguments.
@@ -79,12 +79,12 @@
 
   if (is_function_handle (varargin{1}) || ischar (varargin{1}))
     ntot = max (n(:));
-    func = varargin{1};
+    fcn = varargin{1};
     args = varargin(2:end);
   elseif (isnumeric (varargin{1})
           && (is_function_handle (varargin{2}) || ischar (varargin{2})))
     ntot = varargin{1};
-    func = varargin{2};
+    fcn = varargin{2};
     args = varargin(3:end);
   else
     print_usage ();
@@ -97,7 +97,7 @@
   outargs = cell (1, ntot);
 
   try
-    [outargs{:}] = feval (func, args{:});
+    [outargs{:}] = feval (fcn, args{:});
     if (numel (n) > 1)
       arg = outargs(n);
     else
--- a/scripts/miscellaneous/private/display_info_file.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/miscellaneous/private/display_info_file.m	Mon Apr 04 18:14:56 2022 -0700
@@ -26,10 +26,10 @@
 ## news() and citation() are very much alike.  They both do the same thing,
 ## just for different files.  This function does all the work.
 
-function display_info_file (func, package, file)
+function display_info_file (fcn, package, file)
 
   if (! ischar (package))
-    error ("%s: PACKAGE must be a string", func);
+    error ("%s: PACKAGE must be a string", fcn);
   endif
 
   if (strcmpi (package, "octave"))
@@ -40,16 +40,16 @@
     names     = cellfun (@(x) x.name, installed, "UniformOutput", false);
     pos       = strcmpi (names, package);
     if (! any (pos))
-      error ("%s: package '%s' is not installed", func, package);
+      error ("%s: package '%s' is not installed", fcn, package);
     endif
     filepath = fullfile (installed{pos}.dir, "packinfo", file);
   endif
 
   if (! exist (filepath, "file"))
     if (strcmpi (package, "octave"))
-      error ("%s: broken installation -- unable to locate %s file", func, file);
+      error ("%s: broken installation -- unable to locate %s file", fcn, file);
     else
-      error ("%s: unable to locate %s file for package %s", func, file, package);
+      error ("%s: unable to locate %s file for package %s", fcn, file, package);
     endif
   endif
 
--- a/scripts/ode/decic.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/decic.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{y0_new}, @var{yp0_new}] =} decic (@var{fun}, @var{t0}, @var{y0}, @var{fixed_y0}, @var{yp0}, @var{fixed_yp0})
-## @deftypefnx {} {[@var{y0_new}, @var{yp0_new}] =} decic (@var{fun}, @var{t0}, @var{y0}, @var{fixed_y0}, @var{yp0}, @var{fixed_yp0}, @var{options})
+## @deftypefn  {} {[@var{y0_new}, @var{yp0_new}] =} decic (@var{fcn}, @var{t0}, @var{y0}, @var{fixed_y0}, @var{yp0}, @var{fixed_yp0})
+## @deftypefnx {} {[@var{y0_new}, @var{yp0_new}] =} decic (@var{fcn}, @var{t0}, @var{y0}, @var{fixed_y0}, @var{yp0}, @var{fixed_yp0}, @var{options})
 ## @deftypefnx {} {[@var{y0_new}, @var{yp0_new}, @var{resnorm}] =} decic (@dots{})
 ##
 ## Compute consistent implicit ODE initial conditions @var{y0_new} and
@@ -34,12 +34,12 @@
 ## A maximum of @code{length (@var{y0})} components between @var{fixed_y0} and
 ## @var{fixed_yp0} may be chosen as fixed values.
 ##
-## @var{fun} is a function handle.  The function must accept three inputs where
+## @var{fcn} is a function handle.  The function must accept three inputs where
 ## the first is time @var{t}, the second is a column vector of unknowns
 ## @var{y}, and the third is a column vector of unknowns @var{yp}.
 ##
 ## @var{t0} is the initial time such that
-## @code{@var{fun}(@var{t0}, @var{y0_new}, @var{yp0_new}) = 0}, specified as a
+## @code{@var{fcn}(@var{t0}, @var{y0_new}, @var{yp0_new}) = 0}, specified as a
 ## scalar.
 ##
 ## @var{y0} is a vector used as the initial guess for @var{y}.
@@ -93,7 +93,7 @@
 ## @seealso{ode15i, odeset}
 ## @end deftypefn
 
-function [y0_new, yp0_new, resnorm] = decic (fun, t0,
+function [y0_new, yp0_new, resnorm] = decic (fcn, t0,
                                              y0, fixed_y0, yp0, fixed_yp0,
                                              options)
 
@@ -102,9 +102,9 @@
   endif
 
   ## Validate inputs
-  if (! is_function_handle (fun))
+  if (! is_function_handle (fcn))
     error ("Octave:invalid-input-arg",
-           "decic: FUN must be a valid function handle");
+           "decic: FCN must be a valid function handle");
   endif
 
   if (! (isnumeric (t0) && isscalar (t0)))
@@ -171,7 +171,7 @@
   x0 = [y0(! fixed_y0); yp0(! fixed_yp0)];
   opt = optimset ("tolfun", TolFun, "tolx", TolX, "FinDiffType", "central");
   x = ...
-    fminunc (@(x) objective (x, t0, y0, fixed_y0, yp0, fixed_yp0, nl, nu, fun),
+    fminunc (@(x) objective (x, t0, y0, fixed_y0, yp0, fixed_yp0, nl, nu, fcn),
              x0, opt);
 
   y0_new  = y0;
@@ -180,18 +180,18 @@
   y0_new(! fixed_y0)   = x(1:nl);
   yp0_new(! fixed_yp0) = x(nl+1:nl+nu);
   if (isargout (3))
-    resnorm = fun (t0, y0_new, yp0_new);
+    resnorm = fcn (t0, y0_new, yp0_new);
   endif
 
 endfunction
 
-function res = objective (x, t0, y0, fixed_y0, yp0, fixed_yp0, nl, nu, fun)
+function res = objective (x, t0, y0, fixed_y0, yp0, fixed_yp0, nl, nu, fcn)
 
   y = y0;
   y(! fixed_y0) = x(1:nl);
   yp = yp0;
   yp(! fixed_yp0) = x(nl+1:nl+nu);
-  res = sqrt (sum (fun (t0, y, yp) .^ 2));
+  res = sqrt (sum (fcn (t0, y, yp) .^ 2));
 
 endfunction
 
@@ -221,7 +221,7 @@
 %!error <Invalid call> decic (1,2,3)
 %!error <Invalid call> decic (1,2,3,4)
 %!error <Invalid call> decic (1,2,3,4,5)
-%!error <FUN must be a valid function handle>
+%!error <FCN must be a valid function handle>
 %! decic (1, 0, [1; 0; 0], [1; 1; 0], [-1e-4; 1; 0], [0; 0; 0]);
 %!error <T0 must be a numeric scalar>
 %! decic (@rob, [1, 1], [1; 0; 0], [1; 1; 0], [-1e-4; 1; 0], [0; 0; 0]);
--- a/scripts/ode/ode15i.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/ode15i.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t}, @var{y}] =} ode15i (@var{fun}, @var{trange}, @var{y0}, @var{yp0})
-## @deftypefnx {} {[@var{t}, @var{y}] =} ode15i (@var{fun}, @var{trange}, @var{y0}, @var{yp0}, @var{ode_opt})
+## @deftypefn  {} {[@var{t}, @var{y}] =} ode15i (@var{fcn}, @var{trange}, @var{y0}, @var{yp0})
+## @deftypefnx {} {[@var{t}, @var{y}] =} ode15i (@var{fcn}, @var{trange}, @var{y0}, @var{yp0}, @var{ode_opt})
 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode15i (@dots{})
 ## @deftypefnx {} {@var{solution} =} ode15i (@dots{})
 ## @deftypefnx {} {} ode15i (@dots{})
@@ -35,7 +35,7 @@
 ## @code{ode15i} uses a variable step, variable order BDF (Backward
 ## Differentiation Formula) method that ranges from order 1 to 5.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function that defines the ODE: @code{0 = f(t,y,yp)}.  The
 ## function must accept three inputs where the first is time @var{t}, the
 ## second is the function value @var{y} (a column vector), and the third
@@ -96,7 +96,7 @@
 ## @seealso{decic, odeset, odeget}
 ## @end deftypefn
 
-function varargout = ode15i (fun, trange, y0, yp0, varargin)
+function varargout = ode15i (fcn, trange, y0, yp0, varargin)
 
   if (nargin < 4)
     print_usage ();
@@ -112,8 +112,8 @@
    options = odeset ();
   endif
 
-  ## Check fun, trange, y0, yp0
-  fun = check_default_input (fun, trange, solver, y0, yp0);
+  ## Check fcn, trange, y0, yp0
+  fcn = check_default_input (fcn, trange, solver, y0, yp0);
 
   if (! isempty (options.Jacobian))
     if (ischar (options.Jacobian))
@@ -172,7 +172,7 @@
   ## Jacobian
   options.havejac       = false;
   options.havejacsparse = false;
-  options.havejacfun    = false;
+  options.havejacfcn    = false;
 
   if (! isempty (options.Jacobian))
     options.havejac = true;
@@ -196,7 +196,7 @@
       endif
 
     elseif (is_function_handle (options.Jacobian))
-      options.havejacfun = true;
+      options.havejacfcn = true;
       if (nargin (options.Jacobian) == 3)
         [J1, J2] = options.Jacobian (trange(1), y0, yp0);
 
@@ -208,7 +208,7 @@
                  'ode15i: "Jacobian" function must evaluate to a real square matrix');
         endif
         if (issparse (J1) && issparse (J2))
-          options.havejacsparse = true;  # Jac is sparse fun
+          options.havejacsparse = true;  # Jac is sparse fcn
         endif
       else
         error ("Octave:invalid-input-arg",
@@ -256,7 +256,7 @@
   options.haveeventfunction = ! isempty (options.Events);
 
   ## 3 arguments in the event callback of ode15i
-  [t, y, te, ye, ie] = __ode15__ (fun, trange, y0, yp0, options, 3);
+  [t, y, te, ye, ie] = __ode15__ (fcn, trange, y0, yp0, options, 3);
 
   if (nargout == 2)
     varargout{1} = t;
@@ -286,7 +286,7 @@
 
 %!demo
 %! ## Solve Robertson's equations with ode15i
-%! fun = @(t, y, yp) [-(yp(1) + 0.04*y(1) - 1e4*y(2)*y(3));
+%! fcn = @(t, y, yp) [-(yp(1) + 0.04*y(1) - 1e4*y(2)*y(3));
 %!                    -(yp(2) - 0.04*y(1) + 1e4*y(2)*y(3) + 3e7*y(2)^2);
 %!                    y(1) + y(2) + y(3) - 1];
 %!
@@ -295,7 +295,7 @@
 %! yp0 = [-1e-4; 1e-4; 0];
 %! tspan = [0 4*logspace(-6, 6)];
 %!
-%! [t, y] = ode15i (fun, tspan, y0, yp0, opt);
+%! [t, y] = ode15i (fcn, tspan, y0, yp0, opt);
 %!
 %! y(:,2) = 1e4 * y(:, 2);
 %! figure (2);
@@ -442,19 +442,19 @@
 %! [t, y] = ode15i (@rob, [0, 100], [1; 0; 0], [-1e-4; 1e-4; 0], opt);
 %! assert ([t(end), y(end,:)], fref, 1e-3);
 
-## Jacobian fun dense
+## Jacobian fcn dense
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("Jacobian", @jacfundense);
 %! [t, y] = ode15i (@rob, [0, 100], [1; 0; 0], [-1e-4; 1e-4; 0], opt);
 %! assert ([t(end), y(end,:)], fref, 1e-3);
 
-## Jacobian fun dense as string
+## Jacobian fcn dense as string
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("Jacobian", "jacfundense");
 %! [t, y] = ode15i (@rob, [0, 100], [1; 0; 0], [-1e-4; 1e-4; 0], opt);
 %! assert ([t(end), y(end,:)], fref, 1e-3);
 
-## Jacobian fun sparse
+## Jacobian fcn sparse
 %!testif HAVE_SUNDIALS_SUNLINSOL_KLU
 %! opt = odeset ("Jacobian", @jacfunsparse, "AbsTol", 1e-7, "RelTol", 1e-7);
 %! [t, y] = ode15i (@rob, [0, 100], [1; 0; 0], [-1e-4; 1e-4; 0], opt);
@@ -532,7 +532,7 @@
 %! [0, 1], [1, 1], [1; 1]);
 %! assert (size (yout), [20, 2]);
 
-## Jacobian fun wrong dimension
+## Jacobian fcn wrong dimension
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("Jacobian", @jacwrong);
 %! fail ("[t, y] = ode15i (@rob, [0, 4e6], [1; 0; 0], [-1e-4; 1e-4; 0], opt)",
@@ -585,7 +585,7 @@
 %! fail ("[t, y] = ode15i (@rob, [0, 4e6], [1; 0; 0], [-1e-4; 1e-4; 0], opt)",
 %!       '"Jacobian" function "_5yVNhWVJWJn47RKnzxPsyb_" not found');
 
-%!function ydot = fun (t, y, yp)
+%!function ydot = fcn (t, y, yp)
 %!  ydot = [y - yp];
 %!endfunction
 
@@ -602,54 +602,54 @@
 %! fail ("ode15i (1, 1, 1)", "Invalid call to ode15i");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (1, 1, 1, 1)", "ode15i: fun must be of class:");
+%! fail ("ode15i (1, 1, 1, 1)", "ode15i: fcn must be of class:");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (1, 1, 1, 1, 1)", "ode15i: fun must be of class:");
+%! fail ("ode15i (1, 1, 1, 1, 1)", "ode15i: fcn must be of class:");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (1, 1, 1, 1, 1, 1)", "ode15i: fun must be of class:");
+%! fail ("ode15i (1, 1, 1, 1, 1, 1)", "ode15i: fcn must be of class:");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (@fun, 1, 1, 1)",
+%! fail ("ode15i (@fcn, 1, 1, 1)",
 %!       "ode15i: invalid value assigned to field 'trange'");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (@fun, [1, 1], 1, 1)",
+%! fail ("ode15i (@fcn, [1, 1], 1, 1)",
 %!       "ode15i: invalid value assigned to field 'trange'");
 
 %!testif HAVE_SUNDIALS
-%! fail ("ode15i (@fun, [1, 2], 1, [1, 2])",
+%! fail ("ode15i (@fcn, [1, 2], 1, [1, 2])",
 %!       "ode15i: y0 must have 2 elements");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("RelTol", "_5yVNhWVJWJn47RKnzxPsyb_");
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       "ode15i: RelTol must be of class:");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("RelTol", [1, 2]);
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       "ode15i: RelTol must be scalar");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("RelTol", -2);
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       "ode15i: RelTol must be positive");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("AbsTol", "_5yVNhWVJWJn47RKnzxPsyb_");
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       "ode15i: AbsTol must be of class:");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("AbsTol", -1);
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       "ode15i: AbsTol must be positive");
 
 %!testif HAVE_SUNDIALS
 %! opt = odeset ("AbsTol", [1, 1, 1]);
-%! fail ("[t, y] = ode15i (@fun, [0, 2], 2, 2, opt)",
+%! fail ("[t, y] = ode15i (@fcn, [0, 2], 2, 2, opt)",
 %!       'ode15i: invalid value assigned to field "AbsTol"');
 
 %!testif HAVE_SUNDIALS
--- a/scripts/ode/ode15s.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/ode15s.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t}, @var{y}] =} ode15s (@var{fun}, @var{trange}, @var{y0})
-## @deftypefnx {} {[@var{t}, @var{y}] =} ode15s (@var{fun}, @var{trange}, @var{y0}, @var{ode_opt})
+## @deftypefn  {} {[@var{t}, @var{y}] =} ode15s (@var{fcn}, @var{trange}, @var{y0})
+## @deftypefnx {} {[@var{t}, @var{y}] =} ode15s (@var{fcn}, @var{trange}, @var{y0}, @var{ode_opt})
 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode15s (@dots{})
 ## @deftypefnx {} {@var{solution} =} ode15s (@dots{})
 ## @deftypefnx {} {} ode15s (@dots{})
@@ -35,7 +35,7 @@
 ## @code{ode15s} uses a variable step, variable order BDF (Backward
 ## Differentiation Formula) method that ranges from order 1 to 5.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function that defines the ODE: @code{y' = f(t,y)}.  The function
 ## must accept two inputs where the first is time @var{t} and the second is a
 ## column vector of unknowns @var{y}.
@@ -91,15 +91,15 @@
 ## @seealso{decic, odeset, odeget, ode23, ode45}
 ## @end deftypefn
 
-function varargout = ode15s (fun, trange, y0, varargin)
+function varargout = ode15s (fcn, trange, y0, varargin)
 
   if (nargin < 3)
     print_usage ();
   endif
 
   solver = "ode15s";
-  ## Check fun, trange, y0, yp0
-  fun = check_default_input (fun, trange, solver, y0);
+  ## Check fcn, trange, y0, yp0
+  fcn = check_default_input (fcn, trange, solver, y0);
 
   n = numel (y0);
 
@@ -175,14 +175,14 @@
                           classes, attributes, solver);
 
   ## Mass
-  options.havemassfun    = false;
+  options.havemassfcn    = false;
   options.havestatedep   = false;
   options.havetimedep    = false;
   options.havemasssparse = false;
 
   if (! isempty (options.Mass))
     if (is_function_handle (options.Mass))
-      options.havemassfun = true;
+      options.havemassfcn = true;
       if (nargin (options.Mass) == 2)
         options.havestatedep = true;
         M = options.Mass (trange(1), y0);
@@ -216,19 +216,19 @@
   ## Jacobian
   options.havejac       = false;
   options.havejacsparse = false;
-  options.havejacfun    = false;
+  options.havejacfcn    = false;
 
   if (! isempty (options.Jacobian))
     options.havejac = true;
     if (is_function_handle (options.Jacobian))
-      options.havejacfun = true;
+      options.havejacfcn = true;
       if (nargin (options.Jacobian) == 2)
         A = options.Jacobian (trange(1), y0);
         if (! issquare (A) || rows (A) != n || ! isnumeric (A) || ! isreal (A))
           error ("Octave:invalid-input-arg",
                  'ode15s: "Jacobian" function must evaluate to a real square matrix');
         endif
-        options.havejacsparse = issparse (A);  # Jac is sparse fun
+        options.havejacsparse = issparse (A);  # Jac is sparse fcn
       else
         error ("Octave:invalid-input-arg",
                'ode15s: invalid value assigned to field "Jacobian"');
@@ -261,15 +261,15 @@
     options.havejacsparse = false;
   endif
 
-  ## If Mass or Jacobian is fun, then new Jacobian is fun
+  ## If Mass or Jacobian is fcn, then new Jacobian is fcn
   if (options.havejac)
-    if (options.havejacfun || options.havetimedep)
-      options.Jacobian = @ (t, y, yp) wrapjacfun (t, y, yp,
+    if (options.havejacfcn || options.havetimedep)
+      options.Jacobian = @ (t, y, yp) wrapjacfcn (t, y, yp,
                                                   options.Jacobian,
                                                   options.Mass,
                                                   options.havetimedep,
-                                                  options.havejacfun);
-      options.havejacfun = true;
+                                                  options.havejacfcn);
+      options.havejacfcn = true;
     else   # All matrices are constant
       if (! isempty (options.Mass))
         options.Jacobian = {[- options.Jacobian], [options.Mass]};
@@ -324,7 +324,7 @@
   [t, y, te, ye, ie] = __ode15__ (@ (t, y, yp) wrap (t, y, yp, options.Mass,
                                                      options.havetimedep,
                                                      options.havestatedep,
-                                                     fun),
+                                                     fcn),
                                   trange, y0, yp0, options, 2);
 
   if (nargout == 2)
@@ -352,24 +352,24 @@
 
 endfunction
 
-function res = wrap (t, y, yp, Mass, havetimedep, havestatedep, fun)
+function res = wrap (t, y, yp, Mass, havetimedep, havestatedep, fcn)
 
   if (! isempty (Mass) && havestatedep)
-    res = Mass (t, y) * yp - fun (t, y);
+    res = Mass (t, y) * yp - fcn (t, y);
   elseif (! isempty (Mass) && havetimedep)
-    res = Mass (t) * yp - fun (t, y);
+    res = Mass (t) * yp - fcn (t, y);
   elseif (! isempty (Mass))
-    res = Mass * yp - fun (t, y);
+    res = Mass * yp - fcn (t, y);
   else
-    res = yp - fun (t, y);
+    res = yp - fcn (t, y);
   endif
 
 endfunction
 
-function [jac, jact] = wrapjacfun (t, y, yp, Jac, Mass,
-                                   havetimedep, havejacfun)
+function [jac, jact] = wrapjacfcn (t, y, yp, Jac, Mass,
+                                   havetimedep, havejacfcn)
 
-  if (havejacfun)
+  if (havejacfcn)
     jac = - Jac (t, y);
   else
     jac = - Jac;
@@ -388,7 +388,7 @@
 
 %!demo
 %! ## Solve Robertson's equations with ode15s
-%! fun = @ (t, y) [-0.04*y(1) + 1e4*y(2).*y(3);
+%! fcn = @ (t, y) [-0.04*y(1) + 1e4*y(2).*y(3);
 %!                  0.04*y(1) - 1e4*y(2).*y(3) - 3e7*y(2).^2;
 %!                  y(1) + y(2) + y(3) - 1];
 %!
@@ -399,7 +399,7 @@
 %! options = odeset ("RelTol", 1e-4, "AbsTol", [1e-6, 1e-10, 1e-6],
 %!                   "MStateDependence", "none", "Mass", M);
 %!
-%! [t, y] = ode15s (fun, tspan, y0, options);
+%! [t, y] = ode15s (fcn, tspan, y0, options);
 %!
 %! y(:,2) = 1e4 * y(:,2);
 %! figure (2);
--- a/scripts/ode/ode23.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/ode23.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t}, @var{y}] =} ode23 (@var{fun}, @var{trange}, @var{init})
-## @deftypefnx {} {[@var{t}, @var{y}] =} ode23 (@var{fun}, @var{trange}, @var{init}, @var{ode_opt})
+## @deftypefn  {} {[@var{t}, @var{y}] =} ode23 (@var{fcn}, @var{trange}, @var{init})
+## @deftypefnx {} {[@var{t}, @var{y}] =} ode23 (@var{fcn}, @var{trange}, @var{init}, @var{ode_opt})
 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode23 (@dots{})
 ## @deftypefnx {} {@var{solution} =} ode23 (@dots{})
 ## @deftypefnx {} {} ode23 (@dots{})
@@ -33,7 +33,7 @@
 ## Solve a set of non-stiff Ordinary Differential Equations (non-stiff ODEs)
 ## with the well known explicit @nospell{Bogacki-Shampine} method of order 3.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function that defines the ODE: @code{y' = f(t,y)}.  The function
 ## must accept two inputs where the first is time @var{t} and the second is a
 ## column vector of unknowns @var{y}.
@@ -92,7 +92,7 @@
 ## @seealso{odeset, odeget, ode45, ode15s}
 ## @end deftypefn
 
-function varargout = ode23 (fun, trange, init, varargin)
+function varargout = ode23 (fcn, trange, init, varargin)
 
   if (nargin < 3)
     print_usage ();
@@ -103,7 +103,7 @@
 
   if (nargin >= 4)
     if (! isstruct (varargin{1}))
-      ## varargin{1:len} are parameters for fun
+      ## varargin{1:len} are parameters for fcn
       odeopts = odeset ();
       funarguments = varargin;
     elseif (numel (varargin) > 1)
@@ -142,16 +142,16 @@
   endif
   init = init(:);
 
-  if (ischar (fun))
-    if (! exist (fun))
+  if (ischar (fcn))
+    if (! exist (fcn))
       error ("Octave:invalid-input-arg",
-             ['ode23: function "' fun '" not found']);
+             ['ode23: function "' fcn '" not found']);
     endif
-    fun = str2func (fun);
+    fcn = str2func (fcn);
   endif
-  if (! is_function_handle (fun))
+  if (! is_function_handle (fcn))
     error ("Octave:invalid-input-arg",
-           "ode23: FUN must be a valid function handle");
+           "ode23: FCN must be a valid function handle");
   endif
 
   ## Start preprocessing, have a look which options are set in odeopts,
@@ -194,7 +194,7 @@
 
   if (isempty (odeopts.InitialStep))
     odeopts.InitialStep = odeopts.direction * ...
-                          starting_stepsize (order, fun, trange(1), init,
+                          starting_stepsize (order, fcn, trange(1), init,
                                              odeopts.AbsTol, odeopts.RelTol,
                                              strcmpi (odeopts.NormControl, "on"),
                                              odeopts.funarguments);
@@ -220,12 +220,12 @@
     if (! strcmp (odeopts.MStateDependence, "none"))
       ## constant mass matrices have already been handled
       mass = @(t,x) odeopts.Mass (t, x, odeopts.funarguments{:});
-      fun = @(t,x) mass (t, x, odeopts.funarguments{:}) ...
-                   \ fun (t, x, odeopts.funarguments{:});
+      fcn = @(t,x) mass (t, x, odeopts.funarguments{:}) ...
+                   \ fcn (t, x, odeopts.funarguments{:});
     else
       mass = @(t) odeopts.Mass (t, odeopts.funarguments{:});
-      fun = @(t,x) mass (t, odeopts.funarguments{:}) ...
-                   \ fun (t, x, odeopts.funarguments{:});
+      fcn = @(t,x) mass (t, odeopts.funarguments{:}) ...
+                   \ fcn (t, x, odeopts.funarguments{:});
     endif
   endif
 
@@ -239,7 +239,7 @@
   endif
 
   solution = integrate_adaptive (@runge_kutta_23,
-                                 order, fun, trange, init, odeopts);
+                                 order, fcn, trange, init, odeopts);
 
   ## Postprocessing, do whatever when terminating integration algorithm
   if (odeopts.haveoutputfunction)  # Cleanup plotter
@@ -511,4 +511,4 @@
 %!error <invalid time span>  ode23 (@fpol, [1 1], [3 15 1])
 %!error <INIT must be a numeric> ode23 (@fpol, [0 25], {[3 15 1]})
 %!error <INIT must be a .* vector> ode23 (@fpol, [0 25], [3 15 1; 3 15 1])
-%!error <FUN must be a valid function handle> ode23 (1, [0 25], [3 15 1])
+%!error <FCN must be a valid function handle> ode23 (1, [0 25], [3 15 1])
--- a/scripts/ode/ode23s.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/ode23s.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t}, @var{y}] =} ode23s (@var{fun}, @var{trange}, @var{init})
-## @deftypefnx {} {[@var{t}, @var{y}] =} ode23s (@var{fun}, @var{trange}, @var{init}, @var{ode_opt})
+## @deftypefn  {} {[@var{t}, @var{y}] =} ode23s (@var{fcn}, @var{trange}, @var{init})
+## @deftypefnx {} {[@var{t}, @var{y}] =} ode23s (@var{fcn}, @var{trange}, @var{init}, @var{ode_opt})
 ## @deftypefnx {} {[@var{t}, @var{y}] =} ode23s (@dots{}, @var{par1}, @var{par2}, @dots{})
 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode23s (@dots{})
 ## @deftypefnx {} {@var{solution} =} ode23s (@dots{})
@@ -33,7 +33,7 @@
 ## Solve a set of stiff Ordinary Differential Equations (stiff ODEs) with a
 ## @nospell{Rosenbrock} method of order (2,3).
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function that defines the ODE: @code{M y' = f(t,y)}.  The
 ## function must accept two inputs where the first is time @var{t} and the
 ## second is a column vector of unknowns @var{y}.  @var{M} is a constant mass
@@ -93,7 +93,7 @@
 ## @seealso{odeset, daspk, dassl}
 ## @end deftypefn
 
-function varargout = ode23s (fun, trange, init, varargin)
+function varargout = ode23s (fcn, trange, init, varargin)
 
   if (nargin < 3)
     print_usage ();
@@ -104,7 +104,7 @@
 
   if (nargin >= 4)
     if (! isstruct (varargin{1}))
-      ## varargin{1:len} are parameters for fun
+      ## varargin{1:len} are parameters for fcn
       odeopts = odeset ();
       funarguments = varargin;
     elseif (numel (varargin) > 1)
@@ -143,16 +143,16 @@
   endif
   init = init(:);
 
-  if (ischar (fun))
-    if (! exist (fun))
+  if (ischar (fcn))
+    if (! exist (fcn))
       error ("Octave:invalid-input-arg",
-             ['ode23s: function "' fun '" not found']);
+             ['ode23s: function "' fcn '" not found']);
     endif
-    fun = str2func (fun);
+    fcn = str2func (fcn);
   endif
-  if (! is_function_handle (fun))
+  if (! is_function_handle (fcn))
     error ("Octave:invalid-input-arg",
-           "ode23s: FUN must be a valid function handle");
+           "ode23s: FCN must be a valid function handle");
   endif
 
   ## Start preprocessing, have a look which options are set in odeopts,
@@ -184,7 +184,7 @@
 
   if (isempty (odeopts.InitialStep))
     odeopts.InitialStep = odeopts.direction * ...
-                          starting_stepsize (order, fun, trange(1), init,
+                          starting_stepsize (order, fcn, trange(1), init,
                                              odeopts.AbsTol, odeopts.RelTol,
                                              strcmpi (odeopts.NormControl,
                                                       "on"),
@@ -215,7 +215,7 @@
   endif
 
   solution = integrate_adaptive (@runge_kutta_23s, ...
-                                 order, fun, trange, init, odeopts);
+                                 order, fcn, trange, init, odeopts);
 
   ## Postprocessing, do whatever when terminating integration algorithm
   if (odeopts.haveoutputfunction)  # Cleanup plotter
@@ -277,45 +277,45 @@
 
 %!demo
 %! ## Demo function: stiff Van Der Pol equation
-%! fun = @(t,y) [y(2); 10*(1-y(1)^2)*y(2)-y(1)];
+%! fcn = @(t,y) [y(2); 10*(1-y(1)^2)*y(2)-y(1)];
 %! ## Calling ode23s method
 %! tic ()
-%! [vt, vy] = ode23s (fun, [0 20], [2 0]);
+%! [vt, vy] = ode23s (fcn, [0 20], [2 0]);
 %! toc ()
 %! ## Plotting the result
 %! plot (vt,vy(:,1),'-o');
 
 %!demo
 %! ## Demo function: stiff Van Der Pol equation
-%! fun = @(t,y) [y(2); 10*(1-y(1)^2)*y(2)-y(1)];
+%! fcn = @(t,y) [y(2); 10*(1-y(1)^2)*y(2)-y(1)];
 %! ## Calling ode23s method
 %! odeopts = odeset ("Jacobian", @(t,y) [0 1; -20*y(1)*y(2)-1, 10*(1-y(1)^2)],
 %!                   "InitialStep", 1e-3)
 %! tic ()
-%! [vt, vy] = ode23s (fun, [0 20], [2 0], odeopts);
+%! [vt, vy] = ode23s (fcn, [0 20], [2 0], odeopts);
 %! toc ()
 %! ## Plotting the result
 %! plot (vt,vy(:,1),'-o');
 
 %!demo
 %! ## Demo function: stiff Van Der Pol equation
-%! fun = @(t,y) [y(2); 100*(1-y(1)^2)*y(2)-y(1)];
+%! fcn = @(t,y) [y(2); 100*(1-y(1)^2)*y(2)-y(1)];
 %! ## Calling ode23s method
 %! odeopts = odeset ("InitialStep", 1e-4);
 %! tic ()
-%! [vt, vy] = ode23s (fun, [0 200], [2 0]);
+%! [vt, vy] = ode23s (fcn, [0 200], [2 0]);
 %! toc ()
 %! ## Plotting the result
 %! plot (vt,vy(:,1),'-o');
 
 %!demo
 %! ## Demo function: stiff Van Der Pol equation
-%! fun = @(t,y) [y(2); 100*(1-y(1)^2)*y(2)-y(1)];
+%! fcn = @(t,y) [y(2); 100*(1-y(1)^2)*y(2)-y(1)];
 %! ## Calling ode23s method
 %! odeopts = odeset ("Jacobian", @(t,y) [0 1; -200*y(1)*y(2)-1, 100*(1-y(1)^2)],
 %!                   "InitialStep", 1e-4);
 %! tic ()
-%! [vt, vy] = ode23s (fun, [0 200], [2 0], odeopts);
+%! [vt, vy] = ode23s (fcn, [0 200], [2 0], odeopts);
 %! toc ()
 %! ## Plotting the result
 %! plot (vt,vy(:,1),'-o');
@@ -509,4 +509,4 @@
 %!error <invalid time span>  ode23s (@fpol, [1 1], [3 15 1])
 %!error <INIT must be a numeric> ode23s (@fpol, [0 25], {[3 15 1]})
 %!error <INIT must be a .* vector> ode23s (@fpol, [0 25], [3 15 1; 3 15 1])
-%!error <FUN must be a valid function handle> ode23s (1, [0 25], [3 15 1])
+%!error <FCN must be a valid function handle> ode23s (1, [0 25], [3 15 1])
--- a/scripts/ode/ode45.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/ode45.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,8 +24,8 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t}, @var{y}] =} ode45 (@var{fun}, @var{trange}, @var{init})
-## @deftypefnx {} {[@var{t}, @var{y}] =} ode45 (@var{fun}, @var{trange}, @var{init}, @var{ode_opt})
+## @deftypefn  {} {[@var{t}, @var{y}] =} ode45 (@var{fcn}, @var{trange}, @var{init})
+## @deftypefnx {} {[@var{t}, @var{y}] =} ode45 (@var{fcn}, @var{trange}, @var{init}, @var{ode_opt})
 ## @deftypefnx {} {[@var{t}, @var{y}, @var{te}, @var{ye}, @var{ie}] =} ode45 (@dots{})
 ## @deftypefnx {} {@var{solution} =} ode45 (@dots{})
 ## @deftypefnx {} {} ode45 (@dots{})
@@ -33,7 +33,7 @@
 ## Solve a set of non-stiff Ordinary Differential Equations (non-stiff ODEs)
 ## with the well known explicit @nospell{Dormand-Prince} method of order 4.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function that defines the ODE: @code{y' = f(t,y)}.  The function
 ## must accept two inputs where the first is time @var{t} and the second is a
 ## column vector of unknowns @var{y}.
@@ -89,7 +89,7 @@
 ## @seealso{odeset, odeget, ode23, ode15s}
 ## @end deftypefn
 
-function varargout = ode45 (fun, trange, init, varargin)
+function varargout = ode45 (fcn, trange, init, varargin)
 
   if (nargin < 3)
     print_usage ();
@@ -100,7 +100,7 @@
 
   if (nargin >= 4)
     if (! isstruct (varargin{1}))
-      ## varargin{1:len} are parameters for fun
+      ## varargin{1:len} are parameters for fcn
       odeopts = odeset ();
       funarguments = varargin;
     elseif (numel (varargin) > 1)
@@ -139,16 +139,16 @@
   endif
   init = init(:);
 
-  if (ischar (fun))
-    if (! exist (fun))
+  if (ischar (fcn))
+    if (! exist (fcn))
       error ("Octave:invalid-input-arg",
-             ['ode45: function "' fun '" not found']);
+             ['ode45: function "' fcn '" not found']);
     endif
-    fun = str2func (fun);
+    fcn = str2func (fcn);
   endif
-  if (! is_function_handle (fun))
+  if (! is_function_handle (fcn))
     error ("Octave:invalid-input-arg",
-           "ode45: FUN must be a valid function handle");
+           "ode45: FCN must be a valid function handle");
   endif
 
   ## Start preprocessing, have a look which options are set in odeopts,
@@ -194,7 +194,7 @@
 
   if (isempty (odeopts.InitialStep))
     odeopts.InitialStep = odeopts.direction * ...
-                          starting_stepsize (order, fun, trange(1), init,
+                          starting_stepsize (order, fcn, trange(1), init,
                                              odeopts.AbsTol, odeopts.RelTol,
                                              strcmpi (odeopts.NormControl, "on"),
                                              odeopts.funarguments);
@@ -220,12 +220,12 @@
     if (! strcmp (odeopts.MStateDependence, "none"))
       ### constant mass matrices have already
       mass = @(t,x) odeopts.Mass (t, x, odeopts.funarguments{:});
-      fun = @(t,x) mass (t, x, odeopts.funarguments{:}) ...
-                   \ fun (t, x, odeopts.funarguments{:});
+      fcn = @(t,x) mass (t, x, odeopts.funarguments{:}) ...
+                   \ fcn (t, x, odeopts.funarguments{:});
     else
       mass = @(t) odeopts.Mass (t, odeopts.funarguments{:});
-      fun = @(t,x) mass (t, odeopts.funarguments{:}) ...
-                   \ fun (t, x, odeopts.funarguments{:});
+      fcn = @(t,x) mass (t, odeopts.funarguments{:}) ...
+                   \ fcn (t, x, odeopts.funarguments{:});
     endif
   endif
 
@@ -239,7 +239,7 @@
   endif
 
   solution = integrate_adaptive (@runge_kutta_45_dorpri,
-                                 order, fun, trange, init, odeopts);
+                                 order, fcn, trange, init, odeopts);
 
   ## Postprocessing, do whatever when terminating integration algorithm
   if (odeopts.haveoutputfunction)  # Cleanup plotter
@@ -523,4 +523,4 @@
 %!error <invalid time span> ode45 (@fpol, [1 1], [3 15 1])
 %!error <INIT must be a numeric> ode45 (@fpol, [0 25], {[3 15 1]})
 %!error <INIT must be a .* vector> ode45 (@fpol, [0 25], [3 15 1; 3 15 1])
-%!error <FUN must be a valid function handle> ode45 (1, [0 25], [3 15 1])
+%!error <FCN must be a valid function handle> ode45 (1, [0 25], [3 15 1])
--- a/scripts/ode/odeset.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/odeset.m	Mon Apr 04 18:14:56 2022 -0700
@@ -139,7 +139,7 @@
 ## Print solver statistics after simulation.
 ##
 ## @item @code{Vectorized}: @{@qcode{"off"}@} | @qcode{"on"}
-## Specify whether @code{odefun} can be passed multiple values of the
+## Specify whether @code{odefcn} can be passed multiple values of the
 ## state at once.
 ##
 ## @end table
--- a/scripts/ode/private/check_default_input.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/check_default_input.m	Mon Apr 04 18:14:56 2022 -0700
@@ -23,30 +23,30 @@
 ##
 ########################################################################
 
-function fun = check_default_input (fun, trange, solver, y0, yp0);
+function fcn = check_default_input (fcn, trange, solver, y0, yp0);
 
   if (nargin < 4)
     print_usage ();
   endif
 
-  ## Check fun
-  validateattributes (fun, {"function_handle", "char"}, {}, solver, "fun");
+  ## Check fcn
+  validateattributes (fcn, {"function_handle", "char"}, {}, solver, "fcn");
 
-  if (! (nargin (fun) == nargin - 2))
+  if (! (nargin (fcn) == nargin - 2))
     error ("Octave:invalid-input-arg",
-           [solver ": invalid value assigned to field 'fun'"]);
+           [solver ": invalid value assigned to field 'fcn'"]);
   endif
 
-  if (ischar (fun))
-    if (! exist (fun))
+  if (ischar (fcn))
+    if (! exist (fcn))
       error ("Octave:invalid-input-arg",
-             [solver ": function '" fun "' not found"]);
+             [solver ": function '" fcn "' not found"]);
     endif
-    fun = str2func (fun);
+    fcn = str2func (fcn);
   endif
-  if (! is_function_handle (fun))
+  if (! is_function_handle (fcn))
     error ("Octave:invalid-input-arg",
-           [solver ": invalid value assigned to field '" fun "'"]);
+           [solver ": invalid value assigned to field '" fcn "'"]);
   endif
 
   ## Check trange
@@ -74,10 +74,10 @@
     endif
     yp0 = yp0(:);
 
-    n = numel (feval (fun, trange(1), y0, yp0));
+    n = numel (feval (fcn, trange(1), y0, yp0));
     validateattributes (yp0, {"float"}, {"numel", n}, solver, "yp0");
   else
-    n = numel (feval (fun, trange (1), y0));
+    n = numel (feval (fcn, trange(1), y0));
   endif
 
   validateattributes (y0, {"float"}, {"numel", n}, solver, "y0");
--- a/scripts/ode/private/integrate_adaptive.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/integrate_adaptive.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,7 +24,7 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{solution} =} integrate_adaptive (@var{@@stepper}, @var{order}, @var{@@func}, @var{tspan}, @var{x0}, @var{options})
+## @deftypefn {} {@var{solution} =} integrate_adaptive (@var{@@stepper}, @var{order}, @var{@@fcn}, @var{tspan}, @var{x0}, @var{options})
 ##
 ## This function file can be called by an ODE solver function in order to
 ## integrate the set of ODEs on the interval @var{[t0, t1]} with an adaptive
@@ -68,7 +68,7 @@
 ##
 ## @end deftypefn
 
-function solution = integrate_adaptive (stepper, order, func, tspan, x0,
+function solution = integrate_adaptive (stepper, order, fcn, tspan, x0,
                                         options)
 
   fixed_times = numel (tspan) > 2;
@@ -79,7 +79,7 @@
   ## Get first initial timestep
   dt = options.InitialStep;
   if (isempty (dt))
-    dt = starting_stepsize (order, func, t, x,
+    dt = starting_stepsize (order, fcn, t, x,
                             options.AbsTol, options.RelTol,
                             strcmp (options.NormControl, "on"),
                             options.funarguments);
@@ -133,7 +133,7 @@
     ## Compute integration step from t_old to t_new = t_old + dt
     [t_new, options.comp] = kahan (t_old, options.comp, dt);
     [t_new, x_new, x_est, new_k_vals] = ...
-      stepper (func, t_old, x_old, dt, options, k_vals, t_new);
+      stepper (fcn, t_old, x_old, dt, options, k_vals, t_new);
 
     solution.cntcycles += 1;
 
@@ -163,7 +163,7 @@
           iout = max (t_caught);
           x(:, t_caught) = ...
             runge_kutta_interpolate (order, [t_old t_new], [x_old x_new],
-                                     t(t_caught), new_k_vals, dt, func,
+                                     t(t_caught), new_k_vals, dt, fcn,
                                      options.funarguments);
 
           istep += 1;
--- a/scripts/ode/private/ode_event_handler.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/ode_event_handler.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,9 +24,9 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{retval} =} ode_event_handler (@var{@@evtfun}, @var{t}, @var{y}, @var{flag}, @var{par1}, @var{par2}, @dots{})
+## @deftypefn {} {@var{retval} =} ode_event_handler (@var{@@evtfcn}, @var{t}, @var{y}, @var{flag}, @var{par1}, @var{par2}, @dots{})
 ##
-## Return the solution of the event function (@var{@@evtfun}) which is
+## Return the solution of the event function (@var{@@evtfcn}) which is
 ## specified in the form of a function handle.
 ##
 ## The second input argument @var{t} is a scalar double and specifies the time
@@ -60,7 +60,7 @@
 ## necessary to call it directly.
 ## @end deftypefn
 
-function retval = ode_event_handler (evtfun, t, y, flag = "", varargin)
+function retval = ode_event_handler (evtfcn, t, y, flag = "", varargin)
 
   ## No error handling has been implemented in this function to achieve
   ## the highest performance possible.
@@ -85,9 +85,9 @@
     ## Process the event, i.e.,
     ## find the zero crossings for either a rising or falling edge
     if (! iscell (y))
-      inpargs = {evtfun, t, y};
+      inpargs = {evtfcn, t, y};
     else
-      inpargs = {evtfun, t, y{1}, y{2}};
+      inpargs = {evtfcn, t, y{1}, y{2}};
       y = y{1};  # Delete cell element 2
     endif
     if (nargin > 4)
@@ -141,9 +141,9 @@
     firstrun = true;
 
     if (! iscell (y))
-      inpargs = {evtfun, t, y};
+      inpargs = {evtfcn, t, y};
     else
-      inpargs = {evtfun, t, y{1}, y{2}};
+      inpargs = {evtfcn, t, y{1}, y{2}};
       y = y{1};  # Delete cell element 2
     endif
     if (nargin > 4)
--- a/scripts/ode/private/runge_kutta_23.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/runge_kutta_23.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,10 +24,10 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fun}, @var{t}, @var{x}, @var{dt})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
+## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fcn}, @var{t}, @var{x}, @var{dt})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23 (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}] =} runge_kutta_23 (@dots{})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}, @var{k_vals_out}] =} runge_kutta_23 (@dots{})
 ##
@@ -36,7 +36,7 @@
 ## method of third order.  For the definition of this method see
 ## @url{http://en.wikipedia.org/wiki/List_of_Runge%E2%80%93Kutta_methods}.
 ##
-## @var{fun} is a function handle that defines the ODE: @code{y' = f(tau,y)}.
+## @var{fcn} is a function handle that defines the ODE: @code{y' = f(tau,y)}.
 ## The function must accept two inputs where the first is time @var{tau} and
 ## the second is a column vector of unknowns @var{y}.
 ##
@@ -49,7 +49,7 @@
 ## The optional fourth argument @var{options} specifies options for the ODE
 ## solver.  It is a structure generated by @code{odeset}.  In particular it
 ## contains the field @var{funarguments} with the optional arguments to be used
-## in the evaluation of @var{fun}.
+## in the evaluation of @var{fcn}.
 ##
 ## The optional fifth argument @var{k_vals_in} contains the Runge-Kutta
 ## evaluations of the previous step to use in a FSAL scheme.
@@ -66,7 +66,7 @@
 ## @seealso{runge_kutta_45_dorpri}
 ## @end deftypefn
 
-function [t_next, x_next, x_est, k] = runge_kutta_23 (fun, t, x, dt,
+function [t_next, x_next, x_est, k] = runge_kutta_23 (fcn, t, x, dt,
                                                       options = [],
                                                       k_vals = [],
                                                       t_next = t + dt)
@@ -92,11 +92,11 @@
   if (! isempty (k_vals))    # k values from previous step are passed
     k(:,1) = k_vals(:,end);  # FSAL property
   else
-    k(:,1) = feval (fun, t, x, args{:});
+    k(:,1) = feval (fcn, t, x, args{:});
   endif
 
-  k(:,2) = feval (fun, s(2), x + k(:,1) * aa(2, 1).', args{:});
-  k(:,3) = feval (fun, s(3), x + k(:,2) * aa(3, 2).', args{:});
+  k(:,2) = feval (fcn, s(2), x + k(:,1) * aa(2, 1).', args{:});
+  k(:,3) = feval (fcn, s(3), x + k(:,2) * aa(3, 2).', args{:});
 
   ## compute new time and new values for the unknowns
   ## t_next = t + dt;
@@ -105,7 +105,7 @@
   ## if the estimation of the error is required
   if (nargout >= 3)
     ## new solution to be compared with the previous one
-    k(:,4) = feval (fun, t_next, x_next, args{:});
+    k(:,4) = feval (fcn, t_next, x_next, args{:});
     cc_prime = dt * c_prime;
     x_est = x + k * cc_prime(:);
   endif
--- a/scripts/ode/private/runge_kutta_23s.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/runge_kutta_23s.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,10 +24,10 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fun}, @var{t}, @var{x}, @var{dt})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
+## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fcn}, @var{t}, @var{x}, @var{dt})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_23s (@var{fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}] =} runge_kutta_23s (@dots{})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}, @var{k_vals_out}] =} runge_kutta_23s (@dots{})
 ##
@@ -50,7 +50,7 @@
 ## The optional fourth argument @var{options} specifies options for the ODE
 ## solver.  It is a structure generated by @code{odeset}.  In particular it
 ## contains the field @var{funarguments} with the optional arguments to be used
-## in the evaluation of @var{fun}.
+## in the evaluation of @var{fcn}.
 ##
 ## The optional fifth argument @var{k_vals_in} contains the Runge-Kutta
 ## evaluations of the previous step to use in a FSAL scheme.
@@ -67,7 +67,7 @@
 ## @seealso{runge_kutta_23}
 ## @end deftypefn
 
-function [t_next, x_next, x_est, k] = runge_kutta_23s (fun, t, x, dt,
+function [t_next, x_next, x_est, k] = runge_kutta_23s (fcn, t, x, dt,
                                                        options = [],
                                                        k_vals = [],
                                                        t_next = t + dt)
@@ -83,14 +83,14 @@
     args = {};
   endif
 
-  jacfun = false;
+  jacfcn = false;
   jacmat = false;
   if (! isempty (options.Jacobian))
     if (ischar (options.Jacobian))
-      jacfun = true;
+      jacfcn = true;
       jac = str2fun (options.Jacobian);
     elseif (is_function_handle (options.Jacobian))
-      jacfun = true;
+      jacfcn = true;
       jac = options.Jacobian;
     elseif (ismatrix (options.Jacobian))
       jacmat = true;
@@ -110,15 +110,15 @@
   ## Jacobian matrix, dfxpdp
   if (jacmat)
     J = jac;
-  elseif (jacfun)
+  elseif (jacfcn)
     J = jac (t, x);
   elseif (! jacpat)
-    J = __dfxpdp__ (x, @(a) feval (fun, t, a, args{:}), options.RelTol);
+    J = __dfxpdp__ (x, @(a) feval (fcn, t, a, args{:}), options.RelTol);
   elseif (jacpat)
-    J = __dfxpdp__ (x, @(a) feval (fun, t, a, args{:}), options.RelTol, pattern);
+    J = __dfxpdp__ (x, @(a) feval (fcn, t, a, args{:}), options.RelTol, pattern);
   endif
 
-  T = (feval (fun, t + .1 * dt, x) - feval (fun, t, x, args{:})) / (.1 * dt);
+  T = (feval (fcn, t + .1 * dt, x) - feval (fcn, t, x, args{:})) / (.1 * dt);
 
   ## Wolfbrandt coefficient
   if (isempty (options.Mass))
@@ -135,13 +135,13 @@
   endif
 
   ## compute the slopes
-  F(:,1) = feval (fun, t, x, args{:});
+  F(:,1) = feval (fcn, t, x, args{:});
   if (issparse (W))
     k(:,1) = Qw * (Uw \ (Lw \ (Pw * (Rw \ (F(:,1) + dt*d*T)))));
   else
     k(:,1) = Uw \ (Lw \ (Pw * (F(:,1) + dt*d*T)));
   endif
-  F(:,2) = feval (fun, t+a*dt, x+a*dt*k(:,1), args{:});
+  F(:,2) = feval (fcn, t+a*dt, x+a*dt*k(:,1), args{:});
   if (issparse (W))
     k(:,2) = Uw * (Uw \ (Lw \ (Pw * (Rw \ (F(:,2) - M*k(:,1)))))) + k(:,1);
   else
@@ -153,7 +153,7 @@
 
   if (nargout >= 3)
     ## 3rd order, needed in error formula
-    F(:,3) = feval (fun, t+dt, x_next, args{:});
+    F(:,3) = feval (fcn, t+dt, x_next, args{:});
     if (issparse (W))
       k(:,3) = Qw * (Uw \ (Lw \ (Pw * (Rw \ (F(:,3) - e32 * (M*k(:,2) - F(:,2)) - 2 * (M*k(:,1) - F(:,1)) + dt*d*T)))));
     else
@@ -170,12 +170,12 @@
 endfunction
 
 
-function prt = __dfxpdp__ (p, func, rtol, pattern)
+function prt = __dfxpdp__ (p, fcn, rtol, pattern)
 
   ## This subfunction was copied 2011 from the OF "optim" package
   ## "inst/private/__dfdp__.m".
 
-  f = func (p)(:);
+  f = fcn (p)(:);
   m = numel (f);
   n = numel (p);
 
@@ -193,9 +193,9 @@
     prt = pattern;  # initialize Jacobian
     for j = find (any (pattern, 1))
       ps(j) = p1(j);
-      tp1 = func (ps);
+      tp1 = fcn (ps);
       ps(j) = p2(j);
-      tp2 = func (ps);
+      tp2 = fcn (ps);
       pattern_nnz = find (pattern(:, j));
       prt(pattern_nnz, j) = (tp1(pattern_nnz) - tp2(pattern_nnz)) / absdel(j);
       ps(j) = p(j);
@@ -204,9 +204,9 @@
     prt = zeros (m, n); # initialize Jacobian
     for j = 1:n
       ps(j) = p1(j);
-      tp1 = func (ps);
+      tp1 = fcn (ps);
       ps(j) = p2(j);
-      tp2 = func (ps);
+      tp2 = fcn (ps);
       prt(:, j) = (tp1(:) - tp2(:)) / absdel(j);
       ps(j) = p(j);
     endfor
--- a/scripts/ode/private/runge_kutta_45_dorpri.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/runge_kutta_45_dorpri.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,10 +24,10 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fun}, @var{t}, @var{x}, @var{dt})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fun}, @var{t}, @var{x}, @var{dt}, @var{options})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
-## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fun}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
+## @deftypefn  {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fcn}, @var{t}, @var{x}, @var{dt})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fcn}, @var{t}, @var{x}, @var{dt}, @var{options})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals})
+## @deftypefnx {} {[@var{t_next}, @var{x_next}] =} runge_kutta_45_dorpri (@var{@@fcn}, @var{t}, @var{x}, @var{dt}, @var{options}, @var{k_vals}, @var{t_next})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}] =} runge_kutta_45_dorpri (@dots{})
 ## @deftypefnx {} {[@var{t_next}, @var{x_next}, @var{x_est}, @var{k_vals_out}] =} runge_kutta_45_dorpri (@dots{})
 ##
@@ -67,7 +67,7 @@
 ## to use in an FSAL scheme or for dense output.
 ## @end deftypefn
 
-function [t_next, x_next, x_est, k] = runge_kutta_45_dorpri (fun, t, x, dt,
+function [t_next, x_next, x_est, k] = runge_kutta_45_dorpri (fcn, t, x, dt,
                                                              options = [],
                                                              k_vals = [],
                                                              t_next = t + dt)
@@ -100,14 +100,14 @@
   if (! isempty (k_vals))    # k values from previous step are passed
     k(:,1) = k_vals(:,end);  # FSAL property
   else
-    k(:,1) = feval (fun, t, x, args{:});
+    k(:,1) = feval (fcn, t, x, args{:});
   endif
 
-  k(:,2) = feval (fun, s(2), x + k(:,1)   * aa(2, 1).'  , args{:});
-  k(:,3) = feval (fun, s(3), x + k(:,1:2) * aa(3, 1:2).', args{:});
-  k(:,4) = feval (fun, s(4), x + k(:,1:3) * aa(4, 1:3).', args{:});
-  k(:,5) = feval (fun, s(5), x + k(:,1:4) * aa(5, 1:4).', args{:});
-  k(:,6) = feval (fun, s(6), x + k(:,1:5) * aa(6, 1:5).', args{:});
+  k(:,2) = feval (fcn, s(2), x + k(:,1)   * aa(2, 1).'  , args{:});
+  k(:,3) = feval (fcn, s(3), x + k(:,1:2) * aa(3, 1:2).', args{:});
+  k(:,4) = feval (fcn, s(4), x + k(:,1:3) * aa(4, 1:3).', args{:});
+  k(:,5) = feval (fcn, s(5), x + k(:,1:4) * aa(5, 1:4).', args{:});
+  k(:,6) = feval (fcn, s(6), x + k(:,1:5) * aa(6, 1:5).', args{:});
 
   ## compute new time and new values for the unknowns
   ## t_next = t + dt;
@@ -116,7 +116,7 @@
   ## if the estimation of the error is required
   if (nargout >= 3)
     ## new solution to be compared with the previous one
-    k(:,7) = feval (fun, t_next, x_next, args{:});
+    k(:,7) = feval (fcn, t_next, x_next, args{:});
     cc_prime = dt * c_prime;
     x_est = x + k * cc_prime(:);
   endif
--- a/scripts/ode/private/runge_kutta_interpolate.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/runge_kutta_interpolate.m	Mon Apr 04 18:14:56 2022 -0700
@@ -23,7 +23,7 @@
 ##
 ########################################################################
 
-function u_interp = runge_kutta_interpolate (order, z, u, t, k_vals, dt, func, args)
+function u_interp = runge_kutta_interpolate (order, z, u, t, k_vals, dt, fcn, args)
 
   switch (order)
 
@@ -34,7 +34,7 @@
       if (! isempty (k_vals))
         der = k_vals(:,1);
       else
-        der = feval (func, z(1) , u(:,1), args);
+        der = feval (fcn, z(1) , u(:,1), args);
       endif
       u_interp = quadratic_interpolation (z, u, der, t);
 
@@ -48,8 +48,8 @@
     otherwise
       warning (["High order interpolation not yet implemented: ", ...
                 "using cubic interpolation instead"]);
-      der(:,1) = feval (func, z(1), u(:,1), args);
-      der(:,2) = feval (func, z(2), u(:,2), args);
+      der(:,1) = feval (fcn, z(1), u(:,1), args);
+      der(:,2) = feval (fcn, z(2), u(:,2), args);
       u_interp = hermite_cubic_interpolation (z, u, der, t);
 
   endswitch
--- a/scripts/ode/private/starting_stepsize.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/ode/private/starting_stepsize.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,12 +24,12 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {@var{h} =} starting_stepsize (@var{order}, @var{func}, @var{t0}, @var{x0}, @var{AbsTol}, @var{RelTol}, @var{normcontrol}, @var{args})
+## @deftypefn {} {@var{h} =} starting_stepsize (@var{order}, @var{fcn}, @var{t0}, @var{x0}, @var{AbsTol}, @var{RelTol}, @var{normcontrol}, @var{args})
 ##
 ## Determine a good initial timestep for an ODE solver of order @var{order}
 ## using the algorithm described in reference [1].
 ##
-## The input argument @var{func}, is the function describing the differential
+## The input argument @var{fcn}, is the function describing the differential
 ## equations, @var{t0} is the initial time, and @var{x0} is the initial
 ## condition.  @var{AbsTol} and @var{RelTol} are the absolute and relative
 ## tolerance on the ODE integration taken from an ode options structure.
@@ -42,7 +42,7 @@
 ##
 ## @seealso{odepkg}
 
-function h = starting_stepsize (order, func, t0, x0,
+function h = starting_stepsize (order, fcn, t0, x0,
                                 AbsTol, RelTol, normcontrol,
                                 args = {})
 
@@ -50,7 +50,7 @@
   d0 = AbsRel_norm (x0, x0, AbsTol, RelTol, normcontrol);
 
   ## compute norm of the function evaluated at initial conditions
-  y = func (t0, x0, args{:});
+  y = fcn (t0, x0, args{:});
   if (iscell (y))
     y = y{1};
   endif
@@ -66,7 +66,7 @@
   x1 = x0 + h0 * y;
 
   ## approximate the derivative norm
-  yh = func (t0+h0, x1, args{:});
+  yh = fcn (t0+h0, x1, args{:});
   if (iscell (yh))
     yh = yh{1};
   endif
--- a/scripts/optimization/__all_opts__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/__all_opts__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -55,8 +55,8 @@
     for i = 1:nargin
       try
         opts = optimset (varargin{i});
-        fn = fieldnames (opts).';
-        names = [names, fn];
+        fcn = fieldnames (opts).';
+        names = [names, fcn];
       catch
         ## throw the error as a warning.
         warning (lasterr ());
--- a/scripts/optimization/fminbnd.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/fminbnd.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,12 +24,12 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{x} =} fminbnd (@var{fun}, @var{a}, @var{b})
-## @deftypefnx {} {@var{x} =} fminbnd (@var{fun}, @var{a}, @var{b}, @var{options})
+## @deftypefn  {} {@var{x} =} fminbnd (@var{fcn}, @var{a}, @var{b})
+## @deftypefnx {} {@var{x} =} fminbnd (@var{fcn}, @var{a}, @var{b}, @var{options})
 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{info}, @var{output}] =} fminbnd (@dots{})
 ## Find a minimum point of a univariate function.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.
 ##
 ## The starting interval is specified by @var{a} (left boundary) and @var{b}
@@ -88,10 +88,10 @@
 ## PKG_ADD: ## Discard result to avoid polluting workspace with ans at startup.
 ## PKG_ADD: [~] = __all_opts__ ("fminbnd");
 
-function [x, fval, info, output] = fminbnd (fun, a, b, options = struct ())
+function [x, fval, info, output] = fminbnd (fcn, a, b, options = struct ())
 
   ## Get default options if requested.
-  if (nargin == 1 && ischar (fun) && strcmp (fun, "defaults"))
+  if (nargin == 1 && ischar (fcn) && strcmp (fcn, "defaults"))
     x = struct ("Display", "notify", "FunValCheck", "off",
                 "MaxFunEvals", 500, "MaxIter", 500,
                 "OutputFcn", [], "TolX", 1e-4);
@@ -107,8 +107,8 @@
            "fminbnd: the lower bound cannot be greater than the upper one");
   endif
 
-  if (ischar (fun))
-    fun = str2func (fun);
+  if (ischar (fcn))
+    fcn = str2func (fcn);
   endif
 
   displ = optimget (options, "Display", "notify");
@@ -119,8 +119,8 @@
   maxfev = optimget (options, "MaxFunEvals", 500);
 
   if (funvalchk)
-    ## Replace fun with a guarded version.
-    fun = @(x) guarded_eval (fun, x);
+    ## Replace fcn with a guarded version.
+    fcn = @(x) guarded_eval (fcn, x);
   endif
 
   ## The default exit flag if exceeded number of iterations.
@@ -132,7 +132,7 @@
   v = a + c*(b-a);
   w = x = v;
   e = 0;
-  fv = fw = fval = fun (x);
+  fv = fw = fval = fcn (x);
   nfev += 1;
 
   if (isa (a, "single") || isa (b, "single") || isa (fval, "single"))
@@ -199,7 +199,7 @@
 
     ## f must not be evaluated too close to x.
     u = x + max (abs (d), tol) * (sign (d) + (d == 0));
-    fu = fun (u);
+    fu = fcn (u);
 
     niter += 1;
 
@@ -275,9 +275,9 @@
 endfunction
 
 ## A helper function that evaluates a function and checks for bad results.
-function fx = guarded_eval (fun, x)
+function fx = guarded_eval (fcn, x)
 
-  fx = fun (x);
+  fx = fcn (x);
   fx = fx(1);
   if (! isreal (fx))
     error ("Octave:fmindbnd:notreal", "fminbnd: non-real value encountered");
@@ -289,7 +289,7 @@
 
 ## A hack for printing a formatted table
 function print_formatted_table (table)
-  printf ("\n Func-count     x          f(x)         Procedure\n");
+  printf ("\n Fcn-count     x          f(x)         Procedure\n");
   for row=table
     printf ("%5.5s        %7.7s    %8.8s\t%s\n",
             int2str (row.funccount), num2str (row.x,"%.5f"),
--- a/scripts/optimization/fminsearch.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/fminsearch.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,15 +24,15 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{x} =} fminsearch (@var{fun}, @var{x0})
-## @deftypefnx {} {@var{x} =} fminsearch (@var{fun}, @var{x0}, @var{options})
+## @deftypefn  {} {@var{x} =} fminsearch (@var{fcn}, @var{x0})
+## @deftypefnx {} {@var{x} =} fminsearch (@var{fcn}, @var{x0}, @var{options})
 ## @deftypefnx {} {@var{x} =} fminsearch (@var{problem})
 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{exitflag}, @var{output}] =} fminsearch (@dots{})
 ##
 ## Find a value of @var{x} which minimizes the multi-variable function
-## @var{fun}.
+## @var{fcn}.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.
 ##
 ## The search begins at the point @var{x0} and iterates using the
@@ -103,7 +103,7 @@
 ##
 ## The fourth output is a structure @var{output} containing runtime
 ## about the algorithm.  Fields in the structure are @code{funcCount}
-## containing the number of function calls to @var{fun}, @code{iterations}
+## containing the number of function calls to @var{fcn}, @code{iterations}
 ## containing the number of iteration steps, @code{algorithm} with the name of
 ## the search algorithm (always:
 ## @nospell{@qcode{"Nelder-Mead simplex direct search"}}), and @code{message}
@@ -146,7 +146,7 @@
     if (! isstruct (problem))
       error ("fminsearch: PROBLEM must be a structure");
     endif
-    fun = problem.objective;
+    fcn = problem.objective;
     x0 = problem.x0;
     if (! strcmp (problem.solver, "fminsearch"))
       error ('fminsearch: problem.solver must be set to "fminsearch"');
@@ -157,7 +157,7 @@
       options = [];
     endif
   elseif (nargin > 1)
-    fun = varargin{1};
+    fcn = varargin{1};
     x0 = varargin{2};
     if (nargin > 2)
       options = varargin{3};
@@ -168,25 +168,25 @@
     endif
   endif
 
-  if (ischar (fun))
-    fun = str2func (fun);
+  if (ischar (fcn))
+    fcn = str2func (fcn);
   endif
 
   if (isempty (options))
     options = struct ();
   endif
 
-  [x, exitflag, output] = nmsmax (fun, x0, options, varargin{:});
+  [x, exitflag, output] = nmsmax (fcn, x0, options, varargin{:});
 
   if (isargout (2))
-    fval = feval (fun, x);
+    fval = feval (fcn, x);
   endif
 
 endfunction
 
 ## NMSMAX  Nelder-Mead simplex method for direct search optimization.
-##        [x, fmax, nf] = NMSMAX(FUN, x0, STOPIT, SAVIT) attempts to
-##        maximize the function FUN, using the starting vector x0.
+##        [x, fmax, nf] = NMSMAX(FCN, x0, STOPIT, SAVIT) attempts to
+##        maximize the function FCN, using the starting vector x0.
 ##        The Nelder-Mead direct search method is used.
 ##        Output arguments:
 ##               x    = vector yielding largest function value found,
@@ -210,8 +210,8 @@
 ##        'SAVE SAVIT x fmax nf' is executed after each inner iteration.
 ##        NB: x0 can be a matrix.  In the output argument, in SAVIT saves,
 ##            and in function calls, x has the same shape as x0.
-##        NMSMAX(fun, x0, STOPIT, SAVIT, P1, P2,...) allows additional
-##        arguments to be passed to fun, via feval(fun,x,P1,P2,...).
+##        NMSMAX(FCN, x0, STOPIT, SAVIT, P1, P2,...) allows additional
+##        arguments to be passed to FCN, via feval(FCN,x,P1,P2,...).
 ## References:
 ## N. J. Higham, Optimization by direct search in matrix computations,
 ##    SIAM J. Matrix Anal. Appl, 14(2): 317-333, 1993.
@@ -270,16 +270,16 @@
 
 endfunction
 
-function [x, exitflag, output] = nmsmax (fun, x, options, varargin)
+function [x, exitflag, output] = nmsmax (fcn, x, options, varargin)
 
   [stopit, savit, dirn, trace, tol, maxiter, tol_f, outfcn] = ...
                                                     parse_options (options, x);
 
   if (strcmpi (optimget (options, "FunValCheck", "off"), "on"))
     ## Replace fcn with a guarded version.
-    fun = @(x) guarded_eval (fun, x, varargin{:});
+    fcn = @(x) guarded_eval (fcn, x, varargin{:});
   else
-    fun = @(x) real (fun (x, varargin{:}));
+    fcn = @(x) real (fcn (x, varargin{:}));
   endif
 
   x0 = x(:);  # Work with column vector internally.
@@ -288,7 +288,7 @@
   V = [zeros(n,1) eye(n)];
   f = zeros (n+1,1);
   V(:,1) = x0;
-  f(1) = dirn * fun (x);
+  f(1) = dirn * fcn (x);
   fmax_old = f(1);
   nf = 1;
 
@@ -308,7 +308,7 @@
     for j = 2:n+1
       V(j-1,j) = x0(j-1) + alpha(1);
       x(:) = V(:,j);
-      f(j) = dirn * fun (x);
+      f(j) = dirn * fcn (x);
     endfor
   else
     ## Right-angled simplex based on co-ordinate axes.
@@ -316,7 +316,7 @@
     for j = 2:n+1
       V(:,j) = x0 + alpha(j)*V(:,j);
       x(:) = V(:,j);
-      f(j) = dirn * fun (x);
+      f(j) = dirn * fcn (x);
     endfor
   endif
   nf += n;
@@ -419,14 +419,14 @@
     vbar = (sum (V(:,1:n)')/n)';  # Mean value
     vr = (1 + alpha)*vbar - alpha*V(:,n+1);
     x(:) = vr;
-    fr = dirn * fun (x);
+    fr = dirn * fcn (x);
     nf += 1;
     vk = vr;  fk = fr; how = "reflect";
     if (fr > f(n))
       if (fr > f(1))
         ve = gamma*vr + (1-gamma)*vbar;
         x(:) = ve;
-        fe = dirn * fun (x);
+        fe = dirn * fcn (x);
         nf += 1;
         if (fe > f(1))
           vk = ve;
@@ -443,7 +443,7 @@
       endif
       vc = beta*vt + (1-beta)*vbar;
       x(:) = vc;
-      fc = dirn * fun (x);
+      fc = dirn * fcn (x);
       nf += 1;
       if (fc > f(n))
         vk = vc; fk = fc;
@@ -452,12 +452,12 @@
         for j = 2:n
           V(:,j) = (V(:,1) + V(:,j))/2;
           x(:) = V(:,j);
-          f(j) = dirn * fun (x);
+          f(j) = dirn * fcn (x);
         endfor
         nf += n-1;
         vk = (V(:,1) + V(:,n+1))/2;
         x(:) = vk;
-        fk = dirn * fun (x);
+        fk = dirn * fcn (x);
         nf += 1;
         how = "shrink";
       endif
@@ -497,9 +497,9 @@
 endfunction
 
 ## A helper function that evaluates a function and checks for bad results.
-function y = guarded_eval (fun, x, varargin)
+function y = guarded_eval (fcn, x, varargin)
 
-  y = fun (x, varargin{:});
+  y = fcn (x, varargin{:});
 
   if (! (isreal (y)))
     error ("fminsearch:notreal", "fminsearch: non-real value encountered");
--- a/scripts/optimization/fminunc.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/fminunc.m	Mon Apr 04 18:14:56 2022 -0700
@@ -37,7 +37,7 @@
 ## @code{fminunc} attempts to determine a vector @var{x} such that
 ## @code{@var{fcn} (@var{x})} is a local minimum.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.  @var{fcn} should accept a vector (array)
 ## defining the unknown variables, and return the objective function value,
 ## optionally with gradient.
@@ -413,12 +413,12 @@
 endfunction
 
 ## A helper function that evaluates a function and checks for bad results.
-function [fx, gx] = guarded_eval (fun, x)
+function [fx, gx] = guarded_eval (fcn, x)
 
   if (nargout > 1)
-    [fx, gx] = fun (x);
+    [fx, gx] = fcn (x);
   else
-    fx = fun (x);
+    fx = fcn (x);
     gx = [];
   endif
 
--- a/scripts/optimization/fsolve.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/fsolve.m	Mon Apr 04 18:14:56 2022 -0700
@@ -32,7 +32,7 @@
 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{info}, @var{output}, @var{fjac}] =} fsolve (@dots{})
 ## Solve a system of nonlinear equations defined by the function @var{fcn}.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.  @var{fcn} should accept a vector (array)
 ## defining the unknown variables, and return a vector of left-hand sides of
 ## the equations.  Right-hand sides are defined to be zeros.  In other words,
@@ -152,7 +152,7 @@
 ## recent vector.  A short example how this can be achieved follows:
 ##
 ## @example
-## function [fval, fjac] = user_func (x, optimvalues, state)
+## function [fval, fjac] = user_fcn (x, optimvalues, state)
 ## persistent sav = [], sav0 = [];
 ## if (nargin == 1)
 ##   ## evaluation call
@@ -173,7 +173,7 @@
 ##
 ## ## @dots{}
 ##
-## fsolve (@@user_func, x0, optimset ("OutputFcn", @@user_func, @dots{}))
+## fsolve (@@user_fcn, x0, optimset ("OutputFcn", @@user_fcn, @dots{}))
 ## @end example
 ## @seealso{fzero, optimset}
 ## @end deftypefn
@@ -509,12 +509,12 @@
 endfunction
 
 ## A helper function that evaluates a function and checks for bad results.
-function [fx, jx] = guarded_eval (fun, x, complexeqn)
+function [fx, jx] = guarded_eval (fcn, x, complexeqn)
 
   if (nargout > 1)
-    [fx, jx] = fun (x);
+    [fx, jx] = fcn (x);
   else
-    fx = fun (x);
+    fx = fcn (x);
     jx = [];
   endif
 
@@ -675,7 +675,7 @@
 %! assert (norm (c - c_opt, Inf) < tol);
 %! assert (norm (fval) < norm (noise));
 
-%!function y = cfun (x)
+%!function y = cfcn (x)
 %!  y(1) = (1+i)*x(1)^2 - (1-i)*x(2) - 2;
 %!  y(2) = sqrt (x(1)*x(2)) - (1-2i)*x(3) + (3-4i);
 %!  y(3) = x(1) * x(2) - x(3)^2 + (3+2i);
@@ -685,7 +685,7 @@
 %! x_opt = [-1+i, 1-i, 2+i];
 %! x = [i, 1, 1+i];
 %!
-%! [x, f, info] = fsolve (@cfun, x, optimset ("ComplexEqn", "on"));
+%! [x, f, info] = fsolve (@cfcn, x, optimset ("ComplexEqn", "on"));
 %! tol = 1e-5;
 %! assert (norm (f) < tol);
 %! assert (norm (x - x_opt, Inf) < tol);
--- a/scripts/optimization/fzero.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/fzero.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,14 +24,14 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{x} =} fzero (@var{fun}, @var{x0})
-## @deftypefnx {} {@var{x} =} fzero (@var{fun}, @var{x0}, @var{options})
+## @deftypefn  {} {@var{x} =} fzero (@var{fcn}, @var{x0})
+## @deftypefnx {} {@var{x} =} fzero (@var{fcn}, @var{x0}, @var{options})
 ## @deftypefnx {} {[@var{x}, @var{fval}] =} fzero (@dots{})
 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{info}] =} fzero (@dots{})
 ## @deftypefnx {} {[@var{x}, @var{fval}, @var{info}, @var{output}] =} fzero (@dots{})
 ## Find a zero of a univariate function.
 ##
-## @var{fun} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.
 ##
 ## @var{x0} should be a two-element vector specifying two points which
@@ -40,7 +40,7 @@
 ## following must hold
 ##
 ## @example
-## sign (@var{fun}(@var{x0}(1))) * sign (@var{fun}(@var{x0}(2))) <= 0
+## sign (@var{fcn}(@var{x0}(1))) * sign (@var{fcn}(@var{x0}(2))) <= 0
 ## @end example
 ##
 ## If @var{x0} is a single scalar then several nearby and distant values are
@@ -127,10 +127,10 @@
 ## PKG_ADD: ## Discard result to avoid polluting workspace with ans at startup.
 ## PKG_ADD: [~] = __all_opts__ ("fzero");
 
-function [x, fval, info, output] = fzero (fun, x0, options = struct ())
+function [x, fval, info, output] = fzero (fcn, x0, options = struct ())
 
   ## Get default options if requested.
-  if (nargin == 1 && ischar (fun) && strcmp (fun, "defaults"))
+  if (nargin == 1 && ischar (fcn) && strcmp (fcn, "defaults"))
     x = struct ("Display", "notify", "FunValCheck", "off",
                 "MaxFunEvals", Inf, "MaxIter", Inf,
                 "OutputFcn", [], "TolX", eps);
@@ -141,8 +141,8 @@
     print_usage ();
   endif
 
-  if (ischar (fun))
-    fun = str2func (fun);
+  if (ischar (fcn))
+    fcn = str2func (fcn);
   endif
 
   displev = optimget (options, "Display", "notify");
@@ -166,8 +166,8 @@
   mu = 0.5;
 
   if (funvalchk)
-    ## Replace fun with a guarded version.
-    fun = @(x) guarded_eval (fun, x);
+    ## Replace fcn with a guarded version.
+    fcn = @(x) guarded_eval (fcn, x);
   endif
 
   info = 0;  # The default exit flag if number of iterations exceeded.
@@ -178,20 +178,20 @@
 
   ## Prepare...
   a = x0(1);
-  fa = fun (a);
+  fa = fcn (a);
   nfev = 1;
   if (length (x0) > 1)
     b = x0(2);
-    fb = fun (b);
+    fb = fcn (b);
     nfev += 1;
   else
     ## Try to find a value for b which brackets a zero-crossing
     if (displev == 1)
       printf ( ...
         "\nSearch for an interval around %g containing a sign change:\n", a);
-      printf (" Func-count    a          f(a)             b          ");
+      printf (" Fcn-count    a          f(a)             b          ");
       printf ("f(b)        Procedure\n");
-      fmt_str = " %4d   %13.6g %13.6g %13.6g %13.6g   %s\n";
+      fmt_str = " %4d  %13.6g %13.6g %13.6g %13.6g   %s\n";
     endif
 
     ## For very small values, switch to absolute rather than relative search
@@ -206,7 +206,7 @@
     ## Search in an ever-widening range around the initial point.
     for srch = [-.01 +.025 -.05 +.10 -.25 +.50 -1 +2.5 -5 +10 -50 +100 -500 +1000]
       b = aa + aa*srch;
-      fb = fun (b);
+      fb = fcn (b);
       nfev += 1;
       if (displev == 1)
         printf (fmt_str, nfev, a, fa, b, fb, "search");
@@ -233,8 +233,8 @@
 
   if (displev == 1)
     printf ("\nSearch for a zero in the interval [%g, %g]:\n", a, b);
-    disp (" Func-count    x          f(x)             Procedure");
-    fmt_str = " %4d   %13.6g %13.6g        %s\n";
+    disp (" Fcn-count    x          f(x)             Procedure");
+    fmt_str = " %4d  %13.6g %13.6g        %s\n";
   endif
 
   slope0 = (fb - fa) / (b - a);
@@ -345,7 +345,7 @@
 
     ## Calculate new point.
     x = c;
-    fval = fc = fun (c);
+    fval = fc = fcn (c);
     niter += 1; nfev += 1;
     if (displev == 1)
       printf (fmt_str, nfev, x, fc, "interpolation");
@@ -440,9 +440,9 @@
 endfunction
 
 ## A helper function that evaluates a function and checks for bad results.
-function fx = guarded_eval (fun, x)
+function fx = guarded_eval (fcn, x)
 
-  fx = fun (x);
+  fx = fcn (x);
   fx = fx(1);
   if (! isreal (fx))
     error ("Octave:fzero:notreal", "fzero: non-real value encountered");
--- a/scripts/optimization/sqp.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/optimization/sqp.m	Mon Apr 04 18:14:56 2022 -0700
@@ -213,14 +213,14 @@
   if (iscell (objf))
     switch (numel (objf))
       case 1
-        obj_fun = objf{1};
-        obj_grd = @(x, obj) fd_obj_grd (x, obj, obj_fun);
+        obj_fcn = objf{1};
+        obj_grd = @(x, obj) fd_obj_grd (x, obj, obj_fcn);
       case 2
-        obj_fun = objf{1};
+        obj_fcn = objf{1};
         obj_grd = objf{2};
         have_grd = 1;
       case 3
-        obj_fun = objf{1};
+        obj_fcn = objf{1};
         obj_grd = objf{2};
         obj_hess = objf{3};
         have_grd = 1;
@@ -229,31 +229,31 @@
         error ("sqp: invalid objective function specification");
     endswitch
   else
-    obj_fun = objf;   # No cell array, only obj_fun set
-    obj_grd = @(x, obj) fd_obj_grd (x, obj, obj_fun);
+    obj_fcn = objf;   # No cell array, only obj_fcn set
+    obj_grd = @(x, obj) fd_obj_grd (x, obj, obj_fcn);
   endif
 
-  ce_fun = @empty_cf;
+  ce_fcn = @empty_cf;
   ce_grd = @empty_jac;
   if (nargin > 2)
     if (iscell (cef))
       switch (numel (cef))
         case 1
-          ce_fun = cef{1};
-          ce_grd = @(x) fd_ce_jac (x, ce_fun);
+          ce_fcn = cef{1};
+          ce_grd = @(x) fd_ce_jac (x, ce_fcn);
         case 2
-          ce_fun = cef{1};
+          ce_fcn = cef{1};
           ce_grd = cef{2};
         otherwise
           error ("sqp: invalid equality constraint function specification");
       endswitch
     elseif (! isempty (cef))
-      ce_fun = cef;   # No cell array, only constraint equality function set
-      ce_grd = @(x) fd_ce_jac (x, ce_fun);
+      ce_fcn = cef;   # No cell array, only constraint equality function set
+      ce_grd = @(x) fd_ce_jac (x, ce_fcn);
     endif
   endif
 
-  ci_fun = @empty_cf;
+  ci_fcn = @empty_cf;
   ci_grd = @empty_jac;
   if (nargin > 3)
     ## constraint function given by user with possible gradient
@@ -274,15 +274,15 @@
       if (iscell (cif))
         switch (length (cif))
           case 1
-            ci_fun = cif{1};
+            ci_fcn = cif{1};
           case 2
-            ci_fun = cif{1};
+            ci_fcn = cif{1};
             ci_grd = cif{2};
           otherwise
            error ("sqp: invalid inequality constraint function specification");
         endswitch
       elseif (! isempty (cif))
-        ci_fun = cif;   # No cell array, only constraint inequality function set
+        ci_fcn = cif;   # No cell array, only constraint inequality function set
       endif
     else
       ## constraint inequality function with bounds present
@@ -322,7 +322,7 @@
         error ("sqp: upper bound smaller than lower bound");
       endif
       bounds_grad = [lb_grad; ub_grad];
-      ci_fun = @(x) cf_ub_lb (x, lb_idx, ub_idx, globals);
+      ci_fcn = @(x) cf_ub_lb (x, lb_idx, ub_idx, globals);
       ci_grd = @(x) cigrad_ub_lb (x, bounds_grad, globals);
     endif
 
@@ -350,15 +350,15 @@
   ## Seed x with initial guess and evaluate objective function, constraints,
   ## and gradients at initial value x0.
   ##
-  ## obj_fun   -- objective function
+  ## obj_fcn   -- objective function
   ## obj_grad  -- objective gradient
-  ## ce_fun    -- equality constraint functions
-  ## ci_fun    -- inequality constraint functions
+  ## ce_fcn    -- equality constraint functions
+  ## ci_fcn    -- inequality constraint functions
   ## A == [grad_{x_1} cx_fun, grad_{x_2} cx_fun, ..., grad_{x_n} cx_fun]^T
   x = x0;
 
-  obj = feval (obj_fun, x0);
-  globals.nfun = 1;
+  obj = feval (obj_fcn, x0);
+  globals.nfev = 1;
 
   if (have_grd)
     c = feval (obj_grd, x0);
@@ -374,10 +374,10 @@
     B = eye (n, n);
   endif
 
-  ce = feval (ce_fun, x0);
+  ce = feval (ce_fcn, x0);
   F = feval (ce_grd, x0);
 
-  ci = feval (ci_fun, x0);
+  ci = feval (ci_fcn, x0);
   C = feval (ci_grd, x0);
 
   A = [F; C];
@@ -391,7 +391,7 @@
   info = 0;
   iter = 0;
   ## report ();  # Called with no arguments to initialize reporting
-  ## report (iter, qp_iter, alpha, __sqp_nfun__, obj);
+  ## report (iter, qp_iter, alpha, __sqp_nfev__, obj);
 
   while (++iter < iter_max)
 
@@ -445,7 +445,7 @@
     ## Choose mu such that p is a descent direction for the chosen
     ## merit function phi.
     [x_new, alpha, obj_new, globals] = ...
-        linesearch_L1 (x, p, obj_fun, obj_grd, ce_fun, ci_fun, lambda, ...
+        linesearch_L1 (x, p, obj_fcn, obj_grd, ce_fcn, ci_fcn, lambda, ...
                        obj, c, globals);
 
     delx = x_new - x;
@@ -463,10 +463,10 @@
       c_new = feval (obj_grd, x_new, obj_new);
     endif
 
-    ce_new = feval (ce_fun, x_new);
+    ce_new = feval (ce_fcn, x_new);
     F_new = feval (ce_grd, x_new);
 
-    ci_new = feval (ci_fun, x_new);
+    ci_new = feval (ci_fcn, x_new);
     C_new = feval (ci_grd, x_new);
 
     A_new = [F_new; C_new];
@@ -533,7 +533,7 @@
 
     A = A_new;
 
-    ## report (iter, qp_iter, alpha, __sqp_nfun__, obj);
+    ## report (iter, qp_iter, alpha, __sqp_nfev__, obj);
 
   endwhile
 
@@ -542,24 +542,24 @@
     info = 103;
   endif
 
-  nf = globals.nfun;
+  nf = globals.nfev;
 
 endfunction
 
 
-function [merit, obj, globals] = phi_L1 (obj, obj_fun, ce_fun, ci_fun, ...
+function [merit, obj, globals] = phi_L1 (obj, obj_fcn, ce_fcn, ci_fcn, ...
                                          x, mu, globals)
 
-  ce = feval (ce_fun, x);
-  ci = feval (ci_fun, x);
+  ce = feval (ce_fcn, x);
+  ci = feval (ci_fcn, x);
 
   idx = ci < 0;
 
   con = [ce; ci(idx)];
 
   if (isempty (obj))
-    obj = feval (obj_fun, x);
-    globals.nfun += 1;
+    obj = feval (obj_fcn, x);
+    globals.nfev += 1;
   endif
 
   merit = obj;
@@ -573,7 +573,7 @@
 
 
 function [x_new, alpha, obj, globals] = ...
-  linesearch_L1 (x, p, obj_fun, obj_grd, ce_fun, ci_fun, lambda, obj, c, globals)
+  linesearch_L1 (x, p, obj_fcn, obj_grd, ce_fcn, ci_fcn, lambda, obj, c, globals)
 
   ## Choose parameters
   ##
@@ -593,13 +593,13 @@
 
   alpha = 1;
 
-  ce = feval (ce_fun, x);
+  ce = feval (ce_fcn, x);
 
-  [phi_x_mu, obj, globals] = phi_L1 (obj, obj_fun, ce_fun, ci_fun, x, ...
+  [phi_x_mu, obj, globals] = phi_L1 (obj, obj_fcn, ce_fcn, ci_fcn, x, ...
                                      mu, globals);
 
   D_phi_x_mu = c' * p;
-  d = feval (ci_fun, x);
+  d = feval (ci_fcn, x);
   ## only those elements of d corresponding
   ## to violated constraints should be included.
   idx = d < 0;
@@ -609,7 +609,7 @@
   endif
 
   while (1)
-    [p1, obj, globals] = phi_L1 ([], obj_fun, ce_fun, ci_fun, ...
+    [p1, obj, globals] = phi_L1 ([], obj_fcn, ce_fcn, ci_fcn, ...
                                  x+alpha*p, mu, globals);
     p2 = phi_x_mu+eta*alpha*D_phi_x_mu;
     if (p1 > p2)
@@ -668,9 +668,9 @@
 endfunction
 
 
-function grd = fd_obj_grd (x, obj, obj_fun)
+function grd = fd_obj_grd (x, obj, obj_fcn)
 
-  grd = fdgrd (obj_fun, x, obj);
+  grd = fdgrd (obj_fcn, x, obj);
 
 endfunction
 
@@ -689,9 +689,9 @@
 endfunction
 
 
-function jac = fd_ce_jac (x, ce_fun)
+function jac = fd_ce_jac (x, ce_fcn)
 
-  jac = fdjac (ce_fun, x);
+  jac = fdjac (ce_fcn, x);
 
 endfunction
 
@@ -734,12 +734,12 @@
 endfunction
 
 ## Utility function used to debug sqp
-function report (iter, qp_iter, alpha, nfun, obj)
+function report (iter, qp_iter, alpha, nfev, obj)
 
   if (nargin == 0)
-    printf ("  Itn ItQP     Step  Nfun     Objective\n");
+    printf ("  Itn ItQP     Step  Nfev     Objective\n");
   else
-    printf ("%5d %4d %8.1g %5d %13.6e\n", iter, qp_iter, alpha, nfun, obj);
+    printf ("%5d %4d %8.1g %5d %13.6e\n", iter, qp_iter, alpha, nfev, obj);
   endif
 
 endfunction
--- a/scripts/plot/draw/fplot.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/plot/draw/fplot.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,17 +24,17 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn  {} {} fplot (@var{fn})
-## @deftypefnx {} {} fplot (@var{fn}, @var{limits})
+## @deftypefn  {} {} fplot (@var{fcn})
+## @deftypefnx {} {} fplot (@var{fcn}, @var{limits})
 ## @deftypefnx {} {} fplot (@dots{}, @var{tol})
 ## @deftypefnx {} {} fplot (@dots{}, @var{n})
 ## @deftypefnx {} {} fplot (@dots{}, @var{fmt})
 ## @deftypefnx {} {} fplot (@dots{}, @var{property}, @var{value}, @dots{})
 ## @deftypefnx {} {} fplot (@var{hax}, @dots{})
 ## @deftypefnx {} {[@var{x}, @var{y}] =} fplot (@dots{})
-## Plot a function @var{fn} within the range defined by @var{limits}.
+## Plot a function @var{fcn} within the range defined by @var{limits}.
 ##
-## @var{fn} is a function handle, inline function, or string containing the
+## @var{fcn} is a function handle, inline function, or string containing the
 ## name of the function to evaluate.
 ##
 ## The limits of the plot are of the form @w{@code{[@var{xlo}, @var{xhi}]}} or
@@ -99,19 +99,19 @@
     print_usage ();
   endif
 
-  fn = varargin{1};
-  if (isa (fn, "inline"))
-    fn = vectorize (inline (fn));
-    nam = formula (fn);
-  elseif (is_function_handle (fn))
-    nam = func2str (fn);
-  elseif (all (isalnum (fn)))
-    nam = fn;
-  elseif (ischar (fn))
-    fn = vectorize (inline (fn));
-    nam = formula (fn);
+  fcn = varargin{1};
+  if (isa (fcn, "inline"))
+    fcn = vectorize (inline (fcn));
+    nam = formula (fcn);
+  elseif (is_function_handle (fcn))
+    nam = func2str (fcn);
+  elseif (all (isalnum (fcn)))
+    nam = fcn;
+  elseif (ischar (fcn))
+    fcn = vectorize (inline (fcn));
+    nam = formula (fcn);
   else
-    error ("fplot: FN must be a function handle, inline function, or string");
+    error ("fplot: FCN must be a function handle, inline function, or string");
   endif
 
   if (nargin > 1 && isnumeric (varargin{2}))
@@ -163,21 +163,21 @@
   endif
 
   try
-    y0 = feval (fn, x0);
+    y0 = feval (fcn, x0);
     if (isscalar (y0))
-      warning ("fplot: FN is not a vectorized function which reduces performance");
-      fn = @(x) arrayfun (fn, x);  # Create a new fn that accepts vectors
-      y0 = feval (fn, x0);
+      warning ("fplot: FCN is not a vectorized function which reduces performance");
+      fcn = @(x) arrayfun (fcn, x);  # Create a new fcn that accepts vectors
+      y0 = feval (fcn, x0);
     endif
   catch
     ## feval failed, maybe it is because the function is not vectorized?
-    fn = @(x) arrayfun (fn, x);  # Create a new fn that accepts vectors
-    y0 = feval (fn, x0);
-    warning ("fplot: FN is not a vectorized function which reduces performance");
+    fcn = @(x) arrayfun (fcn, x);  # Create a new fcn that accepts vectors
+    y0 = feval (fcn, x0);
+    warning ("fplot: FCN is not a vectorized function which reduces performance");
   end_try_catch
 
   x = linspace (limits(1), limits(2), n)';
-  y = feval (fn, x);
+  y = feval (fcn, x);
 
   if (rows (x0) == rows (y0))
     fcn_transpose = false;
@@ -186,7 +186,7 @@
     y0 = y0.';
     y = y.';
   else
-    error ("fplot: invalid function FN (# of outputs not equal to inputs)");
+    error ("fplot: invalid function FCN (# of outputs not equal to inputs)");
   endif
 
   err0 = Inf;
@@ -213,7 +213,7 @@
     err0 = err;
     n = 2 * (n - 1) + 1;
     x = linspace (limits(1), limits(2), n)';
-    y = feval (fn, x);
+    y = feval (fcn, x);
     if (fcn_transpose)
       y = y.';
     endif
@@ -274,16 +274,16 @@
 
 %!test
 %! ## Function requiring transpose
-%! fn = @(x) 2 * x(:).';
-%! [x, y] = fplot (fn, [-1, 1]);
+%! fcn = @(x) 2 * x(:).';
+%! [x, y] = fplot (fcn, [-1, 1]);
 %! assert (columns (y) == 1);
 %! assert (rows (x) == rows (y));
 %! assert (y, 2*x);
 
 %!test
 %! ## Constant value function
-%! fn = @(x) 0;
-%! [x, y] = fplot (fn, [-1, 1]);
+%! fcn = @(x) 0;
+%! [x, y] = fplot (fcn, [-1, 1]);
 %! assert (columns (y) == 1);
 %! assert (rows (x) == rows (y));
 %! assert (y, repmat ([0], size (x)));
@@ -313,15 +313,15 @@
 ## Test input validation
 %!error <Invalid call> fplot ()
 %!error <Invalid call> fplot (1,2,3,4,5,6)
-%!error <FN must be a function handle> fplot (1, [0 1])
+%!error <FCN must be a function handle> fplot (1, [0 1])
 %!error <LIMITS must be a real vector> fplot (@cos, [i, 2*i])
 %!error <LIMITS must be a real vector with 2 or 4> fplot (@cos, [1])
 %!error <LIMITS must be a real vector with 2 or 4> fplot (@cos, [1 2 3])
 %!error <bad input in position 2> fplot (@cos, "linewidth")
 %!error <bad input in position 3> fplot (@cos, [-1,1], {1})
-%!warning <FN is not a vectorized function>
-%! fn = @(x) 0;
-%! [x,y] = fplot (fn, [-1,1]);
-%!error <invalid function FN>
-%! fn = @(x) [x;x];
-%! fplot (fn, [-1,1]);
+%!warning <FCN is not a vectorized function>
+%! fcn = @(x) 0;
+%! [x,y] = fplot (fcn, [-1,1]);
+%!error <invalid function FCN>
+%! fcn = @(x) [x;x];
+%! fplot (fcn, [-1,1]);
--- a/scripts/plot/draw/plotyy.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/plot/draw/plotyy.m	Mon Apr 04 18:14:56 2022 -0700
@@ -25,7 +25,7 @@
 
 ## -*- texinfo -*-
 ## @deftypefn  {} {} plotyy (@var{x1}, @var{y1}, @var{x2}, @var{y2})
-## @deftypefnx {} {} plotyy (@dots{}, @var{fun})
+## @deftypefnx {} {} plotyy (@dots{}, @var{fcn})
 ## @deftypefnx {} {} plotyy (@dots{}, @var{fun1}, @var{fun2})
 ## @deftypefnx {} {} plotyy (@var{hax}, @dots{})
 ## @deftypefnx {} {[@var{ax}, @var{h1}, @var{h2}] =} plotyy (@dots{})
@@ -36,8 +36,8 @@
 ##
 ## By default the arguments are evaluated with
 ## @code{feval (@@plot, @var{x}, @var{y})}.  However the type of plot can be
-## modified with the @var{fun} argument, in which case the plots are
-## generated by @code{feval (@var{fun}, @var{x}, @var{y})}.  @var{fun} can be
+## modified with the @var{fcn} argument, in which case the plots are
+## generated by @code{feval (@var{fcn}, @var{x}, @var{y})}.  @var{fcn} can be
 ## a function handle, an inline function, or a string of a function name.
 ##
 ## The function to use for each of the plots can be independently defined
--- a/scripts/plot/draw/private/__bar__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/plot/draw/private/__bar__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,16 +24,16 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {} __bar__ (@var{vertical}, @var{func}, @dots{})
+## @deftypefn {} {} __bar__ (@var{vertical}, @var{fcn}, @dots{})
 ## Undocumented internal function.
 ## @end deftypefn
 
-function varargout = __bar__ (func, vertical, varargin)
+function varargout = __bar__ (fcn, vertical, varargin)
 
-  [hax, varargin, nargin] = __plt_get_axis_arg__ (func, varargin{:});
+  [hax, varargin, nargin] = __plt_get_axis_arg__ (fcn, varargin{:});
 
   if (! isnumeric (varargin{1}))
-    error ("%s: Y must be numeric", func);
+    error ("%s: Y must be numeric", fcn);
   endif
 
   width = 0.8;
@@ -63,9 +63,9 @@
       idx = 2;
     else
       if (! isvector (x))
-        error ("%s: X must be a vector", func);
+        error ("%s: X must be a vector", fcn);
       elseif (numel (unique (x)) != numel (x))
-        error ("%s: X vector values must be unique", func);
+        error ("%s: X vector values must be unique", fcn);
       endif
       idx = 3;
     endif
@@ -99,7 +99,7 @@
     else
       if ((ischar (varargin{idx}) || iscellstr (varargin{idx}))
           && ! have_line_spec)
-        [linespec, valid] = __pltopt__ (func, varargin{idx}, false);
+        [linespec, valid] = __pltopt__ (fcn, varargin{idx}, false);
         if (valid)
           have_line_spec = true;
           ## FIXME: strange parse error requires semicolon to be spaced
@@ -133,7 +133,7 @@
     y = y.';
   endif
   if (ngrp != rows (y))
-    error ("%s: length of X and Y must be equal", func);
+    error ("%s: length of X and Y must be equal", fcn);
   endif
 
   nbars = columns (y);
--- a/scripts/plot/draw/private/__ezplot__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/plot/draw/private/__ezplot__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,7 +24,7 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {[@var{h}, @var{needusage}] =} __ezplot__ (@var{pltfunc}, @var{varargin})
+## @deftypefn {} {[@var{h}, @var{needusage}] =} __ezplot__ (@var{pltfcn}, @var{varargin})
 ## Undocumented internal function.
 ## @end deftypefn
 
@@ -34,11 +34,11 @@
 ##           called was called correctly.  The actual plotting occurs near
 ##           the end in an unwind_protect block.
 
-function [h, needusage] = __ezplot__ (pltfunc, varargin)
+function [h, needusage] = __ezplot__ (pltfcn, varargin)
 
-  ezfunc = ["ez" pltfunc];
+  ezfcn = ["ez" pltfcn];
 
-  [hax, varargin, nargin] = __plt_get_axis_arg__ (ezfunc, varargin{:});
+  [hax, varargin, nargin] = __plt_get_axis_arg__ (ezfcn, varargin{:});
 
   ## Define outputs early in case of shorting out of function with return;
   h = [];
@@ -48,14 +48,14 @@
     return;
   endif
 
-  iscontour = strncmp (pltfunc, "contour", 7);
+  iscontour = strncmp (pltfcn, "contour", 7);
 
   ## Defaults for ezplot
   isplot  = true;
   isplot3 = false;
   ispolar = false;
   nargs = 1;
-  switch (pltfunc)
+  switch (pltfcn)
     case "plot"
       ## defaults already set
 
@@ -75,24 +75,24 @@
   endswitch
 
   parametric = false;
-  fun = varargin{1};
-  if (ischar (fun))
-    if (exist (fun, "file") || exist (fun, "builtin"))
-      fun = str2func (fun);            # convert to function handle
+  fcn = varargin{1};
+  if (ischar (fcn))
+    if (exist (fcn, "file") || exist (fcn, "builtin"))
+      fcn = str2func (fcn);            # convert to function handle
     else
-      fun = vectorize (inline (fun));  # convert to inline function
+      fcn = vectorize (inline (fcn));  # convert to inline function
     endif
   endif
 
-  if (isa (fun, "inline"))
-    argids = argnames (fun);
+  if (isa (fcn, "inline"))
+    argids = argnames (fcn);
     if (isplot && length (argids) == 2)
       nargs = 2;
     elseif (numel (argids) != nargs)
-      error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+      error ("%s: expecting a function of %d arguments", ezfcn, nargs);
     endif
-    fun = vectorize (fun);
-    fstr = formula (fun);
+    fcn = vectorize (fcn);
+    fstr = formula (fcn);
     if (isplot)
       xarg = argids{1};
       if (nargs == 2)
@@ -110,8 +110,8 @@
       xarg = argids{1};
       yarg = argids{2};
     endif
-  elseif (is_function_handle (fun))
-    fstr = func2str (fun);
+  elseif (is_function_handle (fcn))
+    fstr = func2str (fcn);
     idx = index (fstr, ')');
     if (idx != 0)
       args = regexp (fstr(3:(idx-1)), '\w+', 'match');
@@ -119,7 +119,7 @@
     else
       args = {"x"};
       try
-        if (builtin ("nargin", fun) == 2)
+        if (builtin ("nargin", fcn) == 2)
           args{2} = "y";
         endif
       end_try_catch
@@ -127,7 +127,7 @@
     if (isplot && length (args) == 2)
       nargs = 2;
     elseif (numel (args) != nargs)
-      error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+      error ("%s: expecting a function of %d arguments", ezfcn, nargs);
     endif
     if (isplot)
       xarg = args{1};
@@ -147,11 +147,11 @@
       yarg = args{2};
     endif
   else
-    error ("%s: F must be a string or function handle", ezfunc);
+    error ("%s: F must be a string or function handle", ezfcn);
   endif
 
   if (nargin > 2 || (nargin == 2 && isplot))
-    funx = fun;
+    funx = fcn;
     fstrx = fstr;
     funy = varargin{2};
     if (ischar (funy) && ! strcmp (funy, "circ") && ! strcmp (funy, "animate"))
@@ -162,13 +162,13 @@
         funy = vectorize (inline (funy));
       endif
       if (numel (argnames (funy)) != nargs)
-        error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+        error ("%s: expecting a function of %d arguments", ezfcn, nargs);
       endif
       fstry = formula (funy);
     elseif (isa (funy, "inline"))
       parametric = true;
       if (numel (argnames (funy)) != nargs)
-        error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+        error ("%s: expecting a function of %d arguments", ezfcn, nargs);
       endif
       funy = vectorize (funy);
       fstry = formula (funy);
@@ -183,7 +183,7 @@
         args = {"y"};
       endif
       if (numel (args) != nargs)
-        error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+        error ("%s: expecting a function of %d arguments", ezfcn, nargs);
       endif
     endif
 
@@ -192,7 +192,7 @@
       return;
     elseif (parametric && isplot)
       if (nargs == 2)
-        error ("%s: can not define a parametric function in this manner", ezfunc);
+        error ("%s: can not define a parametric function in this manner", ezfcn);
       else
         xarg = "x";
         yarg = "y";
@@ -207,12 +207,12 @@
           funz = vectorize (inline (funz));
         endif
         if (numel (argnames (funz)) > nargs)
-          error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+          error ("%s: expecting a function of %d arguments", ezfcn, nargs);
         endif
         fstrz = formula (funz);
       elseif (isa (funz, "inline"))
         if (numel (argnames (funz)) != nargs)
-          error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+          error ("%s: expecting a function of %d arguments", ezfcn, nargs);
         endif
         funz = vectorize (funz);
         fstrz = formula (funz);
@@ -226,10 +226,10 @@
           args = {"z"};
         endif
         if (numel (args) != nargs)
-          error ("%s: expecting a function of %d arguments", ezfunc, nargs);
+          error ("%s: expecting a function of %d arguments", ezfcn, nargs);
         endif
       else
-        error ("%s: parametric plots require 3 functions", ezfunc);
+        error ("%s: parametric plots require 3 functions", ezfcn);
       endif
     endif
   endif
@@ -264,7 +264,7 @@
     elseif (numel (arg) == 4 && isempty (domain))
       domain = arg(:).';
     else
-      error ("%s: expecting scalar N, or 2-/4-element vector DOM", ezfunc);
+      error ("%s: expecting scalar N, or 2-/4-element vector DOM", ezfcn);
     endif
   endwhile
 
@@ -273,11 +273,11 @@
     return;
   elseif (circ && parametric)
     error ("%s: can not have both circular domain and parametric function",
-           ezfunc);
+           ezfcn);
   endif
 
   if (animate && ! isplot3)
-    error ("%s: animate option only valid for ezplot3", ezfunc);
+    error ("%s: animate option only valid for ezplot3", ezfcn);
   endif
 
   if (parametric)
@@ -366,7 +366,7 @@
       endif
     else  # non-parametric plots
       if (isplot && nargs == 2)
-        Z = feval (fun, X, Y);
+        Z = feval (fcn, X, Y);
 
         ## Matlab returns line objects for this case and so can't call
         ## contour directly as it returns patch objects to allow colormaps
@@ -385,18 +385,18 @@
         endwhile
       else
         if (ispolar)
-          Z = feval (fun, X);
+          Z = feval (fcn, X);
           ## FIXME: Why aren't singularities eliminated for polar plots?
         elseif (isplot)
-          Z = feval (fun, X);
+          Z = feval (fcn, X);
           ## Eliminate the singularities
           Z = __eliminate_sing__ (Z);
           domain = find_valid_domain (X, [], Z);
         elseif (iscontour)
-          Z = feval (fun, X, Y);
+          Z = feval (fcn, X, Y);
           Z = __eliminate_sing__ (Z);
         else  #  mesh, surf plots
-          Z = feval (fun, X, Y);
+          Z = feval (fcn, X, Y);
           Z = __eliminate_sing__ (Z);
           if (circ)
             ## Use domain calculated at the start.
@@ -426,7 +426,7 @@
   unwind_protect
     hax = newplot (hax);
     if (iscontour)
-      [~, h] = feval (pltfunc, hax, X, Y, Z);
+      [~, h] = feval (pltfcn, hax, X, Y, Z);
     elseif (isplot && nargs == 2)
       h = zeros (length (XX), 1);
       hold_state = get (hax, "nextplot");
@@ -442,7 +442,7 @@
       set (hax, "nextplot", hold_state);
       axis (hax, domain);
     elseif (isplot || ispolar)
-      h = feval (pltfunc, hax, X, Z);
+      h = feval (pltfcn, hax, X, Z);
       if (isplot)
         if (! parametric)
           axis (hax, domain);
@@ -454,12 +454,12 @@
       if (animate)
         comet3 (hax, X, Y, Z);
       else
-        h = feval (pltfunc, hax, X, Y, Z);
+        h = feval (pltfcn, hax, X, Y, Z);
       endif
       grid (hax, "on");
       zlabel (hax, "z");
     else  # mesh and surf plots
-      h = feval (pltfunc, hax, X, Y, Z);
+      h = feval (pltfcn, hax, X, Y, Z);
       ## FIXME: surf, mesh should really do a better job of setting zlim
       if (! parametric)
         axis (hax, domain);
--- a/scripts/profiler/html/flat_entry.html	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/profiler/html/flat_entry.html	Mon Apr 04 18:14:56 2022 -0700
@@ -1,5 +1,5 @@
 <tr>
-  <td><a href="func-%num.html">%name</a></td>
+  <td><a href="function-%num.html">%name</a></td>
   <td>%timeabs</td>
   <td>%timerel</td>
   <td>%calls</td>
--- a/scripts/profiler/profexport.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/profiler/profexport.m	Mon Apr 04 18:14:56 2022 -0700
@@ -87,7 +87,7 @@
 
   __writeFlat (fullfile (dir, "index.html"), name, data.FunctionTable);
   for i = 1 : length (data.FunctionTable)
-    __writeFunc (fullfile (dir, sprintf ("func-%d.html", i)), name, ...
+    __writeFunc (fullfile (dir, sprintf ("function-%d.html", i)), name, ...
                  data.FunctionTable, i);
   endfor
 
@@ -175,7 +175,7 @@
     return;
   endif
 
-  template = "<a href='func-%num.html'>%name</a>";
+  template = "<a href='function-%num.html'>%name</a>";
 
   lst = "";
   for i = 1 : length (inds)
--- a/scripts/signal/movfun.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/signal/movfun.m	Mon Apr 04 18:14:56 2022 -0700
@@ -227,29 +227,29 @@
 
   ## Obtain function for boundary conditions
   if (isnumeric (bc))
-    bcfunc = @replaceval_bc;
-    bcfunc (true, bc);  # initialize replaceval function with value
+    bcfcn = @replaceval_bc;
+    bcfcn (true, bc);  # initialize replaceval function with value
   else
     switch (tolower (bc))
       case "shrink"
-        bcfunc = @shrink_bc;
+        bcfcn = @shrink_bc;
 
       case "discard"
-        bcfunc = [];
+        bcfcn = [];
         C -= length (Cpre);
         Cpre = Cpos = [];
         N = length (C);
         szx(dperm(1)) = N;
 
       case "fill"
-        bcfunc = @replaceval_bc;
-        bcfunc (true, NaN);
+        bcfcn = @replaceval_bc;
+        bcfcn (true, NaN);
 
       case "same"
-        bcfunc = @same_bc;
+        bcfcn = @same_bc;
 
       case "periodic"
-        bcfunc = @periodic_bc;
+        bcfcn = @periodic_bc;
 
     endswitch
   endif
@@ -279,7 +279,7 @@
   ## FIXME: Is it faster with cellfun?  Don't think so, but needs testing.
   y = zeros (N, ncols, soutdim);
   parfor i = 1:ncols
-    y(:,i,:) = movfun_oncol (fcn_, x(:,i), wlen, bcfunc,
+    y(:,i,:) = movfun_oncol (fcn_, x(:,i), wlen, bcfcn,
                              slc, C, Cpre, Cpos, win, soutdim);
   endparfor
 
@@ -290,7 +290,7 @@
 
 endfunction
 
-function y = movfun_oncol (fcn, x, wlen, bcfunc, slcidx, C, Cpre, Cpos, win, odim)
+function y = movfun_oncol (fcn, x, wlen, bcfcn, slcidx, C, Cpre, Cpos, win, odim)
 
   N = length (Cpre) + length (C) + length (Cpos);
   y = NA (N, odim);
@@ -300,10 +300,10 @@
 
   ## Process boundaries
   if (! isempty (Cpre))
-    y(Cpre,:) = bcfunc (fcn, x, Cpre, win, wlen, odim);
+    y(Cpre,:) = bcfcn (fcn, x, Cpre, win, wlen, odim);
   endif
   if (! isempty (Cpos))
-    y(Cpos,:) = bcfunc (fcn, x, Cpos, win, wlen, odim);
+    y(Cpos,:) = bcfcn (fcn, x, Cpos, win, wlen, odim);
   endif
 
 endfunction
--- a/scripts/sparse/bicg.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/bicg.m	Mon Apr 04 18:14:56 2022 -0700
@@ -43,9 +43,9 @@
 ##
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline function
-## @code{Afun} such that @w{@code{Afun (x, "notransp") = A * x}} and
-## @w{@code{Afun (x, "transp") = A' * x}}.  Additional parameters to
-## @code{Afun} may be passed after @var{x0}.
+## @code{Afcn} such that @w{@code{Afcn (x, "notransp") = A * x}} and
+## @w{@code{Afcn (x, "transp") = A' * x}}.  Additional parameters to
+## @code{Afcn} may be passed after @var{x0}.
 ##
 ## @item @var{b} is the right-hand side vector.  It must be a column vector
 ## with the same number of rows as @var{A}.
@@ -84,7 +84,7 @@
 ## @end itemize
 ##
 ## Any arguments which follow @var{x0} are treated as parameters, and passed in
-## an appropriate manner to any of the functions (@var{Afun} or @var{Mfun}) or
+## an appropriate manner to any of the functions (@var{Afcn} or @var{Mfcn}) or
 ## that have been given to @code{bicg}.
 ##
 ## The output parameters are:
@@ -140,13 +140,13 @@
 ## restart = 5;
 ## [M1, M2] = ilu (A);  # in this tridiag case, it corresponds to lu (A)
 ## M = M1 * M2;
-## Afun = @@(x, string) strcmp (string, "notransp") * (A * x) + ...
+## Afcn = @@(x, string) strcmp (string, "notransp") * (A * x) + ...
 ##                      strcmp (string, "transp") * (A' * x);
-## Mfun = @@(x, string) strcmp (string, "notransp") * (M \ x) + ...
+## Mfcn = @@(x, string) strcmp (string, "notransp") * (M \ x) + ...
 ##                      strcmp (string, "transp") * (M' \ x);
-## M1fun = @@(x, string) strcmp (string, "notransp") * (M1 \ x) + ...
+## M1fcn = @@(x, string) strcmp (string, "notransp") * (M1 \ x) + ...
 ##                      strcmp (string, "transp") * (M1' \ x);
-## M2fun = @@(x, string) strcmp (string, "notransp") * (M2 \ x) + ...
+## M2fcn = @@(x, string) strcmp (string, "notransp") * (M2 \ x) + ...
 ##                      strcmp (string, "transp") * (M2' \ x);
 ## @end group
 ## @end example
@@ -161,7 +161,7 @@
 ## @code{@var{A}*@var{x}} and @code{@var{A'}*@var{x}}
 ##
 ## @example
-## x = bicg (Afun, b, [], n)
+## x = bicg (Afcn, b, [], n)
 ## @end example
 ##
 ## @sc{Example 3:} @code{bicg} with a preconditioner matrix @var{M}
@@ -173,7 +173,7 @@
 ## @sc{Example 4:} @code{bicg} with a function as preconditioner
 ##
 ## @example
-## x = bicg (Afun, b, 1e-6, n, Mfun)
+## x = bicg (Afcn, b, 1e-6, n, Mfcn)
 ## @end example
 ##
 ## @sc{Example 5:} @code{bicg} with preconditioner matrices @var{M1}
@@ -186,7 +186,7 @@
 ## @sc{Example 6:} @code{bicg} with functions as preconditioners
 ##
 ## @example
-## x = bicg (Afun, b, 1e-6, n, M1fun, M2fun)
+## x = bicg (Afcn, b, 1e-6, n, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 7:} @code{bicg} with as input a function requiring an argument
@@ -207,8 +207,8 @@
 ##   endif
 ## endfunction
 ##
-## Apfun = @@(x, string, p) Ap (A, x, string, p);
-## x = bicg (Apfun, b, [], [], [], [], [], 2);
+## Apfcn = @@(x, string, p) Ap (A, x, string, p);
+## x = bicg (Apfcn, b, [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -223,7 +223,7 @@
 function [x_min, flag, relres, iter_min, resvec] = ...
          bicg (A, b, tol = [], maxit = [], M1 = [], M2 = [], x0 = [], varargin)
 
-  [Afun, M1fun, M2fun] =  __alltohandles__ (A, b, M1, M2, "bicg");
+  [Afcn, M1fcn, M2fcn] =  __alltohandles__ (A, b, M1, M2, "bicg");
 
   [tol, maxit, x0] = __default__input__ ({1e-06, min(rows(b), 20), ...
                                           zeros(rows (b),1)}, tol, maxit, x0);
@@ -253,17 +253,17 @@
   iter = iter_min = 0;
   flag = 1;  # Default flag is "maximum number of iterations reached"
   resvec = zeros (maxit + 1, 1);
-  r0 = b - Afun (x, "notransp", varargin{:});  # Residual of the system
-  s0 = c - Afun (x, "transp", varargin{:});    # Res. of the "dual system"
+  r0 = b - Afcn (x, "notransp", varargin{:});  # Residual of the system
+  s0 = c - Afcn (x, "transp", varargin{:});    # Res. of the "dual system"
   resvec(1) = norm (r0, 2);
 
   try
     warning ("error", "Octave:singular-matrix", "local");
-    prec_r0 = M1fun (r0, "notransp", varargin{:});  # r0 preconditioned
+    prec_r0 = M1fcn (r0, "notransp", varargin{:});  # r0 preconditioned
     prec_s0 = s0;
-    prec_r0 = M2fun (prec_r0, "notransp", varargin{:});
-    prec_s0 = M2fun (prec_s0, "transp", varargin{:});
-    prec_s0 = M1fun (prec_s0, "transp", varargin{:});  # s0 preconditioned
+    prec_r0 = M2fcn (prec_r0, "notransp", varargin{:});
+    prec_s0 = M2fcn (prec_s0, "transp", varargin{:});
+    prec_s0 = M1fcn (prec_s0, "transp", varargin{:});  # s0 preconditioned
     p = prec_r0;  # Direction of the system
     q = prec_s0;  # Direction of the "dual system"
   catch
@@ -271,7 +271,7 @@
   end_try_catch
 
   while ((flag != 2) && (iter < maxit) && (resvec(iter+1) >= norm_b * tol))
-    v = Afun (p, "notransp", varargin{:});
+    v = Afcn (p, "notransp", varargin{:});
     prod_qv = q' * v;
     alpha = (s0' * prec_r0);
     if (abs (prod_qv) <= eps * abs (alpha))
@@ -282,18 +282,18 @@
     x += alpha * p;
     prod_rs = (s0' * prec_r0);  # Product between r0 and s0
     r0 -= alpha * v;
-    s0 -= conj (alpha) * Afun (q, "transp", varargin{:});
-    prec_r0 = M1fun (r0, "notransp", varargin{:});
+    s0 -= conj (alpha) * Afcn (q, "transp", varargin{:});
+    prec_r0 = M1fcn (r0, "notransp", varargin{:});
     prec_s0 = s0;
-    prec_r0 = M2fun (prec_r0, "notransp", varargin{:});
+    prec_r0 = M2fcn (prec_r0, "notransp", varargin{:});
     beta = s0' * prec_r0;
     if (abs (prod_rs) <= abs (beta))
       flag = 4;
       break;
     endif
     beta ./= prod_rs;
-    prec_s0 = M2fun (prec_s0, "transp", varargin{:});
-    prec_s0 = M1fun (prec_s0, "transp", varargin{:});
+    prec_s0 = M2fcn (prec_s0, "transp", varargin{:});
+    prec_s0 = M1fcn (prec_s0, "transp", varargin{:});
     iter += 1;
     resvec(iter+1) = norm (r0);
     if (resvec(iter+1) <= resvec(iter_min+1))
@@ -378,11 +378,11 @@
 %!   endif
 %! endfunction
 %!
-%! Afun = @(x, string) Ap (A, x, string, 1);
-%! x = bicg (Afun, b, [], n);
+%! Afcn = @(x, string) Ap (A, x, string, 1);
+%! x = bicg (Afcn, b, [], n);
 %! x = bicg (A, b, 1e-6, n, M);
 %! x = bicg (A, b, 1e-6, n, M1, M2);
-%! function y = Mfun (M, x, string)
+%! function y = Mfcn (M, x, string)
 %!   if (strcmp (string, "notransp"))
 %!     y = M \ x;
 %!   else
@@ -390,14 +390,14 @@
 %!   endif
 %! endfunction
 %!
-%! M1fun = @(x, string) Mfun (M, x, string);
-%! x = bicg (Afun, b, 1e-6, n, M1fun);
-%! M1fun = @(x, string) Mfun (M1, x, string);
-%! M2fun = @(x, string) Mfun (M2, x, string);
-%! x = bicg (Afun, b, 1e-6, n, M1fun, M2fun);
-%! Afun = @(x, string, p) Ap (A, x, string, p);
+%! M1fcn = @(x, string) Mfcn (M, x, string);
+%! x = bicg (Afcn, b, 1e-6, n, M1fcn);
+%! M1fcn = @(x, string) Mfcn (M1, x, string);
+%! M2fcn = @(x, string) Mfcn (M2, x, string);
+%! x = bicg (Afcn, b, 1e-6, n, M1fcn, M2fcn);
+%! Afcn = @(x, string, p) Ap (A, x, string, p);
 %! ## Solution of A^2 * x = b
-%! x = bicg (Afun, b, [], 2*n, [], [], [], 2);
+%! x = bicg (Afcn, b, [], 2*n, [], [], [], 2);
 
 %!test
 %! ## Check that all type of inputs work
@@ -405,31 +405,31 @@
 %! b = A * ones (5, 1);
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;
-%! Afun = @(z, string) strcmp (string, "notransp") * (A * z) + ...
+%! Afcn = @(z, string) strcmp (string, "notransp") * (A * z) + ...
 %!                     strcmp (string, "transp") * (A' * z);
-%! M1_fun = @(z, string) strcmp (string,"notransp") * (M1 \ z) + ...
+%! M1_fcn = @(z, string) strcmp (string,"notransp") * (M1 \ z) + ...
 %!                         strcmp (string, "transp") * (M1' \ z);
-%! M2_fun = @(z, string) strcmp (string, "notransp") * (M2 \ z) + ...
+%! M2_fcn = @(z, string) strcmp (string, "notransp") * (M2 \ z) + ...
 %!                         strcmp (string, "transp") * (M2' \ z);
 %! [x, flag] = bicg (A, b);
 %! assert (flag, 0);
 %! [x, flag] = bicg (A, b, [], [], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicg (A, b, [], [], M1_fun, M2_fun);
+%! [x, flag] = bicg (A, b, [], [], M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicg (A, b,[], [], M1_fun, M2);
+%! [x, flag] = bicg (A, b,[], [], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicg (A, b,[], [], M1, M2_fun);
+%! [x, flag] = bicg (A, b,[], [], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicg (Afun, b);
+%! [x, flag] = bicg (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = bicg (Afun, b,[], [], M1, M2);
+%! [x, flag] = bicg (Afcn, b,[], [], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicg (Afun, b,[], [], M1_fun, M2);
+%! [x, flag] = bicg (Afcn, b,[], [], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicg (Afun, b,[], [], M1, M2_fun);
+%! [x, flag] = bicg (Afcn, b,[], [], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicg (Afun, b,[], [], M1_fun, M2_fun);
+%! [x, flag] = bicg (Afcn, b,[], [], M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!test
@@ -443,7 +443,7 @@
 %! [x, flag, relres, iter, resvec] = bicg (A, b, tol, maxit, M1, M2);
 %! assert (norm (b - A*x) / norm (b), 0, tol);
 
-%!function y = afun (x, t, a)
+%!function y = afcn (x, t, a)
 %!  switch (t)
 %!    case "notransp"
 %!      y = a * x;
@@ -461,7 +461,7 @@
 %! M1 = spdiags ([ones(n,1)/(-2) ones(n,1)],-1:0, n, n);
 %! M2 = spdiags ([4*ones(n,1) -ones(n,1)], 0:1, n, n);
 %!
-%! [x, flag, relres, iter, resvec] = bicg (@(x, t) afun (x, t, A),
+%! [x, flag, relres, iter, resvec] = bicg (@(x, t) afcn (x, t, A),
 %!                                         b, tol, maxit, M1, M2);
 %! assert (x, ones (size (b)), 1e-7);
 
@@ -516,7 +516,7 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x, trans)
+%!function y = Afcn (x, trans)
 %!  A = sparse (toeplitz ([2, 1, 0, 0], [2, -1, 0, 0]));
 %!  if (strcmp (trans, "notransp"))
 %!     y = A * x;
@@ -525,7 +525,7 @@
 %!  endif
 %!endfunction
 %!
-%! [x, flag] = bicg ("Afun", [1; 2; 2; 3]);
+%! [x, flag] = bicg ("Afcn", [1; 2; 2; 3]);
 %! assert (x, ones (4, 1), 1e-6);
 
 %!test
--- a/scripts/sparse/bicgstab.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/bicgstab.m	Mon Apr 04 18:14:56 2022 -0700
@@ -36,8 +36,8 @@
 ##
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline
-## function @code{Afun} such that @code{Afun(x) = A * x}.  Additional
-## parameters to @code{Afun} are passed after @var{x0}.
+## function @code{Afcn} such that @code{Afcn(x) = A * x}.  Additional
+## parameters to @code{Afcn} are passed after @var{x0}.
 ##
 ## @item @var{b} is the right hand side vector.  It must be a column vector
 ## with the same number of rows as @var{A}.
@@ -118,10 +118,10 @@
 ## restart = 5;
 ## [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A)
 ## M = M1 * M2;
-## Afun = @@(x) A * x;
-## Mfun = @@(x) M \ x;
-## M1fun = @@(x) M1 \ x;
-## M2fun = @@(x) M2 \ x;
+## Afcn = @@(x) A * x;
+## Mfcn = @@(x) M \ x;
+## M1fcn = @@(x) M1 \ x;
+## M2fcn = @@(x) M2 \ x;
 ## @end group
 ## @end example
 ##
@@ -135,7 +135,7 @@
 ## @code{@var{A} * @var{x}}
 ##
 ## @example
-## x = bicgstab (Afun, b, [], n)
+## x = bicgstab (Afcn, b, [], n)
 ## @end example
 ##
 ## @sc{Example 3:} @code{bicgstab} with a preconditioner matrix @var{M}
@@ -147,7 +147,7 @@
 ## @sc{Example 4:} @code{bicgstab} with a function as preconditioner
 ##
 ## @example
-## x = bicgstab (Afun, b, 1e-6, n, Mfun)
+## x = bicgstab (Afcn, b, 1e-6, n, Mfcn)
 ## @end example
 ##
 ## @sc{Example 5:} @code{bicgstab} with preconditioner matrices @var{M1}
@@ -160,7 +160,7 @@
 ## @sc{Example 6:} @code{bicgstab} with functions as preconditioners
 ##
 ## @example
-## x = bicgstab (Afun, b, 1e-6, n, M1fun, M2fun)
+## x = bicgstab (Afcn, b, 1e-6, n, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 7:} @code{bicgstab} with as input a function requiring
@@ -174,8 +174,8 @@
 ##      y = A * y;
 ##    endfor
 ##  endfunction
-## Apfun = @@(x, string, p) Ap (A, x, string, p);
-## x = bicgstab (Apfun, b, [], [], [], [], [], 2);
+## Apfcn = @@(x, string, p) Ap (A, x, string, p);
+## x = bicgstab (Apfcn, b, [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -211,7 +211,7 @@
                    x0 = [], varargin)
 
   ## Check consistency and  type of A, M1, M2
-  [Afun, M1fun, M2fun] =  __alltohandles__ (A, b, M1, M2, "bicgstab");
+  [Afcn, M1fcn, M2fcn] =  __alltohandles__ (A, b, M1, M2, "bicgstab");
 
   ## Check if input tol are empty (set them to default if necessary)
   [tol, maxit, x0] = __default__input__ ({1e-06, min(rows(b), 20), ...
@@ -240,7 +240,7 @@
   ## default setting of flag is 1 (i.e. max number of iterations reached)
   flag = 1;
 
-  res = b - feval (Afun, x, varargin{:});
+  res = b - feval (Afcn, x, varargin{:});
   rr = p = res; # rr is r_star
   rho_1 = rr' * res;
   resvec (1) = norm (res,2);
@@ -249,14 +249,14 @@
   ## To check if the preconditioners are singular or they have some NaN
   try
     warning ("error", "Octave:singular-matrix", "local");
-    p_hat = feval (M1fun, p, varargin{:});
-    p_hat = feval (M2fun, p_hat, varargin{:});
+    p_hat = feval (M1fcn, p, varargin{:});
+    p_hat = feval (M2fcn, p_hat, varargin{:});
   catch
     flag = 2;
   end_try_catch
 
   while (flag !=2) && (iter < d_maxit) && (resvec (iter + 1) >= real_tol)
-    v = feval (Afun, p_hat, varargin{:});
+    v = feval (Afcn, p_hat, varargin{:});
     prod_tmp = (rr' * v);
     if (prod_tmp == 0)
       flag = 4;
@@ -275,9 +275,9 @@
       x_min = x;
       iter_min = iter;
     endif
-    s_hat = feval (M1fun, s, varargin{:});
-    s_hat = feval (M2fun, s_hat, varargin{:});
-    t = feval (Afun, s_hat, varargin{:});
+    s_hat = feval (M1fcn, s, varargin{:});
+    s_hat = feval (M2fcn, s_hat, varargin{:});
+    t = feval (Afcn, s_hat, varargin{:});
     omega = (t' * s) / (t' * t);
     if (omega == 0) # x and residual don't change and the next it will be NaN
       flag = 4;
@@ -304,8 +304,8 @@
     endif
     beta = (rho_1 / rho_2) * (alpha / omega);
     p = res + beta * (p - omega*  v);
-    p_hat = feval (M1fun, p, varargin{:});
-    p_hat = feval (M2fun, p_hat, varargin{:});
+    p_hat = feval (M1fcn, p, varargin{:});
+    p_hat = feval (M2fcn, p_hat, varargin{:});
   endwhile
   resvec = resvec (1:iter+1,1);
 
@@ -362,15 +362,15 @@
 %! [M1, M2] = ilu (A + 0.1 * eye (n));
 %! M = M1 * M2;
 %! x = bicgstab (A, b, [], n);
-%! Afun = @(x) A * x;
-%! x = bicgstab (Afun, b, [], n);
+%! Afcn = @(x) A * x;
+%! x = bicgstab (Afcn, b, [], n);
 %! x = bicgstab (A, b, 1e-6, n, M);
 %! x = bicgstab (A, b, 1e-6, n, M1, M2);
-%! Mfun = @(z) M \ z;
-%! x = bicgstab (Afun, b, 1e-6, n, Mfun);
-%! M1fun = @(z) M1 \ z;
-%! M2fun = @(z) M2 \ z;
-%! x = bicgstab (Afun, b, 1e-6, n, M1fun, M2fun);
+%! Mfcn = @(z) M \ z;
+%! x = bicgstab (Afcn, b, 1e-6, n, Mfcn);
+%! M1fcn = @(z) M1 \ z;
+%! M2fcn = @(z) M2 \ z;
+%! x = bicgstab (Afcn, b, 1e-6, n, M1fcn, M2fcn);
 %! function y = Ap (A, x, z)
 %!   ## compute A^z * x or (A^z)' * x
 %!   y = x;
@@ -378,8 +378,8 @@
 %!     y = A * y;
 %!   endfor
 %! endfunction
-%! Afun = @(x, p) Ap (A, x, p);
-%! x = bicgstab (Afun, b, [], 2 * n, [], [], [], 2); # solution of A^2 * x = b
+%! Afcn = @(x, p) Ap (A, x, p);
+%! x = bicgstab (Afcn, b, [], 2 * n, [], [], [], 2); # solution of A^2 * x = b
 
 %!demo
 %! n = 10;
@@ -406,28 +406,28 @@
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;
 %! maxit = 20;
-%! Afun = @(z) A*z;
-%! M1_fun = @(z) M1 \ z;
-%! M2_fun = @(z) M2 \ z;
+%! Afcn = @(z) A*z;
+%! M1_fcn = @(z) M1 \ z;
+%! M2_fcn = @(z) M2 \ z;
 %! [x, flag] = bicgstab (A,b );
 %! assert (flag, 0);
 %! [x, flag] = bicgstab (A, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (A, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = bicgstab (A, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (A, b, [], maxit, M1_fun, M2);
+%! [x, flag] = bicgstab (A, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (A, b, [], maxit, M1, M2_fun);
+%! [x, flag] = bicgstab (A, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (Afun, b);
+%! [x, flag] = bicgstab (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (Afun, b, [], maxit, M1, M2);
+%! [x, flag] = bicgstab (Afcn, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (Afun, b, [], maxit, M1_fun, M2);
+%! [x, flag] = bicgstab (Afcn, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (Afun, b, [], maxit, M1, M2_fun);
+%! [x, flag] = bicgstab (Afcn, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = bicgstab (Afun, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = bicgstab (Afcn, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!shared n, A, b, tol, maxit, M1, M2
@@ -443,12 +443,12 @@
 %! [x, flag, relres, iter, resvec] = bicgstab (A, b, tol, maxit, M1, M2);
 %! assert (norm (b - A*x) / norm (b), 0, tol);
 
-%!function y = afun (x, a)
+%!function y = afcn (x, a)
 %!  y = a * x;
 %!endfunction
 %!
 %!test
-%! [x, flag, relres, iter, resvec] = bicgstab (@(x) afun (x, A), b,
+%! [x, flag, relres, iter, resvec] = bicgstab (@(x) afcn (x, A), b,
 %!                                             tol, maxit, M1, M2);
 %! assert (norm (b - A*x) / norm (b), 0, tol);
 
@@ -501,13 +501,13 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x)
+%!function y = Afcn (x)
 %!  A = sparse (toeplitz ([2, 1, 0, 0], [2, -1, 0, 0]));
 %!  y = A * x;
 %!endfunction
 %!
 %! b = [1; 2; 2; 3];
-%! [x, flag] = bicgstab ("Afun", b);
+%! [x, flag] = bicgstab ("Afcn", b);
 %! assert (norm (b - A*x) / norm (b), 0, 1e-6);
 
 %!test
--- a/scripts/sparse/cgs.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/cgs.m	Mon Apr 04 18:14:56 2022 -0700
@@ -36,8 +36,8 @@
 ##
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline
-## function @code{Afun} such that @code{Afun(x) = A * x}.  Additional
-## parameters to @code{Afun} are passed after @var{x0}.
+## function @code{Afcn} such that @code{Afcn(x) = A * x}.  Additional
+## parameters to @code{Afcn} are passed after @var{x0}.
 ##
 ## @item @var{b} is the right hand side vector.  It must be a column vector
 ## with same number of rows of @var{A}.
@@ -107,10 +107,10 @@
 ## restart = 5;
 ## [M1, M2] = ilu (A); # in this tridiag case it corresponds to chol (A)'
 ## M = M1 * M2;
-## Afun = @@(x) A * x;
-## Mfun = @@(x) M \ x;
-## M1fun = @@(x) M1 \ x;
-## M2fun = @@(x) M2 \ x;
+## Afcn = @@(x) A * x;
+## Mfcn = @@(x) M \ x;
+## M1fcn = @@(x) M1 \ x;
+## M2fcn = @@(x) M2 \ x;
 ## @end group
 ## @end example
 ##
@@ -124,7 +124,7 @@
 ## @code{@var{A} * @var{x}}
 ##
 ## @example
-## x = cgs (Afun, b, [], n)
+## x = cgs (Afcn, b, [], n)
 ## @end example
 ##
 ## @sc{Example 3:} @code{cgs} with a preconditioner matrix @var{M}
@@ -136,7 +136,7 @@
 ## @sc{Example 4:} @code{cgs} with a function as preconditioner
 ##
 ## @example
-## x = cgs (Afun, b, 1e-6, n, Mfun)
+## x = cgs (Afcn, b, 1e-6, n, Mfcn)
 ## @end example
 ##
 ## @sc{Example 5:} @code{cgs} with preconditioner matrices @var{M1}
@@ -149,7 +149,7 @@
 ## @sc{Example 6:} @code{cgs} with functions as preconditioners
 ##
 ## @example
-## x = cgs (Afun, b, 1e-6, n, M1fun, M2fun)
+## x = cgs (Afcn, b, 1e-6, n, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 7:} @code{cgs} with as input a function requiring an argument
@@ -162,8 +162,8 @@
 ##      y = A * y;
 ##    endfor
 ##  endfunction
-## Apfun = @@(x, string, p) Ap (A, x, string, p);
-## x = cgs (Apfun, b, [], [], [], [], [], 2);
+## Apfcn = @@(x, string, p) Ap (A, x, string, p);
+## x = cgs (Apfcn, b, [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -196,7 +196,7 @@
 function [x_min, flag, relres, iter_min, resvec] = ...
          cgs (A, b, tol = [], maxit = [], M1 = [] , M2 = [], x0 = [], varargin)
 
-  [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, "cgs");
+  [Afcn, M1fcn, M2fcn] = __alltohandles__ (A, b, M1, M2, "cgs");
 
   [tol, maxit, x0] = __default__input__ ({1e-06, min( rows(b), 20), ...
                                           zeros(size (b))}, tol, maxit, x0);
@@ -224,21 +224,21 @@
   ## x_min approximation with the minimum residual
   ## x_pr approximation at the previous iteration (to check stagnation)
 
-  r0 = rr = u = p = b - feval (Afun, x, varargin{:});
+  r0 = rr = u = p = b - feval (Afcn, x, varargin{:});
   resvec (1) = norm (r0, 2);
   rho_1 = rr' * r0;
 
   try
     warning ("error","Octave:singular-matrix","local");
-    p_hat = feval (M1fun, p, varargin{:});
-    p_hat = feval (M2fun, p_hat, varargin {:});
+    p_hat = feval (M1fcn, p, varargin{:});
+    p_hat = feval (M2fcn, p_hat, varargin {:});
   catch
     flag = 2;
   end_try_catch
 
   while ((flag != 2) && (iter < maxit) && ...
          (resvec (iter + 1) >= tol * norm_b))
-    v = feval (Afun, p_hat, varargin{:});
+    v = feval (Afcn, p_hat, varargin{:});
     prod_tmp = (rr' * v);
     if (prod_tmp == 0)
       flag = 4;
@@ -246,10 +246,10 @@
     endif
     alpha = rho_1 / prod_tmp;
     q = u - alpha * v;
-    u_hat = feval(M1fun, u + q, varargin{:});
-    u_hat = feval (M2fun, u_hat, varargin{:});
+    u_hat = feval(M1fcn, u + q, varargin{:});
+    u_hat = feval (M2fcn, u_hat, varargin{:});
     x += alpha*u_hat;
-    r0 -= alpha* feval (Afun, u_hat, varargin{:});
+    r0 -= alpha* feval (Afcn, u_hat, varargin{:});
     iter += 1;
     resvec (iter + 1) = norm (r0, 2);
     if (norm (x - x_pr, 2) <= norm (x, 2) * eps) # Stagnation
@@ -270,8 +270,8 @@
     beta = rho_1 / rho_2;
     u = r0 + beta * q;
     p = u + beta * (q + beta * p);
-    p_hat = feval (M1fun, p, varargin {:});
-    p_hat = feval (M2fun, p_hat, varargin{:});
+    p_hat = feval (M1fcn, p, varargin {:});
+    p_hat = feval (M2fcn, p_hat, varargin{:});
   endwhile
   resvec = resvec (1: (iter + 1));
 
@@ -331,15 +331,15 @@
 %! [M1, M2] = ilu (A + 0.1 * eye (n));
 %! M = M1 * M2;
 %! x = cgs (A, b, [], n);
-%! Afun = @(x) A * x;
-%! x = cgs (Afun, b, [], n);
+%! Afcn = @(x) A * x;
+%! x = cgs (Afcn, b, [], n);
 %! x = cgs (A, b, 1e-6, n, M);
 %! x = cgs (A, b, 1e-6, n, M1, M2);
-%! Mfun = @(z) M \ z;
-%! x = cgs (Afun, b, 1e-6, n, Mfun);
-%! M1fun = @(z) M1 \ z;
-%! M2fun = @(z) M2 \ z;
-%! x = cgs (Afun, b, 1e-6, n, M1fun, M2fun);
+%! Mfcn = @(z) M \ z;
+%! x = cgs (Afcn, b, 1e-6, n, Mfcn);
+%! M1fcn = @(z) M1 \ z;
+%! M2fcn = @(z) M2 \ z;
+%! x = cgs (Afcn, b, 1e-6, n, M1fcn, M2fcn);
 %! function y = Ap (A, x, z)
 %!   ## compute A^z * x or (A^z)' * x
 %!   y = x;
@@ -347,8 +347,8 @@
 %!     y = A * y;
 %!   endfor
 %! endfunction
-%! Afun = @(x, p) Ap (A, x, p);
-%! x = cgs (Afun, b, [], 2*n, [], [], [], 2); # solution of A^2 * x = b
+%! Afcn = @(x, p) Ap (A, x, p);
+%! x = cgs (Afcn, b, [], 2*n, [], [], [], 2); # solution of A^2 * x = b
 
 %!demo
 %! n = 10;
@@ -375,28 +375,28 @@
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;
 %! maxit = 10;
-%! Afun = @(z) A * z;
-%! M1_fun = @(z) M1 \ z;
-%! M2_fun = @(z) M2 \ z;
+%! Afcn = @(z) A * z;
+%! M1_fcn = @(z) M1 \ z;
+%! M2_fcn = @(z) M2 \ z;
 %! [x, flag] = cgs (A,b);
 %! assert (flag, 0);
 %! [x, flag] = cgs (A, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = cgs (A, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = cgs (A, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = cgs (A, b, [], maxit, M1_fun, M2);
+%! [x, flag] = cgs (A, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = cgs (A, b, [], maxit, M1, M2_fun);
+%! [x, flag] = cgs (A, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = cgs (Afun, b);
+%! [x, flag] = cgs (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = cgs (Afun, b, [], maxit, M1, M2);
+%! [x, flag] = cgs (Afcn, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = cgs (Afun, b, [], maxit, M1_fun, M2);
+%! [x, flag] = cgs (Afcn, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = cgs (Afun, b, [], maxit, M1, M2_fun);
+%! [x, flag] = cgs (Afcn, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = cgs (Afun, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = cgs (Afcn, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!shared n, A, b, tol, maxit, M
@@ -451,11 +451,11 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x)
+%!function y = Afcn (x)
 %!  A = toeplitz ([2, 1, 0, 0], [2, -1, 0, 0]);
 %!  y = A * x;
 %!endfunction
-%! [x, flag] = cgs ("Afun", [1; 2; 2; 3]);
+%! [x, flag] = cgs ("Afcn", [1; 2; 2; 3]);
 %! assert (norm (b - A*x) / norm (b), 0, 1e-6);
 
 %!test
--- a/scripts/sparse/eigs.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/eigs.m	Mon Apr 04 18:14:56 2022 -0700
@@ -534,25 +534,25 @@
 %!testif HAVE_ARPACK, HAVE_UMFPACK
 %! assert (eigs (A, k, 4.1), eigs (A, speye (n), k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (d1, d0(end:-1:(end-k+1)), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (d1, d0(k:-1:1), 1e-11);
 %!testif HAVE_ARPACK, HAVE_UMFPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (d1, eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK, HAVE_UMFPACK, HAVE_CHOLMOD
 %! AA = speye (10);
-%! fn = @(x) AA * x;
+%! fcn = @(x) AA * x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! assert (eigs (fn, 10, AA, 3, "lm", opts), [1; 1; 1], 10*eps);
+%! assert (eigs (fcn, 10, AA, 3, "lm", opts), [1; 1; 1], 10*eps);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
 %! d1 = diag (d1);
@@ -662,19 +662,19 @@
 %! assert (sort (imag (eigs (A, k, 4.1))),
 %!         sort (imag (eigs (A, speye (n), k, 4.1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (abs (d1), abs (d0(end:-1:(end-k+1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (abs (d1), d0(1:k), 1e-11);
 %!testif HAVE_ARPACK, HAVE_UMFPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (abs (d1), eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
@@ -794,19 +794,19 @@
 %! assert (sort (imag (eigs (A, k, 4.1))),
 %!         sort (imag (eigs (A, speye (n), k, 4.1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (abs (d1), abs (d0(end:-1:(end-k+1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (abs (d1), d0(1:k), 1e-11);
 %!testif HAVE_ARPACK, HAVE_UMFPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (abs (d1), eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
@@ -1033,19 +1033,19 @@
 %!testif HAVE_ARPACK
 %! assert (eigs (A, k, 4.1), eigs (A, eye (n), k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (d1, d0(end:-1:(end-k+1)), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (d1, d0(k:-1:1), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 1;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (d1, eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
@@ -1156,19 +1156,19 @@
 %! assert (sort (imag (eigs (A, k, 4.1))),
 %!         sort (imag (eigs (A, eye (n), k, 4.1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (abs (d1), abs (d0(end:-1:(end-k+1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (abs (d1), d0(1:k), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 0;  opts.isreal = 1;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (abs (d1), eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
@@ -1287,19 +1287,19 @@
 %! assert (sort (imag (eigs (A, k, 4.1))),
 %!         sort (imag (eigs (A, eye (n), k, 4.1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A * x;
+%! fcn = @(x) A * x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, "lm", opts);
+%! d1 = eigs (fcn, n, k, "lm", opts);
 %! assert (abs (d1), abs (d0(end:-1:(end-k+1))), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) A \ x;
+%! fcn = @(x) A \ x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, "sm", opts);
+%! d1 = eigs (fcn, n, k, "sm", opts);
 %! assert (abs (d1), d0(1:k), 1e-11);
 %!testif HAVE_ARPACK
-%! fn = @(x) (A - 4.1 * eye (n)) \ x;
+%! fcn = @(x) (A - 4.1 * eye (n)) \ x;
 %! opts.issym = 0;  opts.isreal = 0;
-%! d1 = eigs (fn, n, k, 4.1, opts);
+%! d1 = eigs (fcn, n, k, 4.1, opts);
 %! assert (abs (d1), eigs (A, k, 4.1), 1e-11);
 %!testif HAVE_ARPACK
 %! [v1,d1] = eigs (A, k, "lm");
@@ -1452,12 +1452,12 @@
 %!testif HAVE_ARPACK
 %! A = toeplitz ([-2, 1, zeros(1, 8)]);
 %! A = kron (A, eye (10)) + kron (eye (10), A);
-%! Afun = @(x) A * x;
+%! Afcn = @(x) A * x;
 %! opts.v0 = (1:100)';
 %! opts.maxit = 3;
 %! opts.issym = true;
 %! warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
-%! d = eigs (Afun, 100, 4, "sm", opts);
+%! d = eigs (Afcn, 100, 4, "sm", opts);
 %! assert (d(3:4), [NaN; NaN]);
 %!testif HAVE_ARPACK
 %! A = toeplitz ([-2, 1, zeros(1, 8)]);
@@ -1510,11 +1510,11 @@
 %!testif HAVE_ARPACK
 %! A = toeplitz ([0, 1, zeros(1, 8)], [0, -1, zeros(1, 8)]);
 %! A = kron (A, eye (10)) + kron (eye (10), A);
-%! Afun = @(x) A * x;
+%! Afcn = @(x) A * x;
 %! opts.v0 = (1:100)';
 %! opts.maxit = 4;
 %! warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
-%! d = eigs (Afun, 100, 4, "lm", opts);
+%! d = eigs (Afcn, 100, 4, "lm", opts);
 %! assert (d(3:4), [NaN+1i*NaN; NaN+1i*NaN]);
 %!testif HAVE_ARPACK
 %! A = 1i * magic (100);
@@ -1532,12 +1532,12 @@
 %! assert (d(9:10), [NaN+1i*NaN; NaN+1i*NaN]);
 %!testif HAVE_ARPACK
 %! A = 1i * magic (100);
-%! Afun = @(x) A * x;
+%! Afcn = @(x) A * x;
 %! opts.v0 = (1:100)';
 %! opts.maxit = 1;
 %! opts.isreal = false;
 %! warning ("off", "Octave:eigs:UnconvergedEigenvalues", "local");
-%! d = eigs (Afun, 100, 6, "lm", opts);
+%! d = eigs (Afcn, 100, 6, "lm", opts);
 %! assert (d(6), NaN+1i*NaN);
 %!testif HAVE_ARPACK, HAVE_CHOLMOD
 %! A = sparse (magic (10));
@@ -1567,8 +1567,8 @@
 %! B = B * B';
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "LM", opts);
-%! Afun = @(x) A * x;
-%! [Evector_f Evalues_f] = eigs (Afun, 10, B, 4, "LM", opts);
+%! Afcn = @(x) A * x;
+%! [Evector_f Evalues_f] = eigs (Afcn, 10, B, 4, "LM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 %!testif HAVE_ARPACK
@@ -1579,8 +1579,8 @@
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "SM", opts);
 %! [L, U, P] = lu (A);
-%! Afun = @(x) U \ (L \ (P * x));
-%! [Evector_f Evalues_f] = eigs (Afun, 10, B, 4, "SM", opts);
+%! Afcn = @(x) U \ (L \ (P * x));
+%! [Evector_f Evalues_f] = eigs (Afcn, 10, B, 4, "SM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 %!testif HAVE_ARPACK
@@ -1591,9 +1591,9 @@
 %! B = B * B';
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "LM", opts);
-%! Afun = @(x) A * x;
+%! Afcn = @(x) A * x;
 %! opts.issym = true;
-%! [Evector_f Evalues_f] = eigs (Afun, 10, B, 4, "LM", opts);
+%! [Evector_f Evalues_f] = eigs (Afcn, 10, B, 4, "LM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 %!testif HAVE_ARPACK
@@ -1605,9 +1605,9 @@
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "SM", opts);
 %! [L, U, P] = lu (A);
-%! Afun = @(x) U \ (L \ (P * x));
+%! Afcn = @(x) U \ (L \ (P * x));
 %! opts.issym = true;
-%! [Evector_f Evalues_f] = eigs (Afun, 10, B, 4, "SM", opts);
+%! [Evector_f Evalues_f] = eigs (Afcn, 10, B, 4, "SM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 %!testif HAVE_ARPACK
@@ -1617,9 +1617,9 @@
 %! B = B * B';
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "LM", opts);
-%! Afun = @(x) A * x;
+%! Afcn = @(x) A * x;
 %! opts.isreal = false;
-%! [Evector_f Evalues_f] = eigs (Afun, 10, B, 4, "LM", opts);
+%! [Evector_f Evalues_f] = eigs (Afcn, 10, B, 4, "LM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 %!testif HAVE_ARPACK
@@ -1630,9 +1630,9 @@
 %! opts.v0 = (1:10)';
 %! [Evector, Evalues] = eigs (A, B, 4, "SM", opts);
 %! [L, U, P] = lu (A);
-%! Afun = @(x) U \ (L \ (P *x));
+%! Afcn = @(x) U \ (L \ (P *x));
 %! opts.isreal = false;
-%! [Evector_f, Evalues_f] = eigs (Afun, 10, B, 4, "SM", opts);
+%! [Evector_f, Evalues_f] = eigs (Afcn, 10, B, 4, "SM", opts);
 %! assert (Evector, Evector_f);
 %! assert (Evalues, Evalues_f);
 
--- a/scripts/sparse/gmres.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/gmres.m	Mon Apr 04 18:14:56 2022 -0700
@@ -36,8 +36,8 @@
 ##
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline
-## function @code{Afun} such that @code{Afun(x) = A * x}.  Additional
-## parameters to @code{Afun} are passed after @var{x0}.
+## function @code{Afcn} such that @code{Afcn(x) = A * x}.  Additional
+## parameters to @code{Afcn} are passed after @var{x0}.
 ##
 ## @item @var{b} is the right hand side vector.  It must be a column vector
 ## with the same numbers of rows as @var{A}.
@@ -141,10 +141,10 @@
 ## restart = 5;
 ## [M1, M2] = ilu (A); # in this tridiag case, it corresponds to lu (A)
 ## M = M1 * M2;
-## Afun = @@(x) A * x;
-## Mfun = @@(x) M \ x;
-## M1fun = @@(x) M1 \ x;
-## M2fun = @@(x) M2 \ x;
+## Afcn = @@(x) A * x;
+## Mfcn = @@(x) M \ x;
+## M1fcn = @@(x) M1 \ x;
+## M2fcn = @@(x) M2 \ x;
 ## @end group
 ## @end example
 ##
@@ -158,7 +158,7 @@
 ## @code{@var{A} * @var{x}}
 ##
 ## @example
-## x = gmres (Afun, b, [], [], n)
+## x = gmres (Afcn, b, [], [], n)
 ## @end example
 ##
 ## @sc{Example 3:} usage of @code{gmres} with the restart
@@ -180,7 +180,7 @@
 ## @sc{Example 5:} @code{gmres} with a function as preconditioner
 ##
 ## @example
-## x = gmres (Afun, b, [], 1e-6, n, Mfun)
+## x = gmres (Afcn, b, [], 1e-6, n, Mfcn)
 ## @end example
 ##
 ## @sc{Example 6:} @code{gmres} with preconditioner matrices @var{M1}
@@ -193,7 +193,7 @@
 ## @sc{Example 7:} @code{gmres} with functions as preconditioners
 ##
 ## @example
-## x = gmres (Afun, b, 1e-6, n, M1fun, M2fun)
+## x = gmres (Afcn, b, 1e-6, n, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 8:} @code{gmres} with as input a function requiring an argument
@@ -206,8 +206,8 @@
 ##        y = A * y;
 ##      endfor
 ##   endfunction
-## Apfun = @@(x, p) Ap (A, x, p);
-## x = gmres (Apfun, b, [], [], [], [], [], [], 2);
+## Apfcn = @@(x, p) Ap (A, x, p);
+## x = gmres (Apfcn, b, [], [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -248,7 +248,7 @@
     class_name = "double";
   endif
 
-  [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, "gmres");
+  [Afcn, M1fcn, M2fcn] = __alltohandles__ (A, b, M1, M2, "gmres");
 
   ## Check if the inputs are empty, and in case set them
   [tol, x0] = __default__input__ ({1e-06, zeros(size (b))}, tol, x0);
@@ -329,15 +329,15 @@
   flag = 1; # Default flag is maximum # of iterations exceeded
 
   ## begin loop
-  u = feval (Afun, x_old, varargin{:});
+  u = feval (Afcn, x_old, varargin{:});
   try
     warning ("error", "Octave:singular-matrix", "local");
-    prec_res = feval (M1fun, b - u, varargin{:});  # M1*(b-u)
-    prec_res = feval (M2fun, prec_res, varargin{:});
+    prec_res = feval (M1fcn, b - u, varargin{:});  # M1*(b-u)
+    prec_res = feval (M2fcn, prec_res, varargin{:});
     presn = norm (prec_res, 2);
     resvec(1) = presn;
-    z = feval (M1fun, b, varargin{:});
-    z = feval (M2fun, z, varargin{:});
+    z = feval (M1fcn, b, varargin{:});
+    z = feval (M2fcn, z, varargin{:});
     prec_b_norm = norm (z, 2);
     B (1) = presn;
     V(:, 1) = prec_res / presn;
@@ -352,18 +352,18 @@
       restart_it = 1;
       outer_it += 1;
       x_old = x;
-      u = feval (Afun, x_old, varargin{:});
-      prec_res = feval (M1fun, b - u, varargin{:});
-      prec_res = feval (M2fun, prec_res, varargin{:});
+      u = feval (Afcn, x_old, varargin{:});
+      prec_res = feval (M1fcn, b - u, varargin{:});
+      prec_res = feval (M2fcn, prec_res, varargin{:});
       presn = norm (prec_res, 2);
       B(1) = presn;
       H(:) = 0;
       V(:, 1) = prec_res / presn;
     endif
     ## basic iteration
-    u = feval (Afun, V(:, restart_it), varargin{:});
-    tmp = feval (M1fun, u, varargin{:});
-    tmp = feval (M2fun, tmp, varargin{:});
+    u = feval (Afcn, V(:, restart_it), varargin{:});
+    tmp = feval (M1fcn, u, varargin{:});
+    tmp = feval (M2fcn, tmp, varargin{:});
     [V(:,restart_it + 1), H(1:restart_it + 1, restart_it)] = ...
       mgorth (tmp, V(:,1:restart_it));
     Y = (H(1:restart_it + 1, 1:restart_it) \ B(1:restart_it + 1));
@@ -469,23 +469,23 @@
 %! M = M1 * M2;
 %! x = gmres (A, b, [], [], n);
 %! x = gmres (A, b, restart, [], n);  # gmres with restart
-%! Afun = @(x) A * x;
-%! x = gmres (Afun, b, [], [], n);
+%! Afcn = @(x) A * x;
+%! x = gmres (Afcn, b, [], [], n);
 %! x = gmres (A, b, [], 1e-6, n, M);  # gmres without restart
 %! x = gmres (A, b, [], 1e-6, n, M1, M2);
-%! Mfun = @(x) M \ x;
-%! x = gmres (Afun, b, [], 1e-6, n, Mfun);
-%! M1fun = @(x) M1 \ x;
-%! M2fun = @(x) M2 \ x;
-%! x = gmres (Afun, b, [], 1e-6, n, M1fun, M2fun);
+%! Mfcn = @(x) M \ x;
+%! x = gmres (Afcn, b, [], 1e-6, n, Mfcn);
+%! M1fcn = @(x) M1 \ x;
+%! M2fcn = @(x) M2 \ x;
+%! x = gmres (Afcn, b, [], 1e-6, n, M1fcn, M2fcn);
 %! function y = Ap (A, x, p)  # compute A^p * x
 %!    y = x;
 %!    for i = 1:p
 %!      y = A * y;
 %!    endfor
 %!  endfunction
-%! Afun = @(x, p) Ap (A, x, p);
-%! x = gmres (Afun, b, [], [], n, [], [], [], 2);  # solution of A^2 * x = b
+%! Afcn = @(x, p) Ap (A, x, p);
+%! x = gmres (Afcn, b, [], [], n, [], [], [], 2);  # solution of A^2 * x = b
 
 %!demo
 %! n = 10;
@@ -510,28 +510,28 @@
 %! b = sum (A, 2);
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;
-%! Afun = @(z) A * z;
-%! M1_fun = @(z) M1 \ z;
-%! M2_fun = @(z) M2 \ z;
+%! Afcn = @(z) A * z;
+%! M1_fcn = @(z) M1 \ z;
+%! M2_fcn = @(z) M2 \ z;
 %! [x, flag] = gmres (A, b);
 %! assert (flag, 0);
 %! [x, flag] = gmres (A, b, [], [], [], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = gmres (A, b, [], [], [], M1_fun, M2_fun);
+%! [x, flag] = gmres (A, b, [], [], [], M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = gmres (A, b, [], [], [], M1_fun, M2);
+%! [x, flag] = gmres (A, b, [], [], [], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = gmres (A, b, [], [], [], M1, M2_fun);
+%! [x, flag] = gmres (A, b, [], [], [], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = gmres (Afun, b);
+%! [x, flag] = gmres (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = gmres (Afun, b, [],[],[], M1, M2);
+%! [x, flag] = gmres (Afcn, b, [],[],[], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = gmres (Afun, b, [],[],[], M1_fun, M2);
+%! [x, flag] = gmres (Afcn, b, [],[],[], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = gmres (Afun, b, [],[],[], M1, M2_fun);
+%! [x, flag] = gmres (Afcn, b, [],[],[], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = gmres (Afun, b, [],[],[], M1_fun, M2_fun);
+%! [x, flag] = gmres (Afcn, b, [],[],[], M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!test
@@ -634,11 +634,11 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x)
+%!function y = Afcn (x)
 %!   A = toeplitz ([2, 1, 0, 0], [2, -1, 0, 0]);
 %!   y = A * x;
 %!endfunction
-%! [x, flag] = gmres ("Afun", [1; 2; 2; 3]);
+%! [x, flag] = gmres ("Afcn", [1; 2; 2; 3]);
 %! assert (x, ones (4, 1), 1e-6);
 
 %!test # preconditioned residual
--- a/scripts/sparse/pcg.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/pcg.m	Mon Apr 04 18:14:56 2022 -0700
@@ -36,8 +36,8 @@
 ## @itemize
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline function
-## @code{Afun} such that @code{Afun(x) = A * x}.  Additional parameters to
-## @code{Afun} may be passed after @var{x0}.
+## @code{Afcn} such that @code{Afcn(x) = A * x}.  Additional parameters to
+## @code{Afcn} may be passed after @var{x0}.
 ##
 ## @var{A} has to be Hermitian and Positive Definite (@nospell{HPD})@.  If
 ## @code{pcg} detects @var{A} not to be positive definite, a warning is printed
@@ -163,10 +163,10 @@
 ## M1 = ichol (A); # in this tridiagonal case it corresponds to chol (A)'
 ## M2 = M1';
 ## M = M1 * M2;
-## Afun = @@(x) A * x;
-## Mfun = @@(x) M \ x;
-## M1fun = @@(x) M1 \ x;
-## M2fun = @@(x) M2 \ x;
+## Afcn = @@(x) A * x;
+## Mfcn = @@(x) M \ x;
+## M1fcn = @@(x) M1 \ x;
+## M2fcn = @@(x) M2 \ x;
 ## @end group
 ## @end example
 ##
@@ -180,7 +180,7 @@
 ## @code{@var{A} * @var{x}}
 ##
 ## @example
-## x = pcg (Afun, b)
+## x = pcg (Afcn, b)
 ## @end example
 ##
 ## @sc{Example 3:} @code{pcg} with a preconditioner matrix @var{M}
@@ -192,7 +192,7 @@
 ## @sc{Example 4:} @code{pcg} with a function as preconditioner
 ##
 ## @example
-## x = pcg (Afun, b, 1e-6, 100, Mfun)
+## x = pcg (Afcn, b, 1e-6, 100, Mfcn)
 ## @end example
 ##
 ## @sc{Example 5:} @code{pcg} with preconditioner matrices @var{M1}
@@ -205,7 +205,7 @@
 ## @sc{Example 6:} @code{pcg} with functions as preconditioners
 ##
 ## @example
-## x = pcg (Afun, b, 1e-6, 100, M1fun, M2fun)
+## x = pcg (Afcn, b, 1e-6, 100, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 7:} @code{pcg} with as input a function requiring an argument
@@ -218,8 +218,8 @@
 ##        y = A * y;
 ##      endfor
 ##   endfunction
-## Apfun = @@(x, p) Ap (A, x, p);
-## x = pcg (Apfun, b, [], [], [], [], [], 2);
+## Apfcn = @@(x, p) Ap (A, x, p);
+## x = pcg (Apfcn, b, [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -275,7 +275,7 @@
   ## Check if the input data A,b,m1,m2 are consistent (i.e. if they are
   ## matrix or function handle)
 
-  [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, "pcg");
+  [Afcn, M1fcn, M2fcn] = __alltohandles__ (A, b, M1, M2, "pcg");
 
   maxit += 2;
   n_arg_out = nargout;
@@ -301,7 +301,7 @@
   ## x_pr (x previous) needs to check the stagnation
   ## x_min needs to save the iterated with minimum residual
 
-  r = b - feval (Afun, x, varargin{:});
+  r = b - feval (Afcn, x, varargin{:});
   iter = 2;
   iter_min = 0;
   flag = 1;
@@ -320,15 +320,15 @@
     if (iter == 2) # Check whether M1 or M2 are singular
       try
         warning ("error","Octave:singular-matrix","local");
-        z = feval (M1fun, r, varargin{:});
-        z = feval (M2fun, z, varargin{:});
+        z = feval (M1fcn, r, varargin{:});
+        z = feval (M2fcn, z, varargin{:});
       catch
         flag = 2;
         break;
       end_try_catch
     else
-      z = feval (M1fun, r, varargin{:});
-      z = feval (M2fun, z, varargin{:});
+      z = feval (M1fcn, r, varargin{:});
+      z = feval (M2fcn, z, varargin{:});
     endif
 
     tau = z' * r;
@@ -336,7 +336,7 @@
     beta = tau / old_tau;
     old_tau = tau;
     p = z + beta * p;
-    w = feval (Afun, p, varargin{:});
+    w = feval (Afcn, p, varargin{:});
 
     ## Needed only for eigest.
 
@@ -378,8 +378,8 @@
   if (n_arg_out > 5)
   ## Apply the preconditioner once more and finish with the precond
   ## residual.
-    z = feval (M1fun, r, varargin{:});
-    z = feval (M2fun, z, varargin{:});
+    z = feval (M1fcn, r, varargin{:});
+    z = feval (M2fcn, z, varargin{:});
   endif
 
   ## (Eventually) computes the eigenvalue of inv(m2)*inv(m1)*A
@@ -466,24 +466,24 @@
 %! M2 = M1';
 %! M = M1 * M2;
 %! x = pcg (A, b);
-%! Afun = @(x) A * x;
-%! x = pcg (Afun, b);
+%! Afcn = @(x) A * x;
+%! x = pcg (Afcn, b);
 %! x = pcg (A, b, 1e-6, 100, M);
 %! x = pcg (A, b, 1e-6, 100, M1, M2);
-%! Mfun = @(x) M \ x;
-%! x = pcg (Afun, b, 1e-6, 100, Mfun);
-%! M1fun = @(x) M1 \ x;
-%! M2fun = @(x) M2 \ x;
-%! x = pcg (Afun, b, 1e-6, 100, M1fun, M2fun);
+%! Mfcn = @(x) M \ x;
+%! x = pcg (Afcn, b, 1e-6, 100, Mfcn);
+%! M1fcn = @(x) M1 \ x;
+%! M2fcn = @(x) M2 \ x;
+%! x = pcg (Afcn, b, 1e-6, 100, M1fcn, M2fcn);
 %! function y = Ap (A, x, p)  # compute A^p * x
 %!    y = x;
 %!    for i = 1:p
 %!      y = A * y;
 %!    endfor
 %!  endfunction
-%! Afun = @(x, p) Ap (A, x, p);
+%! Afcn = @(x, p) Ap (A, x, p);
 %! ## solution of A^2 * x = b
-%! x = pcg (Afun, b, [], [], [], [], [], 2);
+%! x = pcg (Afcn, b, [], [], [], [], [], 2);
 
 %!demo
 %! n = 10;
@@ -506,32 +506,32 @@
 %! b = A * ones (5, 1);
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;  # M1 * M2 is the Jacobi preconditioner
-%! Afun = @(z) A*z;
-%! M1_fun = @(z) M1 \ z;
-%! M2_fun = @(z) M2 \ z;
+%! Afcn = @(z) A*z;
+%! M1_fcn = @(z) M1 \ z;
+%! M2_fcn = @(z) M2 \ z;
 %! [x, flag, ~, iter] = pcg (A,b);
 %! assert (flag, 0);
 %! [x, flag, ~ , iter] = pcg (A, b, [], [], M1 * M2);
 %! assert (flag, 0);
 %! [x, flag, ~ , iter] = pcg (A, b, [], [], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = pcg (A, b, [], [], M1_fun, M2_fun);
+%! [x, flag] = pcg (A, b, [], [], M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = pcg (A, b,[],[], M1_fun, M2);
+%! [x, flag] = pcg (A, b,[],[], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = pcg (A, b,[],[], M1, M2_fun);
+%! [x, flag] = pcg (A, b,[],[], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b);
+%! [x, flag] = pcg (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b,[],[], M1 * M2);
+%! [x, flag] = pcg (Afcn, b,[],[], M1 * M2);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b,[],[], M1, M2);
+%! [x, flag] = pcg (Afcn, b,[],[], M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b,[],[], M1_fun, M2);
+%! [x, flag] = pcg (Afcn, b,[],[], M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b,[],[], M1, M2_fun);
+%! [x, flag] = pcg (Afcn, b,[],[], M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = pcg (Afun, b,[],[], M1_fun, M2_fun);
+%! [x, flag] = pcg (Afcn, b,[],[], M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!test
@@ -658,11 +658,11 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x)
+%!function y = Afcn (x)
 %!   A = toeplitz ([2, 1, 0, 0]);
 %!   y = A * x;
 %!endfunction
-%! [x, flag] = pcg ("Afun", [3; 4; 4; 3]);
+%! [x, flag] = pcg ("Afcn", [3; 4; 4; 3]);
 %! assert (x, ones (4, 1), 1e-6);
 
 %!test
--- a/scripts/sparse/private/__alltohandles__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/private/__alltohandles__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -24,7 +24,7 @@
 ########################################################################
 
 ## -*- texinfo -*-
-## @deftypefn {} {[@var{Afun}, @var{M1fun}, @var{M2fun}] =} __alltohandles__ (@var{A}, @var{b}, @var{M1}, @var{M2}, @var{solver_name})
+## @deftypefn {} {[@var{Afcn}, @var{M1fcn}, @var{M2fcn}] =} __alltohandles__ (@var{A}, @var{b}, @var{M1}, @var{M2}, @var{solver_name})
 ##
 ## Check if the parameters @var{A} (matrix of our linear system), @var{b}
 ## (right hand side vector), @var{M1}, @var{M2} (preconditioner matrices) are
@@ -47,13 +47,13 @@
 ##
 ## @itemize
 ##
-## @item @var{Afun}, @var{M1fun}, @var{M2fun} are the corresponding
+## @item @var{Afcn}, @var{M1fcn}, @var{M2fcn} are the corresponding
 ## function handles.
 ##
 ## @end itemize
 ## @end deftypefn
 
-function [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, solver_name)
+function [Afcn, M1fcn, M2fcn] = __alltohandles__ (A, b, M1, M2, solver_name)
 
   A_is_numeric = false;
   M1_is_numeric = false;
@@ -61,9 +61,9 @@
 
   ## Check A and set its type
   if (is_function_handle (A))
-     Afun = A;
+     Afcn = A;
   elseif (ischar (A))
-    Afun = str2func (A);
+    Afcn = str2func (A);
   elseif (! isnumeric (A) || ! issquare (A))
     error ([solver_name, ": A must be a square matrix or a function handle"]);
   else
@@ -78,18 +78,18 @@
     switch (solver_name)
       case {"pcg", "gmres", "bicgstab", "cgs", "tfqmr"}
         ## methods which do not require the transpose
-        M1fun = @(x) x;
+        M1fcn = @(x) x;
       case {"bicg"}
         ## methods which do require the transpose
-        M1fun = @(x, ~) x;
+        M1fcn = @(x, ~) x;
       otherwise
         error (["__alltohandles__: unknown method: ", solver_name]);
     endswitch
   else # M1 not empty
     if (is_function_handle (M1))
-      M1fun = M1;
+      M1fcn = M1;
     elseif (ischar (M1))
-      M1fun = str2func (M1);
+      M1fcn = str2func (M1);
     elseif (! isnumeric (M1) || ! issquare (M1))
       error ([solver_name, ": M1 must be a square matrix or a function handle"]);
     else
@@ -101,18 +101,18 @@
     switch (solver_name)
       case {"pcg", "gmres", "bicgstab", "cgs", "tfqmr"}
         ## methods which do not require the transpose
-        M2fun = @(x) x;
+        M2fcn = @(x) x;
       case {"bicg"}
         ## methods which do require the transpose
-        M2fun = @(x, ~) x;
+        M2fcn = @(x, ~) x;
       otherwise
         error (["__alltohandles__: unknown method: ", solver_name]);
     endswitch
   else # M2 not empty
     if (is_function_handle (M2))
-      M2fun = M2;
+      M2fcn = M2;
     elseif (ischar (M2))
-      M2fun = str2func (M2);
+      M2fcn = str2func (M2);
     elseif (! isnumeric (M2) || ! issquare (M2))
       error ([solver_name, ": M2 must be a square matrix or a function handle"]);
     else
@@ -124,24 +124,24 @@
     case {"pcg", "gmres", "bicgstab", "cgs", "tfqmr"}
       ## methods which do not require the transpose
       if (A_is_numeric)
-        Afun = @(x) A * x;
+        Afcn = @(x) A * x;
       endif
       if (M1_is_numeric)
-        M1fun = @(x) M1 \ x;
+        M1fcn = @(x) M1 \ x;
       endif
       if (M2_is_numeric)
-        M2fun = @(x) M2 \ x;
+        M2fcn = @(x) M2 \ x;
       endif
     case {"bicg"}
       ## methods which do require the transpose
       if (A_is_numeric)
-        Afun = @(x, trans) A_sub (A, x, trans);
+        Afcn = @(x, trans) A_sub (A, x, trans);
       endif
       if (M1_is_numeric)
-        M1fun = @(x, trans) M_sub (M1, x, trans);
+        M1fcn = @(x, trans) M_sub (M1, x, trans);
       endif
       if (M2_is_numeric)
-        M2fun = @(x, trans) M_sub (M2, x, trans);
+        M2fcn = @(x, trans) M_sub (M2, x, trans);
       endif
     otherwise
       error (["__alltohandles__: unknown method: ", solver_name]);
--- a/scripts/sparse/private/__sprand__.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/private/__sprand__.m	Mon Apr 04 18:14:56 2022 -0700
@@ -27,9 +27,9 @@
 ## public domain.
 
 ## -*- texinfo -*-
-## @deftypefn  {} {@var{S} =} __sprand__ (@var{s}, @var{randfun})
-## @deftypefnx {} {@var{S} =} __sprand__ (@var{m}, @var{n}, @var{d}, @var{fcnname}, @var{randfun})
-## @deftypefnx {} {@var{S} =} __sprand__ (@var{m}, @var{n}, @var{d}, @var{rc}, @var{fcnname}, @var{randfun})
+## @deftypefn  {} {@var{S} =} __sprand__ (@var{s}, @var{randfcn})
+## @deftypefnx {} {@var{S} =} __sprand__ (@var{m}, @var{n}, @var{d}, @var{fcnname}, @var{randfcn})
+## @deftypefnx {} {@var{S} =} __sprand__ (@var{m}, @var{n}, @var{d}, @var{rc}, @var{fcnname}, @var{randfcn})
 ## Undocumented internal function.
 ## @end deftypefn
 
@@ -38,15 +38,15 @@
 function S = __sprand__ (varargin)
 
   if (nargin == 2)
-    [m, randfun] = deal (varargin{1:2});
+    [m, randfcn] = deal (varargin{1:2});
     [i, j] = find (m);
     [nr, nc] = size (m);
-    S = sparse (i, j, randfun (size (i)), nr, nc);
+    S = sparse (i, j, randfcn (size (i)), nr, nc);
   else
     if (nargin == 5)
-      [m, n, d, fcnname, randfun] = deal (varargin{:});
+      [m, n, d, fcnname, randfcn] = deal (varargin{:});
     else
-      [m, n, d, rc, fcnname, randfun] = deal (varargin{:});
+      [m, n, d, rc, fcnname, randfcn] = deal (varargin{:});
     endif
 
     if (! (isscalar (m) && m == fix (m) && m >= 0))
@@ -87,7 +87,7 @@
         [i, j] = ind2sub ([m, n], idx);
       endif
 
-      S = sparse (i, j, randfun (k, 1), m, n);
+      S = sparse (i, j, randfcn (k, 1), m, n);
 
     elseif (nargin == 6)
       ## Create a matrix with specified reciprocal condition number.
--- a/scripts/sparse/qmr.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/qmr.m	Mon Apr 04 18:14:56 2022 -0700
@@ -97,9 +97,9 @@
   if (nargin >= 2 && isvector (full (b)))
 
     if (ischar (A))
-      fun = str2func (A);
-      Ax  = @(x) feval (fun, x, "notransp");
-      Atx = @(x) feval (fun, x, "transp");
+      fcn = str2func (A);
+      Ax  = @(x) feval (fcn, x, "notransp");
+      Atx = @(x) feval (fcn, x, "transp");
     elseif (is_function_handle (A))
       Ax  = @(x) feval (A, x, "notransp");
       Atx = @(x) feval (A, x, "transp");
@@ -124,9 +124,9 @@
       M1m1x = @(x, ignore) x;
       M1tm1x = M1m1x;
     elseif (ischar (M1))
-      fun = str2func (M1);
-      M1m1x  = @(x) feval (fun, x, "notransp");
-      M1tm1x = @(x) feval (fun, x, "transp");
+      fcn = str2func (M1);
+      M1m1x  = @(x) feval (fcn, x, "notransp");
+      M1tm1x = @(x) feval (fcn, x, "transp");
     elseif (is_function_handle (M1))
       M1m1x  = @(x) feval (M1, x, "notransp");
       M1tm1x = @(x) feval (M1, x, "transp");
@@ -141,9 +141,9 @@
       M2m1x = @(x, ignore) x;
       M2tm1x = M2m1x;
     elseif (ischar (M2))
-      fun = str2func (M2);
-      M2m1x  = @(x) feval (fun, x, "notransp");
-      M2tm1x = @(x) feval (fun, x, "transp");
+      fcn = str2func (M2);
+      M2m1x  = @(x) feval (fcn, x, "notransp");
+      M2tm1x = @(x) feval (fcn, x, "transp");
     elseif (is_function_handle (M2))
       M2m1x  = @(x) feval (M2, x, "notransp");
       M2tm1x = @(x) feval (M2, x, "transp");
@@ -290,7 +290,7 @@
 %! [x, flag, relres, iter, resvec] = qmr (A, b, rtol, maxit, M1, M2);
 %! assert (x, ones (size (b)), 1e-7);
 
-%!function y = afun (x, t, a)
+%!function y = afcn (x, t, a)
 %!  switch (t)
 %!    case "notransp"
 %!      y = a * x;
@@ -308,7 +308,7 @@
 %! M1 = spdiags ([ones(n,1)/(-2) ones(n,1)],-1:0, n, n);
 %! M2 = spdiags ([4*ones(n,1) -ones(n,1)], 0:1, n, n);
 %!
-%! [x, flag, relres, iter, resvec] = qmr (@(x, t) afun (x, t, A),
+%! [x, flag, relres, iter, resvec] = qmr (@(x, t) afcn (x, t, A),
 %!                                         b, rtol, maxit, M1, M2);
 %! assert (x, ones (size (b)), 1e-7);
 
--- a/scripts/sparse/tfqmr.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/sparse/tfqmr.m	Mon Apr 04 18:14:56 2022 -0700
@@ -35,8 +35,8 @@
 ##
 ## @item @var{A} is the matrix of the linear system and it must be square.
 ## @var{A} can be passed as a matrix, function handle, or inline
-## function @code{Afun} such that @code{Afun(x) = A * x}.  Additional
-## parameters to @code{Afun} are passed after @var{x0}.
+## function @code{Afcn} such that @code{Afcn(x) = A * x}.  Additional
+## parameters to @code{Afcn} are passed after @var{x0}.
 ##
 ## @item @var{b} is the right hand side vector.  It must be a column vector
 ## with the same number of rows as @var{A}.
@@ -115,10 +115,10 @@
 ## restart = 5;
 ## [M1, M2] = ilu (A); # in this tridiag case it corresponds to chol (A)'
 ## M = M1 * M2;
-## Afun = @@(x) A * x;
-## Mfun = @@(x) M \ x;
-## M1fun = @@(x) M1 \ x;
-## M2fun = @@(x) M2 \ x;
+## Afcn = @@(x) A * x;
+## Mfcn = @@(x) M \ x;
+## M1fcn = @@(x) M1 \ x;
+## M2fcn = @@(x) M2 \ x;
 ## @end group
 ## @end example
 ##
@@ -132,7 +132,7 @@
 ## @code{@var{A} * @var{x}}
 ##
 ## @example
-## x = tfqmr (Afun, b, [], n)
+## x = tfqmr (Afcn, b, [], n)
 ## @end example
 ##
 ## @sc{Example 3:} @code{tfqmr} with a preconditioner matrix @var{M}
@@ -144,7 +144,7 @@
 ## @sc{Example 4:} @code{tfqmr} with a function as preconditioner
 ##
 ## @example
-## x = tfqmr (Afun, b, 1e-6, n, Mfun)
+## x = tfqmr (Afcn, b, 1e-6, n, Mfcn)
 ## @end example
 ##
 ## @sc{Example 5:} @code{tfqmr} with preconditioner matrices @var{M1}
@@ -157,7 +157,7 @@
 ## @sc{Example 6:} @code{tfmqr} with functions as preconditioners
 ##
 ## @example
-## x = tfqmr (Afun, b, 1e-6, n, M1fun, M2fun)
+## x = tfqmr (Afcn, b, 1e-6, n, M1fcn, M2fcn)
 ## @end example
 ##
 ## @sc{Example 7:} @code{tfqmr} with as input a function requiring an argument
@@ -170,8 +170,8 @@
 ##      y = A * y;
 ##    endfor
 ##  endfunction
-## Apfun = @@(x, string, p) Ap (A, x, string, p);
-## x = tfqmr (Apfun, b, [], [], [], [], [], 2);
+## Apfcn = @@(x, string, p) Ap (A, x, string, p);
+## x = tfqmr (Apfcn, b, [], [], [], [], [], 2);
 ## @end group
 ## @end example
 ##
@@ -206,7 +206,7 @@
          tfqmr (A, b, tol = [], maxit = [], M1 = [], M2 = [], ...
                 x0 = [], varargin)
 
-  [Afun, M1fun, M2fun] = __alltohandles__ (A, b, M1, M2, "tfqmr");
+  [Afcn, M1fcn, M2fcn] = __alltohandles__ (A, b, M1, M2, "tfqmr");
 
   [tol, maxit, x0] = __default__input__ ({1e-06, 2 * min(20, rows (b)), ...
                                           zeros(rows (b), 1)}, tol, ...
@@ -233,7 +233,7 @@
   resvec = zeros (maxit, 1);
   flag = 1;
 
-  w = u = r = r_star = b - feval (Afun, x0, varargin{:});
+  w = u = r = r_star = b - feval (Afcn, x0, varargin{:});
   rho_1 = (r_star' * r);
   d = 0;
   tau = norm (r, 2);
@@ -243,9 +243,9 @@
 
   try
     warning ("error", "Octave:singular-matrix", "local");
-    u_hat = feval (M1fun, u, varargin{:});
-    u_hat = feval (M2fun, u_hat, varargin{:});
-    v = feval (Afun, u_hat, varargin{:});
+    u_hat = feval (M1fcn, u, varargin{:});
+    u_hat = feval (M2fcn, u_hat, varargin{:});
+    v = feval (Afcn, u_hat, varargin{:});
   catch
     flag = 2;
   end_try_catch
@@ -262,16 +262,16 @@
       alpha = rho_1 / v_r;
       u_1 = u - alpha * v;  # u at the after iteration
     endif
-    u_hat = feval (M1fun, u, varargin{:});
-    u_hat = feval (M2fun, u_hat, varargin{:});
-    w -= alpha * feval (Afun, u_hat, varargin{:});
+    u_hat = feval (M1fcn, u, varargin{:});
+    u_hat = feval (M2fcn, u_hat, varargin{:});
+    w -= alpha * feval (Afcn, u_hat, varargin{:});
     d = u_hat + ((theta * theta) / alpha) * eta * d;
     theta = norm (w, 2) / tau;
     c = 1 / sqrt (1 + theta * theta);
     tau *= theta * c;
     eta = (c * c) * alpha;
     x += eta * d;
-    r -= eta * feval (Afun, d, varargin{:});
+    r -= eta * feval (Afcn, d, varargin{:});
     if (it < 0) # iter is odd
       rho_2 = rho_1;
       rho_1 = (r_star' * w);
@@ -283,10 +283,10 @@
       endif
       beta = rho_1 / rho_2;
       u_1 = w + beta * u; # u at the after iteration
-      u1_hat = feval (M1fun, u_1, varargin{:});
-      u1_hat = feval (M2fun, u1_hat, varargin{:});
-      v = feval (Afun, u1_hat, varargin{:}) + ...
-          beta * (feval (Afun, u_hat, varargin{:}) + beta * v);
+      u1_hat = feval (M1fcn, u_1, varargin{:});
+      u1_hat = feval (M2fcn, u1_hat, varargin{:});
+      v = feval (Afcn, u1_hat, varargin{:}) + ...
+          beta * (feval (Afcn, u_hat, varargin{:}) + beta * v);
     endif
     u = u_1;
     iter += 1;
@@ -355,28 +355,28 @@
 %! M1 = diag (sqrt (diag (A)));
 %! M2 = M1;
 %! maxit = 10;
-%! Afun = @(z) A * z;
-%! M1_fun = @(z) M1 \ z;
-%! M2_fun = @(z) M2 \ z;
+%! Afcn = @(z) A * z;
+%! M1_fcn = @(z) M1 \ z;
+%! M2_fcn = @(z) M2 \ z;
 %! [x, flag] = tfqmr (A,b);
 %! assert (flag, 0);
 %! [x, flag] = tfqmr (A, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (A, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = tfqmr (A, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (A, b, [], maxit, M1_fun, M2);
+%! [x, flag] = tfqmr (A, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (A, b, [], maxit, M1, M2_fun);
+%! [x, flag] = tfqmr (A, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (Afun, b);
+%! [x, flag] = tfqmr (Afcn, b);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (Afun, b, [], maxit, M1, M2);
+%! [x, flag] = tfqmr (Afcn, b, [], maxit, M1, M2);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (Afun, b, [], maxit, M1_fun, M2);
+%! [x, flag] = tfqmr (Afcn, b, [], maxit, M1_fcn, M2);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (Afun, b, [], maxit, M1, M2_fun);
+%! [x, flag] = tfqmr (Afcn, b, [], maxit, M1, M2_fcn);
 %! assert (flag, 0);
-%! [x, flag] = tfqmr (Afun, b, [], maxit, M1_fun, M2_fun);
+%! [x, flag] = tfqmr (Afcn, b, [], maxit, M1_fcn, M2_fcn);
 %! assert (flag, 0);
 
 %!shared A, b, n, M1, M2
@@ -393,15 +393,15 @@
 %! assert (x, ones (size (b)), 1e-7);
 %!
 %!test
-%!function y = afun (x, a)
+%!function y = afcn (x, a)
 %!  y = a * x;
 %!endfunction
 %!
 %! tol = 1e-8;
 %! maxit = 15;
 %!
-%! [x, flag, relres, iter, resvec] = tfqmr (@(x) afun (x, A), b,
-%!                                             tol, maxit, M1, M2);
+%! [x, flag, relres, iter, resvec] = tfqmr (@(x) afcn (x, A), b,
+%!                                          tol, maxit, M1, M2);
 %! assert (x, ones (size (b)), 1e-7);
 
 %!test
@@ -458,11 +458,11 @@
 %! assert (class (x), "single");
 
 %!test
-%!function y = Afun (x)
+%!function y = Afcn (x)
 %!  A = toeplitz ([2, 1, 0, 0], [2, -1, 0, 0]);
 %!  y = A * x;
 %!endfunction
-%! [x, flag] = tfqmr ("Afun", [1; 2; 2; 3]);
+%! [x, flag] = tfqmr ("Afcn", [1; 2; 2; 3]);
 %! assert (x, ones (4, 1), 1e-6);
 
 %!test # unpreconditioned residual
@@ -481,23 +481,23 @@
 %! [M1, M2] = ilu (A + 0.1 * eye (n));
 %! M = M1 * M2;
 %! x = tfqmr (A, b, [], n);
-%! Afun = @(x) A * x;
-%! x = tfqmr (Afun, b, [], n);
+%! Afcn = @(x) A * x;
+%! x = tfqmr (Afcn, b, [], n);
 %! x = tfqmr (A, b, 1e-6, n, M);
 %! x = tfqmr (A, b, 1e-6, n, M1, M2);
-%! Mfun = @(z) M \ z;
-%! x = tfqmr (Afun, b, 1e-6, n, Mfun);
-%! M1fun = @(z) M1 \ z;
-%! M2fun = @(z) M2 \ z;
-%! x = tfqmr (Afun, b, 1e-6, n, M1fun, M2fun);
+%! Mfcn = @(z) M \ z;
+%! x = tfqmr (Afcn, b, 1e-6, n, Mfcn);
+%! M1fcn = @(z) M1 \ z;
+%! M2fcn = @(z) M2 \ z;
+%! x = tfqmr (Afcn, b, 1e-6, n, M1fcn, M2fcn);
 %! function y = Ap (A, x, z) # compute A^z * x or (A^z)' * x
 %!    y = x;
 %!    for i = 1:z
 %!      y = A * y;
 %!    endfor
 %!  endfunction
-%! Afun = @(x, p) Ap (A, x, p);
-%! x = tfqmr (Afun, b, [], 2*n, [], [], [], 2); # solution of A^2 * x = b
+%! Afcn = @(x, p) Ap (A, x, p);
+%! x = tfqmr (Afcn, b, [], 2*n, [], [], [], 2); # solution of A^2 * x = b
 
 %!demo
 %! n = 10;
--- a/scripts/testfun/private/dump_demos.m	Mon Apr 04 11:22:26 2022 -0700
+++ b/scripts/testfun/private/dump_demos.m	Mon Apr 04 18:14:56 2022 -0700
@@ -196,11 +196,11 @@
     demos = get_demos (fcn);
     for d = 1:numel (demos)
       idx = sprintf ("%02d", d);
-      base_fn = sprintf ("%s_%s", fcn, idx);
-      fn = sprintf ('%s.%s', base_fn, fmt);
+      base_fcn = sprintf ("%s_%s", fcn, idx);
+      fn = sprintf ('%s.%s', base_fcn, fmt);
       ## Wrap each demo in a function which create a local scope
       ## to prevent that a previous demo overwrites i or pi, for example
-      fprintf (fid, "\nfunction %s ()\n", base_fn);
+      fprintf (fid, "\nfunction %s ()\n", base_fcn);
       fprintf (fid, "    %s\n\n", demos{d});
       fprintf (fid, "end\n\n");
 
@@ -212,7 +212,7 @@
       ## here (https://savannah.gnu.org/bugs/?42557).
       fprintf (fid, "    rand ('seed', 1);\n");
       fprintf (fid, "    tic ();\n");
-      fprintf (fid, "    %s ();\n", base_fn);
+      fprintf (fid, "    %s ();\n", base_fcn);
       fprintf (fid, "    t_plot = toc ();\n");
       fprintf (fid, "    fig = (get (0, 'currentfigure'));\n");
       fprintf (fid, "    if (~ isempty (fig))\n");
@@ -229,8 +229,8 @@
       fprintf (fid, "    fprintf ('File ""%s"" already exists.\\n');\n", fn);
       fprintf (fid, "  end\n");
       fprintf (fid, "catch\n");
-      fprintf (fid, "  fprintf ('ERROR in %s: %%s\\n', lasterr ());\n", base_fn);
-      fprintf (fid, "  err_fid = fopen ('%s.err', 'w');\n", base_fn);
+      fprintf (fid, "  fprintf ('ERROR in %s: %%s\\n', lasterr ());\n", base_fcn);
+      fprintf (fid, "  err_fid = fopen ('%s.err', 'w');\n", base_fcn);
       fprintf (fid, "  fprintf (err_fid, '%%s', lasterr ());\n");
       fprintf (fid, "  fclose (err_fid);\n");
       fprintf (fid, "end\n");