1 #include <petsc/private/vecimpl_kokkos.hpp> 2 #include <petsc/private/dmdaimpl.h> 3 #include <petscdmda_kokkos.hpp> 4 5 /* SUBMANSEC = DMDA */ 6 7 /* Use macro instead of inlined function to avoid annoying warnings like: 'dof' may be used uninitialized in this function [-Wmaybe-uninitialized] */ 8 #define DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof) \ 9 do { \ 10 PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm)); \ 11 PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm)); \ 12 PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL)); \ 13 /* Handle case where user passes in global vector as opposed to local */ \ 14 PetscCall(VecGetLocalSize(vec, &N)); \ 15 if (N == xm * ym * zm * dof) { \ 16 gxm = xm; \ 17 gym = ym; \ 18 gzm = zm; \ 19 gxs = xs; \ 20 gys = ys; \ 21 gzs = zs; \ 22 } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof); \ 23 } while (0) 24 25 /* -------------------- 1D ---------------- */ 26 template <class MemorySpace> 27 PetscErrorCode DMDAVecGetKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView1DType<MemorySpace> *ov, PetscBool overwrite) { 28 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 29 PetscScalarKokkosViewType<MemorySpace> kv; 30 31 PetscFunctionBegin; 32 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 33 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 34 PetscValidPointer(ov, 3); 35 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 36 PetscCheck(dim == 1, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 1D but DMDA is %dD", (int)dim); 37 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 38 else PetscCall(VecGetKokkosView(vec, &kv)); 39 /* Construct the unmanaged OffsetView with {begin0,begin1,begins2},{end0,end1,end2} */ 40 *ov = PetscScalarKokkosOffsetView1DType<MemorySpace>(kv.data(), {gxs * dof}, {(gxs + gxm) * dof}); 41 PetscFunctionReturn(0); 42 } 43 44 template <class MemorySpace> 45 PetscErrorCode DMDAVecRestoreKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView1DType<MemorySpace> *ov, PetscBool overwrite) { 46 PetscScalarKokkosViewType<MemorySpace> kv; 47 48 PetscFunctionBegin; 49 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 50 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 51 PetscValidPointer(ov, 3); 52 kv = ov->view(); /* OffsetView to View */ 53 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 54 else PetscCall(VecRestoreKokkosView(vec, &kv)); 55 PetscFunctionReturn(0); 56 } 57 58 template <class MemorySpace> 59 PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView1DType<MemorySpace> *ov) { 60 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 61 ConstPetscScalarKokkosViewType<MemorySpace> kv; 62 63 PetscFunctionBegin; 64 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 65 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 66 PetscValidPointer(ov, 3); 67 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 68 PetscCheck(dim == 1, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 1D but DMDA is %dD", (int)dim); 69 PetscCall(VecGetKokkosView(vec, &kv)); 70 *ov = ConstPetscScalarKokkosOffsetView1DType<MemorySpace>(kv.data(), {gxs * dof}, {(gxs + gxm) * dof}); 71 PetscFunctionReturn(0); 72 } 73 74 template <class MemorySpace> 75 PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView1DType<MemorySpace> *ov) { 76 ConstPetscScalarKokkosViewType<MemorySpace> kv; 77 78 PetscFunctionBegin; 79 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 80 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 81 PetscValidPointer(ov, 3); 82 kv = ov->view(); 83 PetscCall(VecRestoreKokkosView(vec, &kv)); 84 PetscFunctionReturn(0); 85 } 86 87 /* ============================== 2D ================================= */ 88 template <class MemorySpace> 89 PetscErrorCode DMDAVecGetKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView2DType<MemorySpace> *ov, PetscBool overwrite) { 90 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 91 PetscScalarKokkosViewType<MemorySpace> kv; 92 93 PetscFunctionBegin; 94 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 95 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 96 PetscValidPointer(ov, 3); 97 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 98 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 2D but DMDA is %dD", (int)dim); 99 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 100 else PetscCall(VecGetKokkosView(vec, &kv)); 101 *ov = PetscScalarKokkosOffsetView2DType<MemorySpace>(kv.data(), {gys * dof, gxs * dof}, {(gys + gym) * dof, (gxs + gxm) * dof}); 102 PetscFunctionReturn(0); 103 } 104 105 template <class MemorySpace> 106 PetscErrorCode DMDAVecRestoreKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView2DType<MemorySpace> *ov, PetscBool overwrite) { 107 PetscScalarKokkosViewType<MemorySpace> kv; 108 109 PetscFunctionBegin; 110 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 111 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 112 PetscValidPointer(ov, 3); 113 // kv = ov->view(); /* 2D OffsetView => 2D View => 1D View. Why does it not work? */ 114 kv = PetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1)); 115 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 116 else PetscCall(VecRestoreKokkosView(vec, &kv)); 117 PetscFunctionReturn(0); 118 } 119 120 template <class MemorySpace> 121 PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView2DType<MemorySpace> *ov) { 122 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 123 ConstPetscScalarKokkosViewType<MemorySpace> kv; 124 125 PetscFunctionBegin; 126 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 127 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 128 PetscValidPointer(ov, 3); 129 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 130 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 2D but DMDA is %dD", (int)dim); 131 PetscCall(VecGetKokkosView(vec, &kv)); 132 *ov = ConstPetscScalarKokkosOffsetView2DType<MemorySpace>(kv.data(), {gys * dof, gxs * dof}, {(gys + gym) * dof, (gxs + gxm) * dof}); 133 PetscFunctionReturn(0); 134 } 135 136 template <class MemorySpace> 137 PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView2DType<MemorySpace> *ov) { 138 ConstPetscScalarKokkosViewType<MemorySpace> kv; 139 140 PetscFunctionBegin; 141 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 142 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 143 PetscValidPointer(ov, 3); 144 kv = ConstPetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1)); 145 PetscCall(VecRestoreKokkosView(vec, &kv)); 146 PetscFunctionReturn(0); 147 } 148 149 /* ============================== 3D ================================= */ 150 template <class MemorySpace> 151 PetscErrorCode DMDAVecGetKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView3DType<MemorySpace> *ov, PetscBool overwrite) { 152 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 153 PetscScalarKokkosViewType<MemorySpace> kv; 154 155 PetscFunctionBegin; 156 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 157 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 158 PetscValidPointer(ov, 3); 159 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 160 PetscCheck(dim == 3, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 3D but DMDA is %dD", (int)dim); 161 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 162 else PetscCall(VecGetKokkosView(vec, &kv)); 163 *ov = PetscScalarKokkosOffsetView3DType<MemorySpace>(kv.data(), {gzs * dof, gys * dof, gxs * dof}, {(gzs + gzm) * dof, (gys + gym) * dof, (gxs + gxm) * dof}); 164 PetscFunctionReturn(0); 165 } 166 167 template <class MemorySpace> 168 PetscErrorCode DMDAVecRestoreKokkosOffsetView_Private(DM da, Vec vec, PetscScalarKokkosOffsetView3DType<MemorySpace> *ov, PetscBool overwrite) { 169 PetscScalarKokkosViewType<MemorySpace> kv; 170 171 PetscFunctionBegin; 172 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 173 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 174 PetscValidPointer(ov, 3); 175 kv = PetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2)); 176 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 177 else PetscCall(VecRestoreKokkosView(vec, &kv)); 178 PetscFunctionReturn(0); 179 } 180 181 template <class MemorySpace> 182 PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView3DType<MemorySpace> *ov) { 183 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 184 ConstPetscScalarKokkosViewType<MemorySpace> kv; 185 186 PetscFunctionBegin; 187 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 188 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 189 PetscValidPointer(ov, 3); 190 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 191 PetscCheck(dim == 3, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 3D but DMDA is %dD", (int)dim); 192 PetscCall(VecGetKokkosView(vec, &kv)); 193 *ov = ConstPetscScalarKokkosOffsetView3DType<MemorySpace>(kv.data(), {gzs * dof, gys * dof, gxs * dof}, {(gzs + gzm) * dof, (gys + gym) * dof, (gxs + gxm) * dof}); 194 PetscFunctionReturn(0); 195 } 196 197 template <class MemorySpace> 198 PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, ConstPetscScalarKokkosOffsetView3DType<MemorySpace> *ov) { 199 ConstPetscScalarKokkosViewType<MemorySpace> kv; 200 201 PetscFunctionBegin; 202 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 203 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 204 PetscValidPointer(ov, 3); 205 kv = ConstPetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2)); 206 PetscCall(VecRestoreKokkosView(vec, &kv)); 207 PetscFunctionReturn(0); 208 } 209 210 /* Function template explicit instantiation */ 211 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView1D *); 212 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView1D *); 213 template <> 214 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView1D *ov) { 215 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 216 } 217 template <> 218 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView1D *ov) { 219 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 220 } 221 template <> 222 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView1D *ov) { 223 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 224 } 225 template <> 226 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView1D *ov) { 227 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 228 } 229 230 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView2D *); 231 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView2D *); 232 template <> 233 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 234 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 235 } 236 template <> 237 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 238 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 239 } 240 template <> 241 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 242 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 243 } 244 template <> 245 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 246 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 247 } 248 249 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView3D *); 250 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView3D *); 251 template <> 252 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 253 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 254 } 255 template <> 256 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 257 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 258 } 259 template <> 260 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 261 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 262 } 263 template <> 264 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 265 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 266 } 267 268 #if !defined(KOKKOS_ENABLE_DEFAULT_DEVICE_TYPE_HOST) /* Get host views if the default memory space is not host space */ 269 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView1DHost *); 270 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView1DHost *); 271 template <> 272 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView1DHost *ov) { 273 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 274 } 275 template <> 276 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView1DHost *ov) { 277 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 278 } 279 template <> 280 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView1DHost *ov) { 281 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 282 } 283 template <> 284 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView1DHost *ov) { 285 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 286 } 287 288 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView2DHost *); 289 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView2DHost *); 290 template <> 291 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 292 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 293 } 294 template <> 295 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 296 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 297 } 298 template <> 299 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 300 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 301 } 302 template <> 303 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 304 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 305 } 306 307 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView3DHost *); 308 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM, Vec, ConstPetscScalarKokkosOffsetView3DHost *); 309 template <> 310 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 311 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 312 } 313 template <> 314 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 315 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_FALSE); 316 } 317 template <> 318 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 319 return DMDAVecGetKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 320 } 321 template <> 322 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 323 return DMDAVecRestoreKokkosOffsetView_Private(da, vec, ov, PETSC_TRUE); 324 } 325 #endif 326 327 /* ============================== 2D including DOF ================================= */ 328 template <class MemorySpace> 329 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView2DType<MemorySpace> *ov, PetscBool overwrite) { 330 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 331 PetscScalarKokkosViewType<MemorySpace> kv; 332 333 PetscFunctionBegin; 334 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 335 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 336 PetscValidPointer(ov, 3); 337 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 338 PetscCheck(dim == 1, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 2D but DMDA is %dD", (int)dim); 339 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 340 else PetscCall(VecGetKokkosView(vec, &kv)); 341 *ov = PetscScalarKokkosOffsetView2DType<MemorySpace>(kv.data(), {gxs, 0}, {gxs + gxm, dof}); 342 PetscFunctionReturn(0); 343 } 344 345 template <class MemorySpace> 346 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView2DType<MemorySpace> *ov, PetscBool overwrite) { 347 PetscScalarKokkosViewType<MemorySpace> kv; 348 349 PetscFunctionBegin; 350 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 351 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 352 PetscValidPointer(ov, 3); 353 kv = PetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1)); 354 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 355 else PetscCall(VecRestoreKokkosView(vec, &kv)); 356 PetscFunctionReturn(0); 357 } 358 359 template <class MemorySpace> 360 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView2DType<MemorySpace> *ov) { 361 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 362 ConstPetscScalarKokkosViewType<MemorySpace> kv; 363 364 PetscFunctionBegin; 365 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 366 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 367 PetscValidPointer(ov, 3); 368 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 369 PetscCheck(dim == 1, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 2D but DMDA is %dD", (int)dim); 370 PetscCall(VecGetKokkosView(vec, &kv)); 371 *ov = ConstPetscScalarKokkosOffsetView2DType<MemorySpace>(kv.data(), {gxs, 0}, {gxs + gxm, dof}); 372 PetscFunctionReturn(0); 373 } 374 375 template <class MemorySpace> 376 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView2DType<MemorySpace> *ov) { 377 ConstPetscScalarKokkosViewType<MemorySpace> kv; 378 379 PetscFunctionBegin; 380 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 381 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 382 PetscValidPointer(ov, 3); 383 kv = ConstPetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1)); 384 PetscCall(VecRestoreKokkosView(vec, &kv)); 385 PetscFunctionReturn(0); 386 } 387 388 /* ============================== 3D including DOF ================================= */ 389 template <class MemorySpace> 390 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView3DType<MemorySpace> *ov, PetscBool overwrite) { 391 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 392 PetscScalarKokkosViewType<MemorySpace> kv; 393 394 PetscFunctionBegin; 395 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 396 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 397 PetscValidPointer(ov, 3); 398 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 399 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 3D but DMDA is %dD", (int)dim); 400 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 401 else PetscCall(VecGetKokkosView(vec, &kv)); 402 *ov = PetscScalarKokkosOffsetView3DType<MemorySpace>(kv.data(), {gys, gxs, 0}, {gys + gym, gxs + gxm, dof}); 403 PetscFunctionReturn(0); 404 } 405 406 template <class MemorySpace> 407 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView3DType<MemorySpace> *ov, PetscBool overwrite) { 408 PetscScalarKokkosViewType<MemorySpace> kv; 409 410 PetscFunctionBegin; 411 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 412 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 413 PetscValidPointer(ov, 3); 414 kv = PetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2)); 415 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 416 else PetscCall(VecRestoreKokkosView(vec, &kv)); 417 PetscFunctionReturn(0); 418 } 419 420 template <class MemorySpace> 421 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView3DType<MemorySpace> *ov) { 422 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 423 ConstPetscScalarKokkosViewType<MemorySpace> kv; 424 425 PetscFunctionBegin; 426 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 427 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 428 PetscValidPointer(ov, 3); 429 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 430 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 3D but DMDA is %dD", (int)dim); 431 PetscCall(VecGetKokkosView(vec, &kv)); 432 *ov = ConstPetscScalarKokkosOffsetView3DType<MemorySpace>(kv.data(), {gys, gxs, 0}, {gys + gym, gxs + gxm, dof}); 433 PetscFunctionReturn(0); 434 } 435 436 template <class MemorySpace> 437 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView3DType<MemorySpace> *ov) { 438 ConstPetscScalarKokkosViewType<MemorySpace> kv; 439 440 PetscFunctionBegin; 441 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 442 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 443 PetscValidPointer(ov, 3); 444 kv = ConstPetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2)); 445 PetscCall(VecRestoreKokkosView(vec, &kv)); 446 PetscFunctionReturn(0); 447 } 448 449 /* ============================== 4D including DOF ================================= */ 450 template <class MemorySpace> 451 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView4DType<MemorySpace> *ov, PetscBool overwrite) { 452 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 453 PetscScalarKokkosViewType<MemorySpace> kv; 454 455 PetscFunctionBegin; 456 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 457 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 458 PetscValidPointer(ov, 3); 459 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 460 PetscCheck(dim == 3, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 4D but DMDA is %dD", (int)dim); 461 if (overwrite) PetscCall(VecGetKokkosViewWrite(vec, &kv)); 462 else PetscCall(VecGetKokkosView(vec, &kv)); 463 *ov = PetscScalarKokkosOffsetView4DType<MemorySpace>(kv.data(), {gzs, gys, gxs, 0}, {gzs + gzm, gys + gym, gxs + gxm, dof}); 464 PetscFunctionReturn(0); 465 } 466 467 template <class MemorySpace> 468 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF_Private(DM da, Vec vec, PetscScalarKokkosOffsetView4DType<MemorySpace> *ov, PetscBool overwrite) { 469 PetscScalarKokkosViewType<MemorySpace> kv; 470 471 PetscFunctionBegin; 472 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 473 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 474 PetscValidPointer(ov, 3); 475 kv = PetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2) * ov->extent(3)); 476 if (overwrite) PetscCall(VecRestoreKokkosViewWrite(vec, &kv)); 477 else PetscCall(VecRestoreKokkosView(vec, &kv)); 478 PetscFunctionReturn(0); 479 } 480 481 template <class MemorySpace> 482 PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView4DType<MemorySpace> *ov) { 483 PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof; 484 ConstPetscScalarKokkosViewType<MemorySpace> kv; 485 486 PetscFunctionBegin; 487 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 488 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 489 PetscValidPointer(ov, 3); 490 DMDA_VEC_GET_SHAPE(da, vec, xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof); 491 PetscCheck(dim == 3, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "KokkosOffsetView is 4D but DMDA is %dD", (int)dim); 492 PetscCall(VecGetKokkosView(vec, &kv)); 493 *ov = ConstPetscScalarKokkosOffsetView4DType<MemorySpace>(kv.data(), {gzs, gys, gxs, 0}, {gzs + gzm, gys + gym, gxs + gxm, dof}); 494 PetscFunctionReturn(0); 495 } 496 497 template <class MemorySpace> 498 PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, ConstPetscScalarKokkosOffsetView4DType<MemorySpace> *ov) { 499 ConstPetscScalarKokkosViewType<MemorySpace> kv; 500 501 PetscFunctionBegin; 502 PetscValidHeaderSpecificType(da, DM_CLASSID, 1, DMDA); 503 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 504 PetscValidPointer(ov, 3); 505 kv = ConstPetscScalarKokkosViewType<MemorySpace>(ov->data(), ov->extent(0) * ov->extent(1) * ov->extent(2) * ov->extent(3)); 506 PetscCall(VecRestoreKokkosView(vec, &kv)); 507 PetscFunctionReturn(0); 508 } 509 510 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView2D *); 511 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView2D *); 512 template <> 513 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 514 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 515 } 516 template <> 517 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 518 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 519 } 520 template <> 521 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 522 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 523 } 524 template <> 525 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2D *ov) { 526 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 527 } 528 529 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView3D *); 530 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView3D *); 531 template <> 532 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 533 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 534 } 535 template <> 536 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 537 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 538 } 539 template <> 540 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 541 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 542 } 543 template <> 544 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3D *ov) { 545 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 546 } 547 548 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView4D *); 549 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView4D *); 550 template <> 551 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView4D *ov) { 552 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 553 } 554 template <> 555 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView4D *ov) { 556 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 557 } 558 template <> 559 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView4D *ov) { 560 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 561 } 562 template <> 563 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView4D *ov) { 564 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 565 } 566 567 #if !defined(KOKKOS_ENABLE_DEFAULT_DEVICE_TYPE_HOST) /* Get host views if the default memory space is not host space */ 568 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView2DHost *); 569 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView2DHost *); 570 template <> 571 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 572 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 573 } 574 template <> 575 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 576 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 577 } 578 template <> 579 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 580 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 581 } 582 template <> 583 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView2DHost *ov) { 584 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 585 } 586 587 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView3DHost *); 588 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView3DHost *); 589 template <> 590 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 591 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 592 } 593 template <> 594 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 595 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 596 } 597 template <> 598 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 599 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 600 } 601 template <> 602 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView3DHost *ov) { 603 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 604 } 605 606 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView4DHost *); 607 template PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM, Vec, ConstPetscScalarKokkosOffsetView4DHost *); 608 template <> 609 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView4DHost *ov) { 610 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 611 } 612 template <> 613 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da, Vec vec, PetscScalarKokkosOffsetView4DHost *ov) { 614 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_FALSE); 615 } 616 template <> 617 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView4DHost *ov) { 618 return DMDAVecGetKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 619 } 620 template <> 621 PETSC_VISIBILITY_PUBLIC PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da, Vec vec, PetscScalarKokkosOffsetView4DHost *ov) { 622 return DMDAVecRestoreKokkosOffsetViewDOF_Private(da, vec, ov, PETSC_TRUE); 623 } 624 #endif 625