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:   PetscCallExternal(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:   { \
 60:     hypre_ParVector *par_vector   = (hypre_ParVector *)hypre_IJVectorObject(((hypre_IJVector *)b)); \
 61:     hypre_Vector    *local_vector = hypre_ParVectorLocalVector(par_vector); \
 62:     savedvalue                    = local_vector->data; \
 63:     local_vector->data            = newvalue; \
 64:   }

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

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

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

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

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

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

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

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

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

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

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

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

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