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-2009, Universidad 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

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

 36:   svd->matvecs++;
 37:   if (trans) {
 38:     if (svd->AT) {
 39:       MatMult(svd->AT,x,y);
 40:     } else {
 41:       MatMultTranspose(svd->A,x,y);
 42:     }
 43:   } else {
 44:     if (svd->A) {
 45:       MatMult(svd->A,x,y);
 46:     } else {
 47:       MatMultTranspose(svd->AT,x,y);
 48:     }
 49:   }
 50:   return(0);
 51: }

 55: PetscErrorCode SVDMatGetVecs(SVD svd,Vec *x,Vec *y)
 56: {

 61:   if (svd->A) {
 62:     MatGetVecs(svd->A,x,y);
 63:   } else {
 64:     MatGetVecs(svd->AT,y,x);
 65:   }
 66:   return(0);
 67: }

 71: PetscErrorCode SVDMatGetSize(SVD svd,PetscInt *m,PetscInt *n)
 72: {

 77:   if (svd->A) {
 78:     MatGetSize(svd->A,m,n);
 79:   } else {
 80:     MatGetSize(svd->AT,n,m);
 81:   }
 82:   return(0);
 83: }

 87: PetscErrorCode SVDMatGetLocalSize(SVD svd,PetscInt *m,PetscInt *n)
 88: {

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