On Thu, Apr 3, 2014 at 11:41 PM, Adrian Croucher
<a.crouc...@auckland.ac.nz>wrote:
> hi,
>
> What is the right way to call DMPlexPointLocalRef() from fortran?
>
Right now, there is none. Jed did not write a Fortran interface, and did
not turn off the automatic one which is non-functional.
> I have a test subroutine that I'm basically translating from the
> CreatePartitionVec() function in the TS ex 33 example problem.
>
> Near the end there's a loop which in C goes:
>
> ierr = VecGetArray(*partition, &part);CHKERRQ(ierr);
> for (c = cStart; c < cEnd; ++c) {
> PetscScalar *p;
> ierr = DMPlexPointLocalRef(*dmCell, c, part, &p);CHKERRQ(ierr);
> p[0] = rank;
> }
> ierr = VecRestoreArray(*partition, &part);CHKERRQ(ierr);
>
>
> So in fortran I've done this:
>
> call VecGetArrayF90(partition, part, ierr); CHKERRQ(ierr)
> do c = cStart, cEnd-1
> call DMPlexPointLocalRef(cell_dm, c, part, p, ierr); CHKERRQ(ierr)
>
Instead of this, use
call DMPlexGetDefaultSection(cell_dm, section, ierr); (Move outside of
loop)
call PetscSectionGetOffset(section, c, off, ierr);
part(off) =rank
Thanks,
Matt
p(1) = rank
> end do
> call VecRestoreArrayF90(partition, part, ierr); CHKERRQ(ierr)
>
> where I've previously declared
>
> PetscScalar, pointer :: part(:)
> PetscScalar, pointer:: p(:)
>
> However the compiler complains:
>
> Error: Rank mismatch in argument 'array' at (1) (scalar and rank-1)
>
> It's not liking the 'part' argument, but I'm not sure how it's supposed to
> be passed in. Any clues?
>
> Cheers, Adrian
>
> --
> Dr Adrian Croucher
> Department of Engineering Science
> University of Auckland
> New Zealand
> tel 64-9-373-7599 ext 84611
>
>
--
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
