Actual source code: dvd_schm.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7:
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) \
25: ( (b)->max_size_V + (b)->max_size_auxV + (b)->max_size_auxS + \
26: (b)->own_vecs + (b)->own_scalars + (b)->max_size_oldX )
30: PetscErrorCode dvd_schm_basic_preconf(dvdDashboard *d, dvdBlackboard *b,
31: PetscInt max_size_V, PetscInt mpd, PetscInt min_size_V, PetscInt bs,
32: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk, PC pc,
33: HarmType_t harmMode, KSP ksp, InitType_t init, PetscTruth allResiduals)
34: {
36: PetscInt check_sum0, check_sum1;
40: PetscMemzero(b, sizeof(dvdBlackboard));
41: b->state = DVD_STATE_PRECONF;
43: for(check_sum0=-1,check_sum1=DVD_CHECKSUM(b); check_sum0 != check_sum1;
44: check_sum0 = check_sum1, check_sum1 = DVD_CHECKSUM(b)) {
45: b->own_vecs = b->own_scalars = 0;
47: /* Setup basic management of V */
48: dvd_managementV_basic(d, b, bs, max_size_V, mpd, min_size_V, plusk,
49: harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE,
50: allResiduals);
51:
52:
53: /* Setup the initial subspace for V */
54: if (size_initV) {
55: dvd_initV_user(d, b, size_initV, ini_size_V);
56: } else switch(init) {
57: case DVD_INITV_CLASSIC:
58: dvd_initV_classic(d, b, ini_size_V); break;
59: case DVD_INITV_KRYLOV:
60: dvd_initV_krylov(d, b, ini_size_V); break;
61: }
62:
63: /* Setup the convergence in order to use the SLEPc convergence test */
64: dvd_testconv_slepc(d, b);
65:
66: /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
67: dvd_calcpairs_qz(d, b, PETSC_NULL);
68: if (harmMode != DVD_HARM_NONE) {
69: dvd_harm_conf(d, b, harmMode, PETSC_FALSE, 0.0);
70: }
71:
72: /* Setup the preconditioner */
73: dvd_static_precond_PC(d, b, pc);
75: /* Setup the method for improving the eigenvectors */
76: dvd_improvex_jd(d, b, ksp, bs);
77: dvd_improvex_jd_proj_uv(d, b, DVD_PROJ_KBXX);
78: dvd_improvex_jd_lit_const(d, b, 0, 0.0, 0.0);
79: }
81: return(0);
82: }
87: PetscErrorCode dvd_schm_basic_conf(dvdDashboard *d, dvdBlackboard *b,
88: PetscInt max_size_V, PetscInt mpd, PetscInt min_size_V, PetscInt bs,
89: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk, PC pc,
90: IP ip, HarmType_t harmMode, PetscTruth fixedTarget, PetscScalar t, KSP ksp,
91: PetscReal fix, InitType_t init, PetscTruth allResiduals)
92: {
93: PetscInt check_sum0, check_sum1, maxits;
94: Vec *fv;
95: PetscScalar *fs;
96: PetscReal tol;
97: PetscErrorCode ierr;
101: b->state = DVD_STATE_CONF;
102: check_sum0 = DVD_CHECKSUM(b);
103: b->own_vecs = 0; b->own_scalars = 0;
104: fv = b->free_vecs; fs = b->free_scalars;
106: /* Setup basic management of V */
107: dvd_managementV_basic(d, b, bs, max_size_V, mpd, min_size_V, plusk,
108: harmMode==DVD_HARM_NONE?PETSC_FALSE:PETSC_TRUE,
109: allResiduals);
110:
112: /* Setup the initial subspace for V */
113: if (size_initV) {
114: dvd_initV_user(d, b, size_initV, ini_size_V);
115: } else switch(init) {
116: case DVD_INITV_CLASSIC:
117: dvd_initV_classic(d, b, ini_size_V); break;
118: case DVD_INITV_KRYLOV:
119: dvd_initV_krylov(d, b, ini_size_V); break;
120: }
122: /* Setup the convergence in order to use the SLEPc convergence test */
123: dvd_testconv_slepc(d, b);
125: /* Setup Raileigh-Ritz for selecting the best eigenpairs in V */
126: dvd_calcpairs_qz(d, b, ip);
127: if (harmMode != DVD_HARM_NONE) {
128: dvd_harm_conf(d, b, harmMode, fixedTarget, t);
129: }
131: /* Setup the preconditioner */
132: dvd_static_precond_PC(d, b, pc);
134: /* Setup the method for improving the eigenvectors */
135: dvd_improvex_jd(d, b, ksp, bs);
136: dvd_improvex_jd_proj_uv(d, b, DVD_IS(d->sEP, DVD_EP_HERMITIAN)?
137: DVD_PROJ_KBXZ:DVD_PROJ_KBXY);
138:
139: KSPGetTolerances(ksp, &tol, PETSC_NULL, PETSC_NULL, &maxits);
140:
141: dvd_improvex_jd_lit_const(d, b, maxits, tol, fix);
143: check_sum1 = DVD_CHECKSUM(b);
144: if ((check_sum0 != check_sum1) ||
145: (b->free_vecs - fv > b->own_vecs) ||
146: (b->free_scalars - fs > b->own_scalars))
147: SETERRQ(1, "Something awful happened!");
148:
149: return(0);
150: }