Matt, I don't think there exists such a subroutine:
$ cd petsc-dev/src $ ack "GetArrayDOF" | wc -l 51 $ ack "GetArrayDOFF90" | wc -l 0 Åsmund On 07. okt. 2013 14:39, Matthew Knepley wrote: > On Mon, Oct 7, 2013 at 7:30 AM, Åsmund Ervik <[email protected]> wrote: > >> Hi again, >> >> I've started to delve into the world of DMDA. At this point I want to >> write code that goes as follows, but something tells me this may be a >> bad idea? Synopsis of the code given below the snippet. >> >> ! The "array" array is allocated like this (taken from examples): >> ! call DMDAGetGhostCorners(SolScal,ib1,ibn,jb1,jbn,kb1,kbn,ierr) >> ! allocate(array(0:dof-1,ib1:ibn,jb1:jbn,kb1:kbn)) >> >> subroutine petsc_to_local(da,vec,array,f,dof) >> DM, intent(in) :: da >> Vec,intent(in) :: vec >> PetscScalar, pointer, intent(inout) :: array(:,:,:,:) >> integer,intent(in) :: dof >> real,intent(inout),dimension(1-stw:,1-stw:1-stw:,dof) >> ! >> if (dof == 1) then >> call VecGetArrayF90(da,vec,array(0,:,:,:),ierr) >> else >> call VecGetArrayF90(da,vec,array,ierr) >> endif >> call transform_petsc_us(array,f,dof) >> end subroutine petsc_to_local >> subroutine transform_petsc_us(array,f,dof) >> !Note: assumed shape-array does the "coordinate transformation" >> real, intent(in), dimension(dof,1-stw:,1-stw:,1-stw:) :: array >> real,intent(inout),dimension(1-stw:,1-stw:1-stw:,dof) :: f >> integer :: n >> ! >> do n=1,dof >> f(:,:,:,n) = array(n-1,:,:,:) >> enddo >> end subroutine transform_petsc_us >> >> This is basically a wrapper around VecGetArrayF90 that does some things: >> it handles dof=1 or higher in the same way, it reindexes dof to start at >> 1, moves the dof from the first to the last index (these two are for >> compatibility with existing code), and perhaps most importantly it does >> a "coordinate transformation" such that the end result (in the array f >> that my code subsequently uses) is indexed from 1-stw (stw: stencil >> width) to ilmax+stw (ilmax: local index i max). >> > > You might want VecGetArrayDOFF90(). Its too bad F90 pointers cannot be > easily > manipulated like C or having a 1-offset would not require a copy. > > Matt > > >> This enables me to plug existing code straight in and manipulate the >> array f. Then I do the opposite remap and call VecRestoreArrayF90, and >> then DMLocalToGlobalBegin/End. The code between VecGet... and >> VecRestore... will run for several seconds (or longer), it includes a >> Poisson solve with KSP (I haven't figured that part out yet...) >> >> Am I being stupid here, or is this a good solution? I'm guessing this >> ends up making a copy into f, so it may hurt performance somewhat, but I >> don't think that will be critical? >> >> Best regards, >> Åsmund >> >> >> On 04. okt. 2013 22:40, Jed Brown wrote: >>> Matthew Knepley <[email protected]> writes: >>> >>>> On Fri, Oct 4, 2013 at 2:10 PM, Åsmund Ervik <[email protected]> >> wrote: >>>> >>>>> Is there any reason I should prefer Xdmf over Cgns? I think both use >>>>> hdf5 in the background. >>>>> >>>> >>>> I guess not. It was just easier for me to get Xdmf going. We have trial >>>> CGNS as well. >>> >>> I don't think we have a CGNS writer, but DMPlex can read it. Xdmf gives >>> you a little more flexibility in choosing a schema for the HDF5 file, >>> but you can put other stuff in CGNS files and Xdmf has its share of >>> shortcomings. If you are happy with CGNS, I would say to use it and >>> report any issues you have. We'll enhance the implementation so it can >>> do everything you need. >>> >> > > >
