Mercurial > octave

comparison liboctave/numeric/EIG.cc @ 32660:f53ac65ffba6

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
maint: New method rwdata() as clearer alternative to fortran_vec(). * etc/NEWS.10.md: Announce new method and recommend its use in future code. * doc/interpreter/external.txi: Update Octave manual to use rwdata. * liboctave/array/Array.h (rwdata): New function prototype. * liboctave/array/Array.h (fortran_vec): Make inline function that calls rwdata(). * liboctave/array/DiagArray2.h (rwdata): New function prototype. * liboctave/array/DiagArray2.h (fortran_vec): Make inline function that calls rwdata(). * examples/code/fortrandemo.cc, libgui/graphics/QtHandlesUtils.cc, libinterp/corefcn/Cell.cc, libinterp/corefcn/__dsearchn__.cc, libinterp/corefcn/__ichol__.cc, libinterp/corefcn/__ilu__.cc, libinterp/corefcn/__lin_interpn__.cc, libinterp/corefcn/__magick_read__.cc, libinterp/corefcn/__pchip_deriv__.cc, libinterp/corefcn/amd.cc, libinterp/corefcn/data.cc, libinterp/corefcn/dot.cc, libinterp/corefcn/ellipj.cc, libinterp/corefcn/filter.cc, libinterp/corefcn/gcd.cc, libinterp/corefcn/gl-render.cc, libinterp/corefcn/graphics.cc, libinterp/corefcn/graphics.in.h, libinterp/corefcn/kron.cc, libinterp/corefcn/ls-mat4.cc, libinterp/corefcn/ls-mat5.cc, libinterp/corefcn/mex.cc, libinterp/corefcn/oct-map.cc, libinterp/corefcn/oct-stream.cc, libinterp/corefcn/ordqz.cc, libinterp/corefcn/ordschur.cc, libinterp/corefcn/perms.cc, libinterp/corefcn/psi.cc, libinterp/corefcn/quadcc.cc, libinterp/corefcn/qz.cc, libinterp/corefcn/rand.cc, libinterp/corefcn/sqrtm.cc, libinterp/corefcn/strfind.cc, libinterp/corefcn/symrcm.cc, libinterp/corefcn/tril.cc, libinterp/corefcn/typecast.cc, libinterp/corefcn/xdiv.cc, libinterp/dldfcn/__delaunayn__.cc, libinterp/dldfcn/__glpk__.cc, libinterp/dldfcn/__ode15__.cc, libinterp/dldfcn/__voronoi__.cc, libinterp/dldfcn/audioread.cc, libinterp/dldfcn/convhulln.cc, libinterp/octave-value/ov-base-int.cc, libinterp/octave-value/ov-bool-mat.cc, libinterp/octave-value/ov-ch-mat.cc, libinterp/octave-value/ov-cx-diag.cc, libinterp/octave-value/ov-cx-mat.cc, libinterp/octave-value/ov-flt-cx-diag.cc, libinterp/octave-value/ov-flt-cx-mat.cc, libinterp/octave-value/ov-flt-re-diag.cc, libinterp/octave-value/ov-flt-re-mat.cc, libinterp/octave-value/ov-intx.h, libinterp/octave-value/ov-java.cc, libinterp/octave-value/ov-perm.cc, libinterp/octave-value/ov-re-diag.cc, libinterp/octave-value/ov-re-mat.cc, libinterp/octave-value/ov-str-mat.cc, liboctave/array/Array-base.cc, liboctave/array/Array-util.cc, liboctave/array/CColVector.cc, liboctave/array/CDiagMatrix.cc, liboctave/array/CMatrix.cc, liboctave/array/CNDArray.cc, liboctave/array/CRowVector.cc, liboctave/array/CSparse.cc, liboctave/array/DiagArray2.h, liboctave/array/MArray.cc, liboctave/array/PermMatrix.cc, liboctave/array/Range.h, liboctave/array/Sparse.cc, liboctave/array/boolSparse.cc, liboctave/array/dColVector.cc, liboctave/array/dDiagMatrix.cc, liboctave/array/dMatrix.cc, liboctave/array/dNDArray.cc, liboctave/array/dRowVector.cc, liboctave/array/dSparse.cc, liboctave/array/fCColVector.cc, liboctave/array/fCDiagMatrix.cc, liboctave/array/fCMatrix.cc, liboctave/array/fCNDArray.cc, liboctave/array/fCRowVector.cc, liboctave/array/fColVector.cc, liboctave/array/fMatrix.cc, liboctave/array/fNDArray.cc, liboctave/array/fRowVector.cc, liboctave/array/idx-vector.cc, liboctave/numeric/CollocWt.cc, liboctave/numeric/DASPK.cc, liboctave/numeric/DASRT.cc, liboctave/numeric/DASSL.cc, liboctave/numeric/EIG.cc, liboctave/numeric/LSODE.cc, liboctave/numeric/Quad.cc, liboctave/numeric/aepbalance.cc, liboctave/numeric/bsxfun-defs.cc, liboctave/numeric/chol.cc, liboctave/numeric/eigs-base.cc, liboctave/numeric/fEIG.cc, liboctave/numeric/gepbalance.cc, liboctave/numeric/gsvd.cc, liboctave/numeric/hess.cc, liboctave/numeric/lu.cc, liboctave/numeric/oct-convn.cc, liboctave/numeric/oct-rand.cc, liboctave/numeric/qr.cc, liboctave/numeric/qrp.cc, liboctave/numeric/schur.cc, liboctave/numeric/sparse-dmsolve.cc, liboctave/numeric/sparse-lu.cc, liboctave/numeric/sparse-qr.cc, liboctave/numeric/svd.cc, liboctave/operators/mx-inlines.cc, liboctave/operators/mx-op-defs.h, liboctave/util/oct-base64.cc, liboctave/util/oct-binmap.h, liboctave/util/str-vec.cc: Replace calls to fortran_vec() with calls to rwdata().
author Rik <rik@octave.org>
date Wed, 27 Dec 2023 16:55:14 -0800
parents 2e484f9f1f18
children
comparison
equal deleted inserted replaced
32659:71a23564ceca 32660:f53ac65ffba6
51 (*current_liboctave_error_handler) ("EIG requires square matrix"); 51 (*current_liboctave_error_handler) ("EIG requires square matrix");
52 52
53 F77_INT info = 0; 53 F77_INT info = 0;
54 54
55 Matrix atmp = a; 55 Matrix atmp = a;
56 double *tmp_data = atmp.fortran_vec (); 56 double *tmp_data = atmp.rwdata ();
57 57
58 Array<double> wr (dim_vector (n, 1)); 58 Array<double> wr (dim_vector (n, 1));
59 double *pwr = wr.fortran_vec (); 59 double *pwr = wr.rwdata ();
60 60
61 Array<double> wi (dim_vector (n, 1)); 61 Array<double> wi (dim_vector (n, 1));
62 double *pwi = wi.fortran_vec (); 62 double *pwi = wi.rwdata ();
63 63
64 F77_INT tnvr = (calc_rev ? n : 0); 64 F77_INT tnvr = (calc_rev ? n : 0);
65 Matrix vr (tnvr, tnvr); 65 Matrix vr (tnvr, tnvr);
66 double *pvr = vr.fortran_vec (); 66 double *pvr = vr.rwdata ();
67 67
68 F77_INT tnvl = (calc_lev ? n : 0); 68 F77_INT tnvl = (calc_lev ? n : 0);
69 Matrix vl (tnvl, tnvl); 69 Matrix vl (tnvl, tnvl);
70 double *pvl = vl.fortran_vec (); 70 double *pvl = vl.rwdata ();
71 71
72 F77_INT lwork = -1; 72 F77_INT lwork = -1;
73 double dummy_work; 73 double dummy_work;
74 74
75 F77_INT ilo; 75 F77_INT ilo;
76 F77_INT ihi; 76 F77_INT ihi;
77 77
78 Array<double> scale (dim_vector (n, 1)); 78 Array<double> scale (dim_vector (n, 1));
79 double *pscale = scale.fortran_vec (); 79 double *pscale = scale.rwdata ();
80 80
81 double abnrm; 81 double abnrm;
82 82
83 Array<double> rconde (dim_vector (n, 1)); 83 Array<double> rconde (dim_vector (n, 1));
84 double *prconde = rconde.fortran_vec (); 84 double *prconde = rconde.rwdata ();
85 85
86 Array<double> rcondv (dim_vector (n, 1)); 86 Array<double> rcondv (dim_vector (n, 1));
87 double *prcondv = rcondv.fortran_vec (); 87 double *prcondv = rcondv.rwdata ();
88 88
89 F77_INT dummy_iwork; 89 F77_INT dummy_iwork;
90 90
91 F77_XFCN (dgeevx, DGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1), 91 F77_XFCN (dgeevx, DGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1),
92 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 92 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
104 if (info != 0) 104 if (info != 0)
105 (*current_liboctave_error_handler) ("dgeevx workspace query failed"); 105 (*current_liboctave_error_handler) ("dgeevx workspace query failed");
106 106
107 lwork = static_cast<F77_INT> (dummy_work); 107 lwork = static_cast<F77_INT> (dummy_work);
108 Array<double> work (dim_vector (lwork, 1)); 108 Array<double> work (dim_vector (lwork, 1));
109 double *pwork = work.fortran_vec (); 109 double *pwork = work.rwdata ();
110 110
111 F77_XFCN (dgeevx, DGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1), 111 F77_XFCN (dgeevx, DGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1),
112 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 112 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
113 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 113 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
114 F77_CONST_CHAR_ARG2 ("N", 1), 114 F77_CONST_CHAR_ARG2 ("N", 1),
184 (*current_liboctave_error_handler) ("EIG requires square matrix"); 184 (*current_liboctave_error_handler) ("EIG requires square matrix");
185 185
186 F77_INT info = 0; 186 F77_INT info = 0;
187 187
188 Matrix atmp = a; 188 Matrix atmp = a;
189 double *tmp_data = atmp.fortran_vec (); 189 double *tmp_data = atmp.rwdata ();
190 190
191 ColumnVector wr (n); 191 ColumnVector wr (n);
192 double *pwr = wr.fortran_vec (); 192 double *pwr = wr.rwdata ();
193 193
194 F77_INT lwork = -1; 194 F77_INT lwork = -1;
195 double dummy_work; 195 double dummy_work;
196 196
197 F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 197 F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
203 if (info != 0) 203 if (info != 0)
204 (*current_liboctave_error_handler) ("dsyev workspace query failed"); 204 (*current_liboctave_error_handler) ("dsyev workspace query failed");
205 205
206 lwork = static_cast<F77_INT> (dummy_work); 206 lwork = static_cast<F77_INT> (dummy_work);
207 Array<double> work (dim_vector (lwork, 1)); 207 Array<double> work (dim_vector (lwork, 1));
208 double *pwork = work.fortran_vec (); 208 double *pwork = work.rwdata ();
209 209
210 F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 210 F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
211 F77_CONST_CHAR_ARG2 ("U", 1), 211 F77_CONST_CHAR_ARG2 ("U", 1),
212 n, tmp_data, n, pwr, pwork, lwork, info 212 n, tmp_data, n, pwr, pwork, lwork, info
213 F77_CHAR_ARG_LEN (1) 213 F77_CHAR_ARG_LEN (1)
243 (*current_liboctave_error_handler) ("EIG requires square matrix"); 243 (*current_liboctave_error_handler) ("EIG requires square matrix");
244 244
245 F77_INT info = 0; 245 F77_INT info = 0;
246 246
247 ComplexMatrix atmp = a; 247 ComplexMatrix atmp = a;
248 Complex *tmp_data = atmp.fortran_vec (); 248 Complex *tmp_data = atmp.rwdata ();
249 249
250 ComplexColumnVector wr (n); 250 ComplexColumnVector wr (n);
251 Complex *pw = wr.fortran_vec (); 251 Complex *pw = wr.rwdata ();
252 252
253 F77_INT nvr = (calc_rev ? n : 0); 253 F77_INT nvr = (calc_rev ? n : 0);
254 ComplexMatrix vrtmp (nvr, nvr); 254 ComplexMatrix vrtmp (nvr, nvr);
255 Complex *pvr = vrtmp.fortran_vec (); 255 Complex *pvr = vrtmp.rwdata ();
256 256
257 F77_INT nvl = (calc_lev ? n : 0); 257 F77_INT nvl = (calc_lev ? n : 0);
258 ComplexMatrix vltmp (nvl, nvl); 258 ComplexMatrix vltmp (nvl, nvl);
259 Complex *pvl = vltmp.fortran_vec (); 259 Complex *pvl = vltmp.rwdata ();
260 260
261 F77_INT lwork = -1; 261 F77_INT lwork = -1;
262 Complex dummy_work; 262 Complex dummy_work;
263 263
264 F77_INT lrwork = 2*n; 264 F77_INT lrwork = 2*n;
265 Array<double> rwork (dim_vector (lrwork, 1)); 265 Array<double> rwork (dim_vector (lrwork, 1));
266 double *prwork = rwork.fortran_vec (); 266 double *prwork = rwork.rwdata ();
267 267
268 F77_INT ilo; 268 F77_INT ilo;
269 F77_INT ihi; 269 F77_INT ihi;
270 270
271 Array<double> scale (dim_vector (n, 1)); 271 Array<double> scale (dim_vector (n, 1));
272 double *pscale = scale.fortran_vec (); 272 double *pscale = scale.rwdata ();
273 273
274 double abnrm; 274 double abnrm;
275 275
276 Array<double> rconde (dim_vector (n, 1)); 276 Array<double> rconde (dim_vector (n, 1));
277 double *prconde = rconde.fortran_vec (); 277 double *prconde = rconde.rwdata ();
278 278
279 Array<double> rcondv (dim_vector (n, 1)); 279 Array<double> rcondv (dim_vector (n, 1));
280 double *prcondv = rcondv.fortran_vec (); 280 double *prcondv = rcondv.rwdata ();
281 281
282 F77_XFCN (zgeevx, ZGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1), 282 F77_XFCN (zgeevx, ZGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1),
283 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 283 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
284 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 284 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
285 F77_CONST_CHAR_ARG2 ("N", 1), 285 F77_CONST_CHAR_ARG2 ("N", 1),
297 if (info != 0) 297 if (info != 0)
298 (*current_liboctave_error_handler) ("zgeevx workspace query failed"); 298 (*current_liboctave_error_handler) ("zgeevx workspace query failed");
299 299
300 lwork = static_cast<F77_INT> (dummy_work.real ()); 300 lwork = static_cast<F77_INT> (dummy_work.real ());
301 Array<Complex> work (dim_vector (lwork, 1)); 301 Array<Complex> work (dim_vector (lwork, 1));
302 Complex *pwork = work.fortran_vec (); 302 Complex *pwork = work.rwdata ();
303 303
304 F77_XFCN (zgeevx, ZGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1), 304 F77_XFCN (zgeevx, ZGEEVX, (F77_CONST_CHAR_ARG2 (balance ? "B" : "N", 1),
305 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 305 F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
306 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 306 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
307 F77_CONST_CHAR_ARG2 ("N", 1), 307 F77_CONST_CHAR_ARG2 ("N", 1),
338 (*current_liboctave_error_handler) ("EIG requires square matrix"); 338 (*current_liboctave_error_handler) ("EIG requires square matrix");
339 339
340 F77_INT info = 0; 340 F77_INT info = 0;
341 341
342 ComplexMatrix atmp = a; 342 ComplexMatrix atmp = a;
343 Complex *tmp_data = atmp.fortran_vec (); 343 Complex *tmp_data = atmp.rwdata ();
344 344
345 ColumnVector wr (n); 345 ColumnVector wr (n);
346 double *pwr = wr.fortran_vec (); 346 double *pwr = wr.rwdata ();
347 347
348 F77_INT lwork = -1; 348 F77_INT lwork = -1;
349 Complex dummy_work; 349 Complex dummy_work;
350 350
351 F77_INT lrwork = 3*n; 351 F77_INT lrwork = 3*n;
352 Array<double> rwork (dim_vector (lrwork, 1)); 352 Array<double> rwork (dim_vector (lrwork, 1));
353 double *prwork = rwork.fortran_vec (); 353 double *prwork = rwork.rwdata ();
354 354
355 F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 355 F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
356 F77_CONST_CHAR_ARG2 ("U", 1), 356 F77_CONST_CHAR_ARG2 ("U", 1),
357 n, F77_DBLE_CMPLX_ARG (tmp_data), n, pwr, 357 n, F77_DBLE_CMPLX_ARG (tmp_data), n, pwr,
358 F77_DBLE_CMPLX_ARG (&dummy_work), lwork, 358 F77_DBLE_CMPLX_ARG (&dummy_work), lwork,
363 if (info != 0) 363 if (info != 0)
364 (*current_liboctave_error_handler) ("zheev workspace query failed"); 364 (*current_liboctave_error_handler) ("zheev workspace query failed");
365 365
366 lwork = static_cast<F77_INT> (dummy_work.real ()); 366 lwork = static_cast<F77_INT> (dummy_work.real ());
367 Array<Complex> work (dim_vector (lwork, 1)); 367 Array<Complex> work (dim_vector (lwork, 1));
368 Complex *pwork = work.fortran_vec (); 368 Complex *pwork = work.rwdata ();
369 369
370 F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 370 F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
371 F77_CONST_CHAR_ARG2 ("U", 1), 371 F77_CONST_CHAR_ARG2 ("U", 1),
372 n, F77_DBLE_CMPLX_ARG (tmp_data), n, pwr, 372 n, F77_DBLE_CMPLX_ARG (tmp_data), n, pwr,
373 F77_DBLE_CMPLX_ARG (pwork), lwork, prwork, info 373 F77_DBLE_CMPLX_ARG (pwork), lwork, prwork, info
408 (*current_liboctave_error_handler) ("EIG requires same size matrices"); 408 (*current_liboctave_error_handler) ("EIG requires same size matrices");
409 409
410 F77_INT info = 0; 410 F77_INT info = 0;
411 411
412 Matrix tmp = b; 412 Matrix tmp = b;
413 double *tmp_data = tmp.fortran_vec (); 413 double *tmp_data = tmp.rwdata ();
414 414
415 if (! force_qz) 415 if (! force_qz)
416 { 416 {
417 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), 417 F77_XFCN (dpotrf, DPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1),
418 n, tmp_data, n, 418 n, tmp_data, n,
422 if (a.issymmetric () && b.issymmetric () && info == 0) 422 if (a.issymmetric () && b.issymmetric () && info == 0)
423 return symmetric_init (a, b, calc_rev, calc_lev); 423 return symmetric_init (a, b, calc_rev, calc_lev);
424 } 424 }
425 425
426 Matrix atmp = a; 426 Matrix atmp = a;
427 double *atmp_data = atmp.fortran_vec (); 427 double *atmp_data = atmp.rwdata ();
428 428
429 Matrix btmp = b; 429 Matrix btmp = b;
430 double *btmp_data = btmp.fortran_vec (); 430 double *btmp_data = btmp.rwdata ();
431 431
432 Array<double> ar (dim_vector (n, 1)); 432 Array<double> ar (dim_vector (n, 1));
433 double *par = ar.fortran_vec (); 433 double *par = ar.rwdata ();
434 434
435 Array<double> ai (dim_vector (n, 1)); 435 Array<double> ai (dim_vector (n, 1));
436 double *pai = ai.fortran_vec (); 436 double *pai = ai.rwdata ();
437 437
438 Array<double> beta (dim_vector (n, 1)); 438 Array<double> beta (dim_vector (n, 1));
439 double *pbeta = beta.fortran_vec (); 439 double *pbeta = beta.rwdata ();
440 440
441 F77_INT tnvr = (calc_rev ? n : 0); 441 F77_INT tnvr = (calc_rev ? n : 0);
442 Matrix vr (tnvr, tnvr); 442 Matrix vr (tnvr, tnvr);
443 double *pvr = vr.fortran_vec (); 443 double *pvr = vr.rwdata ();
444 444
445 F77_INT tnvl = (calc_lev ? n : 0); 445 F77_INT tnvl = (calc_lev ? n : 0);
446 Matrix vl (tnvl, tnvl); 446 Matrix vl (tnvl, tnvl);
447 double *pvl = vl.fortran_vec (); 447 double *pvl = vl.rwdata ();
448 448
449 F77_INT lwork = -1; 449 F77_INT lwork = -1;
450 double dummy_work; 450 double dummy_work;
451 451
452 F77_XFCN (dggev, DGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 452 F77_XFCN (dggev, DGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
461 if (info != 0) 461 if (info != 0)
462 (*current_liboctave_error_handler) ("dggev workspace query failed"); 462 (*current_liboctave_error_handler) ("dggev workspace query failed");
463 463
464 lwork = static_cast<F77_INT> (dummy_work); 464 lwork = static_cast<F77_INT> (dummy_work);
465 Array<double> work (dim_vector (lwork, 1)); 465 Array<double> work (dim_vector (lwork, 1));
466 double *pwork = work.fortran_vec (); 466 double *pwork = work.rwdata ();
467 467
468 F77_XFCN (dggev, DGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 468 F77_XFCN (dggev, DGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
469 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 469 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
470 n, atmp_data, n, btmp_data, n, 470 n, atmp_data, n, btmp_data, n,
471 par, pai, pbeta, 471 par, pai, pbeta,
545 (*current_liboctave_error_handler) ("EIG requires same size matrices"); 545 (*current_liboctave_error_handler) ("EIG requires same size matrices");
546 546
547 F77_INT info = 0; 547 F77_INT info = 0;
548 548
549 Matrix atmp = a; 549 Matrix atmp = a;
550 double *atmp_data = atmp.fortran_vec (); 550 double *atmp_data = atmp.rwdata ();
551 551
552 Matrix btmp = b; 552 Matrix btmp = b;
553 double *btmp_data = btmp.fortran_vec (); 553 double *btmp_data = btmp.rwdata ();
554 554
555 ColumnVector wr (n); 555 ColumnVector wr (n);
556 double *pwr = wr.fortran_vec (); 556 double *pwr = wr.rwdata ();
557 557
558 F77_INT lwork = -1; 558 F77_INT lwork = -1;
559 double dummy_work; 559 double dummy_work;
560 560
561 F77_XFCN (dsygv, DSYGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 561 F77_XFCN (dsygv, DSYGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
569 if (info != 0) 569 if (info != 0)
570 (*current_liboctave_error_handler) ("dsygv workspace query failed"); 570 (*current_liboctave_error_handler) ("dsygv workspace query failed");
571 571
572 lwork = static_cast<F77_INT> (dummy_work); 572 lwork = static_cast<F77_INT> (dummy_work);
573 Array<double> work (dim_vector (lwork, 1)); 573 Array<double> work (dim_vector (lwork, 1));
574 double *pwork = work.fortran_vec (); 574 double *pwork = work.rwdata ();
575 575
576 F77_XFCN (dsygv, DSYGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 576 F77_XFCN (dsygv, DSYGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
577 F77_CONST_CHAR_ARG2 ("U", 1), 577 F77_CONST_CHAR_ARG2 ("U", 1),
578 n, atmp_data, n, 578 n, atmp_data, n,
579 btmp_data, n, 579 btmp_data, n,
615 (*current_liboctave_error_handler) ("EIG requires same size matrices"); 615 (*current_liboctave_error_handler) ("EIG requires same size matrices");
616 616
617 F77_INT info = 0; 617 F77_INT info = 0;
618 618
619 ComplexMatrix tmp = b; 619 ComplexMatrix tmp = b;
620 Complex *tmp_data = tmp.fortran_vec (); 620 Complex *tmp_data = tmp.rwdata ();
621 621
622 if (! force_qz) 622 if (! force_qz)
623 { 623 {
624 F77_XFCN (zpotrf, ZPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1), 624 F77_XFCN (zpotrf, ZPOTRF, (F77_CONST_CHAR_ARG2 ("L", 1),
625 n, F77_DBLE_CMPLX_ARG (tmp_data), n, 625 n, F77_DBLE_CMPLX_ARG (tmp_data), n,
629 if (a.ishermitian () && b.ishermitian () && info == 0) 629 if (a.ishermitian () && b.ishermitian () && info == 0)
630 return hermitian_init (a, b, calc_rev, calc_lev); 630 return hermitian_init (a, b, calc_rev, calc_lev);
631 } 631 }
632 632
633 ComplexMatrix atmp = a; 633 ComplexMatrix atmp = a;
634 Complex *atmp_data = atmp.fortran_vec (); 634 Complex *atmp_data = atmp.rwdata ();
635 635
636 ComplexMatrix btmp = b; 636 ComplexMatrix btmp = b;
637 Complex *btmp_data = btmp.fortran_vec (); 637 Complex *btmp_data = btmp.rwdata ();
638 638
639 ComplexColumnVector alpha (n); 639 ComplexColumnVector alpha (n);
640 Complex *palpha = alpha.fortran_vec (); 640 Complex *palpha = alpha.rwdata ();
641 641
642 ComplexColumnVector beta (n); 642 ComplexColumnVector beta (n);
643 Complex *pbeta = beta.fortran_vec (); 643 Complex *pbeta = beta.rwdata ();
644 644
645 F77_INT nvr = (calc_rev ? n : 0); 645 F77_INT nvr = (calc_rev ? n : 0);
646 ComplexMatrix vrtmp (nvr, nvr); 646 ComplexMatrix vrtmp (nvr, nvr);
647 Complex *pvr = vrtmp.fortran_vec (); 647 Complex *pvr = vrtmp.rwdata ();
648 648
649 F77_INT nvl = (calc_lev ? n : 0); 649 F77_INT nvl = (calc_lev ? n : 0);
650 ComplexMatrix vltmp (nvl, nvl); 650 ComplexMatrix vltmp (nvl, nvl);
651 Complex *pvl = vltmp.fortran_vec (); 651 Complex *pvl = vltmp.rwdata ();
652 652
653 F77_INT lwork = -1; 653 F77_INT lwork = -1;
654 Complex dummy_work; 654 Complex dummy_work;
655 655
656 F77_INT lrwork = 8*n; 656 F77_INT lrwork = 8*n;
657 Array<double> rwork (dim_vector (lrwork, 1)); 657 Array<double> rwork (dim_vector (lrwork, 1));
658 double *prwork = rwork.fortran_vec (); 658 double *prwork = rwork.rwdata ();
659 659
660 F77_XFCN (zggev, ZGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 660 F77_XFCN (zggev, ZGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
661 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 661 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
662 n, F77_DBLE_CMPLX_ARG (atmp_data), n, 662 n, F77_DBLE_CMPLX_ARG (atmp_data), n,
663 F77_DBLE_CMPLX_ARG (btmp_data), n, 663 F77_DBLE_CMPLX_ARG (btmp_data), n,
673 if (info != 0) 673 if (info != 0)
674 (*current_liboctave_error_handler) ("zggev workspace query failed"); 674 (*current_liboctave_error_handler) ("zggev workspace query failed");
675 675
676 lwork = static_cast<F77_INT> (dummy_work.real ()); 676 lwork = static_cast<F77_INT> (dummy_work.real ());
677 Array<Complex> work (dim_vector (lwork, 1)); 677 Array<Complex> work (dim_vector (lwork, 1));
678 Complex *pwork = work.fortran_vec (); 678 Complex *pwork = work.rwdata ();
679 679
680 F77_XFCN (zggev, ZGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1), 680 F77_XFCN (zggev, ZGGEV, (F77_CONST_CHAR_ARG2 (calc_lev ? "V" : "N", 1),
681 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 681 F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
682 n, F77_DBLE_CMPLX_ARG (atmp_data), n, 682 n, F77_DBLE_CMPLX_ARG (atmp_data), n,
683 F77_DBLE_CMPLX_ARG (btmp_data), n, 683 F77_DBLE_CMPLX_ARG (btmp_data), n,
723 (*current_liboctave_error_handler) ("EIG requires same size matrices"); 723 (*current_liboctave_error_handler) ("EIG requires same size matrices");
724 724
725 F77_INT info = 0; 725 F77_INT info = 0;
726 726
727 ComplexMatrix atmp = a; 727 ComplexMatrix atmp = a;
728 Complex *atmp_data = atmp.fortran_vec (); 728 Complex *atmp_data = atmp.rwdata ();
729 729
730 ComplexMatrix btmp = b; 730 ComplexMatrix btmp = b;
731 Complex *btmp_data = btmp.fortran_vec (); 731 Complex *btmp_data = btmp.rwdata ();
732 732
733 ColumnVector wr (n); 733 ColumnVector wr (n);
734 double *pwr = wr.fortran_vec (); 734 double *pwr = wr.rwdata ();
735 735
736 F77_INT lwork = -1; 736 F77_INT lwork = -1;
737 Complex dummy_work; 737 Complex dummy_work;
738 738
739 F77_INT lrwork = 3*n; 739 F77_INT lrwork = 3*n;
740 Array<double> rwork (dim_vector (lrwork, 1)); 740 Array<double> rwork (dim_vector (lrwork, 1));
741 double *prwork = rwork.fortran_vec (); 741 double *prwork = rwork.rwdata ();
742 742
743 F77_XFCN (zhegv, ZHEGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 743 F77_XFCN (zhegv, ZHEGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
744 F77_CONST_CHAR_ARG2 ("U", 1), 744 F77_CONST_CHAR_ARG2 ("U", 1),
745 n, F77_DBLE_CMPLX_ARG (atmp_data), n, 745 n, F77_DBLE_CMPLX_ARG (atmp_data), n,
746 F77_DBLE_CMPLX_ARG (btmp_data), n, 746 F77_DBLE_CMPLX_ARG (btmp_data), n,
752 if (info != 0) 752 if (info != 0)
753 (*current_liboctave_error_handler) ("zhegv workspace query failed"); 753 (*current_liboctave_error_handler) ("zhegv workspace query failed");
754 754
755 lwork = static_cast<F77_INT> (dummy_work.real ()); 755 lwork = static_cast<F77_INT> (dummy_work.real ());
756 Array<Complex> work (dim_vector (lwork, 1)); 756 Array<Complex> work (dim_vector (lwork, 1));
757 Complex *pwork = work.fortran_vec (); 757 Complex *pwork = work.rwdata ();
758 758
759 F77_XFCN (zhegv, ZHEGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1), 759 F77_XFCN (zhegv, ZHEGV, (1, F77_CONST_CHAR_ARG2 (calc_rev ? "V" : "N", 1),
760 F77_CONST_CHAR_ARG2 ("U", 1), 760 F77_CONST_CHAR_ARG2 ("U", 1),
761 n, F77_DBLE_CMPLX_ARG (atmp_data), n, 761 n, F77_DBLE_CMPLX_ARG (atmp_data), n,
762 F77_DBLE_CMPLX_ARG (btmp_data), n, 762 F77_DBLE_CMPLX_ARG (btmp_data), n,