Actual source code: dvdschm.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */

 11: #include "davidson.h"

 13: #define DVD_CHECKSUM(b) ((b)->max_size_V + (b)->max_size_oldX)

 15: PetscErrorCode dvd_schm_basic_preconf(dvdDashboard *d,dvdBlackboard *b,PetscInt mpd,PetscInt min_size_V,PetscInt bs,PetscInt ini_size_V,PetscInt size_initV,PetscInt plusk,HarmType_t harmMode,KSP ksp,InitType_t init,PetscBool allResiduals,PetscBool orth,PetscBool doubleexp)
 16: {
 18:   PetscInt       check_sum0,check_sum1;

 21:   PetscMemzero(b,sizeof(dvdBlackboard));
 22:   b->state = DVD_STATE_PRECONF;

 24:   for (check_sum0=-1,check_sum1=DVD_CHECKSUM(b); check_sum0 != check_sum1; check_sum0 = check_sum1, check_sum1 = DVD_CHECKSUM(b)) {

 26:     /* Setup basic management of V */
 27:     dvd_managementV_basic(d,b,bs,mpd,min_size_V,plusk,PetscNot(harmMode==DVD_HARM_NONE),allResiduals);

 29:     /* Setup the initial subspace for V */
 30:     dvd_initV(d,b,ini_size_V,size_initV,(init==DVD_INITV_KRYLOV)?PETSC_TRUE:PETSC_FALSE);

 32:     /* Setup the convergence in order to use the SLEPc convergence test */
 33:     dvd_testconv_slepc(d,b);

 35:     /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
 36:     dvd_calcpairs_qz(d,b,orth,PetscNot(harmMode==DVD_HARM_NONE));
 37:     if (harmMode != DVD_HARM_NONE) {
 38:       dvd_harm_conf(d,b,harmMode,PETSC_FALSE,0.0);
 39:     }

 41:     /* Setup the method for improving the eigenvectors */
 42:     if (doubleexp) {
 43:       dvd_improvex_gd2(d,b,ksp,bs);
 44:     } else {
 45:       dvd_improvex_jd(d,b,ksp,bs,PETSC_FALSE);
 46:       dvd_improvex_jd_proj_uv(d,b);
 47:       dvd_improvex_jd_lit_const(d,b,0,0.0,0.0);
 48:     }
 49:   }
 50:   return(0);
 51: }

 53: PetscErrorCode dvd_schm_basic_conf(dvdDashboard *d,dvdBlackboard *b,PetscInt mpd,PetscInt min_size_V,PetscInt bs,PetscInt ini_size_V,PetscInt size_initV,PetscInt plusk,HarmType_t harmMode,PetscBool fixedTarget,PetscScalar t,KSP ksp,PetscReal fix,InitType_t init,PetscBool allResiduals,PetscBool orth,PetscBool dynamic,PetscBool doubleexp)
 54: {
 55:   PetscInt       check_sum0,check_sum1,maxits;
 56:   PetscReal      tol;

 60:   b->state = DVD_STATE_CONF;
 61:   check_sum0 = DVD_CHECKSUM(b);

 63:   /* Setup basic management of V */
 64:   dvd_managementV_basic(d,b,bs,mpd,min_size_V,plusk,PetscNot(harmMode==DVD_HARM_NONE),allResiduals);

 66:   /* Setup the initial subspace for V */
 67:   dvd_initV(d,b,ini_size_V,size_initV,(init==DVD_INITV_KRYLOV)?PETSC_TRUE:PETSC_FALSE);

 69:   /* Setup the convergence in order to use the SLEPc convergence test */
 70:   dvd_testconv_slepc(d,b);

 72:   /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
 73:   dvd_calcpairs_qz(d,b,orth,PetscNot(harmMode==DVD_HARM_NONE));
 74:   if (harmMode != DVD_HARM_NONE) {
 75:     dvd_harm_conf(d,b,harmMode,fixedTarget,t);
 76:   }

 78:   /* Setup the method for improving the eigenvectors */
 79:   if (doubleexp) {
 80:     dvd_improvex_gd2(d,b,ksp,bs);
 81:   } else {
 82:     dvd_improvex_jd(d,b,ksp,bs,dynamic);
 83:     dvd_improvex_jd_proj_uv(d,b);
 84:     KSPGetTolerances(ksp,&tol,NULL,NULL,&maxits);
 85:     dvd_improvex_jd_lit_const(d,b,maxits,tol,fix);
 86:   }

 88:   check_sum1 = DVD_CHECKSUM(b);
 89:   if (check_sum0 != check_sum1) SETERRQ(PETSC_COMM_SELF,1,"Something awful happened");
 90:   return(0);
 91: }