On Wed 2009-04-08 08:48, Kirk, Benjamin (JSC-EG311) wrote: > As for the performance implications, I think we can get nearly all of it > back by optimizing the underlying PetscVector (but actually do it in the > NumericVector base class) to bypass the PETSc API whenever possible. > > Specifically, we would use the VecCreateGhostWithArray (sp?) method to > provide the storage buffer and then access it directly whenever possible. I
No need to manage the memory yourself unless you really already have the
array. Usually you would just call VecGetArray on the local form (when
you want ghosted values, on the global form when you don't, although
you'll get a pointer to the same memory either way). Note that
VecGetArray is inlined for PETSc native vectors
#define VecGetArray(x,a) ((x)->petscnative ? (*(a) = *((PetscScalar
**)(x)->data),0) : VecGetArray_Private((x),(a)))
so it is very cheap.
I would be interested to see profiling output from Tim's code with
ghosted vectors.
Also, note that there is a bit of insanity in petsc_vector.C, consider
this from PetscVector<T>::add(const T)
for (int i=0; i<n; i++)
{
ierr = VecGetArray (loc_vec, &values);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
PetscScalar value = (values[i] + v);
ierr = VecRestoreArray (loc_vec, &values);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
ierr = VecSetValues (loc_vec, 1, &i, &value, INSERT_VALUES);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
}
which can be changed to
ierr = VecGetArray (loc_vec, &values);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
for (int i=0; i<n; i++) values[i] += v;
ierr = VecRestoreArray (loc_vec, &values);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
or
ierr = VecShift(loc_vec, v);
CHKERRABORT(libMesh::COMM_WORLD,ierr);
Either way, you eliminate a function call with bounds checking for every
entry. The non-ghosted implementation of this method is equally goofy.
Note that you only need to work with the local form when you need to
manipulate the ghosted entries. Many operations don't require this in
which case it's natural to work with the global form and call
VecGhostUpdateBegin/End when you need the ghosted values to be current.
Jed
pgpacx3tGoZKl.pgp
Description: PGP signature
------------------------------------------------------------------------------ This SF.net email is sponsored by: High Quality Requirements in a Collaborative Environment. Download a free trial of Rational Requirements Composer Now! http://p.sf.net/sfu/www-ibm-com
_______________________________________________ Libmesh-devel mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/libmesh-devel
