Actual source code: dvd_testconv.c
slepc-3.5.2 2014-10-10
1: /*
2: SLEPc eigensolver: "davidson"
4: Step: test for convergence
6: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
7: SLEPc - Scalable Library for Eigenvalue Problem Computations
8: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
10: This file is part of SLEPc.
12: SLEPc is free software: you can redistribute it and/or modify it under the
13: terms of version 3 of the GNU Lesser General Public License as published by
14: the Free Software Foundation.
16: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
17: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
19: more details.
21: You should have received a copy of the GNU Lesser General Public License
22: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
23: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
24: */
26: #include davidson.h
28: static PetscBool dvd_testconv_basic_0(dvdDashboard*,PetscScalar,PetscScalar,PetscReal,PetscReal*);
29: static PetscBool dvd_testconv_slepc_0(dvdDashboard*,PetscScalar,PetscScalar,PetscReal,PetscReal*);
33: PetscErrorCode dvd_testconv_basic(dvdDashboard *d, dvdBlackboard *b)
34: {
38: /* Setup the step */
39: if (b->state >= DVD_STATE_CONF) {
40: PetscFree(d->testConv_data);
41: d->testConv = dvd_testconv_basic_0;
42: }
43: return(0);
44: }
48: static PetscBool dvd_testconv_basic_0(dvdDashboard *d,PetscScalar eigvr,PetscScalar eigvi,PetscReal r,PetscReal *err)
49: {
50: PetscBool conv;
51: PetscReal eig_norm, errest;
54: eig_norm = SlepcAbsEigenvalue(eigvr, eigvi);
55: errest = r/eig_norm;
56: conv = (errest <= d->tol) ? PETSC_TRUE : PETSC_FALSE;
57: if (err) *err = errest;
58: PetscFunctionReturn(conv);
59: }
63: PetscErrorCode dvd_testconv_slepc(dvdDashboard *d, dvdBlackboard *b)
64: {
68: /* Setup the step */
69: if (b->state >= DVD_STATE_CONF) {
70: PetscFree(d->testConv_data);
71: d->testConv = dvd_testconv_slepc_0;
72: }
73: return(0);
74: }
78: static PetscBool dvd_testconv_slepc_0(dvdDashboard *d,PetscScalar eigvr,PetscScalar eigvi,PetscReal r,PetscReal *err)
79: {
83: (*d->eps->converged)(d->eps, eigvr, eigvi, r, err, d->eps->convergedctx);
84: CHKERRABORT(PetscObjectComm((PetscObject)d->eps), ierr);
85: PetscFunctionReturn(*err<d->eps->tol ? PETSC_TRUE : PETSC_FALSE);
86: }