[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