Actual source code: petscdmda_kokkos.hpp

  1: #if !defined(PETSCDMDA_KOKKOS_HPP)
  2: #define PETSCDMDA_KOKKOS_HPP

  4: #include <petscvec_kokkos.hpp>
  5: #include <petscdmda.h>

  7: #if defined(PETSC_HAVE_KOKKOS)
  8: #include <Kokkos_Core.hpp>
  9: #include <Kokkos_OffsetView.hpp>

 11: /*@C
 12:    DMDAVecGetKokkosOffsetView - Gets a Kokkos OffsetView that contains up-to-date data of a vector in the given memory space.

 14:    Synopsis:
 15: #include <petscdmda_kokkos.hpp>
 16:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar*,Kokkos::LayoutLeft,MemorySpace>* kv);
 17:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar*,Kokkos::LayoutRight,MemorySpace>* kv);
 18:    PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar*,Kokkos::LayoutRight,MemorySpace>* kv);
 19:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
 20:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
 21:    PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
 22:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
 23:    PetscErrorCode DMDAVecGetKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
 24:    PetscErrorCode DMDAVecGetKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);

 26:    Logically collective on da

 28:    Input Parameters:
 29: +  da - the distributed array
 30: -  v - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()

 32:    Output Parameter:
 33: .  kv - the Kokkos OffsetView with a user-specified template parameter MemorySpace

 35:    Notes:
 36:     Call DMDAVecRestoreKokkosOffsetView() or DMDAVecRestoreKokkosOffsetViewWrite() once you have finished accessing the OffsetView.

 38:     If the vector is not a Kokkos vector, an error will be raised.

 40:     If the vector is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
 41:     a global vector then the ghost points are not accessible. Of course with the local vector you will have to do the
 42:     appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.

 44:     These routines are similar to DMDAVecGetArray() and friends. One can read-only, write-only or read/write access the returned
 45:     Kokkos OffsetView.  Note that passing in a constant OffsetView enables read-only access.
 46:     Currently, only two memory spaces are supported: Kokkos::HostSpace and Kokkos::DefaultExecutionSpace::memory_space.
 47:     If needed, a memory copy will be internally called to copy the latest vector data to the specified memory space.

 49:     In C, to access the returned array of DMDAVecGetArray(), the indexing is "backwards", i.e., array[k][j][i] (instead of array[i][j][k]),
 50:     where i, j, k are loop variables for the x, y, z dimensions respectively specified in DMDACreate3d(), for example.

 52:     To give users the same experience as DMDAVecGetArray(), we mandate the returned OffsetView always has Kokkos::LayoutRight (that is, rightest
 53:     subscript has a stride 1, as in C multi-dimensional arrays), regardless of whether the memory space is host or device. Thus it is important
 54:     to use Iterate::Right as IterateInner if one uses Kokkos::MDRangePolicy to access the OffsetView.

 56:     Note that the OffsetView kv's first dimension (i.e., the leftest, dim 0) corresponds to DA's z direction, and its last dimension
 57:     (rightest) corresponds to DA's x direction.

 59:     If the vector is a global vector, we have
 60: .vb
 61:       kv.extent(0) = zm*dof,  kv.begin(0) = zs*dof, kv.end(0) = (zs+zm)*dof
 62:       kv.extent(1) = ym*dof,  kv.begin(1) = ys*dof, kv.end(1) = (ys+ym)*dof
 63:       kv.extent(2) = xm*dof,  kv.begin(2) = xs*dof, kv.end(2) = (xs+xm)*dof
 64: .ve
 65:     If the vector is a local vector, we have
 66: .vb
 67:       kv.extent(0) = gzm*dof,  kv.begin(0) = gzs*dof, kv.end(0) = (gzs+gzm)*dof
 68:       kv.extent(1) = gym*dof,  kv.begin(1) = gys*dof, kv.end(1) = (gys+gym)*dof
 69:       kv.extent(2) = gxm*dof,  kv.begin(2) = gxs*dof, kv.end(2) = (gxs+gxm)*dof
 70: .ve

 72:     The starts and widths above are obtained by
 73: .vb
 74:      DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
 75:      DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
 76: .ve

 78:     For example, to initialize a grid,

 80: .vb
 81:     typedef Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace> PetscScalarKokkosOffsetView3D;

 83:     PetscScalarKokkosOffsetView3D kv;
 84:     DMDAVecGetKokkosOffsetViewWrite(da,v,&kv); // v is a global vector and we assume dof is 1
 85:     DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);

 87:     parallel_for ("Label",MDRangePolicy <Rank<3, Iterate::Right, Iterate::Right>>(
 88:       {zs,ys,xs},{zs+zm,ys+ym,xs+xm}), KOKKOS_LAMBDA (PetscInt k,PetscInt j,PetscInt i) {
 89:       kv(k,j,i) = ...;
 90:     });
 91:     DMDAVecRestoreKokkosOffsetViewWrite(da,v,&kv);
 92: .ve

 94:     For a multi-component problem, one could cast the returned OffsetView to a user's type. But one has also to shrink
 95:     the OffsetView's extent accordingly. For example,
 96: .vb
 97:     typedef struct {
 98:       PetscScalar omega,temperature;
 99:     } Node;

101:     using NodeKokkosOffsetView3D = Kokkos::Experimental::OffsetView<const Node***,Kokkos::LayoutRight,MemorySpace>;
102:     DMDAVecGetKokkosOffsetViewWrite(da,v,&tv);
103:     NodeKokkosOffsetView3D kv(reinterpret_cast<Node*>(tv.data()),{tv.begin(0)/dof,tv.begin(1)/dof,tv.begin(2)/dof}, {tv.end(0)/dof,tv.end(1)/dof,tv.end(2)/dof});

105:     parallel_for ("Label",MDRangePolicy<Rank<3, Iterate::Right, Iterate::Right>>(
106:       {zs,ys,xs},{zs+zm,ys+ym,xs+xm}), KOKKOS_LAMBDA (PetscInt k,PetscInt j,PetscInt i) {
107:       kv(k,j,i).omega = ...;
108:     });
109:     DMDAVecRestoreKokkosOffsetViewWrite(da,v,&tv);
110: .ve

112:   Level: intermediate

114: .seealso: DMDAVecRestoreKokkosOffsetView(), DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
115:           DMDAVecGetArrayDOF(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
116:           DMStagVecGetArray()
117: @*/
118: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar*,MemorySpace>*);
119: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar*,MemorySpace>*);
120: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar*,MemorySpace>*);

122: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
123: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
124: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);

126: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);
127: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetView         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);
128: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);

130: /*@C
131:    DMDAVecRestoreKokkosOffsetView - Returns the Kokkos OffsetView that gotten from DMDAVecGetKokkosOffsetView()

133:    Synopsis:
134: #include <petscdmda_kokkos.hpp>
135:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar*,Kokkos::LayoutLeft,MemorySpace>* kv);
136:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar*,Kokkos::LayoutRight,MemorySpace>* kv);
137:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar*,Kokkos::LayoutRight,MemorySpace>* kv);
138:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
139:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
140:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
141:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
142:    PetscErrorCode DMDAVecRestoreKokkosOffsetView(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
143:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);

145:    Logically collective on da

147:    Input Parameters:
148: +  da - the distributed array
149: .  v - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
150: -  kv - the Kokkos OffsetView with a user-specified template parameter MemorySpace

152:    Notes:
153:     If the vector is not of type VECKOKKOS, an error will be raised.

155:   Level: intermediate

157: .seealso: DMDAVecGetKokkosOffsetView(), DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
158:           DMDAVecGetArrayDOF(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
159:           DMStagVecGetArray()
160: @*/
161: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar*,MemorySpace>*);
162: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar*,MemorySpace>*);
163: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar*,MemorySpace>*);

165: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
166: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
167: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);

169: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);
170: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetView     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);
171: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***,Kokkos::LayoutRight,MemorySpace>*);

173: /*@C
174:    DMDAVecGetKokkosOffsetViewDOF - Gets a Kokkos OffsetView that contains up-to-date data of a vector in the given memory space, with DOF as the rightest dimension of the OffsetView

176:    Synopsis:
177: #include <petscdmda_kokkos.hpp>
178:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar**,Kokkos::LayoutLeft,MemorySpace>* kv);
179:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
180:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
181:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
182:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
183:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
184:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);
185:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);
186:    PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);

188:    Logically collective on da

190:    Input Parameters:
191: +  da - the distributed array
192: -  v - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()

194:    Output Parameter:
195: .  kv - the Kokkos OffsetView with a user-specified template parameter MemorySpace

197:    Notes:
198:     Call DMDAVecRestoreKokkosOffsetViewDOF() or DMDAVecRestoreKokkosOffsetViewDOFWrite() once you have finished accessing the OffsetView.

200:     If the vector is not a Kokkos vector, an error will be raised.

202:     If the vector is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
203:     a global vector then the ghost points are not accessible. Of course with the local vector you will have to do the
204:     appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.

206:     These routines are similar to DMDAVecGetArrayDOF() and friends. One can read-only, write-only or read/write access the returned
207:     Kokkos OffsetView.  Note that passing in a constant OffsetView enables read-only access.
208:     Currently, only two memory spaces are supported: Kokkos::HostSpace and Kokkos::DefaultExecutionSpace::memory_space.
209:     If needed, a memory copy will be internally called to copy the latest vector data to the given memory space.

211:     In C, to access the returned array of DMDAVecGetArrayDOF(), the indexing is "backwards", i.e., array[k][j][i][c] (instead of array[c][i][j][k]),
212:     where i, j, k are loop variables for the x, y, z dimensions respectively, and c is the loop variable for DOFs, as specified in DMDACreate3d(),
213:     for example.

215:     To give users the same experience as DMDAVecGetArrayDOF(), we mandate the returned OffsetView always has Kokkos::LayoutRight (that is, rightest
216:     subscript has a stride 1, as in C multi-dimensional arrays), regardless of whether the memory space is host or device. Thus it is important
217:     to use Iterate::Right as IterateInner if one uses Kokkos::MDRangePolicy to access the OffsetView.

219:     Note that for a 3D DA, the OffsetView kv's first dimension (i.e., the leftest, dim 0) corresponds to DA's z direction, and its second-to-last dimension
220:     (rightest) corresponds to DA's x direction.

222:     If the vector is a global vector, we have
223: .vb
224:       kv.extent(0) = zm,  kv.begin(0) = zs, kv.end(0) = zs+zm
225:       kv.extent(1) = ym,  kv.begin(1) = ys, kv.end(1) = ys+ym
226:       kv.extent(2) = xm,  kv.begin(2) = xs, kv.end(2) = xs+xm
227:       kv.extent(3) = dof, kv.begin(3) = 0,  kv.end(3) = dof
228: .ve
229:     If the vector is a local vector, we have
230: .vb
231:       kv.extent(0) = gzm, kv.begin(0) = gzs, kv.end(0) = gzs+gzm
232:       kv.extent(1) = gym, kv.begin(1) = gys, kv.end(1) = gys+gym
233:       kv.extent(2) = gxm, kv.begin(2) = gxs, kv.end(2) = gxs+gxm
234:       kv.extent(3) = dof, kv.begin(3) = 0,   kv.end(3) = dof
235: .ve

237:     The starts and widths above are obtained by
238: .vb
239:      DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
240:      DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
241: .ve

243:     For example, to initialize a grid,
244: .vb
245:     typedef Kokkos::Experimental::OffsetView<const PetscScalar****,Kokkos::LayoutRight,MemorySpace> PetscScalarKokkosOffsetView4D;

247:     PetscScalarKokkosOffsetView4D kv;
248:     DMDAVecGetKokkosOffsetViewDOFWrite(da,v,&kv); // v is a global vector
249:     DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);

251:     parallel_for ("Label",MDRangePolicy <Rank<4, Iterate::Right, Iterate::Right>>(
252:       {zs,ys,xs,0},{zs+zm,ys+ym,xs+xm,dof}), KOKKOS_LAMBDA (PetscInt k,PetscInt j,PetscInt i,PetscInt c) {
253:       kv(k,j,i,c) = ...;
254:     });
255:     DMDAVecRestoreKokkosOffsetViewDOFWrite(da,v,&kv);
256: .ve

258:   Level: intermediate

260: .seealso: DMDAVecRestoreKokkosOffsetViewDOF(), DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
261:           DMDAVecGetArrayDOF(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
262:           DMStagVecGetArray()
263: @*/
264: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
265: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
266: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);

268: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);
269: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);
270: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);

272: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);
273: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOF         (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);
274: template<class MemorySpace> PetscErrorCode DMDAVecGetKokkosOffsetViewDOFWrite    (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);

276: /*@C
277:    DMDAVecRestoreKokkosOffsetViewDOF - Returns the Kokkos OffsetView that gotten from DMDAVecGetKokkosOffsetViewDOF()

279:    Synopsis:
280: #include <petscdmda_kokkos.hpp>
281:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
282:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
283:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar**,Kokkos::LayoutRight,MemorySpace>* kv);
284:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
285:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
286:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar***,Kokkos::LayoutRight,MemorySpace>* kv);
287:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<const PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);
288:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);
289:    PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM da,Vec v,Kokkos::Experimental::OffsetView<PetscScalar****,Kokkos::LayoutRight,MemorySpace>* kv);

291:    Logically collective on da

293:    Input Parameters:
294: +  da - the distributed array
295: .  v - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
296: -  kv - the Kokkos OffsetView with a user-specified template parameter MemorySpace

298:    Notes:
299:     If the vector is not of type VECKOKKOS, an error will be raised.

301:   Level: intermediate

303: .seealso: DMDAVecGetKokkosOffsetViewDOF(), DMDAVecGetKokkosOffsetView(), DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
304:           DMDAVecGetArrayDOF(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
305:           DMStagVecGetArray()
306: @*/
307: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
308: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);
309: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar**, Kokkos::LayoutRight,MemorySpace>*);

311: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);
312: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);
313: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar***, Kokkos::LayoutRight,MemorySpace>*);

315: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<const PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);
316: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOF     (DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);
317: template<class MemorySpace> PetscErrorCode DMDAVecRestoreKokkosOffsetViewDOFWrite(DM,Vec,Kokkos::Experimental::OffsetView<      PetscScalar****, Kokkos::LayoutRight,MemorySpace>*);
318: #endif

320: #endif