[petsc-users] Field wise residual norm

Matthew Knepley knepley at gmail.com
Tue Jan 14 21:09:39 CST 2014


On Mon, Jan 13, 2014 at 5:33 PM, Dharmendar Reddy
<dharmareddy84 at gmail.com>wrote:

> Hello,
>         Is the following code on the lines of the approach you suggested:
>
> How do i handle the fortran call to DMCreateFieldIS where i am not
> interested  in fieldNames and numFields as output.
> Do i pass PETSC_NULL_OBJECT ?
>
> subroutine computeResidualNormV2(meshDM,F,fieldId,fNorm)
>   implicit none
> #include "finclude/petsc.h90"
>   ! IO Variables
>   DM, intent(inout) :: meshDM
>   Vec, intent(inout) :: F
>   integer, intent(in ) :: fieldId
>   real(WP), intent(inout) :: fNorm
>
>   ! Local Variables
>   PetscErrorCode :: ierr
>   integer :: numField
>   IS :: subfieldIS
>   DM :: subDM
>   character(len=50), pointer :: fieldName(:)
>   Is, pointer :: fieldIS(:)
>   VecScatter :: vscat
>
>   call DMCreateSubDM(meshDM,1, fieldId, subfieldIS,subdm, ierr)
>

Yes, this looks right.


>   call DMCreateFieldIS(meshDM, numField, fieldName, fieldIS, ierr)
>

You do not need this.


>   call DMGetGlobalVector(subDM, gVecsubField, ierr)
>
>   call
> VecScatterCreate(F,fieldIS(fieldId),gVecsubField,subfieldIS,vscat,ierr)
>

Here you should have

   call
VecScatterCreate(F,subfieldIS,gVecsubField,PETSC_NULL_OBJECT,vscat,ierr)

   Matt


>   call VecScatterBegin(vscat, F, gVecsubField, INSET_VALUES,
> SCATTER_FORWARD, ierr)
>   call VecScatterEnd(vscat, F, gVecsubField, INSET_VALUES,
> SCATTER_FORWARD, ierr)
>
>   call VecNorm( gVecsubField, NORM_2, fNorm, ierr)
>
>   call DMRestoreGlobalVector(subDM, gVecsubField, ierr)
>
> end subroutine computeResidualNormV2
>
> On Mon, Jan 13, 2014 at 4:53 PM, Dharmendar Reddy
> <dharmareddy84 at gmail.com> wrote:
> > I may need to read more about those functions to see how i can get the
> > norms. But a quick look seem to suggest me that i will need lesser
> > function calls than what i am thinking of doing right now.
> >
> > This is how my approach...Am i doing it right ?
> >
> > subroutine computeResidualNorm(meshDM,F,fieldId,fNorm)
> >   implicit none
> > #include "finclude/petsc.h90"
> >   ! IO Variables
> >   DM, intent(inout) :: meshDM
> >   Vec, intent(inout) :: F
> >   integer, intent(in ) :: fieldId
> >   real(WP), intent(inout) :: fNorm
> >
> >   ! Local Variables
> >   PetscErrorCode :: ierr
> >   PetscSection :: Section
> >   PetscScalar, pointer :: FArray(:)
> >   Vec :: gVec, lVec
> >   PetscScalar, pointer :: lArray(:)
> >   integer :: numFieldDof,offset,dofId
> >   integer :: CellId,CellIdStart,CellIdEnd
> >   PetscScalar :: zero
> >
> >   zero = 0.0_WP
> >   call DMGetGlobalVector(meshDM, gVec, ierr)
> >   call DMGetLocalVector(meshDM, lVec, ierr)
> >   call Vecset(lVec, zero, ierr)
> >
> >   ! Get acess to default section
> >   call DMGetDefaultSection(meshDM,section,ierr)
> >   ! Get chart
> >   call  PetscSectionGetChart(section,cellIdstart,CellIdEnd,ierr)
> >   call vecGetArrayF90(F, FArray, ierr)
> >   call VecGetArrayF90(lVec, lArray, ierr)
> >   do cellId=cellIdstart,cellIdEnd-1
> >     call PetscSectionGetFieldDof(section,CellId,fieldId,numFieldDof,ierr)
> >     if(numFieldDof > 0) then
> >       call PetscSectionGetFieldOffset(section,cellId,fieldId,offset,ierr)
> >       do dofId=0,numFieldDof-1
> >         ic = offset+dofId
> >         lArray(ic) = FArray(ic)
> >       end do
> >     end if
> >   end do
> >   call VecRestoreArrayF90(lVec, lArray, ierr)
> >   call vecRestoreArrayF90(F, Farray,ierr)
> >
> >   call DMLocalToGlobalBegin(meshDM, lVec, INSERT_VALUES, gVec, ierr)
> >   call DMLocalToGlobalEnd(meshDM, lVec, INSERT_VALUES, gVec, ierr)
> >
> >   call DMRestoreLocalVector(meshDM, lVec, ierr)
> >
> >   call VecNorm( gVec, NORM_2, fNorm, ierr)
> >
> >   call DMRestoreGlobalVector(meshDM, gVec, ierr)
> > end subroutine computeResidualNorm
> >
> > On Mon, Jan 13, 2014 at 4:15 PM, Matthew Knepley <knepley at gmail.com>
> wrote:
> >> On Mon, Jan 13, 2014 at 4:11 PM, Dharmendar Reddy <
> dharmareddy84 at gmail.com>
> >> wrote:
> >>>
> >>> Hello,
> >>>         Is there a petsc function which can help accomplish the
> following:
> >>> I am solving a coupled system of equations:
> >>>
> >>> F1(x1,x2,x3) = 0
> >>>
> >>> F2(x1,x2,x3) = 0
> >>>
> >>> F3(x1,x2,x3) = 0
> >>>
> >>>
> >>> I am using DMPlex to handle to mesh data. I would like to compute the
> >>> norm of residual of each function.
> >>>
> >>> I assemble the function element wise. I am looking for some interface
> >>> like,
> >>> ComputNorm(dm, F, fieldIndex).
> >>
> >>
> >> Our thinking is generally that you would use
> >>
> >>   DMCreateSubDM()
> >>
> >> and then use VecScatterCreate() with the IS it returns. This will map F
> into
> >> one of the F1/F2/F3.
> >>
> >> Does this make sense?
> >>
> >>    Matt
> >>
> >>>
> >>> Thanks
> >>> Reddy
> >>
> >>
> >>
> >>
> >> --
> >> 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
> >
> >
> >
> > --
> > -----------------------------------------------------
> > Dharmendar Reddy Palle
> > Graduate Student
> > Microelectronics Research center,
> > University of Texas at Austin,
> > 10100 Burnet Road, Bldg. 160
> > MER 2.608F, TX 78758-4445
> > e-mail: dharmareddy84 at gmail.com
> > Phone: +1-512-350-9082
> > United States of America.
> > Homepage: https://webspace.utexas.edu/~dpr342
>
>
>
> --
> -----------------------------------------------------
> Dharmendar Reddy Palle
> Graduate Student
> Microelectronics Research center,
> University of Texas at Austin,
> 10100 Burnet Road, Bldg. 160
> MER 2.608F, TX 78758-4445
> e-mail: dharmareddy84 at gmail.com
> Phone: +1-512-350-9082
> United States of America.
> Homepage: https://webspace.utexas.edu/~dpr342
>



-- 
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/20140114/38405c46/attachment.html>


More information about the petsc-users mailing list