I am trying to pass an array of Vec's in PETSc to a function, modify it
internally and retrieve the results. The idea is to copy a handful of
eigenvectors from the EPS solver to the main routine. A pseudocode is as
follows:
#include <slepceps.h>
PetscErrorCode foo(Vec **y, int n) {
EPS eps; // eigenproblem solver context
// ...
ierr = MatCreateVecs(H, &x, NULL); CHKERRQ(ierr);
ierr = EPSSolve(eps); CHKERRQ(ierr);
// ...
ierr = VecDuplicateVecs(x, n, y); CHKERRQ(ierr);
for (int i = 0; i < n ; i++) { // I can guarantee that n < nconv
ierr = EPSGetEigenvector(eps, i, *y[i], NULL);
CHKERRQ(ierr); // this breaks for i = 1
ierr = VecNorm(*y[i],NORM_2,&norm); CHKERRQ(ierr); // this
prints out fine for i = 0 (norm = 1)
printf("norm = %f\n", norm);
}
ierr = EPSDestroy(&eps); CHKERRQ(ierr);
ierr = VecDestroy(&x); CHKERRQ(ierr);
return ierr;
}
int main(int argc,char **argv)
{
PetscErrorCode ierr;
PetscScalar norm;
Vec *y;
foo(&y, 3);
ierr = VecDestroyVecs(3, &y); CHKERRQ(ierr);
return 0;
}
Am I making a naive mistake here?
