#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
Go to the source code of this file.
Defines | |
#define | ABS(x) (x>0?x:-x) |
#define | EIGENVALUE(re, im) |
Functions | |
void | LAPACK_DGEESX (const char *jobvs, const char *sort, char *sel, const char *sense, const int *n, double *a, const int *lda, int *sdim, double *wr, double *wi, double *vs, const int *ldvs, double *rconde, double *rcondv, double *work, const int *lwork, int *iwork, const int *liwork, int *bwork, int *info) |
static PetscErrorCode | eigenvaluecomp (Mat A) |
static PetscErrorCode | Departure (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MinEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MinEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyRealRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyRealIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyImRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
static PetscErrorCode | MaxEVbyImIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *dummy, PetscTruth *flg) |
PetscErrorCode | RegisterLapackModules () |
PetscErrorCode RegisterLapackModules();
Compute these elements with
ComputeQuantity("lapack",<element>,A,(AnalysisItem*)&res,&flg);
Available elements are:
Definition in file lapack.c.
#define ABS | ( | x | ) | (x>0?x:-x) |
Referenced by eigenvaluecomp().
#define EIGENVALUE | ( | re, | |||
im | ) |
Value:
{ \ PetscReal mag = sqrt(re*re+im*im); \ if (!maxbyreset) {maxbyreset = 1; maxbyrere = re; maxbyreim = im; } \ else if (re>maxbyrere) { maxbyrere = re; maxbyreim = im; } \ if (!maxbyimset) {maxbyimset = 1; maxbyimre = re; maxbyimim = im; } \ else if (im>maxbyimre) { maxbyimre = re; maxbyimim = im; } \ if (!maxbymagset) {maxbymagset = 1; maxmag = mag; \ maxbymagre = re; maxbymagim = im; } \ else if (mag>maxmag) { maxbymagre = re; maxbymagim = im; } \ if (!minbymagset) {minbymagset = 1; minmag = mag; \ minbymagre = re; minbymagim = im; } \ else if (mag<minmag) { minbymagre = re; minbymagim = im; } \ }
Referenced by eigenvaluecomp().
static PetscErrorCode Departure | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 196 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00197 { 00198 Mat A = (Mat)prob; 00199 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00200 PetscFunctionBegin; 00201 ierr = GetDataID("lapack","departure",&id,&has); CHKERRQ(ierr); 00202 HASTOEXIST(has); 00203 ierr = PetscObjectComposedDataGetScalar 00204 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00205 if (*flg) { 00206 rv->r = v; 00207 } else { 00208 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00209 ierr = PetscObjectComposedDataGetScalar 00210 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00211 if (*flg) { 00212 rv->r = v; 00213 } else SETERRQ(1,"Should have computed departure by now"); 00214 } 00215 00216 PetscFunctionReturn(0); 00217 }
static PetscErrorCode eigenvaluecomp | ( | Mat | A | ) | [static] |
This is the computational routine for lapack eigenvalue calculation.
Definition at line 41 of file lapack.c.
References ABS, EIGENVALUE, GetDataID(), HASTOEXIST, id, and LAPACK_DGEESX().
Referenced by Departure(), MaxEVbyImIm(), MaxEVbyImRe(), MaxEVbyMagIm(), MaxEVbyMagRe(), MaxEVbyRealIm(), MaxEVbyRealRe(), MinEVbyMagIm(), and MinEVbyMagRe().
00042 { 00043 MPI_Comm comm; Mat Adense; 00044 int maxbyreset=0,maxbyimset=0,maxbymagset=0,minbymagset=0; 00045 PetscReal dep=0, 00046 maxbymagre,maxbymagim, minbymagre,minbymagim, 00047 maxbyrere,maxbyreim, maxbyimre,maxbyimim, 00048 maxmag,minmag; 00049 PetscErrorCode ierr; 00050 00051 PetscFunctionBegin; 00052 00053 ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); 00054 ierr = MatConvert(A,MATSEQDENSE,MAT_INITIAL_MATRIX,&Adense); CHKERRQ(ierr); 00055 00056 { 00057 char JOBVS = 'N', SORT = 'N', SENSE = 'N', SELECT='X'; 00058 int N,LDA,SDIM=0,LDVS=1,LWORK,LIWORK,*IWORK,INFO; 00059 PetscScalar *AA,*VS=NULL,*RCONDE=NULL,*RCONDV=NULL; 00060 PetscReal *WR,*WI,*WORK,*BWORK=NULL; 00061 00062 ierr = MatGetSize(A,&N,PETSC_NULL); CHKERRQ(ierr); 00063 LDA = N; 00064 ierr = MatGetArray(Adense,&AA); CHKERRQ(ierr); 00065 ierr = PetscMalloc(N*sizeof(PetscReal),&WR); CHKERRQ(ierr); 00066 ierr = PetscMalloc(N*sizeof(PetscReal),&WI); CHKERRQ(ierr); 00067 LWORK = 3*N+1; 00068 ierr = PetscMalloc(LWORK*sizeof(PetscReal),&WORK); CHKERRQ(ierr); 00069 LIWORK = 10; 00070 ierr = PetscMalloc(LIWORK*sizeof(int),&IWORK); CHKERRQ(ierr); 00071 CHKMEMQ; 00072 LAPACK_DGEESX(&JOBVS, &SORT, &SELECT, &SENSE, &N, AA, &LDA, &SDIM, 00073 WR, WI, VS, &LDVS, RCONDE, RCONDV, WORK, &LWORK, 00074 IWORK, &LIWORK, (int*)BWORK, &INFO); 00075 if (INFO<0) { 00076 INFO = -INFO; 00077 SETERRQ1(1,"DGEESX Problem with parameter %d",INFO); 00078 } else if (INFO>0) { 00079 SETERRQ1(1,"DGEESX error %d",INFO); 00080 } 00081 CHKMEMQ; 00082 { 00083 int i,j,rind=0; PetscReal f=0; 00084 for (i=0; i<N; i++) { 00085 #define ABS(x) (x>0?x:-x) 00086 #define EIGENVALUE(re,im) { \ 00087 PetscReal mag = sqrt(re*re+im*im); \ 00088 if (!maxbyreset) {maxbyreset = 1; maxbyrere = re; maxbyreim = im; } \ 00089 else if (re>maxbyrere) { maxbyrere = re; maxbyreim = im; } \ 00090 if (!maxbyimset) {maxbyimset = 1; maxbyimre = re; maxbyimim = im; } \ 00091 else if (im>maxbyimre) { maxbyimre = re; maxbyimim = im; } \ 00092 if (!maxbymagset) {maxbymagset = 1; maxmag = mag; \ 00093 maxbymagre = re; maxbymagim = im; } \ 00094 else if (mag>maxmag) { maxbymagre = re; maxbymagim = im; } \ 00095 if (!minbymagset) {minbymagset = 1; minmag = mag; \ 00096 minbymagre = re; minbymagim = im; } \ 00097 else if (mag<minmag) { minbymagre = re; minbymagim = im; } \ 00098 } 00099 PetscReal dd = AA[rind],db = ABS(AA[rind+1]),dt; 00100 dt = 1.e-14 * ABS(dd); 00101 /*printf("row %d, diagonal %e, below %e, test %e\n",i,dd,db,dt);*/ 00102 if (db<dt) { 00103 /* case of 1x1 eigenvalue block */ 00104 int ind=rind+N; 00105 for (j=i+1; j<N; j++) { 00106 PetscReal aa = AA[ind]; 00107 f += aa*aa; ind += N; 00108 } 00109 EIGENVALUE(dd,0.); 00110 } else { 00111 /* case of 2x2 eigenvalue block */ 00112 int ind; 00113 PetscReal dd = AA[rind],do1=AA[rind+1],do2=AA[rind+N],dn=AA[rind+N+1]; 00114 if (ABS((dd-dn)/dd)>1.e-13) 00115 SETERRQ4(1,"Strange 2x2 block %e,%e,%e,%e\n",dd,do1,do2,dn); 00116 EIGENVALUE(dd,sqrt(do1*do2)); 00117 ind = rind+2*N; 00118 for (j=i+2; j<N; j++) { 00119 PetscReal aa = AA[ind]; 00120 f += aa*aa; ind += N; 00121 } 00122 i++; rind += N+1; ind = rind+N; 00123 for (j=i+1; j<N; j++) { 00124 PetscReal aa = AA[ind]; 00125 f += aa*aa; ind += N; 00126 } 00127 } 00128 rind += N+1; 00129 } 00130 dep = sqrt(f); 00131 } 00132 ierr = PetscFree(WR); CHKERRQ(ierr); 00133 ierr = PetscFree(WI); CHKERRQ(ierr); 00134 } 00135 00136 ierr = MatDestroy(Adense); CHKERRQ(ierr); 00137 00138 { 00139 int id; PetscTruth has; 00140 ierr = GetDataID("lapack","departure",&id,&has); CHKERRQ(ierr); 00141 HASTOEXIST(has); 00142 ierr = PetscObjectComposedDataSetScalar 00143 ((PetscObject)A,id,dep); CHKERRQ(ierr); 00144 00145 ierr = GetDataID 00146 ("lapack","lambda-max-by-magnitude;re",&id,&has); CHKERRQ(ierr); 00147 HASTOEXIST(has); 00148 ierr = PetscObjectComposedDataSetReal 00149 ((PetscObject)A,id,maxbymagre); CHKERRQ(ierr); 00150 ierr = GetDataID 00151 ("lapack","lambda-max-by-magnitude;im",&id,&has); CHKERRQ(ierr); 00152 HASTOEXIST(has); 00153 ierr = PetscObjectComposedDataSetReal 00154 ((PetscObject)A,id,maxbymagim); CHKERRQ(ierr); 00155 00156 ierr = GetDataID 00157 ("lapack","lambda-min-by-magnitude;re",&id,&has); CHKERRQ(ierr); 00158 HASTOEXIST(has); 00159 ierr = PetscObjectComposedDataSetReal 00160 ((PetscObject)A,id,minbymagre); CHKERRQ(ierr); 00161 ierr = GetDataID 00162 ("lapack","lambda-min-by-magnitude;im",&id,&has); CHKERRQ(ierr); 00163 HASTOEXIST(has); 00164 ierr = PetscObjectComposedDataSetReal 00165 ((PetscObject)A,id,minbymagim); CHKERRQ(ierr); 00166 00167 ierr = GetDataID 00168 ("lapack","lambda-max-by-real-part;re",&id,&has); CHKERRQ(ierr); 00169 HASTOEXIST(has); 00170 ierr = PetscObjectComposedDataSetReal 00171 ((PetscObject)A,id,maxbyrere); CHKERRQ(ierr); 00172 ierr = GetDataID 00173 ("lapack","lambda-max-by-real-part;im",&id,&has); CHKERRQ(ierr); 00174 HASTOEXIST(has); 00175 ierr = PetscObjectComposedDataSetReal 00176 ((PetscObject)A,id,maxbyreim); CHKERRQ(ierr); 00177 00178 ierr = GetDataID 00179 ("lapack","lambda-max-by-im-part;re",&id,&has); CHKERRQ(ierr); 00180 HASTOEXIST(has); 00181 ierr = PetscObjectComposedDataSetReal 00182 ((PetscObject)A,id,maxbyimre); CHKERRQ(ierr); 00183 ierr = GetDataID 00184 ("lapack","lambda-max-by-im-part;im",&id,&has); CHKERRQ(ierr); 00185 HASTOEXIST(has); 00186 ierr = PetscObjectComposedDataSetReal 00187 ((PetscObject)A,id,maxbyimim); CHKERRQ(ierr); 00188 } 00189 00190 PetscFunctionReturn(0); 00191 }
void LAPACK_DGEESX | ( | const char * | jobvs, | |
const char * | sort, | |||
char * | sel, | |||
const char * | sense, | |||
const int * | n, | |||
double * | a, | |||
const int * | lda, | |||
int * | sdim, | |||
double * | wr, | |||
double * | wi, | |||
double * | vs, | |||
const int * | ldvs, | |||
double * | rconde, | |||
double * | rcondv, | |||
double * | work, | |||
const int * | lwork, | |||
int * | iwork, | |||
const int * | liwork, | |||
int * | bwork, | |||
int * | info | |||
) |
Referenced by eigenvaluecomp().
static PetscErrorCode MaxEVbyImIm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 397 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00398 { 00399 Mat A = (Mat)prob; 00400 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00401 PetscFunctionBegin; 00402 ierr = GetDataID 00403 ("lapack","lambda-max-by-im-part;im",&id,&has); CHKERRQ(ierr); 00404 HASTOEXIST(has); 00405 ierr = PetscObjectComposedDataGetReal 00406 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00407 if (*flg) rv->r = v; 00408 else { 00409 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00410 ierr = PetscObjectComposedDataGetScalar 00411 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00412 if (*flg) { 00413 rv->r = v; 00414 } else SETERRQ(1,"Should have computed lambda-max-by-im-part;im by now"); 00415 } 00416 PetscFunctionReturn(0); 00417 }
static PetscErrorCode MaxEVbyImRe | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 372 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00373 { 00374 Mat A = (Mat)prob; 00375 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00376 PetscFunctionBegin; 00377 ierr = GetDataID 00378 ("lapack","lambda-max-by-im-part;re",&id,&has); CHKERRQ(ierr); 00379 HASTOEXIST(has); 00380 ierr = PetscObjectComposedDataGetReal 00381 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00382 if (*flg) rv->r = v; 00383 else { 00384 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00385 ierr = PetscObjectComposedDataGetScalar 00386 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00387 if (*flg) { 00388 rv->r = v; 00389 } else SETERRQ(1,"Should have computed lambda-max-by-im-part;re by now"); 00390 } 00391 PetscFunctionReturn(0); 00392 }
static PetscErrorCode MaxEVbyMagIm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 247 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00248 { 00249 Mat A = (Mat)prob; 00250 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00251 PetscFunctionBegin; 00252 ierr = GetDataID 00253 ("lapack","lambda-max-by-magnitude;im",&id,&has); CHKERRQ(ierr); 00254 HASTOEXIST(has); 00255 ierr = PetscObjectComposedDataGetReal 00256 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00257 if (*flg) rv->r = v; 00258 else { 00259 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00260 ierr = PetscObjectComposedDataGetScalar 00261 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00262 if (*flg) { 00263 rv->r = v; 00264 } else SETERRQ(1,"Should have computed lambda-max-by-magnitude;im by now"); 00265 } 00266 PetscFunctionReturn(0); 00267 }
static PetscErrorCode MaxEVbyMagRe | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 222 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00223 { 00224 Mat A = (Mat)prob; 00225 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00226 PetscFunctionBegin; 00227 ierr = GetDataID 00228 ("lapack","lambda-max-by-magnitude;re",&id,&has); CHKERRQ(ierr); 00229 HASTOEXIST(has); 00230 ierr = PetscObjectComposedDataGetReal 00231 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00232 if (*flg) rv->r = v; 00233 else { 00234 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00235 ierr = PetscObjectComposedDataGetScalar 00236 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00237 if (*flg) { 00238 rv->r = v; 00239 } else SETERRQ(1,"Should have computed lambda-max-by-magnitude;re by now"); 00240 } 00241 PetscFunctionReturn(0); 00242 }
static PetscErrorCode MaxEVbyRealIm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 347 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00348 { 00349 Mat A = (Mat)prob; 00350 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00351 PetscFunctionBegin; 00352 ierr = GetDataID 00353 ("lapack","lambda-max-by-real-part;im",&id,&has); CHKERRQ(ierr); 00354 HASTOEXIST(has); 00355 ierr = PetscObjectComposedDataGetReal 00356 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00357 if (*flg) rv->r = v; 00358 else { 00359 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00360 ierr = PetscObjectComposedDataGetScalar 00361 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00362 if (*flg) { 00363 rv->r = v; 00364 } else SETERRQ(1,"Should have computed lambda-max-by-real-part;im by now"); 00365 } 00366 PetscFunctionReturn(0); 00367 }
static PetscErrorCode MaxEVbyRealRe | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 322 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00323 { 00324 Mat A = (Mat)prob; 00325 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00326 PetscFunctionBegin; 00327 ierr = GetDataID 00328 ("lapack","lambda-max-by-real-part;re",&id,&has); CHKERRQ(ierr); 00329 HASTOEXIST(has); 00330 ierr = PetscObjectComposedDataGetReal 00331 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00332 if (*flg) rv->r = v; 00333 else { 00334 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00335 ierr = PetscObjectComposedDataGetScalar 00336 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00337 if (*flg) { 00338 rv->r = v; 00339 } else SETERRQ(1,"Should have computed lambda-max-by-real-part;re by now"); 00340 } 00341 PetscFunctionReturn(0); 00342 }
static PetscErrorCode MinEVbyMagIm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 297 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00298 { 00299 Mat A = (Mat)prob; 00300 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00301 PetscFunctionBegin; 00302 ierr = GetDataID 00303 ("lapack","lambda-min-by-magnitude;im",&id,&has); CHKERRQ(ierr); 00304 HASTOEXIST(has); 00305 ierr = PetscObjectComposedDataGetReal 00306 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00307 if (*flg) rv->r = v; 00308 else { 00309 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00310 ierr = PetscObjectComposedDataGetScalar 00311 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00312 if (*flg) { 00313 rv->r = v; 00314 } else SETERRQ(1,"Should have computed lambda-min-by-magnitude;im by now"); 00315 } 00316 PetscFunctionReturn(0); 00317 }
static PetscErrorCode MinEVbyMagRe | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | dummy, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 272 of file lapack.c.
References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterLapackModules().
00273 { 00274 Mat A = (Mat)prob; 00275 PetscTruth has; int id; PetscErrorCode ierr; PetscReal v = 0; 00276 PetscFunctionBegin; 00277 ierr = GetDataID 00278 ("lapack","lambda-min-by-magnitude;re",&id,&has); CHKERRQ(ierr); 00279 HASTOEXIST(has); 00280 ierr = PetscObjectComposedDataGetReal 00281 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00282 if (*flg) rv->r = v; 00283 else { 00284 ierr = eigenvaluecomp(A); CHKERRQ(ierr); 00285 ierr = PetscObjectComposedDataGetScalar 00286 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00287 if (*flg) { 00288 rv->r = v; 00289 } else SETERRQ(1,"Should have computed lambda-min-by-magnitude;re by now"); 00290 } 00291 PetscFunctionReturn(0); 00292 }
PetscErrorCode RegisterLapackModules | ( | ) |
Declare norm-like modules that can be performed with lapack calculations.
Definition at line 424 of file lapack.c.
References ANALYSISDOUBLE, Departure(), MaxEVbyImIm(), MaxEVbyImRe(), MaxEVbyMagIm(), MaxEVbyMagRe(), MaxEVbyRealIm(), MaxEVbyRealRe(), MinEVbyMagIm(), MinEVbyMagRe(), and RegisterModule().
Referenced by AnaModRegisterStandardModules().
00425 { 00426 PetscErrorCode ierr; 00427 PetscFunctionBegin; 00428 00429 ierr = RegisterModule 00430 ("lapack","departure",ANALYSISDOUBLE,&Departure); CHKERRQ(ierr); 00431 00432 ierr = RegisterModule 00433 ("lapack","lambda-max-by-magnitude;re", 00434 ANALYSISDOUBLE,&MaxEVbyMagRe); CHKERRQ(ierr); 00435 ierr = RegisterModule 00436 ("lapack","lambda-max-by-magnitude;im", 00437 ANALYSISDOUBLE,&MaxEVbyMagIm); CHKERRQ(ierr); 00438 ierr = RegisterModule 00439 ("lapack","lambda-min-by-magnitude;re", 00440 ANALYSISDOUBLE,&MinEVbyMagRe); CHKERRQ(ierr); 00441 ierr = RegisterModule 00442 ("lapack","lambda-min-by-magnitude;im", 00443 ANALYSISDOUBLE,&MinEVbyMagIm); CHKERRQ(ierr); 00444 ierr = RegisterModule 00445 ("lapack","lambda-max-by-real-part;re", 00446 ANALYSISDOUBLE,&MaxEVbyRealRe); CHKERRQ(ierr); 00447 ierr = RegisterModule 00448 ("lapack","lambda-max-by-real-part;im", 00449 ANALYSISDOUBLE,&MaxEVbyRealIm); CHKERRQ(ierr); 00450 ierr = RegisterModule 00451 ("lapack","lambda-max-by-im-part;re", 00452 ANALYSISDOUBLE,&MaxEVbyImRe); CHKERRQ(ierr); 00453 ierr = RegisterModule 00454 ("lapack","lambda-max-by-im-part;im", 00455 ANALYSISDOUBLE,&MaxEVbyImIm); CHKERRQ(ierr); 00456 00457 PetscFunctionReturn(0); 00458 }