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
CmplxQR.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 1994-2013 John W. Eaton
4 Copyright (C) 2008-2009 Jaroslav Hajek
5 Copyright (C) 2009 VZLU Prague
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 "CmplxQR.h"
30 #include "f77-fcn.h"
31 #include "lo-error.h"
32 #include "Range.h"
33 #include "idx-vector.h"
34 #include "oct-locbuf.h"
35 
36 #include "base-qr.cc"
37 
38 template class base_qr<ComplexMatrix>;
39 
40 extern "C"
41 {
42  F77_RET_T
43  F77_FUNC (zgeqrf, ZGEQRF) (const octave_idx_type&, const octave_idx_type&,
44  Complex*, const octave_idx_type&, Complex*,
45  Complex*, const octave_idx_type&,
46  octave_idx_type&);
47 
48  F77_RET_T
49  F77_FUNC (zungqr, ZUNGQR) (const octave_idx_type&, const octave_idx_type&,
50  const octave_idx_type&, Complex*,
51  const octave_idx_type&, Complex*, Complex*,
52  const octave_idx_type&, octave_idx_type&);
53 
54 #ifdef HAVE_QRUPDATE
55 
56  F77_RET_T
57  F77_FUNC (zqr1up, ZQR1UP) (const octave_idx_type&, const octave_idx_type&,
58  const octave_idx_type&, Complex*,
59  const octave_idx_type&, Complex*,
60  const octave_idx_type&, Complex*,
61  Complex*, Complex*, double*);
62 
63  F77_RET_T
64  F77_FUNC (zqrinc, ZQRINC) (const octave_idx_type&, const octave_idx_type&,
65  const octave_idx_type&, Complex*,
66  const octave_idx_type&, Complex*,
67  const octave_idx_type&, const octave_idx_type&,
68  const Complex*, double*);
69 
70  F77_RET_T
71  F77_FUNC (zqrdec, ZQRDEC) (const octave_idx_type&, const octave_idx_type&,
72  const octave_idx_type&, Complex*,
73  const octave_idx_type&, Complex*,
74  const octave_idx_type&, const octave_idx_type&,
75  double*);
76 
77  F77_RET_T
78  F77_FUNC (zqrinr, ZQRINR) (const octave_idx_type&, const octave_idx_type&,
79  Complex*, const octave_idx_type&, Complex*,
80  const octave_idx_type&, const octave_idx_type&,
81  const Complex*, double*);
82 
83  F77_RET_T
84  F77_FUNC (zqrder, ZQRDER) (const octave_idx_type&, const octave_idx_type&,
85  Complex*, const octave_idx_type&, Complex*,
86  const octave_idx_type&, const octave_idx_type&,
87  Complex*, double*);
88 
89  F77_RET_T
90  F77_FUNC (zqrshc, ZQRSHC) (const octave_idx_type&, const octave_idx_type&,
91  const octave_idx_type&, Complex*,
92  const octave_idx_type&, Complex*,
93  const octave_idx_type&, const octave_idx_type&,
94  const octave_idx_type&, Complex*, double*);
95 
96 #endif
97 }
98 
99 ComplexQR::ComplexQR (const ComplexMatrix& a, qr_type_t qr_type)
100 {
101  init (a, qr_type);
102 }
103 
104 void
105 ComplexQR::init (const ComplexMatrix& a, qr_type_t qr_type)
106 {
107  octave_idx_type m = a.rows ();
108  octave_idx_type n = a.cols ();
109 
110  octave_idx_type min_mn = m < n ? m : n;
111  OCTAVE_LOCAL_BUFFER (Complex, tau, min_mn);
112 
113  octave_idx_type info = 0;
114 
115  ComplexMatrix afact = a;
116  if (m > n && qr_type == qr_type_std)
117  afact.resize (m, m);
118 
119  if (m > 0)
120  {
121  // workspace query.
122  Complex clwork;
123  F77_XFCN (zgeqrf, ZGEQRF, (m, n, afact.fortran_vec (), m, tau,
124  &clwork, -1, info));
125 
126  // allocate buffer and do the job.
127  octave_idx_type lwork = clwork.real ();
128  lwork = std::max (lwork, static_cast<octave_idx_type> (1));
129  OCTAVE_LOCAL_BUFFER (Complex, work, lwork);
130  F77_XFCN (zgeqrf, ZGEQRF, (m, n, afact.fortran_vec (), m, tau,
131  work, lwork, info));
132  }
133 
134  form (n, afact, tau, qr_type);
135 }
136 
137 void ComplexQR::form (octave_idx_type n, ComplexMatrix& afact,
138  Complex *tau, qr_type_t qr_type)
139 {
140  octave_idx_type m = afact.rows (), min_mn = std::min (m, n);
141  octave_idx_type info;
142 
143  if (qr_type == qr_type_raw)
144  {
145  for (octave_idx_type j = 0; j < min_mn; j++)
146  {
147  octave_idx_type limit = j < min_mn - 1 ? j : min_mn - 1;
148  for (octave_idx_type i = limit + 1; i < m; i++)
149  afact.elem (i, j) *= tau[j];
150  }
151 
152  r = afact;
153  }
154  else
155  {
156  // Attempt to minimize copying.
157  if (m >= n)
158  {
159  // afact will become q.
160  q = afact;
161  octave_idx_type k = qr_type == qr_type_economy ? n : m;
162  r = ComplexMatrix (k, n);
163  for (octave_idx_type j = 0; j < n; j++)
164  {
165  octave_idx_type i = 0;
166  for (; i <= j; i++)
167  r.xelem (i, j) = afact.xelem (i, j);
168  for (; i < k; i++)
169  r.xelem (i, j) = 0;
170  }
171  afact = ComplexMatrix (); // optimize memory
172  }
173  else
174  {
175  // afact will become r.
176  q = ComplexMatrix (m, m);
177  for (octave_idx_type j = 0; j < m; j++)
178  for (octave_idx_type i = j + 1; i < m; i++)
179  {
180  q.xelem (i, j) = afact.xelem (i, j);
181  afact.xelem (i, j) = 0;
182  }
183  r = afact;
184  }
185 
186 
187  if (m > 0)
188  {
189  octave_idx_type k = q.columns ();
190  // workspace query.
191  Complex clwork;
192  F77_XFCN (zungqr, ZUNGQR, (m, k, min_mn, q.fortran_vec (), m, tau,
193  &clwork, -1, info));
194 
195  // allocate buffer and do the job.
196  octave_idx_type lwork = clwork.real ();
197  lwork = std::max (lwork, static_cast<octave_idx_type> (1));
198  OCTAVE_LOCAL_BUFFER (Complex, work, lwork);
199  F77_XFCN (zungqr, ZUNGQR, (m, k, min_mn, q.fortran_vec (), m, tau,
200  work, lwork, info));
201  }
202  }
203 }
204 
205 #ifdef HAVE_QRUPDATE
206 
207 void
208 ComplexQR::update (const ComplexColumnVector& u, const ComplexColumnVector& v)
209 {
210  octave_idx_type m = q.rows ();
211  octave_idx_type n = r.columns ();
212  octave_idx_type k = q.columns ();
213 
214  if (u.length () == m && v.length () == n)
215  {
216  ComplexColumnVector utmp = u, vtmp = v;
217  OCTAVE_LOCAL_BUFFER (Complex, w, k);
218  OCTAVE_LOCAL_BUFFER (double, rw, k);
219  F77_XFCN (zqr1up, ZQR1UP, (m, n, k, q.fortran_vec (),
220  m, r.fortran_vec (), k,
221  utmp.fortran_vec (), vtmp.fortran_vec (),
222  w, rw));
223  }
224  else
225  (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
226 }
227 
228 void
229 ComplexQR::update (const ComplexMatrix& u, const ComplexMatrix& v)
230 {
231  octave_idx_type m = q.rows ();
232  octave_idx_type n = r.columns ();
233  octave_idx_type k = q.columns ();
234 
235  if (u.rows () == m && v.rows () == n && u.cols () == v.cols ())
236  {
237  OCTAVE_LOCAL_BUFFER (Complex, w, k);
238  OCTAVE_LOCAL_BUFFER (double, rw, k);
239  for (volatile octave_idx_type i = 0; i < u.cols (); i++)
240  {
241  ComplexColumnVector utmp = u.column (i), vtmp = v.column (i);
242  F77_XFCN (zqr1up, ZQR1UP, (m, n, k, q.fortran_vec (),
243  m, r.fortran_vec (), k,
244  utmp.fortran_vec (), vtmp.fortran_vec (),
245  w, rw));
246  }
247  }
248  else
249  (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
250 }
251 
252 void
253 ComplexQR::insert_col (const ComplexColumnVector& u, octave_idx_type j)
254 {
255  octave_idx_type m = q.rows ();
256  octave_idx_type n = r.columns ();
257  octave_idx_type k = q.columns ();
258 
259  if (u.length () != m)
260  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
261  else if (j < 0 || j > n)
262  (*current_liboctave_error_handler) ("qrinsert: index out of range");
263  else
264  {
265  if (k < m)
266  {
267  q.resize (m, k+1);
268  r.resize (k+1, n+1);
269  }
270  else
271  {
272  r.resize (k, n+1);
273  }
274 
275  ComplexColumnVector utmp = u;
276  OCTAVE_LOCAL_BUFFER (double, rw, k);
277  F77_XFCN (zqrinc, ZQRINC, (m, n, k, q.fortran_vec (), q.rows (),
278  r.fortran_vec (), r.rows (), j + 1,
279  utmp.data (), rw));
280  }
281 }
282 
283 void
284 ComplexQR::insert_col (const ComplexMatrix& u, const Array<octave_idx_type>& j)
285 {
286  octave_idx_type m = q.rows ();
287  octave_idx_type n = r.columns ();
288  octave_idx_type k = q.columns ();
289 
290  Array<octave_idx_type> jsi;
291  Array<octave_idx_type> js = j.sort (jsi, 0, ASCENDING);
292  octave_idx_type nj = js.length ();
293  bool dups = false;
294  for (octave_idx_type i = 0; i < nj - 1; i++)
295  dups = dups && js(i) == js(i+1);
296 
297  if (dups)
298  (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
299  else if (u.length () != m || u.columns () != nj)
300  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
301  else if (nj > 0 && (js(0) < 0 || js(nj-1) > n))
302  (*current_liboctave_error_handler) ("qrinsert: index out of range");
303  else if (nj > 0)
304  {
305  octave_idx_type kmax = std::min (k + nj, m);
306  if (k < m)
307  {
308  q.resize (m, kmax);
309  r.resize (kmax, n + nj);
310  }
311  else
312  {
313  r.resize (k, n + nj);
314  }
315 
316  OCTAVE_LOCAL_BUFFER (double, rw, kmax);
317  for (volatile octave_idx_type i = 0; i < js.length (); i++)
318  {
319  octave_idx_type ii = i;
320  ComplexColumnVector utmp = u.column (jsi(i));
321  F77_XFCN (zqrinc, ZQRINC, (m, n + ii, std::min (kmax, k + ii),
322  q.fortran_vec (), q.rows (),
323  r.fortran_vec (), r.rows (), js(ii) + 1,
324  utmp.data (), rw));
325  }
326  }
327 }
328 
329 void
330 ComplexQR::delete_col (octave_idx_type j)
331 {
332  octave_idx_type m = q.rows ();
333  octave_idx_type k = r.rows ();
334  octave_idx_type n = r.columns ();
335 
336  if (j < 0 || j > n-1)
337  (*current_liboctave_error_handler) ("qrdelete: index out of range");
338  else
339  {
340  OCTAVE_LOCAL_BUFFER (double, rw, k);
341  F77_XFCN (zqrdec, ZQRDEC, (m, n, k, q.fortran_vec (), q.rows (),
342  r.fortran_vec (), r.rows (), j + 1, rw));
343 
344  if (k < m)
345  {
346  q.resize (m, k-1);
347  r.resize (k-1, n-1);
348  }
349  else
350  {
351  r.resize (k, n-1);
352  }
353  }
354 }
355 
356 void
357 ComplexQR::delete_col (const Array<octave_idx_type>& j)
358 {
359  octave_idx_type m = q.rows ();
360  octave_idx_type n = r.columns ();
361  octave_idx_type k = q.columns ();
362 
363  Array<octave_idx_type> jsi;
364  Array<octave_idx_type> js = j.sort (jsi, 0, DESCENDING);
365  octave_idx_type nj = js.length ();
366  bool dups = false;
367  for (octave_idx_type i = 0; i < nj - 1; i++)
368  dups = dups && js(i) == js(i+1);
369 
370  if (dups)
371  (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
372  else if (nj > 0 && (js(0) > n-1 || js(nj-1) < 0))
373  (*current_liboctave_error_handler) ("qrinsert: index out of range");
374  else if (nj > 0)
375  {
376  OCTAVE_LOCAL_BUFFER (double, rw, k);
377  for (volatile octave_idx_type i = 0; i < js.length (); i++)
378  {
379  octave_idx_type ii = i;
380  F77_XFCN (zqrdec, ZQRDEC, (m, n - ii, k == m ? k : k - ii,
381  q.fortran_vec (), q.rows (),
382  r.fortran_vec (), r.rows (),
383  js(ii) + 1, rw));
384  }
385  if (k < m)
386  {
387  q.resize (m, k - nj);
388  r.resize (k - nj, n - nj);
389  }
390  else
391  {
392  r.resize (k, n - nj);
393  }
394 
395  }
396 }
397 
398 void
399 ComplexQR::insert_row (const ComplexRowVector& u, octave_idx_type j)
400 {
401  octave_idx_type m = r.rows ();
402  octave_idx_type n = r.columns ();
403  octave_idx_type k = std::min (m, n);
404 
405  if (! q.is_square () || u.length () != n)
406  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
407  else if (j < 0 || j > m)
408  (*current_liboctave_error_handler) ("qrinsert: index out of range");
409  else
410  {
411  q.resize (m + 1, m + 1);
412  r.resize (m + 1, n);
413  ComplexRowVector utmp = u;
414  OCTAVE_LOCAL_BUFFER (double, rw, k);
415  F77_XFCN (zqrinr, ZQRINR, (m, n, q.fortran_vec (), q.rows (),
416  r.fortran_vec (), r.rows (),
417  j + 1, utmp.fortran_vec (), rw));
418 
419  }
420 }
421 
422 void
423 ComplexQR::delete_row (octave_idx_type j)
424 {
425  octave_idx_type m = r.rows ();
426  octave_idx_type n = r.columns ();
427 
428  if (! q.is_square ())
429  (*current_liboctave_error_handler) ("qrdelete: dimensions mismatch");
430  else if (j < 0 || j > m-1)
431  (*current_liboctave_error_handler) ("qrdelete: index out of range");
432  else
433  {
434  OCTAVE_LOCAL_BUFFER (Complex, w, m);
435  OCTAVE_LOCAL_BUFFER (double, rw, m);
436  F77_XFCN (zqrder, ZQRDER, (m, n, q.fortran_vec (), q.rows (),
437  r.fortran_vec (), r.rows (), j + 1,
438  w, rw));
439 
440  q.resize (m - 1, m - 1);
441  r.resize (m - 1, n);
442  }
443 }
444 
445 void
446 ComplexQR::shift_cols (octave_idx_type i, octave_idx_type j)
447 {
448  octave_idx_type m = q.rows ();
449  octave_idx_type k = r.rows ();
450  octave_idx_type n = r.columns ();
451 
452  if (i < 0 || i > n-1 || j < 0 || j > n-1)
453  (*current_liboctave_error_handler) ("qrshift: index out of range");
454  else
455  {
456  OCTAVE_LOCAL_BUFFER (Complex, w, k);
457  OCTAVE_LOCAL_BUFFER (double, rw, k);
458  F77_XFCN (zqrshc, ZQRSHC, (m, n, k,
459  q.fortran_vec (), q.rows (),
460  r.fortran_vec (), r.rows (),
461  i + 1, j + 1, w, rw));
462  }
463 }
464 
465 #else
466 
467 // Replacement update methods.
468 
469 void
470 ComplexQR::update (const ComplexColumnVector& u, const ComplexColumnVector& v)
471 {
472  warn_qrupdate_once ();
473 
474  octave_idx_type m = q.rows ();
475  octave_idx_type n = r.columns ();
476 
477  if (u.length () == m && v.length () == n)
478  {
479  init (q*r + ComplexMatrix (u) * ComplexMatrix (v).hermitian (),
480  get_type ());
481  }
482  else
483  (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
484 }
485 
486 void
487 ComplexQR::update (const ComplexMatrix& u, const ComplexMatrix& v)
488 {
489  warn_qrupdate_once ();
490 
491  octave_idx_type m = q.rows ();
492  octave_idx_type n = r.columns ();
493 
494  if (u.rows () == m && v.rows () == n && u.cols () == v.cols ())
495  {
496  init (q*r + u * v.hermitian (), get_type ());
497  }
498  else
499  (*current_liboctave_error_handler) ("qrupdate: dimensions mismatch");
500 }
501 
502 static
503 ComplexMatrix insert_col (const ComplexMatrix& a, octave_idx_type i,
504  const ComplexColumnVector& x)
505 {
506  ComplexMatrix retval (a.rows (), a.columns () + 1);
507  retval.assign (idx_vector::colon, idx_vector (0, i),
508  a.index (idx_vector::colon, idx_vector (0, i)));
509  retval.assign (idx_vector::colon, idx_vector (i), x);
510  retval.assign (idx_vector::colon, idx_vector (i+1, retval.columns ()),
511  a.index (idx_vector::colon, idx_vector (i, a.columns ())));
512  return retval;
513 }
514 
515 static
516 ComplexMatrix insert_row (const ComplexMatrix& a, octave_idx_type i,
517  const ComplexRowVector& x)
518 {
519  ComplexMatrix retval (a.rows () + 1, a.columns ());
520  retval.assign (idx_vector (0, i), idx_vector::colon,
521  a.index (idx_vector (0, i), idx_vector::colon));
522  retval.assign (idx_vector (i), idx_vector::colon, x);
523  retval.assign (idx_vector (i+1, retval.rows ()), idx_vector::colon,
524  a.index (idx_vector (i, a.rows ()), idx_vector::colon));
525  return retval;
526 }
527 
528 static
529 ComplexMatrix delete_col (const ComplexMatrix& a, octave_idx_type i)
530 {
531  ComplexMatrix retval = a;
532  retval.delete_elements (1, idx_vector (i));
533  return retval;
534 }
535 
536 static
537 ComplexMatrix delete_row (const ComplexMatrix& a, octave_idx_type i)
538 {
539  ComplexMatrix retval = a;
540  retval.delete_elements (0, idx_vector (i));
541  return retval;
542 }
543 
544 static
545 ComplexMatrix shift_cols (const ComplexMatrix& a,
546  octave_idx_type i, octave_idx_type j)
547 {
548  octave_idx_type n = a.columns ();
549  Array<octave_idx_type> p (n);
550  for (octave_idx_type k = 0; k < n; k++) p(k) = k;
551  if (i < j)
552  {
553  for (octave_idx_type k = i; k < j; k++) p(k) = k+1;
554  p(j) = i;
555  }
556  else if (j < i)
557  {
558  p(j) = i;
559  for (octave_idx_type k = j+1; k < i+1; k++) p(k) = k-1;
560  }
561 
562  return a.index (idx_vector::colon, idx_vector (p));
563 }
564 
565 void
566 ComplexQR::insert_col (const ComplexColumnVector& u, octave_idx_type j)
567 {
568  warn_qrupdate_once ();
569 
570  octave_idx_type m = q.rows ();
571  octave_idx_type n = r.columns ();
572 
573  if (u.length () != m)
574  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
575  else if (j < 0 || j > n)
576  (*current_liboctave_error_handler) ("qrinsert: index out of range");
577  else
578  {
579  init (::insert_col (q*r, j, u), get_type ());
580  }
581 }
582 
583 void
584 ComplexQR::insert_col (const ComplexMatrix& u, const Array<octave_idx_type>& j)
585 {
586  warn_qrupdate_once ();
587 
588  octave_idx_type m = q.rows ();
589  octave_idx_type n = r.columns ();
590 
591  Array<octave_idx_type> jsi;
592  Array<octave_idx_type> js = j.sort (jsi, 0, ASCENDING);
593  octave_idx_type nj = js.length ();
594  bool dups = false;
595  for (octave_idx_type i = 0; i < nj - 1; i++)
596  dups = dups && js(i) == js(i+1);
597 
598  if (dups)
599  (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
600  else if (u.length () != m || u.columns () != nj)
601  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
602  else if (nj > 0 && (js(0) < 0 || js(nj-1) > n))
603  (*current_liboctave_error_handler) ("qrinsert: index out of range");
604  else if (nj > 0)
605  {
606  ComplexMatrix a = q*r;
607  for (octave_idx_type i = 0; i < js.length (); i++)
608  a = ::insert_col (a, js(i), u.column (i));
609  init (a, get_type ());
610  }
611 }
612 
613 void
614 ComplexQR::delete_col (octave_idx_type j)
615 {
616  warn_qrupdate_once ();
617 
618  octave_idx_type m = q.rows ();
619  octave_idx_type n = r.columns ();
620 
621  if (j < 0 || j > n-1)
622  (*current_liboctave_error_handler) ("qrdelete: index out of range");
623  else
624  {
625  init (::delete_col (q*r, j), get_type ());
626  }
627 }
628 
629 void
630 ComplexQR::delete_col (const Array<octave_idx_type>& j)
631 {
632  warn_qrupdate_once ();
633 
634  octave_idx_type m = q.rows ();
635  octave_idx_type n = r.columns ();
636 
637  Array<octave_idx_type> jsi;
638  Array<octave_idx_type> js = j.sort (jsi, 0, DESCENDING);
639  octave_idx_type nj = js.length ();
640  bool dups = false;
641  for (octave_idx_type i = 0; i < nj - 1; i++)
642  dups = dups && js(i) == js(i+1);
643 
644  if (dups)
645  (*current_liboctave_error_handler) ("qrinsert: duplicate index detected");
646  else if (nj > 0 && (js(0) > n-1 || js(nj-1) < 0))
647  (*current_liboctave_error_handler) ("qrinsert: index out of range");
648  else if (nj > 0)
649  {
650  ComplexMatrix a = q*r;
651  for (octave_idx_type i = 0; i < js.length (); i++)
652  a = ::delete_col (a, js(i));
653  init (a, get_type ());
654  }
655 }
656 
657 void
658 ComplexQR::insert_row (const ComplexRowVector& u, octave_idx_type j)
659 {
660  warn_qrupdate_once ();
661 
662  octave_idx_type m = r.rows ();
663  octave_idx_type n = r.columns ();
664 
665  if (! q.is_square () || u.length () != n)
666  (*current_liboctave_error_handler) ("qrinsert: dimensions mismatch");
667  else if (j < 0 || j > m)
668  (*current_liboctave_error_handler) ("qrinsert: index out of range");
669  else
670  {
671  init (::insert_row (q*r, j, u), get_type ());
672  }
673 }
674 
675 void
676 ComplexQR::delete_row (octave_idx_type j)
677 {
678  warn_qrupdate_once ();
679 
680  octave_idx_type m = r.rows ();
681  octave_idx_type n = r.columns ();
682 
683  if (! q.is_square ())
684  (*current_liboctave_error_handler) ("qrdelete: dimensions mismatch");
685  else if (j < 0 || j > m-1)
686  (*current_liboctave_error_handler) ("qrdelete: index out of range");
687  else
688  {
689  init (::delete_row (q*r, j), get_type ());
690  }
691 }
692 
693 void
694 ComplexQR::shift_cols (octave_idx_type i, octave_idx_type j)
695 {
696  warn_qrupdate_once ();
697 
698  octave_idx_type m = q.rows ();
699  octave_idx_type n = r.columns ();
700 
701  if (i < 0 || i > n-1 || j < 0 || j > n-1)
702  (*current_liboctave_error_handler) ("qrshift: index out of range");
703  else
704  {
705  init (::shift_cols (q*r, i, j), get_type ());
706  }
707 }
708 
709 #endif
Generated on Mon Dec 30 2013 03:04:48 for GNU Octave by   doxygen 1.8.1.2