[petsc-users] question about arbitrary eigenvector selection in SLEPC

[Jose E. Roman]

El 29/05/2014, a las 17:03, Giacomo Mulas escribió:

> Hi Jose, and list.
> 
> I am in the process of writing the function to use with
> EPSSetArbitrarySelection().
> 
> Inside it, I will need to take some given component (which one is included
> in the info passed via ctx) of the eigenvector and square it. To do this,
> since the eigenvector is not necessarily local, I will need to first do a
> scatter to a local 1-component vector. So this would be like:
> 
> ... some omitted machinery to cast the info from *ctx to more easily
> accessible form...
> 
>  ierr = ISCreateStride(PETXC_COMM_WORLD,1,myindex,1,&is1_from);CHKERRQ(ierr);
>  ierr = ISCreateStride(PETSC_COMM_WORLD,1,0,1,&is1_to);CHKERRQ(ierr);
>  ierr = VecCreateSeq(PETSC_COMM_SELF, 1, &localx1);CHKERRQ(ierr);
>  ierr = VecScatterCreate(xr,is1_from,localx1,is1_to,&scatter1); CHKERRQ(ierr);
>  ierr = VecScatterBegin(scatter1,xr,localx1,INSERT_VALUES,
> 			 SCATTER_FORWARD);
>  ierr = VecScatterEnd(scatter1,xr,localx1,INSERT_VALUES,
> 		       SCATTER_FORWARD);
>  ierr = VecGetArray(localx1,&comp);
>  *rr = comp*comp;
>  ierr = VecRestoreArray(localx1, &comp);
>  ierr = VecDestroy(localx1);
>  ierr = VecScatterDestroy(&scatter1);
>  ierr = ISDestroy(&is1_from);
>  ierr = ISDestroy(&is1_to);
>  *ri = 0;
> 
> ... some internal housekeeping omitted
> 
>  return 0;
> 
> The questions are:
> 
> 1) when the arbitrary function is called, is it called on all nodes
> simultaneously, so that collective functions can be expected to work
> properly, being called on all involved nodes at the same time? Should all
> processes compute the *rr and *ri to be returned, and return the same value?
> would it be more efficient to create a unit vector uv containing only one
> nonzero component, and then use VecDot(xr, uv, &comp), instead of pulling
> the component I need and squaring it as I did above?
> 

Yes, all processes must have the same values.
Use the code snippet proposed by Matt.

> 
> 2) since the stride, the 1-component vector, the scatter are presumably the same through all calls within one EPSSolve, can I take them out of the arbitrary function, and make them available to it through *ctx? For this
> to work, the structure of xr, the eigenvector passed to the arbitrary
> function, must be known outside of EPSSolve.

Yes.
Internally, all vectors are basically cloned from a template vector created with MatGetVecs(A,xr,NULL) so you can do the same outside EPSSolve() to determine local sizes.

Jose


> 
> 
> Thanks, bye
> Giacomo
> 
> -- 
> _________________________________________________________________
> 
> Giacomo Mulas <gmulas at oa-cagliari.inaf.it>
> _________________________________________________________________
> 
> INAF - Osservatorio Astronomico di Cagliari
> via della scienza 5 - 09047 Selargius (CA)
> 
> tel.   +39 070 71180244
> mob. : +39 329  6603810
> _________________________________________________________________
> 
> "When the storms are raging around you, stay right where you are"
>                         (Freddy Mercury)
> _________________________________________________________________