Mercurial > hg > octave-nkf

comparison liboctave/EIG.cc @ 3585:d9803711e047

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2000-02-08 04:35:39 by jwe]
author jwe
date Tue, 08 Feb 2000 04:35:47 +0000
parents 09a3064a3a17
children 7da18459c08b
comparison
equal deleted inserted replaced
3584:aa31644d9779 3585:d9803711e047
158 int info = 0; 158 int info = 0;
159 159
160 Matrix atmp = a; 160 Matrix atmp = a;
161 double *tmp_data = atmp.fortran_vec (); 161 double *tmp_data = atmp.fortran_vec ();
162 162
163 Array<double> wr (n); 163 ColumnVector wr (n);
164 double *pwr = wr.fortran_vec (); 164 double *pwr = wr.fortran_vec ();
165 165
166 // XXX FIXME XXX -- it might be possible to choose a better value of 166 // XXX FIXME XXX -- it might be possible to choose a better value of
167 // lwork that would result in more efficient computations. 167 // lwork that would result in more efficient computations.
168 168
173 F77_XFCN (dsyev, DSYEV, ("V", "U", n, tmp_data, n, pwr, pwork, 173 F77_XFCN (dsyev, DSYEV, ("V", "U", n, tmp_data, n, pwr, pwork,
174 lwork, info, 1L, 1L)); 174 lwork, info, 1L, 1L));
175 175
176 if (f77_exception_encountered || info < 0) 176 if (f77_exception_encountered || info < 0)
177 (*current_liboctave_error_handler) ("unrecoverable error in dsyev"); 177 (*current_liboctave_error_handler) ("unrecoverable error in dsyev");
178 else 178 else if (info > 0)
179 { 179 (*current_liboctave_error_handler) ("dsyev failed to converge");
180 if (info > 0) 180 else
181 (*current_liboctave_error_handler) ("dsyev failed to converge"); 181 {
182 else 182 lambda = ComplexColumnVector (wr);
183 { 183 v = ComplexMatrix (atmp);
184 lambda.resize (n);
185
186 for (int j = 0; j < n; j++)
187 lambda.elem (j) = Complex (wr.elem (j));
188
189 v = atmp;
190 }
191 } 184 }
192 185
193 return info; 186 return info;
194 } 187 }
195 188
262 int info = 0; 255 int info = 0;
263 256
264 ComplexMatrix atmp = a; 257 ComplexMatrix atmp = a;
265 Complex *tmp_data = atmp.fortran_vec (); 258 Complex *tmp_data = atmp.fortran_vec ();
266 259
267 ColumnVector w (n); 260 ColumnVector wr (n);
268 double *pw = w.fortran_vec (); 261 double *pwr = wr.fortran_vec ();
269 262
270 // XXX FIXME XXX -- it might be possible to choose a better value of 263 // XXX FIXME XXX -- it might be possible to choose a better value of
271 // lwork that would result in more efficient computations. 264 // lwork that would result in more efficient computations.
272 265
273 int lwork = 8*n; 266 int lwork = 8*n;
276 269
277 int lrwork = 3*n; 270 int lrwork = 3*n;
278 Array<double> rwork (lrwork); 271 Array<double> rwork (lrwork);
279 double *prwork = rwork.fortran_vec (); 272 double *prwork = rwork.fortran_vec ();
280 273
281 F77_XFCN (zheev, ZHEEV, ("V", "U", n, tmp_data, n, pw, pwork, 274 F77_XFCN (zheev, ZHEEV, ("V", "U", n, tmp_data, n, pwr, pwork,
282 lwork, prwork, info, 1L, 1L)); 275 lwork, prwork, info, 1L, 1L));
283 276
284 if (f77_exception_encountered || info < 0) 277 if (f77_exception_encountered || info < 0)
285 (*current_liboctave_error_handler) ("unrecoverable error in zheev"); 278 (*current_liboctave_error_handler) ("unrecoverable error in zheev");
286 else if (info > 0) 279 else if (info > 0)
287 (*current_liboctave_error_handler) ("zheev failed to converge"); 280 (*current_liboctave_error_handler) ("zheev failed to converge");
288 else 281 else
289 { 282 {
290 lambda = w; 283 lambda = ComplexColumnVector (wr);
291 v = atmp; 284 v = ComplexMatrix (atmp);
292 } 285 }
293 286
294 return info; 287 return info;
295 } 288 }
296 289