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