Actual source code: dvd_schm.c

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

  6:    This file is part of SLEPc.

  8:    SLEPc is free software: you can redistribute it and/or modify it under  the
  9:    terms of version 3 of the GNU Lesser General Public License as published by
 10:    the Free Software Foundation.

 12:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 13:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 14:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 15:    more details.

 17:    You  should have received a copy of the GNU Lesser General  Public  License
 18:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 19:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 20: */

 22:  #include davidson.h

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

 28: 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,PetscInt cX_proj,PetscInt cX_impr,Method_t method)
 29: {
 31:   PetscInt       check_sum0, check_sum1;

 34:   PetscMemzero(b, sizeof(dvdBlackboard));
 35:   b->state = DVD_STATE_PRECONF;

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

 40:     /* Setup basic management of V */
 41:     dvd_managementV_basic(d, b, bs, mpd, min_size_V, plusk,
 42:                                harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE,
 43:                                allResiduals);

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

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

 52:     /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
 53:     dvd_calcpairs_qz(d, b, orth, cX_proj,
 54:                 harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE);
 55:     if (harmMode != DVD_HARM_NONE) {
 56:       dvd_harm_conf(d, b, harmMode, PETSC_FALSE, 0.0);
 57:     }

 59:     /* Setup the method for improving the eigenvectors */
 60:     switch (method) {
 61:       case DVD_METH_GD:
 62:       case DVD_METH_JD:
 63:       dvd_improvex_jd(d, b, ksp, bs, cX_impr, PETSC_FALSE);
 64:       dvd_improvex_jd_proj_uv(d, b, DVD_PROJ_KZX);
 65:       dvd_improvex_jd_lit_const(d, b, 0, 0.0, 0.0);
 66:       break;
 67:       case DVD_METH_GD2:
 68:       dvd_improvex_gd2(d, b, ksp, bs);
 69:       break;
 70:     }
 71:   }
 72:   return(0);
 73: }

 77: 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,PetscInt cX_proj,PetscInt cX_impr,PetscBool dynamic,Method_t method)
 78: {
 79:   PetscInt        check_sum0, check_sum1, maxits;
 80:   PetscReal       tol;
 81:   PetscErrorCode  ierr;

 84:   b->state = DVD_STATE_CONF;
 85:   check_sum0 = DVD_CHECKSUM(b);

 87:   /* Setup basic management of V */
 88:   dvd_managementV_basic(d, b, bs, mpd, min_size_V, plusk,
 89:                         harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE,
 90:                         allResiduals);

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

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

 99:   /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
100:   dvd_calcpairs_qz(d, b, orth, cX_proj,
101:                 harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE);
102:   if (harmMode != DVD_HARM_NONE) {
103:     dvd_harm_conf(d, b, harmMode, fixedTarget, t);
104:   }

106:   /* Setup the method for improving the eigenvectors */
107:   switch (method) {
108:     case DVD_METH_GD:
109:     case DVD_METH_JD:
110:       dvd_improvex_jd(d, b, ksp, bs, cX_impr, dynamic);
111:       dvd_improvex_jd_proj_uv(d, b, DVD_PROJ_KZX);
112:       KSPGetTolerances(ksp, &tol, NULL, NULL, &maxits);
113:       dvd_improvex_jd_lit_const(d, b, maxits, tol, fix);
114:       break;
115:     case DVD_METH_GD2:
116:       dvd_improvex_gd2(d, b, ksp, bs);
117:       break;
118:   }

120:   check_sum1 = DVD_CHECKSUM(b);
121:   if (check_sum0 != check_sum1) SETERRQ(PETSC_COMM_SELF,1, "Something awful happened");
122:   return(0);
123: }