| Line 26... |
Line 26... |
#define __FUNCT__ "SVDSetUp_LANCZOS"
|
#define __FUNCT__ "SVDSetUp_LANCZOS"
|
PetscErrorCode SVDSetUp_LANCZOS(SVD svd)
|
PetscErrorCode SVDSetUp_LANCZOS(SVD svd)
|
{
|
{
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data;
|
SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data;
|
PetscInt N;
|
PetscInt i,N;
|
int i;
|
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
ierr = SVDMatGetSize(svd,PETSC_NULL,&N);CHKERRQ(ierr);
|
ierr = SVDMatGetSize(svd,PETSC_NULL,&N);CHKERRQ(ierr);
|
if (svd->ncv == PETSC_DECIDE)
|
if (svd->ncv == PETSC_DECIDE)
|
svd->ncv = PetscMin(N,PetscMax(2*svd->nsv,10));
|
svd->ncv = PetscMin(N,PetscMax(2*svd->nsv,10));
|
| Line 48... |
Line 47... |
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|
}
|
}
|
|
|
#undef __FUNCT__
|
#undef __FUNCT__
|
#define __FUNCT__ "SVDTwoSideLanczos"
|
#define __FUNCT__ "SVDTwoSideLanczos"
|
PetscErrorCode SVDTwoSideLanczos(SVD svd,PetscReal *alpha,PetscReal *beta,Vec *V,Vec v,Vec *U,int k,int n,PetscScalar* work,Vec wv,Vec wu)
|
PetscErrorCode SVDTwoSideLanczos(SVD svd,PetscReal *alpha,PetscReal *beta,Vec *V,Vec v,Vec *U,PetscInt k,PetscInt n,PetscScalar* work,Vec wv,Vec wu)
|
{
|
{
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
int i;
|
PetscInt i;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
ierr = SVDMatMult(svd,PETSC_FALSE,V[k],U[k]);CHKERRQ(ierr);
|
ierr = SVDMatMult(svd,PETSC_FALSE,V[k],U[k]);CHKERRQ(ierr);
|
ierr = IPOrthogonalize(svd->ip,k,PETSC_NULL,U,U[k],work,alpha,PETSC_NULL,wu);CHKERRQ(ierr);
|
ierr = IPOrthogonalize(svd->ip,k,PETSC_NULL,U,U[k],work,alpha,PETSC_NULL,wu);CHKERRQ(ierr);
|
ierr = VecScale(U[k],1.0/alpha[0]);CHKERRQ(ierr);
|
ierr = VecScale(U[k],1.0/alpha[0]);CHKERRQ(ierr);
|
| Line 73... |
Line 72... |
PetscFunctionReturn(0);
|
PetscFunctionReturn(0);
|
}
|
}
|
|
|
#undef __FUNCT__
|
#undef __FUNCT__
|
#define __FUNCT__ "SVDOneSideLanczos"
|
#define __FUNCT__ "SVDOneSideLanczos"
|
static PetscErrorCode SVDOneSideLanczos(SVD svd,PetscReal *alpha,PetscReal *beta,Vec *V,Vec v,Vec u,Vec u_1,int k,int n,PetscScalar* work,Vec wv)
|
static PetscErrorCode SVDOneSideLanczos(SVD svd,PetscReal *alpha,PetscReal *beta,Vec *V,Vec v,Vec u,Vec u_1,PetscInt k,PetscInt n,PetscScalar* work,Vec wv)
|
{
|
{
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
int i,j;
|
PetscInt i,j;
|
PetscReal a,b;
|
PetscReal a,b;
|
Vec temp;
|
Vec temp;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
ierr = SVDMatMult(svd,PETSC_FALSE,V[k],u);CHKERRQ(ierr);
|
ierr = SVDMatMult(svd,PETSC_FALSE,V[k],u);CHKERRQ(ierr);
|
| Line 136... |
Line 135... |
#else
|
#else
|
PetscErrorCode ierr;
|
PetscErrorCode ierr;
|
SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data;
|
SVD_LANCZOS *lanczos = (SVD_LANCZOS *)svd->data;
|
PetscReal *alpha,*beta,norm,*work,*Q,*PT;
|
PetscReal *alpha,*beta,norm,*work,*Q,*PT;
|
PetscScalar *swork;
|
PetscScalar *swork;
|
PetscInt *perm;
|
PetscInt i,j,k,m,n,info,nwork=0,*iwork,*perm;
|
int i,j,k,m,n,info,nwork=0,*iwork;
|
|
Vec v,u,u_1,wv,wu,*workV,*workU,*permV,*permU;
|
Vec v,u,u_1,wv,wu,*workV,*workU,*permV,*permU;
|
PetscTruth conv;
|
PetscTruth conv;
|
|
|
PetscFunctionBegin;
|
PetscFunctionBegin;
|
/* allocate working space */
|
/* allocate working space */
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n,&alpha);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n,&alpha);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n,&beta);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n,&beta);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n*svd->n,&Q);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n*svd->n,&Q);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n*svd->n,&PT);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*svd->n*svd->n,&PT);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*(3*svd->n+4)*svd->n,&work);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscReal)*(3*svd->n+4)*svd->n,&work);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(int)*8*svd->n,&iwork);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscInt)*8*svd->n,&iwork);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscScalar)*svd->n,&swork);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(PetscScalar)*svd->n,&swork);CHKERRQ(ierr);
|
ierr = VecDuplicate(svd->V[0],&v);CHKERRQ(ierr);
|
ierr = VecDuplicate(svd->V[0],&v);CHKERRQ(ierr);
|
ierr = VecDuplicate(svd->V[0],&wv);CHKERRQ(ierr);
|
ierr = VecDuplicate(svd->V[0],&wv);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(Vec)*svd->n,&workV);CHKERRQ(ierr);
|
ierr = PetscMalloc(sizeof(Vec)*svd->n,&workV);CHKERRQ(ierr);
|
if (lanczos->oneside) {
|
if (lanczos->oneside) {
|