Actual source code: svdmat.c

  1: /*
  2:      SVD routines for accessing the problem matrix.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2013, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

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

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

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

 24: #include <slepc-private/svdimpl.h>      /*I "slepcsvd.h" I*/

 28: PetscErrorCode SVDMatMult(SVD svd,PetscBool trans,Vec x,Vec y)
 29: {

 33:   svd->matvecs++;
 34:   if (trans) {
 35:     if (svd->AT) {
 36:       MatMult(svd->AT,x,y);
 37:     } else {
 38: #if defined(PETSC_USE_COMPLEX)
 39:       MatMultHermitianTranspose(svd->A,x,y);
 40: #else
 41:       MatMultTranspose(svd->A,x,y);
 42: #endif
 43:     }
 44:   } else {
 45:     if (svd->A) {
 46:       MatMult(svd->A,x,y);
 47:     } else {
 48: #if defined(PETSC_USE_COMPLEX)
 49:       MatMultHermitianTranspose(svd->AT,x,y);
 50: #else
 51:       MatMultTranspose(svd->AT,x,y);
 52: #endif
 53:     }
 54:   }
 55:   return(0);
 56: }

 60: PetscErrorCode SVDMatGetVecs(SVD svd,Vec *x,Vec *y)
 61: {

 65:   if (svd->A) {
 66:     MatGetVecs(svd->A,x,y);
 67:   } else {
 68:     MatGetVecs(svd->AT,y,x);
 69:   }
 70:   return(0);
 71: }

 75: PetscErrorCode SVDMatGetSize(SVD svd,PetscInt *m,PetscInt *n)
 76: {

 80:   if (svd->A) {
 81:     MatGetSize(svd->A,m,n);
 82:   } else {
 83:     MatGetSize(svd->AT,n,m);
 84:   }
 85:   return(0);
 86: }

 90: PetscErrorCode SVDMatGetLocalSize(SVD svd,PetscInt *m,PetscInt *n)
 91: {

 95:   if (svd->A) {
 96:     MatGetLocalSize(svd->A,m,n);
 97:   } else {
 98:     MatGetLocalSize(svd->AT,n,m);
 99:   }
100:   return(0);
101: }