Actual source code: dvd_orth.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

 26: PetscErrorCode dvd_orthV(BV V,PetscInt V_new_s,PetscInt V_new_e,PetscRandom rand)
 27: {
 29:   PetscInt       i,j,l,k;
 30:   PetscBool      lindep;
 31:   PetscReal      norm;

 34:   BVGetActiveColumns(V,&l,&k);
 35:   for (i=V_new_s;i<V_new_e;i++) {
 36:     BVSetActiveColumns(V,0,i);
 37:     for (j=0;j<3;j++) {
 38:       if (j>0) {
 39:         BVSetRandomColumn(V,i,rand);
 40:         PetscInfo1(V,"Orthonormalization problems adding the vector %D to the searching subspace\n",i);
 41:       }
 42:       BVOrthogonalizeColumn(V,i,NULL,&norm,&lindep);
 43:       if (!lindep && (PetscAbsReal(norm) > PETSC_SQRT_MACHINE_EPSILON)) break;
 44:     }
 45:     if (lindep || (PetscAbsReal(norm) < PETSC_SQRT_MACHINE_EPSILON)) SETERRQ(PetscObjectComm((PetscObject)V),1, "Error during the orthonormalization of the vectors");
 46:     BVScaleColumn(V,i,1.0/norm);
 47:   }
 48:   BVSetActiveColumns(V,l,k);
 49:   return(0);
 50: }