GNU Octave  4.0.0
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
max.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 1996-2015 John W. Eaton
4 Copyright (C) 2009 VZLU Prague
5 
6 This file is part of Octave.
7 
8 Octave is free software; you can redistribute it and/or modify it
9 under the terms of the GNU General Public License as published by the
10 Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12 
13 Octave is distributed in the hope that it will be useful, but WITHOUT
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
16 for more details.
17 
18 You should have received a copy of the GNU General Public License
19 along with Octave; see the file COPYING. If not, see
20 <http://www.gnu.org/licenses/>.
21 
22 */
23 
24 #ifdef HAVE_CONFIG_H
25 #include <config.h>
26 #endif
27 
28 #include "lo-ieee.h"
29 #include "lo-mappers.h"
30 #include "lo-math.h"
31 #include "dNDArray.h"
32 #include "CNDArray.h"
33 #include "quit.h"
34 
35 #include "defun.h"
36 #include "error.h"
37 #include "gripes.h"
38 #include "oct-obj.h"
39 
40 #include "ov-cx-mat.h"
41 #include "ov-re-sparse.h"
42 #include "ov-cx-sparse.h"
43 
44 template <class ArrayType>
45 static octave_value_list
47  int nargout, int dim, bool ismin)
48 {
49  octave_value_list retval;
50  ArrayType array = octave_value_extract<ArrayType> (arg);
51 
52  if (error_state)
53  return retval;
54 
55  if (nargout == 2)
56  {
57  retval.resize (2);
59  if (ismin)
60  retval(0) = array.min (idx, dim);
61  else
62  retval(0) = array.max (idx, dim);
63 
64  retval(1) = octave_value (idx, true, true);
65  }
66  else
67  {
68  if (ismin)
69  retval(0) = array.min (dim);
70  else
71  retval(0) = array.max (dim);
72  }
73 
74  return retval;
75 }
76 
77 // Matlab returns double arrays for min/max operations on character
78 // arrays, so we specialize here to get that behavior. Other possible
79 // solutions are to convert the argument to double here and call the
80 // code for double, but that could waste memory, or to have the
81 // underlying charNDArray::min/max functions return NDArray instead of
82 // charNDArray, but that is inconsistent with the way other min/max
83 // functions work.
84 
85 template <>
88  int nargout, int dim, bool ismin)
89 {
90  octave_value_list retval;
92 
93  if (error_state)
94  return retval;
95 
96  if (nargout == 2)
97  {
98  retval.resize (2);
100  if (ismin)
101  retval(0) = NDArray (array.min (idx, dim));
102  else
103  retval(0) = NDArray (array.max (idx, dim));
104 
105  retval(1) = octave_value (idx, true, true);
106  }
107  else
108  {
109  if (ismin)
110  retval(0) = NDArray (array.min (dim));
111  else
112  retval(0) = NDArray (array.max (dim));
113  }
114 
115  return retval;
116 }
117 
118 // Specialization for bool arrays.
119 template <>
122  int nargout, int dim, bool ismin)
123 {
124  octave_value_list retval;
125 
126  if (nargout <= 1)
127  {
128  // This case can be handled using any/all.
129  boolNDArray array = arg.bool_array_value ();
130 
131  if (array.is_empty ())
132  retval(0) = array;
133  else if (ismin)
134  retval(0) = array.all (dim);
135  else
136  retval(0) = array.any (dim);
137  }
138  else
139  {
140  // any/all don't have indexed versions, so do it via a conversion.
141  retval = do_minmax_red_op<int8NDArray> (arg, nargout, dim, ismin);
142  if (! error_state)
143  retval(0) = retval(0).bool_array_value ();
144  }
145 
146  return retval;
147 }
148 
149 template <class ArrayType>
150 static octave_value
151 do_minmax_bin_op (const octave_value& argx, const octave_value& argy,
152  bool ismin)
153 {
154  typedef typename ArrayType::element_type ScalarType;
155 
156  octave_value retval;
157 
158  if (argx.is_scalar_type ())
159  {
160  ScalarType x = octave_value_extract<ScalarType> (argx);
161  ArrayType y = octave_value_extract<ArrayType> (argy);
162 
163  if (error_state)
164  ;
165  else if (ismin)
166  retval = min (x, y);
167  else
168  retval = max (x, y);
169  }
170  else if (argy.is_scalar_type ())
171  {
172  ArrayType x = octave_value_extract<ArrayType> (argx);
173  ScalarType y = octave_value_extract<ScalarType> (argy);
174 
175  if (error_state)
176  ;
177  else if (ismin)
178  retval = min (x, y);
179  else
180  retval = max (x, y);
181  }
182  else
183  {
184  ArrayType x = octave_value_extract<ArrayType> (argx);
185  ArrayType y = octave_value_extract<ArrayType> (argy);
186 
187  if (error_state)
188  ;
189  else if (ismin)
190  retval = min (x, y);
191  else
192  retval = max (x, y);
193  }
194 
195  return retval;
196 }
197 
198 // Matlab returns double arrays for min/max operations on character
199 // arrays, so we specialize here to get that behavior. Other possible
200 // solutions are to convert the arguments to double here and call the
201 // code for double, but that could waste a lot of memory, or to have the
202 // underlying charNDArray::min/max functions return NDArray instead of
203 // charNDArray, but that is inconsistent with the way other min/max
204 // functions work.
205 
206 template <>
209  const octave_value& argy, bool ismin)
210 {
211  octave_value retval;
212 
215 
216  if (error_state)
217  ;
218  else if (ismin)
219  {
220  if (x.numel () == 1)
221  retval = NDArray (min (x(0), y));
222  else if (y.numel () == 1)
223  retval = NDArray (min (x, y(0)));
224  else
225  retval = NDArray (min (x, y));
226  }
227  else
228  {
229  if (x.numel () == 1)
230  retval = NDArray (max (x(0), y));
231  else if (y.numel () == 1)
232  retval = NDArray (max (x, y(0)));
233  else
234  retval = NDArray (max (x, y));
235  }
236 
237  return retval;
238 }
239 
240 static octave_value_list
242  int nargout, bool ismin)
243 {
244  octave_value_list retval;
245 
246  const char *func = ismin ? "min" : "max";
247 
248  int nargin = args.length ();
249 
250  if (nargin == 3 || nargin == 1)
251  {
252  octave_value arg = args(0);
253  int dim = -1;
254  if (nargin == 3)
255  {
256  dim = args(2).int_value (true) - 1;
257  if (error_state || dim < 0)
258  {
259  error ("%s: DIM must be a valid dimension", func);
260  return retval;
261  }
262 
263  if (! args(1).is_empty ())
264  warning ("%s: second argument is ignored", func);
265  }
266 
267  switch (arg.builtin_type ())
268  {
269  case btyp_double:
270  {
271  if (arg.is_range () && (dim == -1 || dim == 1))
272  {
273  Range range = arg.range_value ();
274  if (range.nelem () < 1)
275  {
276  retval(0) = arg;
277  if (nargout > 1)
278  retval(1) = arg;
279  }
280  else if (ismin)
281  {
282  retval(0) = range.min ();
283  if (nargout > 1)
284  retval(1) = static_cast<double>
285  (range.inc () < 0 ? range.nelem () : 1);
286  }
287  else
288  {
289  retval(0) = range.max ();
290  if (nargout > 1)
291  retval(1) = static_cast<double>
292  (range.inc () >= 0 ? range.nelem () : 1);
293  }
294  }
295  else if (arg.is_sparse_type ())
296  retval = do_minmax_red_op<SparseMatrix> (arg, nargout, dim,
297  ismin);
298  else
299  retval = do_minmax_red_op<NDArray> (arg, nargout, dim, ismin);
300  break;
301  }
302  case btyp_complex:
303  {
304  if (arg.is_sparse_type ())
305  retval = do_minmax_red_op<SparseComplexMatrix> (arg, nargout, dim,
306  ismin);
307  else
308  retval = do_minmax_red_op<ComplexNDArray> (arg, nargout, dim,
309  ismin);
310  break;
311  }
312  case btyp_float:
313  retval = do_minmax_red_op<FloatNDArray> (arg, nargout, dim, ismin);
314  break;
315  case btyp_float_complex:
316  retval = do_minmax_red_op<FloatComplexNDArray> (arg, nargout, dim,
317  ismin);
318  break;
319  case btyp_char:
320  retval = do_minmax_red_op<charNDArray> (arg, nargout, dim, ismin);
321  break;
322 #define MAKE_INT_BRANCH(X) \
323  case btyp_ ## X: \
324  retval = do_minmax_red_op<X ## NDArray> (arg, nargout, dim, ismin); \
325  break;
326  MAKE_INT_BRANCH (int8);
327  MAKE_INT_BRANCH (int16);
328  MAKE_INT_BRANCH (int32);
329  MAKE_INT_BRANCH (int64);
330  MAKE_INT_BRANCH (uint8);
331  MAKE_INT_BRANCH (uint16);
332  MAKE_INT_BRANCH (uint32);
333  MAKE_INT_BRANCH (uint64);
334 #undef MAKE_INT_BRANCH
335  case btyp_bool:
336  retval = do_minmax_red_op<boolNDArray> (arg, nargout, dim, ismin);
337  break;
338  default:
339  gripe_wrong_type_arg (func, arg);
340  }
341  }
342  else if (nargin == 2)
343  {
344  octave_value argx = args(0);
345  octave_value argy = args(1);
346  builtin_type_t xtyp = argx.builtin_type ();
347  builtin_type_t ytyp = argy.builtin_type ();
348  builtin_type_t rtyp;
349  if (xtyp == btyp_char && ytyp == btyp_char)
350  rtyp = btyp_char;
351  /*
352  FIXME: This is what should happen when boolNDArray has max()
353  else if (xtyp == btyp_bool && ytyp == btyp_bool)
354  rtyp = btyp_bool;
355  */
356  else
357  rtyp = btyp_mixed_numeric (xtyp, ytyp);
358 
359  switch (rtyp)
360  {
361  case btyp_double:
362  {
363  if ((argx.is_sparse_type ()
364  && (argy.is_sparse_type () || argy.is_scalar_type ()))
365  || (argy.is_sparse_type () && argx.is_scalar_type ()))
366  retval = do_minmax_bin_op<SparseMatrix> (argx, argy, ismin);
367  else
368  retval = do_minmax_bin_op<NDArray> (argx, argy, ismin);
369  break;
370  }
371  case btyp_complex:
372  {
373  if ((argx.is_sparse_type ()
374  && (argy.is_sparse_type () || argy.is_scalar_type ()))
375  || (argy.is_sparse_type () && argx.is_scalar_type ()))
376  retval = do_minmax_bin_op<SparseComplexMatrix> (argx, argy,
377  ismin);
378  else
379  retval = do_minmax_bin_op<ComplexNDArray> (argx, argy, ismin);
380  break;
381  }
382  case btyp_float:
383  retval = do_minmax_bin_op<FloatNDArray> (argx, argy, ismin);
384  break;
385  case btyp_float_complex:
386  retval = do_minmax_bin_op<FloatComplexNDArray> (argx, argy, ismin);
387  break;
388  case btyp_char:
389  retval = do_minmax_bin_op<charNDArray> (argx, argy, ismin);
390  break;
391 #define MAKE_INT_BRANCH(X) \
392  case btyp_ ## X: \
393  retval = do_minmax_bin_op<X ## NDArray> (argx, argy, ismin); \
394  break;
395  MAKE_INT_BRANCH (int8);
396  MAKE_INT_BRANCH (int16);
397  MAKE_INT_BRANCH (int32);
398  MAKE_INT_BRANCH (int64);
399  MAKE_INT_BRANCH (uint8);
400  MAKE_INT_BRANCH (uint16);
401  MAKE_INT_BRANCH (uint32);
402  MAKE_INT_BRANCH (uint64);
403 #undef MAKE_INT_BRANCH
404  /*
405  FIXME: This is what should happen when boolNDArray has max()
406  case btyp_bool:
407  retval = do_minmax_bin_op<boolNDArray> (argx, argy, ismin);
408  break;
409  */
410  default:
411  error ("%s: cannot compute %s (%s, %s)", func, func,
412  argx.type_name ().c_str (), argy.type_name ().c_str ());
413  }
414 
415  // FIXME: Delete when boolNDArray has max()
416  if (xtyp == btyp_bool && ytyp == btyp_bool)
417  retval(0) = retval(0).bool_array_value ();
418 
419  }
420  else
421  print_usage ();
422 
423  return retval;
424 }
425 
426 DEFUN (min, args, nargout,
427  "-*- texinfo -*-\n\
428 @deftypefn {Built-in Function} {} min (@var{x})\n\
429 @deftypefnx {Built-in Function} {} min (@var{x}, [], @var{dim})\n\
430 @deftypefnx {Built-in Function} {[@var{w}, @var{iw}] =} min (@var{x})\n\
431 @deftypefnx {Built-in Function} {} min (@var{x}, @var{y})\n\
432 Find minimum values in the array @var{x}.\n\
433 \n\
434 For a vector argument, return the minimum value. For a matrix argument,\n\
435 return a row vector with the minimum value of each column. For a\n\
436 multi-dimensional array, @code{min} operates along the first non-singleton\n\
437 dimension.\n\
438 \n\
439 If the optional third argument @var{dim} is present then operate along\n\
440 this dimension. In this case the second argument is ignored and should be\n\
441 set to the empty matrix.\n\
442 \n\
443 For two matrices (or a matrix and a scalar), return the pairwise minimum.\n\
444 \n\
445 Thus,\n\
446 \n\
447 @example\n\
448 min (min (@var{x}))\n\
449 @end example\n\
450 \n\
451 @noindent\n\
452 returns the smallest element of the 2-D matrix @var{x}, and\n\
453 \n\
454 @example\n\
455 @group\n\
456 min (2:5, pi)\n\
457  @result{} 2.0000 3.0000 3.1416 3.1416\n\
458 @end group\n\
459 @end example\n\
460 \n\
461 @noindent\n\
462 compares each element of the range @code{2:5} with @code{pi}, and returns a\n\
463 row vector of the minimum values.\n\
464 \n\
465 For complex arguments, the magnitude of the elements are used for\n\
466 comparison. If the magnitudes are identical, then the results are ordered\n\
467 by phase angle in the range (-pi, pi]. Hence,\n\
468 \n\
469 @example\n\
470 @group\n\
471 min ([-1 i 1 -i])\n\
472  @result{} -i\n\
473 @end group\n\
474 @end example\n\
475 \n\
476 @noindent\n\
477 because all entries have magnitude 1, but -i has the smallest phase angle\n\
478 with value -pi/2.\n\
479 \n\
480 If called with one input and two output arguments, @code{min} also returns\n\
481 the first index of the minimum value(s). Thus,\n\
482 \n\
483 @example\n\
484 @group\n\
485 [x, ix] = min ([1, 3, 0, 2, 0])\n\
486  @result{} x = 0\n\
487  ix = 3\n\
488 @end group\n\
489 @end example\n\
490 @seealso{max, cummin, cummax}\n\
491 @end deftypefn")
492 {
493  return do_minmax_body (args, nargout, true);
494 }
495 
496 /*
497 ## Test generic double class
498 %!assert (min ([1, 4, 2, 3]), 1)
499 %!assert (min ([1; -10; 5; -2]), -10)
500 %!assert (min ([4, 2i 4.999; -2, 2, 3+4i]), [-2, 2, 4.999])
501 ## Special routines for char arrays
502 %!assert (min (["abc", "ABC"]), 65)
503 %!assert (min (["abc"; "CBA"]), [67 66 65])
504 ## Special routines for logical arrays
505 %!assert (min (logical ([])), logical ([]))
506 %!assert (min (logical ([0 0 1 0])), false)
507 %!assert (min (logical ([0 0 1 0; 0 1 1 0])), logical ([0 0 1 0]))
508 ## Single values
509 %!assert (min (single ([1, 4, 2, 3])), single (1))
510 %!assert (min (single ([1; -10; 5; -2])), single (-10))
511 %!assert (min (single ([4, 2i 4.999; -2, 2, 3+4i])), single ([-2, 2, 4.999]))
512 ## Integer values
513 %!assert (min (uint8 ([1, 4, 2, 3])), uint8 (1))
514 %!assert (min (uint8 ([1; -10; 5; -2])), uint8 (-10))
515 %!assert (min (int8 ([1, 4, 2, 3])), int8 (1))
516 %!assert (min (int8 ([1; -10; 5; -2])), int8 (-10))
517 %!assert (min (uint16 ([1, 4, 2, 3])), uint16 (1))
518 %!assert (min (uint16 ([1; -10; 5; -2])), uint16 (-10))
519 %!assert (min (int16 ([1, 4, 2, 3])), int16 (1))
520 %!assert (min (int16 ([1; -10; 5; -2])), int16 (-10))
521 %!assert (min (uint32 ([1, 4, 2, 3])), uint32 (1))
522 %!assert (min (uint32 ([1; -10; 5; -2])), uint32 (-10))
523 %!assert (min (int32 ([1, 4, 2, 3])), int32 (1))
524 %!assert (min (int32 ([1; -10; 5; -2])), int32 (-10))
525 %!assert (min (uint64 ([1, 4, 2, 3])), uint64 (1))
526 %!assert (min (uint64 ([1; -10; 5; -2])), uint64 (-10))
527 %!assert (min (int64 ([1, 4, 2, 3])), int64 (1))
528 %!assert (min (int64 ([1; -10; 5; -2])), int64 (-10))
529 ## Sparse double values
530 %!assert (min (sparse ([1, 4, 2, 3])), sparse (1))
531 %!assert (min (sparse ([1; -10; 5; -2])), sparse(-10))
532 ## FIXME: sparse doesn't order complex values by phase angle
533 %!xtest
534 %! assert (min (sparse ([4, 2i 4.999; -2, 2, 3+4i])), sparse ([-2, 2, 4.999]));
535 
536 ## Test dimension argument
537 %!test
538 %! x = reshape (1:8, [2,2,2]);
539 %! assert (min (x, [], 1), reshape ([1, 3, 5, 7], [1,2,2]));
540 %! assert (min (x, [], 2), reshape ([1, 2, 5, 6], [2,1,2]));
541 %! [y, i] = min (x, [], 3);
542 %! assert (ndims (y), 2);
543 %! assert (y, [1, 3; 2, 4]);
544 %! assert (ndims (i), 2);
545 %! assert (i, [1, 1; 1, 1]);
546 
547 ## Test 2-output forms for various arg types
548 ## Special routines for char arrays
549 %!test
550 %! [y, i] = min (["abc", "ABC"]);
551 %! assert (y, 65);
552 %! assert (i, 4);
553 ## Special routines for logical arrays
554 %!test
555 %! x = logical ([0 0 1 0]);
556 %! [y, i] = min (x);
557 %! assert (y, false);
558 %! assert (i, 1);
559 ## Special handling of ranges
560 %!test
561 %! rng = 1:2:10;
562 %! [y, i] = min (rng);
563 %! assert (y, 1);
564 %! assert (i, 1);
565 %! rng = 10:-2:1;
566 %! [y, i] = min (rng);
567 %! assert (y, 2);
568 %! assert (i, 5);
569 
570 ## Test 2-input calling form for various arg types
571 ## Test generic double class
572 %!test
573 %! x = [1, 2, 3, 4]; y = fliplr (x);
574 %! assert (min (x, y), [1 2 2 1]);
575 %! assert (min (x, 3), [1 2 3 3]);
576 %! assert (min (2, x), [1 2 2 2]);
577 %! assert (min (x, 2.1i), [1 2 2.1i 2.1i]);
578 ## FIXME: Ordering of complex results with equal magnitude is not by phase
579 ## angle in the 2-input form. Instead, it is in the order in which it
580 ## appears in the argument list.
581 %!xtest
582 %! x = [1, 2, 3, 4]; y = fliplr (x);
583 %! assert (min (x, 2i), [2i 2i 3 4]);
584 ## Special routines for char arrays
585 %!assert (min ("abc", "b"), [97 98 98])
586 %!assert (min ("b", "cba"), [98 98 97])
587 ## Special handling for logical arrays
588 %!assert (min ([true false], false), [false false])
589 %!assert (min (true, [true false]), [true false])
590 ## Single values
591 %!test
592 %! x = single ([1, 2, 3, 4]); y = fliplr (x);
593 %! assert (min (x, y), single ([1 2 2 1]));
594 %! assert (min (x, 3), single ([1 2 3 3]));
595 %! assert (min (2, x), single ([1 2 2 2]));
596 %! assert (min (x, 2.1i), single ([1 2 2.1i 2.1i]));
597 ## Integer values
598 %!test
599 %! x = uint8 ([1, 2, 3, 4]); y = fliplr (x);
600 %! assert (min (x, y), uint8 ([1 2 2 1]));
601 %! assert (min (x, 3), uint8 ([1 2 3 3]));
602 %! assert (min (2, x), uint8 ([1 2 2 2]));
603 %! x = int8 ([1, 2, 3, 4]); y = fliplr (x);
604 %! assert (min (x, y), int8 ([1 2 2 1]));
605 %! assert (min (x, 3), int8 ([1 2 3 3]));
606 %! assert (min (2, x), int8 ([1 2 2 2]));
607 %! x = uint16 ([1, 2, 3, 4]); y = fliplr (x);
608 %! assert (min (x, y), uint16 ([1 2 2 1]));
609 %! assert (min (x, 3), uint16 ([1 2 3 3]));
610 %! assert (min (2, x), uint16 ([1 2 2 2]));
611 %! x = int16 ([1, 2, 3, 4]); y = fliplr (x);
612 %! assert (min (x, y), int16 ([1 2 2 1]));
613 %! assert (min (x, 3), int16 ([1 2 3 3]));
614 %! assert (min (2, x), int16 ([1 2 2 2]));
615 %! x = uint32 ([1, 2, 3, 4]); y = fliplr (x);
616 %! assert (min (x, y), uint32 ([1 2 2 1]));
617 %! assert (min (x, 3), uint32 ([1 2 3 3]));
618 %! assert (min (2, x), uint32 ([1 2 2 2]));
619 %! x = int32 ([1, 2, 3, 4]); y = fliplr (x);
620 %! assert (min (x, y), int32 ([1 2 2 1]));
621 %! assert (min (x, 3), int32 ([1 2 3 3]));
622 %! assert (min (2, x), int32 ([1 2 2 2]));
623 %! x = uint64 ([1, 2, 3, 4]); y = fliplr (x);
624 %! assert (min (x, y), uint64 ([1 2 2 1]));
625 %! assert (min (x, 3), uint64 ([1 2 3 3]));
626 %! assert (min (2, x), uint64 ([1 2 2 2]));
627 %! x = int64 ([1, 2, 3, 4]); y = fliplr (x);
628 %! assert (min (x, y), int64 ([1 2 2 1]));
629 %! assert (min (x, 3), int64 ([1 2 3 3]));
630 %! assert (min (2, x), int64 ([1 2 2 2]));
631 ## Sparse double values
632 %!test
633 %! x = sparse ([1, 2, 3, 4]); y = fliplr (x);
634 %! assert (min (x, y), sparse ([1 2 2 1]));
635 %! assert (min (x, 3), sparse ([1 2 3 3]));
636 %! assert (min (2, x), sparse ([1 2 2 2]));
637 %! assert (min (x, 2.1i), sparse ([1 2 2.1i 2.1i]));
638 
639 %!error min ()
640 %!error min (1, 2, 3, 4)
641 %!error <DIM must be a valid dimension> min ([1 2; 3 4], [], -3)
642 %!warning <second argument is ignored> min ([1 2 3 4], 2, 2);
643 %!error <wrong type argument 'cell'> min ({1 2 3 4})
644 %!error <cannot compute min \(cell, scalar\)> min ({1, 2, 3}, 2)
645 */
646 
647 DEFUN (max, args, nargout,
648  "-*- texinfo -*-\n\
649 @deftypefn {Built-in Function} {} max (@var{x})\n\
650 @deftypefnx {Built-in Function} {} max (@var{x}, [], @var{dim})\n\
651 @deftypefnx {Built-in Function} {[@var{w}, @var{iw}] =} max (@var{x})\n\
652 @deftypefnx {Built-in Function} {} max (@var{x}, @var{y})\n\
653 Find maximum values in the array @var{x}.\n\
654 \n\
655 For a vector argument, return the maximum value. For a matrix argument,\n\
656 return a row vector with the maximum value of each column. For a\n\
657 multi-dimensional array, @code{max} operates along the first non-singleton\n\
658 dimension.\n\
659 \n\
660 If the optional third argument @var{dim} is present then operate along\n\
661 this dimension. In this case the second argument is ignored and should be\n\
662 set to the empty matrix.\n\
663 \n\
664 For two matrices (or a matrix and a scalar), return the pairwise maximum.\n\
665 \n\
666 Thus,\n\
667 \n\
668 @example\n\
669 max (max (@var{x}))\n\
670 @end example\n\
671 \n\
672 @noindent\n\
673 returns the largest element of the 2-D matrix @var{x}, and\n\
674 \n\
675 @example\n\
676 @group\n\
677 max (2:5, pi)\n\
678  @result{} 3.1416 3.1416 4.0000 5.0000\n\
679 @end group\n\
680 @end example\n\
681 \n\
682 @noindent\n\
683 compares each element of the range @code{2:5} with @code{pi}, and returns a\n\
684 row vector of the maximum values.\n\
685 \n\
686 For complex arguments, the magnitude of the elements are used for\n\
687 comparison. If the magnitudes are identical, then the results are ordered\n\
688 by phase angle in the range (-pi, pi]. Hence,\n\
689 \n\
690 @example\n\
691 @group\n\
692 max ([-1 i 1 -i])\n\
693  @result{} -1\n\
694 @end group\n\
695 @end example\n\
696 \n\
697 @noindent\n\
698 because all entries have magnitude 1, but -1 has the largest phase angle\n\
699 with value pi.\n\
700 \n\
701 If called with one input and two output arguments, @code{max} also returns\n\
702 the first index of the maximum value(s). Thus,\n\
703 \n\
704 @example\n\
705 @group\n\
706 [x, ix] = max ([1, 3, 5, 2, 5])\n\
707  @result{} x = 5\n\
708  ix = 3\n\
709 @end group\n\
710 @end example\n\
711 @seealso{min, cummax, cummin}\n\
712 @end deftypefn")
713 {
714  return do_minmax_body (args, nargout, false);
715 }
716 
717 /*
718 ## Test generic double class
719 %!assert (max ([1, 4, 2, 3]), 4)
720 %!assert (max ([1; -10; 5; -2]), 5)
721 %!assert (max ([4, 2i 4.999; -2, 2, 3+4i]), [4, 2i, 3+4i])
722 ## Special routines for char arrays
723 %!assert (max (["abc", "ABC"]), 99)
724 %!assert (max (["abc"; "CBA"]), [97 98 99])
725 ## Special routines for logical arrays
726 %!assert (max (logical ([])), logical ([]))
727 %!assert (max (logical ([0 0 1 0])), true)
728 %!assert (max (logical ([0 0 1 0; 0 1 0 0])), logical ([0 1 1 0]))
729 ## Single values
730 %!assert (max (single ([1, 4, 2, 3])), single (4))
731 %!assert (max (single ([1; -10; 5; -2])), single (5))
732 %!assert (max (single ([4, 2i 4.999; -2, 2, 3+4i])), single ([4, 2i, 3+4i]))
733 ## Integer values
734 %!assert (max (uint8 ([1, 4, 2, 3])), uint8 (4))
735 %!assert (max (uint8 ([1; -10; 5; -2])), uint8 (5))
736 %!assert (max (int8 ([1, 4, 2, 3])), int8 (4))
737 %!assert (max (int8 ([1; -10; 5; -2])), int8 (5))
738 %!assert (max (uint16 ([1, 4, 2, 3])), uint16 (4))
739 %!assert (max (uint16 ([1; -10; 5; -2])), uint16 (5))
740 %!assert (max (int16 ([1, 4, 2, 3])), int16 (4))
741 %!assert (max (int16 ([1; -10; 5; -2])), int16 (5))
742 %!assert (max (uint32 ([1, 4, 2, 3])), uint32 (4))
743 %!assert (max (uint32 ([1; -10; 5; -2])), uint32 (5))
744 %!assert (max (int32 ([1, 4, 2, 3])), int32 (4))
745 %!assert (max (int32 ([1; -10; 5; -2])), int32 (5))
746 %!assert (max (uint64 ([1, 4, 2, 3])), uint64 (4))
747 %!assert (max (uint64 ([1; -10; 5; -2])), uint64 (5))
748 %!assert (max (int64 ([1, 4, 2, 3])), int64 (4))
749 %!assert (max (int64 ([1; -10; 5; -2])), int64 (5))
750 ## Sparse double values
751 %!assert (max (sparse ([1, 4, 2, 3])), sparse (4))
752 %!assert (max (sparse ([1; -10; 5; -2])), sparse(5))
753 %!assert (max (sparse ([4, 2i 4.999; -2, 2, 3+4i])), sparse ([4, 2i, 3+4i]))
754 
755 ## Test dimension argument
756 %!test
757 %! x = reshape (1:8, [2,2,2]);
758 %! assert (min (x, [], 1), reshape ([1, 3, 5, 7], [1,2,2]));
759 %! assert (min (x, [], 2), reshape ([1, 2, 5, 6], [2,1,2]));
760 %! [y, i] = min (x, [], 3);
761 %! assert (ndims (y), 2);
762 %! assert (y, [1, 3; 2, 4]);
763 %! assert (ndims (i), 2);
764 %! assert (i, [1, 1; 1, 1]);
765 
766 ## Test 2-output forms for various arg types
767 ## Special routines for char arrays
768 %!test
769 %! [y, i] = max (["abc", "ABC"]);
770 %! assert (y, 99);
771 %! assert (i, 3);
772 ## Special routines for logical arrays
773 %!test
774 %! x = logical ([0 0 1 0]);
775 %! [y, i] = max (x);
776 %! assert (y, true);
777 %! assert (i, 3);
778 ## Special handling of ranges
779 %!test
780 %! rng = 1:2:10;
781 %! [y, i] = max (rng);
782 %! assert (y, 9);
783 %! assert (i, 5);
784 %! rng = 10:-2:1;
785 %! [y, i] = max (rng);
786 %! assert (y, 10);
787 %! assert (i, 1);
788 
789 ## Test 2-input calling form for various arg types
790 ## Test generic double class
791 %!test
792 %! x = [1, 2, 3, 4]; y = fliplr (x);
793 %! assert (max (x, y), [4 3 3 4]);
794 %! assert (max (x, 3), [3 3 3 4]);
795 %! assert (max (2, x), [2 2 3 4]);
796 %! assert (max (x, 2.1i), [2.1i 2.1i 3 4]);
797 ## FIXME: Ordering of complex results with equal magnitude is not by phase
798 ## angle in the 2-input form. Instead, it is in the order in which it
799 ## appears in the argument list.
800 %!xtest
801 %! x = [1, 2, 3, 4]; y = fliplr (x);
802 %! assert (max (x, 2i), [2i 2i 3 4]);
803 ## Special routines for char arrays
804 %!assert (max ("abc", "b"), [98 98 99])
805 %!assert (max ("b", "cba"), [99 98 98])
806 ## Special handling for logical arrays
807 %!assert (max ([true false], false), [true false])
808 %!assert (max (true, [false false]), [true true])
809 ## Single values
810 %!test
811 %! x = single ([1, 2, 3, 4]); y = fliplr (x);
812 %! assert (max (x, y), single ([4 3 3 4]));
813 %! assert (max (x, 3), single ([3 3 3 4]));
814 %! assert (max (2, x), single ([2 2 3 4]));
815 %! assert (max (x, 2.1i), single ([2.1i 2.1i 3 4]));
816 ## Integer values
817 %!test
818 %! x = uint8 ([1, 2, 3, 4]); y = fliplr (x);
819 %! assert (max (x, y), uint8 ([4 3 3 4]));
820 %! assert (max (x, 3), uint8 ([3 3 3 4]));
821 %! assert (max (2, x), uint8 ([2 2 3 4]));
822 %! x = int8 ([1, 2, 3, 4]); y = fliplr (x);
823 %! assert (max (x, y), int8 ([4 3 3 4]));
824 %! assert (max (x, 3), int8 ([3 3 3 4]));
825 %! assert (max (2, x), int8 ([2 2 3 4]));
826 %! x = uint16 ([1, 2, 3, 4]); y = fliplr (x);
827 %! assert (max (x, y), uint16 ([4 3 3 4]));
828 %! assert (max (x, 3), uint16 ([3 3 3 4]));
829 %! assert (max (2, x), uint16 ([2 2 3 4]));
830 %! x = int16 ([1, 2, 3, 4]); y = fliplr (x);
831 %! assert (max (x, y), int16 ([4 3 3 4]));
832 %! assert (max (x, 3), int16 ([3 3 3 4]));
833 %! assert (max (2, x), int16 ([2 2 3 4]));
834 %! x = uint32 ([1, 2, 3, 4]); y = fliplr (x);
835 %! assert (max (x, y), uint32 ([4 3 3 4]));
836 %! assert (max (x, 3), uint32 ([3 3 3 4]));
837 %! assert (max (2, x), uint32 ([2 2 3 4]));
838 %! x = int32 ([1, 2, 3, 4]); y = fliplr (x);
839 %! assert (max (x, y), int32 ([4 3 3 4]));
840 %! assert (max (x, 3), int32 ([3 3 3 4]));
841 %! assert (max (2, x), int32 ([2 2 3 4]));
842 %! x = uint64 ([1, 2, 3, 4]); y = fliplr (x);
843 %! assert (max (x, y), uint64 ([4 3 3 4]));
844 %! assert (max (x, 3), uint64 ([3 3 3 4]));
845 %! assert (max (2, x), uint64 ([2 2 3 4]));
846 %! x = int64 ([1, 2, 3, 4]); y = fliplr (x);
847 %! assert (max (x, y), int64 ([4 3 3 4]));
848 %! assert (max (x, 3), int64 ([3 3 3 4]));
849 %! assert (max (2, x), int64 ([2 2 3 4]));
850 ## Sparse double values
851 %!test
852 %! x = sparse ([1, 2, 3, 4]); y = fliplr (x);
853 %! assert (max (x, y), sparse ([4 3 3 4]));
854 %! assert (max (x, 3), sparse ([3 3 3 4]));
855 %! assert (max (2, x), sparse ([2 2 3 4]));
856 %! assert (max (x, 2.1i), sparse ([2.1i 2.1i 3 4]));
857 
858 ## Test for bug #40743
859 %!assert (max (zeros (1,0), ones (1,1)), zeros (1,0))
860 %!assert (max (sparse (zeros (1,0)), sparse (ones (1,1))), sparse (zeros (1,0)))
861 
862 %!error max ()
863 %!error max (1, 2, 3, 4)
864 %!error <DIM must be a valid dimension> max ([1 2; 3 4], [], -3)
865 %!warning <second argument is ignored> max ([1 2 3 4], 2, 2);
866 %!error <wrong type argument 'cell'> max ({1 2 3 4})
867 %!error <cannot compute max \(cell, scalar\)> max ({1, 2, 3}, 2)
868 
869 */
870 
871 template <class ArrayType>
872 static octave_value_list
874  int nargout, int dim, bool ismin)
875 {
876  octave_value_list retval;
877  ArrayType array = octave_value_extract<ArrayType> (arg);
878 
879  if (error_state)
880  return retval;
881 
882  if (nargout == 2)
883  {
884  retval.resize (2);
886  if (ismin)
887  retval(0) = array.cummin (idx, dim);
888  else
889  retval(0) = array.cummax (idx, dim);
890 
891  retval(1) = octave_value (idx, true, true);
892  }
893  else
894  {
895  if (ismin)
896  retval(0) = array.cummin (dim);
897  else
898  retval(0) = array.cummax (dim);
899  }
900 
901  return retval;
902 }
903 
904 static octave_value_list
906  int nargout, bool ismin)
907 {
908  octave_value_list retval;
909 
910  const char *func = ismin ? "cummin" : "cummax";
911 
912  int nargin = args.length ();
913 
914  if (nargin == 1 || nargin == 2)
915  {
916  octave_value arg = args(0);
917  int dim = -1;
918  if (nargin == 2)
919  {
920  dim = args(1).int_value (true) - 1;
921  if (error_state || dim < 0)
922  {
923  error ("%s: DIM must be a valid dimension", func);
924  return retval;
925  }
926  }
927 
928  switch (arg.builtin_type ())
929  {
930  case btyp_double:
931  retval = do_cumminmax_red_op<NDArray> (arg, nargout, dim, ismin);
932  break;
933  case btyp_complex:
934  retval = do_cumminmax_red_op<ComplexNDArray> (arg, nargout, dim,
935  ismin);
936  break;
937  case btyp_float:
938  retval = do_cumminmax_red_op<FloatNDArray> (arg, nargout, dim, ismin);
939  break;
940  case btyp_float_complex:
941  retval = do_cumminmax_red_op<FloatComplexNDArray> (arg, nargout, dim,
942  ismin);
943  break;
944 #define MAKE_INT_BRANCH(X) \
945  case btyp_ ## X: \
946  retval = do_cumminmax_red_op<X ## NDArray> (arg, nargout, dim, \
947  ismin); \
948  break;
949  MAKE_INT_BRANCH (int8);
950  MAKE_INT_BRANCH (int16);
951  MAKE_INT_BRANCH (int32);
952  MAKE_INT_BRANCH (int64);
953  MAKE_INT_BRANCH (uint8);
954  MAKE_INT_BRANCH (uint16);
955  MAKE_INT_BRANCH (uint32);
956  MAKE_INT_BRANCH (uint64);
957 #undef MAKE_INT_BRANCH
958  case btyp_bool:
959  {
960  retval = do_cumminmax_red_op<int8NDArray> (arg, nargout, dim,
961  ismin);
962  if (retval.length () > 0)
963  retval(0) = retval(0).bool_array_value ();
964  break;
965  }
966  default:
967  gripe_wrong_type_arg (func, arg);
968  }
969  }
970  else
971  print_usage ();
972 
973  return retval;
974 }
975 
976 DEFUN (cummin, args, nargout,
977  "-*- texinfo -*-\n\
978 @deftypefn {Built-in Function} {} cummin (@var{x})\n\
979 @deftypefnx {Built-in Function} {} cummin (@var{x}, @var{dim})\n\
980 @deftypefnx {Built-in Function} {[@var{w}, @var{iw}] =} cummin (@var{x})\n\
981 Return the cumulative minimum values along dimension @var{dim}.\n\
982 \n\
983 If @var{dim} is unspecified it defaults to column-wise operation. For\n\
984 example:\n\
985 \n\
986 @example\n\
987 @group\n\
988 cummin ([5 4 6 2 3 1])\n\
989  @result{} 5 4 4 2 2 1\n\
990 @end group\n\
991 @end example\n\
992 \n\
993 If called with two output arguments the index of the minimum value is also\n\
994 returned.\n\
995 \n\
996 @example\n\
997 @group\n\
998 [w, iw] = cummin ([5 4 6 2 3 1])\n\
999 @result{}\n\
1000 w = 5 4 4 2 2 1\n\
1001 iw = 1 2 2 4 4 6\n\
1002 @end group\n\
1003 @end example\n\
1004 \n\
1005 @seealso{cummax, min, max}\n\
1006 @end deftypefn")
1007 {
1008  return do_cumminmax_body (args, nargout, true);
1009 }
1010 
1011 /*
1012 %!assert (cummin ([1, 4, 2, 3]), [1 1 1 1])
1013 %!assert (cummin ([1; -10; 5; -2]), [1; -10; -10; -10])
1014 %!assert (cummin ([4, i; -2, 2]), [4, i; -2, i])
1015 %!assert (cummin ([1 2; NaN 1], 2), [1 1; NaN 1])
1016 
1017 %!test
1018 %! x = reshape (1:8, [2,2,2]);
1019 %! assert (cummin (x, 1), reshape ([1 1 3 3 5 5 7 7], [2,2,2]));
1020 %! assert (cummin (x, 2), reshape ([1 2 1 2 5 6 5 6], [2,2,2]));
1021 %! [w, iw] = cummin (x, 3);
1022 %! assert (ndims (w), 3);
1023 %! assert (w, repmat ([1 3; 2 4], [1 1 2]));
1024 %! assert (ndims (iw), 3);
1025 %! assert (iw, ones (2,2,2));
1026 
1027 
1028 %!error cummin ()
1029 %!error cummin (1, 2, 3)
1030 */
1031 
1032 DEFUN (cummax, args, nargout,
1033  "-*- texinfo -*-\n\
1034 @deftypefn {Built-in Function} {} cummax (@var{x})\n\
1035 @deftypefnx {Built-in Function} {} cummax (@var{x}, @var{dim})\n\
1036 @deftypefnx {Built-in Function} {[@var{w}, @var{iw}] =} cummax (@dots{})\n\
1037 Return the cumulative maximum values along dimension @var{dim}.\n\
1038 \n\
1039 If @var{dim} is unspecified it defaults to column-wise operation. For\n\
1040 example:\n\
1041 \n\
1042 @example\n\
1043 @group\n\
1044 cummax ([1 3 2 6 4 5])\n\
1045  @result{} 1 3 3 6 6 6\n\
1046 @end group\n\
1047 @end example\n\
1048 \n\
1049 If called with two output arguments the index of the maximum value is also\n\
1050 returned.\n\
1051 \n\
1052 @example\n\
1053 @group\n\
1054 [w, iw] = cummax ([1 3 2 6 4 5])\n\
1055 @result{}\n\
1056 w = 1 3 3 6 6 6\n\
1057 iw = 1 2 2 4 4 4\n\
1058 @end group\n\
1059 @end example\n\
1060 \n\
1061 @seealso{cummin, max, min}\n\
1062 @end deftypefn")
1063 {
1064  return do_cumminmax_body (args, nargout, false);
1065 }
1066 
1067 /*
1068 %!assert (cummax ([1, 4, 2, 3]), [1 4 4 4])
1069 %!assert (cummax ([1; -10; 5; -2]), [1; 1; 5; 5])
1070 %!assert (cummax ([4, i 4.9, -2, 2, 3+4i]), [4, 4, 4.9, 4.9, 4.9, 3+4i])
1071 %!assert (cummax ([1 NaN 0; NaN NaN 1], 2), [1 1 1; NaN NaN 1])
1072 
1073 %!test
1074 %! x = reshape (8:-1:1, [2,2,2]);
1075 %! assert (cummax (x, 1), reshape ([8 8 6 6 4 4 2 2], [2,2,2]));
1076 %! assert (cummax (x, 2), reshape ([8 7 8 7 4 3 4 3], [2,2,2]));
1077 %! [w, iw] = cummax (x, 3);
1078 %! assert (ndims (w), 3);
1079 %! assert (w, repmat ([8 6; 7 5], [1 1 2]));
1080 %! assert (ndims (iw), 3);
1081 %! assert (iw, ones (2,2,2));
1082 
1083 %!error cummax ()
1084 %!error cummax (1, 2, 3)
1085 */
1086 
boolNDArray all(int dim=-1) const
Definition: boolNDArray.cc:61
charNDArray max(int dim=-1) const
Definition: chNDArray.cc:145
bool is_empty(void) const
Definition: Array.h:472
bool is_range(void) const
Definition: ov.h:571
void gripe_wrong_type_arg(const char *name, const char *s, bool is_error)
Definition: gripes.cc:135
OCTINTERP_API void print_usage(void)
Definition: defun.cc:51
octave_idx_type numel(void) const
Number of elements in the array.
Definition: Array.h:275
octave_idx_type length(void) const
Definition: oct-obj.h:89
octave_idx_type nelem(void) const
Definition: Range.h:64
bool is_scalar_type(void) const
Definition: ov.h:657
Definition: Range.h:31
#define DEFUN(name, args_name, nargout_name, doc)
Definition: defun.h:44
void error(const char *fmt,...)
Definition: error.cc:476
boolNDArray any(int dim=-1) const
Definition: boolNDArray.cc:67
double min(void) const
Definition: Range.cc:194
static octave_value_list do_cumminmax_body(const octave_value_list &args, int nargout, bool ismin)
Definition: max.cc:905
builtin_type_t
Definition: ov-base.h:59
octave_value_list do_minmax_red_op< boolNDArray >(const octave_value &arg, int nargout, int dim, bool ismin)
Definition: max.cc:121
builtin_type_t btyp_mixed_numeric(builtin_type_t x, builtin_type_t y)
Definition: ov-base.cc:59
double max(void) const
Definition: Range.cc:215
bool is_sparse_type(void) const
Definition: ov.h:666
Range range_value(void) const
Definition: ov.h:903
static octave_value_list do_minmax_body(const octave_value_list &args, int nargout, bool ismin)
Definition: max.cc:241
double inc(void) const
Definition: Range.h:63
static octave_value_list do_cumminmax_red_op(const octave_value &arg, int nargout, int dim, bool ismin)
Definition: max.cc:873
static octave_value do_minmax_bin_op(const octave_value &argx, const octave_value &argy, bool ismin)
Definition: max.cc:151
int error_state
Definition: error.cc:101
charNDArray octave_value_extract< charNDArray >(const octave_value &v)
Definition: ov.h:1401
static octave_value_list do_minmax_red_op(const octave_value &arg, int nargout, int dim, bool ismin)
Definition: max.cc:46
#define MAKE_INT_BRANCH(X)
double arg(double x)
Definition: lo-mappers.h:37
void warning(const char *fmt,...)
Definition: error.cc:681
charNDArray max(char d, const charNDArray &m)
Definition: chNDArray.cc:233
charNDArray min(int dim=-1) const
Definition: chNDArray.cc:157
void resize(octave_idx_type n, const octave_value &rfv=octave_value())
Definition: oct-obj.h:93
octave_value do_minmax_bin_op< charNDArray >(const octave_value &argx, const octave_value &argy, bool ismin)
Definition: max.cc:208
builtin_type_t builtin_type(void) const
Definition: ov.h:603
octave_value_list do_minmax_red_op< charNDArray >(const octave_value &arg, int nargout, int dim, bool ismin)
Definition: max.cc:87
return octave_value(v1.char_array_value().concat(v2.char_array_value(), ra_idx),((a1.is_sq_string()||a2.is_sq_string())? '\'': '"'))
F77_RET_T const double * x
charNDArray min(char d, const charNDArray &m)
Definition: chNDArray.cc:210