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:   PetscErrorCode    ierr;
 13:   VecHYPRE_IJVector nij;

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

 30: PetscErrorCode VecHYPRE_IJVectorDestroy(VecHYPRE_IJVector *ij)
 31: {

 35:   if (!*ij) return(0);
 36:   if ((*ij)->pvec) SETERRQ(PetscObjectComm((PetscObject)((*ij)->pvec)),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
 37:   PetscStackCallStandard(HYPRE_IJVectorDestroy,((*ij)->ij));
 38:   PetscFree(*ij);
 39:   return(0);
 40: }

 42: PetscErrorCode VecHYPRE_IJVectorCopy(Vec v,VecHYPRE_IJVector ij)
 43: {
 44:   PetscErrorCode    ierr;
 45:   const PetscScalar *array;

 48: #if defined(PETSC_HAVE_HYPRE_DEVICE)
 49:   HYPRE_IJVectorInitialize_v2(ij->ij,HYPRE_MEMORY_DEVICE);
 50: #else
 51:   HYPRE_IJVectorInitialize(ij->ij);
 52: #endif
 53:   VecGetArrayRead(v,&array);
 54:   HYPRE_IJVectorSetValues(ij->ij,v->map->n,NULL,(HYPRE_Complex*)array);
 55:   VecRestoreArrayRead(v,&array);
 56:   HYPRE_IJVectorAssemble(ij->ij);
 57:   return(0);
 58: }

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

 72: /*
 73:   This routine access the pointer to the raw data of the "v" to be passed to HYPRE
 74:    - rw values indicate the type of access : 0 -> read, 1 -> write, 2 -> read-write
 75:    - hmem is the location HYPRE is expecting
 76:    - the function returns a pointer to the data (ptr) and the corresponding restore
 77:   Could be extended to VECKOKKOS if we had a way to access the raw pointer to device data.
 78: */
 79: PETSC_STATIC_INLINE PetscErrorCode VecGetArrayForHYPRE(Vec v, int rw, HYPRE_MemoryLocation hmem, PetscScalar **ptr, PetscErrorCode(**res)(Vec,PetscScalar**))
 80: {
 81:   PetscBool      usehip = PETSC_FALSE,usecuda = PETSC_FALSE;

 85: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
 86:   hmem = HYPRE_MEMORY_HOST; /* this is just a convenience because HYPRE_MEMORY_HOST and HYPRE_MEMORY_DEVICE are the same in this case */
 87: #else
 88: #if defined(HYPRE_USING_HIP)
 89:   usehip = PETSC_TRUE;
 90: #elif defined(HYPRE_USING_CUDA)
 91:   usecuda = PETSC_TRUE;
 92: #else
 93: #error Not yet coded!
 94: #endif
 95: #endif
 96:   *ptr = NULL;
 97:   *res = NULL;
 98:   switch (rw) {
 99:   case 0: /* read */
100:     if (hmem == HYPRE_MEMORY_HOST) {
101:       VecGetArrayRead(v,(const PetscScalar**)ptr);
102:       *res = (PetscErrorCode(*)(Vec,PetscScalar**))VecRestoreArrayRead;
103:     } else if (usehip) {
104:       VecHIPGetArrayRead(v,(const PetscScalar**)ptr);
105:       *res = (PetscErrorCode(*)(Vec,PetscScalar**))VecHIPRestoreArrayRead;
106:     } else if (usecuda) {
107:       VecCUDAGetArrayRead(v,(const PetscScalar**)ptr);
108:       *res = (PetscErrorCode(*)(Vec,PetscScalar**))VecCUDARestoreArrayRead;
109:     }
110:     break;
111:   case 1: /* write */
112:     if (hmem == HYPRE_MEMORY_HOST) {
113:       VecGetArrayWrite(v,ptr);
114:       *res = VecRestoreArrayWrite;
115:     } else if (usehip) {
116:       VecHIPGetArrayWrite(v,(PetscScalar**)ptr);
117:       *res = VecHIPRestoreArrayWrite;
118:     } else if (usecuda) {
119:       VecCUDAGetArrayWrite(v,(PetscScalar**)ptr);
120:       *res = VecCUDARestoreArrayWrite;
121:     }
122:     break;
123:   case 2: /* read/write */
124:     if (hmem == HYPRE_MEMORY_HOST) {
125:       VecGetArray(v,ptr);
126:       *res = VecRestoreArray;
127:     } else if (usehip) {
128:       VecHIPGetArray(v,(PetscScalar**)ptr);
129:       *res = VecHIPRestoreArray;
130:     } else if (usecuda) {
131:       VecCUDAGetArray(v,(PetscScalar**)ptr);
132:       *res = VecCUDARestoreArray;
133:     }
134:     break;
135:   default:
136:     SETERRQ1(PetscObjectComm((PetscObject)v),PETSC_ERR_SUP,"Unhandled case %d",rw);
137:   }
138:   return(0);
139: }

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

143: /* Temporarily pushes the array of the data in v to ij (read access)
144:    depending on the value of the ij memory location
145:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
146: PetscErrorCode VecHYPRE_IJVectorPushVecRead(VecHYPRE_IJVector ij, Vec v)
147: {
148:   HYPRE_Complex  *pv;

153:   if (ij->pvec) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
154:   if (ij->hv) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
155:   VecGetArrayForHYPRE(v,0,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
156:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
157:   ij->pvec = v;
158:   return(0);
159: }

161: /* Temporarily pushes the array of the data in v to ij (write access)
162:    depending on the value of the ij memory location
163:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
164: PetscErrorCode VecHYPRE_IJVectorPushVecWrite(VecHYPRE_IJVector ij, Vec v)
165: {
166:   HYPRE_Complex  *pv;

171:   if (ij->pvec) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
172:   if (ij->hv) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
173:   VecGetArrayForHYPRE(v,1,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
174:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
175:   ij->pvec = v;
176:   return(0);
177: }

179: /* Temporarily pushes the array of the data in v to ij (read/write access)
180:    depending on the value of the ij memory location
181:    Must be completed with a call to VecHYPRE_IJVectorPopVec */
182: PetscErrorCode VecHYPRE_IJVectorPushVec(VecHYPRE_IJVector ij, Vec v)
183: {
184:   HYPRE_Complex  *pv;

189:   if (ij->pvec) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
190:   if (ij->hv) SETERRQ(PetscObjectComm((PetscObject)v),PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
191:   VecGetArrayForHYPRE(v,2,VecHYPRE_IJVectorMemoryLocation(ij->ij),(PetscScalar**)&pv,&ij->restore);
192:   VecHYPRE_ParVectorReplacePointer(ij->ij,pv,ij->hv);
193:   ij->pvec = v;
194:   return(0);
195: }

197: /* Restores the pointer data to v */
198: PetscErrorCode VecHYPRE_IJVectorPopVec(VecHYPRE_IJVector ij)
199: {
200:   HYPRE_Complex  *pv;

204:   if (!ij->pvec) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPushVec()");
205:   if (!ij->restore) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPushVec()");
206:   VecHYPRE_ParVectorReplacePointer(ij->ij,ij->hv,pv);
207:   (*ij->restore)(ij->pvec,(PetscScalar**)&pv);
208:   ij->hv = NULL;
209:   ij->pvec = NULL;
210:   ij->restore = NULL;
211:   return(0);
212: }

214: PetscErrorCode VecHYPRE_IJBindToCPU(VecHYPRE_IJVector ij,PetscBool bind)
215: {
216:   HYPRE_MemoryLocation hmem = bind ? HYPRE_MEMORY_HOST : HYPRE_MEMORY_DEVICE;
217:   hypre_ParVector      *hij;

220:   if (ij->pvec) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
221:   if (ij->hv) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Forgot to call VecHYPRE_IJVectorPopVec()");
222: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
223:   hmem = HYPRE_MEMORY_HOST;
224: #endif
225: #if PETSC_PKG_HYPRE_VERSION_GT(2,19,0)
226:   if (hmem != VecHYPRE_IJVectorMemoryLocation(ij->ij)) {
227:     PetscStackCallStandard(HYPRE_IJVectorGetObject,(ij->ij,(void**)&hij));
228:     PetscStackCallStandard(hypre_ParVectorMigrate,(hij,hmem));
229:   }
230: #endif
231:   return(0);
232: }