| /* Remove the initial subspace */ |
| qep->nini = 0; |
| PetscFunctionReturn(0); |
| } |
| if (Vi) { ierr = VecSet(Vi,0.0);CHKERRQ(ierr); } |
| } |
| #endif |
| PetscFunctionReturn(0); |
| } |
| PetscFunctionBegin; |
| ierr = QEPComputeResidualNorm_Private(qep,kr,ki,xr,xi,&norm);CHKERRQ(ierr); |
| #ifndef PETSC_USE_COMPLEX |
| if (ki == 0 || |
| PetscAbsScalar(ki) < PetscAbsScalar(kr*PETSC_MACHINE_EPSILON)) { |
| } |
| } |
| #endif |
| PetscFunctionReturn(0); |
| } |
| @*/ |
| PetscErrorCode QEPCompareEigenvalues(QEP qep,PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *result) |
| { |
| PetscReal a,b; |
| PetscReal a,b; |
| PetscFunctionBegin; |
| switch(qep->which) { |
| PetscReal x,y; |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| if (!its) { |
| PetscDrawLG lg; |
| PetscErrorCode ierr; |
| PetscReal *x,*y; |
| PetscInt i; |
| int n = PetscMin(qep->nev,255); |
| PetscInt i,n = PetscMin(qep->nev,255); |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| if (!its) { |
| #undef __FUNCT__ |
| #define __FUNCT__ "QEPRegisterAll" |
| /*@C |
| QEPRegisterAll - Registers all the solvers in the QEP package. |
| QEPRegisterAll - Registers all the solvers in the QEP package. |
| Not Collective |
| Not Collective |
| Level: advanced |
| Level: advanced |
| .seealso: QEPRegisterDynamic() |
| @*/ |
| #endif |
| ierr = PetscLogEventEnd(QEP_Dense,0,0,0,0);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| #endif |
| } |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| if (qep->setupcalled) PetscFunctionReturn(0); |
| ierr = PetscLogEventBegin(QEP_SetUp,qep,0,0,0);CHKERRQ(ierr); |
| /* Set default solver type */ |
| ierr = PetscFree(qep->ISL);CHKERRQ(ierr); |
| } |
| } |
| ierr = PetscLogEventEnd(QEP_SetUp,qep,0,0,0);CHKERRQ(ierr); |
| qep->setupcalled = 1; |
| PetscFunctionReturn(0); |
| #include <private/qepimpl.h> /*I "slepcqep.h" I*/ |
| PetscFList QEPList = 0; |
| PetscClassId QEP_CLASSID = 0; |
| PetscLogEvent QEP_SetUp = 0, QEP_Solve = 0, QEP_Dense = 0; |
| PetscFList QEPList = 0; |
| PetscClassId QEP_CLASSID = 0; |
| PetscLogEvent QEP_SetUp = 0, QEP_Solve = 0, QEP_Dense = 0; |
| static PetscBool QEPPackageInitialized = PETSC_FALSE; |
| #undef __FUNCT__ |
| #define __FUNCT__ "QEPFinalizePackage" |
| /*@C |
| QEPFinalizePackage - This function destroys everything in the Slepc interface to the QEP package. It is |
| called from SlepcFinalize(). |
| QEPFinalizePackage - This function destroys everything in the Slepc interface |
| to the QEP package. It is called from SlepcFinalize(). |
| Level: developer |
| Level: developer |
| .seealso: SlepcFinalize() |
| @*/ |
| #undef __FUNCT__ |
| #define __FUNCT__ "QEPInitializePackage" |
| /*@C |
| QEPInitializePackage - This function initializes everything in the QEP package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to QEPCreate() |
| when using static libraries. |
| QEPInitializePackage - This function initializes everything in the QEP package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to QEPCreate() |
| when using static libraries. |
| Input Parameter: |
| path - The dynamic library path, or PETSC_NULL |
| Input Parameter: |
| . path - The dynamic library path, or PETSC_NULL |
| Level: developer |
| Level: developer |
| .seealso: SlepcInitialize() |
| @*/ |
| PetscErrorCode QEPInitializePackage(const char *path) { |
| PetscErrorCode QEPInitializePackage(const char *path) |
| { |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| PetscValidCharPointer(type,2); |
| ierr = PetscTypeCompare((PetscObject)qep,type,&match);CHKERRQ(ierr); |
| if (match) PetscFunctionReturn(0); |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| ierr = PetscFree(qep->data);CHKERRQ(ierr); |
| /* free work vectors */ |
| ierr = QEPDefaultFreeWork(qep);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| PetscInt i; |
| PetscFunctionBegin; |
| if (qep->nwork != nw) { |
| if (qep->nwork > 0) { |
| ierr = VecDestroyVecs(qep->nwork,&qep->work); CHKERRQ(ierr); |
| } |
| ierr = PetscLogObjectParents(qep,nw,qep->work); |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode QEPDefaultConverged(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx) |
| { |
| PetscReal w; |
| PetscFunctionBegin; |
| w = SlepcAbsEigenvalue(eigr,eigi); |
| *errest = res; |
| PetscFunctionBegin; |
| ncv = PetscBLASIntCast(qep->ncv); |
| nconv = PetscBLASIntCast(qep->nconv); |
| ierr = PetscMalloc(nconv*nconv*sizeof(PetscScalar),&Z);CHKERRQ(ierr); |
| ierr = PetscMalloc(3*nconv*sizeof(PetscScalar),&work);CHKERRQ(ierr); |
| #if defined(PETSC_USE_COMPLEX) |
| } |
| if (marker!=-1) k = marker; |
| *kout = k; |
| PetscFunctionReturn(0); |
| } |
| @*/ |
| PetscErrorCode QEPSetFromOptions(QEP qep) |
| { |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg,val; |
| PetscReal r; |
| PetscInt i,j,k; |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg,val; |
| PetscReal r; |
| PetscInt i,j,k; |
| PetscViewerASCIIMonitor monviewer; |
| SlepcConvMonitor ctx; |
| SlepcConvMonitor ctx; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| if( nev != PETSC_IGNORE ) { |
| if (nev<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0"); |
| qep->nev = nev; |
| #undef __FUNCT__ |
| #define __FUNCT__ "QEPSetTrackAll" |
| /*@ |
| QEPSetTrackAll - Specifies if the solver must compute the residual of all |
| approximate eigenpairs or not. |
| QEPSetTrackAll - Specifies if the solver must compute the residual of all |
| approximate eigenpairs or not. |
| Collective on QEP |
| Collective on QEP |
| Input Parameters: |
| + qep - the eigensolver context |
| - trackall - whether compute all residuals or not |
| Input Parameters: |
| + qep - the eigensolver context |
| - trackall - whether compute all residuals or not |
| Notes: |
| If the user sets trackall=PETSC_TRUE then the solver explicitly computes |
| the residual for each eigenpair approximation. Computing the residual is |
| usually an expensive operation and solvers commonly compute the associated |
| residual to the first unconverged eigenpair. |
| The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically |
| activates this option. |
| Notes: |
| If the user sets trackall=PETSC_TRUE then the solver explicitly computes |
| the residual for each eigenpair approximation. Computing the residual is |
| usually an expensive operation and solvers commonly compute the associated |
| residual to the first unconverged eigenpair. |
| The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically |
| activates this option. |
| Level: intermediate |
| Level: intermediate |
| .seealso: QEPGetTrackAll() |
| @*/ |
| #undef __FUNCT__ |
| #define __FUNCT__ "QEPGetTrackAll" |
| /*@ |
| QEPGetTrackAll - Returns the flag indicating whether all residual norms must |
| be computed or not. |
| QEPGetTrackAll - Returns the flag indicating whether all residual norms must |
| be computed or not. |
| Not Collective |
| Not Collective |
| Input Parameter: |
| . qep - the eigensolver context |
| Input Parameter: |
| . qep - the eigensolver context |
| Output Parameter: |
| . trackall - the returned flag |
| Output Parameter: |
| . trackall - the returned flag |
| Level: intermediate |
| Level: intermediate |
| .seealso: QEPSetTrackAll() |
| @*/ |
| PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| ierr = PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);CHKERRQ(ierr); |
| PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| ierr = PetscObjectAppendOptionsPrefix((PetscObject)qep, prefix);CHKERRQ(ierr); |
| PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[]) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(qep,QEP_CLASSID,1); |
| PetscValidPointer(prefix,2); |
| QEP_QARNOLDI *ctx = (QEP_QARNOLDI *)qep->data; |
| PetscFunctionBegin; |
| if (qep->ncv) { /* ncv set */ |
| if (qep->ncv<qep->nev) SETERRQ(((PetscObject)qep)->comm,1,"The value of ncv must be at least nev"); |
| } |
| ierr = KSPSetOperators(ctx->ksp,qep->M,qep->M,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); |
| ierr = KSPSetUp(ctx->ksp);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscBool breakdown; |
| PetscFunctionBegin; |
| ierr = PetscMemzero(S,qep->ncv*qep->ncv*sizeof(PetscScalar));CHKERRQ(ierr); |
| ierr = PetscMalloc(qep->ncv*qep->ncv*sizeof(PetscScalar),&Q);CHKERRQ(ierr); |
| lwork = 7*qep->ncv; |
| ierr = SlepcUpdateVectors(nv,qep->V,qep->nconv,k+l,Q,nv,PETSC_FALSE);CHKERRQ(ierr); |
| qep->nconv = k; |
| ierr = QEPMonitor(qep,qep->its,qep->nconv,qep->eigr,qep->eigi,qep->errest,nv);CHKERRQ(ierr); |
| } |
| for (j=0;j<qep->nconv;j++) { |
| #define __FUNCT__ "QEPSetUp_LINEAR" |
| PetscErrorCode QEPSetUp_LINEAR(QEP qep) |
| { |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscInt i=0; |
| EPSWhich which; |
| PetscBool trackall; |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscInt i=0; |
| EPSWhich which; |
| PetscBool trackall; |
| /* function tables */ |
| PetscErrorCode (*fcreate[][2])(MPI_Comm,QEP_LINEAR*,Mat*) = { |
| { MatCreateExplicit_QEPLINEAR_N1A, MatCreateExplicit_QEPLINEAR_N1B }, /* N1 */ |
| ierr = EPSGetDimensions(ctx->eps,PETSC_NULL,&qep->ncv,&qep->mpd);CHKERRQ(ierr); |
| ierr = EPSGetTolerances(ctx->eps,&qep->tol,&qep->max_it);CHKERRQ(ierr); |
| if (qep->nini>0 || qep->ninil>0) PetscInfo(qep,"Ignoring initial vectors\n"); |
| PetscFunctionReturn(0); |
| } |
| ierr = EPSGetOperators(eps,&A,PETSC_NULL);CHKERRQ(ierr); |
| ierr = MatGetVecs(A,&xr,PETSC_NULL);CHKERRQ(ierr); |
| ierr = VecDuplicate(xr,&xi);CHKERRQ(ierr); |
| ierr = VecCreateMPIWithArray(((PetscObject)qep)->comm,qep->nloc,qep->n,PETSC_NULL,&wr);CHKERRQ(ierr); |
| ierr = VecCreateMPIWithArray(((PetscObject)qep)->comm,qep->nloc,qep->n,PETSC_NULL,&wi);CHKERRQ(ierr); |
| for (i=0;i<qep->nconv;i++) { |
| ierr = VecDestroy(&wi);CHKERRQ(ierr); |
| ierr = VecDestroy(&xr);CHKERRQ(ierr); |
| ierr = VecDestroy(&xi);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| ierr = EPSGetOperators(eps,&A,PETSC_NULL);CHKERRQ(ierr); |
| ierr = MatGetVecs(A,&xr,PETSC_NULL);CHKERRQ(ierr); |
| ierr = VecDuplicate(xr,&xi);CHKERRQ(ierr); |
| ierr = VecCreateMPIWithArray(((PetscObject)qep)->comm,qep->nloc,qep->n,PETSC_NULL,&w);CHKERRQ(ierr); |
| for (i=0;i<qep->nconv;i++) { |
| ierr = EPSGetEigenpair(eps,i,&qep->eigr[i],&qep->eigi[i],xr,xi);CHKERRQ(ierr); |
| ierr = VecDestroy(&w);CHKERRQ(ierr); |
| ierr = VecDestroy(&xr);CHKERRQ(ierr); |
| ierr = VecDestroy(&xi);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscFunctionBegin; |
| ierr = PetscOptionsBegin(((PetscObject)qep)->comm,((PetscObject)qep)->prefix,"LINEAR Quadratic Eigenvalue Problem solver Options","QEP");CHKERRQ(ierr); |
| ierr = PetscOptionsInt("-qep_linear_cform","Number of the companion form","QEPLinearSetCompanionForm",ctx->cform,&i,&set);CHKERRQ(ierr); |
| if (set) { |
| ierr = QEPLinearSetCompanionForm(qep,i);CHKERRQ(ierr); |
| } |
| ierr = PetscOptionsBool("-qep_linear_explicitmatrix","Use explicit matrix in linearization","QEPLinearSetExplicitMatrix",ctx->explicitmatrix,&val,&set);CHKERRQ(ierr); |
| if (set) { |
| ierr = QEPLinearSetExplicitMatrix(qep,val);CHKERRQ(ierr); |
| ierr = EPSGetST(ctx->eps,&st);CHKERRQ(ierr); |
| ierr = STSetMatMode(st,ST_MATMODE_SHELL);CHKERRQ(ierr); |
| } |
| ierr = PetscOptionsEnd();CHKERRQ(ierr); |
| ctx->setfromoptionscalled = PETSC_TRUE; |
| PetscFunctionReturn(0); |
| #define __FUNCT__ "QEPLinearSetEPS_LINEAR" |
| PetscErrorCode QEPLinearSetEPS_LINEAR(QEP qep,EPS eps) |
| { |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,2); |
| #define __FUNCT__ "QEPView_LINEAR" |
| PetscErrorCode QEPView_LINEAR(QEP qep,PetscViewer viewer) |
| { |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscFunctionBegin; |
| if (ctx->explicitmatrix) { |
| #define __FUNCT__ "QEPDestroy_LINEAR" |
| PetscErrorCode QEPDestroy_LINEAR(QEP qep) |
| { |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscErrorCode ierr; |
| QEP_LINEAR *ctx = (QEP_LINEAR *)qep->data; |
| PetscFunctionBegin; |
| ierr = EPSDestroy(&ctx->eps);CHKERRQ(ierr); |
| #undef __FUNCT__ |
| #define __FUNCT__ "STComputeExplicitOperator" |
| /*@ |
| STComputeExplicitOperator - Computes the explicit operator associated |
| to the eigenvalue problem with the specified spectral transformation. |
| STComputeExplicitOperator - Computes the explicit operator associated |
| to the eigenvalue problem with the specified spectral transformation. |
| Collective on ST |
| Collective on ST |
| Input Parameter: |
| . st - the spectral transform context |
| Input Parameter: |
| . st - the spectral transform context |
| Output Parameter: |
| . mat - the explicit operator |
| Output Parameter: |
| . mat - the explicit operator |
| Notes: |
| This routine builds a matrix containing the explicit operator. For |
| example, in generalized problems with shift-and-invert spectral |
| transformation the result would be matrix (A - s B)^-1 B. |
| This computation is done by applying the operator to columns of the |
| identity matrix. Note that the result is a dense matrix. |
| Notes: |
| This routine builds a matrix containing the explicit operator. For |
| example, in generalized problems with shift-and-invert spectral |
| transformation the result would be matrix (A - s B)^-1 B. |
| Level: advanced |
| This computation is done by applying the operator to columns of the |
| identity matrix. Note that the result is a dense matrix. |
| Level: advanced |
| .seealso: STApply() |
| @*/ |
| PetscErrorCode STComputeExplicitOperator(ST st,Mat *mat) |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| PetscInfo(st,"Setting up new ST\n"); |
| if (st->setupcalled) PetscFunctionReturn(0); |
| ierr = PetscLogEventBegin(ST_SetUp,st,0,0,0);CHKERRQ(ierr); |
| if (st->ops->postsolve) { |
| ierr = (*st->ops->postsolve)(st);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| if (st->ops->backtr) { |
| ierr = (*st->ops->backtr)(st,n,eigr,eigi);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #include <private/stimpl.h> /*I "slepcst.h" I*/ |
| PetscClassId ST_CLASSID = 0; |
| PetscLogEvent ST_SetUp = 0, ST_Apply = 0, ST_ApplyTranspose = 0; |
| PetscClassId ST_CLASSID = 0; |
| PetscLogEvent ST_SetUp = 0, ST_Apply = 0, ST_ApplyTranspose = 0; |
| static PetscBool STPackageInitialized = PETSC_FALSE; |
| #undef __FUNCT__ |
| #define __FUNCT__ "STFinalizePackage" |
| /*@C |
| STFinalizePackage - This function destroys everything in the Slepc interface to the ST package. It is |
| called from SlepcFinalize(). |
| STFinalizePackage - This function destroys everything in the Slepc interface |
| to the ST package. It is called from SlepcFinalize(). |
| Level: developer |
| Level: developer |
| .seealso: SlepcFinalize() |
| @*/ |
| #undef __FUNCT__ |
| #define __FUNCT__ "STInitializePackage" |
| /*@C |
| STInitializePackage - This function initializes everything in the ST package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to STCreate() |
| when using static libraries. |
| STInitializePackage - This function initializes everything in the ST package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to STCreate() |
| when using static libraries. |
| Input Parameter: |
| path - The dynamic library path, or PETSC_NULL |
| Input Parameter: |
| . path - The dynamic library path, or PETSC_NULL |
| Level: developer |
| Level: developer |
| .seealso: SlepcInitialize() |
| @*/ |
| PetscErrorCode STInitializePackage(const char *path) { |
| PetscErrorCode STInitializePackage(const char *path) |
| { |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| *newst = st; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| PetscErrorCode STSetOperators(ST st,Mat A,Mat B) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| PetscValidHeaderSpecific(A,MAT_CLASSID,2); |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| if (shift) *shift = st->sigma; |
| if (shift) *shift = st->sigma; |
| PetscFunctionReturn(0); |
| } |
| ierr = VecDestroy(&st->D); CHKERRQ(ierr); |
| st->D = D; |
| if (!st->wb) { |
| ierr = VecDuplicate(st->D,&st->wb); CHKERRQ(ierr); |
| ierr = VecDuplicate(st->D,&st->wb);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "STRegister" |
| /*@C |
| STRegister - See STRegisterDynamic() |
| STRegister - See STRegisterDynamic() |
| Level: advanced |
| Level: advanced |
| @*/ |
| PetscErrorCode STRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(ST)) |
| { |
| } |
| ierr = PetscOptionsEnd();CHKERRQ(ierr); |
| ierr = KSPSetFromOptions(st->ksp);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "STMatShellMult" |
| PetscErrorCode STMatShellMult(Mat A,Vec x,Vec y) |
| static PetscErrorCode STMatShellMult(Mat A,Vec x,Vec y) |
| { |
| PetscErrorCode ierr; |
| ST ctx; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext(A,(void**)&ctx);CHKERRQ(ierr); |
| ierr = MatMult(ctx->A,x,y);CHKERRQ(ierr); |
| if (ctx->sigma != 0.0) { |
| if (ctx->B) { /* y = (A - sB) x */ |
| #undef __FUNCT__ |
| #define __FUNCT__ "STMatShellMultTranspose" |
| PetscErrorCode STMatShellMultTranspose(Mat A,Vec x,Vec y) |
| static PetscErrorCode STMatShellMultTranspose(Mat A,Vec x,Vec y) |
| { |
| PetscErrorCode ierr; |
| ST ctx; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext(A,(void**)&ctx);CHKERRQ(ierr); |
| ierr = MatMultTranspose(ctx->A,x,y);CHKERRQ(ierr); |
| if (ctx->sigma != 0.0) { |
| if (ctx->B) { /* y = (A - sB) x */ |
| #undef __FUNCT__ |
| #define __FUNCT__ "STMatShellGetDiagonal" |
| PetscErrorCode STMatShellGetDiagonal(Mat A,Vec diag) |
| static PetscErrorCode STMatShellGetDiagonal(Mat A,Vec diag) |
| { |
| PetscErrorCode ierr; |
| ST ctx; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext(A,(void**)&ctx);CHKERRQ(ierr); |
| ierr = MatGetDiagonal(ctx->A,diag);CHKERRQ(ierr); |
| if (ctx->sigma != 0.0) { |
| if (ctx->B) { |
| if ( (hasA && !st->B) || (hasA && hasB) ) { |
| ierr = MatShellSetOperation(*mat,MATOP_GET_DIAGONAL,(void(*)(void))STMatShellGetDiagonal);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| */ |
| PetscErrorCode STAssociatedKSPSolve(ST st,Vec b,Vec x) |
| { |
| PetscErrorCode ierr; |
| PetscInt its; |
| PetscErrorCode ierr; |
| PetscInt its; |
| KSPConvergedReason reason; |
| PetscFunctionBegin; |
| PetscErrorCode STBackTransform_Fold(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi) |
| { |
| PetscInt j; |
| PetscFunctionBegin; |
| PetscValidScalarPointer(eigr,3); |
| PetscValidScalarPointer(eigi,4); |
| const KSPType ksptype; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp,&pc);CHKERRQ(ierr); |
| ierr = KSPGetType(st->ksp,&ksptype);CHKERRQ(ierr); |
| ierr = PCGetType(pc,&pctype);CHKERRQ(ierr); |
| ierr = PCSetType(pc,PCREDUNDANT);CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| ST_FOLD *ctx; |
| PetscFunctionBegin; |
| ierr = PetscNew(ST_FOLD,&ctx); CHKERRQ(ierr); |
| PetscLogObjectMemory(st,sizeof(ST_FOLD)); |
| st->data = (void *) ctx; |
| st->ops->setfromoptions = STSetFromOptions_Fold; |
| st->ops->destroy = STDestroy_Fold; |
| st->checknullspace = 0; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscFunctionBegin; |
| ierr = MatGetLocalSize(st->B,&n,&m);CHKERRQ(ierr); |
| ierr = MatCreateShell(((PetscObject)st)->comm,n,m,PETSC_DETERMINE,PETSC_DETERMINE,st,B);CHKERRQ(ierr); |
| ierr = MatShellSetOperation(*B,MATOP_MULT,(void(*)(void))STBilinearMatMult_Cayley);CHKERRQ(ierr); |
| ierr = MatShellSetOperation(*B,MATOP_MULT,(void(*)(void))STBilinearMatMult_Cayley);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscInt j; |
| #ifndef PETSC_USE_COMPLEX |
| PetscScalar t,i,r; |
| #endif |
| PetscFunctionBegin; |
| PetscValidPointer(eigr,3); |
| #ifndef PETSC_USE_COMPLEX |
| PetscValidPointer(eigi,4); |
| for (j=0;j<n;j++) { |
| if (eigi[j] == 0.0) eigr[j] = (ctx->nu + eigr[j] * st->sigma) / (eigr[j] - 1.0); |
| } |
| } |
| #else |
| PetscFunctionBegin; |
| PetscValidPointer(eigr,2); |
| for (j=0;j<n;j++) { |
| eigr[j] = (ctx->nu + eigr[j] * st->sigma) / (eigr[j] - 1.0); |
| } |
| ST_CAYLEY *ctx = (ST_CAYLEY *) st->data; |
| PetscFunctionBegin; |
| ierr = MatDestroy(&st->mat);CHKERRQ(ierr); |
| /* if the user did not set the shift, use the target value */ |
| const KSPType ksptype; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp,&pc);CHKERRQ(ierr); |
| ierr = KSPGetType(st->ksp,&ksptype);CHKERRQ(ierr); |
| ierr = PCGetType(pc,&pctype);CHKERRQ(ierr); |
| ctx->nu_set = PETSC_FALSE; |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STCayleySetAntishift_C","STCayleySetAntishift_Cayley",STCayleySetAntishift_Cayley);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "STShellGetContext" |
| /*@C |
| STShellGetContext - Returns the user-provided context associated with a shell ST |
| STShellGetContext - Returns the user-provided context associated with a shell ST |
| Not Collective |
| Not Collective |
| Input Parameter: |
| . st - spectral transformation context |
| Input Parameter: |
| . st - spectral transformation context |
| Output Parameter: |
| . ctx - the user provided context |
| Output Parameter: |
| . ctx - the user provided context |
| Level: advanced |
| Level: advanced |
| Notes: |
| This routine is intended for use within various shell routines |
| Notes: |
| This routine is intended for use within various shell routines |
| .seealso: STShellSetContext() |
| @*/ |
| @*/ |
| PetscErrorCode STShellSetContext(ST st,void *ctx) |
| { |
| ST_Shell *shell = (ST_Shell*)st->data; |
| ST_Shell *shell = (ST_Shell*)st->data; |
| PetscErrorCode ierr; |
| PetscBool flg; |
| if (!shell->apply) SETERRQ(((PetscObject)st)->comm,PETSC_ERR_USER,"No apply() routine provided to Shell ST"); |
| PetscStackPush("STSHELL apply() user function"); |
| CHKMEMQ; |
| ierr = (*shell->apply)(st,x,y);CHKERRQ(ierr); |
| ierr = (*shell->apply)(st,x,y);CHKERRQ(ierr); |
| CHKMEMQ; |
| PetscStackPop; |
| PetscFunctionReturn(0); |
| if (!shell->applytrans) SETERRQ(((PetscObject)st)->comm,PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST"); |
| PetscStackPush("STSHELL applytranspose() user function"); |
| CHKMEMQ; |
| ierr = (*shell->applytrans)(st,x,y);CHKERRQ(ierr); |
| ierr = (*shell->applytrans)(st,x,y);CHKERRQ(ierr); |
| CHKMEMQ; |
| PetscStackPop; |
| PetscFunctionReturn(0); |
| if (shell->backtr) { |
| PetscStackPush("STSHELL backtransform() user function"); |
| CHKMEMQ; |
| ierr = (*shell->backtr)(st,n,eigr,eigi);CHKERRQ(ierr); |
| ierr = (*shell->backtr)(st,n,eigr,eigi);CHKERRQ(ierr); |
| CHKMEMQ; |
| PetscStackPop; |
| } |
| const KSPType ksptype; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp,&pc);CHKERRQ(ierr); |
| ierr = KSPGetType(st->ksp,&ksptype);CHKERRQ(ierr); |
| ierr = PCGetType(pc,&pctype);CHKERRQ(ierr); |
| ierr = PCSetType(pc,PCREDUNDANT);CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","STShellSetApply_Shell",STShellSetApply_Shell);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","STShellSetApplyTranspose_Shell",STShellSetApplyTranspose_Shell);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","STShellSetBackTransform_Shell",STShellSetBackTransform_Shell);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "SLEPcNotImplemented_Precond" |
| PetscErrorCode SLEPcNotImplemented_Precond(ST st, Vec x, Vec y) { |
| PetscErrorCode SLEPcNotImplemented_Precond(ST st, Vec x, Vec y) |
| { |
| SETERRQ(((PetscObject)st)->comm,PETSC_ERR_SUP,"Operation not implemented in STPRECOND."); |
| } |
| PetscBool t0, t1; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp, &pc); CHKERRQ(ierr); |
| ierr = PetscObjectGetType((PetscObject)pc, &pctype); CHKERRQ(ierr); |
| ierr = STPrecondGetMatForPC(st, &P); CHKERRQ(ierr); |
| ierr = PCSetType(pc, (t0 && t1)? PCJACOBI:PCNONE); CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* if the user did not set the shift, use the target value */ |
| if (!st->sigma_set) st->sigma = st->defsigma; |
| } |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* Nothing to be done if STSetUp has not been called yet */ |
| if (!st->setupcalled) PetscFunctionReturn(0); |
| st->sigma = newshift; |
| if (st->shift_matrix != ST_MATMODE_SHELL) { |
| ierr = STSetUp_Precond(st); CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| ST_PRECOND *data; |
| PetscFunctionBegin; |
| ierr = PetscNew(ST_PRECOND, &data); CHKERRQ(ierr); |
| st->data = data; |
| ierr = STPrecondSetKSPHasMat_Precond(st, PETSC_TRUE); CHKERRQ(ierr); |
| ierr = KSPSetType(st->ksp, KSPPREONLY); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STPrecondGetMatForPC_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STPrecondSetMatForPC_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STPrecondGetKSPHasMat_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)st,"STPrecondSetKSPHasMat_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscFree(st->data); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| ierr = KSPGetPC(st->ksp, &pc); CHKERRQ(ierr); |
| ierr = PCGetOperatorsSet(pc, PETSC_NULL, &flag); CHKERRQ(ierr); |
| if (flag) { |
| ierr = PCGetOperators(pc, PETSC_NULL, mat, PETSC_NULL); CHKERRQ(ierr); |
| } else |
| *mat = PETSC_NULL; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| PetscValidHeaderSpecific(mat,MAT_CLASSID,2); |
| ierr = KSPGetPC(st->ksp, &pc); CHKERRQ(ierr); |
| /* Yes, all these lines are needed to safely set mat as the preconditioner |
| matrix in pc */ |
| ierr = MatDestroy(&A);CHKERRQ(ierr); |
| ierr = MatDestroy(&mat);CHKERRQ(ierr); |
| ierr = STPrecondSetKSPHasMat(st, PETSC_TRUE); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #define __FUNCT__ "STPrecondSetKSPHasMat_Precond" |
| PetscErrorCode STPrecondSetKSPHasMat_Precond(ST st,PetscBool setmat) |
| { |
| ST_PRECOND *data = (ST_PRECOND*)st->data; |
| ST_PRECOND *data = (ST_PRECOND*)st->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| data->setmat = setmat; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #define __FUNCT__ "STPrecondGetKSPHasMat_Precond" |
| PetscErrorCode STPrecondGetKSPHasMat_Precond(ST st,PetscBool *setmat) |
| { |
| ST_PRECOND *data = (ST_PRECOND*)st->data; |
| ST_PRECOND *data = (ST_PRECOND*)st->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(st,ST_CLASSID,1); |
| *setmat = data->setmat; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #define __FUNCT__ "STBackTransform_Sinvert" |
| PetscErrorCode STBackTransform_Sinvert(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi) |
| { |
| PetscInt j; |
| PetscInt j; |
| #ifndef PETSC_USE_COMPLEX |
| PetscScalar t; |
| #endif |
| PetscFunctionBegin; |
| PetscValidPointer(eigr,2); |
| PetscValidPointer(eigi,3); |
| PetscValidPointer(eigr,3); |
| #ifndef PETSC_USE_COMPLEX |
| PetscValidPointer(eigi,4); |
| for (j=0;j<n;j++) { |
| if (eigi[j] == 0) eigr[j] = 1.0 / eigr[j] + st->sigma; |
| else { |
| } |
| } |
| #else |
| PetscFunctionBegin; |
| PetscValidPointer(eigr,2); |
| for (j=0;j<n;j++) { |
| eigr[j] = 1.0 / eigr[j] + st->sigma; |
| } |
| MatStructure flg; |
| PetscFunctionBegin; |
| /* Nothing to be done if STSetUp has not been called yet */ |
| if (!st->setupcalled) PetscFunctionReturn(0); |
| const KSPType ksptype; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp,&pc);CHKERRQ(ierr); |
| ierr = KSPGetType(st->ksp,&ksptype);CHKERRQ(ierr); |
| ierr = PCGetType(pc,&pctype);CHKERRQ(ierr); |
| ierr = PCSetType(pc,PCREDUNDANT);CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| { |
| PetscFunctionBegin; |
| st->data = 0; |
| st->ops->apply = STApply_Sinvert; |
| st->ops->getbilinearform = STGetBilinearForm_Default; |
| st->ops->applytrans = STApplyTranspose_Sinvert; |
| st->ops->setshift = STSetShift_Sinvert; |
| st->ops->view = STView_Default; |
| st->ops->setfromoptions = STSetFromOptions_Sinvert; |
| st->checknullspace = STCheckNullSpace_Default; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscErrorCode STBackTransform_Shift(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi) |
| { |
| PetscInt j; |
| PetscFunctionBegin; |
| PetscValidPointer(eigr,3); |
| for (j=0;j<n;j++) { |
| const KSPType ksptype; |
| PetscFunctionBegin; |
| ierr = KSPGetPC(st->ksp,&pc);CHKERRQ(ierr); |
| ierr = KSPGetType(st->ksp,&ksptype);CHKERRQ(ierr); |
| ierr = PCGetType(pc,&pctype);CHKERRQ(ierr); |
| ierr = PCSetType(pc,PCREDUNDANT);CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| #else |
| PetscErrorCode ierr; |
| PetscScalar *work; |
| PetscBLASInt itype = 1,*iwork,info,n, |
| liwork; |
| PetscBLASInt itype = 1,*iwork,info,n,liwork; |
| const char *jobz; |
| #if defined(PETSC_USE_COMPLEX) |
| PetscReal *rwork; |
| norm = BLASnrm2_(&mm,Y,&inc); |
| tmp = 1.0 / norm; |
| BLASscal_(&mm,&tmp,Y,&inc); |
| PetscFunctionReturn(0); |
| #endif |
| } |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if (eps->setupcalled) PetscFunctionReturn(0); |
| ierr = PetscLogEventBegin(EPS_SetUp,eps,0,0,0);CHKERRQ(ierr); |
| /* Set default solver type */ |
| ierr = STSetOperators(eps->OP,A,B);CHKERRQ(ierr); |
| eps->setupcalled = 0; /* so that next solve call will call setup */ |
| ierr = VecDestroy(&eps->D);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| .seealso: SlepcInitialize() |
| @*/ |
| PetscErrorCode EPSInitializePackage(const char *path) { |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if (target) { |
| *target = eps->target; |
| } |
| if (target) *target = eps->target; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| if (eps->nwork != nw) { |
| if (eps->nwork > 0) { |
| ierr = VecDestroyVecs(eps->nwork,&eps->work); CHKERRQ(ierr); |
| ierr = VecDuplicateVecs(eps->V[0],nw,&eps->work); CHKERRQ(ierr); |
| ierr = PetscLogObjectParents(eps,nw,eps->work); |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode EPSEigRelativeConverged(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx) |
| { |
| PetscReal w; |
| PetscFunctionBegin; |
| w = SlepcAbsEigenvalue(eigr,eigi); |
| *errest = res; |
| PetscErrorCode EPSNormRelativeConverged(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx) |
| { |
| PetscReal w; |
| PetscFunctionBegin; |
| w = SlepcAbsEigenvalue(eigr,eigi); |
| *errest = res / (eps->nrma + w*eps->nrmb); |
| PetscReal norm; |
| PetscFunctionBegin; |
| /* allocate workspace */ |
| ierr = VecDuplicate(V[0],&x);CHKERRQ(ierr); |
| ierr = VecDuplicate(V[0],&y);CHKERRQ(ierr); |
| @*/ |
| PetscErrorCode EPSSetFromOptions(EPS eps) |
| { |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg,val; |
| PetscReal r,nrma,nrmb; |
| PetscScalar s; |
| PetscInt i,j,k; |
| const char *bal_list[4] = { "none", "oneside", "twoside","user" }; |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg,val; |
| PetscReal r,nrma,nrmb; |
| PetscScalar s; |
| PetscInt i,j,k; |
| const char *bal_list[4] = { "none", "oneside", "twoside","user" }; |
| PetscViewerASCIIMonitor monviewer; |
| SlepcConvMonitor ctx; |
| SlepcConvMonitor ctx; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if( nev != PETSC_IGNORE ) { |
| if (nev<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0"); |
| eps->nev = nev; |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| switch (type) { |
| case EPS_HEP: |
| eps->isgeneralized = PETSC_FALSE; |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type"); |
| } |
| eps->problem_type = type; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| ierr = PetscObjectSetOptionsPrefix((PetscObject)eps, prefix);CHKERRQ(ierr); |
| PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| ierr = PetscObjectAppendOptionsPrefix((PetscObject)eps, prefix);CHKERRQ(ierr); |
| PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[]) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| PetscValidPointer(prefix,2); |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = EPSRegisterDynamic(EPSPOWER,path,"EPSCreate_POWER",EPSCreate_POWER);CHKERRQ(ierr); |
| ierr = EPSRegisterDynamic(EPSSUBSPACE,path,"EPSCreate_SUBSPACE",EPSCreate_SUBSPACE);CHKERRQ(ierr); |
| ierr = EPSRegisterDynamic(EPSARNOLDI,path,"EPSCreate_ARNOLDI",EPSCreate_ARNOLDI);CHKERRQ(ierr); |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *ctx = (eps->monitorcontext[0]); |
| *ctx = eps->monitorcontext[0]; |
| PetscFunctionReturn(0); |
| } |
| PetscScalar er,ei; |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| if (!its) { |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSMonitorLGAll" |
| PetscErrorCode EPSMonitorLGAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx) |
| PetscScalar er,ei; |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| if (!its) { |
| typedef struct { |
| /* old values of eps */ |
| EPSWhich |
| old_which; |
| PetscErrorCode |
| (*old_which_func)(EPS,PetscScalar,PetscScalar,PetscScalar,PetscScalar, |
| PetscInt*,void*); |
| void |
| *old_which_ctx; |
| EPSWhich old_which; |
| PetscErrorCode (*old_which_func)(EPS,PetscScalar,PetscScalar,PetscScalar,PetscScalar, PetscInt*,void*); |
| void *old_which_ctx; |
| } EPSSortForSTData; |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSortForSTFunc" |
| PetscErrorCode EPSSortForSTFunc(EPS eps, PetscScalar ar, PetscScalar ai, |
| PetscScalar br, PetscScalar bi, PetscInt *r, |
| void *ctx) |
| PetscScalar br, PetscScalar bi, PetscInt *r, void *ctx) |
| { |
| EPSSortForSTData *data = (EPSSortForSTData*)ctx; |
| PetscErrorCode ierr; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* Back-transform the harmonic values */ |
| ierr = STBackTransform(eps->OP,1,&ar,&ai);CHKERRQ(ierr); |
| ierr = STBackTransform(eps->OP,1,&br,&bi);CHKERRQ(ierr); |
| eps->which = EPS_WHICH_USER; |
| eps->which_func = EPSSortForSTFunc; |
| eps->which_ctx = data; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSolve" |
| /*@ |
| /* Remove the initial subspaces */ |
| eps->nini = 0; |
| eps->ninil = 0; |
| PetscFunctionReturn(0); |
| } |
| } |
| ierr = EPSGetEigenvalue(eps,i,eigr,eigi);CHKERRQ(ierr); |
| ierr = EPSGetEigenvector(eps,i,Vr,Vi);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| if (i<0 || i>=eps->nconv) { |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE, "Argument 2 out of range"); |
| } |
| if (!eps->perm) k = i; |
| else k = eps->perm[i]; |
| #ifdef PETSC_USE_COMPLEX |
| if (eigr) *eigr = eps->eigr[k]; |
| if (eigi) *eigi = eps->eigi[k]; |
| #endif |
| PetscFunctionReturn(0); |
| } |
| if (!eps->evecsavailable && (Vr || Vi) ) { |
| ierr = (*eps->ops->computevectors)(eps);CHKERRQ(ierr); |
| } |
| if (!eps->perm) k = i; |
| else k = eps->perm[i]; |
| #ifdef PETSC_USE_COMPLEX |
| if (Vi) { ierr = VecSet(Vi,0.0); CHKERRQ(ierr); } |
| } |
| #endif |
| PetscFunctionReturn(0); |
| } |
| if (!eps->evecsavailable && (Wr || Wi) ) { |
| ierr = (*eps->ops->computevectors)(eps);CHKERRQ(ierr); |
| } |
| if (!eps->perm) k = i; |
| else k = eps->perm[i]; |
| #ifdef PETSC_USE_COMPLEX |
| if (Wi) { ierr = VecSet(Wi,0.0); CHKERRQ(ierr); } |
| } |
| #endif |
| PetscFunctionReturn(0); |
| } |
| { |
| PetscErrorCode ierr; |
| Vec x,y; |
| PetscBLASInt i; |
| PetscBLASInt i; |
| PetscFunctionBegin; |
| ierr = VecCreateMPIWithArray(((PetscObject)globaleps)->comm,*n,PETSC_DECIDE,PETSC_NULL,&x);CHKERRQ(ierr); |
| PetscScalar *pV; |
| PetscFunctionBegin; |
| ncv = PetscBLASIntCast(eps->ncv); |
| n = PetscBLASIntCast(eps->nloc); |
| eps->reason = EPS_CONVERGED_TOL; |
| if (stat!=0) { SETERRQ1(((PetscObject)eps)->comm,PETSC_ERR_LIB,"Error in TRLAN (code=%d)",stat);} |
| PetscFunctionReturn(0); |
| } |
| EPS_ARPACK *ar = (EPS_ARPACK *)eps->data; |
| char bmat[1], howmny[] = "A"; |
| const char *which; |
| PetscBLASInt n, iparam[11], ipntr[14], ido, info, nev, ncv; |
| PetscBLASInt n, iparam[11], ipntr[14], ido, info, nev, ncv, fcomm; |
| PetscScalar sigmar, *pV, *resid; |
| Vec x, y, w = eps->work[0]; |
| Mat A; |
| PetscBool isSinv, isShift, rvec; |
| PetscBLASInt fcomm; |
| #if !defined(PETSC_USE_COMPLEX) |
| PetscScalar sigmai = 0.0; |
| #endif |
| ierr = VecDestroy(&x);CHKERRQ(ierr); |
| ierr = VecDestroy(&y);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PC pc; |
| PetscFunctionBegin; |
| ierr = VecCreateMPIWithArray(((PetscObject)eps)->comm,eps->nloc,PETSC_DECIDE,PETSC_NULL,&x);CHKERRQ(ierr); |
| ierr = VecCreateMPIWithArray(((PetscObject)eps)->comm,eps->nloc,PETSC_DECIDE,PETSC_NULL,&y);CHKERRQ(ierr); |
| ierr = VecGetArray(eps->V[0],&pV);CHKERRQ(ierr); |
| lflag = 0; /* reverse communication interface flag */ |
| do { |
| BLZpack_( blz->istor, blz->rstor, &sigma, &nneig, blz->u, blz->v, |
| &lflag, &nvopu, blz->eig, pV ); |
| } |
| ierr = VecDestroy(&x);CHKERRQ(ierr); |
| ierr = VecDestroy(&y);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| KSP ksp; |
| } EPS_BLOPEX; |
| #undef __FUNCT__ |
| #define __FUNCT__ "Precond_FnSingleVector" |
| static void Precond_FnSingleVector(void *data,void *x,void *y) |
| { |
| PetscErrorCode ierr; |
| EPS eps = (EPS)data; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX*)eps->data; |
| PetscInt lits; |
| PetscErrorCode ierr; |
| EPS eps = (EPS)data; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX*)eps->data; |
| PetscInt lits; |
| PetscFunctionBegin; |
| ierr = KSPSolve(blopex->ksp,(Vec)x,(Vec)y); CHKERRABORT(PETSC_COMM_WORLD,ierr); |
| ierr = KSPGetIterationNumber(blopex->ksp, &lits); CHKERRABORT(PETSC_COMM_WORLD,ierr); |
| eps->OP->lineariterations+= lits; |
| PetscFunctionReturnVoid(); |
| } |
| #define __FUNCT__ "Precond_FnMultiVector" |
| static void Precond_FnMultiVector(void *data,void *x,void *y) |
| { |
| EPS eps = (EPS)data; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX*)eps->data; |
| EPS eps = (EPS)data; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX*)eps->data; |
| PetscFunctionBegin; |
| blopex->ii.Eval(Precond_FnSingleVector,data,x,y); |
| PetscFunctionReturnVoid(); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSetUp_BLOPEX" |
| PetscErrorCode EPSSetUp_BLOPEX(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX *)eps->data; |
| PetscBool isPrecond, isPreonly; |
| PetscErrorCode ierr; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX *)eps->data; |
| PetscBool isPrecond, isPreonly; |
| PetscFunctionBegin; |
| if (!eps->ishermitian) { |
| #define __FUNCT__ "EPSSolve_BLOPEX" |
| PetscErrorCode EPSSolve_BLOPEX(EPS eps) |
| { |
| EPS_BLOPEX *blopex = (EPS_BLOPEX *)eps->data; |
| int info,its; |
| EPS_BLOPEX *blopex = (EPS_BLOPEX *)eps->data; |
| int info,its; |
| PetscFunctionBegin; |
| eps->nconv = eps->ncv; |
| if (info==-1) eps->reason = EPS_DIVERGED_ITS; |
| else eps->reason = EPS_CONVERGED_TOL; |
| PetscFunctionReturn(0); |
| } |
| const char* prefix; |
| PetscFunctionBegin; |
| ierr = STSetType(eps->OP, STPRECOND); CHKERRQ(ierr); |
| ierr = STPrecondSetKSPHasMat(eps->OP, PETSC_TRUE); CHKERRQ(ierr); |
| static void multMatvec_PRIMME(void *in, void *out, int *blockSize, primme_params *primme); |
| static void applyPreconditioner_PRIMME(void *in, void *out, int *blockSize, struct primme_params *primme); |
| void par_GlobalSumDouble(void *sendBuf, void *recvBuf, int *count, primme_params *primme) { |
| void par_GlobalSumDouble(void *sendBuf, void *recvBuf, int *count, primme_params *primme) |
| { |
| PetscErrorCode ierr; |
| ierr = MPI_Allreduce((double*)sendBuf, (double*)recvBuf, *count, MPI_DOUBLE, MPI_SUM, ((PetscObject)(primme->commInfo))->comm);CHKERRABORT(((PetscObject)(primme->commInfo))->comm,ierr); |
| } |
| PetscBool t; |
| PetscFunctionBegin; |
| ierr = MPI_Comm_size(((PetscObject)eps)->comm,&numProcs);CHKERRQ(ierr); |
| ierr = MPI_Comm_rank(((PetscObject)eps)->comm,&procID);CHKERRQ(ierr); |
| #endif |
| PetscFunctionBegin; |
| /* Reset some parameters left from previous runs */ |
| ops->primme.aNorm = 0.0; |
| ops->primme.initSize = eps->nini; |
| eps->reason = eps->ncv >= eps->nev ? EPS_CONVERGED_TOL : EPS_DIVERGED_ITS; |
| eps->its = ops->primme.stats.numOuterIterations; |
| eps->OP->applys = ops->primme.stats.numMatvecs; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscInt i, N = primme->n; |
| EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix; |
| Vec x = ops->x; |
| Vec y = ops->y; |
| Vec x = ops->x, y = ops->y; |
| Mat A = ops->A; |
| PetscFunctionBegin; |
| PetscErrorCode ierr; |
| PetscInt i, N = primme->n, lits; |
| EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix; |
| Vec x = ops->x; |
| Vec y = ops->y; |
| Vec x = ops->x, y = ops->y; |
| PetscFunctionBegin; |
| for (i=0;i<*blockSize;i++) { |
| PetscErrorCode EPSDestroy_PRIMME(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_PRIMME *ops = (EPS_PRIMME *)eps->data; |
| EPS_PRIMME *ops = (EPS_PRIMME *)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| primme_Free(&ops->primme); |
| ierr = VecDestroy(&ops->x);CHKERRQ(ierr); |
| ierr = VecDestroy(&ops->y);CHKERRQ(ierr); |
| ierr = PetscFree(eps->data);CHKERRQ(ierr); |
| ierr = EPSFreeSolution(eps);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetBlockSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetMethod_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetBlockSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetMethod_C","",PETSC_NULL);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| /* Display PRIMME params */ |
| primme_display_params(*primme); |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode EPSSetFromOptions_PRIMME(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_PRIMME *ops = (EPS_PRIMME *)eps->data; |
| EPS_PRIMME *ops = (EPS_PRIMME *)eps->data; |
| PetscInt op; |
| PetscBool flg; |
| PetscFunctionBegin; |
| ierr = PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"PRIMME Options","EPS");CHKERRQ(ierr); |
| op = ops->primme.maxBlockSize; |
| op = 0; |
| ierr = PetscOptionsEList("-eps_primme_method","set method for solving the eigenproblem", |
| "EPSPRIMMESetMethod",methodList,15,methodList[1],&op,&flg); CHKERRQ(ierr); |
| if (flg) {ierr = EPSPRIMMESetMethod(eps, methodN[op]);CHKERRQ(ierr);} |
| if (flg) { ierr = EPSPRIMMESetMethod(eps, methodN[op]);CHKERRQ(ierr); } |
| ierr = PetscOptionsEnd();CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EPS_PRIMME *ops = (EPS_PRIMME *) eps->data; |
| PetscFunctionBegin; |
| if (bs == PETSC_DEFAULT) ops->primme.maxBlockSize = 1; |
| else if (bs <= 0) { |
| SETERRQ(((PetscObject)eps)->comm,1, "PRIMME: wrong block size"); |
| } else ops->primme.maxBlockSize = bs; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSPRIMMESetBlockSize" |
| /*@ |
| EPSPRIMMESetBlockSize - The maximum block size the code will try to use. |
| The user should set |
| this based on the architecture specifics of the target computer, |
| as well as any a priori knowledge of multiplicities. The code does |
| NOT require BlockSize > 1 to find multiple eigenvalues. For some |
| methods, keeping BlockSize = 1 yields the best overall performance. |
| EPSPRIMMESetBlockSize - The maximum block size the code will try to use. |
| The user should set |
| this based on the architecture specifics of the target computer, |
| as well as any a priori knowledge of multiplicities. The code does |
| NOT require BlockSize > 1 to find multiple eigenvalues. For some |
| methods, keeping BlockSize = 1 yields the best overall performance. |
| Collective on EPS |
| EPS_PRIMME *ops = (EPS_PRIMME *) eps->data; |
| PetscFunctionBegin; |
| if (bs) *bs = ops->primme.maxBlockSize; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSPRIMMEGetBlockSize" |
| /*@ |
| EPSPRIMMEGetBlockSize - Get the maximum block size the code will try to use. |
| Collective on EPS |
| EPSPRIMMEGetBlockSize - Get the maximum block size the code will try to use. |
| Collective on EPS |
| Input Parameters: |
| . eps - the eigenproblem solver context |
| EPS_PRIMME *ops = (EPS_PRIMME *) eps->data; |
| PetscFunctionBegin; |
| if (method == PETSC_DEFAULT) ops->method = DEFAULT_MIN_TIME; |
| else ops->method = (primme_preset_method)method; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| EPS_PRIMME *ops = (EPS_PRIMME *) eps->data; |
| PetscFunctionBegin; |
| if (method) |
| *method = (EPSPRIMMEMethod)ops->method; |
| if (method) *method = (EPSPRIMMEMethod)ops->method; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| EPS_PRIMME *primme; |
| PetscFunctionBegin; |
| ierr = STSetType(eps->OP, STPRECOND); CHKERRQ(ierr); |
| ierr = STPrecondSetKSPHasMat(eps->OP, PETSC_TRUE); CHKERRQ(ierr); |
| if( eps->nconv == eps->nev ) eps->reason = EPS_CONVERGED_TOL; |
| else eps->reason = EPS_DIVERGED_ITS; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode EPSSolve_DSITRLANCZOS(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscInt i,k,l,lds,lt,nv,m; |
| PetscInt i,k,l,lds,lt,nv,m,*iwork; |
| Vec u=eps->work[0]; |
| PetscScalar *Q, sigma, lambda, zero = 0.0; |
| PetscReal *a,*b,*work,beta,distance = 1e-3; |
| PetscInt *iwork; |
| PetscBool breakdown; |
| PetscFunctionBegin; |
| ierr = EPSMonitor(eps,eps->its,k,eps->eigr,eps->eigi,eps->errest,nv+eps->nconv);CHKERRQ(ierr); |
| eps->nconv = k; |
| } |
| ierr = PetscFree(Q);CHKERRQ(ierr); |
| ierr = PetscFree(a);CHKERRQ(ierr); |
| ierr = PetscFree(b);CHKERRQ(ierr); |
| PetscErrorCode EPSCreate_DSITRLANCZOS(EPS eps) |
| { |
| PetscFunctionBegin; |
| eps->data = PETSC_NULL; |
| eps->ops->setup = EPSSetUp_DSITRLANCZOS; |
| eps->ops->setfromoptions = PETSC_NULL; |
| eps->ops->destroy = EPSDestroy_Default; |
| eps->ops->view = PETSC_NULL; |
| eps->ops->computevectors = EPSComputeVectors_Default; |
| eps->data = PETSC_NULL; |
| eps->ops->setup = EPSSetUp_DSITRLANCZOS; |
| eps->ops->setfromoptions = PETSC_NULL; |
| eps->ops->destroy = EPSDestroy_Default; |
| eps->ops->view = PETSC_NULL; |
| eps->ops->computevectors = EPSComputeVectors_Default; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscReal norm; |
| PetscFunctionBegin; |
| for (j=k;j<m-1;j++) { |
| if (trans) { ierr = STApplyTranspose(eps->OP,V[j],V[j+1]);CHKERRQ(ierr); } |
| else { ierr = STApply(eps->OP,V[j],V[j+1]);CHKERRQ(ierr); } |
| if (trans) { ierr = STApplyTranspose(eps->OP,V[m-1],f);CHKERRQ(ierr); } |
| else { ierr = STApply(eps->OP,V[m-1],f);CHKERRQ(ierr); } |
| ierr = IPOrthogonalize(eps->ip,eps->nds,eps->DS,m,PETSC_NULL,V,f,H+ldh*(m-1),beta,PETSC_NULL);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| } |
| if (marker!=-1) k = marker; |
| *kout = k; |
| PetscFunctionReturn(0); |
| } |
| PetscReal norm=0.0; |
| PetscFunctionBegin; |
| for (j=k;j<m;j++) { |
| ierr = STApply(eps->OP,V[j],f);CHKERRQ(ierr); |
| ierr = IPOrthogonalize(eps->ip,0,PETSC_NULL,eps->nds,PETSC_NULL,eps->DS,f,PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); |
| ierr = IPNorm(eps->ip,f,beta);CHKERRQ(ierr); |
| ierr = VecScale(f,1.0 / *beta);CHKERRQ(ierr); |
| *breakdown = PETSC_FALSE; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSArnoldiSetDelayed_ARNOLDI" |
| PetscErrorCode EPSArnoldiSetDelayed_ARNOLDI(EPS eps,PetscBool delayed) |
| { |
| EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data; |
| EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data; |
| PetscFunctionBegin; |
| arnoldi->delayed = delayed; |
| #define __FUNCT__ "EPSArnoldiGetDelayed_ARNOLDI" |
| PetscErrorCode EPSArnoldiGetDelayed_ARNOLDI(EPS eps,PetscBool *delayed) |
| { |
| EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data; |
| EPS_ARNOLDI *arnoldi = (EPS_ARNOLDI *)eps->data; |
| PetscFunctionBegin; |
| *delayed = arnoldi->delayed; |
| PetscInt i,j,m = *M; |
| PetscReal norm,*omega,lomega[100],*omega_old,lomega_old[100],eps1,delta,eta; |
| PetscBool *which,lwhich[100],*which2,lwhich2[100], |
| reorth = PETSC_FALSE,force_reorth = PETSC_FALSE,fro = PETSC_FALSE,estimate_anorm = PETSC_FALSE; |
| reorth = PETSC_FALSE,force_reorth = PETSC_FALSE, |
| fro = PETSC_FALSE,estimate_anorm = PETSC_FALSE; |
| PetscScalar *hwork,lhwork[100]; |
| PetscFunctionBegin; |
| */ |
| static PetscErrorCode EPSBasicLanczos(EPS eps,PetscReal *alpha,PetscReal *beta,Vec *V,PetscInt k,PetscInt *m,Vec f,PetscBool *breakdown,PetscReal anorm) |
| { |
| EPS_LANCZOS *lanczos = (EPS_LANCZOS *)eps->data; |
| EPS_LANCZOS *lanczos = (EPS_LANCZOS *)eps->data; |
| PetscScalar *T; |
| PetscInt i,n=*m; |
| PetscReal betam; |
| PetscErrorCode ierr; |
| IPOrthogonalizationRefinementType orthog_ref; |
| PetscScalar *T; |
| PetscInt i,n=*m; |
| PetscReal betam; |
| PetscFunctionBegin; |
| switch (lanczos->reorthog) { |
| #define __FUNCT__ "EPSSolve_LANCZOS" |
| PetscErrorCode EPSSolve_LANCZOS(EPS eps) |
| { |
| EPS_LANCZOS *lanczos = (EPS_LANCZOS *)eps->data; |
| EPS_LANCZOS *lanczos = (EPS_LANCZOS *)eps->data; |
| PetscErrorCode ierr; |
| PetscInt nconv,i,j,k,l,x,n,m,*perm,restart,ncv=eps->ncv,r; |
| Vec w=eps->work[1],f=eps->work[0]; |
| PetscErrorCode EPSLanczosGetReorthog_LANCZOS(EPS eps,EPSLanczosReorthogType *reorthog) |
| { |
| EPS_LANCZOS *lanczos = (EPS_LANCZOS *)eps->data; |
| PetscFunctionBegin; |
| *reorthog = lanczos->reorthog; |
| PetscFunctionReturn(0); |
| PetscFunctionBegin; |
| ierr = PetscLogEventBegin(EPS_Dense,0,0,0,0);CHKERRQ(ierr); |
| n = PetscBLASIntCast(n_); |
| /* quick return */ |
| if (n == 1) { |
| PetscErrorCode EPSSolve_KRYLOVSCHUR_SYMM(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscInt i,k,l,lds,lt,nv,m; |
| PetscInt i,k,l,lds,lt,nv,m,*iwork; |
| Vec u=eps->work[0]; |
| PetscScalar *Q; |
| PetscReal *a,*b,*work,beta; |
| PetscInt *iwork; |
| PetscBool breakdown; |
| PetscFunctionBegin; |
| ierr = EPSMonitor(eps,eps->its,k,eps->eigr,eps->eigi,eps->errest,nv+eps->nconv);CHKERRQ(ierr); |
| eps->nconv = k; |
| } |
| ierr = PetscFree(Q);CHKERRQ(ierr); |
| ierr = PetscFree(a);CHKERRQ(ierr); |
| ierr = PetscFree(b);CHKERRQ(ierr); |
| /* S = S + g*b' */ |
| for (i=0;i<m;i++) |
| S[i+(m-1)*lds] = S[i+(m-1)*lds] + g[i]*beta; |
| PetscFunctionReturn(0); |
| #endif |
| } |
| ierr = VecCopy(u,eps->V[k+l]);CHKERRQ(ierr); |
| } |
| eps->nconv = k; |
| ierr = EPSMonitor(eps,eps->its,eps->nconv,eps->eigr,eps->eigi,eps->errest,nv);CHKERRQ(ierr); |
| } |
| ierr = PetscFree(Q);CHKERRQ(ierr); |
| eps->nconv = k; |
| ierr = EPSMonitor(eps,eps->its,eps->nconv,eps->eigr,eps->eigi,eps->errest,nv);CHKERRQ(ierr); |
| } |
| ierr = PetscFree(Q);CHKERRQ(ierr); |
| ierr = PetscFree(work);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| PetscErrorCode EPSCreate_KRYLOVSCHUR(EPS eps) |
| { |
| PetscFunctionBegin; |
| eps->data = PETSC_NULL; |
| eps->ops->setup = EPSSetUp_KRYLOVSCHUR; |
| eps->ops->setfromoptions = PETSC_NULL; |
| eps->ops->destroy = EPSDestroy_Default; |
| eps->ops->view = PETSC_NULL; |
| eps->ops->backtransform = EPSBackTransform_Default; |
| eps->ops->computevectors = EPSComputeVectors_Schur; |
| eps->data = PETSC_NULL; |
| eps->ops->setup = EPSSetUp_KRYLOVSCHUR; |
| eps->ops->setfromoptions = PETSC_NULL; |
| eps->ops->destroy = EPSDestroy_Default; |
| eps->ops->view = PETSC_NULL; |
| eps->ops->backtransform = EPSBackTransform_Default; |
| eps->ops->computevectors = EPSComputeVectors_Schur; |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSCreate_DAVIDSON" |
| PetscErrorCode EPSCreate_DAVIDSON(EPS eps) { |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data; |
| PetscErrorCode EPSCreate_DAVIDSON(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data; |
| PetscFunctionBegin; |
| ierr = STSetType(eps->OP, STPRECOND); CHKERRQ(ierr); |
| ierr = STPrecondSetKSPHasMat(eps->OP, PETSC_FALSE); CHKERRQ(ierr); |
| ierr = EPSDAVIDSONSetRestart_DAVIDSON(eps, 6, 0); CHKERRQ(ierr); |
| ierr = EPSDAVIDSONSetInitialSize_DAVIDSON(eps, 5); CHKERRQ(ierr); |
| ierr = EPSDAVIDSONSetFix_DAVIDSON(eps, 0.01); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSetUp_DAVIDSON" |
| PetscErrorCode EPSSetUp_DAVIDSON(EPS eps) { |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *dvd = &data->ddb; |
| dvdBlackboard b; |
| PetscInt i,nvecs,nscalars,min_size_V,plusk,bs,initv; |
| Mat A,B; |
| KSP ksp; |
| PetscBool t,ipB,ispositive; |
| HarmType_t harm; |
| InitType_t init; |
| PetscReal fix; |
| PetscErrorCode EPSSetUp_DAVIDSON(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *dvd = &data->ddb; |
| dvdBlackboard b; |
| PetscInt i,nvecs,nscalars,min_size_V,plusk,bs,initv; |
| Mat A,B; |
| KSP ksp; |
| PetscBool t,ipB,ispositive; |
| HarmType_t harm; |
| InitType_t init; |
| PetscReal fix; |
| PetscFunctionBegin; |
| /* Setup EPS options and get the problem specification */ |
| ierr = EPSDAVIDSONGetBlockSize_DAVIDSON(eps, &bs); CHKERRQ(ierr); |
| if (bs <= 0) bs = 1; |
| initv, |
| PetscAbs(eps->nini), |
| plusk, harm, |
| PETSC_NULL, init, eps->trackall); |
| CHKERRQ(ierr); |
| PETSC_NULL, init, eps->trackall);CHKERRQ(ierr); |
| /* Reserve memory */ |
| nvecs = b.max_size_auxV + b.own_vecs; |
| data->size_wV = nvecs; |
| for (i=0; i<nvecs; i++) { |
| ierr = VecCreateMPIWithArray(((PetscObject)eps)->comm, eps->nloc, PETSC_DECIDE, |
| data->wS+i*eps->nloc, &data->wV[i]); |
| CHKERRQ(ierr); |
| data->wS+i*eps->nloc, &data->wV[i]);CHKERRQ(ierr); |
| } |
| b.free_vecs = data->wV; |
| b.free_scalars = data->wS + nvecs*eps->nloc; |
| PetscAbs(eps->nini), plusk, |
| eps->ip, harm, dvd->withTarget, |
| eps->target, ksp, |
| fix, init, eps->trackall); |
| CHKERRQ(ierr); |
| fix, init, eps->trackall);CHKERRQ(ierr); |
| /* Associate the eigenvalues to the EPS */ |
| eps->eigr = dvd->eigr; |
| eps->eigi = dvd->eigi; |
| eps->errest = dvd->errest; |
| eps->V = dvd->V; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSolve_DAVIDSON" |
| PetscErrorCode EPSSolve_DAVIDSON(EPS eps) { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *d = &data->ddb; |
| PetscErrorCode ierr; |
| PetscErrorCode EPSSolve_DAVIDSON(EPS eps) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *d = &data->ddb; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* Call the starting routines */ |
| DVD_FL_CALL(d->startList, d); |
| if( eps->nconv >= eps->nev ) eps->reason = EPS_CONVERGED_TOL; |
| else eps->reason = EPS_DIVERGED_ITS; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSDestroy_DAVIDSON" |
| PetscErrorCode EPSDestroy_DAVIDSON(EPS eps) { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *dvd = &data->ddb; |
| PetscErrorCode ierr; |
| PetscInt i; |
| PetscErrorCode EPSDestroy_DAVIDSON(EPS eps) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *dvd = &data->ddb; |
| PetscErrorCode ierr; |
| PetscInt i; |
| PetscFunctionBegin; |
| /* Call step destructors and destroys the list */ |
| DVD_FL_CALL(dvd->destroyList, dvd); |
| DVD_FL_DEL(dvd->destroyList); |
| ierr = PetscFree(data->wV);CHKERRQ(ierr); |
| ierr = PetscFree(data->wS);CHKERRQ(ierr); |
| ierr = PetscFree(data);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSView_DAVIDSON" |
| PetscErrorCode EPSView_DAVIDSON(EPS eps,PetscViewer viewer) |
| { |
| PetscErrorCode ierr; |
| PetscBool isascii; |
| PetscInt opi, opi0; |
| PetscBool opb; |
| const char* name; |
| PetscErrorCode ierr; |
| PetscBool isascii, opb; |
| PetscInt opi, opi0; |
| const char* name; |
| PetscFunctionBegin; |
| name = ((PetscObject)eps)->type_name; |
| ierr = PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); |
| if (!isascii) { |
| ierr = EPSDAVIDSONGetRestart_DAVIDSON(eps, &opi, &opi0); CHKERRQ(ierr); |
| ierr = PetscViewerASCIIPrintf(viewer,"size of the subspace after restarting: %d\n", opi);CHKERRQ(ierr); |
| ierr = PetscViewerASCIIPrintf(viewer,"number of vectors after restarting from the previous iteration: %d\n", opi0);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "SLEPcNotImplemented" |
| PetscErrorCode SLEPcNotImplemented() { |
| PetscErrorCode SLEPcNotImplemented() |
| { |
| SETERRQ(PETSC_COMM_WORLD,1, "Do not call this function!"); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSDAVIDSONSetKrylovStart_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONSetKrylovStart_DAVIDSON(EPS eps,PetscBool krylovstart) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| data->krylovstart = krylovstart; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONGetKrylovStart_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONGetKrylovStart_DAVIDSON(EPS eps,PetscBool *krylovstart) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *krylovstart = data->krylovstart; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONSetBlockSize_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONSetBlockSize_DAVIDSON(EPS eps,PetscInt blocksize) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if(blocksize == PETSC_DEFAULT || blocksize == PETSC_DECIDE) blocksize = 1; |
| if(blocksize <= 0) |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid blocksize value"); |
| data->blocksize = blocksize; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONGetBlockSize_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONGetBlockSize_DAVIDSON(EPS eps,PetscInt *blocksize) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *blocksize = data->blocksize; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONSetRestart_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONSetRestart_DAVIDSON(EPS eps,PetscInt minv,PetscInt plusk) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if(minv == PETSC_DEFAULT || minv == PETSC_DECIDE) minv = 5; |
| if(minv <= 0) |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid minv value"); |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid plusk value"); |
| data->minv = minv; |
| data->plusk = plusk; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONGetRestart_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONGetRestart_DAVIDSON(EPS eps,PetscInt *minv,PetscInt *plusk) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *minv = data->minv; |
| *plusk = data->plusk; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONGetInitialSize_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONGetInitialSize_DAVIDSON(EPS eps,PetscInt *initialsize) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *initialsize = data->initialsize; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONSetInitialSize_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONSetInitialSize_DAVIDSON(EPS eps,PetscInt initialsize) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if(initialsize == PETSC_DEFAULT || initialsize == PETSC_DECIDE) initialsize = 5; |
| if(initialsize <= 0) |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid initial size value"); |
| data->initialsize = initialsize; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONGetFix_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONGetFix_DAVIDSON(EPS eps,PetscReal *fix) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| *fix = data->fix; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSDAVIDSONSetFix_DAVIDSON" |
| PetscErrorCode EPSDAVIDSONSetFix_DAVIDSON(EPS eps,PetscReal fix) |
| { |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,1); |
| if(fix == PETSC_DEFAULT || fix == PETSC_DECIDE) fix = 0.01; |
| if(fix < 0.0) |
| SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid fix value"); |
| data->fix = fix; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSComputeVectors_QZ" |
| /* |
| */ |
| PetscErrorCode EPSComputeVectors_QZ(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *d = &data->ddb; |
| PetscScalar *pX, *auxS; |
| PetscInt size_auxS; |
| PetscErrorCode ierr; |
| EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data; |
| dvdDashboard *d = &data->ddb; |
| PetscScalar *pX, *auxS; |
| PetscInt size_auxS; |
| PetscFunctionBegin; |
| /* Compute the eigenvectors associated to (cS, cT) */ |
| ierr = PetscMalloc(sizeof(PetscScalar)*d->nconv*d->nconv, &pX); CHKERRQ(ierr); |
| size_auxS = 11*d->nconv + 4*d->nconv*d->nconv; |
| size_auxS, PETSC_FALSE); CHKERRQ(ierr); |
| /* pX[i] <- pX[i] / ||pX[i]|| */ |
| ierr = SlepcDenseNorm(pX, d->nconv, d->nconv, d->nconv, d->ceigi); |
| CHKERRQ(ierr); |
| ierr = SlepcDenseNorm(pX, d->nconv, d->nconv, d->nconv, d->ceigi);CHKERRQ(ierr); |
| /* V <- cX * pX */ |
| ierr = SlepcUpdateVectorsZ(eps->V, 0.0, 1.0, d->cX, d->size_cX, pX, |
| ierr = PetscFree(auxS); CHKERRQ(ierr); |
| eps->evecsavailable = PETSC_TRUE; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSSetFromOptions_GD" |
| PetscErrorCode EPSSetFromOptions_GD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscBool flg,op; |
| PetscInt opi,opi0; |
| PetscErrorCode ierr; |
| PetscBool flg,op; |
| PetscInt opi,opi0; |
| PetscFunctionBegin; |
| ierr = PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"GD Options","EPS");CHKERRQ(ierr); |
| ierr = EPSGDGetKrylovStart(eps, &op); CHKERRQ(ierr); |
| if(flg) { ierr = EPSGDSetInitialSize(eps, opi); CHKERRQ(ierr); } |
| ierr = PetscOptionsEnd();CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSetUp_GD" |
| PetscErrorCode EPSSetUp_GD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscBool t; |
| KSP ksp; |
| PetscErrorCode ierr; |
| PetscBool t; |
| KSP ksp; |
| PetscFunctionBegin; |
| /* Setup common for all davidson solvers */ |
| ierr = EPSSetUp_DAVIDSON(eps); CHKERRQ(ierr); |
| ierr = STGetKSP(eps->OP, &ksp); CHKERRQ(ierr); |
| ierr = PetscTypeCompare((PetscObject)ksp, KSPPREONLY, &t); CHKERRQ(ierr); |
| if (!t) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP, "EPSGD only works with KSPPREONLY"); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_BEGIN |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSCreate_GD" |
| PetscErrorCode EPSCreate_GD(EPS eps) { |
| PetscErrorCode EPSCreate_GD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* Load the DAVIDSON solver */ |
| ierr = EPSCreate_DAVIDSON(eps); CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDGetRestart_C","EPSDAVIDSONGetRestart_DAVIDSON",EPSDAVIDSONGetRestart_DAVIDSON);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDSetInitialSize_C","EPSDAVIDSONSetInitialSize_DAVIDSON",EPSDAVIDSONSetInitialSize_DAVIDSON);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDGetInitialSize_C","EPSDAVIDSONGetInitialSize_DAVIDSON",EPSDAVIDSONGetInitialSize_DAVIDSON);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDSetKrylovStart_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDGetKrylovStart_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDSetBlockSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDGetRestart_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDSetInitialSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSGDGetInitialSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = EPSDestroy_DAVIDSON(eps); |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "EPSSetFromOptions_JD" |
| PetscErrorCode EPSSetFromOptions_JD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscBool flg,op; |
| PetscInt opi,opi0; |
| PetscReal opf; |
| PetscErrorCode ierr; |
| PetscBool flg,op; |
| PetscInt opi,opi0; |
| PetscReal opf; |
| PetscFunctionBegin; |
| ierr = PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"JD Options","EPS");CHKERRQ(ierr); |
| ierr = EPSJDGetKrylovStart(eps, &op); CHKERRQ(ierr); |
| if(flg) { ierr = EPSJDSetFix(eps, opf); CHKERRQ(ierr); } |
| ierr = PetscOptionsEnd();CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSSetUp_JD" |
| PetscErrorCode EPSSetUp_JD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| PetscBool t; |
| KSP ksp; |
| PetscErrorCode ierr; |
| PetscBool t; |
| KSP ksp; |
| PetscFunctionBegin; |
| /* Setup common for all davidson solvers */ |
| ierr = EPSSetUp_DAVIDSON(eps); CHKERRQ(ierr); |
| ierr = STGetKSP(eps->OP, &ksp); CHKERRQ(ierr); |
| ierr = PetscTypeCompare((PetscObject)ksp, KSPPREONLY, &t); CHKERRQ(ierr); |
| if (t) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP, "EPSJD does not work with KSPPREONLY"); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_BEGIN |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSCreate_JD" |
| PetscErrorCode EPSCreate_JD(EPS eps) { |
| PetscErrorCode ierr; |
| KSP ksp; |
| PetscErrorCode EPSCreate_JD(EPS eps) |
| { |
| PetscErrorCode ierr; |
| KSP ksp; |
| PetscFunctionBegin; |
| /* Load the DAVIDSON solver */ |
| ierr = EPSCreate_DAVIDSON(eps); CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetInitialSize_C","EPSDAVIDSONGetInitialSize_DAVIDSON",EPSDAVIDSONGetInitialSize_DAVIDSON);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetFix_C","EPSDAVIDSONSetFix_DAVIDSON",EPSDAVIDSONSetFix_DAVIDSON);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetFix_C","EPSDAVIDSONGetFix_DAVIDSON",EPSDAVIDSONGetFix_DAVIDSON);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetKrylovStart_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetKrylovStart_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetBlockSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetInitialSize_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetFix_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetFix_C","",PETSC_NULL);CHKERRQ(ierr); |
| ierr = EPSDestroy_DAVIDSON(eps); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "EPSJDSetKrylovStart" |
| /*@ |
| ierr = PetscFree(pW);CHKERRQ(ierr); |
| } |
| } |
| eps->nconv = eps->ncv; |
| eps->its = 1; |
| eps->reason = EPS_CONVERGED_TOL; |
| PetscFunctionReturn(0); |
| } |
| } |
| } |
| } |
| if (eps->its >= eps->max_it) eps->reason = EPS_DIVERGED_ITS; |
| } |
| ierr = PetscFree(select);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| ierr = EPSGetStartVectorLeft(eps,eps->nconv,w,PETSC_NULL);CHKERRQ(ierr); |
| } |
| } |
| if( eps->nconv == eps->nev ) eps->reason = EPS_CONVERGED_TOL; |
| else eps->reason = EPS_DIVERGED_ITS; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode EPSBackTransform_POWER(EPS eps) |
| { |
| PetscErrorCode ierr; |
| EPS_POWER *power = (EPS_POWER *)eps->data; |
| EPS_POWER *power = (EPS_POWER *)eps->data; |
| PetscFunctionBegin; |
| if (power->shift_type == EPS_POWER_SHIFT_CONSTANT) { |
| PetscErrorCode EPSPowerGetShiftType_POWER(EPS eps,EPSPowerShiftType *shift) |
| { |
| EPS_POWER *power = (EPS_POWER *)eps->data; |
| PetscFunctionBegin; |
| *shift = power->shift_type; |
| PetscFunctionReturn(0); |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidPointer(norm,3); |
| if (!ip->matrix && ip->bilinear_form == IP_INNER_HERMITIAN) { |
| ierr = VecNorm(x,NORM_2,norm);CHKERRQ(ierr); |
| } else { |
| *norm = PetscSqrtScalar(p); |
| #endif |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidPointer(norm,3); |
| if (!ip->matrix && ip->bilinear_form == IP_INNER_HERMITIAN) { |
| ierr = VecNormBegin(x,NORM_2,norm);CHKERRQ(ierr); |
| } else { |
| ierr = IPInnerProductBegin(ip,x,x,&p);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidPointer(norm,3); |
| if (!ip->matrix && ip->bilinear_form == IP_INNER_HERMITIAN) { |
| ierr = VecNormEnd(x,NORM_2,norm);CHKERRQ(ierr); |
| } else { |
| *norm = PetscSqrtScalar(p); |
| #endif |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidHeaderSpecific(y,VEC_CLASSID,3); |
| PetscValidScalarPointer(p,4); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| ip->innerproducts++; |
| if (ip->matrix) { |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidHeaderSpecific(y,VEC_CLASSID,3); |
| PetscValidScalarPointer(p,4); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| ip->innerproducts++; |
| if (ip->matrix) { |
| PetscValidHeaderSpecific(x,VEC_CLASSID,2); |
| PetscValidHeaderSpecific(y,VEC_CLASSID,3); |
| PetscValidScalarPointer(p,4); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| if (ip->matrix) { |
| if (ip->bilinear_form == IP_INNER_HERMITIAN) { |
| PetscValidPointer(y,4); |
| PetscValidHeaderSpecific(*y,VEC_CLASSID,4); |
| PetscValidScalarPointer(p,5); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| ip->innerproducts += n; |
| if (ip->matrix) { |
| PetscValidPointer(y,4); |
| PetscValidHeaderSpecific(*y,VEC_CLASSID,4); |
| PetscValidScalarPointer(p,5); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| ip->innerproducts += n; |
| if (ip->matrix) { |
| PetscValidPointer(y,4); |
| PetscValidHeaderSpecific(*y,VEC_CLASSID,4); |
| PetscValidScalarPointer(p,5); |
| ierr = PetscLogEventBegin(IP_InnerProduct,ip,x,0,0);CHKERRQ(ierr); |
| if (ip->matrix) { |
| if (ip->bilinear_form == IP_INNER_HERMITIAN) { |
| #include <private/ipimpl.h> /*I "slepcip.h" I*/ |
| PetscClassId IP_CLASSID = 0; |
| PetscClassId IP_CLASSID = 0; |
| PetscLogEvent IP_InnerProduct = 0, IP_Orthogonalize = 0, IP_ApplyMatrix = 0; |
| #undef __FUNCT__ |
| PetscErrorCode IPInitializePackage(const char *path) |
| { |
| static PetscBool initialized = PETSC_FALSE; |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscErrorCode ierr; |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| if (initialized) PetscFunctionReturn(0); |
| @*/ |
| PetscErrorCode IPCreate(MPI_Comm comm,IP *newip) |
| { |
| IP ip; |
| PetscErrorCode ierr; |
| IP ip; |
| PetscFunctionBegin; |
| PetscValidPointer(newip,2); |
| ip->Bx = PETSC_NULL; |
| ip->xid = 0; |
| ip->xstate = 0; |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode IPSetOptionsPrefix(IP ip,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| ierr = PetscObjectSetOptionsPrefix((PetscObject)ip,prefix);CHKERRQ(ierr); |
| PetscErrorCode IPAppendOptionsPrefix(IP ip,const char *prefix) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| ierr = PetscObjectAppendOptionsPrefix((PetscObject)ip,prefix);CHKERRQ(ierr); |
| @*/ |
| PetscErrorCode IPGetOptionsPrefix(IP ip,const char *prefix[]) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| PetscValidPointer(prefix,2); |
| ierr = PetscObjectGetOptionsPrefix((PetscObject)ip, prefix);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| PetscValidPointer(prefix,2); |
| ierr = PetscObjectGetOptionsPrefix((PetscObject)ip, prefix);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| @*/ |
| PetscErrorCode IPSetFromOptions(IP ip) |
| { |
| PetscErrorCode ierr; |
| const char *orth_list[2] = { "mgs" , "cgs" }; |
| const char *ref_list[3] = { "never" , "ifneeded", "always" }; |
| PetscReal r; |
| PetscInt i,j; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| @*/ |
| PetscErrorCode IPView(IP ip,PetscViewer viewer) |
| { |
| PetscBool isascii; |
| PetscErrorCode ierr; |
| PetscBool isascii; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)ip)->comm); |
| PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2); |
| PetscCheckSameComm(ip,1,viewer,2); |
| ierr = PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr); |
| if (isascii) { |
| ierr = PetscObjectPrintClassNamePrefixType((PetscObject)ip,viewer,"IP Object");CHKERRQ(ierr); |
| PetscReal onrm,nrm; |
| PetscFunctionBegin; |
| if (H) { |
| for (i=0;i<n;i++) |
| H[i] = 0; |
| default: |
| SETERRQ(((PetscObject)ip)->comm,PETSC_ERR_ARG_WRONG,"Unknown orthogonalization refinement"); |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscFunctionBegin; |
| ierr = PetscLogEventBegin(IP_Orthogonalize,ip,0,0,0);CHKERRQ(ierr); |
| if (nds==0 && n==0) { |
| if (norm) { ierr = IPNorm(ip,v,norm);CHKERRQ(ierr); } |
| if (lindep) *lindep = PETSC_FALSE; |
| SETERRQ(((PetscObject)ip)->comm,PETSC_ERR_ARG_WRONG,"Unknown orthogonalization type"); |
| } |
| } |
| ierr = PetscLogEventEnd(IP_Orthogonalize,ip,0,0,0);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscRandom rctx=PETSC_NULL; |
| PetscFunctionBegin; |
| for (k=m; k<n; k++) { |
| /* orthogonalize v_k with respect to v_0, ..., v_{k-1} */ |
| } |
| ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscScalar shh[100],*lhh,*vw,*tau,one=1.0,*work; |
| PetscFunctionBegin; |
| /* Don't allocate small arrays */ |
| if (n<=100) lhh = shh; |
| else { ierr = PetscMalloc(n*sizeof(PetscScalar),&lhh);CHKERRQ(ierr); } |
| PetscScalar lh[100],*h; |
| PetscBool allocated = PETSC_FALSE; |
| PetscReal lhnrm,*hnrm,lnrm,*nrm; |
| PetscFunctionBegin; |
| if (n==0) { |
| if (norm) { ierr = IPNorm(ip,v,norm);CHKERRQ(ierr); } |
| } else { |
| ierr = PetscLogEventBegin(IP_Orthogonalize,ip,0,0,0);CHKERRQ(ierr); |
| /* allocate H if needed */ |
| if (!H) { |
| if (n<=100) h = lh; |
| } |
| if (allocated) { ierr = PetscFree(h);CHKERRQ(ierr); } |
| ierr = PetscLogEventEnd(IP_Orthogonalize,ip,0,0,0);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(ip,IP_CLASSID,1); |
| if (mat) *mat = ip->matrix; |
| if (mat) *mat = ip->matrix; |
| if (form) *form = ip->bilinear_form; |
| PetscFunctionReturn(0); |
| } |
| #include <private/svdimpl.h> /*I "slepcsvd.h" I*/ |
| PetscFList SVDList = 0; |
| PetscClassId SVD_CLASSID = 0; |
| PetscLogEvent SVD_SetUp = 0, SVD_Solve = 0, SVD_Dense = 0; |
| PetscFList SVDList = 0; |
| PetscClassId SVD_CLASSID = 0; |
| PetscLogEvent SVD_SetUp = 0, SVD_Solve = 0, SVD_Dense = 0; |
| static PetscBool SVDPackageInitialized = PETSC_FALSE; |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDFinalizePackage" |
| /*@C |
| SVDFinalizePackage - This function destroys everything in the Slepc interface to the SVD package. It is |
| called from SlepcFinalize(). |
| SVDFinalizePackage - This function destroys everything in the Slepc interface |
| to the SVD package. It is called from SlepcFinalize(). |
| Level: developer |
| Level: developer |
| .seealso: SlepcFinalize() |
| @*/ |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDInitializePackage" |
| /*@C |
| SVDInitializePackage - This function initializes everything in the SVD package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to SVDCreate() |
| when using static libraries. |
| SVDInitializePackage - This function initializes everything in the SVD package. It is called |
| from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to SVDCreate() |
| when using static libraries. |
| Input Parameter: |
| path - The dynamic library path, or PETSC_NULL |
| Input Parameter: |
| . path - The dynamic library path, or PETSC_NULL |
| Level: developer |
| Level: developer |
| .seealso: SlepcInitialize() |
| @*/ |
| PetscErrorCode SVDInitializePackage(const char *path) |
| { |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscErrorCode ierr; |
| char logList[256]; |
| char *className; |
| PetscBool opt; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| if (SVDPackageInitialized) PetscFunctionReturn(0); |
| PetscFunctionBegin; |
| PetscValidPointer(outsvd,2); |
| ierr = PetscHeaderCreate(svd,_p_SVD,struct _SVDOps,SVD_CLASSID,-1,"SVD",comm,SVDDestroy,SVDView);CHKERRQ(ierr); |
| *outsvd = svd; |
| ierr = IPSetOptionsPrefix(svd->ip,((PetscObject)svd)->prefix); |
| ierr = IPAppendOptionsPrefix(svd->ip,"svd_"); |
| PetscLogObjectParent(svd,svd->ip); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDRegister" |
| /*@C |
| SVDRegister - See SVDRegisterDynamic() |
| SVDRegister - See SVDRegisterDynamic() |
| Level: advanced |
| Level: advanced |
| @*/ |
| PetscErrorCode SVDRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(SVD)) |
| { |
| PetscValidHeaderSpecific(ip,IP_CLASSID,2); |
| PetscCheckSameComm(svd,1,ip,2); |
| ierr = PetscObjectReference((PetscObject)ip);CHKERRQ(ierr); |
| ierr = IPDestroy(&svd->ip); CHKERRQ(ierr); |
| ierr = IPDestroy(&svd->ip);CHKERRQ(ierr); |
| svd->ip = ip; |
| PetscFunctionReturn(0); |
| } |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(svd,SVD_CLASSID,1); |
| if (nsv != PETSC_IGNORE) { |
| if (nsv<1) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0"); |
| svd->nsv = nsv; |
| @*/ |
| PetscErrorCode SVDSetFromOptions(SVD svd) |
| { |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg; |
| const char *mode_list[2] = { "explicit", "implicit" }; |
| PetscInt i,j,k; |
| PetscReal r; |
| PetscErrorCode ierr; |
| char type[256],monfilename[PETSC_MAX_PATH_LEN]; |
| PetscBool flg; |
| const char *mode_list[2] = { "explicit", "implicit" }; |
| PetscInt i,j,k; |
| PetscReal r; |
| PetscViewerASCIIMonitor monviewer; |
| SlepcConvMonitor ctx; |
| SlepcConvMonitor ctx; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(svd,SVD_CLASSID,1); |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDSetTrackAll" |
| /*@ |
| SVDSetTrackAll - Specifies if the solver must compute the residual norm of all |
| approximate singular value or not. |
| SVDSetTrackAll - Specifies if the solver must compute the residual norm of all |
| approximate singular value or not. |
| Collective on SVD |
| Collective on SVD |
| Input Parameters: |
| + svd - the singular value solver context |
| - trackall - whether to compute all residuals or not |
| Input Parameters: |
| + svd - the singular value solver context |
| - trackall - whether to compute all residuals or not |
| Notes: |
| If the user sets trackall=PETSC_TRUE then the solver computes (or estimates) |
| the residual norm for each singular value approximation. Computing the residual is |
| usually an expensive operation and solvers commonly compute only the residual |
| associated to the first unconverged singular value. |
| Notes: |
| If the user sets trackall=PETSC_TRUE then the solver computes (or estimates) |
| the residual norm for each singular value approximation. Computing the residual is |
| usually an expensive operation and solvers commonly compute only the residual |
| associated to the first unconverged singular value. |
| The options '-svd_monitor_all' and '-svd_monitor_draw_all' automatically |
| activate this option. |
| The options '-svd_monitor_all' and '-svd_monitor_draw_all' automatically |
| activate this option. |
| Level: intermediate |
| Level: intermediate |
| .seealso: SVDGetTrackAll() |
| @*/ |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDGetTrackAll" |
| /*@ |
| SVDGetTrackAll - Returns the flag indicating whether all residual norms must |
| be computed or not. |
| SVDGetTrackAll - Returns the flag indicating whether all residual norms must |
| be computed or not. |
| Not Collective |
| Not Collective |
| Input Parameter: |
| . svd - the singular value solver context |
| Input Parameter: |
| . svd - the singular value solver context |
| Output Parameter: |
| . trackall - the returned flag |
| Output Parameter: |
| . trackall - the returned flag |
| Level: intermediate |
| Level: intermediate |
| .seealso: SVDSetTrackAll() |
| @*/ |
| PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[]) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(svd,SVD_CLASSID,1); |
| PetscValidPointer(prefix,2); |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(svd,SVD_CLASSID,1); |
| if (!svd->setupcalled) { ierr = SVDSetUp(svd);CHKERRQ(ierr); } |
| svd->its = 0; |
| svd->matvecs = 0; |
| /* Remove the initial subspace */ |
| svd->nini = 0; |
| PetscFunctionReturn(0); |
| } |
| { |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(svd,SVD_CLASSID,1); |
| *ctx = (svd->monitorcontext[0]); |
| *ctx = (svd->monitorcontext[0]); |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx) |
| { |
| PetscViewer viewer = (PetscViewer) monctx; |
| PetscDraw draw; |
| PetscDrawLG lg; |
| PetscDraw draw,draw1; |
| PetscDrawLG lg,lg1; |
| PetscErrorCode ierr; |
| PetscReal x,y,p; |
| PetscDraw draw1; |
| PetscDrawLG lg1; |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDraw(viewer,1,&draw1);CHKERRQ(ierr); |
| PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx) |
| { |
| PetscViewer viewer = (PetscViewer) monctx; |
| PetscDraw draw; |
| PetscDrawLG lg; |
| PetscDraw draw,draw1; |
| PetscDrawLG lg,lg1; |
| PetscErrorCode ierr; |
| PetscReal *x,*y,p; |
| int n = PetscMin(svd->nsv,255); |
| PetscInt i; |
| PetscDraw draw1; |
| PetscDrawLG lg1; |
| PetscInt i,n = PetscMin(svd->nsv,255); |
| PetscFunctionBegin; |
| if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); } |
| ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDrawLG(viewer,0,&lg);CHKERRQ(ierr); |
| ierr = PetscViewerDrawGetDraw(viewer,1,&draw1);CHKERRQ(ierr); |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDRegisterAll" |
| /*@C |
| SVDRegisterAll - Registers all the singular value solvers in the SVD package. |
| SVDRegisterAll - Registers all the singular value solvers in the SVD package. |
| Not Collective |
| Not Collective |
| Level: advanced |
| Level: advanced |
| .seealso: SVDRegisterDynamic() |
| @*/ |
| PetscScalar qwork,*work; |
| PetscBLASInt n,info,lwork,*iwork,M,N; |
| #if defined(PETSC_USE_COMPLEX) |
| PetscReal *rwork; |
| PetscReal *rwork; |
| #endif |
| PetscFunctionBegin; |
| PetscFunctionReturn(0); |
| #endif |
| } |
| #define __FUNCT__ "SVDSetUp_TRLANCZOS" |
| PetscErrorCode SVDSetUp_TRLANCZOS(SVD svd) |
| { |
| PetscErrorCode ierr; |
| PetscInt i,N,nloc; |
| PetscScalar *pU; |
| PetscErrorCode ierr; |
| PetscInt i,N,nloc; |
| PetscScalar *pU; |
| PetscFunctionBegin; |
| ierr = SVDMatGetSize(svd,PETSC_NULL,&N);CHKERRQ(ierr); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDDestroy_TRLANCZOS" |
| PetscErrorCode SVDDestroy_TRLANCZOS(SVD svd) |
| #undef __FUNCT__ |
| #define __FUNCT__ "ShellMatMult_CYCLIC" |
| PetscErrorCode ShellMatMult_CYCLIC(Mat B,Vec x, Vec y) |
| static PetscErrorCode ShellMatMult_CYCLIC(Mat B,Vec x, Vec y) |
| { |
| PetscErrorCode ierr; |
| SVD svd; |
| #undef __FUNCT__ |
| #define __FUNCT__ "ShellMatGetDiagonal_CYCLIC" |
| PetscErrorCode ShellMatGetDiagonal_CYCLIC(Mat B,Vec diag) |
| static PetscErrorCode ShellMatGetDiagonal_CYCLIC(Mat B,Vec diag) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "SVDSetUp_CYCLIC" |
| PetscErrorCode SVDSetUp_CYCLIC(SVD svd) |
| { |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscInt M,N,m,n,i,nloc,isl; |
| PetscScalar *pU; |
| PetscBool trackall; |
| Vec v; |
| Mat Zm,Zn; |
| PetscScalar *isa,*va; |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscInt M,N,m,n,i,nloc,isl; |
| PetscScalar *pU,*isa,*va; |
| PetscBool trackall; |
| Vec v; |
| Mat Zm,Zn; |
| PetscFunctionBegin; |
| ierr = MatDestroy(&cyclic->mat);CHKERRQ(ierr); |
| ierr = VecDestroy(&cyclic->x1);CHKERRQ(ierr); |
| ierr = VecDestroy(&cyclic->x2);CHKERRQ(ierr); |
| ierr = VecCreateMPIWithArray(((PetscObject)svd)->comm,nloc,PETSC_DECIDE,pU+i*nloc,&svd->U[i]);CHKERRQ(ierr); |
| } |
| } |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode SVDCyclicSetEPS_CYCLIC(SVD svd,EPS eps) |
| { |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(eps,EPS_CLASSID,2); |
| #define __FUNCT__ "SVDView_CYCLIC" |
| PetscErrorCode SVDView_CYCLIC(SVD svd,PetscViewer viewer) |
| { |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscFunctionBegin; |
| if (cyclic->explicitmatrix) { |
| #define __FUNCT__ "SVDDestroy_CYCLIC" |
| PetscErrorCode SVDDestroy_CYCLIC(SVD svd) |
| { |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic = (SVD_CYCLIC *)svd->data; |
| PetscFunctionBegin; |
| ierr = EPSDestroy(&cyclic->eps);CHKERRQ(ierr); |
| #define __FUNCT__ "SVDCreate_CYCLIC" |
| PetscErrorCode SVDCreate_CYCLIC(SVD svd) |
| { |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic; |
| PetscErrorCode ierr; |
| SVD_CYCLIC *cyclic; |
| PetscFunctionBegin; |
| ierr = PetscNew(SVD_CYCLIC,&cyclic);CHKERRQ(ierr); |
| #define __FUNCT__ "SVDSetUp_CROSS" |
| PetscErrorCode SVDSetUp_CROSS(SVD svd) |
| { |
| PetscErrorCode ierr; |
| SVD_CROSS *cross = (SVD_CROSS *)svd->data; |
| PetscInt n,i; |
| PetscBool trackall; |
| PetscErrorCode ierr; |
| SVD_CROSS *cross = (SVD_CROSS *)svd->data; |
| PetscInt n,i; |
| PetscBool trackall; |
| PetscFunctionBegin; |
| ierr = MatDestroy(&cross->mat);CHKERRQ(ierr); |
| #define __FUNCT__ "SVDSetFromOptions_CROSS" |
| PetscErrorCode SVDSetFromOptions_CROSS(SVD svd) |
| { |
| SVD_CROSS *cross = (SVD_CROSS *)svd->data; |
| SVD_CROSS *cross = (SVD_CROSS *)svd->data; |
| PetscFunctionBegin; |
| cross->setfromoptionscalled = PETSC_TRUE; |
| #define __FUNCT__ "SVDSetUp_LANCZOS" |
| PetscErrorCode SVDSetUp_LANCZOS(SVD svd) |
| { |
| PetscErrorCode ierr; |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| PetscInt i,N,nloc; |
| PetscScalar *pU; |
| PetscErrorCode ierr; |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| PetscInt i,N,nloc; |
| PetscScalar *pU; |
| PetscFunctionBegin; |
| ierr = SVDMatGetSize(svd,PETSC_NULL,&N);CHKERRQ(ierr); |
| #define __FUNCT__ "SVDLanczosSetOneSide_LANCZOS" |
| PetscErrorCode SVDLanczosSetOneSide_LANCZOS(SVD svd,PetscBool oneside) |
| { |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| PetscFunctionBegin; |
| if (lanczos->oneside != oneside) { |
| #define __FUNCT__ "SVDLanczosGetOneSide_LANCZOS" |
| PetscErrorCode SVDLanczosGetOneSide_LANCZOS(SVD svd,PetscBool *oneside) |
| { |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data; |
| PetscFunctionBegin; |
| PetscValidPointer(oneside,2); |
| #define __FUNCT__ "SVDSetUp_LAPACK" |
| PetscErrorCode SVDSetUp_LAPACK(SVD svd) |
| { |
| PetscErrorCode ierr; |
| PetscInt N,i,nloc; |
| PetscScalar *pU; |
| PetscErrorCode ierr; |
| PetscInt N,i,nloc; |
| PetscScalar *pU; |
| PetscFunctionBegin; |
| ierr = SVDMatGetSize(svd,PETSC_NULL,&N);CHKERRQ(ierr); |
| #define __FUNCT__ "SVDSolve_LAPACK" |
| PetscErrorCode SVDSolve_LAPACK(SVD svd) |
| { |
| PetscErrorCode ierr; |
| PetscInt M,N,n,i,j,k; |
| Mat mat; |
| PetscScalar *pU,*pVT,*pmat,*pu,*pv; |
| PetscReal *sigma; |
| PetscErrorCode ierr; |
| PetscInt M,N,n,i,j,k; |
| Mat mat; |
| PetscScalar *pU,*pVT,*pmat,*pu,*pv; |
| PetscReal *sigma; |
| PetscFunctionBegin; |
| ierr = MatConvert(svd->OP,MATSEQDENSE,MAT_INITIAL_MATRIX,&mat);CHKERRQ(ierr); |
| #define PetscValidVecComp(y) |
| #endif |
| static PetscErrorCode VecCreate_Comp_Private(Vec v, Vec *x, PetscInt nx, |
| PetscBool x_to_me, Vec_Comp_N* n); |
| static PetscErrorCode VecCreate_Comp_Private(Vec v, Vec *x, PetscInt nx, PetscBool x_to_me, Vec_Comp_N* n); |
| #include "veccomp0.h" |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = VecRegisterDynamic(VECCOMP, path, "VecCreate_Comp", VecCreate_Comp); |
| CHKERRQ(ierr); |
| ierr = VecRegisterDynamic(VECCOMP, path, "VecCreate_Comp", VecCreate_Comp);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| /* if memory was published with AMS then destroy it */ |
| ierr = PetscObjectDepublish(v);CHKERRQ(ierr); |
| } |
| ierr = PetscFree(vs->x); CHKERRQ(ierr); |
| ierr = PetscFree(vs); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| static PetscErrorCode VecCreate_Comp_Private(Vec v, Vec *x, PetscInt nx, |
| PetscBool x_to_me, Vec_Comp_N *n) |
| { |
| Vec_Comp *s; |
| PetscErrorCode ierr; |
| PetscInt N=0, lN=0, i, k; |
| Vec_Comp *s; |
| PetscErrorCode ierr; |
| PetscInt N=0, lN=0, i, k; |
| PetscFunctionBegin; |
| /* Allocate a new Vec_Comp */ |
| if (v->data) { ierr = PetscFree(v->data); CHKERRQ(ierr); } |
| ierr = PetscNewLog(v, Vec_Comp, &s); CHKERRQ(ierr); |
| ierr = VecSetSizes(v, s->n->lN, s->n->N); CHKERRQ(ierr); |
| ierr = PetscObjectChangeTypeName((PetscObject)v,VECCOMP); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_BEGIN |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecCreate_Comp" |
| PetscErrorCode VecCreate_Comp(Vec V) |
| { |
| PetscErrorCode ierr; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = VecCreate_Comp_Private(V, PETSC_NULL, 0, PETSC_FALSE, PETSC_NULL); |
| CHKERRQ(ierr); |
| ierr = VecCreate_Comp_Private(V, PETSC_NULL, 0, PETSC_FALSE, PETSC_NULL);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| EXTERN_C_END |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecCreateComp" |
| PetscErrorCode VecCreateComp(MPI_Comm comm, PetscInt *Nx, |
| PetscInt n, const VecType t, |
| Vec Vparent, Vec *V) |
| { |
| PetscErrorCode ierr; |
| Vec *x; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec *x; |
| PetscInt i; |
| PetscFunctionBegin; |
| ierr = VecCreate(comm, V); CHKERRQ(ierr); |
| ierr = PetscMalloc(sizeof(Vec)*n, &x); CHKERRQ(ierr); |
| for(i=0; i<n; i++) { |
| ierr = VecSetType(x[i], t); CHKERRQ(ierr); |
| } |
| ierr = VecCreate_Comp_Private(*V, x, n, PETSC_TRUE, |
| Vparent?((Vec_Comp*)Vparent->data)->n:PETSC_NULL); |
| CHKERRQ(ierr); |
| Vparent?((Vec_Comp*)Vparent->data)->n:PETSC_NULL);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecCreateCompWithVecs" |
| PetscErrorCode VecCreateCompWithVecs(Vec *x, PetscInt n, |
| Vec Vparent, Vec *V) |
| PetscErrorCode VecCreateCompWithVecs(Vec *x, PetscInt n, Vec Vparent, Vec *V) |
| { |
| PetscErrorCode ierr; |
| PetscInt i; |
| PetscErrorCode ierr; |
| PetscInt i; |
| PetscFunctionBegin; |
| ierr = VecCreate(((PetscObject)x[0])->comm, V); CHKERRQ(ierr); |
| for(i=0; i<n; i++) { |
| ierr = PetscObjectReference((PetscObject)x[i]); CHKERRQ(ierr); |
| } |
| ierr = VecCreate_Comp_Private(*V, x, n, PETSC_FALSE, |
| Vparent?((Vec_Comp*)Vparent->data)->n:PETSC_NULL); |
| CHKERRQ(ierr); |
| Vparent?((Vec_Comp*)Vparent->data)->n:PETSC_NULL);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecDuplicate_Comp" |
| PetscErrorCode VecDuplicate_Comp(Vec win, Vec *V) |
| { |
| PetscErrorCode ierr; |
| Vec *x; |
| PetscInt i; |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| PetscErrorCode ierr; |
| Vec *x; |
| PetscInt i; |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| PetscFunctionBegin; |
| for(i=0; i<s->nx; i++) { |
| ierr = VecDuplicate(s->x[i], &x[i]); CHKERRQ(ierr); |
| } |
| ierr = VecCreate_Comp_Private(*V, x, s->nx, PETSC_TRUE, s->n); CHKERRQ(ierr); |
| ierr = VecCreate_Comp_Private(*V, x, s->nx, PETSC_TRUE, s->n);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecCompGetVecs" |
| PetscErrorCode VecCompGetVecs(Vec win, const Vec **x, PetscInt *n) |
| { |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| PetscFunctionBegin; |
| PetscValidVecComp(win); |
| if(x) *x = s->x; |
| if(n) *n = s->nx; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecCompSetVecs" |
| PetscErrorCode VecCompSetVecs(Vec win, Vec *x, PetscInt n) |
| { |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| PetscErrorCode ierr; |
| Vec_Comp *s = (Vec_Comp*)win->data; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| PetscValidVecComp(win); |
| if(x) { |
| s->nx = n; |
| } |
| s->n->n = n; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecAXPY_Comp" |
| PetscErrorCode VecAXPY_Comp(Vec v, PetscScalar alpha, Vec w) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| PetscValidVecComp(w); |
| for(i=0; i<vs->n->n; i++) { |
| ierr = VecAXPY(vs->x[i], alpha, ws->x[i]); CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecAYPX_Comp" |
| PetscErrorCode VecAYPX_Comp(Vec v, PetscScalar alpha, Vec w) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| PetscValidVecComp(w); |
| for(i=0; i<vs->n->n; i++) { |
| ierr = VecAYPX(vs->x[i], alpha, ws->x[i]); CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecAXPBY_Comp" |
| PetscErrorCode VecAXPBY_Comp(Vec v, PetscScalar alpha, PetscScalar beta, Vec w) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, *ws = (Vec_Comp*)w->data; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| PetscValidVecComp(w); |
| for(i=0; i<vs->n->n; i++) { |
| ierr = VecAXPBY(vs->x[i], alpha, beta, ws->x[i]); CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecMAXPY_Comp" |
| PetscErrorCode VecMAXPY_Comp(Vec v, PetscInt n, const PetscScalar *alpha, |
| Vec *w) |
| PetscErrorCode VecMAXPY_Comp(Vec v, PetscInt n, const PetscScalar *alpha, Vec *w) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| Vec *wx; |
| PetscInt i, j; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| Vec *wx; |
| PetscInt i, j; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| for(i=0; i<n; i++) PetscValidVecComp(w[i]); |
| } |
| ierr = PetscFree(wx); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecWAXPY_Comp" |
| PetscErrorCode VecWAXPY_Comp(Vec v, PetscScalar alpha, Vec w, Vec z) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data, |
| *zs = (Vec_Comp*)z->data; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data, |
| *zs = (Vec_Comp*)z->data; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| PetscValidVecComp(w); |
| PetscValidVecComp(z); |
| for(i=0; i<vs->n->n; i++) { |
| ierr = VecWAXPY(vs->x[i], alpha, ws->x[i], zs->x[i]); CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecAXPBYPCZ_Comp" |
| PetscErrorCode VecAXPBYPCZ_Comp(Vec v, PetscScalar alpha, PetscScalar beta, |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecGetSize_Comp" |
| PetscErrorCode VecGetSize_Comp(Vec v, PetscInt *size) |
| { |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| if (size) *size = vs->n->N; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecGetLocalSize_Comp" |
| PetscErrorCode VecGetLocalSize_Comp(Vec v, PetscInt *size) |
| { |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| if (size) *size = vs->n->lN; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecMax_Comp" |
| PetscErrorCode VecMax_Comp(Vec v, PetscInt *idx, PetscReal *z) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscInt idxp, s=0, s0; |
| PetscReal zp, z0; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscInt idxp, s=0, s0; |
| PetscReal zp, z0; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| if(!idx && !z) |
| PetscFunctionReturn(0); |
| } |
| } |
| if (z) *z = zp; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecMin_Comp" |
| PetscErrorCode VecMin_Comp(Vec v, PetscInt *idx, PetscReal *z) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscInt idxp, s=0, s0; |
| PetscReal zp, z0; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data; |
| PetscInt idxp, s=0, s0; |
| PetscReal zp, z0; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| if(!idx && !z) |
| PetscFunctionReturn(0); |
| } |
| } |
| if (z) *z = zp; |
| PetscFunctionReturn(0); |
| } |
| #undef __FUNCT__ |
| #define __FUNCT__ "VecMaxPointwiseDivide_Comp" |
| PetscErrorCode VecMaxPointwiseDivide_Comp(Vec v, Vec w, PetscReal *m) |
| { |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, |
| *ws = (Vec_Comp*)w->data; |
| PetscReal work; |
| PetscInt i; |
| PetscErrorCode ierr; |
| Vec_Comp *vs = (Vec_Comp*)v->data, *ws = (Vec_Comp*)w->data; |
| PetscReal work; |
| PetscInt i; |
| PetscFunctionBegin; |
| PetscValidVecComp(v); |
| PetscValidVecComp(w); |
| ierr = VecMaxPointwiseDivide(vs->x[i], ws->x[i], &work); CHKERRQ(ierr); |
| *m = PetscMax(*m, work); |
| } |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* operator matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| const EPSType type; |
| PetscScalar kr, ki; |
| PetscInt N, n=10, m, Istart, Iend, II, nev, maxit, i, j, its, nconv; |
| PetscBool flag; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| PetscReal error, tol, re, im; |
| PetscScalar kr, ki; |
| PetscMPIInt size; |
| PetscErrorCode ierr; |
| PetscInt N, n=10, nev, maxit, i, its, nconv; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); |
| PetscErrorCode MatLaplacian2D_Mult( Mat A, Vec x, Vec y ) |
| { |
| void *ctx; |
| PetscErrorCode ierr; |
| int nx, lo, j, one=1; |
| PetscScalar *px, *py, dmone=-1.0; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext( A, &ctx ); CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* operator matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| const EPSType type; |
| char filename[256]; |
| PetscViewer viewer; |
| PetscBool flg; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Vec v0; /* initial vector */ |
| Mat A; /* operator matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| PetscReal error, tol, re, im; |
| PetscScalar kr, ki; |
| PetscInt N, m=15, nev, maxit, i, its, nconv; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| { |
| const PetscReal cst = 0.5/(PetscReal)(m-1); |
| PetscReal pd, pu; |
| PetscInt Istart, Iend, i, j, jmax, ix=0; |
| PetscErrorCode ierr; |
| PetscInt Istart, Iend, i, j, jmax, ix=0; |
| PetscFunctionBegin; |
| ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A,B; /* matrices */ |
| EPS eps; /* eigenproblem solver context */ |
| const EPSType type; |
| char filename[256]; |
| PetscViewer viewer; |
| PetscBool flg; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* Grcar matrix */ |
| SVD svd; /* singular value solver context */ |
| PetscInt N=30, Istart, Iend, i, col[5], nconv1, nconv2; |
| PetscScalar value[] = { -1, 1, 1, 1, 1 }; |
| PetscReal sigma_1, sigma_n; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* eigenvalue problem matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| const EPSType type; |
| PetscInt N=30, n, i, col[3], Istart, Iend, nev, maxit, its, nconv; |
| PetscBool FirstBlock=PETSC_FALSE, LastBlock=PETSC_FALSE; |
| CTX_BRUSSEL *ctx; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "MatBrussel_Mult" |
| PetscErrorCode MatBrussel_Mult(Mat A,Vec x,Vec y) |
| { |
| PetscErrorCode ierr; |
| PetscInt n; |
| PetscScalar *px, *py; |
| CTX_BRUSSEL *ctx; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext(A,(void**)&ctx);CHKERRQ(ierr); |
| PetscErrorCode MatBrussel_GetDiagonal(Mat A,Vec diag) |
| { |
| Vec d1, d2; |
| PetscErrorCode ierr; |
| PetscInt n; |
| PetscScalar *pd; |
| MPI_Comm comm; |
| CTX_BRUSSEL *ctx; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = MatShellGetContext(A,(void**)&ctx);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* operator matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| ST st; /* spectral transformation context */ |
| SampleShellST *shell; /* user-defined spectral transform context */ |
| const EPSType type; |
| PetscReal error, tol, re, im; |
| PetscScalar kr, ki; |
| PetscScalar kr, ki, value[3]; |
| PetscInt n=30, i, col[3], Istart, Iend, FirstBlock=0, LastBlock=0, nev, maxit, its, nconv; |
| PetscScalar value[3]; |
| PetscBool isShell; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| PetscErrorCode SampleShellSTBackTransform(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi) |
| { |
| PetscInt j; |
| PetscFunctionBegin; |
| for (j=0;j<n;j++) { |
| eigr[j] = 1.0 / eigr[j]; |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| EPS eps; /* eigenproblem solver context */ |
| Mat A; /* operator matrix */ |
| Vec x; |
| PetscInt N, n=10, m, i, j, II, Istart, Iend, nev, maxit, its, nconv; |
| PetscScalar w; |
| PetscBool flag; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Vec v0,w0; /* initial vector */ |
| Vec *X,*Y; /* right and left eigenvectors */ |
| Mat A; /* operator matrix */ |
| PetscReal error1, error2, tol, re, im; |
| PetscScalar kr, ki; |
| PetscInt nev, maxit, i, its, nconv, N, m=15; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| { |
| const PetscReal cst = 0.5/(PetscReal)(m-1); |
| PetscReal pd, pu; |
| PetscInt i, j, jmax, ix=0, Istart, Iend; |
| PetscErrorCode ierr; |
| PetscInt i, j, jmax, ix=0, Istart, Iend; |
| PetscFunctionBegin; |
| ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A, B; /* matrices */ |
| EPS eps; /* eigenproblem solver context */ |
| ST st; /* spectral transformation context */ |
| PetscScalar kr, ki; |
| PetscInt N, n=10, m, Istart, Iend, II, nev, maxit, i, j, its, nconv, nulldim=0; |
| PetscBool flag; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* operator matrix */ |
| SVD svd; /* singular value problem solver context */ |
| const SVDType type; |
| char filename[256]; |
| PetscViewer viewer; |
| PetscBool flg; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat A; /* operator matrix */ |
| Vec u,v; /* left and right singular vectors */ |
| SVD svd; /* singular value problem solver context */ |
| const SVDType type; |
| PetscReal error, tol, sigma, mu=PETSC_SQRT_MACHINE_EPSILON; |
| PetscInt n=100, i, j, Istart, Iend, nsv, maxit, its, nconv; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat M, C, K; /* problem matrices */ |
| QEP qep; /* quadratic eigenproblem solver context */ |
| const QEPType type; |
| PetscScalar kr, ki; |
| PetscInt N, n=10, m, Istart, Iend, II, nev, maxit, i, j, its, nconv; |
| PetscBool flag; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Mat M, C, K; /* problem matrices */ |
| QEP qep; /* quadratic eigenproblem solver context */ |
| const QEPType type; |
| char filename[256]; |
| PetscViewer viewer; |
| PetscBool flg; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| */ |
| PetscErrorCode MyEigenSort(EPS eps, PetscScalar ar, PetscScalar ai, PetscScalar br, PetscScalar bi, PetscInt *r, void *ctx); |
| PetscErrorCode MatMarkovModel( PetscInt m, Mat A ); |
| #undef __FUNCT__ |
| #define __FUNCT__ "main" |
| int main( int argc, char **argv ) |
| { |
| PetscErrorCode ierr; |
| Vec v0; /* initial vector */ |
| Mat A; /* operator matrix */ |
| EPS eps; /* eigenproblem solver context */ |
| PetscReal error, tol, re, im; |
| PetscScalar kr, ki, target=0.5; |
| PetscInt N, m=15, nev, maxit, i, its, nconv; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| { |
| const PetscReal cst = 0.5/(PetscReal)(m-1); |
| PetscReal pd, pu; |
| PetscErrorCode ierr; |
| PetscInt Istart, Iend, i, j, jmax, ix=0; |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); |
| */ |
| PetscErrorCode MyEigenSort(EPS eps, PetscScalar ar, PetscScalar ai, PetscScalar br, PetscScalar bi, PetscInt *r, void *ctx) |
| { |
| PetscScalar target = *(PetscScalar*)ctx; |
| PetscReal da,db; |
| PetscBool aisright, bisright; |
| PetscScalar target = *(PetscScalar*)ctx; |
| PetscReal da,db; |
| PetscBool aisright, bisright; |
| PetscFunctionBegin; |
| if (PetscRealPart(target) < PetscRealPart(ar)) aisright = PETSC_TRUE; |
| else aisright = PETSC_FALSE; |
| if (PetscRealPart(target) < PetscRealPart(br)) bisright = PETSC_TRUE; |
| *r = -1; /* 'a' is on the right */ |
| else |
| *r = 1; /* 'b' is on the right */ |
| PetscFunctionReturn(0); |
| } |
| EPS eps; /* eigenproblem solver context */ |
| const EPSType type; |
| PetscReal error, tol, re, im; |
| PetscScalar kr, ki; |
| PetscScalar kr, ki, value[3]; |
| Vec xr, xi; |
| PetscErrorCode ierr; |
| PetscInt n=30, i, Istart, Iend, col[3], nev, maxit, its, nconv; |
| PetscBool FirstBlock=PETSC_FALSE, LastBlock=PETSC_FALSE; |
| PetscScalar value[3]; |
| PetscErrorCode ierr; |
| SlepcInitialize(&argc,&argv,(char*)0,help); |
| PetscBool has; |
| PetscFunctionBegin; |
| #if !defined(PETSC_USE_COMPLEX) |
| ierr = PetscTypeCompare((PetscObject)A,MATSEQSBAIJ,is);CHKERRQ(ierr); |
| if (*is) PetscFunctionReturn(0); |
| ierr = VecDestroy(&x);CHKERRQ(ierr); |
| ierr = VecDestroy(&w1);CHKERRQ(ierr); |
| ierr = VecDestroy(&w2);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscReal SlepcAbsEigenvalue(PetscScalar x,PetscScalar y) |
| { |
| PetscReal xabs,yabs,w,z,t; |
| PetscFunctionBegin; |
| xabs = PetscAbsReal(x); |
| yabs = PetscAbsReal(y); |
| PetscFunctionBegin; |
| PetscValidHeaderSpecific(mat,MAT_CLASSID,1); |
| PetscValidPointer(newmat,2); |
| ierr = PetscObjectGetComm((PetscObject)mat,&comm);CHKERRQ(ierr); |
| ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); |
| /* convert matrix to MatSeqDense */ |
| ierr = MatConvert(mat,MATSEQDENSE,MAT_INITIAL_MATRIX,newmat);CHKERRQ(ierr); |
| } |
| PetscFunctionReturn(0); |
| } |
| #define __FUNCT__ "SlepcMatTile_SEQAIJ" |
| static PetscErrorCode SlepcMatTile_SEQAIJ(PetscScalar a,Mat A,PetscScalar b,Mat B,PetscScalar c,Mat C,PetscScalar d,Mat D,Mat G) |
| { |
| PetscErrorCode ierr; |
| PetscInt i,j,M1,M2,N1,N2; |
| PetscInt *nnz,ncols,*scols; |
| PetscScalar *svals,*buf; |
| PetscErrorCode ierr; |
| PetscInt i,j,M1,M2,N1,N2,*nnz,ncols,*scols; |
| PetscScalar *svals,*buf; |
| const PetscInt *cols; |
| const PetscScalar *vals; |
| static PetscErrorCode SlepcMatTile_MPIAIJ(PetscScalar a,Mat A,PetscScalar b,Mat B,PetscScalar c,Mat C,PetscScalar d,Mat D,Mat G) |
| { |
| PetscErrorCode ierr; |
| PetscInt p,i,j,N1,N2,m1,m2,n1,n2; |
| PetscInt *dnz,*onz,ncols,*scols,start,gstart; |
| PetscInt *map1,*map2,np; |
| PetscScalar *svals,*buf; |
| PetscInt p,i,j,N1,N2,m1,m2,n1,n2,*map1,*map2,np; |
| PetscInt *dnz,*onz,ncols,*scols,start,gstart; |
| PetscScalar *svals,*buf; |
| const PetscInt *cols,*mapptr1,*mapptr2; |
| const PetscScalar *vals; |
| } |
| else SETERRQ(((PetscObject)A)->comm,PETSC_ERR_SUP,"Not implemented for this matrix type"); |
| } |
| ierr = MatAssemblyBegin(*G,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); |
| ierr = MatAssemblyEnd(*G,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = SlepcUpdateStrideVectors(n_,V,s,1,e,Q,ldq_,qtrans);CHKERRQ(ierr); |
| PetscFunctionReturn(0); |
| } |
| PetscErrorCode SlepcConvMonitorDestroy(SlepcConvMonitor *ctx) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| if (!*ctx) PetscFunctionReturn(0); |
| ierr = PetscViewerASCIIMonitorDestroy(&(*ctx)->viewer);CHKERRQ(ierr); |
| PetscErrorCode SlepcInitialize_Packages(void) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = STInitializePackage(PETSC_NULL);CHKERRQ(ierr); |
| ierr = EPSInitializePackage(PETSC_NULL);CHKERRQ(ierr); |
| PetscErrorCode SlepcInitialize_LogEvents(void) |
| { |
| PetscErrorCode ierr; |
| PetscFunctionBegin; |
| ierr = PetscLogEventRegister("UpdateVectors",0,&SLEPC_UpdateVectors);CHKERRQ(ierr); |
| ierr = PetscLogEventRegister("VecMAXPBY",0,&SLEPC_VecMAXPBY);CHKERRQ(ierr); |
| @*/ |
| PetscErrorCode SlepcInitialize(int *argc,char ***args,char file[],const char help[]) |
| { |
| PetscErrorCode ierr; |
| PetscErrorCode info=0; |
| PetscErrorCode ierr,info=0; |
| PetscBool flg; |
| PetscFunctionBegin; |
| PetscFunctionBegin; |
| if (SlepcInitializeCalled) { |
| PetscFunctionReturn(0); |
| } |
| info = PetscSetHelpVersionFunctions(SlepcPrintHelpIntro,SlepcPrintVersion);CHKERRQ(info); |
| ierr = PetscInitialized(&flg);CHKERRQ(ierr); |
| if (!flg) { |
| info = PetscInitialize(argc,args,file,help);CHKERRQ(info); |
| PetscFunctionBegin; |
| PetscInfo(0,"SLEPc successfully ended!\n"); |
| if (SlepcBeganPetsc) { |
| info = PetscFinalize();CHKERRQ(info); |
| } |
| SlepcInitializeCalled = PETSC_FALSE; |
| PetscFunctionReturn(info); |
| } |