GNU Octave  3.8.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
fCRowVector.cc
Go to the documentation of this file.
1 // RowVector manipulations.
2 /*
3 
4 Copyright (C) 1994-2013 John W. Eaton
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 <iostream>
29 
30 #include "Array-util.h"
31 #include "f77-fcn.h"
32 #include "functor.h"
33 #include "lo-error.h"
34 #include "mx-base.h"
35 #include "mx-inlines.cc"
36 #include "oct-cmplx.h"
37 
38 // Fortran functions we call.
39 
40 extern "C"
41 {
42  F77_RET_T
43  F77_FUNC (cgemv, CGEMV) (F77_CONST_CHAR_ARG_DECL,
45  const FloatComplex&, const FloatComplex*,
46  const octave_idx_type&, const FloatComplex*,
47  const octave_idx_type&, const FloatComplex&,
48  FloatComplex*, const octave_idx_type&
50 
51  F77_RET_T
52  F77_FUNC (xcdotu, XCDOTU) (const octave_idx_type&, const FloatComplex*,
53  const octave_idx_type&, const FloatComplex*,
54  const octave_idx_type&, FloatComplex&);
55 }
56 
57 // FloatComplex Row Vector class
58 
59 bool
61 {
62  octave_idx_type len = length ();
63  if (len != a.length ())
64  return 0;
65  return mx_inline_equal (len, data (), a.data ());
66 }
67 
68 bool
70 {
71  return !(*this == a);
72 }
73 
74 // destructive insert/delete/reorder operations
75 
78 {
79  octave_idx_type a_len = a.length ();
80 
81  if (c < 0 || c + a_len > length ())
82  {
83  (*current_liboctave_error_handler) ("range error for insert");
84  return *this;
85  }
86 
87  if (a_len > 0)
88  {
89  make_unique ();
90 
91  for (octave_idx_type i = 0; i < a_len; i++)
92  xelem (c+i) = a.elem (i);
93  }
94 
95  return *this;
96 }
97 
100  octave_idx_type c)
101 {
102  octave_idx_type a_len = a.length ();
103 
104  if (c < 0 || c + a_len > length ())
105  {
106  (*current_liboctave_error_handler) ("range error for insert");
107  return *this;
108  }
109 
110  if (a_len > 0)
111  {
112  make_unique ();
113 
114  for (octave_idx_type i = 0; i < a_len; i++)
115  xelem (c+i) = a.elem (i);
116  }
117 
118  return *this;
119 }
120 
123 {
124  octave_idx_type len = length ();
125 
126  if (len > 0)
127  {
128  make_unique ();
129 
130  for (octave_idx_type i = 0; i < len; i++)
131  xelem (i) = val;
132  }
133 
134  return *this;
135 }
136 
139 {
140  octave_idx_type len = length ();
141 
142  if (len > 0)
143  {
144  make_unique ();
145 
146  for (octave_idx_type i = 0; i < len; i++)
147  xelem (i) = val;
148  }
149 
150  return *this;
151 }
152 
155 {
156  octave_idx_type len = length ();
157 
158  if (c1 < 0 || c2 < 0 || c1 >= len || c2 >= len)
159  {
160  (*current_liboctave_error_handler) ("range error for fill");
161  return *this;
162  }
163 
164  if (c1 > c2) { std::swap (c1, c2); }
165 
166  if (c2 >= c1)
167  {
168  make_unique ();
169 
170  for (octave_idx_type i = c1; i <= c2; i++)
171  xelem (i) = val;
172  }
173 
174  return *this;
175 }
176 
180 {
181  octave_idx_type len = length ();
182 
183  if (c1 < 0 || c2 < 0 || c1 >= len || c2 >= len)
184  {
185  (*current_liboctave_error_handler) ("range error for fill");
186  return *this;
187  }
188 
189  if (c1 > c2) { std::swap (c1, c2); }
190 
191  if (c2 >= c1)
192  {
193  make_unique ();
194 
195  for (octave_idx_type i = c1; i <= c2; i++)
196  xelem (i) = val;
197  }
198 
199  return *this;
200 }
201 
204 {
205  octave_idx_type len = length ();
206  octave_idx_type nc_insert = len;
207  FloatComplexRowVector retval (len + a.length ());
208  retval.insert (*this, 0);
209  retval.insert (a, nc_insert);
210  return retval;
211 }
212 
215 {
216  octave_idx_type len = length ();
217  octave_idx_type nc_insert = len;
218  FloatComplexRowVector retval (len + a.length ());
219  retval.insert (*this, 0);
220  retval.insert (a, nc_insert);
221  return retval;
222 }
223 
226 {
228 }
229 
232 {
234 }
235 
238 {
239  return do_mx_unary_map<FloatComplex, FloatComplex, std::conj<float> > (a);
240 }
241 
242 // resize is the destructive equivalent for this one
243 
246 {
247  if (c1 > c2) { std::swap (c1, c2); }
248 
249  octave_idx_type new_c = c2 - c1 + 1;
250 
251  FloatComplexRowVector result (new_c);
252 
253  for (octave_idx_type i = 0; i < new_c; i++)
254  result.elem (i) = elem (c1+i);
255 
256  return result;
257 }
258 
261 {
262  FloatComplexRowVector result (n);
263 
264  for (octave_idx_type i = 0; i < n; i++)
265  result.elem (i) = elem (r1+i);
266 
267  return result;
268 }
269 
270 // row vector by row vector -> row vector operations
271 
274 {
275  octave_idx_type len = length ();
276 
277  octave_idx_type a_len = a.length ();
278 
279  if (len != a_len)
280  {
281  gripe_nonconformant ("operator +=", len, a_len);
282  return *this;
283  }
284 
285  if (len == 0)
286  return *this;
287 
288  FloatComplex *d = fortran_vec (); // Ensures only 1 reference to my privates!
289 
290  mx_inline_add2 (len, d, a.data ());
291  return *this;
292 }
293 
296 {
297  octave_idx_type len = length ();
298 
299  octave_idx_type a_len = a.length ();
300 
301  if (len != a_len)
302  {
303  gripe_nonconformant ("operator -=", len, a_len);
304  return *this;
305  }
306 
307  if (len == 0)
308  return *this;
309 
310  FloatComplex *d = fortran_vec (); // Ensures only 1 reference to my privates!
311 
312  mx_inline_sub2 (len, d, a.data ());
313  return *this;
314 }
315 
316 // row vector by matrix -> row vector
317 
320 {
321  FloatComplexRowVector retval;
322 
323  octave_idx_type len = v.length ();
324 
325  octave_idx_type a_nr = a.rows ();
326  octave_idx_type a_nc = a.cols ();
327 
328  if (a_nr != len)
329  gripe_nonconformant ("operator *", 1, len, a_nr, a_nc);
330  else
331  {
332  if (len == 0)
333  retval.resize (a_nc, 0.0);
334  else
335  {
336  // Transpose A to form A'*x == (x'*A)'
337 
338  octave_idx_type ld = a_nr;
339 
340  retval.resize (a_nc);
341  FloatComplex *y = retval.fortran_vec ();
342 
343  F77_XFCN (cgemv, CGEMV, (F77_CONST_CHAR_ARG2 ("T", 1),
344  a_nr, a_nc, 1.0, a.data (),
345  ld, v.data (), 1, 0.0, y, 1
346  F77_CHAR_ARG_LEN (1)));
347  }
348  }
349 
350  return retval;
351 }
352 
355 {
356  FloatComplexRowVector tmp (v);
357  return tmp * a;
358 }
359 
360 // other operations
361 
364 {
365  octave_idx_type len = length ();
366  if (len == 0)
367  return FloatComplex (0.0);
368 
369  FloatComplex res = elem (0);
370  float absres = std::abs (res);
371 
372  for (octave_idx_type i = 1; i < len; i++)
373  if (std::abs (elem (i)) < absres)
374  {
375  res = elem (i);
376  absres = std::abs (res);
377  }
378 
379  return res;
380 }
381 
384 {
385  octave_idx_type len = length ();
386  if (len == 0)
387  return FloatComplex (0.0);
388 
389  FloatComplex res = elem (0);
390  float absres = std::abs (res);
391 
392  for (octave_idx_type i = 1; i < len; i++)
393  if (std::abs (elem (i)) > absres)
394  {
395  res = elem (i);
396  absres = std::abs (res);
397  }
398 
399  return res;
400 }
401 
402 // i/o
403 
404 std::ostream&
405 operator << (std::ostream& os, const FloatComplexRowVector& a)
406 {
407 // int field_width = os.precision () + 7;
408  for (octave_idx_type i = 0; i < a.length (); i++)
409  os << " " /* setw (field_width) */ << a.elem (i);
410  return os;
411 }
412 
413 std::istream&
414 operator >> (std::istream& is, FloatComplexRowVector& a)
415 {
416  octave_idx_type len = a.length ();
417 
418  if (len > 0)
419  {
420  FloatComplex tmp;
421  for (octave_idx_type i = 0; i < len; i++)
422  {
423  is >> tmp;
424  if (is)
425  a.elem (i) = tmp;
426  else
427  break;
428  }
429  }
430  return is;
431 }
432 
433 // row vector by column vector -> scalar
434 
435 // row vector by column vector -> scalar
436 
439 {
440  FloatComplexColumnVector tmp (a);
441  return v * tmp;
442 }
443 
446 {
447  FloatComplex retval (0.0, 0.0);
448 
449  octave_idx_type len = v.length ();
450 
451  octave_idx_type a_len = a.length ();
452 
453  if (len != a_len)
454  gripe_nonconformant ("operator *", len, a_len);
455  else if (len != 0)
456  F77_FUNC (xcdotu, XCDOTU) (len, v.data (), 1, a.data (), 1, retval);
457 
458  return retval;
459 }
460 
461 // other operations
462 
465 {
466  if (n < 1) n = 1;
467 
469 
470  FloatComplex delta = (x2 - x1) / (n - 1.0f);
471  retval(0) = x1;
472  for (octave_idx_type i = 1; i < n-1; i++)
473  retval(i) = x1 + static_cast<float> (i)*delta;
474  retval(n-1) = x2;
475 
476  return retval;
477 }