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
dNDArray.cc
Go to the documentation of this file.
1 // N-D Array manipulations.
2 /*
3 
4 Copyright (C) 1996-2015 John W. Eaton
5 Copyright (C) 2009 VZLU Prague, a.s.
6 
7 This file is part of Octave.
8 
9 Octave is free software; you can redistribute it and/or modify it
10 under the terms of the GNU General Public License as published by the
11 Free Software Foundation; either version 3 of the License, or (at your
12 option) any later version.
13 
14 Octave is distributed in the hope that it will be useful, but WITHOUT
15 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
17 for more details.
18 
19 You should have received a copy of the GNU General Public License
20 along with Octave; see the file COPYING. If not, see
21 <http://www.gnu.org/licenses/>.
22 
23 */
24 
25 #ifdef HAVE_CONFIG_H
26 #include <config.h>
27 #endif
28 
29 #include <cfloat>
30 
31 #include <vector>
32 
33 #include "Array-util.h"
34 #include "dNDArray.h"
35 #include "f77-fcn.h"
36 #include "functor.h"
37 #include "lo-error.h"
38 #include "lo-ieee.h"
39 #include "lo-mappers.h"
40 #include "mx-base.h"
41 #include "mx-op-defs.h"
42 #include "oct-fftw.h"
43 #include "oct-locbuf.h"
44 
45 #include "bsxfun-defs.cc"
46 
47 NDArray::NDArray (const Array<octave_idx_type>& a, bool zero_based,
48  bool negative_to_nan)
49 {
50  const octave_idx_type *pa = a.fortran_vec ();
51  resize (a.dims ());
52  double *ptmp = fortran_vec ();
53  if (negative_to_nan)
54  {
55  double nan_val = lo_ieee_nan_value ();
56 
57  if (zero_based)
58  for (octave_idx_type i = 0; i < a.numel (); i++)
59  {
60  double val = static_cast<double>
61  (pa[i] + static_cast<octave_idx_type> (1));
62  if (val <= 0)
63  ptmp[i] = nan_val;
64  else
65  ptmp[i] = val;
66  }
67  else
68  for (octave_idx_type i = 0; i < a.numel (); i++)
69  {
70  double val = static_cast<double> (pa[i]);
71  if (val <= 0)
72  ptmp[i] = nan_val;
73  else
74  ptmp[i] = val;
75  }
76  }
77  else
78  {
79  if (zero_based)
80  for (octave_idx_type i = 0; i < a.numel (); i++)
81  ptmp[i] = static_cast<double>
82  (pa[i] + static_cast<octave_idx_type> (1));
83  else
84  for (octave_idx_type i = 0; i < a.numel (); i++)
85  ptmp[i] = static_cast<double> (pa[i]);
86  }
87 }
88 
90  : MArray<double> (a.dims ())
91 {
92  octave_idx_type n = a.numel ();
93  for (octave_idx_type i = 0; i < n; i++)
94  xelem (i) = static_cast<unsigned char> (a(i));
95 }
96 
97 #if defined (HAVE_FFTW)
98 
100 NDArray::fourier (int dim) const
101 {
102  dim_vector dv = dims ();
103 
104  if (dim > dv.length () || dim < 0)
105  return ComplexNDArray ();
106 
107  octave_idx_type stride = 1;
108  octave_idx_type n = dv(dim);
109 
110  for (int i = 0; i < dim; i++)
111  stride *= dv(i);
112 
113  octave_idx_type howmany = numel () / dv (dim);
114  howmany = (stride == 1 ? howmany : (howmany > stride ? stride : howmany));
115  octave_idx_type nloop = (stride == 1 ? 1 : numel () / dv (dim) / stride);
116  octave_idx_type dist = (stride == 1 ? n : 1);
117 
118  const double *in (fortran_vec ());
119  ComplexNDArray retval (dv);
120  Complex *out (retval.fortran_vec ());
121 
122  // Need to be careful here about the distance between fft's
123  for (octave_idx_type k = 0; k < nloop; k++)
124  octave_fftw::fft (in + k * stride * n, out + k * stride * n,
125  n, howmany, stride, dist);
126 
127  return retval;
128 }
129 
131 NDArray::ifourier (int dim) const
132 {
133  dim_vector dv = dims ();
134 
135  if (dim > dv.length () || dim < 0)
136  return ComplexNDArray ();
137 
138  octave_idx_type stride = 1;
139  octave_idx_type n = dv(dim);
140 
141  for (int i = 0; i < dim; i++)
142  stride *= dv(i);
143 
144  octave_idx_type howmany = numel () / dv (dim);
145  howmany = (stride == 1 ? howmany : (howmany > stride ? stride : howmany));
146  octave_idx_type nloop = (stride == 1 ? 1 : numel () / dv (dim) / stride);
147  octave_idx_type dist = (stride == 1 ? n : 1);
148 
149  ComplexNDArray retval (*this);
150  Complex *out (retval.fortran_vec ());
151 
152  // Need to be careful here about the distance between fft's
153  for (octave_idx_type k = 0; k < nloop; k++)
154  octave_fftw::ifft (out + k * stride * n, out + k * stride * n,
155  n, howmany, stride, dist);
156 
157  return retval;
158 }
159 
161 NDArray::fourier2d (void) const
162 {
163  dim_vector dv = dims ();
164  if (dv.length () < 2)
165  return ComplexNDArray ();
166 
167  dim_vector dv2(dv(0), dv(1));
168  const double *in = fortran_vec ();
169  ComplexNDArray retval (dv);
170  Complex *out = retval.fortran_vec ();
171  octave_idx_type howmany = numel () / dv(0) / dv(1);
172  octave_idx_type dist = dv(0) * dv(1);
173 
174  for (octave_idx_type i=0; i < howmany; i++)
175  octave_fftw::fftNd (in + i*dist, out + i*dist, 2, dv2);
176 
177  return retval;
178 }
179 
182 {
183  dim_vector dv = dims ();
184  if (dv.length () < 2)
185  return ComplexNDArray ();
186 
187  dim_vector dv2(dv(0), dv(1));
188  ComplexNDArray retval (*this);
189  Complex *out = retval.fortran_vec ();
190  octave_idx_type howmany = numel () / dv(0) / dv(1);
191  octave_idx_type dist = dv(0) * dv(1);
192 
193  for (octave_idx_type i=0; i < howmany; i++)
194  octave_fftw::ifftNd (out + i*dist, out + i*dist, 2, dv2);
195 
196  return retval;
197 }
198 
200 NDArray::fourierNd (void) const
201 {
202  dim_vector dv = dims ();
203  int rank = dv.length ();
204 
205  const double *in (fortran_vec ());
206  ComplexNDArray retval (dv);
207  Complex *out (retval.fortran_vec ());
208 
209  octave_fftw::fftNd (in, out, rank, dv);
210 
211  return retval;
212 }
213 
216 {
217  dim_vector dv = dims ();
218  int rank = dv.length ();
219 
220  ComplexNDArray tmp (*this);
221  Complex *in (tmp.fortran_vec ());
222  ComplexNDArray retval (dv);
223  Complex *out (retval.fortran_vec ());
224 
225  octave_fftw::ifftNd (in, out, rank, dv);
226 
227  return retval;
228 }
229 
230 #else
231 
232 extern "C"
233 {
234  // Note that the original complex fft routines were not written for
235  // double complex arguments. They have been modified by adding an
236  // implicit double precision (a-h,o-z) statement at the beginning of
237  // each subroutine.
238 
239  F77_RET_T
240  F77_FUNC (zffti, ZFFTI) (const octave_idx_type&, Complex*);
241 
242  F77_RET_T
243  F77_FUNC (zfftf, ZFFTF) (const octave_idx_type&, Complex*, Complex*);
244 
245  F77_RET_T
246  F77_FUNC (zfftb, ZFFTB) (const octave_idx_type&, Complex*, Complex*);
247 }
248 
250 NDArray::fourier (int dim) const
251 {
252  dim_vector dv = dims ();
253 
254  if (dim > dv.length () || dim < 0)
255  return ComplexNDArray ();
256 
257  ComplexNDArray retval (dv);
258  octave_idx_type npts = dv(dim);
259  octave_idx_type nn = 4*npts+15;
260  Array<Complex> wsave (dim_vector (nn, 1));
261  Complex *pwsave = wsave.fortran_vec ();
262 
263  OCTAVE_LOCAL_BUFFER (Complex, tmp, npts);
264 
265  octave_idx_type stride = 1;
266 
267  for (int i = 0; i < dim; i++)
268  stride *= dv(i);
269 
270  octave_idx_type howmany = numel () / npts;
271  howmany = (stride == 1 ? howmany : (howmany > stride ? stride : howmany));
272  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
273  octave_idx_type dist = (stride == 1 ? npts : 1);
274 
275  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
276 
277  for (octave_idx_type k = 0; k < nloop; k++)
278  {
279  for (octave_idx_type j = 0; j < howmany; j++)
280  {
281  octave_quit ();
282 
283  for (octave_idx_type i = 0; i < npts; i++)
284  tmp[i] = elem ((i + k*npts)*stride + j*dist);
285 
286  F77_FUNC (zfftf, ZFFTF) (npts, tmp, pwsave);
287 
288  for (octave_idx_type i = 0; i < npts; i++)
289  retval((i + k*npts)*stride + j*dist) = tmp[i];
290  }
291  }
292 
293  return retval;
294 }
295 
297 NDArray::ifourier (int dim) const
298 {
299  dim_vector dv = dims ();
300 
301  if (dim > dv.length () || dim < 0)
302  return ComplexNDArray ();
303 
304  ComplexNDArray retval (dv);
305  octave_idx_type npts = dv(dim);
306  octave_idx_type nn = 4*npts+15;
307  Array<Complex> wsave (dim_vector (nn, 1));
308  Complex *pwsave = wsave.fortran_vec ();
309 
310  OCTAVE_LOCAL_BUFFER (Complex, tmp, npts);
311 
312  octave_idx_type stride = 1;
313 
314  for (int i = 0; i < dim; i++)
315  stride *= dv(i);
316 
317  octave_idx_type howmany = numel () / npts;
318  howmany = (stride == 1 ? howmany : (howmany > stride ? stride : howmany));
319  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
320  octave_idx_type dist = (stride == 1 ? npts : 1);
321 
322  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
323 
324  for (octave_idx_type k = 0; k < nloop; k++)
325  {
326  for (octave_idx_type j = 0; j < howmany; j++)
327  {
328  octave_quit ();
329 
330  for (octave_idx_type i = 0; i < npts; i++)
331  tmp[i] = elem ((i + k*npts)*stride + j*dist);
332 
333  F77_FUNC (zfftb, ZFFTB) (npts, tmp, pwsave);
334 
335  for (octave_idx_type i = 0; i < npts; i++)
336  retval((i + k*npts)*stride + j*dist) = tmp[i] /
337  static_cast<double> (npts);
338  }
339  }
340 
341  return retval;
342 }
343 
345 NDArray::fourier2d (void) const
346 {
347  dim_vector dv = dims ();
348  dim_vector dv2 (dv(0), dv(1));
349  int rank = 2;
350  ComplexNDArray retval (*this);
351  octave_idx_type stride = 1;
352 
353  for (int i = 0; i < rank; i++)
354  {
355  octave_idx_type npts = dv2(i);
356  octave_idx_type nn = 4*npts+15;
357  Array<Complex> wsave (dim_vector (nn, 1));
358  Complex *pwsave = wsave.fortran_vec ();
359  Array<Complex> row (dim_vector (npts, 1));
360  Complex *prow = row.fortran_vec ();
361 
362  octave_idx_type howmany = numel () / npts;
363  howmany = (stride == 1 ? howmany :
364  (howmany > stride ? stride : howmany));
365  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
366  octave_idx_type dist = (stride == 1 ? npts : 1);
367 
368  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
369 
370  for (octave_idx_type k = 0; k < nloop; k++)
371  {
372  for (octave_idx_type j = 0; j < howmany; j++)
373  {
374  octave_quit ();
375 
376  for (octave_idx_type l = 0; l < npts; l++)
377  prow[l] = retval((l + k*npts)*stride + j*dist);
378 
379  F77_FUNC (zfftf, ZFFTF) (npts, prow, pwsave);
380 
381  for (octave_idx_type l = 0; l < npts; l++)
382  retval((l + k*npts)*stride + j*dist) = prow[l];
383  }
384  }
385 
386  stride *= dv2(i);
387  }
388 
389  return retval;
390 }
391 
393 NDArray::ifourier2d (void) const
394 {
395  dim_vector dv = dims ();
396  dim_vector dv2 (dv(0), dv(1));
397  int rank = 2;
398  ComplexNDArray retval (*this);
399  octave_idx_type stride = 1;
400 
401  for (int i = 0; i < rank; i++)
402  {
403  octave_idx_type npts = dv2(i);
404  octave_idx_type nn = 4*npts+15;
405  Array<Complex> wsave (dim_vector (nn, 1));
406  Complex *pwsave = wsave.fortran_vec ();
407  Array<Complex> row (dim_vector (npts, 1));
408  Complex *prow = row.fortran_vec ();
409 
410  octave_idx_type howmany = numel () / npts;
411  howmany = (stride == 1 ? howmany :
412  (howmany > stride ? stride : howmany));
413  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
414  octave_idx_type dist = (stride == 1 ? npts : 1);
415 
416  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
417 
418  for (octave_idx_type k = 0; k < nloop; k++)
419  {
420  for (octave_idx_type j = 0; j < howmany; j++)
421  {
422  octave_quit ();
423 
424  for (octave_idx_type l = 0; l < npts; l++)
425  prow[l] = retval((l + k*npts)*stride + j*dist);
426 
427  F77_FUNC (zfftb, ZFFTB) (npts, prow, pwsave);
428 
429  for (octave_idx_type l = 0; l < npts; l++)
430  retval((l + k*npts)*stride + j*dist) =
431  prow[l] / static_cast<double> (npts);
432  }
433  }
434 
435  stride *= dv2(i);
436  }
437 
438  return retval;
439 }
440 
442 NDArray::fourierNd (void) const
443 {
444  dim_vector dv = dims ();
445  int rank = dv.length ();
446  ComplexNDArray retval (*this);
447  octave_idx_type stride = 1;
448 
449  for (int i = 0; i < rank; i++)
450  {
451  octave_idx_type npts = dv(i);
452  octave_idx_type nn = 4*npts+15;
453  Array<Complex> wsave (dim_vector (nn, 1));
454  Complex *pwsave = wsave.fortran_vec ();
455  Array<Complex> row (dim_vector (npts, 1));
456  Complex *prow = row.fortran_vec ();
457 
458  octave_idx_type howmany = numel () / npts;
459  howmany = (stride == 1 ? howmany :
460  (howmany > stride ? stride : howmany));
461  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
462  octave_idx_type dist = (stride == 1 ? npts : 1);
463 
464  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
465 
466  for (octave_idx_type k = 0; k < nloop; k++)
467  {
468  for (octave_idx_type j = 0; j < howmany; j++)
469  {
470  octave_quit ();
471 
472  for (octave_idx_type l = 0; l < npts; l++)
473  prow[l] = retval((l + k*npts)*stride + j*dist);
474 
475  F77_FUNC (zfftf, ZFFTF) (npts, prow, pwsave);
476 
477  for (octave_idx_type l = 0; l < npts; l++)
478  retval((l + k*npts)*stride + j*dist) = prow[l];
479  }
480  }
481 
482  stride *= dv(i);
483  }
484 
485  return retval;
486 }
487 
489 NDArray::ifourierNd (void) const
490 {
491  dim_vector dv = dims ();
492  int rank = dv.length ();
493  ComplexNDArray retval (*this);
494  octave_idx_type stride = 1;
495 
496  for (int i = 0; i < rank; i++)
497  {
498  octave_idx_type npts = dv(i);
499  octave_idx_type nn = 4*npts+15;
500  Array<Complex> wsave (dim_vector (nn, 1));
501  Complex *pwsave = wsave.fortran_vec ();
502  Array<Complex> row (dim_vector (npts, 1));
503  Complex *prow = row.fortran_vec ();
504 
505  octave_idx_type howmany = numel () / npts;
506  howmany = (stride == 1 ? howmany :
507  (howmany > stride ? stride : howmany));
508  octave_idx_type nloop = (stride == 1 ? 1 : numel () / npts / stride);
509  octave_idx_type dist = (stride == 1 ? npts : 1);
510 
511  F77_FUNC (zffti, ZFFTI) (npts, pwsave);
512 
513  for (octave_idx_type k = 0; k < nloop; k++)
514  {
515  for (octave_idx_type j = 0; j < howmany; j++)
516  {
517  octave_quit ();
518 
519  for (octave_idx_type l = 0; l < npts; l++)
520  prow[l] = retval((l + k*npts)*stride + j*dist);
521 
522  F77_FUNC (zfftb, ZFFTB) (npts, prow, pwsave);
523 
524  for (octave_idx_type l = 0; l < npts; l++)
525  retval((l + k*npts)*stride + j*dist) =
526  prow[l] / static_cast<double> (npts);
527  }
528  }
529 
530  stride *= dv(i);
531  }
532 
533  return retval;
534 }
535 
536 #endif
537 
538 // unary operations
539 
542 {
543  if (any_element_is_nan ())
545 
546  return do_mx_unary_op<bool, double> (*this, mx_inline_not);
547 }
548 
549 bool
550 NDArray::any_element_is_negative (bool neg_zero) const
551 {
552  return (neg_zero ? test_all (xnegative_sign)
553  : do_mx_check<double> (*this, mx_inline_any_negative));
554 }
555 
556 bool
557 NDArray::any_element_is_positive (bool neg_zero) const
558 {
559  return (neg_zero ? test_all (xpositive_sign)
560  : do_mx_check<double> (*this, mx_inline_any_positive));
561 }
562 
563 bool
565 {
566  return do_mx_check<double> (*this, mx_inline_any_nan);
567 }
568 
569 bool
571 {
572  return ! do_mx_check<double> (*this, mx_inline_all_finite);
573 }
574 
575 bool
577 {
578  return ! test_all (xis_one_or_zero);
579 }
580 
581 bool
583 {
584  return test_all (xis_zero);
585 }
586 
587 bool
589 {
591 }
592 
593 // Return nonzero if any element of M is not an integer. Also extract
594 // the largest and smallest values and return them in MAX_VAL and MIN_VAL.
595 
596 bool
597 NDArray::all_integers (double& max_val, double& min_val) const
598 {
599  octave_idx_type nel = nelem ();
600 
601  if (nel > 0)
602  {
603  max_val = elem (0);
604  min_val = elem (0);
605  }
606  else
607  return false;
608 
609  for (octave_idx_type i = 0; i < nel; i++)
610  {
611  double val = elem (i);
612 
613  if (val > max_val)
614  max_val = val;
615 
616  if (val < min_val)
617  min_val = val;
618 
619  if (! xisinteger (val))
620  return false;
621  }
622 
623  return true;
624 }
625 
626 bool
628 {
629  return test_all (xisinteger);
630 }
631 
632 bool
634 {
636 }
637 
638 // FIXME: this is not quite the right thing.
639 
641 NDArray::all (int dim) const
642 {
643  return do_mx_red_op<bool, double> (*this, dim, mx_inline_all);
644 }
645 
647 NDArray::any (int dim) const
648 {
649  return do_mx_red_op<bool, double> (*this, dim, mx_inline_any);
650 }
651 
652 NDArray
653 NDArray::cumprod (int dim) const
654 {
655  return do_mx_cum_op<double, double> (*this, dim, mx_inline_cumprod);
656 }
657 
658 NDArray
659 NDArray::cumsum (int dim) const
660 {
661  return do_mx_cum_op<double, double> (*this, dim, mx_inline_cumsum);
662 }
663 
664 NDArray
665 NDArray::prod (int dim) const
666 {
667  return do_mx_red_op<double, double> (*this, dim, mx_inline_prod);
668 }
669 
670 NDArray
671 NDArray::sum (int dim) const
672 {
673  return do_mx_red_op<double, double> (*this, dim, mx_inline_sum);
674 }
675 
676 NDArray
677 NDArray::xsum (int dim) const
678 {
679  return do_mx_red_op<double, double> (*this, dim, mx_inline_xsum);
680 }
681 
682 NDArray
683 NDArray::sumsq (int dim) const
684 {
685  return do_mx_red_op<double, double> (*this, dim, mx_inline_sumsq);
686 }
687 
688 NDArray
689 NDArray::max (int dim) const
690 {
691  return do_mx_minmax_op<double> (*this, dim, mx_inline_max);
692 }
693 
694 NDArray
695 NDArray::max (Array<octave_idx_type>& idx_arg, int dim) const
696 {
697  return do_mx_minmax_op<double> (*this, idx_arg, dim, mx_inline_max);
698 }
699 
700 NDArray
701 NDArray::min (int dim) const
702 {
703  return do_mx_minmax_op<double> (*this, dim, mx_inline_min);
704 }
705 
706 NDArray
707 NDArray::min (Array<octave_idx_type>& idx_arg, int dim) const
708 {
709  return do_mx_minmax_op<double> (*this, idx_arg, dim, mx_inline_min);
710 }
711 
712 NDArray
713 NDArray::cummax (int dim) const
714 {
715  return do_mx_cumminmax_op<double> (*this, dim, mx_inline_cummax);
716 }
717 
718 NDArray
719 NDArray::cummax (Array<octave_idx_type>& idx_arg, int dim) const
720 {
721  return do_mx_cumminmax_op<double> (*this, idx_arg, dim, mx_inline_cummax);
722 }
723 
724 NDArray
725 NDArray::cummin (int dim) const
726 {
727  return do_mx_cumminmax_op<double> (*this, dim, mx_inline_cummin);
728 }
729 
730 NDArray
731 NDArray::cummin (Array<octave_idx_type>& idx_arg, int dim) const
732 {
733  return do_mx_cumminmax_op<double> (*this, idx_arg, dim, mx_inline_cummin);
734 }
735 
736 NDArray
737 NDArray::diff (octave_idx_type order, int dim) const
738 {
739  return do_mx_diff_op<double> (*this, dim, order, mx_inline_diff);
740 }
741 
742 NDArray
744 {
745  if (rb.numel () > 0)
746  insert (rb, ra_idx);
747  return *this;
748 }
749 
752 {
753  ComplexNDArray retval (*this);
754  if (rb.numel () > 0)
755  retval.insert (rb, ra_idx);
756  return retval;
757 }
758 
761 {
762  charNDArray retval (dims ());
763  octave_idx_type nel = numel ();
764 
765  for (octave_idx_type i = 0; i < nel; i++)
766  {
767  double d = elem (i);
768 
769  if (xisnan (d))
770  {
771  (*current_liboctave_error_handler)
772  ("invalid conversion from NaN to character");
773  return retval;
774  }
775  else
776  {
777  octave_idx_type ival = NINTbig (d);
778 
779  if (ival < 0 || ival > std::numeric_limits<unsigned char>::max ())
780  // FIXME: is there something better to do? Should we warn the user?
781  ival = 0;
782 
783  retval.elem (i) = static_cast<char>(ival);
784  }
785  }
786 
787  if (rb.numel () == 0)
788  return retval;
789 
790  retval.insert (rb, ra_idx);
791  return retval;
792 }
793 
794 NDArray
796 {
797  return do_mx_unary_op<double, Complex> (a, mx_inline_real);
798 }
799 
800 NDArray
802 {
803  return do_mx_unary_op<double, Complex> (a, mx_inline_imag);
804 }
805 
806 NDArray&
808 {
809  Array<double>::insert (a, r, c);
810  return *this;
811 }
812 
813 NDArray&
815 {
816  Array<double>::insert (a, ra_idx);
817  return *this;
818 }
819 
820 NDArray
821 NDArray::abs (void) const
822 {
823  return do_mx_unary_map<double, double, std::abs> (*this);
824 }
825 
827 NDArray::isnan (void) const
828 {
829  return do_mx_unary_map<bool, double, xisnan> (*this);
830 }
831 
833 NDArray::isinf (void) const
834 {
835  return do_mx_unary_map<bool, double, xisinf> (*this);
836 }
837 
839 NDArray::isfinite (void) const
840 {
841  return do_mx_unary_map<bool, double, xfinite> (*this);
842 }
843 
844 void
846  const dim_vector& dimensions,
847  int start_dimension)
848 {
849  ::increment_index (ra_idx, dimensions, start_dimension);
850 }
851 
854  const dim_vector& dimensions)
855 {
856  return ::compute_index (ra_idx, dimensions);
857 }
858 
859 NDArray
861 {
862  return MArray<double>::diag (k);
863 }
864 
865 NDArray
867 {
868  return MArray<double>::diag (m, n);
869 }
870 
871 // This contains no information on the array structure !!!
872 std::ostream&
873 operator << (std::ostream& os, const NDArray& a)
874 {
875  octave_idx_type nel = a.nelem ();
876 
877  for (octave_idx_type i = 0; i < nel; i++)
878  {
879  os << " ";
880  octave_write_double (os, a.elem (i));
881  os << "\n";
882  }
883  return os;
884 }
885 
886 std::istream&
887 operator >> (std::istream& is, NDArray& a)
888 {
889  octave_idx_type nel = a.nelem ();
890 
891  if (nel > 0)
892  {
893  double tmp;
894  for (octave_idx_type i = 0; i < nel; i++)
895  {
896  tmp = octave_read_value<double> (is);
897  if (is)
898  a.elem (i) = tmp;
899  else
900  goto done;
901  }
902  }
903 
904 done:
905 
906  return is;
907 }
908 
910 
912 NDS_BOOL_OPS (NDArray, double)
913 
914 SND_CMP_OPS (double, NDArray)
915 SND_BOOL_OPS (double, NDArray)
916 
917 NDND_CMP_OPS (NDArray, NDArray)
918 NDND_BOOL_OPS (NDArray, NDArray)
919 
922 
924 BSXFUN_OP2_DEF_MXLOOP (pow, ComplexNDArray, ComplexNDArray,
925  NDArray, mx_inline_pow)
926 BSXFUN_OP2_DEF_MXLOOP (pow, ComplexNDArray, NDArray,
927  ComplexNDArray, mx_inline_pow)
dim_vector dimensions
Definition: Array.h:127
T mx_inline_xsum(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:1344
void octave_write_double(std::ostream &os, double d)
Definition: lo-utils.cc:386
octave_idx_type compute_index(octave_idx_type n, const dim_vector &dims)
Definition: Array-util.cc:178
bool any_element_is_negative(bool=false) const
Definition: dNDArray.cc:550
NDArray cumsum(int dim=-1) const
Definition: dNDArray.cc:659
#define NDS_BOOL_OPS(ND, S)
Definition: mx-op-defs.h:253
#define NDS_CMP_OPS(ND, S)
Definition: mx-op-defs.h:236
void mx_inline_cummax(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:894
std::istream & operator>>(std::istream &is, NDArray &a)
Definition: dNDArray.cc:887
static int fft(const double *in, Complex *out, size_t npts, size_t nsamples=1, octave_idx_type stride=1, octave_idx_type dist=-1)
Definition: oct-fftw.cc:827
const octave_base_value const Array< octave_idx_type > & ra_idx
#define SND_CMP_OPS(S, ND)
Definition: mx-op-defs.h:283
bool xisnan(double x)
Definition: lo-mappers.cc:144
#define BSXFUN_OP_DEF_MXLOOP(OP, ARRAY, LOOP)
Definition: bsxfun-defs.cc:221
octave_idx_type numel(void) const
Number of elements in the array.
Definition: Array.h:275
subroutine zffti(n, wsave)
Definition: zffti.f:1
static octave_idx_type nn
Definition: DASPK.cc:71
NDArray max(double d, const NDArray &m)
Definition: dNDArray.cc:909
NDArray sumsq(int dim=-1) const
Definition: dNDArray.cc:683
void mx_inline_real(size_t n, T *r, const std::complex< T > *x)
Definition: mx-inlines.cc:251
bool mx_inline_all_finite(size_t n, const T *x)
Definition: mx-inlines.cc:195
NDArray prod(int dim=-1) const
Definition: dNDArray.cc:665
#define MINMAX_FCNS(T, S)
Definition: mx-op-defs.h:590
bool all_elements_are_zero(void) const
Definition: dNDArray.cc:582
bool xnegative_sign(double x)
Definition: lo-mappers.cc:618
bool xis_zero(double x)
Definition: lo-utils.cc:53
Array< T > diag(octave_idx_type k=0) const
Get the kth super or subdiagonal.
Definition: Array.cc:2526
NDArray(void)
Definition: dNDArray.h:39
ComplexRowVector row(octave_idx_type i) const
Definition: CMatrix.cc:990
double & elem(octave_idx_type n)
Definition: Array.h:380
double lo_ieee_nan_value(void)
Definition: lo-ieee.cc:126
Definition: MArray.h:36
T mx_inline_sumsq(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:524
std::ostream & operator<<(std::ostream &os, const NDArray &a)
Definition: dNDArray.cc:873
bool mx_inline_any_negative(size_t n, const T *x)
Definition: mx-inlines.cc:208
boolNDArray isnan(void) const
Definition: dNDArray.cc:827
NDArray cumprod(int dim=-1) const
Definition: dNDArray.cc:653
bool any_element_is_nan(void) const
Definition: dNDArray.cc:564
void mx_inline_max(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:719
bool mx_inline_any(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:524
F77_RET_T const double const double double * d
NDArray imag(const ComplexNDArray &a)
Definition: dNDArray.cc:801
subroutine zfftb(n, c, wsave)
Definition: zfftb.f:1
static int fftNd(const double *, Complex *, const int, const dim_vector &)
Definition: oct-fftw.cc:889
NDArray concat(const NDArray &rb, const Array< octave_idx_type > &ra_idx)
Definition: dNDArray.cc:743
void mx_inline_imag(size_t n, T *r, const std::complex< T > *x)
Definition: mx-inlines.cc:254
Array< T > & insert(const Array< T > &a, const Array< octave_idx_type > &idx)
Insert an array into another at a specified position.
Definition: Array.cc:1591
octave_idx_type nelem(void) const
Number of elements in the array.
Definition: Array.h:271
bool any_element_not_one_or_zero(void) const
Definition: dNDArray.cc:576
const dim_vector & dims(void) const
Return a const-reference so that dims ()(i) works efficiently.
Definition: Array.h:337
NDArray min(int dim=-1) const
Definition: dNDArray.cc:701
ComplexNDArray fourierNd(void) const
Definition: dNDArray.cc:200
#define NDND_CMP_OPS(ND1, ND2)
Definition: mx-op-defs.h:330
NDArray real(const ComplexNDArray &a)
Definition: dNDArray.cc:795
bool test_any(F fcn) const
Simpler calls.
Definition: Array.h:710
ComplexNDArray & insert(const NDArray &a, octave_idx_type r, octave_idx_type c)
Definition: CNDArray.cc:753
#define NDND_BOOL_OPS(ND1, ND2)
Definition: mx-op-defs.h:347
ComplexNDArray ifourier(int dim=1) const
Definition: dNDArray.cc:131
NDArray & insert(const NDArray &a, octave_idx_type r, octave_idx_type c)
Definition: dNDArray.cc:807
void gripe_nan_to_logical_conversion(void)
bool any_element_is_positive(bool=false) const
Definition: dNDArray.cc:557
ComplexNDArray ifourier2d(void) const
Definition: dNDArray.cc:181
static int ifft(const Complex *in, Complex *out, size_t npts, size_t nsamples=1, octave_idx_type stride=1, octave_idx_type dist=-1)
Definition: oct-fftw.cc:865
NDArray cummax(int dim=-1) const
Definition: dNDArray.cc:713
void resize(const dim_vector &dv, const double &rfv)
static octave_idx_type compute_index(Array< octave_idx_type > &ra_idx, const dim_vector &dimensions)
Definition: dNDArray.cc:853
octave_int< T > pow(const octave_int< T > &a, const octave_int< T > &b)
#define SND_BOOL_OPS(S, ND)
Definition: mx-op-defs.h:300
bool mx_inline_any_positive(size_t n, const T *x)
Definition: mx-inlines.cc:221
NDArray max(int dim=-1) const
Definition: dNDArray.cc:689
#define F77_RET_T
Definition: f77-fcn.h:264
void mx_inline_cummin(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:893
void mx_inline_pow(size_t n, R *r, const X *x, const Y *y)
Definition: mx-inlines.cc:311
NDArray abs(void) const
Definition: dNDArray.cc:821
bool mx_inline_any_nan(size_t n, const T *x)
Definition: mx-inlines.cc:182
static int ifftNd(const Complex *, Complex *, const int, const dim_vector &)
Definition: oct-fftw.cc:933
#define BSXFUN_STDREL_DEFS_MXLOOP(ARRAY)
Definition: bsxfun-defs.cc:244
double & xelem(octave_idx_type n)
Definition: Array.h:353
T mx_inline_sum(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:522
bool all_elements_are_int_or_inf_or_nan(void) const
Definition: dNDArray.cc:588
charNDArray & insert(const charNDArray &a, octave_idx_type r, octave_idx_type c)
Definition: chNDArray.cc:169
bool test_all(F fcn) const
Definition: Array.h:714
boolNDArray all(int dim=-1) const
Definition: dNDArray.cc:641
NDArray cummin(int dim=-1) const
Definition: dNDArray.cc:725
static void increment_index(Array< octave_idx_type > &ra_idx, const dim_vector &dimensions, int start_dimension=0)
Definition: dNDArray.cc:845
void mx_inline_cumprod(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:628
NDArray sum(int dim=-1) const
Definition: dNDArray.cc:671
boolNDArray any(int dim=-1) const
Definition: dNDArray.cc:647
NDArray diag(octave_idx_type k=0) const
Definition: dNDArray.cc:860
bool xisinteger(double x)
Definition: lo-mappers.h:206
void mx_inline_cumsum(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:627
boolNDArray isfinite(void) const
Definition: dNDArray.cc:839
bool xpositive_sign(double x)
Definition: lo-mappers.h:222
ComplexNDArray fourier(int dim=1) const
Definition: dNDArray.cc:100
bool all_integers(void) const
Definition: dNDArray.cc:627
F77_RET_T F77_FUNC(zgemv, ZGEMV)(F77_CONST_CHAR_ARG_DECL
ComplexNDArray fourier2d(void) const
Definition: dNDArray.cc:161
bool any_element_is_inf_or_nan(void) const
Definition: dNDArray.cc:570
#define OCTAVE_LOCAL_BUFFER(T, buf, size)
Definition: oct-locbuf.h:197
#define BSXFUN_OP2_DEF_MXLOOP(OP, ARRAY, ARRAY1, ARRAY2, LOOP)
Definition: bsxfun-defs.cc:226
bool xtoo_large_for_float(double x)
Definition: lo-utils.cc:56
ComplexNDArray ifourierNd(void) const
Definition: dNDArray.cc:215
bool mx_inline_all(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:524
#define BSXFUN_STDOP_DEFS_MXLOOP(ARRAY)
Definition: bsxfun-defs.cc:236
std::complex< double > Complex
Definition: oct-cmplx.h:29
const T * fortran_vec(void) const
Definition: Array.h:481
octave_idx_type NINTbig(double x)
Definition: lo-mappers.cc:635
void mx_inline_not(size_t n, bool *r, const X *x)
Definition: mx-inlines.cc:126
int length(void) const
Definition: dim-vector.h:281
T mx_inline_prod(const T *v, octave_idx_type n)
Definition: mx-inlines.cc:523
bool xis_one_or_zero(double x)
Definition: lo-utils.cc:50
NDArray diff(octave_idx_type order=1, int dim=-1) const
Definition: dNDArray.cc:737
subroutine zfftf(n, c, wsave)
Definition: zfftf.f:1
NDArray xsum(int dim=-1) const
Definition: dNDArray.cc:677
bool too_large_for_float(void) const
Definition: dNDArray.cc:633
boolNDArray isinf(void) const
Definition: dNDArray.cc:833
void mx_inline_diff(const T *v, T *r, octave_idx_type n, octave_idx_type order)
Definition: mx-inlines.cc:1036
bool xis_int_or_inf_or_nan(double x)
Definition: lo-utils.cc:47
boolNDArray operator!(void) const
Definition: dNDArray.cc:541
void mx_inline_min(const T *v, T *r, octave_idx_type n)
Definition: mx-inlines.cc:718