Actual source code: vhyp.c


  2: /*
  3:     Creates hypre ijvector from PETSc vector
  4: */

  6: #include <petsc/private/vecimpl.h>
  7: #include <../src/vec/vec/impls/hypre/vhyp.h>
  8: #include <HYPRE.h>

 10: PetscErrorCode VecHYPRE_IJVectorCreate(PetscLayout map,VecHYPRE_IJVector *ij)
 11: {
 12:   VecHYPRE_IJVector nij;

 14:   PetscNew(&nij);
 15:   PetscLayoutSetUp(map);
 16:   HYPRE_IJVectorCreate(map->comm,map->rstart,map->rend-1,&nij->ij);
 17:   HYPRE_IJVectorSetObjectType(nij->ij,HYPRE_PARCSR);
 18: #if defined(PETSC_HAVE_HYPRE_DEVICE)
 19:   HYPRE_IJVectorInitialize_v2(nij->ij,HYPRE_MEMORY_DEVICE);
 20: #else
 21:   HYPRE_IJVectorInitialize(nij->ij);
 22: #endif
 23:   HYPRE_IJVectorAssemble(nij->ij);
 24:   *ij  = nij;
 25:   return 0;
 26: }

 28: PetscErrorCode VecHYPRE_IJVectorDestroy(VecHYPRE_IJVector *ij)
 29: {
 30:   if (!*ij) return 0;
 32:   PetscStackCallStandard(HYPRE_IJVectorDestroy,(*ij)->ij);
 33:   PetscFree(*ij);
 34:   return 0;
 35: }

 37: PetscErrorCode VecHYPRE_IJVectorCopy(Vec v,VecHYPRE_IJVector ij)
 38: {
 39:   const PetscScalar *array;

 41: #if defined(PETSC_HAVE_HYPRE_DEVICE)
 42:   HYPRE_IJVectorInitialize_v2(ij->ij,HYPRE_MEMORY_DEVICE);
 43: #else
 44:   HYPRE_IJVectorInitialize(ij->ij);
 45: #endif
 46:   VecGetArrayRead(v,&array);
 47:   HYPRE_IJVectorSetValues(ij->ij,v->map->n,NULL,(HYPRE_Complex*)array);
 48:   VecRestoreArrayRead(v,&array);
 49:   HYPRE_IJVectorAssemble(ij->ij);
 50:   return 0;
 51: }

 53: /*
 54:     Replaces the address where the HYPRE vector points to its data with the address of
 55:   PETSc's data. Saves the old address so it can be reset when we are finished with it.
 56:   Allows use to get the data into a HYPRE vector without the cost of memcopies
 57: */
 58: #define VecHYPRE_ParVectorReplacePointer(b,newvalue,savedvalue) {                               \
 59:   hypre_ParVector *par_vector   = (hypre_ParVector*)hypre_IJVectorObject(((hypre_IJVector*)b)); \
 60:   hypre_Vector    *local_vector = hypre_ParVectorLocalVector(par_vector);                       \
 61:   savedvalue         = local_vector->data;                                               \
 62:   local_vector->data = newvalue;                                                         \
 63: }

 65: /*
 66:   This routine access the pointer to the raw data of the "v" to be passed to HYPRE
 67:    - rw values indicate the type of access : 0 -> read, 1 -> write, 2 -> read-write
 68:    - hmem is the location HYPRE is expecting
 69:    - the function returns a pointer to the data (ptr) and the corresponding restore
 70:   Could be extended to VECKOKKOS if we had a way to access the raw pointer to device data.
 71: */
 72: static inline PetscErrorCode VecGetArrayForHYPRE(Vec v, int rw, HYPRE_MemoryLocation hmem, PetscScalar **ptr, PetscErrorCode(**res)(Vec,PetscScalar**))
 73: {
 74:   PetscMemType   mtype;
 75:   MPI_Comm       comm;

 77: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
 78:   hmem = HYPRE_MEMORY_HOST; /* this is just a convenience because HYPRE_MEMORY_HOST and HYPRE_MEMORY_DEVICE are the same in this case */
 79: #endif
 80:   *ptr = NULL;
 81:   *res = NULL;
 82:   PetscObjectGetComm((PetscObject)v,&comm);
 83:   switch (rw) {
 84:   case 0: /* read */
 85:     if (hmem == HYPRE_MEMORY_HOST) {
 86:       VecGetArrayRead(v,(const PetscScalar**)ptr);
 87:       *res = (PetscErrorCode(*)(Vec,PetscScalar**))VecRestoreArrayRead;
 88:     } else {
 89:       VecGetArrayReadAndMemType(v,(const PetscScalar**)ptr,&mtype);
 91:       *res = (PetscErrorCode(*)(Vec,PetscScalar**))VecRestoreArrayReadAndMemType;
 92:     }
 93:     break;
 94:   case 1: /* write */
 95:     if (hmem == HYPRE_MEMORY_HOST) {
 96:       VecGetArrayWrite(v,ptr);
 97:       *res = VecRestoreArrayWrite;
 98:     } else {
 99:       VecGetArrayWriteAndMemType(v,(PetscScalar**)ptr,&mtype);
101:       *res = VecRestoreArrayWriteAndMemType;
102:     }
103:     break;
104:   case 2: /* read/write */
105:     if (hmem == HYPRE_MEMORY_HOST) {
106:       VecGetArray(v,ptr);
107:       *res = VecRestoreArray;
108:     } else {
109:       VecGetArrayAndMemType(v,(PetscScalar**)ptr,&mtype);
111:       *res = VecRestoreArrayAndMemType;
112:     }
113:     break;
114:   default:
115:     SETERRQ(comm,PETSC_ERR_SUP,"Unhandled case %d",rw);
116:   }
117:   return 0;
118: }

120: #define VecHYPRE_IJVectorMemoryLocation(v) hypre_IJVectorMemoryLocation((hypre_IJVector*)(v))

122: /* Temporarily pushes the array of the data in v to ij (read access)
123:    depending on the value of the ij memory location
124:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
125: PetscErrorCode VecHYPRE_IJVectorPushVecRead(VecHYPRE_IJVector ij, Vec v)
126: {
127:   HYPRE_Complex  *pv;

132:   VecGetArrayForHYPRE(v,0,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
133:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
134:   ij->pvec = v;
135:   return 0;
136: }

138: /* Temporarily pushes the array of the data in v to ij (write access)
139:    depending on the value of the ij memory location
140:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
141: PetscErrorCode VecHYPRE_IJVectorPushVecWrite(VecHYPRE_IJVector ij, Vec v)
142: {
143:   HYPRE_Complex  *pv;

148:   VecGetArrayForHYPRE(v,1,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
149:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
150:   ij->pvec = v;
151:   return 0;
152: }

154: /* Temporarily pushes the array of the data in v to ij (read/write access)
155:    depending on the value of the ij memory location
156:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
157: PetscErrorCode VecHYPRE_IJVectorPushVec(VecHYPRE_IJVector ij, Vec v)
158: {
159:   HYPRE_Complex  *pv;

164:   VecGetArrayForHYPRE(v,2,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
165:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
166:   ij->pvec = v;
167:   return 0;
168: }

170: /* Restores the pointer data to v */
171: PetscErrorCode VecHYPRE_IJVectorPopVec(VecHYPRE_IJVector ij)
172: {
173:   HYPRE_Complex  *pv;

177:   VecHYPRE_ParVectorReplacePointer(ij->ij,ij->hv,pv);
178:   (*ij->restore)(ij->pvec,(PetscScalar**)&pv);
179:   ij->hv = NULL;
180:   ij->pvec = NULL;
181:   ij->restore = NULL;
182:   return 0;
183: }

185: PetscErrorCode VecHYPRE_IJBindToCPU(VecHYPRE_IJVector ij,PetscBool bind)
186: {
187:   HYPRE_MemoryLocation hmem = bind ? HYPRE_MEMORY_HOST : HYPRE_MEMORY_DEVICE;
188:   hypre_ParVector      *hij;

192: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
193:   hmem = HYPRE_MEMORY_HOST;
194: #endif
195: #if PETSC_PKG_HYPRE_VERSION_GT(2,19,0)
196:   if (hmem != VecHYPRE_IJVectorMemoryLocation(ij->ij)) {
197:     PetscStackCallStandard(HYPRE_IJVectorGetObject,ij->ij,(void**)&hij);
198:     PetscStackCallStandard(hypre_ParVectorMigrate,hij,hmem);
199:   }
200: #endif
201:   return 0;
202: }