/*
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
SLEPc - Scalable Library for Eigenvalue Problem Computations
Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
This file is part of SLEPc.
SLEPc is free software: you can redistribute it and/or modify it under the
terms of version 3 of the GNU Lesser General Public License as published by
the Free Software Foundation.
SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
more details.
You should have received a copy of the GNU Lesser General Public License
along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
*/
static char help[] = "Tests matdense interface.\n\n";
#include <slepcmatdense.h>
#undef __FUNCT__
#define __FUNCT__ "main"
int main( int argc, char **argv )
{
MatDense A,B,C,D; /* matrices */
PetscScalar *v;
PetscInt n=45,m,i,j,its,M,reps=1;
PetscBool flag;
Vec *vecs0,*vecs1,t;
PetscLogStage saxpy,smult,supdate;
PetscErrorCode ierr;
SlepcInitialize(&argc,&argv,(char*)0,help);
ierr = PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);CHKERRQ(ierr);
ierr = PetscOptionsGetInt(PETSC_NULL,"-m",&m,&flag);CHKERRQ(ierr);
if(!flag) m=n;
M = PetscMax(m,n);
ierr = PetscOptionsGetInt(PETSC_NULL,"-r",&reps,&flag);CHKERRQ(ierr);
ierr = PetscPrintf(PETSC_COMM_WORLD,"\nDense matrices of size %Dx%D\n\n",m,n);CHKERRQ(ierr);
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
Create matrices
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
ierr = MatDenseCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
ierr = MatDenseSetMaxSizes(A,m,n);CHKERRQ(ierr);
ierr = MatDenseSetSizes(A,0,0,m,n);CHKERRQ(ierr);
ierr = MatDenseSetFromOptions(A);CHKERRQ(ierr);
ierr = MatDenseSetUpPreallocation(A);CHKERRQ(ierr);
ierr = MatDenseDuplicate(A,MATDENSE_DO_NOT_COPY_VALUES,&B);CHKERRQ(ierr);
ierr = MatDenseSetFromOptions(B);CHKERRQ(ierr);
ierr = MatDenseSetUpPreallocation(B);CHKERRQ(ierr);
ierr = MatDenseCreate(PETSC_COMM_WORLD,&C);CHKERRQ(ierr);
ierr = MatDenseSetMaxSizes(C,M,M);CHKERRQ(ierr);
ierr = MatDenseSetFromOptions(C);CHKERRQ(ierr);
ierr = MatDenseSetUpPreallocation(C);CHKERRQ(ierr);
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
Initialize matrices
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
ierr = MatDenseGetArray(A,&v);CHKERRQ(ierr);
for (j=0; j<n; j++) {
for (i=0; i<m; i++) {
v[m*j+i] = 1.0;
}
}
ierr = MatDenseRestoreArray(A,&v);CHKERRQ(ierr);
ierr = MatDenseGetArray(B,&v);CHKERRQ(ierr);
for (j=0; j<n; j++) {
for (i=0; i<m; i++) {
v[m*j+i] = i & 1 ? 1.0 : -1.0;
}
}
ierr = MatDenseRestoreArray(B,&v);CHKERRQ(ierr);
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
Initialize vectors
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
ierr = VecCreate(PETSC_COMM_WORLD,&t);CHKERRQ(ierr);
ierr = VecSetSizes(t,m,PETSC_DECIDE);CHKERRQ(ierr);
ierr = VecSetFromOptions(t);CHKERRQ(ierr);
ierr = SlepcVecSetTemplate(t);CHKERRQ(ierr);
ierr = VecDuplicateVecs(t,n,&vecs0);CHKERRQ(ierr);
ierr = VecDuplicateVecs(t,n,&vecs1);CHKERRQ(ierr);
ierr = VecDestroy(&t);CHKERRQ(ierr);
for (i=0;i<n;i++) { ierr = VecSetRandom(vecs0[i],PETSC_NULL);CHKERRQ(ierr); }
for (i=0;i<n;i++) { ierr = VecSetRandom(vecs1[i],PETSC_NULL);CHKERRQ(ierr); }
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
Test operations
- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
ierr = PetscLogStageRegister("AXPY",&saxpy);CHKERRQ(ierr);
ierr = PetscLogStagePush(saxpy);CHKERRQ(ierr);
for (i=0; i<reps; i++) { ierr = MatDenseAXPY(B,1.0,A);CHKERRQ(ierr); }
for (i=0; i<reps; i++) {
for (j=0; j<n; j++) { ierr = VecAXPY(vecs0[j],1.0,vecs1[j]);CHKERRQ(ierr); }
}
PetscLogStagePop();
ierr = PetscLogStageRegister("MatMult",&smult);CHKERRQ(ierr);
ierr = PetscLogStagePush(smult);CHKERRQ(ierr);
ierr = MatDenseSetSizes(C,0,0,n,n);CHKERRQ(ierr);
for (i=0; i<reps; i++) { ierr = MatDenseMatMult(C,1.0,0.0,A,PETSC_TRUE,B,PETSC_FALSE);CHKERRQ(ierr); }
ierr = MatDenseGetArray(C,&v);CHKERRQ(ierr);
for (i=0; i<reps; i++) {
for (j=0; j<n; j++) { ierr = VecMDot(vecs0[j],n,vecs1,&v[n*j]);CHKERRQ(ierr); }
}
ierr = MatDenseRestoreArray(C,&v);CHKERRQ(ierr);
PetscLogStagePop();
ierr = PetscLogStageRegister("Update",&supdate);CHKERRQ(ierr);
ierr = PetscLogStagePush(supdate);CHKERRQ(ierr);
ierr = MatDenseSetSizes(A,0,0,m,n);CHKERRQ(ierr);
ierr = MatDenseDuplicate(A,MATDENSE_MAKE_SIBLING,&D);CHKERRQ(ierr);
ierr = MatDenseSetSizes(D,0,0,m,n);CHKERRQ(ierr);
ierr = MatDenseSetSizes(C,0,0,n,n);CHKERRQ(ierr);
for (i=0; i<reps; i++) { ierr = MatDenseMatMult(D,1.0,0.0,A,PETSC_FALSE,C,PETSC_FALSE);CHKERRQ(ierr); }
ierr = MatDenseDestroy(&D);CHKERRQ(ierr);
ierr = MatDenseGetArray(C,&v);CHKERRQ(ierr);
for (i=0; i<reps; i++) {
ierr = SlepcUpdateVectors(n,vecs0,0,n,v,m,PETSC_FALSE);CHKERRQ(ierr);
}
ierr = MatDenseRestoreArray(C,&v);CHKERRQ(ierr);
PetscLogStagePop();
ierr = MatDenseDestroy(&A);CHKERRQ(ierr);
ierr = MatDenseDestroy(&B);CHKERRQ(ierr);
ierr = MatDenseDestroy(&C);CHKERRQ(ierr);
ierr = VecDestroyVecs(n,&vecs0);CHKERRQ(ierr);
ierr = VecDestroyVecs(n,&vecs1);CHKERRQ(ierr);
ierr = SlepcFinalize();CHKERRQ(ierr);
return 0;
}