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-2011, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.
  9:       
 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 <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:       MatMultTranspose(svd->A,x,y);
 39:     }
 40:   } else {
 41:     if (svd->A) {
 42:       MatMult(svd->A,x,y);
 43:     } else {
 44:       MatMultTranspose(svd->AT,x,y);
 45:     }
 46:   }
 47:   return(0);
 48: }

 52: PetscErrorCode SVDMatGetVecs(SVD svd,Vec *x,Vec *y)
 53: {

 57:   if (svd->A) {
 58:     MatGetVecs(svd->A,x,y);
 59:   } else {
 60:     MatGetVecs(svd->AT,y,x);
 61:   }
 62:   return(0);
 63: }

 67: PetscErrorCode SVDMatGetSize(SVD svd,PetscInt *m,PetscInt *n)
 68: {

 72:   if (svd->A) {
 73:     MatGetSize(svd->A,m,n);
 74:   } else {
 75:     MatGetSize(svd->AT,n,m);
 76:   }
 77:   return(0);
 78: }

 82: PetscErrorCode SVDMatGetLocalSize(SVD svd,PetscInt *m,PetscInt *n)
 83: {

 87:   if (svd->A) {
 88:     MatGetLocalSize(svd->A,m,n);
 89:   } else {
 90:     MatGetLocalSize(svd->AT,n,m);
 91:   }
 92:   return(0);
 93: }