[petsc-users] updating values in a DA Global array
Barry Smith
bsmith at mcs.anl.gov
Thu Jun 24 15:46:05 CDT 2010
It turns out there was a bug in gfortran versions prior to 2.5 that made DAVecGetArrayF90() unworkable. After Satish told me this, ...
I have implemented DAVecGetArrayF90() for 1, 2, and 3 dimensions with dof = 1 and dof > 1 in petsc-dev. The test example is src/dm/da/examples/tutorials/ex11f90.F
Please report any problems with it to petsc-maint at mcs.anl.gov and if you use gfortran make sure it is version 2.5. Not that if your dof passed to DACreate is greater than one then you need another dimension in the array pointer you declare. Also, in keeping with the Fortran tradition indices for it start at 1 not zero.
This can lead to much simplier structured grid Fortran codes.
Barry
On Jun 22, 2010, at 6:53 PM, Matthew Knepley wrote:
> On Tue, Jun 22, 2010 at 11:30 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> Matt,
>
> Actually DAVecRestoreArrayF90() only works for 1 dimensional DAs. The F90 interface has to be written for 2 and 3d and the creation of Fortran 3d arrays (and 4d for dof > 1) stuff written.
>
> Ah, I see. Mark, if you can bear C for 1 function, it will work. The problem for the F90 is that we
> are cutting out a section of a big array, and making it look multidimensional. It is possible with F90,
> but a pain to get everything right in the array descriptor and no one has had to stamina to do it yet.
>
> Sorry about that,
>
> Matt
>
> Barry
>
>
> On Jun 22, 2010, at 6:21 PM, Matthew Knepley wrote:
>
>> On Tue, Jun 22, 2010 at 12:59 PM, Mark Cheeseman <mark.cheeseman at kaust.edu.sa> wrote:
>> Hi,
>>
>> I am trying to write a PETSc program in FORTRAN90 where I need to update a single value in a global distributed array. I know the global coordinates of the position that needs to be updated in the global array but I cannot get the mapping from the local vector correct. In this case, I am working on a domain with global dimensions [arraysize(1),arraysize(2),arraysize(3)] and I want to alter a single point in the global distributed array, uGLOBAL, at the global position [arraysize(1)/2-1,arraysize(2)-1,3]. I cannot seem to be able to do this... what am I doing wrong?
>>
>> ...
>> DA da
>> Vec uGLOBAL, uLOCAL, tmp
>> PetscErrorCode ierr
>> PetscScalar, pointer :: xx
>> PetscInt rank, source_rank, i,j,k, row
>>
>> ....
>>
>> call MPI_Comm_rank( PETSC_COMM_WORLD, rank, ierr )
>> call DACreate3d( PETSC_COMM_WORLD, DA_NONPERIODIC, DA_STENCIL_BOX, &
>> arraysize(1), arraysize(2), arraysize(3), PETSC_DECIDE, &
>> PETSC_DECIDE, PETSC_DECIDE, 1, 5, PETSC_NULL_INTEGER, &
>> PETSC_NULL_INTEGER, PETSC_NULL_INTEGER, da, ierr)
>> call DACreateGlobalVector( da, pNOW, ierr )
>> call DAGetCorners( da, xs, ys, zs, xl, yl, zl, ierr )
>>
>> call DAVecGetArrayF90()
>>
>> do i = xs,xs+xl-1
>> if ( i.eq.arraysize(1)/2-1 ) then
>> do j = ys,ys+yl-1
>> if ( j.eq.arraysize(2)/2-1 ) then
>> do k = zs,zs+zl-1
>> if ( k.eq.3 ) then
>>
>> array(k,j,i) = pressure
>>
>> endif
>> enddo
>> endif
>> enddo
>> endif
>> enddo
>>
>>
>> call
>
>>
>> That should work.
>>
>> Matt
>>
>> Thank you,
>> Mark
>>
>> --
>> Mark Patrick Cheeseman
>>
>> Computational Scientist
>> KSL (KAUST Supercomputing Laboratory)
>> Building 1, Office #126
>> King Abdullah University of Science & Technology
>> Thuwal 23955-6900
>> Kingdom of Saudi Arabia
>>
>> EMAIL : mark.cheeseman at kaust.edu.sa
>> PHONE : +966 (2) 808 0221 (office)
>> +966 (54) 470 1082 (mobile)
>> SKYPE : mark.patrick.cheeseman
>>
>>
>>
>> --
>> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.
>> -- Norbert Wiener
>
>
>
>
> --
> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.
> -- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20100624/047265ad/attachment.htm>
More information about the petsc-users
mailing list