SALSA Analysis Modules
|
Various estimates of spectrum related quantities. More...
#include <stdlib.h>
#include "string.h"
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscksp.h"
#include "petscmat.h"
Go to the source code of this file.
Functions | |
PetscErrorCode | SpectrumComputePreconditionedSpectrum () |
PetscErrorCode | SpectrumComputeUnpreconditionedSpectrum () |
static PetscErrorCode | fivestepplan (KSP solver, PetscInt it, PetscReal err, KSPConvergedReason *reason, void *ctx) |
static PetscErrorCode | compute_eigenvalues_petsc (Mat A, PetscReal **r, PetscReal **c, int *rneig, PetscReal *rx, PetscReal *rm, PetscBool *success, const char **reason) |
static PetscErrorCode | compute_eigenvalues (Mat A, PetscReal **r, PetscReal **c, int *neig, PetscBool *success, const char **reason) |
static PetscErrorCode | compute_singularvalues (Mat A, PetscBool *success, const char **reason) |
static PetscErrorCode | compute_ellipse_from_Ritz_values (Mat A, PetscReal *r, PetscReal *c, int neig) |
static PetscErrorCode | NRitzValues (AnaModNumericalProblem prob, AnalysisItem *iv, int *lv, PetscBool *flg) |
static PetscErrorCode | RitzValuesR (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | RitzValuesC (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SpectrumAX (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SpectrumAY (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SpectrumCX (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SpectrumCY (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | PosFraction (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SigmaMax (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SigmaMin (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | Kappa (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MinEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MinEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyRealRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyRealIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyImRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | MaxEVbyImIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SpectrumOptions (char *opt) |
PetscErrorCode | RegisterSpectrumModules () |
PetscErrorCode | DeregisterSpectrumModules () |
Variables | |
static PetscBool | trace_spectrum = PETSC_FALSE |
static PetscBool | trace_hessenberg = PETSC_FALSE |
static int | use_prec = 0 |
Various estimates of spectrum related quantities.
The spectrum (eigenvalues and singular values) are both very informative about a matrix, and hard to compute. The spectrum module provides various estimates that are derived from running the system through a GMRES iteration for a few dozen iterations.
Activate this module with
PetscErrorCode RegisterSpectrumModules();
Compute these elements with
ComputeQuantity("spectrum",<element>,A,(AnalysisItem*)&res,&reslen,&flg);
These elements can be computed for both a preconditioned and unpreconditioned matrix. Switch between the two with
PetscErrorCode SpectrumComputePreconditionedSpectrum(); PetscErrorCode SpectrumComputeUnpreconditionedSpectrum();
In the preconditioned case a preconditioner has to have been defined in the Petsc options database with prefix "-ana"; for instance "-ana_pc_type jacobi".
Available elements are:
Jacobi methods, such as GMRES, give a reduction $AV=VH$ of $A$ to upper Hessenberg form $H$. The reduction is defined by the Krylov basis $V$. While a full reduction gives a Hessenberg matrix that has the exact same eigenvalues as $A$, using $V$ with a limited number of columns $n<N$) gives an $H$ with eigenvalues that approximate those of $A$. In particular, the outer eigenvalues are known to be fairly accurate even for modest values of $n$.
This process is used in the Petsc calls KSPSetComputeEigenvalues() and KSPSetComputeSingularValues(), which we use in the compute_eigenvalues() function. Alternatively, if the SLEPc library is available, that can be used for a similar but probably more accurate computation.
See SpectrumOptions().
Definition in file spectrum.c.
static PetscErrorCode compute_eigenvalues | ( | Mat | A, |
PetscReal ** | r, | ||
PetscReal ** | c, | ||
int * | neig, | ||
PetscBool * | success, | ||
const char ** | reason | ||
) | [static] |
Compute a number of eigenvalues (returned as separate real and imaginary parts) of a matrix.
Parameters:
r
and c
arrays. This is zero if success is false. Definition at line 466 of file spectrum.c.
References compute_eigenvalues_petsc(), GetDataID(), HASTOEXIST, and id.
Referenced by NRitzValues(), PosFraction(), RitzValuesC(), RitzValuesR(), SpectrumAX(), SpectrumAY(), SpectrumCX(), and SpectrumCY().
{ #if !defined(HAVE_SLEPC) PetscReal emax,emin; #endif int id; PetscBool has; PetscErrorCode ierr; PetscFunctionBegin; #if defined(HAVE_SLEPC) ierr = compute_eigenvalues_slepc(A,r,c,neig,success,(char**)reason); CHKERRQ(ierr); #else ierr = compute_eigenvalues_petsc (A,r,c,neig,&emax,&emin,success,reason); CHKERRQ(ierr); #endif if (*success) { ierr = GetDataID("spectrum","n-ritz-values",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,*neig); CHKERRQ(ierr); ierr = GetDataID("spectrum","ritz-values-r",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetRealstar ((PetscObject)A,id,*r); CHKERRQ(ierr); ierr = GetDataID("spectrum","ritz-values-c",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetRealstar ((PetscObject)A,id,*c); CHKERRQ(ierr); } PetscFunctionReturn(0); }
static PetscErrorCode compute_eigenvalues_petsc | ( | Mat | A, |
PetscReal ** | r, | ||
PetscReal ** | c, | ||
int * | rneig, | ||
PetscReal * | rx, | ||
PetscReal * | rm, | ||
PetscBool * | success, | ||
const char ** | reason | ||
) | [static] |
Definition at line 296 of file spectrum.c.
References AnaModTraceMessage(), fivestepplan(), trace_hessenberg, trace_spectrum, and use_prec.
Referenced by compute_eigenvalues(), and compute_singularvalues().
{ KSP solver; Mat B; PC pc; Vec x=NULL,b=NULL; MPI_Comm comm; int maxit=60,local_size,neig; const PCType pctype; PetscErrorCode ierr; PetscFunctionBegin; *success = PETSC_TRUE; if (use_prec) { ierr = AnaModTraceMessage ("computing preconditioned spectrum\n"); CHKERRQ(ierr); } else { ierr = AnaModTraceMessage ("computing plain spectrum\n"); CHKERRQ(ierr); } /* * Create solver around A */ ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); ierr = KSPCreate(comm,&solver); CHKERRQ(ierr); ierr = KSPAppendOptionsPrefix(solver,"ana_"); CHKERRQ(ierr); ierr = KSPGetPC(solver,&pc); CHKERRQ(ierr); ierr = PetscObjectQuery ((PetscObject)A,"approximation",(PetscObject*)&B); CHKERRQ(ierr); if (!B) B = A; ierr = PCSetOperators(pc,A,B,DIFFERENT_NONZERO_PATTERN); CHKERRQ(ierr); /* * No pc, or get the current one from the options database */ if (!use_prec) { ierr = PetscOptionsClearValue("-ana_pc_type"); CHKERRQ(ierr); ierr = PCSetType(pc,PCNONE); CHKERRQ(ierr); pctype = "none"; } ierr = KSPSetFromOptions(solver); CHKERRQ(ierr); { /* we can not compute the spectrum for direct methods */ PC pc; const PCType type; PetscBool flg; ierr = KSPGetPC(solver,&pc); CHKERRQ(ierr); ierr = PCGetType(pc,&type); CHKERRQ(ierr); ierr = PetscStrcmp(type,PCLU,&flg); CHKERRQ(ierr); if (flg) { *success = PETSC_FALSE; *reason = "direct method"; goto exit;} } ierr = MatGetLocalSize(A,&local_size,&ierr); CHKERRQ(ierr); ierr = VecCreateMPI(comm,local_size,PETSC_DECIDE,&x); CHKERRQ(ierr); ierr = VecDuplicate(x,&b); CHKERRQ(ierr); /* * GMRES method for eigenvalues */ ierr = KSPSetType(solver,KSPGMRES); CHKERRQ(ierr); ierr = KSPGMRESSetRestart(solver,maxit); CHKERRQ(ierr); ierr = KSPSetTolerances(solver,1.e-10,1.e-30,1.e+6,maxit); CHKERRQ(ierr); ierr = KSPSetComputeEigenvalues(solver,PETSC_TRUE); CHKERRQ(ierr); ierr = KSPSetComputeSingularValues(solver,PETSC_TRUE); CHKERRQ(ierr); ierr = PetscOptionsSetValue ("-ana_ksp_gmres_modifiedgramschmidt",PETSC_NULL); CHKERRQ(ierr); { void *ctx; ierr = KSPDefaultConvergedCreate(&ctx); CHKERRQ(ierr); ierr = KSPSetConvergenceTest (solver,&fivestepplan,ctx,PETSC_NULL); CHKERRQ(ierr); } /* * run one system solution */ { PetscScalar one=1.; KSPConvergedReason convergereason; int its; ierr = VecSet(x,one); CHKERRQ(ierr); ierr = PetscPushErrorHandler(&PetscReturnErrorHandler,NULL); CHKERRQ(ierr); ierr = KSPSetFromOptions(solver); CHKERRQ(ierr); ierr = KSPSetUp(solver); if (!ierr) {ierr = KSPSolve(solver,x,b);} ierr = PetscPopErrorHandler(); CHKERRQ(ierr); if (ierr) { *success = PETSC_FALSE; *reason = "KSPSetup/Solve breakdown"; goto exit;} ierr = KSPGetConvergedReason(solver,&convergereason); CHKERRQ(ierr); ierr = KSPGetIterationNumber(solver,&its); CHKERRQ(ierr); printf("%d\n",its); if (convergereason < 0 && convergereason != KSP_DIVERGED_ITS) { *success = PETSC_FALSE; *reason = "KSPSolve divergence"; goto exit; } } /* * analyse the returns */ { PetscReal emax,emin, *real_part,*im_part; emax = 10; emin = 1; ierr = KSPComputeExtremeSingularValues(solver,&emax,&emin); CHKERRQ(ierr); *rx = emax; *rm = emin; ierr = PetscMalloc((maxit+1)*sizeof(PetscReal),&real_part); CHKERRQ(ierr); ierr = PetscMalloc((maxit+1)*sizeof(PetscReal),&im_part); CHKERRQ(ierr); neig = 0; ierr = KSPComputeEigenvalues (solver,maxit,real_part,im_part,&neig); CHKERRQ(ierr); *rneig = neig; if (trace_spectrum) { int i; printf("Spectrum calculated:\n"); for (i=0; i<neig; i++) printf("%d: %e+%e i\n",i,real_part[i],im_part[i]); } *r = real_part; *c = im_part; #ifdef PETSC_INTERNALS if (trace_hessenberg) { Mat H; void *solverdata = solver->data; KSP_GMRES *gmres = (KSP_GMRES*)solverdata; PetscInt /* n = gmres->it + 1,i,*/ N = gmres->max_k + 2; // PetscScalar *hess; ierr = MatCreateSeqDense (MPI_COMM_SELF,N,N,gmres->hh_origin,&H); CHKERRQ(ierr); if (trace_hessenberg) { ierr = MatView(H,0); CHKERRQ(ierr); } ierr = MatDestroy(H); CHKERRQ(ierr); /* ierr = PetscMalloc(N*N*sizeof(PetscScalar),&hess); CHKERRQ(ierr); ierr = PetscMemcpy (hess,gmres->hh_origin,N*N*sizeof(PetscScalar)); CHKERRQ(ierr); ierr = GetDataID("spectrum","hessenberg",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetRealstar ((PetscObject)A,id,hess); CHKERRQ(ierr); *h = hess; */ } #endif } exit: if (x) {ierr = VecDestroy(&x); CHKERRQ(ierr);} if (b) {ierr = VecDestroy(&b); CHKERRQ(ierr);} ierr = KSPDestroy(&solver); CHKERRQ(ierr); PetscFunctionReturn(0); }
static PetscErrorCode compute_ellipse_from_Ritz_values | ( | Mat | A, |
PetscReal * | r, | ||
PetscReal * | c, | ||
int | neig | ||
) | [static] |
Definition at line 544 of file spectrum.c.
References GetDataID(), HASTOEXIST, and id.
Referenced by NRitzValues(), PosFraction(), RitzValuesC(), RitzValuesR(), SpectrumAX(), SpectrumAY(), SpectrumCX(), and SpectrumCY().
{ PetscReal np,ax,ay,cx,cy; PetscReal mx,Mx,my,My, lambda_min,lambda_max; int npos,nneg, i,id; PetscBool has; PetscErrorCode ierr; PetscFunctionBegin; /* #positive vs #negative */ nneg = 0; for (i=0; i<neig; i++) if (r[i]<0) nneg++; npos = neig-nneg; if (nneg==0) np = 1; else np = npos/(1.*neig); ierr = GetDataID("spectrum","positive-fraction",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,np); CHKERRQ(ierr); /* shape of enclosing ellipse */ Mx = mx = My = my = 0.; lambda_max = lambda_min = 0; for (i=0; i<neig; i++) { double lambda = sqrt(r[i]*r[i]+c[i]*c[i]); if (lambda_min==0 || lambda<lambda_min) lambda_min = lambda; if (lambda>lambda_max) lambda_max = lambda; if (r[i]>Mx) { if (mx==Mx && mx==0) mx=Mx=r[i]; else Mx=r[i]; } if (r[i]<mx) { if (mx==Mx && mx==0) Mx=mx=r[i]; else mx=r[i]; } if (c[i]>My) { if (my==My && my==0) my=My=c[i]; else My=c[i]; } if (c[i]<my) { if (my==My && my==0) My=my=c[i]; else my=c[i]; } } ax = (Mx-mx)/2; ierr = GetDataID("spectrum","ellipse-ax",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,ax); CHKERRQ(ierr); ay = (My-my)/2; ierr = GetDataID("spectrum","ellipse-ay",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,ay); CHKERRQ(ierr); cx = (Mx+mx)/2; ierr = GetDataID("spectrum","ellipse-cx",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); /*printf("set cx=%e @ %d from %d\n",cx,id,(int)A);*/ ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,cx); CHKERRQ(ierr); cy = (My+my)/2; ierr = GetDataID("spectrum","ellipse-cy",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,cy); CHKERRQ(ierr); PetscFunctionReturn(0); }
static PetscErrorCode compute_singularvalues | ( | Mat | A, |
PetscBool * | success, | ||
const char ** | reason | ||
) | [static] |
Compute a number of singular values of a matrix.
Parameters:
Definition at line 512 of file spectrum.c.
References compute_eigenvalues_petsc(), GetDataID(), HASTOEXIST, and id.
Referenced by SigmaMax(), and SigmaMin().
{ PetscReal emax,emin; int id; PetscBool has; PetscErrorCode ierr; PetscFunctionBegin; #if defined(HAVE_SLEPC) ierr = compute_sigmas_slepc(A,&emax,&emin,success,(char**)reason); CHKERRQ(ierr); #else { PetscReal *r,*c; int neig; ierr = compute_eigenvalues_petsc (A,&r,&c,&neig,&emax,&emin,success,reason); CHKERRQ(ierr); } #endif if (*success) { ierr = GetDataID("spectrum","sigma-max",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,emax); CHKERRQ(ierr); ierr = GetDataID("spectrum","sigma-min",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,emin); CHKERRQ(ierr); } PetscFunctionReturn(0); }
PetscErrorCode DeregisterSpectrumModules | ( | void | ) |
Definition at line 1414 of file spectrum.c.
Referenced by AnaModDeregisterSalsaModules().
{ #if defined(HAVE_SLEPC) PetscErrorCode ierr; PetscFunctionBegin; ierr = SlepcFinalize(); CHKERRQ(ierr); PetscFunctionReturn(0); #else PetscFunctionBegin; PetscFunctionReturn(0); #endif }
static PetscErrorCode fivestepplan | ( | KSP | solver, |
PetscInt | it, | ||
PetscReal | err, | ||
KSPConvergedReason * | reason, | ||
void * | ctx | ||
) | [static] |
Definition at line 114 of file spectrum.c.
Referenced by compute_eigenvalues_petsc().
{ PetscErrorCode ierr; PetscFunctionBegin; if (it<5) { if (!it) { ierr = KSPDefaultConverged(solver,it,err,reason,ctx); CHKERRQ(ierr); } *reason = KSP_CONVERGED_ITERATING; } else { ierr = KSPDefaultConverged(solver,it,err,reason,ctx); CHKERRQ(ierr); } PetscFunctionReturn(0); }
static PetscErrorCode Kappa | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 970 of file spectrum.c.
References GetDataID(), HASTOEXIST, id, AnalysisItem::r, SigmaMax(), and SigmaMin().
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","kappa",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem ret; PetscReal sigmin,sigmax; ierr = SigmaMin(prob,&ret,lv,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); sigmin = ret.r; ierr = SigmaMax(prob,&ret,lv,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); sigmax = ret.r; rv->r = sigmax/sigmin; /* PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } */ } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyImIm | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1282 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-im-part-im",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { if (PetscAbsScalar(im[i])>mag) {mag = PetscAbsScalar(im[i]); N = i;} } v = im[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyImRe | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1244 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-im-part-re",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { if (PetscAbsScalar(im[i])>mag) {mag = PetscAbsScalar(im[i]); N = i;} } v = re[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyMagIm | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1051 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-magnitude-im",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { PetscReal m = sqrt(re[i]*re[i]+im[i]*im[i]); if (m>mag) {mag = m; N = i;} } v = im[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyMagRe | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1012 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-magnitude-re",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { PetscReal m = sqrt(re[i]*re[i]+im[i]*im[i]); if (m>mag) {mag = m; N = i;} } v = re[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyRealIm | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1206 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-real-part-im",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { if (PetscAbsScalar(re[i])>mag) {mag = PetscAbsScalar(re[i]); N = i;} } v = im[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MaxEVbyRealRe | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1168 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-max-by-real-part-re",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { if (PetscAbsScalar(re[i])>mag) {mag = PetscAbsScalar(re[i]); N = i;} } v = re[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MinEVbyMagIm | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1129 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-min-by-magnitude-im",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { PetscReal m = sqrt(re[i]*re[i]+im[i]*im[i]); if (mag==0 || m<mag) {mag = m; N = i;} } v = im[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode MinEVbyMagRe | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 1090 of file spectrum.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, id, AnalysisItem::r, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0; PetscFunctionBegin; ierr = GetDataID ("spectrum","lambda-min-by-magnitude-re",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { AnalysisItem q; int n,i,N; PetscReal *re,*im,mag; ierr = ComputeQuantity (prob,"spectrum","ritz-values-r",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); re = q.rr; ierr = ComputeQuantity (prob,"spectrum","ritz-values-c",&q,&n,flg); CHKERRQ(ierr); if (!*flg) PetscFunctionReturn(0); im = q.rr; *lv = n; N = 0; mag = 0; for (i=0; i<n; i++) { PetscReal m = sqrt(re[i]*re[i]+im[i]*im[i]); if (mag==0 || m<mag) {mag = m; N = i;} } v = re[N]; ierr = PetscObjectComposedDataSetReal((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode NRitzValues | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | iv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 612 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; int v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","n-ritz-values",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) iv->i = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) iv->i = v; } } PetscFunctionReturn(0); }
static PetscErrorCode PosFraction | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 886 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","positive-fraction",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } } PetscFunctionReturn(0); }
PetscErrorCode RegisterSpectrumModules | ( | void | ) |
Definition at line 1343 of file spectrum.c.
References ANALYSISDBLARRAY, ANALYSISDOUBLE, ANALYSISINTEGER, DeclareCategoryOptionFunction(), Kappa(), MaxEVbyImIm(), MaxEVbyImRe(), MaxEVbyMagIm(), MaxEVbyMagRe(), MaxEVbyRealIm(), MaxEVbyRealRe(), MinEVbyMagIm(), MinEVbyMagRe(), NRitzValues(), PosFraction(), RegisterModule(), RitzValuesC(), RitzValuesR(), SigmaMax(), SigmaMin(), SpectrumAX(), SpectrumAY(), SpectrumCX(), SpectrumCY(), and SpectrumOptions().
Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().
{ PetscErrorCode ierr; PetscFunctionBegin; #if defined(HAVE_SLEPC) ierr = SlepcInitialize(0,0,0,0); CHKERRQ(ierr); #endif ierr = RegisterModule ("spectrum","n-ritz-values",ANALYSISINTEGER,&NRitzValues); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","ritz-values-r",ANALYSISDBLARRAY,&RitzValuesR); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","ritz-values-c",ANALYSISDBLARRAY,&RitzValuesC); CHKERRQ(ierr); /* ierr = RegisterModule ("spectrum","hessenberg",ANALYSISDBLARRAY,&Hessenberg); CHKERRQ(ierr); */ ierr = RegisterModule ("spectrum","ellipse-ax",ANALYSISDOUBLE,&SpectrumAX); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","ellipse-ay",ANALYSISDOUBLE,&SpectrumAY); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","ellipse-cx",ANALYSISDOUBLE,&SpectrumCX); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","ellipse-cy",ANALYSISDOUBLE,&SpectrumCY); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","kappa",ANALYSISDOUBLE,&Kappa); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","positive-fraction",ANALYSISDOUBLE,&PosFraction); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","sigma-max",ANALYSISDOUBLE,&SigmaMax); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","sigma-min",ANALYSISDOUBLE,&SigmaMin); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-magnitude-re", ANALYSISDOUBLE,&MaxEVbyMagRe); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-magnitude-im", ANALYSISDOUBLE,&MaxEVbyMagIm); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-min-by-magnitude-re", ANALYSISDOUBLE,&MinEVbyMagRe); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-min-by-magnitude-im", ANALYSISDOUBLE,&MinEVbyMagIm); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-real-part-re", ANALYSISDOUBLE,&MaxEVbyRealRe); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-real-part-im", ANALYSISDOUBLE,&MaxEVbyRealIm); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-im-part-re", ANALYSISDOUBLE,&MaxEVbyImRe); CHKERRQ(ierr); ierr = RegisterModule ("spectrum","lambda-max-by-im-part-im", ANALYSISDOUBLE,&MaxEVbyImIm); CHKERRQ(ierr); ierr = DeclareCategoryOptionFunction ("spectrum",&SpectrumOptions); CHKERRQ(ierr); PetscFunctionReturn(0); }
static PetscErrorCode RitzValuesC | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 685 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal *v = NULL; PetscBool has; int id,id2; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ritz-values-c",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetRealstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); ierr = GetDataID("spectrum","n-ritz-values",&id2,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,*lv,*flg); CHKERRQ(ierr); if (*flg) rv->rr = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetRealstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,*lv,*flg); CHKERRQ(ierr); if (*flg) rv->rr = v; } } PetscFunctionReturn(0); }
static PetscErrorCode RitzValuesR | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 645 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::rr.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal *v = NULL; PetscBool has; int id,id2; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ritz-values-r",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetRealstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); ierr = GetDataID("spectrum","n-ritz-values",&id2,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,*lv,*flg); CHKERRQ(ierr); if (*flg) rv->rr = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r+1,c+1,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetRealstar ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id2,*lv,*flg); CHKERRQ(ierr); if (*flg) rv->rr = v; } } PetscFunctionReturn(0); }
static PetscErrorCode SigmaMax | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 918 of file spectrum.c.
References compute_singularvalues(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by Kappa(), and RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","sigma-max",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscBool success; const char *reason; ierr = compute_singularvalues(A,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode SigmaMin | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 944 of file spectrum.c.
References compute_singularvalues(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by Kappa(), and RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","sigma-min",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscBool success; const char *reason; ierr = compute_singularvalues(A,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode SpectrumAX | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 757 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ellipse-ax",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } } PetscFunctionReturn(0); }
static PetscErrorCode SpectrumAY | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 789 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ellipse-ay",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } } PetscFunctionReturn(0); }
PetscErrorCode SpectrumComputePreconditionedSpectrum | ( | void | ) |
Definition at line 95 of file spectrum.c.
References use_prec.
{ PetscFunctionBegin; use_prec = 1; PetscFunctionReturn(0); }
PetscErrorCode SpectrumComputeUnpreconditionedSpectrum | ( | void | ) |
Definition at line 104 of file spectrum.c.
References use_prec.
Referenced by AnaModRegisterSalsaModules().
{ PetscFunctionBegin; use_prec = 0; PetscFunctionReturn(0); }
static PetscErrorCode SpectrumCX | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 821 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ellipse-cx",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); /*printf("finding cx @ %d from %d: %d\n",id,(int)A,*flg);*/ if (*flg) rv->r = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } } PetscFunctionReturn(0); }
static PetscErrorCode SpectrumCY | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 854 of file spectrum.c.
References compute_eigenvalues(), compute_ellipse_from_Ritz_values(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSpectrumModules().
{ Mat A = (Mat)prob; PetscReal v = 0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("spectrum","ellipse-cy",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { PetscReal *r,*c; int neig; PetscBool success; const char *reason; ierr = compute_eigenvalues(A,&r,&c,&neig,&success,&reason); CHKERRQ(ierr); if (!success) { printf("Spectrum computation failure in <%s>: <%s>\n",__FUNCT__,reason); *flg = PETSC_FALSE; PetscFunctionReturn(0);} if (neig==0) { printf("Spectrum computation failed, reason unknown\n"); *flg = PETSC_FALSE; } else { ierr = compute_ellipse_from_Ritz_values(A,r,c,neig); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } } PetscFunctionReturn(0); }
static PetscErrorCode SpectrumOptions | ( | char * | opt | ) | [static] |
Spectrum options can be set with "-anamod_spectrum <someoption>"
. Here we process the options. The possibilities are:
"trace_spectrum"
: list the computed spectrum values on the screen"dump_hessenberg"
: generate a Matlab file of the Hessenberg matrixSee optionsfile.
Definition at line 1327 of file spectrum.c.
References trace_hessenberg, and trace_spectrum.
Referenced by RegisterSpectrumModules().
{ PetscBool flg; PetscErrorCode ierr; PetscFunctionBegin; ierr = PetscStrcmp(opt,"trace_spectrum",&flg); CHKERRQ(ierr); if (flg) trace_spectrum = PETSC_TRUE; ierr = PetscStrcmp(opt,"trace_hessenberg",&flg); CHKERRQ(ierr); if (flg) trace_hessenberg = PETSC_TRUE; PetscFunctionReturn(0); }
PetscBool trace_hessenberg = PETSC_FALSE [static] |
Definition at line 90 of file spectrum.c.
Referenced by compute_eigenvalues_petsc(), and SpectrumOptions().
PetscBool trace_spectrum = PETSC_FALSE [static] |
Definition at line 89 of file spectrum.c.
Referenced by compute_eigenvalues_petsc(), and SpectrumOptions().
int use_prec = 0 [static] |
Definition at line 92 of file spectrum.c.
Referenced by compute_eigenvalues_petsc(), SpectrumComputePreconditionedSpectrum(), and SpectrumComputeUnpreconditionedSpectrum().