[petsc-users] Getting DMDA local coordinates

[Åsmund Ervik]

Thanks Jed,

The following works fine (for future reference):

  call DMDASetUniformCoordinates(SolScal,0.1,1.0,0.0,1.2,0.0,1.0,ierr)

  call DMGetCoordinatesLocal(SolScal,coordVec,ierr)
  call DMGetCoordinateDM(SolScal,CoordDM,ierr)
  call DMDAVecGetArrayF90(CoordDM,coordVec,coords,ierr)

  write(*,*) coords(:,0,0,0)
  write(*,*) coords(:,imax-1,jmax-1,kmax-1)
gives
  0.1  0.0  0.0
  1.0  1.2  1.0

Perfect. Not very logical though, that we need this new DM, so perhaps
add an example/something in the manual?


I'm by the way planning to send you a Fortran example solving van der
Pol with a 3D DMDA (if you want it). I could include the above code in
that example if desired.  The example is written with some utility
functions (discussed in emails a few weeks ago) such that it is easy to
reuse 95% of a legacy code solving van der Pol (or other PDE). The main
motivation for adding this example (IMO) is that none of the current DM
examples show how DMDAs are used for solving a PDE in a very simple
case, particularly not for Fortran.

One question in this regard: can you have examples split into two files?
That makes it easier to use Fortran modules in this case.


Åsmund


On 21. okt. 2013 13:55, Jed Brown wrote:
> Åsmund Ervik <asmund.ervik at ntnu.no> writes:
> 
>> Dear all,
>>
>> I'm using a DMDA to manage my grid, and I want to have uniform
>> coordinates, so I use DMDASetUniformCoordinates. Since I let PETSc
>> decide the grid partitioning, I need to get the local coordinates, i.e.
>> xmin,xmax for the local grid. I found the function
>> DMGetCoordinatesLocal, but I can't seem to figure out how it works (not
>> in the manual, examples are confusing). What I do is
>>
>>   call DMDASetUniformCoordinates(SolScal,0.1,1.0,0.0,1.2,0.0,1.0,ierr)
>>   call DMGetCoordinatesLocal(SolScal,coordVec,ierr)
>>   call DMDAVecGetArrayF90(SolScal,coordVec,coords,ierr)
>>
>> Then I have the "coords" F90 array, which has the same dimensions as my
>> local part of the DA, i.e. dof:imax:jmax:kmax (all starting at 0). 
> 
> Actually, the dimensions are dim:imax:jmax:kmax.  Use
> 
>   call DMGetCoordinateDM(SolScal,CoordDM,ierr)
>   call DMDAVecGetArrayF90(CoordDM,coordVec,coords,ierr)
> 
>> Where in this array can I find the values xmin,xmax etc? It looks
>> really strange, e.g. a slice in i-direction gives 0.10 0.00 1.20 0.55
>> 0.00
>>
>> I do realize I can compute these on my own, but I think using
>> DMGetCoordinatesLocal has less potential for errors.
>>
>> Best regards,
>> Åsmund Ervik