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: PetscFunctionBegin;
15: PetscCall(PetscNew(&nij));
16: PetscCall(PetscLayoutSetUp(map));
17: PetscCallExternal(HYPRE_IJVectorCreate, map->comm, map->rstart, map->rend - 1, &nij->ij);
18: PetscCallExternal(HYPRE_IJVectorSetObjectType, nij->ij, HYPRE_PARCSR);
19: #if defined(PETSC_HAVE_HYPRE_DEVICE)
20: PetscCallExternal(HYPRE_IJVectorInitialize_v2, nij->ij, HYPRE_MEMORY_DEVICE);
21: #else
22: PetscCallExternal(HYPRE_IJVectorInitialize, nij->ij);
23: #endif
24: PetscCallExternal(HYPRE_IJVectorAssemble, nij->ij);
25: *ij = nij;
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: PetscErrorCode VecHYPRE_IJVectorDestroy(VecHYPRE_IJVector *ij)
30: {
31: PetscFunctionBegin;
32: if (!*ij) PetscFunctionReturn(PETSC_SUCCESS);
33: PetscCheck(!(*ij)->pvec, PetscObjectComm((PetscObject)((*ij)->pvec)), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
34: PetscCallExternal(HYPRE_IJVectorDestroy, (*ij)->ij);
35: PetscCall(PetscFree(*ij));
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: PetscErrorCode VecHYPRE_IJVectorCopy(Vec v, VecHYPRE_IJVector ij)
40: {
41: const PetscScalar *array;
43: PetscFunctionBegin;
44: #if defined(PETSC_HAVE_HYPRE_DEVICE)
45: PetscCallExternal(HYPRE_IJVectorInitialize_v2, ij->ij, HYPRE_MEMORY_DEVICE);
46: #else
47: PetscCallExternal(HYPRE_IJVectorInitialize, ij->ij);
48: #endif
49: PetscCall(VecGetArrayRead(v, &array));
50: PetscCallExternal(HYPRE_IJVectorSetValues, ij->ij, v->map->n, NULL, (HYPRE_Complex *)array);
51: PetscCall(VecRestoreArrayRead(v, &array));
52: PetscCallExternal(HYPRE_IJVectorAssemble, ij->ij);
53: PetscFunctionReturn(PETSC_SUCCESS);
54: }
56: /*
57: Replaces the address where the HYPRE vector points to its data with the address of
58: PETSc's data. Saves the old address so it can be reset when we are finished with it.
59: Allows use to get the data into a HYPRE vector without the cost of memcopies
60: */
61: #define VecHYPRE_ParVectorReplacePointer(b, newvalue, savedvalue) \
62: { \
63: hypre_ParVector *par_vector = (hypre_ParVector *)hypre_IJVectorObject(((hypre_IJVector *)b)); \
64: hypre_Vector *local_vector = hypre_ParVectorLocalVector(par_vector); \
65: savedvalue = local_vector->data; \
66: local_vector->data = newvalue; \
67: }
69: /*
70: This routine access the pointer to the raw data of the "v" to be passed to HYPRE
71: - rw values indicate the type of access : 0 -> read, 1 -> write, 2 -> read-write
72: - hmem is the location HYPRE is expecting
73: - the function returns a pointer to the data (ptr) and the corresponding restore
74: Could be extended to VECKOKKOS if we had a way to access the raw pointer to device data.
75: */
76: static inline PetscErrorCode VecGetArrayForHYPRE(Vec v, int rw, HYPRE_MemoryLocation hmem, PetscScalar **ptr, PetscErrorCode (**res)(Vec, PetscScalar **))
77: {
78: PetscMemType mtype;
79: MPI_Comm comm;
81: PetscFunctionBegin;
82: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
83: hmem = HYPRE_MEMORY_HOST; /* this is just a convenience because HYPRE_MEMORY_HOST and HYPRE_MEMORY_DEVICE are the same in this case */
84: #endif
85: *ptr = NULL;
86: *res = NULL;
87: PetscCall(PetscObjectGetComm((PetscObject)v, &comm));
88: switch (rw) {
89: case 0: /* read */
90: if (hmem == HYPRE_MEMORY_HOST) {
91: PetscCall(VecGetArrayRead(v, (const PetscScalar **)ptr));
92: *res = (PetscErrorCode(*)(Vec, PetscScalar **))VecRestoreArrayRead;
93: } else {
94: PetscCall(VecGetArrayReadAndMemType(v, (const PetscScalar **)ptr, &mtype));
95: PetscCheck(PetscMemTypeDevice(mtype), comm, PETSC_ERR_ARG_WRONG, "HYPRE_MEMORY_DEVICE expects a device vector. You need to enable PETSc device support, for example, in some cases, -vec_type cuda");
96: *res = (PetscErrorCode(*)(Vec, PetscScalar **))VecRestoreArrayReadAndMemType;
97: }
98: break;
99: case 1: /* write */
100: if (hmem == HYPRE_MEMORY_HOST) {
101: PetscCall(VecGetArrayWrite(v, ptr));
102: *res = VecRestoreArrayWrite;
103: } else {
104: PetscCall(VecGetArrayWriteAndMemType(v, (PetscScalar **)ptr, &mtype));
105: PetscCheck(PetscMemTypeDevice(mtype), comm, PETSC_ERR_ARG_WRONG, "HYPRE_MEMORY_DEVICE expects a device vector. You need to enable PETSc device support, for example, in some cases, -vec_type cuda");
106: *res = VecRestoreArrayWriteAndMemType;
107: }
108: break;
109: case 2: /* read/write */
110: if (hmem == HYPRE_MEMORY_HOST) {
111: PetscCall(VecGetArray(v, ptr));
112: *res = VecRestoreArray;
113: } else {
114: PetscCall(VecGetArrayAndMemType(v, (PetscScalar **)ptr, &mtype));
115: PetscCheck(PetscMemTypeDevice(mtype), comm, PETSC_ERR_ARG_WRONG, "HYPRE_MEMORY_DEVICE expects a device vector. You need to enable PETSc device support, for example, in some cases, -vec_type cuda");
116: *res = VecRestoreArrayAndMemType;
117: }
118: break;
119: default:
120: SETERRQ(comm, PETSC_ERR_SUP, "Unhandled case %d", rw);
121: }
122: PetscFunctionReturn(PETSC_SUCCESS);
123: }
125: #define VecHYPRE_IJVectorMemoryLocation(v) hypre_IJVectorMemoryLocation((hypre_IJVector *)(v))
127: /* Temporarily pushes the array of the data in v to ij (read access)
128: depending on the value of the ij memory location
129: Must be completed with a call to VecHYPRE_IJVectorPopVec */
130: PetscErrorCode VecHYPRE_IJVectorPushVecRead(VecHYPRE_IJVector ij, Vec v)
131: {
132: HYPRE_Complex *pv;
134: PetscFunctionBegin;
136: PetscCheck(!ij->pvec, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
137: PetscCheck(!ij->hv, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
138: PetscCall(VecGetArrayForHYPRE(v, 0, VecHYPRE_IJVectorMemoryLocation(ij->ij), (PetscScalar **)&pv, &ij->restore));
139: VecHYPRE_ParVectorReplacePointer(ij->ij, pv, ij->hv);
140: ij->pvec = v;
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: /* Temporarily pushes the array of the data in v to ij (write access)
145: depending on the value of the ij memory location
146: Must be completed with a call to VecHYPRE_IJVectorPopVec */
147: PetscErrorCode VecHYPRE_IJVectorPushVecWrite(VecHYPRE_IJVector ij, Vec v)
148: {
149: HYPRE_Complex *pv;
151: PetscFunctionBegin;
153: PetscCheck(!ij->pvec, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
154: PetscCheck(!ij->hv, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
155: PetscCall(VecGetArrayForHYPRE(v, 1, VecHYPRE_IJVectorMemoryLocation(ij->ij), (PetscScalar **)&pv, &ij->restore));
156: VecHYPRE_ParVectorReplacePointer(ij->ij, pv, ij->hv);
157: ij->pvec = v;
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: /* Temporarily pushes the array of the data in v to ij (read/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_IJVectorPushVec(VecHYPRE_IJVector ij, Vec v)
165: {
166: HYPRE_Complex *pv;
168: PetscFunctionBegin;
170: PetscCheck(!ij->pvec, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
171: PetscCheck(!ij->hv, PetscObjectComm((PetscObject)v), PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
172: PetscCall(VecGetArrayForHYPRE(v, 2, VecHYPRE_IJVectorMemoryLocation(ij->ij), (PetscScalar **)&pv, &ij->restore));
173: VecHYPRE_ParVectorReplacePointer(ij->ij, pv, ij->hv);
174: ij->pvec = v;
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }
178: /* Restores the pointer data to v */
179: PetscErrorCode VecHYPRE_IJVectorPopVec(VecHYPRE_IJVector ij)
180: {
181: HYPRE_Complex *pv;
183: PetscFunctionBegin;
184: PetscCheck(ij->pvec, PETSC_COMM_SELF, PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPushVec()");
185: PetscCheck(ij->restore, PETSC_COMM_SELF, PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPushVec()");
186: VecHYPRE_ParVectorReplacePointer(ij->ij, ij->hv, pv);
187: PetscCall((*ij->restore)(ij->pvec, (PetscScalar **)&pv));
188: ij->hv = NULL;
189: ij->pvec = NULL;
190: ij->restore = NULL;
191: PetscFunctionReturn(PETSC_SUCCESS);
192: }
194: PetscErrorCode VecHYPRE_IJBindToCPU(VecHYPRE_IJVector ij, PetscBool bind)
195: {
196: HYPRE_MemoryLocation hmem = bind ? HYPRE_MEMORY_HOST : HYPRE_MEMORY_DEVICE;
197: hypre_ParVector *hij;
199: PetscFunctionBegin;
200: PetscCheck(!ij->pvec, PETSC_COMM_SELF, PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
201: PetscCheck(!ij->hv, PETSC_COMM_SELF, PETSC_ERR_ORDER, "Forgot to call VecHYPRE_IJVectorPopVec()");
202: #if !defined(PETSC_HAVE_HYPRE_DEVICE)
203: hmem = HYPRE_MEMORY_HOST;
204: #endif
205: #if PETSC_PKG_HYPRE_VERSION_GT(2, 19, 0)
206: if (hmem != VecHYPRE_IJVectorMemoryLocation(ij->ij)) {
207: PetscCallExternal(HYPRE_IJVectorGetObject, ij->ij, (void **)&hij);
208: PetscCallExternal(hypre_ParVectorMigrate, hij, hmem);
209: }
210: #endif
211: PetscFunctionReturn(PETSC_SUCCESS);
212: }