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 <[email protected]> > _________________________________________________________________ > > 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) > _________________________________________________________________
