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: }