GNU Octave  4.2.1
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
data.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 1994-2017 John W. Eaton
4 Copyright (C) 2009 Jaroslav Hajek
5 Copyright (C) 2009-2010 VZLU Prague
6 Copyright (C) 2012 Carlo de Falco
7 
8 This file is part of Octave.
9 
10 Octave is free software; you can redistribute it and/or modify it
11 under the terms of the GNU General Public License as published by the
12 Free Software Foundation; either version 3 of the License, or (at your
13 option) any later version.
14 
15 Octave is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
18 for more details.
19 
20 You should have received a copy of the GNU General Public License
21 along with Octave; see the file COPYING. If not, see
22 <http://www.gnu.org/licenses/>.
23 
24 */
25 
26 #if defined (HAVE_CONFIG_H)
27 # include "config.h"
28 #endif
29 
30 #include <sys/types.h>
31 
32 #include <cfloat>
33 #include <ctime>
34 
35 #include <string>
36 
37 #include "lo-ieee.h"
38 #include "lo-math.h"
39 #include "oct-base64.h"
40 #include "oct-time.h"
41 #include "str-vec.h"
42 #include "quit.h"
43 #include "mx-base.h"
44 #include "oct-binmap.h"
45 
46 #include "Cell.h"
47 #include "defun.h"
48 #include "error.h"
49 #include "errwarn.h"
50 #include "oct-map.h"
51 #include "ovl.h"
52 #include "ov.h"
53 #include "ov-class.h"
54 #include "ov-float.h"
55 #include "ov-complex.h"
56 #include "ov-flt-complex.h"
57 #include "ov-cx-mat.h"
58 #include "ov-flt-cx-mat.h"
59 #include "ov-cx-sparse.h"
60 #include "parse.h"
61 #include "pt-mat.h"
62 #include "utils.h"
63 #include "variables.h"
64 #include "pager.h"
65 #include "xnorm.h"
66 
67 #if ! defined (HAVE_HYPOTF) && defined (HAVE__HYPOTF)
68 # define hypotf _hypotf
69 # define HAVE_HYPOTF 1
70 #endif
71 
72 static void
73 index_error (const char *fmt, const std::string& idx, const std::string& msg)
74 {
75  error (fmt, idx.c_str (), msg.c_str ());
76 }
77 
78 DEFUN (all, args, ,
79  doc: /* -*- texinfo -*-
80 @deftypefn {} {} all (@var{x})
81 @deftypefnx {} {} all (@var{x}, @var{dim})
82 For a vector argument, return true (logical 1) if all elements of the vector
83 are nonzero.
84 
85 For a matrix argument, return a row vector of logical ones and
86 zeros with each element indicating whether all of the elements of the
87 corresponding column of the matrix are nonzero. For example:
88 
89 @example
90 @group
91 all ([2, 3; 1, 0])
92  @result{} [ 1, 0 ]
93 @end group
94 @end example
95 
96 If the optional argument @var{dim} is supplied, work along dimension
97 @var{dim}.
98 @seealso{any}
99 @end deftypefn */)
100 {
101  int nargin = args.length ();
102 
103  if (nargin < 1 || nargin > 2)
104  print_usage ();
105 
106  int dim = (nargin == 1 ? -1
107  : args(1).int_value ("all: DIM must be an integer")-1);
108 
109  if (dim < -1)
110  error ("all: invalid dimension argument = %d", dim + 1);
111 
112  return ovl (args(0).all (dim));
113 }
114 
115 /*
116 %!test
117 %! x = ones (3);
118 %! x(1,1) = 0;
119 %! assert (all (all (rand (3) + 1) == [1, 1, 1]) == 1);
120 %! assert (all (all (x) == [0, 1, 1]) == 1);
121 %! assert (all (x, 1) == [0, 1, 1]);
122 %! assert (all (x, 2) == [0; 1; 1]);
123 
124 %!test
125 %! x = ones (3, "single");
126 %! x(1,1) = 0;
127 %! assert (all (all (single (rand (3) + 1)) == [1, 1, 1]) == 1);
128 %! assert (all (all (x) == [0, 1, 1]) == 1);
129 %! assert (all (x, 1) == [0, 1, 1]);
130 %! assert (all (x, 2) == [0; 1; 1]);
131 
132 %!error all ()
133 %!error all (1, 2, 3)
134 */
135 
136 DEFUN (any, args, ,
137  doc: /* -*- texinfo -*-
138 @deftypefn {} {} any (@var{x})
139 @deftypefnx {} {} any (@var{x}, @var{dim})
140 For a vector argument, return true (logical 1) if any element of the vector
141 is nonzero.
142 
143 For a matrix argument, return a row vector of logical ones and
144 zeros with each element indicating whether any of the elements of the
145 corresponding column of the matrix are nonzero. For example:
146 
147 @example
148 @group
149 any (eye (2, 4))
150  @result{} [ 1, 1, 0, 0 ]
151 @end group
152 @end example
153 
154 If the optional argument @var{dim} is supplied, work along dimension
155 @var{dim}. For example:
156 
157 @example
158 @group
159 any (eye (2, 4), 2)
160  @result{} [ 1; 1 ]
161 @end group
162 @end example
163 @seealso{all}
164 @end deftypefn */)
165 {
166  int nargin = args.length ();
167 
168  if (nargin < 1 || nargin > 2)
169  print_usage ();
170 
171  int dim = (nargin == 1 ? -1
172  : args(1).int_value ("any: DIM must be an integer")-1);
173 
174  if (dim < -1)
175  error ("any: invalid dimension argument = %d", dim + 1);
176 
177  return ovl (args(0).any (dim));
178 }
179 
180 /*
181 %!test
182 %! x = zeros (3);
183 %! x(3,3) = 1;
184 %! assert (all (any (x) == [0, 0, 1]) == 1);
185 %! assert (all (any (ones (3)) == [1, 1, 1]) == 1);
186 %! assert (any (x, 1) == [0, 0, 1]);
187 %! assert (any (x, 2) == [0; 0; 1]);
188 
189 %!test
190 %! x = zeros (3, "single");
191 %! x(3,3) = 1;
192 %! assert (all (any (x) == [0, 0, 1]) == 1);
193 %! assert (all (any (ones (3, "single")) == [1, 1, 1]) == 1);
194 %! assert (any (x, 1) == [0, 0, 1]);
195 %! assert (any (x, 2) == [0; 0; 1]);
196 
197 %!error any ()
198 %!error any (1, 2, 3)
199 */
200 
201 // These mapping functions may also be useful in other places, eh?
202 
203 DEFUN (atan2, args, ,
204  doc: /* -*- texinfo -*-
205 @deftypefn {} {} atan2 (@var{y}, @var{x})
206 Compute atan (@var{y} / @var{x}) for corresponding elements of @var{y}
207 and @var{x}.
208 
209 @var{y} and @var{x} must match in size and orientation. The signs of
210 elements of @var{y} and @var{x} are used to determine the quadrants of each
211 resulting value.
212 
213 This function is equivalent to @code{arg (complex (@var{x}, @var{y}))}.
214 @seealso{tan, tand, tanh, atanh}
215 @end deftypefn */)
216 {
217  if (args.length () != 2)
218  print_usage ();
219 
221 
222  if (! args(0).is_numeric_type ())
223  err_wrong_type_arg ("atan2", args(0));
224 
225  if (! args(1).is_numeric_type ())
226  err_wrong_type_arg ("atan2", args(1));
227 
228  if (args(0).is_complex_type () || args(1).is_complex_type ())
229  error ("atan2: not defined for complex numbers");
230 
231  if (args(0).is_single_type () || args(1).is_single_type ())
232  {
233  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
234  retval = atan2f (args(0).float_value (), args(1).float_value ());
235  else
236  {
237  FloatNDArray a0 = args(0).float_array_value ();
238  FloatNDArray a1 = args(1).float_array_value ();
239  retval = binmap<float> (a0, a1, ::atan2f, "atan2");
240  }
241  }
242  else
243  {
244  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
245  retval = atan2 (args(0).scalar_value (), args(1).scalar_value ());
246  else if (args(0).is_sparse_type ())
247  {
248  SparseMatrix m0 = args(0).sparse_matrix_value ();
249  SparseMatrix m1 = args(1).sparse_matrix_value ();
250  retval = binmap<double> (m0, m1, ::atan2, "atan2");
251  }
252  else
253  {
254  NDArray a0 = args(0).array_value ();
255  NDArray a1 = args(1).array_value ();
256  retval = binmap<double> (a0, a1, ::atan2, "atan2");
257  }
258  }
259 
260  return retval;
261 }
262 
263 /*
264 %!assert (size (atan2 (zeros (0, 2), zeros (0, 2))), [0, 2])
265 %!assert (size (atan2 (rand (2, 3, 4), zeros (2, 3, 4))), [2, 3, 4])
266 %!assert (size (atan2 (rand (2, 3, 4), 1)), [2, 3, 4])
267 %!assert (size (atan2 (1, rand (2, 3, 4))), [2, 3, 4])
268 %!assert (size (atan2 (1, 2)), [1, 1])
269 
270 %!test
271 %! rt2 = sqrt (2);
272 %! rt3 = sqrt (3);
273 %! v = [0, pi/6, pi/4, pi/3, -pi/3, -pi/4, -pi/6, 0];
274 %! y = [0, rt3, 1, rt3, -rt3, -1, -rt3, 0];
275 %! x = [1, 3, 1, 1, 1, 1, 3, 1];
276 %! assert (atan2 (y, x), v, sqrt (eps));
277 
278 %!test
279 %! rt2 = sqrt (2);
280 %! rt3 = sqrt (3);
281 %! v = single ([0, pi/6, pi/4, pi/3, -pi/3, -pi/4, -pi/6, 0]);
282 %! y = single ([0, rt3, 1, rt3, -rt3, -1, -rt3, 0]);
283 %! x = single ([1, 3, 1, 1, 1, 1, 3, 1]);
284 %! assert (atan2 (y, x), v, sqrt (eps ("single")));
285 
286 ## Test sparse implementations
287 %!shared xs
288 %! xs = sparse (0:3);
289 %!test
290 %! y = atan2 (1, xs);
291 %! assert (issparse (y), false);
292 %! assert (nnz (y), 4);
293 %! assert (y, atan2 (1, 0:3));
294 %!test
295 %! y = atan2 (0, xs);
296 %! assert (issparse (y), false);
297 %! assert (nnz (y), 0);
298 %! assert (y, zeros (1,4));
299 %!test
300 %! y = atan2 (xs, 1);
301 %! assert (issparse (y));
302 %! assert (nnz (y), 3);
303 %! assert (y, sparse (atan2 (0:3, 1)));
304 %!test
305 %! y = atan2 (xs, 0);
306 %! assert (issparse (y));
307 %! assert (nnz (y), 3);
308 %! assert (y, sparse (atan2 (0:3, 0)));
309 %!test
310 %! y = atan2 (xs, sparse (ones (1, 4)));
311 %! assert (issparse (y));
312 %! assert (nnz (y), 3);
313 %! assert (y, sparse (atan2 (0:3, ones (1,4))));
314 %!test
315 %! y = atan2 (xs, sparse (zeros (1,4)));
316 %! assert (issparse (y));
317 %! assert (nnz (y), 3);
318 %! assert (y, sparse (atan2 (0:3, zeros (1,4))));
319 
320 %!error atan2 ()
321 %!error atan2 (1, 2, 3)
322 */
323 
324 static octave_value
325 do_hypot (const octave_value& x, const octave_value& y)
326 {
328 
329  octave_value arg0 = x;
330  octave_value arg1 = y;
331  if (! arg0.is_numeric_type ())
332  err_wrong_type_arg ("hypot", arg0);
333  if (! arg1.is_numeric_type ())
334  err_wrong_type_arg ("hypot", arg1);
335 
336  if (arg0.is_complex_type ())
337  arg0 = arg0.abs ();
338  if (arg1.is_complex_type ())
339  arg1 = arg1.abs ();
340 
341  if (arg0.is_single_type () || arg1.is_single_type ())
342  {
343  if (arg0.is_scalar_type () && arg1.is_scalar_type ())
344  retval = hypotf (arg0.float_value (), arg1.float_value ());
345  else
346  {
347  FloatNDArray a0 = arg0.float_array_value ();
348  FloatNDArray a1 = arg1.float_array_value ();
349  retval = binmap<float> (a0, a1, ::hypotf, "hypot");
350  }
351  }
352  else
353  {
354  if (arg0.is_scalar_type () && arg1.is_scalar_type ())
355  retval = hypot (arg0.scalar_value (), arg1.scalar_value ());
356  else if (arg0.is_sparse_type () || arg1.is_sparse_type ())
357  {
358  SparseMatrix m0 = arg0.sparse_matrix_value ();
359  SparseMatrix m1 = arg1.sparse_matrix_value ();
360  retval = binmap<double> (m0, m1, ::hypot, "hypot");
361  }
362  else
363  {
364  NDArray a0 = arg0.array_value ();
365  NDArray a1 = arg1.array_value ();
366  retval = binmap<double> (a0, a1, ::hypot, "hypot");
367  }
368  }
369 
370  return retval;
371 }
372 
373 DEFUN (hypot, args, ,
374  doc: /* -*- texinfo -*-
375 @deftypefn {} {} hypot (@var{x}, @var{y})
376 @deftypefnx {} {} hypot (@var{x}, @var{y}, @var{z}, @dots{})
377 Compute the element-by-element square root of the sum of the squares of
378 @var{x} and @var{y}.
379 
380 This is equivalent to
381 @code{sqrt (@var{x}.^2 + @var{y}.^2)}, but is calculated in a manner that
382 avoids overflows for large values of @var{x} or @var{y}.
383 
384 @code{hypot} can also be called with more than 2 arguments; in this case,
385 the arguments are accumulated from left to right:
386 
387 @example
388 @group
389 hypot (hypot (@var{x}, @var{y}), @var{z})
390 hypot (hypot (hypot (@var{x}, @var{y}), @var{z}), @var{w}), etc.
391 @end group
392 @end example
393 @end deftypefn */)
394 {
395  int nargin = args.length ();
396 
397  if (nargin < 2)
398  print_usage ();
399 
401 
402  if (nargin == 2)
403  retval = do_hypot (args(0), args(1));
404  else
405  {
406  retval = args(0);
407 
408  for (int i = 1; i < nargin; i++)
409  retval = do_hypot (retval, args(i));
410  }
411 
412  return retval;
413 }
414 
415 /*
416 %!assert (size (hypot (zeros (0, 2), zeros (0, 2))), [0, 2])
417 %!assert (size (hypot (rand (2, 3, 4), zeros (2, 3, 4))), [2, 3, 4])
418 %!assert (size (hypot (rand (2, 3, 4), 1)), [2, 3, 4])
419 %!assert (size (hypot (1, rand (2, 3, 4))), [2, 3, 4])
420 %!assert (size (hypot (1, 2)), [1, 1])
421 %!assert (hypot (1:10, 1:10), sqrt (2) * [1:10], 16*eps)
422 %!assert (hypot (single (1:10), single (1:10)), single (sqrt (2) * [1:10]))
423 
424 ## Test sparse implementations
425 %!shared xs
426 %! xs = sparse (0:3);
427 %!test
428 %! y = hypot (1, xs);
429 %! assert (nnz (y), 4);
430 %! assert (y, sparse (hypot (1, 0:3)));
431 %!test
432 %! y = hypot (0, xs);
433 %! assert (nnz (y), 3);
434 %! assert (y, xs);
435 %!test
436 %! y = hypot (xs, 1);
437 %! assert (nnz (y), 4);
438 %! assert (y, sparse (hypot (0:3, 1)));
439 %!test
440 %! y = hypot (xs, 0);
441 %! assert (nnz (y), 3);
442 %! assert (y, xs);
443 %!test
444 %! y = hypot (sparse ([0 0]), sparse ([0 1]));
445 %! assert (nnz (y), 1);
446 %! assert (y, sparse ([0 1]));
447 %!test
448 %! y = hypot (sparse ([0 1]), sparse ([0 0]));
449 %! assert (nnz (y), 1);
450 %! assert (y, sparse ([0 1]));
451 
452 */
453 
454 template <typename T, typename ET>
455 void
456 map_2_xlog2 (const Array<T>& x, Array<T>& f, Array<ET>& e)
457 {
458  f = Array<T>(x.dims ());
459  e = Array<ET>(x.dims ());
460  for (octave_idx_type i = 0; i < x.numel (); i++)
461  {
462  int exp;
463  f.xelem (i) = octave::math::log2 (x(i), exp);
464  e.xelem (i) = exp;
465  }
466 }
467 
469  doc: /* -*- texinfo -*-
470 @deftypefn {} {} log2 (@var{x})
471 @deftypefnx {} {[@var{f}, @var{e}] =} log2 (@var{x})
472 Compute the base-2 logarithm of each element of @var{x}.
473 
474 If called with two output arguments, split @var{x} into
475 binary mantissa and exponent so that
476 @tex
477 ${1 \over 2} \le \left| f \right| < 1$
478 @end tex
479 @ifnottex
480 @w{@code{1/2 <= abs(f) < 1}}
481 @end ifnottex
482 and @var{e} is an integer. If
483 @tex
484 $x = 0$, $f = e = 0$.
485 @end tex
486 @ifnottex
487 @w{@code{x = 0}}, @w{@code{f = e = 0}}.
488 @end ifnottex
489 @seealso{pow2, log, log10, exp}
490 @end deftypefn */)
491 {
492  if (args.length () != 1)
493  print_usage ();
494 
496 
497  if (nargout < 2)
498  retval = ovl (args(0).log2 ());
499  else if (args(0).is_single_type ())
500  {
501  if (args(0).is_real_type ())
502  {
503  FloatNDArray f;
504  FloatNDArray x = args(0).float_array_value ();
505  // FIXME: should E be an int value?
506  FloatMatrix e;
507  map_2_xlog2 (x, f, e);
508  retval = ovl (f, e);
509  }
510  else if (args(0).is_complex_type ())
511  {
513  FloatComplexNDArray x = args(0).float_complex_array_value ();
514  // FIXME: should E be an int value?
515  FloatNDArray e;
516  map_2_xlog2 (x, f, e);
517  retval = ovl (f, e);
518  }
519  }
520  else if (args(0).is_real_type ())
521  {
522  NDArray f;
523  NDArray x = args(0).array_value ();
524  // FIXME: should E be an int value?
525  Matrix e;
526  map_2_xlog2 (x, f, e);
527  retval = ovl (f, e);
528  }
529  else if (args(0).is_complex_type ())
530  {
532  ComplexNDArray x = args(0).complex_array_value ();
533  // FIXME: should E be an int value?
534  NDArray e;
535  map_2_xlog2 (x, f, e);
536  retval = ovl (f, e);
537  }
538  else
539  err_wrong_type_arg ("log2", args(0));
540 
541  return retval;
542 }
543 
544 /*
545 %!assert (log2 ([1/4, 1/2, 1, 2, 4]), [-2, -1, 0, 1, 2])
546 %!assert (log2 (Inf), Inf)
547 %!assert (isnan (log2 (NaN)))
548 %!assert (log2 (4*i), 2 + log2 (1*i))
549 %!assert (log2 (complex (0,Inf)), Inf + log2 (i))
550 
551 %!test
552 %! [f, e] = log2 ([0,-1; 2,-4; Inf,-Inf]);
553 %! assert (f, [0,-0.5; 0.5,-0.5; Inf,-Inf]);
554 %! assert (e(1:2,:), [0,1;2,3]);
555 
556 %!test
557 %! [f, e] = log2 (complex (zeros (3, 2), [0,-1; 2,-4; Inf,-Inf]));
558 %! assert (f, complex (zeros (3, 2), [0,-0.5; 0.5,-0.5; Inf,-Inf]));
559 %! assert (e(1:2,:), [0,1; 2,3]);
560 
561 %!assert <42583> (all (log2 (pow2 (-1074:1023)) == -1074:1023))
562 */
563 
564 DEFUN (rem, args, ,
565  doc: /* -*- texinfo -*-
566 @deftypefn {} {} rem (@var{x}, @var{y})
567 Return the remainder of the division @code{@var{x} / @var{y}}.
568 
569 The remainder is computed using the expression
570 
571 @example
572 x - y .* fix (x ./ y)
573 @end example
574 
575 An error message is printed if the dimensions of the arguments do not agree,
576 or if either argument is complex.
577 
578 Programming Notes: Floating point numbers within a few eps of an integer
579 will be rounded to an integer before computation for compatibility with
580 @sc{matlab}.
581 
582 By convention,
583 
584 @example
585 @group
586 rem (@var{x}, 0) = NaN if @var{x} is a floating point variable
587 rem (@var{x}, 0) = 0 if @var{x} is an integer variable
588 rem (@var{x}, @var{y}) returns a value with the signbit from @var{x}
589 @end group
590 @end example
591 
592 For the opposite conventions see the @code{mod} function. In general,
593 @code{rem} is best when computing the remainder after division of two
594 @emph{positive} numbers. For negative numbers, or when the values are
595 periodic, @code{mod} is a better choice.
596 @seealso{mod}
597 @end deftypefn */)
598 {
599  if (args.length () != 2)
600  print_usage ();
601 
603 
604  if (! args(0).is_numeric_type ())
605  err_wrong_type_arg ("rem", args(0));
606 
607  if (! args(1).is_numeric_type ())
608  err_wrong_type_arg ("rem", args(1));
609 
610  if (args(0).is_complex_type () || args(1).is_complex_type ())
611  error ("rem: not defined for complex numbers");
612 
613  if (args(0).is_integer_type () || args(1).is_integer_type ())
614  {
615  builtin_type_t btyp0 = args(0).builtin_type ();
616  builtin_type_t btyp1 = args(1).builtin_type ();
617  if (btyp0 == btyp_double || btyp0 == btyp_float)
618  btyp0 = btyp1;
619  if (btyp1 == btyp_double || btyp1 == btyp_float)
620  btyp1 = btyp0;
621 
622  if (btyp0 != btyp1)
623  error ("rem: cannot combine %s and %d",
624  args(0).class_name ().c_str (),
625  args(1).class_name ().c_str ());
626 
627  switch (btyp0)
628  {
629 #define MAKE_INT_BRANCH(X) \
630  case btyp_ ## X: \
631  { \
632  X##NDArray a0 = args(0).X##_array_value (); \
633  X##NDArray a1 = args(1).X##_array_value (); \
634  retval = binmap<octave_##X,octave_##X,octave_##X> (a0, a1, rem, "rem"); \
635  } \
636  break
637 
638  MAKE_INT_BRANCH (int8);
642  MAKE_INT_BRANCH (uint8);
645  MAKE_INT_BRANCH (uint64);
646 
647 #undef MAKE_INT_BRANCH
648 
649  default:
650  panic_impossible ();
651  }
652  }
653  else if (args(0).is_single_type () || args(1).is_single_type ())
654  {
655  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
656  retval = octave::math::rem (args(0).float_value (), args(1).float_value ());
657  else
658  {
659  FloatNDArray a0 = args(0).float_array_value ();
660  FloatNDArray a1 = args(1).float_array_value ();
661  retval = binmap<float> (a0, a1, octave::math::rem<float>, "rem");
662  }
663  }
664  else
665  {
666  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
667  retval = octave::math::rem (args(0).scalar_value (), args(1).scalar_value ());
668  else if (args(0).is_sparse_type () || args(1).is_sparse_type ())
669  {
670  SparseMatrix m0 = args(0).sparse_matrix_value ();
671  SparseMatrix m1 = args(1).sparse_matrix_value ();
672  retval = binmap<double> (m0, m1, octave::math::rem<double>, "rem");
673  }
674  else
675  {
676  NDArray a0 = args(0).array_value ();
677  NDArray a1 = args(1).array_value ();
678  retval = binmap<double> (a0, a1, octave::math::rem<double>, "rem");
679  }
680  }
681 
682  return retval;
683 }
684 
685 /*
686 %!assert (size (rem (zeros (0, 2), zeros (0, 2))), [0, 2])
687 %!assert (size (rem (rand (2, 3, 4), zeros (2, 3, 4))), [2, 3, 4])
688 %!assert (size (rem (rand (2, 3, 4), 1)), [2, 3, 4])
689 %!assert (size (rem (1, rand (2, 3, 4))), [2, 3, 4])
690 %!assert (size (rem (1, 2)), [1, 1])
691 
692 %!assert (rem ([1, 2, 3; -1, -2, -3], 2), [1, 0, 1; -1, 0, -1])
693 %!assert (rem ([1, 2, 3; -1, -2, -3], 2 * ones (2, 3)),[1, 0, 1; -1, 0, -1])
694 %!assert (rem ([0, 1, 2], [0, 0, 1]), [NaN, NaN, 0])
695 %!assert (rem (uint8 ([1, 2, 3; -1, -2, -3]), uint8 (2)), uint8 ([1, 0, 1; -1, 0, -1]))
696 %!assert (uint8 (rem ([1, 2, 3; -1, -2, -3], 2 * ones (2, 3))),uint8 ([1, 0, 1; -1, 0, -1]))
697 %!assert (rem (uint8 ([0, 1, 2]), [0, 0, 1]), uint8 ([0, 0, 0]))
698 
699 ## Test sparse implementations
700 %!shared xs
701 %! xs = sparse (0:3);
702 %!test
703 %! y = rem (11, xs);
704 %! assert (isnan (y(1)));
705 %! assert (y, sparse (rem (11, 0:3)));
706 %!test
707 %! y = rem (0, xs);
708 %! assert (nnz (y), 1);
709 %! assert (y, sparse ([NaN 0 0 0]));
710 %!test
711 %! y = rem (xs, 2);
712 %! assert (nnz (y), 2);
713 %! assert (y, sparse (rem (0:3, 2)));
714 %!test
715 %! y = rem (xs, 1);
716 %! assert (nnz (y), 0);
717 %! assert (y, sparse (rem (0:3, 1)));
718 %!test
719 %! y = rem (sparse ([11 11 11 11]), xs);
720 %! assert (nnz (y), 3);
721 %! assert (y, sparse (rem (11, 0:3)));
722 %!test
723 %! y = rem (sparse ([0 0 0 0]), xs);
724 %! assert (nnz (y), 1);
725 %! assert (y, sparse ([NaN 0 0 0]));
726 
727 %!assert <45587> (signbit (rem (-0, 1)))
728 %!assert <45587> (! signbit (rem (0, 1)))
729 
730 %!assert <42627> (rem (0.94, 0.01), 0.0)
731 
732 %!error rem (uint (8), int8 (5))
733 %!error rem (uint8 ([1, 2]), uint8 ([3, 4, 5]))
734 %!error rem ()
735 %!error rem (1, 2, 3)
736 %!error rem ([1, 2], [3, 4, 5])
737 %!error rem (i, 1)
738 */
739 
740 DEFUN (mod, args, ,
741  doc: /* -*- texinfo -*-
742 @deftypefn {} {} mod (@var{x}, @var{y})
743 Compute the modulo of @var{x} and @var{y}.
744 
745 Conceptually this is given by
746 
747 @example
748 x - y .* floor (x ./ y)
749 @end example
750 
751 @noindent
752 and is written such that the correct modulus is returned for integer types.
753 This function handles negative values correctly. That is,
754 @w{@code{mod (-1, 3)}} is 2, not -1, as @w{@code{rem (-1, 3)}} returns.
755 
756 An error results if the dimensions of the arguments do not agree, or if
757 either of the arguments is complex.
758 
759 Programming Notes: Floating point numbers within a few eps of an integer
760 will be rounded to an integer before computation for compatibility with
761 @sc{matlab}.
762 
763 By convention,
764 
765 @example
766 @group
767 mod (@var{x}, 0) = @var{x}
768 mod (@var{x}, @var{y}) returns a value with the signbit from @var{y}
769 @end group
770 @end example
771 
772 For the opposite conventions see the @code{rem} function. In general,
773 @code{mod} is a better choice than @code{rem} when any of the inputs are
774 negative numbers or when the values are periodic.
775 @seealso{rem}
776 @end deftypefn */)
777 {
778  if (args.length () != 2)
779  print_usage ();
780 
782 
783  if (! args(0).is_numeric_type ())
784  err_wrong_type_arg ("mod", args(0));
785 
786  if (! args(1).is_numeric_type ())
787  err_wrong_type_arg ("mod", args(1));
788 
789  if (args(0).is_complex_type () || args(1).is_complex_type ())
790  error ("mod: not defined for complex numbers");
791 
792  if (args(0).is_integer_type () || args(1).is_integer_type ())
793  {
794  builtin_type_t btyp0 = args(0).builtin_type ();
795  builtin_type_t btyp1 = args(1).builtin_type ();
796  if (btyp0 == btyp_double || btyp0 == btyp_float)
797  btyp0 = btyp1;
798  if (btyp1 == btyp_double || btyp1 == btyp_float)
799  btyp1 = btyp0;
800 
801  if (btyp0 != btyp1)
802  error ("mod: cannot combine %s and %d",
803  args(0).class_name ().c_str (),
804  args(1).class_name ().c_str ());
805 
806  switch (btyp0)
807  {
808 #define MAKE_INT_BRANCH(X) \
809  case btyp_ ## X: \
810  { \
811  X##NDArray a0 = args(0).X##_array_value (); \
812  X##NDArray a1 = args(1).X##_array_value (); \
813  retval = binmap<octave_##X,octave_##X,octave_##X> (a0, a1, mod, "mod"); \
814  } \
815  break
816 
817  MAKE_INT_BRANCH (int8);
821  MAKE_INT_BRANCH (uint8);
824  MAKE_INT_BRANCH (uint64);
825 
826 #undef MAKE_INT_BRANCH
827 
828  default:
829  panic_impossible ();
830  }
831  }
832  else if (args(0).is_single_type () || args(1).is_single_type ())
833  {
834  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
835  retval = octave::math::mod (args(0).float_value (), args(1).float_value ());
836  else
837  {
838  FloatNDArray a0 = args(0).float_array_value ();
839  FloatNDArray a1 = args(1).float_array_value ();
840  retval = binmap<float> (a0, a1, octave::math::mod<float>, "mod");
841  }
842  }
843  else
844  {
845  if (args(0).is_scalar_type () && args(1).is_scalar_type ())
846  retval = octave::math::mod (args(0).scalar_value (), args(1).scalar_value ());
847  else if (args(0).is_sparse_type () || args(1).is_sparse_type ())
848  {
849  SparseMatrix m0 = args(0).sparse_matrix_value ();
850  SparseMatrix m1 = args(1).sparse_matrix_value ();
851  retval = binmap<double> (m0, m1, octave::math::mod<double>, "mod");
852  }
853  else
854  {
855  NDArray a0 = args(0).array_value ();
856  NDArray a1 = args(1).array_value ();
857  retval = binmap<double> (a0, a1, octave::math::mod<double>, "mod");
858  }
859  }
860 
861  return retval;
862 }
863 
864 /*
865 ## empty input test
866 %!assert (isempty (mod ([], [])))
867 
868 ## x mod y, y != 0 tests
869 %!assert (mod (5, 3), 2)
870 %!assert (mod (-5, 3), 1)
871 %!assert (mod (0, 3), 0)
872 %!assert (mod ([-5, 5, 0], [3, 3, 3]), [1, 2, 0])
873 %!assert (mod ([-5; 5; 0], [3; 3; 3]), [1; 2; 0])
874 %!assert (mod ([-5, 5; 0, 3], [3, 3 ; 3, 1]), [1, 2 ; 0, 0])
875 
876 ## x mod 0 tests
877 %!assert (mod (5, 0), 5)
878 %!assert (mod (-5, 0), -5)
879 %!assert (mod ([-5, 5, 0], [3, 0, 3]), [1, 5, 0])
880 %!assert (mod ([-5; 5; 0], [3; 0; 3]), [1; 5; 0])
881 %!assert (mod ([-5, 5; 0, 3], [3, 0 ; 3, 1]), [1, 5 ; 0, 0])
882 %!assert (mod ([-5, 5; 0, 3], [0, 0 ; 0, 0]), [-5, 5; 0, 3])
883 
884 ## mixed scalar/matrix tests
885 %!assert (mod ([-5, 5; 0, 3], 0), [-5, 5; 0, 3])
886 %!assert (mod ([-5, 5; 0, 3], 3), [1, 2; 0, 0])
887 %!assert (mod (-5, [0,0; 0,0]), [-5, -5; -5, -5])
888 %!assert (mod (-5, [3,0; 3,1]), [1, -5; 1, 0])
889 %!assert (mod (-5, [3,2; 3,1]), [1, 1; 1, 0])
890 
891 ## integer types
892 %!assert (mod (uint8 (5), uint8 (4)), uint8 (1))
893 %!assert (mod (uint8 ([1:5]), uint8 (4)), uint8 ([1,2,3,0,1]))
894 %!assert (mod (uint8 ([1:5]), uint8 (0)), uint8 ([1:5]))
895 %!error (mod (uint8 (5), int8 (4)))
896 
897 ## mixed integer/real types
898 %!assert (mod (uint8 (5), 4), uint8 (1))
899 %!assert (mod (5, uint8 (4)), uint8 (1))
900 %!assert (mod (uint8 ([1:5]), 4), uint8 ([1,2,3,0,1]))
901 
902 ## non-integer real numbers
903 %!assert (mod (2.1, 0.1), 0)
904 %!assert (mod (2.1, 0.2), 0.1, eps)
905 
906 %!assert <45587> (signbit (mod (-0, 0)))
907 %!assert <45587> (! signbit (mod (0, -0)))
908 
909 %!assert <42627> (mod (0.94, 0.01), 0.0)
910 */
911 
912 // FIXME: Macros NATIVE_REDUCTION_1 and NATIVE_REDUCTION seem to be unused.
913 // Checked 1/23/2016. They should probably be removed for clarity.
914 // FIXME: Need to convert reduction functions of this file for single precision
915 
916 #define NATIVE_REDUCTION_1(FCN, TYPE, DIM) \
917  (arg.is_ ## TYPE ## _type ()) \
918  { \
919  TYPE ## NDArray tmp = arg. TYPE ##_array_value (); \
920  \
921  retval = tmp.FCN (DIM); \
922  }
923 
924 #define NATIVE_REDUCTION(FCN, BOOL_FCN) \
925  \
926  int nargin = args.length (); \
927  \
928  bool isnative = false; \
929  bool isdouble = false; \
930  \
931  if (nargin > 1 && args(nargin - 1).is_string ()) \
932  { \
933  std::string str = args(nargin - 1).string_value (); \
934  \
935  if (str == "native") \
936  isnative = true; \
937  else if (str == "double") \
938  isdouble = true; \
939  else \
940  error ("sum: unrecognized string argument"); \
941  \
942  nargin--; \
943  } \
944  \
945  if (nargin < 1 || nargin > 2) \
946  print_usage (); \
947  \
948  octave_value retval; \
949  \
950  octave_value arg = args(0); \
951  \
952  int dim = (nargin == 1 ? -1 : args(1).int_value (true) - 1); \
953  \
954  if (dim < -1) \
955  error (#FCN ": invalid dimension argument = %d", dim + 1); \
956  \
957  if (arg.is_sparse_type ()) \
958  { \
959  if (arg.is_real_type ()) \
960  { \
961  SparseMatrix tmp = arg.sparse_matrix_value (); \
962  \
963  retval = tmp.FCN (dim); \
964  } \
965  else \
966  { \
967  SparseComplexMatrix tmp \
968  = arg.sparse_complex_matrix_value (); \
969  \
970  retval = tmp.FCN (dim); \
971  } \
972  } \
973  else \
974  { \
975  if (isnative) \
976  { \
977  if NATIVE_REDUCTION_1 (FCN, uint8, dim) \
978  else if NATIVE_REDUCTION_1 (FCN, uint16, dim) \
979  else if NATIVE_REDUCTION_1 (FCN, uint32, dim) \
980  else if NATIVE_REDUCTION_1 (FCN, uint64, dim) \
981  else if NATIVE_REDUCTION_1 (FCN, int8, dim) \
982  else if NATIVE_REDUCTION_1 (FCN, int16, dim) \
983  else if NATIVE_REDUCTION_1 (FCN, int32, dim) \
984  else if NATIVE_REDUCTION_1 (FCN, int64, dim) \
985  else if (arg.is_bool_type ()) \
986  { \
987  boolNDArray tmp = arg.bool_array_value (); \
988  \
989  retval = boolNDArray (tmp.BOOL_FCN (dim)); \
990  } \
991  else if (arg.is_char_matrix ()) \
992  { \
993  error (#FCN, ": invalid char type"); \
994  } \
995  else if (! isdouble && arg.is_single_type ()) \
996  { \
997  if (arg.is_complex_type ()) \
998  { \
999  FloatComplexNDArray tmp = \
1000  arg.float_complex_array_value (); \
1001  \
1002  retval = tmp.FCN (dim); \
1003  } \
1004  else if (arg.is_real_type ()) \
1005  { \
1006  FloatNDArray tmp = arg.float_array_value (); \
1007  \
1008  retval = tmp.FCN (dim); \
1009  } \
1010  } \
1011  else if (arg.is_complex_type ()) \
1012  { \
1013  ComplexNDArray tmp = arg.complex_array_value (); \
1014  \
1015  retval = tmp.FCN (dim); \
1016  } \
1017  else if (arg.is_real_type ()) \
1018  { \
1019  NDArray tmp = arg.array_value (); \
1020  \
1021  retval = tmp.FCN (dim); \
1022  } \
1023  else \
1024  err_wrong_type_arg (#FCN, arg); \
1025  } \
1026  else if (arg.is_bool_type ()) \
1027  { \
1028  boolNDArray tmp = arg.bool_array_value (); \
1029  \
1030  retval = tmp.FCN (dim); \
1031  } \
1032  else if (! isdouble && arg.is_single_type ()) \
1033  { \
1034  if (arg.is_real_type ()) \
1035  { \
1036  FloatNDArray tmp = arg.float_array_value (); \
1037  \
1038  retval = tmp.FCN (dim); \
1039  } \
1040  else if (arg.is_complex_type ()) \
1041  { \
1042  FloatComplexNDArray tmp = \
1043  arg.float_complex_array_value (); \
1044  \
1045  retval = tmp.FCN (dim); \
1046  } \
1047  } \
1048  else if (arg.is_real_type ()) \
1049  { \
1050  NDArray tmp = arg.array_value (); \
1051  \
1052  retval = tmp.FCN (dim); \
1053  } \
1054  else if (arg.is_complex_type ()) \
1055  { \
1056  ComplexNDArray tmp = arg.complex_array_value (); \
1057  \
1058  retval = tmp.FCN (dim); \
1059  } \
1060  else \
1061  err_wrong_type_arg (#FCN, arg); \
1062  } \
1063  \
1064  return retval
1065 
1066 #define DATA_REDUCTION(FCN) \
1067  \
1068  int nargin = args.length (); \
1069  \
1070  if (nargin < 1 || nargin > 2) \
1071  print_usage (); \
1072  \
1073  octave_value retval; \
1074  \
1075  octave_value arg = args(0); \
1076  \
1077  int dim = (nargin == 1 ? -1 : args(1).int_value (true) - 1); \
1078  \
1079  if (dim < -1) \
1080  error (#FCN ": invalid dimension argument = %d", dim + 1); \
1081  \
1082  if (arg.is_real_type ()) \
1083  { \
1084  if (arg.is_sparse_type ()) \
1085  { \
1086  SparseMatrix tmp = arg.sparse_matrix_value (); \
1087  \
1088  retval = tmp.FCN (dim); \
1089  } \
1090  else if (arg.is_single_type ()) \
1091  { \
1092  FloatNDArray tmp = arg.float_array_value (); \
1093  \
1094  retval = tmp.FCN (dim); \
1095  } \
1096  else \
1097  { \
1098  NDArray tmp = arg.array_value (); \
1099  \
1100  retval = tmp.FCN (dim); \
1101  } \
1102  } \
1103  else if (arg.is_complex_type ()) \
1104  { \
1105  if (arg.is_sparse_type ()) \
1106  { \
1107  SparseComplexMatrix tmp = arg.sparse_complex_matrix_value (); \
1108  \
1109  retval = tmp.FCN (dim); \
1110  } \
1111  else if (arg.is_single_type ()) \
1112  { \
1113  FloatComplexNDArray tmp \
1114  = arg.float_complex_array_value (); \
1115  \
1116  retval = tmp.FCN (dim); \
1117  } \
1118  else \
1119  { \
1120  ComplexNDArray tmp = arg.complex_array_value (); \
1121  \
1122  retval = tmp.FCN (dim); \
1123  } \
1124  } \
1125  else \
1126  err_wrong_type_arg (#FCN, arg); \
1127  \
1128  return retval
1129 
1130 DEFUN (cumprod, args, ,
1131  doc: /* -*- texinfo -*-
1132 @deftypefn {} {} cumprod (@var{x})
1133 @deftypefnx {} {} cumprod (@var{x}, @var{dim})
1134 Cumulative product of elements along dimension @var{dim}.
1135 
1136 If @var{dim} is omitted, it defaults to the first non-singleton dimension.
1137 For example:
1138 
1139 @example
1140 @group
1141 cumprod ([1, 2; 3, 4; 5, 6])
1142  @result{} 1 2
1143  3 8
1144  15 48
1145 @end group
1146 @end example
1147 @seealso{prod, cumsum}
1148 @end deftypefn */)
1149 {
1150  DATA_REDUCTION (cumprod);
1151 }
1152 
1153 /*
1154 %!assert (cumprod ([1, 2, 3]), [1, 2, 6])
1155 %!assert (cumprod ([-1; -2; -3]), [-1; 2; -6])
1156 %!assert (cumprod ([i, 2+i, -3+2i, 4]), [i, -1+2i, -1-8i, -4-32i])
1157 %!assert (cumprod ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i]), [1, 2, 3; i, 4i, 9i; -1+i, -8+8i, -27+27i])
1158 
1159 %!assert (cumprod (single ([1, 2, 3])), single ([1, 2, 6]))
1160 %!assert (cumprod (single ([-1; -2; -3])), single ([-1; 2; -6]))
1161 %!assert (cumprod (single ([i, 2+i, -3+2i, 4])), single ([i, -1+2i, -1-8i, -4-32i]))
1162 %!assert (cumprod (single ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i])), single ([1, 2, 3; i, 4i, 9i; -1+i, -8+8i, -27+27i]))
1163 
1164 %!assert (cumprod ([2, 3; 4, 5], 1), [2, 3; 8, 15])
1165 %!assert (cumprod ([2, 3; 4, 5], 2), [2, 6; 4, 20])
1166 
1167 %!assert (cumprod (single ([2, 3; 4, 5]), 1), single ([2, 3; 8, 15]))
1168 %!assert (cumprod (single ([2, 3; 4, 5]), 2), single ([2, 6; 4, 20]))
1169 
1170 %!error cumprod ()
1171 */
1172 
1173 DEFUN (cumsum, args, ,
1174  doc: /* -*- texinfo -*-
1175 @deftypefn {} {} cumsum (@var{x})
1176 @deftypefnx {} {} cumsum (@var{x}, @var{dim})
1177 @deftypefnx {} {} cumsum (@dots{}, "native")
1178 @deftypefnx {} {} cumsum (@dots{}, "double")
1179 @deftypefnx {} {} cumsum (@dots{}, "extra")
1180 Cumulative sum of elements along dimension @var{dim}.
1181 
1182 If @var{dim} is omitted, it defaults to the first non-singleton dimension.
1183 For example:
1184 
1185 @example
1186 @group
1187 cumsum ([1, 2; 3, 4; 5, 6])
1188  @result{} 1 2
1189  4 6
1190  9 12
1191 @end group
1192 @end example
1193 
1194 See @code{sum} for an explanation of the optional parameters
1195 @qcode{"native"}, @qcode{"double"}, and @qcode{"extra"}.
1196 @seealso{sum, cumprod}
1197 @end deftypefn */)
1198 {
1199  int nargin = args.length ();
1200 
1201  bool isnative = false;
1202  bool isdouble = false;
1203 
1204  if (nargin > 1 && args(nargin - 1).is_string ())
1205  {
1206  std::string str = args(nargin - 1).string_value ();
1207 
1208  if (str == "native")
1209  isnative = true;
1210  else if (str == "double")
1211  isdouble = true;
1212  else
1213  error ("cumsum: unrecognized string argument");
1214 
1215  nargin--;
1216  }
1217 
1218  if (nargin < 1 || nargin > 2)
1219  print_usage ();
1220 
1221  int dim = -1;
1222  if (nargin == 2)
1223  {
1224  dim = args(1).int_value () - 1;
1225  if (dim < 0)
1226  error ("cumsum: invalid dimension argument = %d", dim + 1);
1227  }
1228 
1230  octave_value arg = args(0);
1231 
1232  switch (arg.builtin_type ())
1233  {
1234  case btyp_double:
1235  if (arg.is_sparse_type ())
1236  retval = arg.sparse_matrix_value ().cumsum (dim);
1237  else
1238  retval = arg.array_value ().cumsum (dim);
1239  break;
1240  case btyp_complex:
1241  if (arg.is_sparse_type ())
1242  retval = arg.sparse_complex_matrix_value ().cumsum (dim);
1243  else
1244  retval = arg.complex_array_value ().cumsum (dim);
1245  break;
1246  case btyp_float:
1247  if (isdouble)
1248  retval = arg.array_value ().cumsum (dim);
1249  else
1250  retval = arg.float_array_value ().cumsum (dim);
1251  break;
1252  case btyp_float_complex:
1253  if (isdouble)
1254  retval = arg.complex_array_value ().cumsum (dim);
1255  else
1256  retval = arg.float_complex_array_value ().cumsum (dim);
1257  break;
1258 
1259 #define MAKE_INT_BRANCH(X) \
1260  case btyp_ ## X: \
1261  if (isnative) \
1262  retval = arg.X ## _array_value ().cumsum (dim); \
1263  else \
1264  retval = arg.array_value ().cumsum (dim); \
1265  break;
1266 
1267  MAKE_INT_BRANCH (int8);
1271  MAKE_INT_BRANCH (uint8);
1274  MAKE_INT_BRANCH (uint64);
1275 
1276 #undef MAKE_INT_BRANCH
1277 
1278  case btyp_bool:
1279  if (arg.is_sparse_type ())
1280  {
1281  SparseMatrix cs = arg.sparse_matrix_value ().cumsum (dim);
1282  if (isnative)
1283  retval = cs != 0.0;
1284  else
1285  retval = cs;
1286  }
1287  else
1288  {
1289  NDArray cs = arg.bool_array_value ().cumsum (dim);
1290  if (isnative)
1291  retval = cs != 0.0;
1292  else
1293  retval = cs;
1294  }
1295  break;
1296 
1297  default:
1298  err_wrong_type_arg ("cumsum", arg);
1299  }
1300 
1301  return retval;
1302 }
1303 
1304 /*
1305 %!assert (cumsum ([1, 2, 3]), [1, 3, 6])
1306 %!assert (cumsum ([-1; -2; -3]), [-1; -3; -6])
1307 %!assert (cumsum ([i, 2+i, -3+2i, 4]), [i, 2+2i, -1+4i, 3+4i])
1308 %!assert (cumsum ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i]), [1, 2, 3; 1+i, 2+2i, 3+3i; 2+2i, 4+4i, 6+6i])
1309 
1310 %!assert (cumsum (single ([1, 2, 3])), single ([1, 3, 6]))
1311 %!assert (cumsum (single ([-1; -2; -3])), single ([-1; -3; -6]))
1312 %!assert (cumsum (single ([i, 2+i, -3+2i, 4])), single ([i, 2+2i, -1+4i, 3+4i]))
1313 %!assert (cumsum (single ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i])), single ([1, 2, 3; 1+i, 2+2i, 3+3i; 2+2i, 4+4i, 6+6i]))
1314 
1315 %!assert (cumsum ([1, 2; 3, 4], 1), [1, 2; 4, 6])
1316 %!assert (cumsum ([1, 2; 3, 4], 2), [1, 3; 3, 7])
1317 
1318 %!assert (cumsum (single ([1, 2; 3, 4]), 1), single ([1, 2; 4, 6]))
1319 %!assert (cumsum (single ([1, 2; 3, 4]), 2), single ([1, 3; 3, 7]))
1320 
1321 %!error cumsum ()
1322 */
1323 
1324 DEFUN (diag, args, ,
1325  doc: /* -*- texinfo -*-
1326 @deftypefn {} {@var{M} =} diag (@var{v})
1327 @deftypefnx {} {@var{M} =} diag (@var{v}, @var{k})
1328 @deftypefnx {} {@var{M} =} diag (@var{v}, @var{m}, @var{n})
1329 @deftypefnx {} {@var{v} =} diag (@var{M})
1330 @deftypefnx {} {@var{v} =} diag (@var{M}, @var{k})
1331 Return a diagonal matrix with vector @var{v} on diagonal @var{k}.
1332 
1333 The second argument is optional. If it is positive, the vector is placed on
1334 the @var{k}-th superdiagonal. If it is negative, it is placed on the
1335 @var{-k}-th subdiagonal. The default value of @var{k} is 0, and the vector
1336 is placed on the main diagonal. For example:
1337 
1338 @example
1339 @group
1340 diag ([1, 2, 3], 1)
1341  @result{} 0 1 0 0
1342  0 0 2 0
1343  0 0 0 3
1344  0 0 0 0
1345 @end group
1346 @end example
1347 
1348 @noindent
1349 The 3-input form returns a diagonal matrix with vector @var{v} on the main
1350 diagonal and the resulting matrix being of size @var{m} rows x @var{n}
1351 columns.
1352 
1353 Given a matrix argument, instead of a vector, @code{diag} extracts the
1354 @var{k}-th diagonal of the matrix.
1355 @end deftypefn */)
1356 {
1357  int nargin = args.length ();
1358 
1359  if (nargin < 1 || nargin > 3)
1360  print_usage ();
1361 
1363 
1364  if (nargin == 1)
1365  retval = args(0).diag ();
1366  else if (nargin == 2)
1367  {
1368  octave_idx_type k = args(1).xidx_type_value ("diag: invalid argument K");
1369 
1370  retval = args(0).diag (k);
1371  }
1372  else
1373  {
1374  octave_value arg0 = args(0);
1375 
1376  if (arg0.ndims () != 2 || (arg0.rows () != 1 && arg0.columns () != 1))
1377  error ("diag: V must be a vector");
1378 
1379  octave_idx_type m = args(1).xidx_type_value ("diag: invalid dimensions");
1380  octave_idx_type n = args(2).xidx_type_value ("diag: invalid dimensions");
1381 
1382  retval = arg0.diag (m, n);
1383  }
1384 
1385  return retval;
1386 }
1387 
1388 /*
1389 
1390 %!assert (full (diag ([1; 2; 3])), [1, 0, 0; 0, 2, 0; 0, 0, 3])
1391 %!assert (diag ([1; 2; 3], 1), [0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0])
1392 %!assert (diag ([1; 2; 3], 2), [0, 0, 1, 0, 0; 0, 0, 0, 2, 0; 0, 0, 0, 0, 3; 0, 0, 0, 0, 0; 0, 0, 0, 0, 0])
1393 %!assert (diag ([1; 2; 3],-1), [0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0])
1394 %!assert (diag ([1; 2; 3],-2), [0, 0, 0, 0, 0; 0, 0, 0, 0, 0; 1, 0, 0, 0, 0; 0, 2, 0, 0, 0; 0, 0, 3, 0, 0])
1395 
1396 %!assert (diag ([1, 0, 0; 0, 2, 0; 0, 0, 3]), [1; 2; 3])
1397 %!assert (diag ([0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0], 1), [1; 2; 3])
1398 %!assert (diag ([0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0], -1), [1; 2; 3])
1399 %!assert (diag (ones (1, 0), 2), zeros (2))
1400 %!assert (diag (1:3, 4, 2), [1, 0; 0, 2; 0, 0; 0, 0])
1401 
1402 %!assert (full (diag (single ([1; 2; 3]))), single ([1, 0, 0; 0, 2, 0; 0, 0, 3]))
1403 %!assert (diag (single ([1; 2; 3]), 1), single ([0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0]))
1404 %!assert (diag (single ([1; 2; 3]), 2), single ([0, 0, 1, 0, 0; 0, 0, 0, 2, 0; 0, 0, 0, 0, 3; 0, 0, 0, 0, 0; 0, 0, 0, 0, 0]))
1405 %!assert (diag (single ([1; 2; 3]),-1), single ([0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0]))
1406 %!assert (diag (single ([1; 2; 3]),-2), single ([0, 0, 0, 0, 0; 0, 0, 0, 0, 0; 1, 0, 0, 0, 0; 0, 2, 0, 0, 0; 0, 0, 3, 0, 0]))
1407 
1408 %!assert (diag (single ([1, 0, 0; 0, 2, 0; 0, 0, 3])), single ([1; 2; 3]))
1409 %!assert (diag (single ([0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0]), 1), single ([1; 2; 3]))
1410 %!assert (diag (single ([0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0]), -1), single ([1; 2; 3]))
1411 
1412 %!assert (diag (int8 ([1; 2; 3])), int8 ([1, 0, 0; 0, 2, 0; 0, 0, 3]))
1413 %!assert (diag (int8 ([1; 2; 3]), 1), int8 ([0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0]))
1414 %!assert (diag (int8 ([1; 2; 3]), 2), int8 ([0, 0, 1, 0, 0; 0, 0, 0, 2, 0; 0, 0, 0, 0, 3; 0, 0, 0, 0, 0; 0, 0, 0, 0, 0]))
1415 %!assert (diag (int8 ([1; 2; 3]),-1), int8 ([0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0]))
1416 %!assert (diag (int8 ([1; 2; 3]),-2), int8 ([0, 0, 0, 0, 0; 0, 0, 0, 0, 0; 1, 0, 0, 0, 0; 0, 2, 0, 0, 0; 0, 0, 3, 0, 0]))
1417 
1418 %!assert (diag (int8 ([1, 0, 0; 0, 2, 0; 0, 0, 3])), int8 ([1; 2; 3]))
1419 %!assert (diag (int8 ([0, 1, 0, 0; 0, 0, 2, 0; 0, 0, 0, 3; 0, 0, 0, 0]), 1), int8 ([1; 2; 3]))
1420 %!assert (diag (int8 ([0, 0, 0, 0; 1, 0, 0, 0; 0, 2, 0, 0; 0, 0, 3, 0]), -1), int8 ([1; 2; 3]))
1421 
1422 %!assert <37411> (diag (diag ([5, 2, 3])(:,1)), diag([5 0 0 ]))
1423 %!assert <37411> (diag (diag ([5, 2, 3])(:,1), 2), [0 0 5 0 0; zeros(4, 5)])
1424 %!assert <37411> (diag (diag ([5, 2, 3])(:,1), -2), [[0 0 5 0 0]', zeros(5, 4)])
1425 
1426 ## Test non-square size
1427 %!assert (diag ([1,2,3], 6, 3), [1 0 0; 0 2 0; 0 0 3; 0 0 0; 0 0 0; 0 0 0])
1428 %!assert (diag (1, 2, 3), [1,0,0; 0,0,0])
1429 %!assert (diag ({1}, 2, 3), {1,[],[]; [],[],[]})
1430 %!assert (diag ({1,2}, 3, 4), {1,[],[],[]; [],2,[],[]; [],[],[],[]})
1431 
1432 ## Test out-of-range diagonals
1433 %!assert (diag (ones (3,3), 4), zeros (0, 1))
1434 %!assert (diag (cell (3,3), 4), cell (0, 1))
1435 %!assert (diag (sparse (ones (3,3)), 4), sparse (zeros (0, 1)))
1436 
1437 %% Test input validation
1438 %!error <Invalid call to diag> diag ()
1439 %!error <Invalid call to diag> diag (1,2,3,4)
1440 %!error diag (ones (2), 3, 3)
1441 %!error diag (1:3, -4, 3)
1442 
1443 %!assert (diag (1, 3, 3), diag ([1, 0, 0]))
1444 %!assert (diag (i, 3, 3), diag ([i, 0, 0]))
1445 %!assert (diag (single (1), 3, 3), diag ([single(1), 0, 0]))
1446 %!assert (diag (single (i), 3, 3), diag ([single(i), 0, 0]))
1447 %!assert (diag ([1, 2], 3, 3), diag ([1, 2, 0]))
1448 %!assert (diag ([1, 2]*i, 3, 3), diag ([1, 2, 0]*i))
1449 %!assert (diag (single ([1, 2]), 3, 3), diag (single ([1, 2, 0])))
1450 %!assert (diag (single ([1, 2]*i), 3, 3), diag (single ([1, 2, 0]*i)))
1451 */
1452 
1453 DEFUN (prod, args, ,
1454  doc: /* -*- texinfo -*-
1455 @deftypefn {} {} prod (@var{x})
1456 @deftypefnx {} {} prod (@var{x}, @var{dim})
1457 @deftypefnx {} {} prod (@dots{}, "native")
1458 @deftypefnx {} {} prod (@dots{}, "double")
1459 Product of elements along dimension @var{dim}.
1460 
1461 If @var{dim} is omitted, it defaults to the first non-singleton dimension.
1462 
1463 The optional @qcode{"type"} input determines the class of the variable
1464 used for calculations. If the argument @qcode{"native"} is given, then
1465 the operation is performed in the same type as the original argument, rather
1466 than the default double type.
1467 
1468 For example:
1469 
1470 @example
1471 @group
1472 prod ([true, true])
1473  @result{} 1
1474 prod ([true, true], "native")
1475  @result{} true
1476 @end group
1477 @end example
1478 
1479 On the contrary, if @qcode{"double"} is given, the operation is performed
1480 in double precision even for single precision inputs.
1481 @seealso{cumprod, sum}
1482 @end deftypefn */)
1483 {
1484  int nargin = args.length ();
1485 
1486  bool isnative = false;
1487  bool isdouble = false;
1488 
1489  if (nargin > 1 && args(nargin - 1).is_string ())
1490  {
1491  std::string str = args(nargin - 1).string_value ();
1492 
1493  if (str == "native")
1494  isnative = true;
1495  else if (str == "double")
1496  isdouble = true;
1497  else
1498  error ("prod: unrecognized type argument '%s'", str.c_str ());
1499 
1500  nargin--;
1501  }
1502 
1503  if (nargin < 1 || nargin > 2)
1504  print_usage ();
1505 
1507 
1508  octave_value arg = args(0);
1509 
1510  int dim = -1;
1511  if (nargin == 2)
1512  {
1513  dim = args(1).int_value () - 1;
1514  if (dim < 0)
1515  error ("prod: invalid dimension DIM = %d", dim + 1);
1516  }
1517 
1518  switch (arg.builtin_type ())
1519  {
1520  case btyp_double:
1521  if (arg.is_sparse_type ())
1522  retval = arg.sparse_matrix_value ().prod (dim);
1523  else
1524  retval = arg.array_value ().prod (dim);
1525  break;
1526  case btyp_complex:
1527  if (arg.is_sparse_type ())
1528  retval = arg.sparse_complex_matrix_value ().prod (dim);
1529  else
1530  retval = arg.complex_array_value ().prod (dim);
1531  break;
1532  case btyp_float:
1533  if (isdouble)
1534  retval = arg.float_array_value ().dprod (dim);
1535  else
1536  retval = arg.float_array_value ().prod (dim);
1537  break;
1538  case btyp_float_complex:
1539  if (isdouble)
1540  retval = arg.float_complex_array_value ().dprod (dim);
1541  else
1542  retval = arg.float_complex_array_value ().prod (dim);
1543  break;
1544 
1545 #define MAKE_INT_BRANCH(X) \
1546  case btyp_ ## X: \
1547  if (isnative) \
1548  retval = arg.X ## _array_value ().prod (dim); \
1549  else \
1550  retval = arg.array_value ().prod (dim); \
1551  break;
1552 
1553  MAKE_INT_BRANCH (int8);
1557  MAKE_INT_BRANCH (uint8);
1560  MAKE_INT_BRANCH (uint64);
1561 
1562 #undef MAKE_INT_BRANCH
1563 
1564  // GAGME: Accursed Matlab compatibility...
1565  case btyp_char:
1566  retval = arg.array_value (true).prod (dim);
1567  break;
1568 
1569  case btyp_bool:
1570  if (arg.is_sparse_type ())
1571  {
1572  if (isnative)
1573  retval = arg.sparse_bool_matrix_value ().all (dim);
1574  else
1575  retval = arg.sparse_matrix_value ().prod (dim);
1576  }
1577  else if (isnative)
1578  retval = arg.bool_array_value ().all (dim);
1579  else
1580  retval = NDArray (arg.bool_array_value ().all (dim));
1581  break;
1582 
1583  default:
1584  err_wrong_type_arg ("prod", arg);
1585  }
1586 
1587  return retval;
1588 }
1589 
1590 /*
1591 %!assert (prod ([1, 2, 3]), 6)
1592 %!assert (prod ([-1; -2; -3]), -6)
1593 %!assert (prod ([i, 2+i, -3+2i, 4]), -4 - 32i)
1594 %!assert (prod ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i]), [-1+i, -8+8i, -27+27i])
1595 
1596 %!assert (prod (single ([1, 2, 3])), single (6))
1597 %!assert (prod (single ([-1; -2; -3])), single (-6))
1598 %!assert (prod (single ([i, 2+i, -3+2i, 4])), single (-4 - 32i))
1599 %!assert (prod (single ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i])), single ([-1+i, -8+8i, -27+27i]))
1600 
1601 %% Test sparse
1602 %!assert (prod (sparse ([1, 2, 3])), sparse (6))
1603 %!assert (prod (sparse ([-1; -2; -3])), sparse (-6))
1604 ## Commented out until bug #42290 is fixed
1605 #%!assert (prod (sparse ([i, 2+i, -3+2i, 4])), sparse (-4 - 32i))
1606 #%!assert (prod (sparse ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i])), sparse ([-1+i, -8+8i, -27+27i]))
1607 
1608 %!assert (prod ([1, 2; 3, 4], 1), [3, 8])
1609 %!assert (prod ([1, 2; 3, 4], 2), [2; 12])
1610 %!assert (prod (zeros (1, 0)), 1)
1611 %!assert (prod (zeros (1, 0), 1), zeros (1, 0))
1612 %!assert (prod (zeros (1, 0), 2), 1)
1613 %!assert (prod (zeros (0, 1)), 1)
1614 %!assert (prod (zeros (0, 1), 1), 1)
1615 %!assert (prod (zeros (0, 1), 2), zeros (0, 1))
1616 %!assert (prod (zeros (2, 0)), zeros (1, 0))
1617 %!assert (prod (zeros (2, 0), 1), zeros (1, 0))
1618 %!assert (prod (zeros (2, 0), 2), [1; 1])
1619 %!assert (prod (zeros (0, 2)), [1, 1])
1620 %!assert (prod (zeros (0, 2), 1), [1, 1])
1621 %!assert (prod (zeros (0, 2), 2), zeros (0, 1))
1622 
1623 %!assert (prod (single ([1, 2; 3, 4]), 1), single ([3, 8]))
1624 %!assert (prod (single ([1, 2; 3, 4]), 2), single ([2; 12]))
1625 %!assert (prod (zeros (1, 0, "single")), single (1))
1626 %!assert (prod (zeros (1, 0, "single"), 1), zeros (1, 0, "single"))
1627 %!assert (prod (zeros (1, 0, "single"), 2), single (1))
1628 %!assert (prod (zeros (0, 1, "single")), single (1))
1629 %!assert (prod (zeros (0, 1, "single"), 1), single (1))
1630 %!assert (prod (zeros (0, 1, "single"), 2), zeros (0, 1, "single"))
1631 %!assert (prod (zeros (2, 0, "single")), zeros (1, 0, "single"))
1632 %!assert (prod (zeros (2, 0, "single"), 1), zeros (1, 0, "single"))
1633 %!assert (prod (zeros (2, 0, "single"), 2), single ([1; 1]))
1634 %!assert (prod (zeros (0, 2, "single")), single ([1, 1]))
1635 %!assert (prod (zeros (0, 2, "single"), 1), single ([1, 1]))
1636 %!assert (prod (zeros (0, 2, "single"), 2), zeros (0, 1, "single"))
1637 
1638 %% Test "double" type argument
1639 %!assert (prod (single ([1, 2, 3]), "double"), 6)
1640 %!assert (prod (single ([-1; -2; -3]), "double"), -6)
1641 %!assert (prod (single ([i, 2+i, -3+2i, 4]), "double"), -4 - 32i)
1642 %!assert (prod (single ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i]), "double"), [-1+i, -8+8i, -27+27i])
1643 
1644 %% Test "native" type argument
1645 %!assert (prod (uint8 ([1, 2, 3]), "native"), uint8 (6))
1646 %!assert (prod (uint8 ([-1; -2; -3]), "native"), uint8 (0))
1647 %!assert (prod (int8 ([1, 2, 3]), "native"), int8 (6))
1648 %!assert (prod (int8 ([-1; -2; -3]), "native"), int8 (-6))
1649 %!assert (prod ([true false; true true], "native"), [true false])
1650 %!assert (prod ([true false; true true], 2, "native"), [false; true])
1651 
1652 %% Test input validation
1653 %!error prod ()
1654 %!error prod (1,2,3)
1655 %!error <unrecognized type argument 'foobar'> prod (1, "foobar")
1656 */
1657 
1658 static bool
1659 all_scalar_1x1 (const octave_value_list& args)
1660 {
1661  int n_args = args.length ();
1662  for (int i = 0; i < n_args; i++)
1663  if (args(i).numel () != 1)
1664  return false;
1665 
1666  return true;
1667 }
1668 
1669 template <typename TYPE, typename T>
1670 static void
1672  const octave_value_list& args,
1673  int dim)
1674 {
1675  int n_args = args.length ();
1678  && all_scalar_1x1 (args))
1679  {
1680  // Optimize all scalars case.
1681  dim_vector dv (1, 1);
1682  if (dim == -1 || dim == -2)
1683  dim = -dim - 1;
1684  else if (dim >= 2)
1685  dv.resize (dim+1, 1);
1686  dv(dim) = n_args;
1687 
1688  result.clear (dv);
1689 
1690  for (int j = 0; j < n_args; j++)
1691  {
1692  octave_quit ();
1693 
1694  result(j) = octave_value_extract<T> (args(j));
1695  }
1696  }
1697  else
1698  {
1699  OCTAVE_LOCAL_BUFFER (Array<T>, array_list, n_args);
1700 
1701  for (int j = 0; j < n_args; j++)
1702  {
1703  octave_quit ();
1704 
1705  array_list[j] = octave_value_extract<TYPE> (args(j));
1706  }
1707 
1708  result = Array<T>::cat (dim, n_args, array_list);
1709  }
1710 }
1711 
1712 template <typename TYPE, typename T>
1713 static void
1714 single_type_concat (Sparse<T>& result,
1715  const octave_value_list& args,
1716  int dim)
1717 {
1718  int n_args = args.length ();
1719  OCTAVE_LOCAL_BUFFER (Sparse<T>, sparse_list, n_args);
1720 
1721  for (int j = 0; j < n_args; j++)
1722  {
1723  octave_quit ();
1724 
1725  sparse_list[j] = octave_value_extract<TYPE> (args(j));
1726  }
1727 
1728  result = Sparse<T>::cat (dim, n_args, sparse_list);
1729 }
1730 
1731 // Dispatcher.
1732 template <typename TYPE>
1733 static TYPE
1734 do_single_type_concat (const octave_value_list& args, int dim)
1735 {
1736  TYPE result;
1737 
1738  single_type_concat<TYPE, typename TYPE::element_type> (result, args, dim);
1739 
1740  return result;
1741 }
1742 
1743 template <typename MAP>
1744 static void
1745 single_type_concat_map (octave_map& result,
1746  const octave_value_list& args,
1747  int dim)
1748 {
1749  int n_args = args.length ();
1750  OCTAVE_LOCAL_BUFFER (MAP, map_list, n_args);
1751 
1752  for (int j = 0; j < n_args; j++)
1753  {
1754  octave_quit ();
1755 
1756  map_list[j] = octave_value_extract<MAP> (args(j));
1757  }
1758 
1759  result = octave_map::cat (dim, n_args, map_list);
1760 }
1761 
1762 static octave_map
1763 do_single_type_concat_map (const octave_value_list& args,
1764  int dim)
1765 {
1767  if (all_scalar_1x1 (args)) // optimize all scalars case.
1768  single_type_concat_map<octave_scalar_map> (result, args, dim);
1769  else
1770  single_type_concat_map<octave_map> (result, args, dim);
1771 
1772  return result;
1773 }
1774 
1775 static octave_value
1776 attempt_type_conversion (const octave_value& ov, std::string dtype)
1777 {
1779 
1780  // First try to find function in the class of OV that can convert to
1781  // the dispatch type dtype. It will have the name of the dispatch
1782  // type.
1783 
1784  std::string cname = ov.class_name ();
1785 
1786  octave_value fcn = symbol_table::find_method (dtype, cname);
1787 
1788  if (fcn.is_defined ())
1789  {
1791 
1792  try
1793  {
1794  result = fcn.do_multi_index_op (1, octave_value_list (1, ov));
1795  }
1796  catch (octave::execution_exception& e)
1797  {
1798  error (e, "conversion from %s to %s failed", dtype.c_str (),
1799  cname.c_str ());
1800  }
1801 
1802  if (result.empty ())
1803  error ("conversion from %s to %s failed", dtype.c_str (),
1804  cname.c_str ());
1805 
1806  retval = result(0);
1807  }
1808  else
1809  {
1810  // No conversion function available. Try the constructor for the
1811  // dispatch type.
1812 
1813  fcn = symbol_table::find_method (dtype, dtype);
1814 
1815  if (! fcn.is_defined ())
1816  error ("no constructor for %s!", dtype.c_str ());
1817 
1819 
1820  try
1821  {
1822  result = fcn.do_multi_index_op (1, octave_value_list (1, ov));
1823  }
1824  catch (octave::execution_exception& e)
1825  {
1826  error (e, "%s constructor failed for %s argument", dtype.c_str (),
1827  cname.c_str ());
1828  }
1829 
1830  if (result.empty ())
1831  error ("%s constructor failed for %s argument", dtype.c_str (),
1832  cname.c_str ());
1833 
1834  retval = result(0);
1835  }
1836 
1837  return retval;
1838 }
1839 
1841 do_class_concat (const octave_value_list& ovl, std::string cattype, int dim)
1842 {
1844 
1845  // Get dominant type for list
1846 
1847  std::string dtype = get_dispatch_type (ovl);
1848 
1849  octave_value fcn = symbol_table::find_method (cattype, dtype);
1850 
1851  if (fcn.is_defined ())
1852  {
1853  // Have method for dominant type. Call it and let it handle conversions.
1854 
1855  octave_value_list tmp2;
1856 
1857  try
1858  {
1859  tmp2 = fcn.do_multi_index_op (1, ovl);
1860  }
1861  catch (octave::execution_exception& e)
1862  {
1863  error (e, "%s/%s method failed", dtype.c_str (), cattype.c_str ());
1864  }
1865 
1866  if (tmp2.empty ())
1867  error ("%s/%s method did not return a value", dtype.c_str (),
1868  cattype.c_str ());
1869 
1870  retval = tmp2(0);
1871  }
1872  else
1873  {
1874  // No method for dominant type, so attempt type conversions for
1875  // all elements that are not of the dominant type, then do the
1876  // default operation for octave_class values.
1877 
1878  octave_idx_type j = 0;
1879  octave_idx_type len = ovl.length ();
1881  for (octave_idx_type k = 0; k < len; k++)
1882  {
1883  octave_value elt = ovl(k);
1884 
1885  std::string t1_type = elt.class_name ();
1886 
1887  if (t1_type == dtype)
1888  tmp(j++) = elt;
1889  else if (elt.is_object () || ! elt.is_empty ())
1890  tmp(j++) = attempt_type_conversion (elt, dtype);
1891  }
1892 
1893  tmp.resize (j);
1894 
1895  octave_map m = do_single_type_concat_map (tmp, dim);
1896 
1897  std::string cname = tmp(0).class_name ();
1898  std::list<std::string> parents = tmp(0).parent_class_name_list ();
1899 
1900  retval = octave_value (new octave_class (m, cname, parents));
1901  }
1902 
1903  return retval;
1904 }
1905 
1906 static octave_value
1907 do_cat (const octave_value_list& xargs, int dim, std::string fname)
1908 {
1910 
1911  // We may need to convert elements of the list to cells, so make a copy.
1912  // This should be efficient, it is done mostly by incrementing reference
1913  // counts.
1914  octave_value_list args = xargs;
1915 
1916  int n_args = args.length ();
1917 
1918  if (n_args == 0)
1919  retval = Matrix ();
1920  else if (n_args == 1)
1921  retval = args(0);
1922  else if (n_args > 1)
1923  {
1924  std::string result_type;
1925 
1926  bool all_strings_p = true;
1927  bool all_sq_strings_p = true;
1928  bool all_dq_strings_p = true;
1929  bool all_real_p = true;
1930  bool all_cmplx_p = true;
1931  bool any_sparse_p = false;
1932  bool any_cell_p = false;
1933  bool any_class_p = false;
1934 
1935  bool first_elem_is_struct = false;
1936 
1937  for (int i = 0; i < n_args; i++)
1938  {
1939  if (i == 0)
1940  {
1941  result_type = args(i).class_name ();
1942 
1943  first_elem_is_struct = args(i).is_map ();
1944  }
1945  else
1946  result_type = get_concat_class (result_type, args(i).class_name ());
1947 
1948  if (all_strings_p && ! args(i).is_string ())
1949  all_strings_p = false;
1950  if (all_sq_strings_p && ! args(i).is_sq_string ())
1951  all_sq_strings_p = false;
1952  if (all_dq_strings_p && ! args(i).is_dq_string ())
1953  all_dq_strings_p = false;
1954  if (all_real_p && ! args(i).is_real_type ())
1955  all_real_p = false;
1956  if (all_cmplx_p && ! (args(i).is_complex_type ()
1957  || args(i).is_real_type ()))
1958  all_cmplx_p = false;
1959  if (! any_sparse_p && args(i).is_sparse_type ())
1960  any_sparse_p = true;
1961  if (! any_cell_p && args(i).is_cell ())
1962  any_cell_p = true;
1963  if (! any_class_p && args(i).is_object ())
1964  any_class_p = true;
1965  }
1966 
1967  if (any_cell_p && ! any_class_p && ! first_elem_is_struct)
1968  {
1969  for (int i = 0; i < n_args; i++)
1970  {
1971  if (! args(i).is_cell ())
1972  args(i) = Cell (args(i));
1973  }
1974  }
1975 
1976  if (any_class_p)
1977  {
1978  retval = do_class_concat (args, fname, dim);
1979  }
1980  else if (result_type == "double")
1981  {
1982  if (any_sparse_p)
1983  {
1984  if (all_real_p)
1985  retval = do_single_type_concat<SparseMatrix> (args, dim);
1986  else
1987  retval = do_single_type_concat<SparseComplexMatrix> (args, dim);
1988  }
1989  else
1990  {
1991  if (all_real_p)
1992  retval = do_single_type_concat<NDArray> (args, dim);
1993  else
1994  retval = do_single_type_concat<ComplexNDArray> (args, dim);
1995  }
1996  }
1997  else if (result_type == "single")
1998  {
1999  if (all_real_p)
2000  retval = do_single_type_concat<FloatNDArray> (args, dim);
2001  else
2002  retval = do_single_type_concat<FloatComplexNDArray> (args, dim);
2003  }
2004  else if (result_type == "char")
2005  {
2006  char type = all_dq_strings_p ? '"' : '\'';
2007 
2008  if (! all_strings_p)
2009  warn_implicit_conversion ("Octave:num-to-str",
2010  "numeric", result_type);
2011  else
2012  maybe_warn_string_concat (all_dq_strings_p, all_sq_strings_p);
2013 
2014  charNDArray result = do_single_type_concat<charNDArray> (args, dim);
2015 
2016  retval = octave_value (result, type);
2017  }
2018  else if (result_type == "logical")
2019  {
2020  if (any_sparse_p)
2021  retval = do_single_type_concat<SparseBoolMatrix> (args, dim);
2022  else
2023  retval = do_single_type_concat<boolNDArray> (args, dim);
2024  }
2025  else if (result_type == "int8")
2026  retval = do_single_type_concat<int8NDArray> (args, dim);
2027  else if (result_type == "int16")
2028  retval = do_single_type_concat<int16NDArray> (args, dim);
2029  else if (result_type == "int32")
2030  retval = do_single_type_concat<int32NDArray> (args, dim);
2031  else if (result_type == "int64")
2032  retval = do_single_type_concat<int64NDArray> (args, dim);
2033  else if (result_type == "uint8")
2034  retval = do_single_type_concat<uint8NDArray> (args, dim);
2035  else if (result_type == "uint16")
2036  retval = do_single_type_concat<uint16NDArray> (args, dim);
2037  else if (result_type == "uint32")
2038  retval = do_single_type_concat<uint32NDArray> (args, dim);
2039  else if (result_type == "uint64")
2040  retval = do_single_type_concat<uint64NDArray> (args, dim);
2041  else if (result_type == "cell")
2042  retval = do_single_type_concat<Cell> (args, dim);
2043  else if (result_type == "struct")
2044  retval = do_single_type_concat_map (args, dim);
2045  else
2046  {
2047  dim_vector dv = args(0).dims ();
2048 
2049  // Default concatenation.
2050  bool (dim_vector::*concat_rule) (const dim_vector&, int)
2051  = &dim_vector::concat;
2052 
2053  if (dim == -1 || dim == -2)
2054  {
2055  concat_rule = &dim_vector::hvcat;
2056  dim = -dim - 1;
2057  }
2058 
2059  for (int i = 1; i < args.length (); i++)
2060  {
2061  if (! (dv.*concat_rule) (args(i).dims (), dim))
2062  error ("cat: dimension mismatch");
2063  }
2064 
2065  // The lines below might seem crazy, since we take a copy
2066  // of the first argument, resize it to be empty and then resize
2067  // it to be full. This is done since it means that there is no
2068  // recopying of data, as would happen if we used a single resize.
2069  // It should be noted that resize operation is also significantly
2070  // slower than the do_cat_op function, so it makes sense to have
2071  // an empty matrix and copy all data.
2072  //
2073  // We might also start with a empty octave_value using
2074  // tmp = octave_value_typeinfo::lookup_type
2075  // (args(1).type_name());
2076  // and then directly resize. However, for some types there might
2077  // be some additional setup needed, and so this should be avoided.
2078 
2079  octave_value tmp = args(0);
2080  tmp = tmp.resize (dim_vector (0,0)).resize (dv);
2081 
2082  int dv_len = dv.ndims ();
2083  Array<octave_idx_type> ra_idx (dim_vector (dv_len, 1), 0);
2084 
2085  for (int j = 0; j < n_args; j++)
2086  {
2087  // Can't fast return here to skip empty matrices as something
2088  // like cat (1,[],single ([])) must return an empty matrix of
2089  // the right type.
2090  tmp = do_cat_op (tmp, args(j), ra_idx);
2091 
2092  dim_vector dv_tmp = args(j).dims ();
2093 
2094  if (dim >= dv_len)
2095  {
2096  if (j > 1)
2097  error ("%s: indexing error", fname.c_str ());
2098 
2099  break;
2100  }
2101  else
2102  ra_idx(dim) += (dim < dv_tmp.ndims () ? dv_tmp(dim) : 1);
2103  }
2104  retval = tmp;
2105  }
2106  }
2107  else
2108  print_usage ();
2109 
2110  return retval;
2111 }
2112 
2113 DEFUN (horzcat, args, ,
2114  doc: /* -*- texinfo -*-
2115 @deftypefn {} {} horzcat (@var{array1}, @var{array2}, @dots{}, @var{arrayN})
2116 Return the horizontal concatenation of N-D array objects, @var{array1},
2117 @var{array2}, @dots{}, @var{arrayN} along dimension 2.
2118 
2119 Arrays may also be concatenated horizontally using the syntax for creating
2120 new matrices. For example:
2121 
2122 @example
2123 @var{hcat} = [ @var{array1}, @var{array2}, @dots{} ]
2124 @end example
2125 @seealso{cat, vertcat}
2126 @end deftypefn */)
2127 {
2128  return do_cat (args, -2, "horzcat");
2129 }
2130 
2131 /*
2132 ## Test concatenation with all zero matrices
2133 %!test
2134 %! warning ("off", "Octave:num-to-str", "local");
2135 %! assert (horzcat ("", 65*ones (1,10)), "AAAAAAAAAA");
2136 %! assert (horzcat (65*ones (1,10), ""), "AAAAAAAAAA");
2137 
2138 %!assert (class (horzcat (int64 (1), int64 (1))), "int64")
2139 %!assert (class (horzcat (int64 (1), int32 (1))), "int64")
2140 %!assert (class (horzcat (int64 (1), int16 (1))), "int64")
2141 %!assert (class (horzcat (int64 (1), int8 (1))), "int64")
2142 %!assert (class (horzcat (int64 (1), uint64 (1))), "int64")
2143 %!assert (class (horzcat (int64 (1), uint32 (1))), "int64")
2144 %!assert (class (horzcat (int64 (1), uint16 (1))), "int64")
2145 %!assert (class (horzcat (int64 (1), uint8 (1))), "int64")
2146 %!assert (class (horzcat (int64 (1), single (1))), "int64")
2147 %!assert (class (horzcat (int64 (1), double (1))), "int64")
2148 %!assert (class (horzcat (int64 (1), cell (1))), "cell")
2149 %!assert (class (horzcat (int64 (1), true)), "int64")
2150 %!test
2151 %! warning ("off", "Octave:num-to-str", "local");
2152 %! assert (class (horzcat (int64 (1), "a")), "char");
2153 
2154 %!assert (class (horzcat (int32 (1), int64 (1))), "int32")
2155 %!assert (class (horzcat (int32 (1), int32 (1))), "int32")
2156 %!assert (class (horzcat (int32 (1), int16 (1))), "int32")
2157 %!assert (class (horzcat (int32 (1), int8 (1))), "int32")
2158 %!assert (class (horzcat (int32 (1), uint64 (1))), "int32")
2159 %!assert (class (horzcat (int32 (1), uint32 (1))), "int32")
2160 %!assert (class (horzcat (int32 (1), uint16 (1))), "int32")
2161 %!assert (class (horzcat (int32 (1), uint8 (1))), "int32")
2162 %!assert (class (horzcat (int32 (1), single (1))), "int32")
2163 %!assert (class (horzcat (int32 (1), double (1))), "int32")
2164 %!assert (class (horzcat (int32 (1), cell (1))), "cell")
2165 %!assert (class (horzcat (int32 (1), true)), "int32")
2166 %!test
2167 %! warning ("off", "Octave:num-to-str", "local");
2168 %! assert (class (horzcat (int32 (1), "a")), "char");
2169 
2170 %!assert (class (horzcat (int16 (1), int64 (1))), "int16")
2171 %!assert (class (horzcat (int16 (1), int32 (1))), "int16")
2172 %!assert (class (horzcat (int16 (1), int16 (1))), "int16")
2173 %!assert (class (horzcat (int16 (1), int8 (1))), "int16")
2174 %!assert (class (horzcat (int16 (1), uint64 (1))), "int16")
2175 %!assert (class (horzcat (int16 (1), uint32 (1))), "int16")
2176 %!assert (class (horzcat (int16 (1), uint16 (1))), "int16")
2177 %!assert (class (horzcat (int16 (1), uint8 (1))), "int16")
2178 %!assert (class (horzcat (int16 (1), single (1))), "int16")
2179 %!assert (class (horzcat (int16 (1), double (1))), "int16")
2180 %!assert (class (horzcat (int16 (1), cell (1))), "cell")
2181 %!assert (class (horzcat (int16 (1), true)), "int16")
2182 %!test
2183 %! warning ("off", "Octave:num-to-str", "local");
2184 %! assert (class (horzcat (int16 (1), "a")), "char");
2185 
2186 %!assert (class (horzcat (int8 (1), int64 (1))), "int8")
2187 %!assert (class (horzcat (int8 (1), int32 (1))), "int8")
2188 %!assert (class (horzcat (int8 (1), int16 (1))), "int8")
2189 %!assert (class (horzcat (int8 (1), int8 (1))), "int8")
2190 %!assert (class (horzcat (int8 (1), uint64 (1))), "int8")
2191 %!assert (class (horzcat (int8 (1), uint32 (1))), "int8")
2192 %!assert (class (horzcat (int8 (1), uint16 (1))), "int8")
2193 %!assert (class (horzcat (int8 (1), uint8 (1))), "int8")
2194 %!assert (class (horzcat (int8 (1), single (1))), "int8")
2195 %!assert (class (horzcat (int8 (1), double (1))), "int8")
2196 %!assert (class (horzcat (int8 (1), cell (1))), "cell")
2197 %!assert (class (horzcat (int8 (1), true)), "int8")
2198 %!test
2199 %! warning ("off", "Octave:num-to-str", "local");
2200 %! assert (class (horzcat (int8 (1), "a")), "char");
2201 
2202 %!assert (class (horzcat (uint64 (1), int64 (1))), "uint64")
2203 %!assert (class (horzcat (uint64 (1), int32 (1))), "uint64")
2204 %!assert (class (horzcat (uint64 (1), int16 (1))), "uint64")
2205 %!assert (class (horzcat (uint64 (1), int8 (1))), "uint64")
2206 %!assert (class (horzcat (uint64 (1), uint64 (1))), "uint64")
2207 %!assert (class (horzcat (uint64 (1), uint32 (1))), "uint64")
2208 %!assert (class (horzcat (uint64 (1), uint16 (1))), "uint64")
2209 %!assert (class (horzcat (uint64 (1), uint8 (1))), "uint64")
2210 %!assert (class (horzcat (uint64 (1), single (1))), "uint64")
2211 %!assert (class (horzcat (uint64 (1), double (1))), "uint64")
2212 %!assert (class (horzcat (uint64 (1), cell (1))), "cell")
2213 %!assert (class (horzcat (uint64 (1), true)), "uint64")
2214 %!test
2215 %! warning ("off", "Octave:num-to-str", "local");
2216 %! assert (class (horzcat (uint64 (1), "a")), "char");
2217 
2218 %!assert (class (horzcat (uint32 (1), int64 (1))), "uint32")
2219 %!assert (class (horzcat (uint32 (1), int32 (1))), "uint32")
2220 %!assert (class (horzcat (uint32 (1), int16 (1))), "uint32")
2221 %!assert (class (horzcat (uint32 (1), int8 (1))), "uint32")
2222 %!assert (class (horzcat (uint32 (1), uint64 (1))), "uint32")
2223 %!assert (class (horzcat (uint32 (1), uint32 (1))), "uint32")
2224 %!assert (class (horzcat (uint32 (1), uint16 (1))), "uint32")
2225 %!assert (class (horzcat (uint32 (1), uint8 (1))), "uint32")
2226 %!assert (class (horzcat (uint32 (1), single (1))), "uint32")
2227 %!assert (class (horzcat (uint32 (1), double (1))), "uint32")
2228 %!assert (class (horzcat (uint32 (1), cell (1))), "cell")
2229 %!assert (class (horzcat (uint32 (1), true)), "uint32")
2230 %!test
2231 %! warning ("off", "Octave:num-to-str", "local");
2232 %! assert (class (horzcat (uint32 (1), "a")), "char");
2233 
2234 %!assert (class (horzcat (uint16 (1), int64 (1))), "uint16")
2235 %!assert (class (horzcat (uint16 (1), int32 (1))), "uint16")
2236 %!assert (class (horzcat (uint16 (1), int16 (1))), "uint16")
2237 %!assert (class (horzcat (uint16 (1), int8 (1))), "uint16")
2238 %!assert (class (horzcat (uint16 (1), uint64 (1))), "uint16")
2239 %!assert (class (horzcat (uint16 (1), uint32 (1))), "uint16")
2240 %!assert (class (horzcat (uint16 (1), uint16 (1))), "uint16")
2241 %!assert (class (horzcat (uint16 (1), uint8 (1))), "uint16")
2242 %!assert (class (horzcat (uint16 (1), single (1))), "uint16")
2243 %!assert (class (horzcat (uint16 (1), double (1))), "uint16")
2244 %!assert (class (horzcat (uint16 (1), cell (1))), "cell")
2245 %!assert (class (horzcat (uint16 (1), true)), "uint16")
2246 %!test
2247 %! warning ("off", "Octave:num-to-str", "local");
2248 %! assert (class (horzcat (uint16 (1), "a")), "char");
2249 
2250 %!assert (class (horzcat (uint8 (1), int64 (1))), "uint8")
2251 %!assert (class (horzcat (uint8 (1), int32 (1))), "uint8")
2252 %!assert (class (horzcat (uint8 (1), int16 (1))), "uint8")
2253 %!assert (class (horzcat (uint8 (1), int8 (1))), "uint8")
2254 %!assert (class (horzcat (uint8 (1), uint64 (1))), "uint8")
2255 %!assert (class (horzcat (uint8 (1), uint32 (1))), "uint8")
2256 %!assert (class (horzcat (uint8 (1), uint16 (1))), "uint8")
2257 %!assert (class (horzcat (uint8 (1), uint8 (1))), "uint8")
2258 %!assert (class (horzcat (uint8 (1), single (1))), "uint8")
2259 %!assert (class (horzcat (uint8 (1), double (1))), "uint8")
2260 %!assert (class (horzcat (uint8 (1), cell (1))), "cell")
2261 %!assert (class (horzcat (uint8 (1), true)), "uint8")
2262 %!test
2263 %! warning ("off", "Octave:num-to-str", "local");
2264 %! assert (class (horzcat (uint8 (1), "a")), "char");
2265 
2266 %!assert (class (horzcat (single (1), int64 (1))), "int64")
2267 %!assert (class (horzcat (single (1), int32 (1))), "int32")
2268 %!assert (class (horzcat (single (1), int16 (1))), "int16")
2269 %!assert (class (horzcat (single (1), int8 (1))), "int8")
2270 %!assert (class (horzcat (single (1), uint64 (1))), "uint64")
2271 %!assert (class (horzcat (single (1), uint32 (1))), "uint32")
2272 %!assert (class (horzcat (single (1), uint16 (1))), "uint16")
2273 %!assert (class (horzcat (single (1), uint8 (1))), "uint8")
2274 %!assert (class (horzcat (single (1), single (1))), "single")
2275 %!assert (class (horzcat (single (1), double (1))), "single")
2276 %!assert (class (horzcat (single (1), cell (1))), "cell")
2277 %!assert (class (horzcat (single (1), true)), "single")
2278 %!test
2279 %! warning ("off", "Octave:num-to-str", "local");
2280 %! assert (class (horzcat (single (1), "a")), "char");
2281 
2282 %!assert (class (horzcat (double (1), int64 (1))), "int64")
2283 %!assert (class (horzcat (double (1), int32 (1))), "int32")
2284 %!assert (class (horzcat (double (1), int16 (1))), "int16")
2285 %!assert (class (horzcat (double (1), int8 (1))), "int8")
2286 %!assert (class (horzcat (double (1), uint64 (1))), "uint64")
2287 %!assert (class (horzcat (double (1), uint32 (1))), "uint32")
2288 %!assert (class (horzcat (double (1), uint16 (1))), "uint16")
2289 %!assert (class (horzcat (double (1), uint8 (1))), "uint8")
2290 %!assert (class (horzcat (double (1), single (1))), "single")
2291 %!assert (class (horzcat (double (1), double (1))), "double")
2292 %!assert (class (horzcat (double (1), cell (1))), "cell")
2293 %!assert (class (horzcat (double (1), true)), "double")
2294 %!test
2295 %! warning ("off", "Octave:num-to-str", "local");
2296 %! assert (class (horzcat (double (1), "a")), "char");
2297 
2298 %!assert (class (horzcat (cell (1), int64 (1))), "cell")
2299 %!assert (class (horzcat (cell (1), int32 (1))), "cell")
2300 %!assert (class (horzcat (cell (1), int16 (1))), "cell")
2301 %!assert (class (horzcat (cell (1), int8 (1))), "cell")
2302 %!assert (class (horzcat (cell (1), uint64 (1))), "cell")
2303 %!assert (class (horzcat (cell (1), uint32 (1))), "cell")
2304 %!assert (class (horzcat (cell (1), uint16 (1))), "cell")
2305 %!assert (class (horzcat (cell (1), uint8 (1))), "cell")
2306 %!assert (class (horzcat (cell (1), single (1))), "cell")
2307 %!assert (class (horzcat (cell (1), double (1))), "cell")
2308 %!assert (class (horzcat (cell (1), cell (1))), "cell")
2309 %!assert (class (horzcat (cell (1), true)), "cell")
2310 %!assert (class (horzcat (cell (1), "a")), "cell")
2311 
2312 %!assert (class (horzcat (true, int64 (1))), "int64")
2313 %!assert (class (horzcat (true, int32 (1))), "int32")
2314 %!assert (class (horzcat (true, int16 (1))), "int16")
2315 %!assert (class (horzcat (true, int8 (1))), "int8")
2316 %!assert (class (horzcat (true, uint64 (1))), "uint64")
2317 %!assert (class (horzcat (true, uint32 (1))), "uint32")
2318 %!assert (class (horzcat (true, uint16 (1))), "uint16")
2319 %!assert (class (horzcat (true, uint8 (1))), "uint8")
2320 %!assert (class (horzcat (true, single (1))), "single")
2321 %!assert (class (horzcat (true, double (1))), "double")
2322 %!assert (class (horzcat (true, cell (1))), "cell")
2323 %!assert (class (horzcat (true, true)), "logical")
2324 %!test
2325 %! warning ("off", "Octave:num-to-str", "local");
2326 %! assert (class (horzcat (true, "a")), "char");
2327 
2328 %!test
2329 %! warning ("off", "Octave:num-to-str", "local");
2330 %! assert (class (horzcat ("a", int64 (1))), "char");
2331 %! assert (class (horzcat ("a", int32 (1))), "char");
2332 %! assert (class (horzcat ("a", int16 (1))), "char");
2333 %! assert (class (horzcat ("a", int8 (1))), "char");
2334 %! assert (class (horzcat ("a", int64 (1))), "char");
2335 %! assert (class (horzcat ("a", int32 (1))), "char");
2336 %! assert (class (horzcat ("a", int16 (1))), "char");
2337 %! assert (class (horzcat ("a", int8 (1))), "char");
2338 %! assert (class (horzcat ("a", single (1))), "char");
2339 %! assert (class (horzcat ("a", double (1))), "char");
2340 %! assert (class (horzcat ("a", cell (1))), "cell");
2341 %! assert (class (horzcat ("a", true)), "char");
2342 %! assert (class (horzcat ("a", "a")), "char");
2343 
2344 %!assert (class (horzcat (cell (1), struct ("foo", "bar"))), "cell")
2345 
2346 %!error horzcat (struct ("foo", "bar"), cell (1))
2347 */
2348 
2349 DEFUN (vertcat, args, ,
2350  doc: /* -*- texinfo -*-
2351 @deftypefn {} {} vertcat (@var{array1}, @var{array2}, @dots{}, @var{arrayN})
2352 Return the vertical concatenation of N-D array objects, @var{array1},
2353 @var{array2}, @dots{}, @var{arrayN} along dimension 1.
2354 
2355 Arrays may also be concatenated vertically using the syntax for creating
2356 new matrices. For example:
2357 
2358 @example
2359 @var{vcat} = [ @var{array1}; @var{array2}; @dots{} ]
2360 @end example
2361 @seealso{cat, horzcat}
2362 @end deftypefn */)
2363 {
2364  return do_cat (args, -1, "vertcat");
2365 }
2366 
2367 /*
2368 %!test
2369 %! c = {"foo"; "bar"; "bazoloa"};
2370 %! assert (vertcat (c, "a", "bc", "def"), {"foo"; "bar"; "bazoloa"; "a"; "bc"; "def"});
2371 */
2372 
2373 DEFUN (cat, args, ,
2374  doc: /* -*- texinfo -*-
2375 @deftypefn {} {} cat (@var{dim}, @var{array1}, @var{array2}, @dots{}, @var{arrayN})
2376 Return the concatenation of N-D array objects, @var{array1},
2377 @var{array2}, @dots{}, @var{arrayN} along dimension @var{dim}.
2378 
2379 @example
2380 @group
2381 A = ones (2, 2);
2382 B = zeros (2, 2);
2383 cat (2, A, B)
2384  @result{} 1 1 0 0
2385  1 1 0 0
2386 @end group
2387 @end example
2388 
2389 Alternatively, we can concatenate @var{A} and @var{B} along the
2390 second dimension in the following way:
2391 
2392 @example
2393 @group
2394 [A, B]
2395 @end group
2396 @end example
2397 
2398 @var{dim} can be larger than the dimensions of the N-D array objects
2399 and the result will thus have @var{dim} dimensions as the
2400 following example shows:
2401 
2402 @example
2403 @group
2404 cat (4, ones (2, 2), zeros (2, 2))
2405  @result{} ans(:,:,1,1) =
2406 
2407  1 1
2408  1 1
2409 
2410  ans(:,:,1,2) =
2411 
2412  0 0
2413  0 0
2414 @end group
2415 @end example
2416 @seealso{horzcat, vertcat}
2417 @end deftypefn */)
2418 {
2419  if (args.length () == 0)
2420  print_usage ();
2421 
2422  int dim = args(0).xint_value ("cat: DIM must be an integer") - 1;
2423 
2424  if (dim < 0)
2425  error ("cat: DIM must be a valid dimension");
2426 
2427  return ovl (do_cat (args.slice (1, args.length () - 1), dim, "cat"));
2428 }
2429 
2430 /*
2431 %!function ret = __testcat (t1, t2, tr, cmplx)
2432 %! assert (cat (1, cast ([], t1), cast ([], t2)), cast ([], tr));
2433 %!
2434 %! assert (cat (1, cast (1, t1), cast (2, t2)), cast ([1; 2], tr));
2435 %! assert (cat (1, cast (1, t1), cast ([2; 3], t2)), cast ([1; 2; 3], tr));
2436 %! assert (cat (1, cast ([1; 2], t1), cast (3, t2)), cast ([1; 2; 3], tr));
2437 %! assert (cat (1, cast ([1; 2], t1), cast ([3; 4], t2)), cast ([1; 2; 3; 4], tr));
2438 %! assert (cat (2, cast (1, t1), cast (2, t2)), cast ([1, 2], tr));
2439 %! assert (cat (2, cast (1, t1), cast ([2, 3], t2)), cast ([1, 2, 3], tr));
2440 %! assert (cat (2, cast ([1, 2], t1), cast (3, t2)), cast ([1, 2, 3], tr));
2441 %! assert (cat (2, cast ([1, 2], t1), cast ([3, 4], t2)), cast ([1, 2, 3, 4], tr));
2442 %!
2443 %! assert ([cast(1, t1); cast(2, t2)], cast ([1; 2], tr));
2444 %! assert ([cast(1, t1); cast([2; 3], t2)], cast ([1; 2; 3], tr));
2445 %! assert ([cast([1; 2], t1); cast(3, t2)], cast ([1; 2; 3], tr));
2446 %! assert ([cast([1; 2], t1); cast([3; 4], t2)], cast ([1; 2; 3; 4], tr));
2447 %! assert ([cast(1, t1), cast(2, t2)], cast ([1, 2], tr));
2448 %! assert ([cast(1, t1), cast([2, 3], t2)], cast ([1, 2, 3], tr));
2449 %! assert ([cast([1, 2], t1), cast(3, t2)], cast ([1, 2, 3], tr));
2450 %! assert ([cast([1, 2], t1), cast([3, 4], t2)], cast ([1, 2, 3, 4], tr));
2451 %!
2452 %! if (nargin == 3 || cmplx)
2453 %! assert (cat (1, cast (1i, t1), cast (2, t2)), cast ([1i; 2], tr));
2454 %! assert (cat (1, cast (1i, t1), cast ([2; 3], t2)), cast ([1i; 2; 3], tr));
2455 %! assert (cat (1, cast ([1i; 2], t1), cast (3, t2)), cast ([1i; 2; 3], tr));
2456 %! assert (cat (1, cast ([1i; 2], t1), cast ([3; 4], t2)), cast ([1i; 2; 3; 4], tr));
2457 %! assert (cat (2, cast (1i, t1), cast (2, t2)), cast ([1i, 2], tr));
2458 %! assert (cat (2, cast (1i, t1), cast ([2, 3], t2)), cast ([1i, 2, 3], tr));
2459 %! assert (cat (2, cast ([1i, 2], t1), cast (3, t2)), cast ([1i, 2, 3], tr));
2460 %! assert (cat (2, cast ([1i, 2], t1), cast ([3, 4], t2)), cast ([1i, 2, 3, 4], tr));
2461 %!
2462 %! assert ([cast(1i, t1); cast(2, t2)], cast ([1i; 2], tr));
2463 %! assert ([cast(1i, t1); cast([2; 3], t2)], cast ([1i; 2; 3], tr));
2464 %! assert ([cast([1i; 2], t1); cast(3, t2)], cast ([1i; 2; 3], tr));
2465 %! assert ([cast([1i; 2], t1); cast([3; 4], t2)], cast ([1i; 2; 3; 4], tr));
2466 %! assert ([cast(1i, t1), cast(2, t2)], cast ([1i, 2], tr));
2467 %! assert ([cast(1i, t1), cast([2, 3], t2)], cast ([1i, 2, 3], tr));
2468 %! assert ([cast([1i, 2], t1), cast(3, t2)], cast ([1i, 2, 3], tr));
2469 %! assert ([cast([1i, 2], t1), cast([3, 4], t2)], cast ([1i, 2, 3, 4], tr));
2470 %!
2471 %! assert (cat (1, cast (1, t1), cast (2i, t2)), cast ([1; 2i], tr));
2472 %! assert (cat (1, cast (1, t1), cast ([2i; 3], t2)), cast ([1; 2i; 3], tr));
2473 %! assert (cat (1, cast ([1; 2], t1), cast (3i, t2)), cast ([1; 2; 3i], tr));
2474 %! assert (cat (1, cast ([1; 2], t1), cast ([3i; 4], t2)), cast ([1; 2; 3i; 4], tr));
2475 %! assert (cat (2, cast (1, t1), cast (2i, t2)), cast ([1, 2i], tr));
2476 %! assert (cat (2, cast (1, t1), cast ([2i, 3], t2)), cast ([1, 2i, 3], tr));
2477 %! assert (cat (2, cast ([1, 2], t1), cast (3i, t2)), cast ([1, 2, 3i], tr));
2478 %! assert (cat (2, cast ([1, 2], t1), cast ([3i, 4], t2)), cast ([1, 2, 3i, 4], tr));
2479 %!
2480 %! assert ([cast(1, t1); cast(2i, t2)], cast ([1; 2i], tr));
2481 %! assert ([cast(1, t1); cast([2i; 3], t2)], cast ([1; 2i; 3], tr));
2482 %! assert ([cast([1; 2], t1); cast(3i, t2)], cast ([1; 2; 3i], tr));
2483 %! assert ([cast([1; 2], t1); cast([3i; 4], t2)], cast ([1; 2; 3i; 4], tr));
2484 %! assert ([cast(1, t1), cast(2i, t2)], cast ([1, 2i], tr));
2485 %! assert ([cast(1, t1), cast([2i, 3], t2)], cast ([1, 2i, 3], tr));
2486 %! assert ([cast([1, 2], t1), cast(3i, t2)], cast ([1, 2, 3i], tr));
2487 %! assert ([cast([1, 2], t1), cast([3i, 4], t2)], cast ([1, 2, 3i, 4], tr));
2488 %!
2489 %! assert (cat (1, cast (1i, t1), cast (2i, t2)), cast ([1i; 2i], tr));
2490 %! assert (cat (1, cast (1i, t1), cast ([2i; 3], t2)), cast ([1i; 2i; 3], tr));
2491 %! assert (cat (1, cast ([1i; 2], t1), cast (3i, t2)), cast ([1i; 2; 3i], tr));
2492 %! assert (cat (1, cast ([1i; 2], t1), cast ([3i; 4], t2)), cast ([1i; 2; 3i; 4], tr));
2493 %! assert (cat (2, cast (1i, t1), cast (2i, t2)), cast ([1i, 2i], tr));
2494 %! assert (cat (2, cast (1i, t1), cast ([2i, 3], t2)), cast ([1i, 2i, 3], tr));
2495 %! assert (cat (2, cast ([1i, 2], t1), cast (3i, t2)), cast ([1i, 2, 3i], tr));
2496 %! assert (cat (2, cast ([1i, 2], t1), cast ([3i, 4], t2)), cast ([1i, 2, 3i, 4], tr));
2497 %!
2498 %! assert ([cast(1i, t1); cast(2i, t2)], cast ([1i; 2i], tr));
2499 %! assert ([cast(1i, t1); cast([2i; 3], t2)], cast ([1i; 2i; 3], tr));
2500 %! assert ([cast([1i; 2], t1); cast(3i, t2)], cast ([1i; 2; 3i], tr));
2501 %! assert ([cast([1i; 2], t1); cast([3i; 4], t2)], cast ([1i; 2; 3i; 4], tr));
2502 %! assert ([cast(1i, t1), cast(2i, t2)], cast ([1i, 2i], tr));
2503 %! assert ([cast(1i, t1), cast([2i, 3], t2)], cast ([1i, 2i, 3], tr));
2504 %! assert ([cast([1i, 2], t1), cast(3i, t2)], cast ([1i, 2, 3i], tr));
2505 %! assert ([cast([1i, 2], t1), cast([3i, 4], t2)], cast ([1i, 2, 3i, 4], tr));
2506 %! endif
2507 %! ret = true;
2508 %!endfunction
2509 
2510 %!assert (__testcat ("double", "double", "double"))
2511 %!assert (__testcat ("single", "double", "single"))
2512 %!assert (__testcat ("double", "single", "single"))
2513 %!assert (__testcat ("single", "single", "single"))
2514 
2515 %!assert (__testcat ("double", "int8", "int8", false))
2516 %!assert (__testcat ("int8", "double", "int8", false))
2517 %!assert (__testcat ("single", "int8", "int8", false))
2518 %!assert (__testcat ("int8", "single", "int8", false))
2519 %!assert (__testcat ("int8", "int8", "int8", false))
2520 %!assert (__testcat ("double", "int16", "int16", false))
2521 %!assert (__testcat ("int16", "double", "int16", false))
2522 %!assert (__testcat ("single", "int16", "int16", false))
2523 %!assert (__testcat ("int16", "single", "int16", false))
2524 %!assert (__testcat ("int16", "int16", "int16", false))
2525 %!assert (__testcat ("double", "int32", "int32", false))
2526 %!assert (__testcat ("int32", "double", "int32", false))
2527 %!assert (__testcat ("single", "int32", "int32", false))
2528 %!assert (__testcat ("int32", "single", "int32", false))
2529 %!assert (__testcat ("int32", "int32", "int32", false))
2530 %!assert (__testcat ("double", "int64", "int64", false))
2531 %!assert (__testcat ("int64", "double", "int64", false))
2532 %!assert (__testcat ("single", "int64", "int64", false))
2533 %!assert (__testcat ("int64", "single", "int64", false))
2534 %!assert (__testcat ("int64", "int64", "int64", false))
2535 
2536 %!assert (__testcat ("double", "uint8", "uint8", false))
2537 %!assert (__testcat ("uint8", "double", "uint8", false))
2538 %!assert (__testcat ("single", "uint8", "uint8", false))
2539 %!assert (__testcat ("uint8", "single", "uint8", false))
2540 %!assert (__testcat ("uint8", "uint8", "uint8", false))
2541 %!assert (__testcat ("double", "uint16", "uint16", false))
2542 %!assert (__testcat ("uint16", "double", "uint16", false))
2543 %!assert (__testcat ("single", "uint16", "uint16", false))
2544 %!assert (__testcat ("uint16", "single", "uint16", false))
2545 %!assert (__testcat ("uint16", "uint16", "uint16", false))
2546 %!assert (__testcat ("double", "uint32", "uint32", false))
2547 %!assert (__testcat ("uint32", "double", "uint32", false))
2548 %!assert (__testcat ("single", "uint32", "uint32", false))
2549 %!assert (__testcat ("uint32", "single", "uint32", false))
2550 %!assert (__testcat ("uint32", "uint32", "uint32", false))
2551 %!assert (__testcat ("double", "uint64", "uint64", false))
2552 %!assert (__testcat ("uint64", "double", "uint64", false))
2553 %!assert (__testcat ("single", "uint64", "uint64", false))
2554 %!assert (__testcat ("uint64", "single", "uint64", false))
2555 %!assert (__testcat ("uint64", "uint64", "uint64", false))
2556 
2557 %!assert (cat (3, [], [1,2;3,4]), [1,2;3,4])
2558 %!assert (cat (3, [1,2;3,4], []), [1,2;3,4])
2559 %!assert (cat (3, [], [1,2;3,4], []), [1,2;3,4])
2560 %!assert (cat (3, [], [], []), zeros (0, 0, 3))
2561 
2562 %!assert (cat (3, [], [], 1, 2), cat (3, 1, 2))
2563 %!assert (cat (3, [], [], [1,2;3,4]), [1,2;3,4])
2564 %!assert (cat (4, [], [], [1,2;3,4]), [1,2;3,4])
2565 
2566 %!assert ([zeros(3,2,2); ones(1,2,2)], repmat ([0;0;0;1],[1,2,2]))
2567 %!assert ([zeros(3,2,2); ones(1,2,2)], vertcat (zeros (3,2,2), ones (1,2,2)))
2568 
2569 %!error <dimension mismatch> cat (3, cat (3, [], []), [1,2;3,4])
2570 %!error <dimension mismatch> cat (3, zeros (0, 0, 2), [1,2;3,4])
2571 */
2572 
2573 static octave_value
2574 do_permute (const octave_value_list& args, bool inv)
2575 {
2576  if (args.length () != 2 || args(1).length () < args(1).ndims ())
2577  print_usage ();
2578 
2579  Array<int> vec = args(1).int_vector_value ();
2580 
2581  // FIXME: maybe we should create an idx_vector object here
2582  // and pass that to permute?
2583  int n = vec.numel ();
2584  for (int i = 0; i < n; i++)
2585  vec(i)--;
2586 
2587  return octave_value (args(0).permute (vec, inv));
2588 }
2589 
2590 DEFUN (permute, args, ,
2591  doc: /* -*- texinfo -*-
2592 @deftypefn {} {} permute (@var{A}, @var{perm})
2593 Return the generalized transpose for an N-D array object @var{A}.
2594 
2595 The permutation vector @var{perm} must contain the elements
2596 @w{@code{1:ndims (A)}} (in any order, but each element must appear only
2597 once). The @var{N}th dimension of @var{A} gets remapped to dimension
2598 @code{@var{PERM}(@var{N})}. For example:
2599 
2600 @example
2601 @group
2602 @var{x} = zeros ([2, 3, 5, 7]);
2603 size (@var{x})
2604  @result{} 2 3 5 7
2605 
2606 size (permute (@var{x}, [2, 1, 3, 4]))
2607  @result{} 3 2 5 7
2608 
2609 size (permute (@var{x}, [1, 3, 4, 2]))
2610  @result{} 2 5 7 3
2611 
2612 ## The identity permutation
2613 size (permute (@var{x}, [1, 2, 3, 4]))
2614  @result{} 2 3 5 7
2615 @end group
2616 @end example
2617 @seealso{ipermute}
2618 @end deftypefn */)
2619 {
2620  return do_permute (args, false);
2621 }
2622 
2623 DEFUN (ipermute, args, ,
2624  doc: /* -*- texinfo -*-
2625 @deftypefn {} {} ipermute (@var{A}, @var{iperm})
2626 The inverse of the @code{permute} function.
2627 
2628 The expression
2629 
2630 @example
2631 ipermute (permute (A, perm), perm)
2632 @end example
2633 
2634 @noindent
2635 returns the original array @var{A}.
2636 @seealso{permute}
2637 @end deftypefn */)
2638 {
2639  return do_permute (args, true);
2640 }
2641 
2642 DEFUN (length, args, ,
2643  doc: /* -*- texinfo -*-
2644 @deftypefn {} {} length (@var{a})
2645 Return the length of the object @var{a}.
2646 
2647 The length is 0 for empty objects, 1 for scalars, and the number of elements
2648 for vectors. For matrix or N-dimensional objects, the length is the number
2649 of elements along the largest dimension
2650 (equivalent to @w{@code{max (size (@var{a}))}}).
2651 @seealso{numel, size}
2652 @end deftypefn */)
2653 {
2654  if (args.length () != 1)
2655  print_usage ();
2656 
2657  return ovl (args(0).length ());
2658 }
2659 
2660 DEFUN (ndims, args, ,
2661  doc: /* -*- texinfo -*-
2662 @deftypefn {} {} ndims (@var{a})
2663 Return the number of dimensions of @var{a}.
2664 
2665 For any array, the result will always be greater than or equal to 2.
2666 Trailing singleton dimensions are not counted.
2667 
2668 @example
2669 @group
2670 ndims (ones (4, 1, 2, 1))
2671  @result{} 3
2672 @end group
2673 @end example
2674 @seealso{size}
2675 @end deftypefn */)
2676 {
2677  if (args.length () != 1)
2678  print_usage ();
2679 
2680  return ovl (args(0).ndims ());
2681 }
2682 
2683 DEFUN (numel, args, ,
2684  doc: /* -*- texinfo -*-
2685 @deftypefn {} {} numel (@var{a})
2686 @deftypefnx {} {} numel (@var{a}, @var{idx1}, @var{idx2}, @dots{})
2687 Return the number of elements in the object @var{a}.
2688 
2689 Optionally, if indices @var{idx1}, @var{idx2}, @dots{} are supplied,
2690 return the number of elements that would result from the indexing
2691 
2692 @example
2693 @var{a}(@var{idx1}, @var{idx2}, @dots{})
2694 @end example
2695 
2696 Note that the indices do not have to be scalar numbers. For example,
2697 
2698 @example
2699 @group
2700 @var{a} = 1;
2701 @var{b} = ones (2, 3);
2702 numel (@var{a}, @var{b})
2703 @end group
2704 @end example
2705 
2706 @noindent
2707 will return 6, as this is the number of ways to index with @var{b}.
2708 Or the index could be the string @qcode{":"} which represents the colon
2709 operator. For example,
2710 
2711 @example
2712 @group
2713 @var{a} = ones (5, 3);
2714 numel (@var{a}, 2, ":")
2715 @end group
2716 @end example
2717 
2718 @noindent
2719 will return 3 as the second row has three column entries.
2720 
2721 This method is also called when an object appears as lvalue with cs-list
2722 indexing, i.e., @code{object@{@dots{}@}} or @code{object(@dots{}).field}.
2723 @seealso{size, length, ndims}
2724 @end deftypefn */)
2725 {
2726  int nargin = args.length ();
2727 
2728  if (nargin == 0)
2729  print_usage ();
2730 
2732 
2733  if (nargin == 1)
2734  retval = args(0).numel ();
2735  else if (nargin > 1)
2736  {
2737  // Don't use numel (const octave_value_list&) here as that corresponds to
2738  // an overloaded call, not to builtin!
2739  retval = dims_to_numel (args(0).dims (), args.slice (1, nargin-1));
2740  }
2741 
2742  return retval;
2743 }
2744 
2745 DEFUN (size, args, nargout,
2746  doc: /* -*- texinfo -*-
2747 @deftypefn {} {@var{sz} =} size (@var{a})
2748 @deftypefnx {} {@var{dim_sz} =} size (@var{a}, @var{dim})
2749 @deftypefnx {} {[@var{rows}, @var{cols}, @dots{}, @var{dim_N_sz}] =} size (@dots{})
2750 Return a row vector with the size (number of elements) of each dimension for
2751 the object @var{a}.
2752 
2753 When given a second argument, @var{dim}, return the size of the corresponding
2754 dimension.
2755 
2756 With a single output argument, @code{size} returns a row vector. When called
2757 with multiple output arguments, @code{size} returns the size of dimension N
2758 in the Nth argument. The number of rows, dimension 1, is returned in the
2759 first argument, the number of columns, dimension 2, is returned in the
2760 second argument, etc. If there are more dimensions in @var{a} then there are
2761 output arguments, @code{size} returns the total number of elements in the
2762 remaining dimensions in the final output argument.
2763 
2764 Example 1: single row vector output
2765 
2766 @example
2767 @group
2768 size ([1, 2; 3, 4; 5, 6])
2769  @result{} [ 3, 2 ]
2770 @end group
2771 @end example
2772 
2773 Example 2: number of elements in 2nd dimension (columns)
2774 
2775 @example
2776 @group
2777 size ([1, 2; 3, 4; 5, 6], 2)
2778  @result{} 2
2779 @end group
2780 @end example
2781 
2782 Example 3: number of output arguments == number of dimensions
2783 
2784 @example
2785 @group
2786 [nr, nc] = size ([1, 2; 3, 4; 5, 6])
2787  @result{} nr = 3
2788  @result{} nc = 2
2789 @end group
2790 @end example
2791 
2792 Example 4: number of output arguments != number of dimensions
2793 
2794 @example
2795 @group
2796 [nr, remainder] = size (ones (2, 3, 4, 5))
2797  @result{} nr = 2
2798  @result{} remainder = 60
2799 @end group
2800 @end example
2801 
2802 @seealso{numel, ndims, length, rows, columns, size_equal, common_size}
2803 @end deftypefn */)
2804 {
2806 
2807  int nargin = args.length ();
2808 
2809  if (nargin == 1)
2810  {
2811  const dim_vector dimensions = args(0).dims ();
2812 
2813  if (nargout > 1)
2814  {
2815  const dim_vector rdims = dimensions.redim (nargout);
2816  retval.resize (nargout);
2817  for (int i = 0; i < nargout; i++)
2818  retval(i) = rdims(i);
2819  }
2820  else
2821  {
2822  int ndims = dimensions.ndims ();
2823 
2824  NoAlias<Matrix> m (1, ndims);
2825 
2826  for (int i = 0; i < ndims; i++)
2827  m(i) = dimensions(i);
2828 
2829  retval(0) = m;
2830  }
2831  }
2832  else if (nargin == 2 && nargout < 2)
2833  {
2834  if (! args(1).is_real_scalar ())
2835  error ("size: DIM must be a positive integer");
2836 
2837  octave_idx_type nd = args(1).idx_type_value ();
2838 
2839  const dim_vector dv = args(0).dims ();
2840 
2841  if (nd < 1)
2842  error ("size: requested dimension DIM (= %d) out of range", nd);
2843 
2844  if (nd <= dv.ndims ())
2845  retval(0) = dv(nd-1);
2846  else
2847  retval(0) = 1;
2848  }
2849  else
2850  print_usage ();
2851 
2852  return retval;
2853 }
2854 
2855 DEFUN (size_equal, args, ,
2856  doc: /* -*- texinfo -*-
2857 @deftypefn {} {} size_equal (@var{a}, @var{b}, @dots{})
2858 Return true if the dimensions of all arguments agree.
2859 
2860 Trailing singleton dimensions are ignored. When called with a single argument,
2861 or no argument, @code{size_equal} returns true.
2862 @seealso{size, numel, ndims, common_size}
2863 @end deftypefn */)
2864 {
2865  int nargin = args.length ();
2866 
2867  if (nargin >= 1)
2868  {
2869  dim_vector a_dims = args(0).dims ();
2870 
2871  for (int i = 1; i < nargin; ++i)
2872  {
2873  dim_vector b_dims = args(i).dims ();
2874 
2875  if (a_dims != b_dims)
2876  return ovl (false);
2877  }
2878  }
2879 
2880  return ovl (true);
2881 }
2882 
2883 DEFUN (nnz, args, ,
2884  doc: /* -*- texinfo -*-
2885 @deftypefn {} {@var{n} =} nnz (@var{a})
2886 Return the number of nonzero elements in @var{a}.
2887 @seealso{nzmax, nonzeros, find}
2888 @end deftypefn */)
2889 {
2890  if (args.length () != 1)
2891  print_usage ();
2892 
2893  return ovl (args(0).nnz ());
2894 }
2895 
2896 DEFUN (nzmax, args, ,
2897  doc: /* -*- texinfo -*-
2898 @deftypefn {} {@var{n} =} nzmax (@var{SM})
2899 Return the amount of storage allocated to the sparse matrix @var{SM}.
2900 
2901 Note that Octave tends to crop unused memory at the first opportunity
2902 for sparse objects. Thus, in general the value of @code{nzmax} will be the
2903 same as @code{nnz} except for some cases of user-created sparse objects.
2904 @seealso{nnz, spalloc, sparse}
2905 @end deftypefn */)
2906 {
2907  if (args.length () != 1)
2908  print_usage ();
2909 
2910  return ovl (args(0).nzmax ());
2911 }
2912 
2913 DEFUN (rows, args, ,
2914  doc: /* -*- texinfo -*-
2915 @deftypefn {} {} rows (@var{a})
2916 Return the number of rows of @var{a}.
2917 @seealso{columns, size, length, numel, isscalar, isvector, ismatrix}
2918 @end deftypefn */)
2919 {
2920  if (args.length () != 1)
2921  print_usage ();
2922 
2923  return ovl (args(0).rows ());
2924 }
2925 
2926 /*
2927 %!assert (rows (ones (2,5)), 2)
2928 %!assert (rows (ones (5,2)), 5)
2929 %!assert (rows (ones (5,4,3,2)), 5)
2930 %!assert (rows (ones (3,4,5,2)), 3)
2931 
2932 %!assert (rows (cell (2,5)), 2)
2933 %!assert (rows (cell (5,2)), 5)
2934 %!assert (rows (cell (5,4,3,2)), 5)
2935 %!assert (rows (cell (3,4,5,2)), 3)
2936 
2937 %!test
2938 %! x(2,5,3).a = 1;
2939 %! assert (rows (x), 2);
2940 %! y(5,4,3).b = 2;
2941 %! assert (rows (y), 5);
2942 
2943 %!assert (rows ("Hello World"), 1)
2944 
2945 %!assert (rows ([]), 0)
2946 %!assert (rows (zeros (2,0)), 2)
2947 
2948 ## Test input validation
2949 %!error rows ()
2950 %!error rows (1,2)
2951 */
2952 
2953 DEFUN (columns, args, ,
2954  doc: /* -*- texinfo -*-
2955 @deftypefn {} {} columns (@var{a})
2956 Return the number of columns of @var{a}.
2957 @seealso{rows, size, length, numel, isscalar, isvector, ismatrix}
2958 @end deftypefn */)
2959 {
2960  if (args.length () != 1)
2961  print_usage ();
2962 
2963  return ovl (args(0).columns ());
2964 }
2965 
2966 DEFUN (sum, args, ,
2967  doc: /* -*- texinfo -*-
2968 @deftypefn {} {} sum (@var{x})
2969 @deftypefnx {} {} sum (@var{x}, @var{dim})
2970 @deftypefnx {} {} sum (@dots{}, "native")
2971 @deftypefnx {} {} sum (@dots{}, "double")
2972 @deftypefnx {} {} sum (@dots{}, "extra")
2973 Sum of elements along dimension @var{dim}.
2974 
2975 If @var{dim} is omitted, it defaults to the first non-singleton dimension.
2976 
2977 The optional @qcode{"type"} input determines the class of the variable
2978 used for calculations. If the argument @qcode{"native"} is given, then
2979 the operation is performed in the same type as the original argument, rather
2980 than the default double type.
2981 
2982 For example:
2983 
2984 @example
2985 @group
2986 sum ([true, true])
2987  @result{} 2
2988 sum ([true, true], "native")
2989  @result{} true
2990 @end group
2991 @end example
2992 
2993 On the contrary, if @qcode{"double"} is given, the sum is performed in
2994 double precision even for single precision inputs.
2995 
2996 For double precision inputs, the @qcode{"extra"} option will use a more
2997 accurate algorithm than straightforward summation. For single precision
2998 inputs, @qcode{"extra"} is the same as @qcode{"double"}. Otherwise,
2999 @qcode{"extra"} has no effect.
3000 @seealso{cumsum, sumsq, prod}
3001 @end deftypefn */)
3002 {
3003  int nargin = args.length ();
3004 
3005  bool isnative = false;
3006  bool isdouble = false;
3007  bool isextra = false;
3008 
3009  if (nargin > 1 && args(nargin - 1).is_string ())
3010  {
3011  std::string str = args(nargin - 1).string_value ();
3012 
3013  if (str == "native")
3014  isnative = true;
3015  else if (str == "double")
3016  isdouble = true;
3017  else if (str == "extra")
3018  isextra = true;
3019  else
3020  error ("sum: unrecognized type argument '%s'", str.c_str ());
3021 
3022  nargin--;
3023  }
3024 
3025  if (nargin < 1 || nargin > 2)
3026  print_usage ();
3027 
3028  int dim = -1;
3029  if (nargin == 2)
3030  {
3031  dim = args(1).int_value () - 1;
3032  if (dim < 0)
3033  error ("sum: invalid dimension DIM = %d", dim + 1);
3034  }
3035 
3037  octave_value arg = args(0);
3038 
3039  switch (arg.builtin_type ())
3040  {
3041  case btyp_double:
3042  if (arg.is_sparse_type ())
3043  {
3044  if (isextra)
3045  warning ("sum: 'extra' not yet implemented for sparse matrices");
3046  retval = arg.sparse_matrix_value ().sum (dim);
3047  }
3048  else if (isextra)
3049  retval = arg.array_value ().xsum (dim);
3050  else
3051  retval = arg.array_value ().sum (dim);
3052  break;
3053 
3054  case btyp_complex:
3055  if (arg.is_sparse_type ())
3056  {
3057  if (isextra)
3058  warning ("sum: 'extra' not yet implemented for sparse matrices");
3059  retval = arg.sparse_complex_matrix_value ().sum (dim);
3060  }
3061  else if (isextra)
3062  retval = arg.complex_array_value ().xsum (dim);
3063  else
3064  retval = arg.complex_array_value ().sum (dim);
3065  break;
3066 
3067  case btyp_float:
3068  if (isdouble || isextra)
3069  retval = arg.float_array_value ().dsum (dim);
3070  else
3071  retval = arg.float_array_value ().sum (dim);
3072  break;
3073 
3074  case btyp_float_complex:
3075  if (isdouble || isextra)
3076  retval = arg.float_complex_array_value ().dsum (dim);
3077  else
3078  retval = arg.float_complex_array_value ().sum (dim);
3079  break;
3080 
3081 #define MAKE_INT_BRANCH(X) \
3082  case btyp_ ## X: \
3083  if (isnative) \
3084  retval = arg.X ## _array_value ().sum (dim); \
3085  else \
3086  retval = arg.X ## _array_value ().dsum (dim); \
3087  break;
3088 
3089  MAKE_INT_BRANCH (int8);
3093  MAKE_INT_BRANCH (uint8);
3096  MAKE_INT_BRANCH (uint64);
3097 
3098 #undef MAKE_INT_BRANCH
3099 
3100  // GAGME: Accursed Matlab compatibility...
3101  case btyp_char:
3102  if (isextra)
3103  retval = arg.array_value (true).xsum (dim);
3104  else
3105  retval = arg.array_value (true).sum (dim);
3106  break;
3107 
3108  case btyp_bool:
3109  if (arg.is_sparse_type ())
3110  {
3111  if (isnative)
3112  retval = arg.sparse_bool_matrix_value ().any (dim);
3113  else
3114  retval = arg.sparse_bool_matrix_value ().sum (dim);
3115  }
3116  else if (isnative)
3117  retval = arg.bool_array_value ().any (dim);
3118  else
3119  retval = arg.bool_array_value ().sum (dim);
3120  break;
3121 
3122  default:
3123  err_wrong_type_arg ("sum", arg);
3124  }
3125 
3126  return retval;
3127 }
3128 
3129 /*
3130 %!assert (sum ([1, 2, 3]), 6)
3131 %!assert (sum ([-1; -2; -3]), -6)
3132 %!assert (sum ([i, 2+i, -3+2i, 4]), 3+4i)
3133 %!assert (sum ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i]), [2+2i, 4+4i, 6+6i])
3134 
3135 %!assert (sum (single ([1, 2, 3])), single (6))
3136 %!assert (sum (single ([-1; -2; -3])), single (-6))
3137 %!assert (sum (single ([i, 2+i, -3+2i, 4])), single (3+4i))
3138 %!assert (sum (single ([1, 2, 3; i, 2i, 3i; 1+i, 2+2i, 3+3i])), single ([2+2i, 4+4i, 6+6i]))
3139 
3140 %!assert (sum ([1, 2; 3, 4], 1), [4, 6])
3141 %!assert (sum ([1, 2; 3, 4], 2), [3; 7])
3142 %!assert (sum (zeros (1, 0)), 0)
3143 %!assert (sum (zeros (1, 0), 1), zeros (1, 0))
3144 %!assert (sum (zeros (1, 0), 2), 0)
3145 %!assert (sum (zeros (0, 1)), 0)
3146 %!assert (sum (zeros (0, 1), 1), 0)
3147 %!assert (sum (zeros (0, 1), 2), zeros (0, 1))
3148 %!assert (sum (zeros (2, 0)), zeros (1, 0))
3149 %!assert (sum (zeros (2, 0), 1), zeros (1, 0))
3150 %!assert (sum (zeros (2, 0), 2), [0; 0])
3151 %!assert (sum (zeros (0, 2)), [0, 0])
3152 %!assert (sum (zeros (0, 2), 1), [0, 0])
3153 %!assert (sum (zeros (0, 2), 2), zeros (0, 1))
3154 %!assert (sum (zeros (2, 2, 0, 3)), zeros (1, 2, 0, 3))
3155 %!assert (sum (zeros (2, 2, 0, 3), 2), zeros (2, 1, 0, 3))
3156 %!assert (sum (zeros (2, 2, 0, 3), 3), zeros (2, 2, 1, 3))
3157 %!assert (sum (zeros (2, 2, 0, 3), 4), zeros (2, 2, 0))
3158 %!assert (sum (zeros (2, 2, 0, 3), 7), zeros (2, 2, 0, 3))
3159 
3160 %!assert (sum (single ([1, 2; 3, 4]), 1), single ([4, 6]))
3161 %!assert (sum (single ([1, 2; 3, 4]), 2), single ([3; 7]))
3162 %!assert (sum (zeros (1, 0, "single")), single (0))
3163 %!assert (sum (zeros (1, 0, "single"), 1), zeros (1, 0, "single"))
3164 %!assert (sum (zeros (1, 0, "single"), 2), single (0))
3165 %!assert (sum (zeros (0, 1, "single")), single (0))
3166 %!assert (sum (zeros (0, 1, "single"), 1), single (0))
3167 %!assert (sum (zeros (0, 1, "single"), 2), zeros (0, 1, "single"))
3168 %!assert (sum (zeros (2, 0, "single")), zeros (1, 0, "single"))
3169 %!assert (sum (zeros (2, 0, "single"), 1), zeros (1, 0, "single"))
3170 %!assert (sum (zeros (2, 0, "single"), 2), single ([0; 0]))
3171 %!assert (sum (zeros (0, 2, "single")), single ([0, 0]))
3172 %!assert (sum (zeros (0, 2, "single"), 1), single ([0, 0]))
3173 %!assert (sum (zeros (0, 2, "single"), 2), zeros (0, 1, "single"))
3174 %!assert (sum (zeros (2, 2, 0, 3, "single")), zeros (1, 2, 0, 3, "single"))
3175 %!assert (sum (zeros (2, 2, 0, 3, "single"), 2), zeros (2, 1, 0, 3, "single"))
3176 %!assert (sum (zeros (2, 2, 0, 3, "single"), 3), zeros (2, 2, 1, 3, "single"))
3177 %!assert (sum (zeros (2, 2, 0, 3, "single"), 4), zeros (2, 2, 0, "single"))
3178 %!assert (sum (zeros (2, 2, 0, 3, "single"), 7), zeros (2, 2, 0, 3, "single"))
3179 
3180 ## Test "native"
3181 %!assert (sum ([true,true]), 2)
3182 %!assert (sum ([true,true], "native"), true)
3183 %!assert (sum (int8 ([127,10,-20])), 117)
3184 %!assert (sum (int8 ([127,10,-20]), "native"), int8 (107))
3185 
3186 ;-)
3187 %!assert (sum ("Octave") + "8", sumsq (primes (17)))
3188 
3189 %!error sum ()
3190 %!error sum (1,2,3)
3191 %!error <unrecognized type argument 'foobar'> sum (1, "foobar")
3192 */
3193 
3194 DEFUN (sumsq, args, ,
3195  doc: /* -*- texinfo -*-
3196 @deftypefn {} {} sumsq (@var{x})
3197 @deftypefnx {} {} sumsq (@var{x}, @var{dim})
3198 Sum of squares of elements along dimension @var{dim}.
3199 
3200 If @var{dim} is omitted, it defaults to the first non-singleton dimension.
3201 
3202 This function is conceptually equivalent to computing
3203 
3204 @example
3205 sum (x .* conj (x), dim)
3206 @end example
3207 
3208 @noindent
3209 but it uses less memory and avoids calling @code{conj} if @var{x} is real.
3210 @seealso{sum, prod}
3211 @end deftypefn */)
3212 {
3213  DATA_REDUCTION (sumsq);
3214 }
3215 
3216 /*
3217 %!assert (sumsq ([1, 2, 3]), 14)
3218 %!assert (sumsq ([-1; -2; 4i]), 21)
3219 %!assert (sumsq ([1, 2, 3; 2, 3, 4; 4i, 6i, 2]), [21, 49, 29])
3220 
3221 %!assert (sumsq (single ([1, 2, 3])), single (14))
3222 %!assert (sumsq (single ([-1; -2; 4i])), single (21))
3223 %!assert (sumsq (single ([1, 2, 3; 2, 3, 4; 4i, 6i, 2])), single ([21, 49, 29]))
3224 
3225 %!assert (sumsq ([1, 2; 3, 4], 1), [10, 20])
3226 %!assert (sumsq ([1, 2; 3, 4], 2), [5; 25])
3227 
3228 %!assert (sumsq (single ([1, 2; 3, 4]), 1), single ([10, 20]))
3229 %!assert (sumsq (single ([1, 2; 3, 4]), 2), single ([5; 25]))
3230 
3231 %!error sumsq ()
3232 */
3233 
3234 DEFUN (islogical, args, ,
3235  doc: /* -*- texinfo -*-
3236 @deftypefn {} {} islogical (@var{x})
3237 @deftypefnx {} {} isbool (@var{x})
3238 Return true if @var{x} is a logical object.
3239 @seealso{isfloat, isinteger, ischar, isnumeric, isa}
3240 @end deftypefn */)
3241 {
3242  if (args.length () != 1)
3243  print_usage ();
3244 
3245  return ovl (args(0).is_bool_type ());
3246 }
3247 
3248 DEFALIAS (isbool, islogical);
3249 
3250 /*
3251 %!assert (islogical (true), true)
3252 %!assert (islogical (false), true)
3253 %!assert (islogical ([true, false]), true)
3254 %!assert (islogical (1), false)
3255 %!assert (islogical (1i), false)
3256 %!assert (islogical ([1,1]), false)
3257 %!assert (islogical (single (1)), false)
3258 %!assert (islogical (single (1i)), false)
3259 %!assert (islogical (single ([1,1])), false)
3260 %!assert (islogical (sparse ([true, false])), true)
3261 %!assert (islogical (sparse ([1, 0])), false)
3262 */
3263 
3264 DEFUN (isinteger, args, ,
3265  doc: /* -*- texinfo -*-
3266 @deftypefn {} {} isinteger (@var{x})
3267 Return true if @var{x} is an integer object (int8, uint8, int16, etc.).
3268 
3269 Note that @w{@code{isinteger (14)}} is false because numeric constants in
3270 Octave are double precision floating point values.
3271 @seealso{isfloat, ischar, islogical, isnumeric, isa}
3272 @end deftypefn */)
3273 {
3274  if (args.length () != 1)
3275  print_usage ();
3276 
3277  return ovl (args(0).is_integer_type ());
3278 }
3279 
3280 /*
3281 %!assert (isinteger (int8 (16)))
3282 %!assert (isinteger (int16 (16)))
3283 %!assert (isinteger (int32 (16)))
3284 %!assert (isinteger (int64 (16)))
3285 
3286 %!assert (isinteger (uint8 (16)))
3287 %!assert (isinteger (uint16 (16)))
3288 %!assert (isinteger (uint32 (16)))
3289 %!assert (isinteger (uint64 (16)))
3290 
3291 %!assert (isinteger (intmax ("int8")))
3292 %!assert (isinteger (intmax ("int16")))
3293 %!assert (isinteger (intmax ("int32")))
3294 %!assert (isinteger (intmax ("int64")))
3295 
3296 %!assert (isinteger (intmax ("uint8")))
3297 %!assert (isinteger (intmax ("uint16")))
3298 %!assert (isinteger (intmax ("uint32")))
3299 %!assert (isinteger (intmax ("uint64")))
3300 
3301 %!assert (isinteger (intmin ("int8")))
3302 %!assert (isinteger (intmin ("int16")))
3303 %!assert (isinteger (intmin ("int32")))
3304 %!assert (isinteger (intmin ("int64")))
3305 
3306 %!assert (isinteger (intmin ("uint8")))
3307 %!assert (isinteger (intmin ("uint16")))
3308 %!assert (isinteger (intmin ("uint32")))
3309 %!assert (isinteger (intmin ("uint64")))
3310 
3311 %!assert (isinteger (uint8 ([1:10])))
3312 %!assert (isinteger (uint8 ([1:10; 1:10])))
3313 
3314 %!assert (! isinteger (16))
3315 %!assert (! isinteger ("parrot"))
3316 %!assert (! isinteger ([1, 2, 3]))
3317 
3318 %!error isinteger ()
3319 %!error isinteger ("multiple", "parameters")
3320 */
3321 
3322 DEFUN (iscomplex, args, ,
3323  doc: /* -*- texinfo -*-
3324 @deftypefn {} {} iscomplex (@var{x})
3325 Return true if @var{x} is a complex-valued numeric object.
3326 @seealso{isreal, isnumeric, islogical, ischar, isfloat, isa}
3327 @end deftypefn */)
3328 {
3329  if (args.length () != 1)
3330  print_usage ();
3331 
3332  return ovl (args(0).is_complex_type ());
3333 }
3334 
3335 DEFUN (isfloat, args, ,
3336  doc: /* -*- texinfo -*-
3337 @deftypefn {} {} isfloat (@var{x})
3338 Return true if @var{x} is a floating-point numeric object.
3339 
3340 Objects of class double or single are floating-point objects.
3341 @seealso{isinteger, ischar, islogical, isnumeric, isa}
3342 @end deftypefn */)
3343 {
3344  if (args.length () != 1)
3345  print_usage ();
3346 
3347  return ovl (args(0).is_float_type ());
3348 }
3349 
3350 // FIXME: perhaps this should be implemented with an
3351 // octave_value member function?
3352 
3353 DEFUN (complex, args, ,
3354  doc: /* -*- texinfo -*-
3355 @deftypefn {} {} complex (@var{x})
3356 @deftypefnx {} {} complex (@var{re}, @var{im})
3357 Return a complex value from real arguments.
3358 
3359 With 1 real argument @var{x}, return the complex result
3360 @w{@code{@var{x} + 0i}}.
3361 
3362 With 2 real arguments, return the complex result
3363 @w{@code{@var{re} + @var{im}}}.
3364 @code{complex} can often be more convenient than expressions such as
3365 @w{@code{a + i*b}}.
3366 For example:
3367 
3368 @example
3369 @group
3370 complex ([1, 2], [3, 4])
3371  @result{} [ 1 + 3i 2 + 4i ]
3372 @end group
3373 @end example
3374 @seealso{real, imag, iscomplex, abs, arg}
3375 @end deftypefn */)
3376 {
3377  int nargin = args.length ();
3378 
3379  if (nargin < 1 || nargin > 2)
3380  print_usage ();
3381 
3383 
3384  if (nargin == 1)
3385  {
3386  octave_value arg = args(0);
3387 
3388  if (arg.is_complex_type ())
3389  retval = arg;
3390  else
3391  {
3392  if (arg.is_sparse_type ())
3393  {
3394  SparseComplexMatrix val = arg.xsparse_complex_matrix_value ("complex: invalid conversion");
3395 
3396  retval = octave_value (new octave_sparse_complex_matrix (val));
3397  }
3398  else if (arg.is_single_type ())
3399  {
3400  if (arg.numel () == 1)
3401  {
3402  FloatComplex val = arg.xfloat_complex_value ("complex: invalid conversion");
3403 
3404  retval = octave_value (new octave_float_complex (val));
3405  }
3406  else
3407  {
3408  FloatComplexNDArray val = arg.xfloat_complex_array_value ("complex: invalid conversion");
3409 
3410  retval = octave_value (new octave_float_complex_matrix (val));
3411  }
3412  }
3413  else
3414  {
3415  if (arg.numel () == 1)
3416  {
3417  Complex val = arg.xcomplex_value ("complex: invalid conversion");
3418 
3419  retval = octave_value (new octave_complex (val));
3420  }
3421  else
3422  {
3423  ComplexNDArray val = arg.xcomplex_array_value ("complex: invalid conversion");
3424 
3425  retval = octave_value (new octave_complex_matrix (val));
3426  }
3427  }
3428  }
3429  }
3430  else
3431  {
3432  octave_value re = args(0);
3433  octave_value im = args(1);
3434 
3435  if (re.is_sparse_type () && im.is_sparse_type ())
3436  {
3437  const SparseMatrix re_val = re.sparse_matrix_value ();
3438  const SparseMatrix im_val = im.sparse_matrix_value ();
3439 
3440  if (re.numel () == 1)
3441  {
3443  if (re_val.nnz () == 0)
3444  result = Complex (0, 1) * SparseComplexMatrix (im_val);
3445  else
3446  {
3447  octave_idx_type nr = im_val.rows ();
3448  octave_idx_type nc = im_val.cols ();
3449  result = SparseComplexMatrix (nr, nc, re_val(0));
3450 
3451  for (octave_idx_type j = 0; j < nc; j++)
3452  {
3453  octave_idx_type off = j * nr;
3454  for (octave_idx_type i = im_val.cidx (j);
3455  i < im_val.cidx (j + 1); i++)
3456  result.data (im_val.ridx (i) + off) +=
3457  Complex (0, im_val.data (i));
3458  }
3459  }
3460  retval = octave_value (new octave_sparse_complex_matrix (result));
3461  }
3462  else if (im.numel () == 1)
3463  {
3465  if (im_val.nnz () == 0)
3466  result = SparseComplexMatrix (re_val);
3467  else
3468  {
3469  octave_idx_type nr = re_val.rows ();
3470  octave_idx_type nc = re_val.cols ();
3471  result = SparseComplexMatrix (nr, nc,
3472  Complex (0, im_val(0)));
3473 
3474  for (octave_idx_type j = 0; j < nc; j++)
3475  {
3476  octave_idx_type off = j * nr;
3477  for (octave_idx_type i = re_val.cidx (j);
3478  i < re_val.cidx (j + 1); i++)
3479  result.data (re_val.ridx (i) + off) +=
3480  re_val.data (i);
3481  }
3482  }
3483  retval = octave_value (new octave_sparse_complex_matrix (result));
3484  }
3485  else
3486  {
3487  if (re_val.dims () != im_val.dims ())
3488  error ("complex: dimension mismatch");
3489 
3491  result = SparseComplexMatrix (re_val)
3492  + Complex (0, 1) * SparseComplexMatrix (im_val);
3493  retval = octave_value (new octave_sparse_complex_matrix (result));
3494  }
3495  }
3496  else if (re.is_single_type () || im.is_single_type ())
3497  {
3498  if (re.numel () == 1)
3499  {
3500  float re_val = re.float_value ();
3501 
3502  if (im.numel () == 1)
3503  {
3504  float im_val = im.double_value ();
3505 
3506  retval = octave_value (new octave_float_complex
3507  (FloatComplex (re_val, im_val)));
3508  }
3509  else
3510  {
3511  const FloatNDArray im_val = im.float_array_value ();
3512 
3513  FloatComplexNDArray result (im_val.dims (),
3514  FloatComplex ());
3515 
3516  for (octave_idx_type i = 0; i < im_val.numel (); i++)
3517  result.xelem (i) = FloatComplex (re_val, im_val(i));
3518 
3520  (result));
3521  }
3522  }
3523  else
3524  {
3525  const FloatNDArray re_val = re.float_array_value ();
3526 
3527  if (im.numel () == 1)
3528  {
3529  float im_val = im.float_value ();
3530 
3531  FloatComplexNDArray result (re_val.dims (),
3532  FloatComplex ());
3533 
3534  for (octave_idx_type i = 0; i < re_val.numel (); i++)
3535  result.xelem (i) = FloatComplex (re_val(i), im_val);
3536 
3538  (result));
3539  }
3540  else
3541  {
3542  const FloatNDArray im_val = im.float_array_value ();
3543 
3544  if (re_val.dims () != im_val.dims ())
3545  error ("complex: dimension mismatch");
3546 
3547  FloatComplexNDArray result (re_val.dims (),
3548  FloatComplex ());
3549 
3550  for (octave_idx_type i = 0; i < re_val.numel (); i++)
3551  result.xelem (i) = FloatComplex (re_val(i),
3552  im_val(i));
3553 
3555  (result));
3556  }
3557  }
3558  }
3559  else if (re.numel () == 1)
3560  {
3561  double re_val = re.double_value ();
3562 
3563  if (im.numel () == 1)
3564  {
3565  double im_val = im.double_value ();
3566 
3567  retval = octave_value (new octave_complex
3568  (Complex (re_val, im_val)));
3569  }
3570  else
3571  {
3572  const NDArray im_val = im.array_value ();
3573 
3574  ComplexNDArray result (im_val.dims (), Complex ());
3575 
3576  for (octave_idx_type i = 0; i < im_val.numel (); i++)
3577  result.xelem (i) = Complex (re_val, im_val(i));
3578 
3579  retval = octave_value (new octave_complex_matrix (result));
3580  }
3581  }
3582  else
3583  {
3584  const NDArray re_val = re.array_value ();
3585 
3586  if (im.numel () == 1)
3587  {
3588  double im_val = im.double_value ();
3589 
3590  ComplexNDArray result (re_val.dims (), Complex ());
3591 
3592  for (octave_idx_type i = 0; i < re_val.numel (); i++)
3593  result.xelem (i) = Complex (re_val(i), im_val);
3594 
3595  retval = octave_value (new octave_complex_matrix (result));
3596  }
3597  else
3598  {
3599  const NDArray im_val = im.array_value ();
3600 
3601  if (re_val.dims () != im_val.dims ())
3602  error ("complex: dimension mismatch");
3603 
3604  ComplexNDArray result (re_val.dims (), Complex ());
3605 
3606  for (octave_idx_type i = 0; i < re_val.numel (); i++)
3607  result.xelem (i) = Complex (re_val(i), im_val(i));
3608 
3609  retval = octave_value (new octave_complex_matrix (result));
3610  }
3611  }
3612  }
3613 
3614  return retval;
3615 }
3616 
3617 DEFUN (isreal, args, ,
3618  doc: /* -*- texinfo -*-
3619 @deftypefn {} {} isreal (@var{x})
3620 Return true if @var{x} is a non-complex matrix or scalar.
3621 
3622 For compatibility with @sc{matlab}, this includes logical and character
3623 matrices.
3624 @seealso{iscomplex, isnumeric, isa}
3625 @end deftypefn */)
3626 {
3627  if (args.length () != 1)
3628  print_usage ();
3629 
3630  return ovl (args(0).is_real_type ());
3631 }
3632 
3633 DEFUN (isempty, args, ,
3634  doc: /* -*- texinfo -*-
3635 @deftypefn {} {} isempty (@var{a})
3636 Return true if @var{a} is an empty matrix (any one of its dimensions is
3637 zero).
3638 @seealso{isnull, isa}
3639 @end deftypefn */)
3640 {
3641  if (args.length () != 1)
3642  print_usage ();
3643 
3644  return ovl (args(0).is_empty ());
3645 }
3646 
3647 /*
3648 %% Debian bug #706376
3649 %!assert (isempty (speye(2^16)), false)
3650 */
3651 
3652 DEFUN (isnumeric, args, ,
3653  doc: /* -*- texinfo -*-
3654 @deftypefn {} {} isnumeric (@var{x})
3655 Return true if @var{x} is a numeric object, i.e., an integer, real, or
3656 complex array.
3657 
3658 Logical and character arrays are not considered to be numeric.
3659 @seealso{isinteger, isfloat, isreal, iscomplex, islogical, ischar, iscell, isstruct, isa}
3660 @end deftypefn */)
3661 {
3662  if (args.length () != 1)
3663  print_usage ();
3664 
3665  return ovl (args(0).is_numeric_type ());
3666 }
3667 
3668 /*
3669 %!assert (isnumeric (1), true)
3670 %!assert (isnumeric (1i), true)
3671 %!assert (isnumeric ([1,1]), true)
3672 %!assert (isnumeric (single (1)), true)
3673 %!assert (isnumeric (single (1i)), true)
3674 %!assert (isnumeric (single ([1,1])), true)
3675 %!assert (isnumeric (int8 (1)), true)
3676 %!assert (isnumeric (uint8 ([1,1])), true)
3677 %!assert (isnumeric ("Hello World"), false)
3678 %!assert (isnumeric (true), false)
3679 %!assert (isnumeric (false), false)
3680 %!assert (isnumeric ([true, false]), false)
3681 %!assert (isnumeric (sparse ([true, false])), false)
3682 */
3683 
3684 DEFUN (isscalar, args, ,
3685  doc: /* -*- texinfo -*-
3686 @deftypefn {} {} isscalar (@var{x})
3687 Return true if @var{x} is a scalar.
3688 @seealso{isvector, ismatrix}
3689 @end deftypefn */)
3690 {
3691  if (args.length () != 1)
3692  print_usage ();
3693 
3694  return ovl (args(0).numel () == 1);
3695 }
3696 
3697 /*
3698 %!assert (isscalar (1))
3699 %!assert (isscalar ([1, 2]), false)
3700 %!assert (isscalar ([]), false)
3701 %!assert (isscalar ([1, 2; 3, 4]), false)
3702 
3703 %!assert (isscalar ("t"))
3704 %!assert (isscalar ("test"), false)
3705 %!assert (isscalar (["test"; "ing"]), false)
3706 
3707 %!test
3708 %! s.a = 1;
3709 %! assert (isscalar (s));
3710 
3711 %% Test input validation
3712 %!error isscalar ()
3713 %!error isscalar (1, 2)
3714 */
3715 
3716 DEFUN (isvector, args, ,
3717  doc: /* -*- texinfo -*-
3718 @deftypefn {} {} isvector (@var{x})
3719 Return true if @var{x} is a vector.
3720 
3721 A vector is a 2-D array where one of the dimensions is equal to 1. As a
3722 consequence a 1x1 array, or scalar, is also a vector.
3723 @seealso{isscalar, ismatrix, size, rows, columns, length}
3724 @end deftypefn */)
3725 {
3726  if (args.length () != 1)
3727  print_usage ();
3728 
3729  dim_vector sz = args(0).dims ();
3730 
3731  return ovl (sz.ndims () == 2 && (sz(0) == 1 || sz(1) == 1));
3732 }
3733 
3734 /*
3735 %!assert (isvector (1), true)
3736 %!assert (isvector ([1; 2; 3]), true)
3737 %!assert (isvector ([1, 2, 3]), true)
3738 %!assert (isvector ([]), false)
3739 %!assert (isvector ([1, 2; 3, 4]), false)
3740 
3741 %!assert (isvector ("t"), true)
3742 %!assert (isvector ("test"), true)
3743 %!assert (isvector (["test"; "ing"]), false)
3744 
3745 %!test
3746 %! s.a = 1;
3747 %! assert (isvector (s), true);
3748 
3749 %% Test input validation
3750 %!error isvector ()
3751 %!error isvector ([1, 2], 2)
3752 */
3753 
3754 DEFUN (isrow, args, ,
3755  doc: /* -*- texinfo -*-
3756 @deftypefn {} {} isrow (@var{x})
3757 Return true if @var{x} is a row vector 1xN with non-negative N.
3758 @seealso{iscolumn, isscalar, isvector, ismatrix}
3759 @end deftypefn */)
3760 {
3761  if (args.length () != 1)
3762  print_usage ();
3763 
3764  dim_vector sz = args(0).dims ();
3765 
3766  return ovl (sz.ndims () == 2 && sz(0) == 1);
3767 }
3768 
3769 /*
3770 %!assert (isrow ([1, 2, 3]))
3771 %!assert (isrow ([1; 2; 3]), false)
3772 %!assert (isrow (1))
3773 %!assert (isrow ([]), false)
3774 %!assert (isrow ([1, 2; 3, 4]), false)
3775 
3776 %!assert (isrow (ones (1, 0)), true)
3777 %!assert (isrow (ones (1, 1)), true)
3778 %!assert (isrow (ones (1, 2)), true)
3779 %!assert (isrow (ones (1, 1, 1)), true)
3780 %!assert (isrow (ones (1, 1, 1, 1)), true)
3781 
3782 %!assert (isrow (ones (0, 0)), false)
3783 %!assert (isrow (ones (1, 1, 0)), false)
3784 
3785 %!assert (isrow ("t"), true)
3786 %!assert (isrow ("test"), true)
3787 %!assert (isrow (["test"; "ing"]), false)
3788 
3789 %!test
3790 %! s.a = 1;
3791 %! assert (isrow (s), true);
3792 
3793 %% Test input validation
3794 %!error isrow ()
3795 %!error isrow ([1, 2], 2)
3796 */
3797 
3798 DEFUN (iscolumn, args, ,
3799  doc: /* -*- texinfo -*-
3800 @deftypefn {} {} iscolumn (@var{x})
3801 Return true if @var{x} is a column vector Nx1 with non-negative N.
3802 @seealso{isrow, isscalar, isvector, ismatrix}
3803 @end deftypefn */)
3804 {
3805  if (args.length () != 1)
3806  print_usage ();
3807 
3808  dim_vector sz = args(0).dims ();
3809 
3810  return ovl (sz.ndims () == 2 && sz(1) == 1);
3811 }
3812 
3813 /*
3814 %!assert (iscolumn ([1, 2, 3]), false)
3815 %!assert (iscolumn ([1; 2; 3]), true)
3816 %!assert (iscolumn (1), true)
3817 %!assert (iscolumn ([]), false)
3818 %!assert (iscolumn ([1, 2; 3, 4]), false)
3819 
3820 %!assert (iscolumn ("t"), true)
3821 %!assert (iscolumn ("test"), false)
3822 %!assert (iscolumn (["test"; "ing"]), false)
3823 
3824 %!assert (iscolumn (ones (0, 1)), true)
3825 %!assert (iscolumn (ones (1, 1)), true)
3826 %!assert (iscolumn (ones (2, 1)), true)
3827 %!assert (iscolumn (ones (1, 1, 1)), true)
3828 %!assert (iscolumn (ones (1, 1, 1, 1)), true)
3829 
3830 %!assert (iscolumn (ones (0, 0)), false)
3831 %!assert (iscolumn (ones (0, 1, 0)), false)
3832 
3833 %!test
3834 %! s.a = 1;
3835 %! assert (iscolumn (s));
3836 
3837 %% Test input validation
3838 %!error iscolumn ()
3839 %!error iscolumn ([1, 2], 2)
3840 */
3841 
3842 DEFUN (ismatrix, args, ,
3843  doc: /* -*- texinfo -*-
3844 @deftypefn {} {} ismatrix (@var{a})
3845 Return true if @var{a} is a 2-D array.
3846 @seealso{isscalar, isvector, iscell, isstruct, issparse, isa}
3847 @end deftypefn */)
3848 {
3849  if (args.length () != 1)
3850  print_usage ();
3851 
3852  dim_vector sz = args(0).dims ();
3853 
3854  return ovl (sz.ndims () == 2 && sz(0) >= 0 && sz(1) >= 0);
3855 }
3856 
3857 /*
3858 %!assert (ismatrix ([]), true)
3859 %!assert (ismatrix (1), true)
3860 %!assert (ismatrix ([1, 2, 3]), true)
3861 %!assert (ismatrix ([1, 2; 3, 4]), true)
3862 
3863 %!assert (ismatrix (zeros (0)), true)
3864 %!assert (ismatrix (zeros (0, 0)), true)
3865 %!assert (ismatrix (zeros (0, 0, 0)), false)
3866 %!assert (ismatrix (zeros (3, 2, 4)), false)
3867 
3868 %!assert (ismatrix (single ([])), true)
3869 %!assert (ismatrix (single (1)), true)
3870 %!assert (ismatrix (single ([1, 2, 3])), true)
3871 %!assert (ismatrix (single ([1, 2; 3, 4])), true)
3872 
3873 %!assert (ismatrix ("t"), true)
3874 %!assert (ismatrix ("test"), true)
3875 %!assert (ismatrix (["test"; "ing"]), true)
3876 
3877 %!test
3878 %! s.a = 1;
3879 %! assert (ismatrix (s), true);
3880 
3881 %!error ismatrix ()
3882 %!error ismatrix ([1, 2; 3, 4], 2)
3883 */
3884 
3885 DEFUN (issquare, args, ,
3886  doc: /* -*- texinfo -*-
3887 @deftypefn {} {} issquare (@var{x})
3888 Return true if @var{x} is a square matrix.
3889 @seealso{isscalar, isvector, ismatrix, size}
3890 @end deftypefn */)
3891 {
3892  if (args.length () != 1)
3893  print_usage ();
3894 
3895  dim_vector sz = args(0).dims ();
3896 
3897  return ovl (sz.ndims () == 2 && sz(0) == sz(1));
3898 }
3899 
3900 /*
3901 %!assert (issquare ([]))
3902 %!assert (issquare (1))
3903 %!assert (! issquare ([1, 2]))
3904 %!assert (issquare ([1, 2; 3, 4]))
3905 %!assert (! issquare ([1, 2; 3, 4; 5, 6]))
3906 %!assert (! issquare (ones (3,3,3)))
3907 %!assert (issquare ("t"))
3908 %!assert (! issquare ("test"))
3909 %!assert (issquare (["test"; "ing"; "1"; "2"]))
3910 %!test
3911 %! s.a = 1;
3912 %! assert (issquare (s));
3913 %!assert (issquare ({1, 2; 3, 4}))
3914 %!assert (sparse (([1, 2; 3, 4])))
3915 
3916 %% Test input validation
3917 %!error issquare ()
3918 %!error issquare ([1, 2; 3, 4], 2)
3919 */
3920 
3921 static octave_value
3922 fill_matrix (const octave_value_list& args, int val, const char *fcn)
3923 {
3925 
3926  int nargin = args.length ();
3927 
3929 
3930  dim_vector dims (1, 1);
3931 
3932  if (nargin > 0 && args(nargin-1).is_string ())
3933  {
3934  std::string nm = args(nargin-1).string_value ();
3935  nargin--;
3936 
3938  }
3939 
3940  switch (nargin)
3941  {
3942  case 0:
3943  break;
3944 
3945  case 1:
3946  get_dimensions (args(0), fcn, dims);
3947  break;
3948 
3949  default:
3950  {
3951  dims.resize (nargin);
3952 
3953  for (int i = 0; i < nargin; i++)
3954  dims(i) = (args(i).is_empty ()
3955  ? 0 : args(i).xidx_type_value ("%s: dimension arguments must be scalar integers", fcn));
3956  }
3957  break;
3958  }
3959 
3960  dims.chop_trailing_singletons ();
3961 
3962  check_dimensions (dims, fcn);
3963 
3964  // FIXME: perhaps this should be made extensible by
3965  // using the class name to lookup a function to call to create
3966  // the new value.
3967 
3968  // Note that automatic narrowing will handle conversion from
3969  // NDArray to scalar.
3970 
3971  switch (dt)
3972  {
3974  retval = int8NDArray (dims, val);
3975  break;
3976 
3978  retval = uint8NDArray (dims, val);
3979  break;
3980 
3982  retval = int16NDArray (dims, val);
3983  break;
3984 
3986  retval = uint16NDArray (dims, val);
3987  break;
3988 
3990  retval = int32NDArray (dims, val);
3991  break;
3992 
3994  retval = uint32NDArray (dims, val);
3995  break;
3996 
3998  retval = int64NDArray (dims, val);
3999  break;
4000 
4002  retval = uint64NDArray (dims, val);
4003  break;
4004 
4006  retval = FloatNDArray (dims, val);
4007  break;
4008 
4010  {
4011  if (val == 1 && dims.ndims () == 2 && dims(0) == 1)
4012  retval = Range (1.0, 0.0, dims(1)); // packed form
4013  else
4014  retval = NDArray (dims, val);
4015  }
4016  break;
4017 
4019  retval = boolNDArray (dims, val);
4020  break;
4021 
4022  default:
4023  error ("%s: invalid class name", fcn);
4024  break;
4025  }
4026 
4027  return retval;
4028 }
4029 
4030 static octave_value
4031 fill_matrix (const octave_value_list& args, double val, float fval,
4032  const char *fcn)
4033 {
4035 
4036  int nargin = args.length ();
4037 
4039 
4040  dim_vector dims (1, 1);
4041 
4042  if (nargin > 0 && args(nargin-1).is_string ())
4043  {
4044  std::string nm = args(nargin-1).string_value ();
4045  nargin--;
4046 
4048  }
4049 
4050  switch (nargin)
4051  {
4052  case 0:
4053  break;
4054 
4055  case 1:
4056  get_dimensions (args(0), fcn, dims);
4057  break;
4058 
4059  default:
4060  {
4061  dims.resize (nargin);
4062 
4063  for (int i = 0; i < nargin; i++)
4064  dims(i) = (args(i).is_empty ()
4065  ? 0 : args(i).xidx_type_value ("%s: dimension arguments must be scalar integers", fcn));
4066  }
4067  break;
4068  }
4069 
4070  dims.chop_trailing_singletons ();
4071 
4072  check_dimensions (dims, fcn);
4073 
4074  // Note that automatic narrowing will handle conversion from
4075  // NDArray to scalar.
4076 
4077  switch (dt)
4078  {
4080  retval = FloatNDArray (dims, fval);
4081  break;
4082 
4084  retval = NDArray (dims, val);
4085  break;
4086 
4087  default:
4088  error ("%s: invalid class name", fcn);
4089  break;
4090  }
4091 
4092  return retval;
4093 }
4094 
4095 static octave_value
4096 fill_matrix (const octave_value_list& args, double val, const char *fcn)
4097 {
4099 
4100  int nargin = args.length ();
4101 
4103 
4104  dim_vector dims (1, 1);
4105 
4106  if (nargin > 0 && args(nargin-1).is_string ())
4107  {
4108  std::string nm = args(nargin-1).string_value ();
4109  nargin--;
4110 
4112  }
4113 
4114  switch (nargin)
4115  {
4116  case 0:
4117  break;
4118 
4119  case 1:
4120  get_dimensions (args(0), fcn, dims);
4121  break;
4122 
4123  default:
4124  {
4125  dims.resize (nargin);
4126 
4127  for (int i = 0; i < nargin; i++)
4128  dims(i) = (args(i).is_empty ()
4129  ? 0 : args(i).xidx_type_value ("%s: dimension arguments must be scalar integers", fcn));
4130  }
4131  break;
4132  }
4133 
4134  dims.chop_trailing_singletons ();
4135 
4136  check_dimensions (dims, fcn);
4137 
4138  // Note that automatic narrowing will handle conversion from
4139  // NDArray to scalar.
4140 
4141  switch (dt)
4142  {
4144  retval = FloatNDArray (dims, static_cast<float> (val));
4145  break;
4146 
4148  retval = NDArray (dims, val);
4149  break;
4150 
4151  default:
4152  error ("%s: invalid class name", fcn);
4153  break;
4154  }
4155 
4156  return retval;
4157 }
4158 
4159 static octave_value
4160 fill_matrix (const octave_value_list& args, const Complex& val,
4161  const char *fcn)
4162 {
4164 
4165  int nargin = args.length ();
4166 
4168 
4169  dim_vector dims (1, 1);
4170 
4171  if (nargin > 0 && args(nargin-1).is_string ())
4172  {
4173  std::string nm = args(nargin-1).string_value ();
4174  nargin--;
4175 
4177  }
4178 
4179  switch (nargin)
4180  {
4181  case 0:
4182  break;
4183 
4184  case 1:
4185  get_dimensions (args(0), fcn, dims);
4186  break;
4187 
4188  default:
4189  {
4190  dims.resize (nargin);
4191 
4192  for (int i = 0; i < nargin; i++)
4193  dims(i) = (args(i).is_empty ()
4194  ? 0 : args(i).xidx_type_value ("%s: dimension arguments must be scalar integers", fcn));
4195  }
4196  break;
4197  }
4198 
4199  dims.chop_trailing_singletons ();
4200 
4201  check_dimensions (dims, fcn);
4202 
4203  // Note that automatic narrowing will handle conversion from
4204  // NDArray to scalar.
4205 
4206  switch (dt)
4207  {
4209  retval = FloatComplexNDArray (dims,
4210  static_cast<FloatComplex> (val));
4211  break;
4212 
4214  retval = ComplexNDArray (dims, val);
4215  break;
4216 
4217  default:
4218  error ("%s: invalid class name", fcn);
4219  break;
4220  }
4221 
4222  return retval;
4223 }
4224 
4225 static octave_value
4226 fill_matrix (const octave_value_list& args, bool val, const char *fcn)
4227 {
4229 
4230  int nargin = args.length ();
4231 
4232  dim_vector dims (1, 1);
4233 
4234  switch (nargin)
4235  {
4236  case 0:
4237  break;
4238 
4239  case 1:
4240  get_dimensions (args(0), fcn, dims);
4241  break;
4242 
4243  default:
4244  {
4245  dims.resize (nargin);
4246 
4247  for (int i = 0; i < nargin; i++)
4248  dims(i) = (args(i).is_empty ()
4249  ? 0 : args(i).xidx_type_value ("%s: dimension arguments must be scalar integers", fcn));
4250  }
4251  break;
4252  }
4253 
4254  dims.chop_trailing_singletons ();
4255 
4256  check_dimensions (dims, fcn);
4257 
4258  // Note that automatic narrowing will handle conversion from
4259  // NDArray to scalar.
4260 
4261  retval = boolNDArray (dims, val);
4262 
4263  return retval;
4264 }
4265 
4266 DEFUN (ones, args, ,
4267  doc: /* -*- texinfo -*-
4268 @deftypefn {} {} ones (@var{n})
4269 @deftypefnx {} {} ones (@var{m}, @var{n})
4270 @deftypefnx {} {} ones (@var{m}, @var{n}, @var{k}, @dots{})
4271 @deftypefnx {} {} ones ([@var{m} @var{n} @dots{}])
4272 @deftypefnx {} {} ones (@dots{}, @var{class})
4273 Return a matrix or N-dimensional array whose elements are all 1.
4274 
4275 If invoked with a single scalar integer argument @var{n}, return a square
4276 @nospell{NxN} matrix.
4277 
4278 If invoked with two or more scalar integer arguments, or a vector of integer
4279 values, return an array with the given dimensions.
4280 
4281 To create a constant matrix whose values are all the same use an expression
4282 such as
4283 
4284 @example
4285 val_matrix = val * ones (m, n)
4286 @end example
4287 
4288 The optional argument @var{class} specifies the class of the return array
4289 and defaults to double. For example:
4290 
4291 @example
4292 val = ones (m,n, "uint8")
4293 @end example
4294 @seealso{zeros}
4295 @end deftypefn */)
4296 {
4297  return fill_matrix (args, 1, "ones");
4298 }
4299 
4300 /*
4301 %!assert (ones (3), [1, 1, 1; 1, 1, 1; 1, 1, 1])
4302 %!assert (ones (2, 3), [1, 1, 1; 1, 1, 1])
4303 %!assert (ones (3, 2), [1, 1; 1, 1; 1, 1])
4304 %!assert (size (ones (3, 4, 5)), [3, 4, 5])
4305 
4306 %!assert (ones (3, "single"), single ([1, 1, 1; 1, 1, 1; 1, 1, 1]))
4307 %!assert (ones (2, 3, "single"), single ([1, 1, 1; 1, 1, 1]))
4308 %!assert (ones (3, 2, "single"), single ([1, 1; 1, 1; 1, 1]))
4309 %!assert (size (ones (3, 4, 5, "single")), [3, 4, 5])
4310 
4311 %!assert (ones (3, "int8"), int8 ([1, 1, 1; 1, 1, 1; 1, 1, 1]))
4312 %!assert (ones (2, 3, "int8"), int8 ([1, 1, 1; 1, 1, 1]))
4313 %!assert (ones (3, 2, "int8"), int8 ([1, 1; 1, 1; 1, 1]))
4314 %!assert (size (ones (3, 4, 5, "int8")), [3, 4, 5])
4315 */
4316 
4317 /*
4318 ## Tests for bug #47298
4319 ## Matlab requires the size to be a row vector. In that logic, it supports
4320 ## n to be a 1x0 vector (returns 0x0) but not a 0x1 vector. Octave supports
4321 ## any vector and therefore must support 0x1, 1x0, and 0x0x1 (but not 0x1x1).
4322 %!test <47298>
4323 %! funcs = {@zeros, @ones, @inf, @nan, @NA, @i, @pi, @e};
4324 %! for idx = 1:numel (funcs)
4325 %! func = funcs{idx};
4326 %! assert (func (zeros (1, 0)), zeros (0, 0));
4327 %! assert (func (zeros (0, 1)), zeros (0, 0));
4328 %! assert (func (zeros (0, 1, 1)), zeros (0, 0));
4329 %! fail ([func2str(func) " ([])"]);
4330 %! fail ([func2str(func) " (zeros (0, 0, 1))"]);
4331 %! endfor
4332 */
4333 
4334 DEFUN (zeros, args, ,
4335  doc: /* -*- texinfo -*-
4336 @deftypefn {} {} zeros (@var{n})
4337 @deftypefnx {} {} zeros (@var{m}, @var{n})
4338 @deftypefnx {} {} zeros (@var{m}, @var{n}, @var{k}, @dots{})
4339 @deftypefnx {} {} zeros ([@var{m} @var{n} @dots{}])
4340 @deftypefnx {} {} zeros (@dots{}, @var{class})
4341 Return a matrix or N-dimensional array whose elements are all 0.
4342 
4343 If invoked with a single scalar integer argument, return a square
4344 @nospell{NxN} matrix.
4345 
4346 If invoked with two or more scalar integer arguments, or a vector of integer
4347 values, return an array with the given dimensions.
4348 
4349 The optional argument @var{class} specifies the class of the return array
4350 and defaults to double. For example:
4351 
4352 @example
4353 val = zeros (m,n, "uint8")
4354 @end example
4355 @seealso{ones}
4356 @end deftypefn */)
4357 {
4358  return fill_matrix (args, 0, "zeros");
4359 }
4360 
4361 /*
4362 %!assert (zeros (3), [0, 0, 0; 0, 0, 0; 0, 0, 0])
4363 %!assert (zeros (2, 3), [0, 0, 0; 0, 0, 0])
4364 %!assert (zeros (3, 2), [0, 0; 0, 0; 0, 0])
4365 %!assert (size (zeros (3, 4, 5)), [3, 4, 5])
4366 
4367 %!assert (zeros (3, "single"), single ([0, 0, 0; 0, 0, 0; 0, 0, 0]))
4368 %!assert (zeros (2, 3, "single"), single ([0, 0, 0; 0, 0, 0]))
4369 %!assert (zeros (3, 2, "single"), single ([0, 0; 0, 0; 0, 0]))
4370 %!assert (size (zeros (3, 4, 5, "single")), [3, 4, 5])
4371 
4372 %!assert (zeros (3, "int8"), int8 ([0, 0, 0; 0, 0, 0; 0, 0, 0]))
4373 %!assert (zeros (2, 3, "int8"), int8 ([0, 0, 0; 0, 0, 0]))
4374 %!assert (zeros (3, 2, "int8"), int8 ([0, 0; 0, 0; 0, 0]))
4375 %!assert (size (zeros (3, 4, 5, "int8")), [3, 4, 5])
4376 */
4377 
4378 DEFUN (Inf, args, ,
4379  doc: /* -*- texinfo -*-
4380 @c List other form of function in documentation index
4381 @findex inf
4382 
4383 @deftypefn {} {} Inf
4384 @deftypefnx {} {} Inf (@var{n})
4385 @deftypefnx {} {} Inf (@var{n}, @var{m})
4386 @deftypefnx {} {} Inf (@var{n}, @var{m}, @var{k}, @dots{})
4387 @deftypefnx {} {} Inf (@dots{}, @var{class})
4388 Return a scalar, matrix or N-dimensional array whose elements are all equal
4389 to the IEEE representation for positive infinity.
4390 
4391 Infinity is produced when results are too large to be represented using the
4392 IEEE floating point format for numbers. Two common examples which produce
4393 infinity are division by zero and overflow.
4394 
4395 @example
4396 @group
4397 [ 1/0 e^800 ]
4398 @result{} Inf Inf
4399 @end group
4400 @end example
4401 
4402 When called with no arguments, return a scalar with the value @samp{Inf}.
4403 
4404 When called with a single argument, return a square matrix with the
4405 dimension specified.
4406 
4407 When called with more than one scalar argument the first two arguments are
4408 taken as the number of rows and columns and any further arguments specify
4409 additional matrix dimensions.
4410 
4411 The optional argument @var{class} specifies the return type and may be
4412 either @qcode{"double"} or @qcode{"single"}.
4413 @seealso{isinf, NaN}
4414 @end deftypefn */)
4415 {
4416  return fill_matrix (args, lo_ieee_inf_value (),
4417  lo_ieee_float_inf_value (), "Inf");
4418 }
4419 
4420 DEFALIAS (inf, Inf);
4421 
4422 /*
4423 %!assert (inf (3), [Inf, Inf, Inf; Inf, Inf, Inf; Inf, Inf, Inf])
4424 %!assert (inf (2, 3), [Inf, Inf, Inf; Inf, Inf, Inf])
4425 %!assert (inf (3, 2), [Inf, Inf; Inf, Inf; Inf, Inf])
4426 %!assert (size (inf (3, 4, 5)), [3, 4, 5])
4427 
4428 %!assert (inf (3, "single"), single ([Inf, Inf, Inf; Inf, Inf, Inf; Inf, Inf, Inf]))
4429 %!assert (inf (2, 3, "single"), single ([Inf, Inf, Inf; Inf, Inf, Inf]))
4430 %!assert (inf (3, 2, "single"), single ([Inf, Inf; Inf, Inf; Inf, Inf]))
4431 %!assert (size (inf (3, 4, 5, "single")), [3, 4, 5])
4432 
4433 %!error (inf (3, "int8"))
4434 %!error (inf (2, 3, "int8"))
4435 %!error (inf (3, 2, "int8"))
4436 %!error (inf (3, 4, 5, "int8"))
4437 */
4438 
4439 DEFUN (NaN, args, ,
4440  doc: /* -*- texinfo -*-
4441 @c List other form of function in documentation index
4442 @findex nan
4443 
4444 @deftypefn {} {} NaN
4445 @deftypefnx {} {} NaN (@var{n})
4446 @deftypefnx {} {} NaN (@var{n}, @var{m})
4447 @deftypefnx {} {} NaN (@var{n}, @var{m}, @var{k}, @dots{})
4448 @deftypefnx {} {} NaN (@dots{}, @var{class})
4449 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4450 to the IEEE symbol NaN (Not a Number).
4451 
4452 NaN is the result of operations which do not produce a well defined
4453 numerical result. Common operations which produce a NaN are arithmetic
4454 with infinity
4455 @tex
4456 ($\infty - \infty$), zero divided by zero ($0/0$),
4457 @end tex
4458 @ifnottex
4459 (Inf - Inf), zero divided by zero (0/0),
4460 @end ifnottex
4461 and any operation involving another NaN value (5 + NaN).
4462 
4463 Note that NaN always compares not equal to NaN (NaN != NaN). This behavior
4464 is specified by the IEEE standard for floating point arithmetic. To find
4465 NaN values, use the @code{isnan} function.
4466 
4467 When called with no arguments, return a scalar with the value @samp{NaN}.
4468 
4469 When called with a single argument, return a square matrix with the
4470 dimension specified.
4471 
4472 When called with more than one scalar argument the first two arguments are
4473 taken as the number of rows and columns and any further arguments specify
4474 additional matrix dimensions.
4475 
4476 The optional argument @var{class} specifies the return type and may be
4477 either @qcode{"double"} or @qcode{"single"}.
4478 @seealso{isnan, Inf}
4479 @end deftypefn */)
4480 {
4481  return fill_matrix (args, lo_ieee_nan_value (),
4482  lo_ieee_float_nan_value (), "NaN");
4483 }
4484 
4485 DEFALIAS (nan, NaN);
4486 
4487 /*
4488 %!assert (NaN (3), [NaN, NaN, NaN; NaN, NaN, NaN; NaN, NaN, NaN])
4489 %!assert (NaN (2, 3), [NaN, NaN, NaN; NaN, NaN, NaN])
4490 %!assert (NaN (3, 2), [NaN, NaN; NaN, NaN; NaN, NaN])
4491 %!assert (size (NaN (3, 4, 5)), [3, 4, 5])
4492 
4493 %!assert (NaN (3, "single"), single ([NaN, NaN, NaN; NaN, NaN, NaN; NaN, NaN, NaN]))
4494 %!assert (NaN (2, 3, "single"), single ([NaN, NaN, NaN; NaN, NaN, NaN]))
4495 %!assert (NaN (3, 2, "single"), single ([NaN, NaN; NaN, NaN; NaN, NaN]))
4496 %!assert (size (NaN (3, 4, 5, "single")), [3, 4, 5])
4497 
4498 %!error (NaN (3, "int8"))
4499 %!error (NaN (2, 3, "int8"))
4500 %!error (NaN (3, 2, "int8"))
4501 %!error (NaN (3, 4, 5, "int8"))
4502 */
4503 
4504 DEFUN (e, args, ,
4505  doc: /* -*- texinfo -*-
4506 @deftypefn {} {} e
4507 @deftypefnx {} {} e (@var{n})
4508 @deftypefnx {} {} e (@var{n}, @var{m})
4509 @deftypefnx {} {} e (@var{n}, @var{m}, @var{k}, @dots{})
4510 @deftypefnx {} {} e (@dots{}, @var{class})
4511 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4512 to the base of natural logarithms.
4513 
4514 The constant
4515 @tex
4516 $e$ satisfies the equation $\log (e) = 1$.
4517 @end tex
4518 @ifnottex
4519 @samp{e} satisfies the equation @code{log} (e) = 1.
4520 @end ifnottex
4521 
4522 When called with no arguments, return a scalar with the value @math{e}.
4523 
4524 When called with a single argument, return a square matrix with the
4525 dimension specified.
4526 
4527 When called with more than one scalar argument the first two arguments are
4528 taken as the number of rows and columns and any further arguments specify
4529 additional matrix dimensions.
4530 
4531 The optional argument @var{class} specifies the return type and may be
4532 either @qcode{"double"} or @qcode{"single"}.
4533 @seealso{log, exp, pi, I}
4534 @end deftypefn */)
4535 {
4536 #if defined (M_E)
4537  double e_val = M_E;
4538 #else
4539  double e_val = exp (1.0);
4540 #endif
4541 
4542  return fill_matrix (args, e_val, "e");
4543 }
4544 
4545 DEFUN (eps, args, ,
4546  doc: /* -*- texinfo -*-
4547 @deftypefn {} {} eps
4548 @deftypefnx {} {} eps (@var{x})
4549 @deftypefnx {} {} eps (@var{n}, @var{m})
4550 @deftypefnx {} {} eps (@var{n}, @var{m}, @var{k}, @dots{})
4551 @deftypefnx {} {} eps (@dots{}, @var{class})
4552 Return a scalar, matrix or N-dimensional array whose elements are all eps,
4553 the machine precision.
4554 
4555 More precisely, @code{eps} is the relative spacing between any two adjacent
4556 numbers in the machine's floating point system. This number is obviously
4557 system dependent. On machines that support IEEE floating point arithmetic,
4558 @code{eps} is approximately
4559 @tex
4560 $2.2204\times10^{-16}$ for double precision and $1.1921\times10^{-7}$
4561 @end tex
4562 @ifnottex
4563 2.2204e-16 for double precision and 1.1921e-07
4564 @end ifnottex
4565 for single precision.
4566 
4567 When called with no arguments, return a scalar with the value
4568 @code{eps (1.0)}.
4569 
4570 Given a single argument @var{x}, return the distance between @var{x} and the
4571 next largest value.
4572 
4573 When called with more than one argument the first two arguments are taken as
4574 the number of rows and columns and any further arguments specify additional
4575 matrix dimensions. The optional argument @var{class} specifies the return
4576 type and may be either @qcode{"double"} or @qcode{"single"}.
4577 @seealso{realmax, realmin, intmax, flintmax}
4578 @end deftypefn */)
4579 {
4581 
4582  if (args.length () == 1 && ! args(0).is_string ())
4583  {
4584  if (args(0).is_single_type ())
4585  {
4586  Array<float> x = args(0).float_array_value ();
4587 
4588  Array<float> epsval (x.dims ());
4589 
4590  for (octave_idx_type i = 0; i < x.numel (); i++)
4591  {
4592  float val = ::fabsf (x(i));
4593  if (octave::math::isnan (val) || octave::math::isinf (val))
4594  epsval(i) = lo_ieee_nan_value ();
4595  else if (val < std::numeric_limits<float>::min ())
4596  epsval(i) = powf (2.0, -149e0);
4597  else
4598  {
4599  int expon;
4600  octave::math::frexp (val, &expon);
4601  epsval(i) = std::pow (2.0f,
4602  static_cast<float> (expon - 24));
4603  }
4604  }
4605 
4606  retval = epsval;
4607  }
4608  else
4609  {
4610  Array<double> x = args(0).array_value ();
4611 
4612  Array<double> epsval (x.dims ());
4613 
4614  for (octave_idx_type i = 0; i < x.numel (); i++)
4615  {
4616  double val = ::fabs (x(i));
4617  if (octave::math::isnan (val) || octave::math::isinf (val))
4618  epsval(i) = lo_ieee_nan_value ();
4619  else if (val < std::numeric_limits<double>::min ())
4620  epsval(i) = pow (2.0, -1074e0);
4621  else
4622  {
4623  int expon;
4624  octave::math::frexp (val, &expon);
4625  epsval(i) = std::pow (2.0,
4626  static_cast<double> (expon - 53));
4627  }
4628 
4629  retval = epsval;
4630  }
4631  }
4632  }
4633  else
4634  retval = fill_matrix (args, std::numeric_limits<double>::epsilon (),
4635  std::numeric_limits<float>::epsilon (), "eps");
4636 
4637  return retval;
4638 }
4639 
4640 /*
4641 %!assert (eps (1/2), 2^(-53))
4642 %!assert (eps (1), 2^(-52))
4643 %!assert (eps (2), 2^(-51))
4644 %!assert (eps (realmax), 2^971)
4645 %!assert (eps (0), 2^(-1074))
4646 %!assert (eps (realmin/2), 2^(-1074))
4647 %!assert (eps (realmin/16), 2^(-1074))
4648 %!assert (eps (Inf), NaN)
4649 %!assert (eps (NaN), NaN)
4650 %!assert (eps ([1/2 1 2 realmax 0 realmin/2 realmin/16 Inf NaN]),
4651 %! [2^(-53) 2^(-52) 2^(-51) 2^971 2^(-1074) 2^(-1074) 2^(-1074) NaN NaN])
4652 %!assert (eps (single (1/2)), single (2^(-24)))
4653 %!assert (eps (single (1)), single (2^(-23)))
4654 %!assert (eps (single (2)), single (2^(-22)))
4655 %!assert (eps (realmax ("single")), single (2^104))
4656 %!assert (eps (single (0)), single (2^(-149)))
4657 %!assert (eps (realmin ("single")/2), single (2^(-149)))
4658 %!assert (eps (realmin ("single")/16), single (2^(-149)))
4659 %!assert (eps (single (Inf)), single (NaN))
4660 %!assert (eps (single (NaN)), single (NaN))
4661 %!assert (eps (single ([1/2 1 2 realmax("single") 0 realmin("single")/2 realmin("single")/16 Inf NaN])),
4662 %! single ([2^(-24) 2^(-23) 2^(-22) 2^104 2^(-149) 2^(-149) 2^(-149) NaN NaN]))
4663 
4664 */
4665 
4666 DEFUN (pi, args, ,
4667  doc: /* -*- texinfo -*-
4668 @deftypefn {} {} pi
4669 @deftypefnx {} {} pi (@var{n})
4670 @deftypefnx {} {} pi (@var{n}, @var{m})
4671 @deftypefnx {} {} pi (@var{n}, @var{m}, @var{k}, @dots{})
4672 @deftypefnx {} {} pi (@dots{}, @var{class})
4673 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4674 to the ratio of the circumference of a circle to its
4675 @tex
4676 diameter($\pi$).
4677 @end tex
4678 @ifnottex
4679 diameter.
4680 @end ifnottex
4681 
4682 Internally, @code{pi} is computed as @samp{4.0 * atan (1.0)}.
4683 
4684 When called with no arguments, return a scalar with the value of
4685 @tex
4686 $\pi$.
4687 @end tex
4688 @ifnottex
4689 pi.
4690 @end ifnottex
4691 
4692 When called with a single argument, return a square matrix with the
4693 dimension specified.
4694 
4695 When called with more than one scalar argument the first two arguments are
4696 taken as the number of rows and columns and any further arguments specify
4697 additional matrix dimensions.
4698 
4699 The optional argument @var{class} specifies the return type and may be
4700 either @qcode{"double"} or @qcode{"single"}.
4701 @seealso{e, I}
4702 @end deftypefn */)
4703 {
4704 #if defined (M_PI)
4705  double pi_val = M_PI;
4706 #else
4707  double pi_val = 4.0 * atan (1.0);
4708 #endif
4709 
4710  return fill_matrix (args, pi_val, "pi");
4711 }
4712 
4713 DEFUN (realmax, args, ,
4714  doc: /* -*- texinfo -*-
4715 @deftypefn {} {} realmax
4716 @deftypefnx {} {} realmax (@var{n})
4717 @deftypefnx {} {} realmax (@var{n}, @var{m})
4718 @deftypefnx {} {} realmax (@var{n}, @var{m}, @var{k}, @dots{})
4719 @deftypefnx {} {} realmax (@dots{}, @var{class})
4720 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4721 to the largest floating point number that is representable.
4722 
4723 The actual value is system dependent. On machines that support IEEE
4724 floating point arithmetic, @code{realmax} is approximately
4725 @tex
4726 $1.7977\times10^{308}$ for double precision and $3.4028\times10^{38}$
4727 @end tex
4728 @ifnottex
4729 1.7977e+308 for double precision and 3.4028e+38
4730 @end ifnottex
4731 for single precision.
4732 
4733 When called with no arguments, return a scalar with the value
4734 @code{realmax (@qcode{"double"})}.
4735 
4736 When called with a single argument, return a square matrix with the
4737 dimension specified.
4738 
4739 When called with more than one scalar argument the first two arguments are
4740 taken as the number of rows and columns and any further arguments specify
4741 additional matrix dimensions.
4742 
4743 The optional argument @var{class} specifies the return type and may be
4744 either @qcode{"double"} or @qcode{"single"}.
4745 @seealso{realmin, intmax, flintmax, eps}
4746 @end deftypefn */)
4747 {
4748  return fill_matrix (args, std::numeric_limits<double>::max (),
4749  std::numeric_limits<float>::max (), "realmax");
4750 }
4751 
4752 DEFUN (realmin, args, ,
4753  doc: /* -*- texinfo -*-
4754 @deftypefn {} {} realmin
4755 @deftypefnx {} {} realmin (@var{n})
4756 @deftypefnx {} {} realmin (@var{n}, @var{m})
4757 @deftypefnx {} {} realmin (@var{n}, @var{m}, @var{k}, @dots{})
4758 @deftypefnx {} {} realmin (@dots{}, @var{class})
4759 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4760 to the smallest normalized floating point number that is representable.
4761 
4762 The actual value is system dependent. On machines that support
4763 IEEE floating point arithmetic, @code{realmin} is approximately
4764 @tex
4765 $2.2251\times10^{-308}$ for double precision and $1.1755\times10^{-38}$
4766 @end tex
4767 @ifnottex
4768 2.2251e-308 for double precision and 1.1755e-38
4769 @end ifnottex
4770 for single precision.
4771 
4772 When called with no arguments, return a scalar with the value
4773 @code{realmin (@qcode{"double"})}.
4774 
4775 When called with a single argument, return a square matrix with the
4776 dimension specified.
4777 
4778 When called with more than one scalar argument the first two arguments are
4779 taken as the number of rows and columns and any further arguments specify
4780 additional matrix dimensions.
4781 
4782 The optional argument @var{class} specifies the return type and may be
4783 either @qcode{"double"} or @qcode{"single"}.
4784 @seealso{realmax, intmin, eps}
4785 @end deftypefn */)
4786 {
4787  return fill_matrix (args, std::numeric_limits<double>::min (),
4788  std::numeric_limits<float>::min (), "realmin");
4789 }
4790 
4791 DEFUN (I, args, ,
4792  doc: /* -*- texinfo -*-
4793 @c List other forms of function in documentation index
4794 @findex i
4795 @findex j
4796 @findex J
4797 
4798 @deftypefn {} {} I
4799 @deftypefnx {} {} I (@var{n})
4800 @deftypefnx {} {} I (@var{n}, @var{m})
4801 @deftypefnx {} {} I (@var{n}, @var{m}, @var{k}, @dots{})
4802 @deftypefnx {} {} I (@dots{}, @var{class})
4803 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4804 to the pure imaginary unit, defined as
4805 @tex
4806 $\sqrt{-1}$.
4807 @end tex
4808 @ifnottex
4809 @w{@code{sqrt (-1)}}.
4810 @end ifnottex
4811 
4812 I, and its equivalents i, j, and J, are functions so any of the names may
4813 be reused for other purposes (such as i for a counter variable).
4814 
4815 When called with no arguments, return a scalar with the value @math{i}.
4816 
4817 When called with a single argument, return a square matrix with the
4818 dimension specified.
4819 
4820 When called with more than one scalar argument the first two arguments are
4821 taken as the number of rows and columns and any further arguments specify
4822 additional matrix dimensions.
4823 
4824 The optional argument @var{class} specifies the return type and may be
4825 either @qcode{"double"} or @qcode{"single"}.
4826 @seealso{e, pi, log, exp}
4827 @end deftypefn */)
4828 {
4829  return fill_matrix (args, Complex (0.0, 1.0), "I");
4830 }
4831 
4832 DEFALIAS (i, I);
4833 DEFALIAS (J, I);
4834 DEFALIAS (j, I);
4835 
4836 DEFUN (NA, args, ,
4837  doc: /* -*- texinfo -*-
4838 @deftypefn {} {} NA
4839 @deftypefnx {} {} NA (@var{n})
4840 @deftypefnx {} {} NA (@var{n}, @var{m})
4841 @deftypefnx {} {} NA (@var{n}, @var{m}, @var{k}, @dots{})
4842 @deftypefnx {} {} NA (@dots{}, @var{class})
4843 Return a scalar, matrix, or N-dimensional array whose elements are all equal
4844 to the special constant used to designate missing values.
4845 
4846 Note that NA always compares not equal to NA (NA != NA).
4847 To find NA values, use the @code{isna} function.
4848 
4849 When called with no arguments, return a scalar with the value @samp{NA}.
4850 
4851 When called with a single argument, return a square matrix with the
4852 dimension specified.
4853 
4854 When called with more than one scalar argument the first two arguments are
4855 taken as the number of rows and columns and any further arguments specify
4856 additional matrix dimensions.
4857 
4858 The optional argument @var{class} specifies the return type and may be
4859 either @qcode{"double"} or @qcode{"single"}.
4860 @seealso{isna}
4861 @end deftypefn */)
4862 {
4863  return fill_matrix (args, lo_ieee_na_value (),
4864  lo_ieee_float_na_value (), "NA");
4865 }
4866 
4867 /*
4868 %!assert (single (NA ("double")), NA ("single"))
4869 %!assert (double (NA ("single")), NA ("double"))
4870 */
4871 
4872 DEFUN (false, args, ,
4873  doc: /* -*- texinfo -*-
4874 @deftypefn {} {} false (@var{x})
4875 @deftypefnx {} {} false (@var{n}, @var{m})
4876 @deftypefnx {} {} false (@var{n}, @var{m}, @var{k}, @dots{})
4877 Return a matrix or N-dimensional array whose elements are all logical 0.
4878 
4879 If invoked with a single scalar integer argument, return a square
4880 matrix of the specified size.
4881 
4882 If invoked with two or more scalar integer arguments, or a vector of integer
4883 values, return an array with given dimensions.
4884 @seealso{true}
4885 @end deftypefn */)
4886 {
4887  return fill_matrix (args, false, "false");
4888 }
4889 
4890 DEFUN (true, args, ,
4891  doc: /* -*- texinfo -*-
4892 @deftypefn {} {} true (@var{x})
4893 @deftypefnx {} {} true (@var{n}, @var{m})
4894 @deftypefnx {} {} true (@var{n}, @var{m}, @var{k}, @dots{})
4895 Return a matrix or N-dimensional array whose elements are all logical 1.
4896 
4897 If invoked with a single scalar integer argument, return a square
4898 matrix of the specified size.
4899 
4900 If invoked with two or more scalar integer arguments, or a vector of integer
4901 values, return an array with given dimensions.
4902 @seealso{false}
4903 @end deftypefn */)
4904 {
4905  return fill_matrix (args, true, "true");
4906 }
4907 
4908 template <typename MT>
4910 identity_matrix (int nr, int nc)
4911 {
4913 
4914  typename MT::element_type one (1);
4915 
4916  if (nr == 1 && nc == 1)
4917  retval = one;
4918  else
4919  {
4920  dim_vector dims (nr, nc);
4921 
4922  typename MT::element_type zero (0);
4923 
4924  MT m (dims, zero);
4925 
4926  if (nr > 0 && nc > 0)
4927  {
4928  int n = std::min (nr, nc);
4929 
4930  for (int i = 0; i < n; i++)
4931  m(i,i) = one;
4932  }
4933 
4934  retval = m;
4935  }
4936 
4937  return retval;
4938 }
4939 
4940 #define INSTANTIATE_EYE(T) \
4941  template octave_value identity_matrix<T> (int, int)
4942 
4954 
4955 static octave_value
4957 {
4959 
4960  // FIXME: perhaps this should be made extensible by using
4961  // the class name to lookup a function to call to create the new
4962  // value.
4963 
4964  switch (dt)
4965  {
4967  retval = identity_matrix<int8NDArray> (nr, nc);
4968  break;
4969 
4971  retval = identity_matrix<uint8NDArray> (nr, nc);
4972  break;
4973 
4975  retval = identity_matrix<int16NDArray> (nr, nc);
4976  break;
4977 
4979  retval = identity_matrix<uint16NDArray> (nr, nc);
4980  break;
4981 
4983  retval = identity_matrix<int32NDArray> (nr, nc);
4984  break;
4985 
4987  retval = identity_matrix<uint32NDArray> (nr, nc);
4988  break;
4989 
4991  retval = identity_matrix<int64NDArray> (nr, nc);
4992  break;
4993 
4995  retval = identity_matrix<uint64NDArray> (nr, nc);
4996  break;
4997 
4999  retval = FloatDiagMatrix (nr, nc, 1.0f);
5000  break;
5001 
5003  retval = DiagMatrix (nr, nc, 1.0);
5004  break;
5005 
5007  retval = identity_matrix<boolNDArray> (nr, nc);
5008  break;
5009 
5010  default:
5011  error ("eye: invalid class name");
5012  break;
5013  }
5014 
5015  return retval;
5016 }
5017 
5018 #undef INT_EYE_MATRIX
5019 
5020 DEFUN (eye, args, ,
5021  doc: /* -*- texinfo -*-
5022 @deftypefn {} {} eye (@var{n})
5023 @deftypefnx {} {} eye (@var{m}, @var{n})
5024 @deftypefnx {} {} eye ([@var{m} @var{n}])
5025 @deftypefnx {} {} eye (@dots{}, @var{class})
5026 Return an identity matrix.
5027 
5028 If invoked with a single scalar argument @var{n}, return a square
5029 @nospell{NxN} identity matrix.
5030 
5031 If supplied two scalar arguments (@var{m}, @var{n}), @code{eye} takes them
5032 to be the number of rows and columns. If given a vector with two elements,
5033 @code{eye} uses the values of the elements as the number of rows and
5034 columns, respectively. For example:
5035 
5036 @example
5037 @group
5038 eye (3)
5039  @result{} 1 0 0
5040  0 1 0
5041  0 0 1
5042 @end group
5043 @end example
5044 
5045 The following expressions all produce the same result:
5046 
5047 @example
5048 @group
5049 eye (2)
5050 @equiv{}
5051 eye (2, 2)
5052 @equiv{}
5053 eye (size ([1, 2; 3, 4]))
5054 @end group
5055 @end example
5056 
5057 The optional argument @var{class}, allows @code{eye} to return an array of
5058 the specified type, like
5059 
5060 @example
5061 val = zeros (n,m, "uint8")
5062 @end example
5063 
5064 Calling @code{eye} with no arguments is equivalent to calling it with an
5065 argument of 1. Any negative dimensions are treated as zero. These odd
5066 definitions are for compatibility with @sc{matlab}.
5067 @seealso{speye, ones, zeros}
5068 @end deftypefn */)
5069 {
5070  int nargin = args.length ();
5071 
5073 
5074  // Check for type information.
5075 
5076  if (nargin > 0 && args(nargin-1).is_string ())
5077  {
5078  std::string nm = args(nargin-1).string_value ();
5079  nargin--;
5080 
5082  }
5083 
5084  if (nargin > 2)
5085  print_usage ();
5086 
5088 
5089  if (nargin == 0)
5090  retval = identity_matrix (1, 1, dt);
5091  else if (nargin == 1)
5092  {
5093  octave_idx_type nr, nc;
5094  get_dimensions (args(0), "eye", nr, nc);
5095 
5096  retval = identity_matrix (nr, nc, dt);
5097  }
5098  else
5099  {
5100  octave_idx_type nr, nc;
5101  get_dimensions (args(0), args(1), "eye", nr, nc);
5102 
5103  retval = identity_matrix (nr, nc, dt);
5104  }
5105 
5106  return retval;
5107 }
5108 
5109 /*
5110 %!assert (full (eye (3)), [1, 0, 0; 0, 1, 0; 0, 0, 1])
5111 %!assert (full (eye (2, 3)), [1, 0, 0; 0, 1, 0])
5112 
5113 %!assert (full (eye (3,"single")), single ([1, 0, 0; 0, 1, 0; 0, 0, 1]))
5114 %!assert (full (eye (2, 3,"single")), single ([1, 0, 0; 0, 1, 0]))
5115 
5116 %!assert (eye (3, "int8"), int8 ([1, 0, 0; 0, 1, 0; 0, 0, 1]))
5117 %!assert (eye (2, 3, "int8"), int8 ([1, 0, 0; 0, 1, 0]))
5118 
5119 %!error eye (1, 2, 3)
5120 */
5121 
5122 template <typename MT>
5123 static octave_value
5124 do_linspace (const octave_value& base, const octave_value& limit,
5125  octave_idx_type n)
5126 {
5127  typedef typename MT::column_vector_type CVT;
5128  typedef typename MT::element_type T;
5129 
5131 
5132  if (base.is_scalar_type ())
5133  {
5134  T bs = octave_value_extract<T> (base);
5135  if (limit.is_scalar_type ())
5136  {
5137  T ls = octave_value_extract<T> (limit);
5138  retval = linspace (bs, ls, n);
5139  }
5140  else
5141  {
5142  CVT lv = octave_value_extract<CVT> (limit);
5143  CVT bv (lv.numel (), bs);
5144  retval = linspace (bv, lv, n);
5145  }
5146  }
5147  else
5148  {
5149  CVT bv = octave_value_extract<CVT> (base);
5150  if (limit.is_scalar_type ())
5151  {
5152  T ls = octave_value_extract<T> (limit);
5153  CVT lv (bv.numel (), ls);
5154  retval = linspace (bv, lv, n);
5155  }
5156  else
5157  {
5158  CVT lv = octave_value_extract<CVT> (limit);
5159  retval = linspace (bv, lv, n);
5160  }
5161  }
5162 
5163  return retval;
5164 }
5165 
5166 DEFUN (linspace, args, ,
5167  doc: /* -*- texinfo -*-
5168 @deftypefn {} {} linspace (@var{base}, @var{limit})
5169 @deftypefnx {} {} linspace (@var{base}, @var{limit}, @var{n})
5170 Return a row vector with @var{n} linearly spaced elements between
5171 @var{base} and @var{limit}.
5172 
5173 If the number of elements is greater than one, then the endpoints @var{base}
5174 and @var{limit} are always included in the range. If @var{base} is greater
5175 than @var{limit}, the elements are stored in decreasing order. If the
5176 number of points is not specified, a value of 100 is used.
5177 
5178 The @code{linspace} function returns a row vector when both @var{base}
5179 and @var{limit} are scalars. If one, or both, inputs are vectors, then
5180 @code{linspace} transforms them to column vectors and returns a matrix where
5181 each row is an independent sequence between
5182 @w{@code{@var{base}(@var{row_n}), @var{limit}(@var{row_n})}}.
5183 
5184 For compatibility with @sc{matlab}, return the second argument (@var{limit})
5185 if fewer than two values are requested.
5186 @seealso{logspace}
5187 @end deftypefn */)
5188 {
5189  int nargin = args.length ();
5190 
5191  if (nargin != 2 && nargin != 3)
5192  print_usage ();
5193 
5194  octave_idx_type npoints = 100;
5195  if (nargin == 3)
5196  {
5197  // Apparently undocumented Matlab. If the third arg is an empty
5198  // numeric value, the number of points defaults to 1.
5199  octave_value arg_3 = args(2);
5200 
5201  if (arg_3.is_numeric_type () && arg_3.is_empty ())
5202  npoints = 1;
5203  else if (! arg_3.is_scalar_type ())
5204  error ("linspace: N must be a scalar");
5205  else
5206  // Even if third arg is not an integer, it must be cast to int
5207  npoints = arg_3.idx_type_value ();
5208  }
5209 
5210  octave_value arg_1 = args(0);
5211  octave_value arg_2 = args(1);
5212 
5213  dim_vector sz1 = arg_1.dims ();
5214  bool isvector1 = sz1.ndims () == 2 && (sz1(0) == 1 || sz1(1) == 1);
5215  dim_vector sz2 = arg_2.dims ();
5216  bool isvector2 = sz2.ndims () == 2 && (sz2(0) == 1 || sz2(1) == 1);
5217 
5218  if (! isvector1 || ! isvector2)
5219  error ("linspace: A, B must be scalars or vectors");
5220 
5222 
5223  if (arg_1.is_single_type () || arg_2.is_single_type ())
5224  {
5225  if (arg_1.is_complex_type () || arg_2.is_complex_type ())
5226  retval = do_linspace<FloatComplexMatrix> (arg_1, arg_2, npoints);
5227  else
5228  retval = do_linspace<FloatMatrix> (arg_1, arg_2, npoints);
5229  }
5230  else
5231  {
5232  if (arg_1.is_complex_type () || arg_2.is_complex_type ())
5233  retval = do_linspace<ComplexMatrix> (arg_1, arg_2, npoints);
5234  else
5235  retval = do_linspace<Matrix> (arg_1, arg_2, npoints);
5236  }
5237 
5238  return retval;
5239 }
5240 
5241 /*
5242 %!test
5243 %! x1 = linspace (1, 2);
5244 %! x2 = linspace (1, 2, 10);
5245 %! x3 = linspace (1, -2, 10);
5246 %! assert (size (x1) == [1, 100] && x1(1) == 1 && x1(100) == 2);
5247 %! assert (size (x2) == [1, 10] && x2(1) == 1 && x2(10) == 2);
5248 %! assert (size (x3) == [1, 10] && x3(1) == 1 && x3(10) == -2);
5249 
5250 ## Test complex values
5251 %!test
5252 %! exp = [1+0i, 2-1.25i, 3-2.5i, 4-3.75i, 5-5i];
5253 %! obs = linspace (1, 5-5i, 5);
5254 %! assert (obs, exp);
5255 
5256 ## Test support for vectors in BASE and LIMIT
5257 %!assert (linspace ([1 2 3], [7 8 9]),
5258 %! [linspace(1, 7); linspace(2, 8); linspace(3, 9)], 10*eps)
5259 %!assert (linspace ([1 2 3]', [7 8 9]'),
5260 %! [linspace(1, 7); linspace(2, 8); linspace(3, 9)], 10*eps)
5261 %!assert (linspace ([1 2 3], 9),
5262 %! [linspace(1, 9); linspace(2, 9); linspace(3, 9)], 10*eps)
5263 %!assert (linspace ([1 2 3]', 9),
5264 %! [linspace(1, 9); linspace(2, 9); linspace(3, 9)], 10*eps)
5265 %!assert (linspace (1, [7 8 9]),
5266 %! [linspace(1, 7); linspace(1, 8); linspace(1, 9)], 10*eps)
5267 %!assert (linspace (1, [7 8 9]'),
5268 %! [linspace(1, 7); linspace(1, 8); linspace(1, 9)], 10*eps)
5269 
5270 ## Test class of output
5271 %!assert (class (linspace (1, 2)), "double")
5272 %!assert (class (linspace (single (1), 2)), "single")
5273 %!assert (class (linspace (1, single (2))), "single")
5274 
5275 ## Test obscure Matlab compatibility options
5276 %!assert (linspace (0, 1, []), 1)
5277 %!assert (linspace (10, 20, 2), [10 20])
5278 %!assert (linspace (10, 20, 1), [20])
5279 %!assert (linspace (10, 20, 0), zeros (1, 0))
5280 %!assert (linspace (10, 20, -1), zeros (1, 0))
5281 %!assert (numel (linspace (0, 1, 2+eps)), 2)
5282 %!assert (numel (linspace (0, 1, 2-eps)), 1)
5283 %!assert (linspace (10, 20, 2.1), [10 20])
5284 %!assert (linspace (10, 20, 2.9), [10 20])
5285 
5286 %!error linspace ()
5287 %!error linspace (1, 2, 3, 4)
5288 %!error <N must be a scalar> linspace (1, 2, [3, 4])
5289 %!error <must be scalars or vectors> linspace (ones (2,2), 2, 3)
5290 %!error <must be scalars or vectors> linspace (2, ones (2,2), 3)
5291 %!error <must be scalars or vectors> linspace (1, [], 3)
5292 */
5293 
5294 // FIXME: should accept dimensions as separate args for N-D
5295 // arrays as well as 1-D and 2-D arrays.
5296 
5297 DEFUN (resize, args, ,
5298  doc: /* -*- texinfo -*-
5299 @deftypefn {} {} resize (@var{x}, @var{m})
5300 @deftypefnx {} {} resize (@var{x}, @var{m}, @var{n}, @dots{})
5301 @deftypefnx {} {} resize (@var{x}, [@var{m} @var{n} @dots{}])
5302 Resize @var{x} cutting off elements as necessary.
5303 
5304 In the result, element with certain indices is equal to the corresponding
5305 element of @var{x} if the indices are within the bounds of @var{x};
5306 otherwise, the element is set to zero.
5307 
5308 In other words, the statement
5309 
5310 @example
5311 y = resize (x, dv)
5312 @end example
5313 
5314 @noindent
5315 is equivalent to the following code:
5316 
5317 @example
5318 @group
5319 y = zeros (dv, class (x));
5320 sz = min (dv, size (x));
5321 for i = 1:length (sz)
5322  idx@{i@} = 1:sz(i);
5323 endfor
5324 y(idx@{:@}) = x(idx@{:@});
5325 @end group
5326 @end example
5327 
5328 @noindent
5329 but is performed more efficiently.
5330 
5331 If only @var{m} is supplied, and it is a scalar, the dimension of the
5332 result is @var{m}-by-@var{m}.
5333 If @var{m}, @var{n}, @dots{} are all scalars, then the dimensions of
5334 the result are @var{m}-by-@var{n}-by-@dots{}.
5335 If given a vector as input, then the
5336 dimensions of the result are given by the elements of that vector.
5337 
5338 An object can be resized to more dimensions than it has;
5339 in such case the missing dimensions are assumed to be 1.
5340 Resizing an object to fewer dimensions is not possible.
5341 @seealso{reshape, postpad, prepad, cat}
5342 @end deftypefn */)
5343 {
5344  int nargin = args.length ();
5345 
5346  if (nargin < 2)
5347  print_usage ();
5348 
5350 
5351  if (nargin == 2)
5352  {
5353  Array<double> vec = args(1).vector_value ();
5354  int ndim = vec.numel ();
5355  if (ndim == 1)
5356  {
5357  octave_idx_type m = static_cast<octave_idx_type> (vec(0));
5358  retval = args(0);
5359  retval = retval.resize (dim_vector (m, m), true);
5360  }
5361  else
5362  {
5363  dim_vector dv;
5364  dv.resize (ndim);
5365  for (int i = 0; i < ndim; i++)
5366  dv(i) = static_cast<octave_idx_type> (vec(i));
5367  retval = args(0);
5368  retval = retval.resize (dv, true);
5369  }
5370  }
5371  else
5372  {
5373  dim_vector dv;
5374  dv.resize (nargin - 1);
5375  for (octave_idx_type i = 1; i < nargin; i++)
5376  dv(i-1) = static_cast<octave_idx_type> (args(i).scalar_value ());
5377 
5378  retval = args(0);
5379  retval = retval.resize (dv, true);
5380  }
5381 
5382  return retval;
5383 }
5384 
5385 // FIXME: should use octave_idx_type for dimensions.
5386 
5387 DEFUN (reshape, args, ,
5388  doc: /* -*- texinfo -*-
5389 @deftypefn {} {} reshape (@var{A}, @var{m}, @var{n}, @dots{})
5390 @deftypefnx {} {} reshape (@var{A}, [@var{m} @var{n} @dots{}])
5391 @deftypefnx {} {} reshape (@var{A}, @dots{}, [], @dots{})
5392 @deftypefnx {} {} reshape (@var{A}, @var{size})
5393 Return a matrix with the specified dimensions (@var{m}, @var{n}, @dots{})
5394 whose elements are taken from the matrix @var{A}.
5395 
5396 The elements of the matrix are accessed in column-major order (like Fortran
5397 arrays are stored).
5398 
5399 The following code demonstrates reshaping a 1x4 row vector into a 2x2 square
5400 matrix.
5401 
5402 @example
5403 @group
5404 reshape ([1, 2, 3, 4], 2, 2)
5405  @result{} 1 3
5406  2 4
5407 @end group
5408 @end example
5409 
5410 @noindent
5411 Note that the total number of elements in the original matrix
5412 (@code{prod (size (@var{A}))}) must match the total number of elements
5413 in the new matrix (@code{prod ([@var{m} @var{n} @dots{}])}).
5414 
5415 A single dimension of the return matrix may be left unspecified and Octave
5416 will determine its size automatically. An empty matrix ([]) is used to flag
5417 the unspecified dimension.
5418 @seealso{resize, vec, postpad, cat, squeeze}
5419 @end deftypefn */)
5420 {
5421  int nargin = args.length ();
5422 
5423  if (nargin < 2)
5424  print_usage ();
5425 
5427 
5428  dim_vector new_dims;
5429 
5430  if (nargin == 2)
5431  {
5432  Array<octave_idx_type> new_size = args(1).octave_idx_type_vector_value ();
5433 
5434  if (new_size.numel () < 2)
5435  error ("reshape: SIZE must have 2 or more dimensions");
5436 
5437  new_dims = dim_vector::alloc (new_size.numel ());
5438 
5439  for (octave_idx_type i = 0; i < new_size.numel (); i++)
5440  {
5441  if (new_size(i) < 0)
5442  error ("reshape: SIZE must be non-negative");
5443 
5444  new_dims(i) = new_size(i);
5445  }
5446  }
5447  else
5448  {
5449  new_dims = dim_vector::alloc (nargin-1);
5450  int empty_dim = -1;
5451 
5452  for (int i = 1; i < nargin; i++)
5453  {
5454  if (args(i).is_empty ())
5455  {
5456  if (empty_dim > 0)
5457  error ("reshape: only a single dimension can be unknown");
5458 
5459  empty_dim = i;
5460  new_dims(i-1) = 1;
5461  }
5462  else
5463  {
5464  new_dims(i-1) = args(i).idx_type_value ();
5465 
5466  if (new_dims(i-1) < 0)
5467  error ("reshape: SIZE must be non-negative");
5468  }
5469  }
5470 
5471  if (empty_dim > 0)
5472  {
5473  octave_idx_type nel = new_dims.numel ();
5474 
5475  if (nel == 0)
5476  new_dims(empty_dim-1) = 0;
5477  else
5478  {
5479  octave_idx_type a_nel = args(0).numel ();
5480  octave_idx_type size_empty_dim = a_nel / nel;
5481 
5482  if (a_nel != size_empty_dim * nel)
5483  error ("reshape: SIZE is not divisible by the product of known dimensions (= %d)",
5484  nel);
5485 
5486  new_dims(empty_dim-1) = size_empty_dim;
5487  }
5488  }
5489  }
5490 
5491  retval = args(0).reshape (new_dims);
5492 
5493  return retval;
5494 }
5495 
5496 /*
5497 %!assert (size (reshape (ones (4, 4), 2, 8)), [2, 8])
5498 %!assert (size (reshape (ones (4, 4), 8, 2)), [8, 2])
5499 %!assert (size (reshape (ones (15, 4), 1, 60)), [1, 60])
5500 %!assert (size (reshape (ones (15, 4), 60, 1)), [60, 1])
5501 
5502 %!assert (size (reshape (ones (4, 4, "single"), 2, 8)), [2, 8])
5503 %!assert (size (reshape (ones (4, 4, "single"), 8, 2)), [8, 2])
5504 %!assert (size (reshape (ones (15, 4, "single"), 1, 60)), [1, 60])
5505 %!assert (size (reshape (ones (15, 4, "single"), 60, 1)), [60, 1])
5506 
5507 %!test
5508 %! s.a = 1;
5509 %! fail ("reshape (s, 2, 3)", "can't reshape 1x1 array to 2x3 array");
5510 
5511 %!error reshape ()
5512 %!error reshape (1, 2, 3, 4)
5513 %!error <SIZE must have 2 or more dimensions> reshape (1:3, 3)
5514 %!error <SIZE must be non-negative> reshape (1:3, [3 -1])
5515 %!error <only a single dimension can be unknown> reshape (1:3, 1,[],[],3)
5516 %!error <SIZE must be non-negative> reshape (1:3, 3, -1)
5517 %!error <SIZE is not divisible> reshape (1:3, 3, [], 2)
5518 */
5519 
5520 DEFUN (vec, args, ,
5521  doc: /* -*- texinfo -*-
5522 @deftypefn {} {@var{v} =} vec (@var{x})
5523 @deftypefnx {} {@var{v} =} vec (@var{x}, @var{dim})
5524 Return the vector obtained by stacking the columns of the matrix @var{x}
5525 one above the other.
5526 
5527 Without @var{dim} this is equivalent to @code{@var{x}(:)}.
5528 
5529 If @var{dim} is supplied, the dimensions of @var{v} are set to @var{dim}
5530 with all elements along the last dimension. This is equivalent to
5531 @code{shiftdim (@var{x}(:), 1-@var{dim})}.
5532 @seealso{vech, resize, cat}
5533 @end deftypefn */)
5534 {
5535  int nargin = args.length ();
5536 
5537  if (nargin < 1 || nargin > 2)
5538  print_usage ();
5539 
5540  int dim = 1;
5541  if (nargin == 2)
5542  {
5543  dim = args(1).idx_type_value ();
5544 
5545  if (dim < 1)
5546  error ("vec: DIM must be greater than zero");
5547  }
5548 
5550  octave_value arg = args(0);
5551 
5552  octave_value retval = arg.single_subsref ("(", colon);
5553 
5554  if (dim > 1)
5555  {
5556  dim_vector new_dims = dim_vector::alloc (dim);
5557 
5558  for (int i = 0; i < dim-1; i++)
5559  new_dims(i) = 1;
5560 
5561  new_dims(dim-1) = retval.numel ();
5562 
5563  retval = retval.reshape (new_dims);
5564  }
5565 
5566  return retval;
5567 }
5568 
5569 /*
5570 %!assert (vec ([1, 2; 3, 4]), [1; 3; 2; 4])
5571 %!assert (vec ([1, 3, 2, 4]), [1; 3; 2; 4])
5572 %!assert (vec ([1, 2, 3, 4], 2), [1, 2, 3, 4])
5573 %!assert (vec ([1, 2; 3, 4]), vec ([1, 2; 3, 4], 1))
5574 %!assert (vec ([1, 2; 3, 4], 1), [1; 3; 2; 4])
5575 %!assert (vec ([1, 2; 3, 4], 2), [1, 3, 2, 4])
5576 %!assert (vec ([1, 3; 2, 4], 3), reshape ([1, 2, 3, 4], 1, 1, 4))
5577 %!assert (vec ([1, 3; 2, 4], 3), shiftdim (vec ([1, 3; 2, 4]), -2))
5578 
5579 %!error vec ()
5580 %!error vec (1, 2, 3)
5581 %!error vec ([1, 2; 3, 4], 0)
5582 */
5583 
5584 DEFUN (squeeze, args, ,
5585  doc: /* -*- texinfo -*-
5586 @deftypefn {} {} squeeze (@var{x})
5587 Remove singleton dimensions from @var{x} and return the result.
5588 
5589 Note that for compatibility with @sc{matlab}, all objects have
5590 a minimum of two dimensions and row vectors are left unchanged.
5591 @seealso{reshape}
5592 @end deftypefn */)
5593 {
5594  if (args.length () != 1)
5595  print_usage ();
5596 
5597  return ovl (args(0).squeeze ());
5598 }
5599 
5600 DEFUN (full, args, ,
5601  doc: /* -*- texinfo -*-
5602 @deftypefn {} {@var{FM} =} full (@var{SM})
5603 Return a full storage matrix from a sparse, diagonal, or permutation matrix,
5604 or a range.
5605 @seealso{sparse, issparse}
5606 @end deftypefn */)
5607 {
5608  if (args.length () != 1)
5609  print_usage ();
5610 
5611  return ovl (args(0).full_value ());
5612 }
5613 
5614 // Compute various norms of the vector X.
5615 
5616 DEFUN (norm, args, ,
5617  doc: /* -*- texinfo -*-
5618 @deftypefn {} {} norm (@var{A})
5619 @deftypefnx {} {} norm (@var{A}, @var{p})
5620 @deftypefnx {} {} norm (@var{A}, @var{p}, @var{opt})
5621 Compute the p-norm of the matrix @var{A}.
5622 
5623 If the second argument is not given, @w{@code{p = 2}} is used.
5624 
5625 If @var{A} is a matrix (or sparse matrix):
5626 
5627 @table @asis
5628 @item @var{p} = @code{1}
5629 1-norm, the largest column sum of the absolute values of @var{A}.
5630 
5631 @item @var{p} = @code{2}
5632 Largest singular value of @var{A}.
5633 
5634 @item @var{p} = @code{Inf} or @qcode{"inf"}
5635 @cindex infinity norm
5636 Infinity norm, the largest row sum of the absolute values of @var{A}.
5637 
5638 @item @var{p} = @qcode{"fro"}
5639 @cindex Frobenius norm
5640 Frobenius norm of @var{A}, @code{sqrt (sum (diag (@var{A}' * @var{A})))}.
5641 
5642 @item other @var{p}, @code{@var{p} > 1}
5643 @cindex general p-norm
5644 maximum @code{norm (A*x, p)} such that @code{norm (x, p) == 1}
5645 @end table
5646 
5647 If @var{A} is a vector or a scalar:
5648 
5649 @table @asis
5650 @item @var{p} = @code{Inf} or @qcode{"inf"}
5651 @code{max (abs (@var{A}))}.
5652 
5653 @item @var{p} = @code{-Inf}
5654 @code{min (abs (@var{A}))}.
5655 
5656 @item @var{p} = @qcode{"fro"}
5657 Frobenius norm of @var{A}, @code{sqrt (sumsq (abs (A)))}.
5658 
5659 @item @var{p} = 0
5660 Hamming norm---the number of nonzero elements.
5661 
5662 @item other @var{p}, @code{@var{p} > 1}
5663 p-norm of @var{A}, @code{(sum (abs (@var{A}) .^ @var{p})) ^ (1/@var{p})}.
5664 
5665 @item other @var{p} @code{@var{p} < 1}
5666 the p-pseudonorm defined as above.
5667 @end table
5668 
5669 If @var{opt} is the value @qcode{"rows"}, treat each row as a vector and
5670 compute its norm. The result is returned as a column vector.
5671 Similarly, if @var{opt} is @qcode{"columns"} or @qcode{"cols"} then
5672 compute the norms of each column and return a row vector.
5673 @seealso{normest, normest1, cond, svd}
5674 @end deftypefn */)
5675 {
5676  int nargin = args.length ();
5677 
5678  if (nargin < 1 || nargin > 3)
5679  print_usage ();
5680 
5681  octave_value x_arg = args(0);
5682 
5683  if (x_arg.ndims () != 2)
5684  error ("norm: only valid for 2-D objects");
5685 
5686  enum { sfmatrix, sfcols, sfrows, sffrob, sfinf } strflag = sfmatrix;
5687  if (nargin > 1 && args(nargin-1).is_string ())
5688  {
5689  std::string str = args(nargin-1).string_value ();
5690  if (str == "cols" || str == "columns")
5691  strflag = sfcols;
5692  else if (str == "rows")
5693  strflag = sfrows;
5694  else if (str == "fro")
5695  strflag = sffrob;
5696  else if (str == "inf")
5697  strflag = sfinf;
5698  else
5699  error ("norm: unrecognized option: %s", str.c_str ());
5700 
5701  // we've handled the last parameter, so act as if it was removed
5702  nargin--;
5703  }
5704 
5705  octave_value p_arg = (nargin > 1) ? args(1) : octave_value (2);
5706 
5707  if (p_arg.is_empty ())
5708  p_arg = octave_value (2);
5709  else if (p_arg.is_string ())
5710  {
5711  std::string str = p_arg.string_value ();
5712  if (strflag != sfcols && strflag != sfrows)
5713  error ("norm: invalid combination of options");
5714 
5715  if (str == "cols" || str == "columns" || str == "rows")
5716  error ("norm: invalid combination of options");
5717 
5718  if (str == "fro")
5719  p_arg = octave_value (2);
5720  else if (str == "inf")
5722  else
5723  error ("norm: unrecognized option: %s", str.c_str ());
5724  }
5725  else if (! p_arg.is_scalar_type ())
5726  err_wrong_type_arg ("norm", p_arg);
5727 
5729 
5730  switch (strflag)
5731  {
5732  case sfmatrix:
5733  retval = xnorm (x_arg, p_arg);
5734  break;
5735 
5736  case sfcols:
5737  retval = xcolnorms (x_arg, p_arg);
5738  break;
5739 
5740  case sfrows:
5741  retval = xrownorms (x_arg, p_arg);
5742  break;
5743 
5744  case sffrob:
5745  retval = xfrobnorm (x_arg);
5746  break;
5747 
5748  case sfinf:
5749  retval = xnorm (x_arg, octave::numeric_limits<double>::Inf ());
5750  break;
5751  }
5752 
5753  return retval;
5754 }
5755 
5756 /*
5757 %!shared x
5758 %! x = [1, -3, 4, 5, -7];
5759 %!assert (norm (x,1), 20)
5760 %!assert (norm (x,2), 10)
5761 %!assert (norm (x,3), 8.24257059961711, -4*eps)
5762 %!assert (norm (x,Inf), 7)
5763 %!assert (norm (x,-Inf), 1)
5764 %!assert (norm (x,"inf"), 7)
5765 %!assert (norm (x,"fro"), 10, -eps)
5766 %!assert (norm (x), 10)
5767 %!assert (norm ([1e200, 1]), 1e200)
5768 %!assert (norm ([3+4i, 3-4i, sqrt(31)]), 9, -4*eps)
5769 %!shared m
5770 %! m = magic (4);
5771 %!assert (norm (m,1), 34)
5772 %!assert (norm (m,2), 34, -eps)
5773 %!assert (norm (m,Inf), 34)
5774 %!assert (norm (m,"inf"), 34)
5775 %!shared m2, flo, fhi
5776 %! m2 = [1,2;3,4];
5777 %! flo = 1e-300;
5778 %! fhi = 1e+300;
5779 %!assert (norm (flo*m2,"fro"), sqrt (30)*flo, -eps)
5780 %!assert (norm (fhi*m2,"fro"), sqrt (30)*fhi, -eps)
5781 
5782 %!shared x
5783 %! x = single ([1, -3, 4, 5, -7]);
5784 %!assert (norm (x,1), single (20))
5785 %!assert (norm (x,2), single (10))
5786 %!assert (norm (x,3), single (8.24257059961711), -4*eps ("single"))
5787 %!assert (norm (x,Inf), single (7))
5788 %!assert (norm (x,-Inf), single (1))
5789 %!assert (norm (x,"inf"), single (7))
5790 %!assert (norm (x,"fro"), single (10), -eps ("single"))
5791 %!assert (norm (x), single (10))
5792 %!assert (norm (single ([1e38, 1])), single (1e38))
5793 %!assert (norm (single ([3+4i, 3-4i, sqrt(31)])), single (9), -4*eps ("single"))
5794 %!shared m
5795 %! m = single (magic (4));
5796 %!assert (norm (m,1), single (34))
5797 %!assert (norm (m,2), single (34), -eps ("single"))
5798 %!assert (norm (m,Inf), single (34))
5799 %!assert (norm (m,"inf"), single (34))
5800 %!shared m2, flo, fhi
5801 %! m2 = single ([1,2;3,4]);
5802 %! flo = single (1e-300);
5803 %! fhi = single (1e+300);
5804 %!assert (norm (flo*m2,"fro"), single (sqrt (30)*flo), -eps ("single"))
5805 %!assert (norm (fhi*m2,"fro"), single (sqrt (30)*fhi), -eps ("single"))
5806 
5807 %!shared q
5808 %! q = rand (1e3, 3);
5809 %!assert (norm (q, 3, "rows"), sum (q.^3, 2).^(1/3), sqrt (eps))
5810 %!assert (norm (q, "fro", "rows"), sum (q.^2, 2).^(1/2), sqrt (eps))
5811 %!assert (norm (q, "fro", "rows"), sqrt (sumsq (q, 2)), sqrt (eps))
5812 %!assert (norm (q, "fro", "cols"), sqrt (sumsq (q, 1)), sqrt (eps))
5813 %!assert (norm (q, 3, "cols"), sum (q.^3, 1).^(1/3), sqrt (eps))
5814 %!assert (norm (q, "inf", "rows"), norm (q, Inf, "rows"))
5815 %!assert (norm (q, "inf", "cols"), norm (q, Inf, "cols"))
5816 %!assert (norm (q, [], "rows"), norm (q, 2, "rows"))
5817 %!assert (norm (q, [], "cols"), norm (q, 2, "cols"))
5818 
5819 %!test <30631>
5820 %! ## Test for norm returning NaN on sparse matrix
5821 %! A = sparse (2,2);
5822 %! A(2,1) = 1;
5823 %! assert (norm (A), 1);
5824 
5825 ## Test input validation
5826 %!error norm ()
5827 %!error norm (1,2,3,4)
5828 %!error <unrecognized option> norm (1, "invalid")
5829 %!error <unrecognized option> norm (1, "rows", "invalid")
5830 %!error <unrecognized option> norm (1, "invalid", "rows")
5831 %!error <invalid combination of options> norm (1, "cols", "rows")
5832 %!error <invalid combination of options> norm (1, "rows", "rows")
5833 */
5834 
5835 static octave_value
5836 unary_op_defun_body (octave_value::unary_op op,
5837  const octave_value_list& args)
5838 {
5839  if (args.length () != 1)
5840  print_usage ();
5841 
5842  return do_unary_op (op, args(0));
5843 }
5844 
5845 DEFUN (not, args, ,
5846  doc: /* -*- texinfo -*-
5847 @deftypefn {} {@var{z} =} not (@var{x})
5848 Return the logical NOT of @var{x}.
5849 
5850 This function is equivalent to the operator syntax @w{@code{! x}}.
5851 @seealso{and, or, xor}
5852 @end deftypefn */)
5853 {
5854  return unary_op_defun_body (octave_value::op_not, args);
5855 }
5856 
5857 DEFUN (uplus, args, ,
5858  doc: /* -*- texinfo -*-
5859 @deftypefn {} {} uplus (@var{x})
5860 This function and @w{@tcode{+ x}} are equivalent.
5861 @seealso{uminus, plus, minus}
5862 @end deftypefn */)
5863 {
5864  return unary_op_defun_body (octave_value::op_uplus, args);
5865 }
5866 
5867 DEFUN (uminus, args, ,
5868  doc: /* -*- texinfo -*-
5869 @deftypefn {} {} uminus (@var{x})
5870 This function and @w{@tcode{- x}} are equivalent.
5871 @seealso{uplus, minus}
5872 @end deftypefn */)
5873 {
5874  return unary_op_defun_body (octave_value::op_uminus, args);
5875 }
5876 
5877 DEFUN (transpose, args, ,
5878  doc: /* -*- texinfo -*-
5879 @deftypefn {} {} transpose (@var{x})
5880 Return the transpose of @var{x}.
5881 
5882 This function and @tcode{x.'} are equivalent.
5883 @seealso{ctranspose}
5884 @end deftypefn */)
5885 {
5886  return unary_op_defun_body (octave_value::op_transpose, args);
5887 }
5888 
5889 /*
5890 %!assert (2.', 2)
5891 %!assert (2i.', 2i)
5892 %!assert ([1:4].', [1;2;3;4])
5893 %!assert ([1;2;3;4].', [1:4])
5894 %!assert ([1,2;3,4].', [1,3;2,4])
5895 %!assert ([1,2i;3,4].', [1,3;2i,4])
5896 
5897 %!assert (transpose ([1,2;3,4]), [1,3;2,4])
5898 
5899 %!assert (single (2).', single (2))
5900 %!assert (single (2i).', single (2i))
5901 %!assert (single ([1:4]).', single ([1;2;3;4]))
5902 %!assert (single ([1;2;3;4]).', single ([1:4]))
5903 %!assert (single ([1,2;3,4]).', single ([1,3;2,4]))
5904 %!assert (single ([1,2i;3,4]).', single ([1,3;2i,4]))
5905 
5906 %!assert (transpose (single ([1,2;3,4])), single ([1,3;2,4]))
5907 */
5908 
5909 DEFUN (ctranspose, args, ,
5910  doc: /* -*- texinfo -*-
5911 @deftypefn {} {} ctranspose (@var{x})
5912 Return the complex conjugate transpose of @var{x}.
5913 
5914 This function and @tcode{x'} are equivalent.
5915 @seealso{transpose}
5916 @end deftypefn */)
5917 {
5918  return unary_op_defun_body (octave_value::op_hermitian, args);
5919 }
5920 
5921 /*
5922 %!assert (2', 2)
5923 %!assert (2i', -2i)
5924 %!assert ([1:4]', [1;2;3;4])
5925 %!assert ([1;2;3;4]', [1:4])
5926 %!assert ([1,2;3,4]', [1,3;2,4])
5927 %!assert ([1,2i;3,4]', [1,3;-2i,4])
5928 
5929 %!assert (ctranspose ([1,2i;3,4]), [1,3;-2i,4])
5930 
5931 %!assert (single (2)', single (2))
5932 %!assert (single (2i)', single (-2i))
5933 %!assert (single ([1:4])', single ([1;2;3;4]))
5934 %!assert (single ([1;2;3;4])', single ([1:4]))
5935 %!assert (single ([1,2;3,4])', single ([1,3;2,4]))
5936 %!assert (single ([1,2i;3,4])', single ([1,3;-2i,4]))
5937 
5938 %!assert (ctranspose (single ([1,2i;3,4])), single ([1,3;-2i,4]))
5939 */
5940 
5941 static octave_value
5942 binary_op_defun_body (octave_value::binary_op op,
5943  const octave_value_list& args)
5944 {
5945  if (args.length () != 2)
5946  print_usage ();
5947 
5948  return do_binary_op (op, args(0), args(1));
5949 }
5950 
5951 static octave_value
5952 binary_assoc_op_defun_body (octave_value::binary_op op,
5954  const octave_value_list& args)
5955 {
5956  int nargin = args.length ();
5957 
5958  if (nargin == 0)
5959  print_usage ();
5960 
5962 
5963  if (nargin == 1)
5964  retval = args(0);
5965  else if (nargin == 2)
5966  retval = do_binary_op (op, args(0), args(1));
5967  else
5968  {
5969  retval = do_binary_op (op, args(0), args(1));
5970 
5971  for (int i = 2; i < nargin; i++)
5972  retval.assign (aop, args(i));
5973  }
5974 
5975  return retval;
5976 }
5977 
5978 DEFUN (plus, args, ,
5979  doc: /* -*- texinfo -*-
5980 @deftypefn {} {} plus (@var{x}, @var{y})
5981 @deftypefnx {} {} plus (@var{x1}, @var{x2}, @dots{})
5982 This function and @w{@tcode{x + y}} are equivalent.
5983 
5984 If more arguments are given, the summation is applied
5985 cumulatively from left to right:
5986 
5987 @example
5988 (@dots{}((x1 + x2) + x3) + @dots{})
5989 @end example
5990 
5991 At least one argument is required.
5992 @seealso{minus, uplus}
5993 @end deftypefn */)
5994 {
5995  return binary_assoc_op_defun_body (octave_value::op_add,
5996  octave_value::op_add_eq, args);
5997 }
5998 
5999 DEFUN (minus, args, ,
6000  doc: /* -*- texinfo -*-
6001 @deftypefn {} {} minus (@var{x}, @var{y})
6002 This function and @w{@tcode{x - y}} are equivalent.
6003 @seealso{plus, uminus}
6004 @end deftypefn */)
6005 {
6006  return binary_op_defun_body (octave_value::op_sub, args);
6007 }
6008 
6009 DEFUN (mtimes, args, ,
6010  doc: /* -*- texinfo -*-
6011 @deftypefn {} {} mtimes (@var{x}, @var{y})
6012 @deftypefnx {} {} mtimes (@var{x1}, @var{x2}, @dots{})
6013 Return the matrix multiplication product of inputs.
6014 
6015 This function and @w{@tcode{x * y}} are equivalent.
6016 If more arguments are given, the multiplication is applied
6017 cumulatively from left to right:
6018 
6019 @example
6020 (@dots{}((x1 * x2) * x3) * @dots{})
6021 @end example
6022 
6023 At least one argument is required.
6024 @seealso{times, plus, minus, rdivide, mrdivide, mldivide, mpower}
6025 @end deftypefn */)
6026 {
6027  return binary_assoc_op_defun_body (octave_value::op_mul,
6028  octave_value::op_mul_eq, args);
6029 }
6030 
6031 DEFUN (mrdivide, args, ,
6032  doc: /* -*- texinfo -*-
6033 @deftypefn {} {} mrdivide (@var{x}, @var{y})
6034 Return the matrix right division of @var{x} and @var{y}.
6035 
6036 This function and @w{@tcode{x / y}} are equivalent.
6037 @seealso{mldivide, rdivide, plus, minus}
6038 @end deftypefn */)
6039 {
6040  return binary_op_defun_body (octave_value::op_div, args);
6041 }
6042 
6043 DEFUN (mpower, args, ,
6044  doc: /* -*- texinfo -*-
6045 @deftypefn {} {} mpower (@var{x}, @var{y})
6046 Return the matrix power operation of @var{x} raised to the @var{y} power.
6047 
6048 This function and @w{@tcode{x ^ y}} are equivalent.
6049 @seealso{power, mtimes, plus, minus}
6050 @end deftypefn */)
6051 {
6052  return binary_op_defun_body (octave_value::op_pow, args);
6053 }
6054 
6055 DEFUN (mldivide, args, ,
6056  doc: /* -*- texinfo -*-
6057 @deftypefn {} {} mldivide (@var{x}, @var{y})
6058 Return the matrix left division of @var{x} and @var{y}.
6059 
6060 This function and @w{@tcode{x @xbackslashchar{} y}} are equivalent.
6061 @seealso{mrdivide, ldivide, rdivide}
6062 @end deftypefn */)
6063 {
6064  return binary_op_defun_body (octave_value::op_ldiv, args);
6065 }
6066 
6067 DEFUN (lt, args, ,
6068  doc: /* -*- texinfo -*-
6069 @deftypefn {} {} lt (@var{x}, @var{y})
6070 This function is equivalent to @w{@code{x < y}}.
6071 @seealso{le, eq, ge, gt, ne}
6072 @end deftypefn */)
6073 {
6074  return binary_op_defun_body (octave_value::op_lt, args);
6075 }
6076 
6077 DEFUN (le, args, ,
6078  doc: /* -*- texinfo -*-
6079 @deftypefn {} {} le (@var{x}, @var{y})
6080 This function is equivalent to @w{@code{x <= y}}.
6081 @seealso{eq, ge, gt, ne, lt}
6082 @end deftypefn */)
6083 {
6084  return binary_op_defun_body (octave_value::op_le, args);
6085 }
6086 
6087 DEFUN (eq, args, ,
6088  doc: /* -*- texinfo -*-
6089 @deftypefn {} {} eq (@var{x}, @var{y})
6090 Return true if the two inputs are equal.
6091 
6092 This function is equivalent to @w{@code{x == y}}.
6093 @seealso{ne, isequal, le, ge, gt, ne, lt}
6094 @end deftypefn */)
6095 {
6096  return binary_op_defun_body (octave_value::op_eq, args);
6097 }
6098 
6099 DEFUN (ge, args, ,
6100  doc: /* -*- texinfo -*-
6101 @deftypefn {} {} ge (@var{x}, @var{y})
6102 This function is equivalent to @w{@code{x >= y}}.
6103 @seealso{le, eq, gt, ne, lt}
6104 @end deftypefn */)
6105 {
6106  return binary_op_defun_body (octave_value::op_ge, args);
6107 }
6108 
6109 DEFUN (gt, args, ,
6110  doc: /* -*- texinfo -*-
6111 @deftypefn {} {} gt (@var{x}, @var{y})
6112 This function is equivalent to @w{@code{x > y}}.
6113 @seealso{le, eq, ge, ne, lt}
6114 @end deftypefn */)
6115 {
6116  return binary_op_defun_body (octave_value::op_gt, args);
6117 }
6118 
6119 DEFUN (ne, args, ,
6120  doc: /* -*- texinfo -*-
6121 @deftypefn {} {} ne (@var{x}, @var{y})
6122 Return true if the two inputs are not equal.
6123 
6124 This function is equivalent to @w{@code{x != y}}.
6125 @seealso{eq, isequal, le, ge, lt}
6126 @end deftypefn */)
6127 {
6128  return binary_op_defun_body (octave_value::op_ne, args);
6129 }
6130 
6131 DEFUN (times, args, ,
6132  doc: /* -*- texinfo -*-
6133 @deftypefn {} {} times (@var{x}, @var{y})
6134 @deftypefnx {} {} times (@var{x1}, @var{x2}, @dots{})
6135 Return the element-by-element multiplication product of inputs.
6136 
6137 This function and @w{@tcode{x .* y}} are equivalent.
6138 If more arguments are given, the multiplication is applied
6139 cumulatively from left to right:
6140 
6141 @example
6142 (@dots{}((x1 .* x2) .* x3) .* @dots{})
6143 @end example
6144 
6145 At least one argument is required.
6146 @seealso{mtimes, rdivide}
6147 @end deftypefn */)
6148 {
6149  return binary_assoc_op_defun_body (octave_value::op_el_mul,
6151 }
6152 
6153 DEFUN (rdivide, args, ,
6154  doc: /* -*- texinfo -*-
6155 @deftypefn {} {} rdivide (@var{x}, @var{y})
6156 Return the element-by-element right division of @var{x} and @var{y}.
6157 
6158 This function and @w{@tcode{x ./ y}} are equivalent.
6159 @seealso{ldivide, mrdivide, times, plus}
6160 @end deftypefn */)
6161 {
6162  return binary_op_defun_body (octave_value::op_el_div, args);
6163 }
6164 
6165 DEFUN (power, args, ,
6166  doc: /* -*- texinfo -*-
6167 @deftypefn {} {} power (@var{x}, @var{y})
6168 Return the element-by-element operation of @var{x} raised to the
6169 @var{y} power.
6170 
6171 This function and @w{@tcode{x .^ y}} are equivalent.
6172 
6173 If several complex results are possible, returns the one with smallest
6174 non-negative argument (angle). Use @code{realpow}, @code{realsqrt},
6175 @code{cbrt}, or @code{nthroot} if a real result is preferred.
6176 
6177 @seealso{mpower, realpow, realsqrt, cbrt, nthroot}
6178 @end deftypefn */)
6179 {
6180  return binary_op_defun_body (octave_value::op_el_pow, args);
6181 }
6182 
6183 DEFUN (ldivide, args, ,
6184  doc: /* -*- texinfo -*-
6185 @deftypefn {} {} ldivide (@var{x}, @var{y})
6186 Return the element-by-element left division of @var{x} and @var{y}.
6187 
6188 This function and @w{@tcode{x .@xbackslashchar{} y}} are equivalent.
6189 @seealso{rdivide, mldivide, times, plus}
6190 @end deftypefn */)
6191 {
6192  return binary_op_defun_body (octave_value::op_el_ldiv, args);
6193 }
6194 
6195 DEFUN (and, args, ,
6196  doc: /* -*- texinfo -*-
6197 @deftypefn {} {@var{z} =} and (@var{x}, @var{y})
6198 @deftypefnx {} {@var{z} =} and (@var{x1}, @var{x2}, @dots{})
6199 Return the logical AND of @var{x} and @var{y}.
6200 
6201 This function is equivalent to the operator syntax @w{@code{x & y}}. If
6202 more than two arguments are given, the logical AND is applied cumulatively
6203 from left to right:
6204 
6205 @example
6206 (@dots{}((x1 & x2) & x3) & @dots{})
6207 @end example
6208 
6209 At least one argument is required.
6210 @seealso{or, not, xor}
6211 @end deftypefn */)
6212 {
6213  return binary_assoc_op_defun_body (octave_value::op_el_and,
6215 }
6216 
6217 DEFUN (or, args, ,
6218  doc: /* -*- texinfo -*-
6219 @deftypefn {} {@var{z} =} or (@var{x}, @var{y})
6220 @deftypefnx {} {@var{z} =} or (@var{x1}, @var{x2}, @dots{})
6221 Return the logical OR of @var{x} and @var{y}.
6222 
6223 This function is equivalent to the operator syntax @w{@code{x | y}}. If
6224 more than two arguments are given, the logical OR is applied cumulatively
6225 from left to right:
6226 
6227 @example
6228 (@dots{}((x1 | x2) | x3) | @dots{})
6229 @end example
6230 
6231 At least one argument is required.
6232 @seealso{and, not, xor}
6233 @end deftypefn */)
6234