<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Mon, Jan 13, 2014 at 4:53 PM, Dharmendar Reddy <span dir="ltr"><<a href="mailto:dharmareddy84@gmail.com" target="_blank">dharmareddy84@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">I may need to read more about those functions to see how i can get the<br>
norms. But a quick look seem to suggest me that i will need lesser<br>
function calls than what i am thinking of doing right now.<br>
<br>
This is how my approach...Am i doing it right ?<br></blockquote><div><br></div><div>This looks like it should work.</div><div><br></div><div>  Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

subroutine computeResidualNorm(meshDM,F,fieldId,fNorm)<br>
  implicit none<br>
#include "finclude/petsc.h90"<br>
  ! IO Variables<br>
  DM, intent(inout) :: meshDM<br>
  Vec, intent(inout) :: F<br>
  integer, intent(in ) :: fieldId<br>
  real(WP), intent(inout) :: fNorm<br>
<br>
  ! Local Variables<br>
  PetscErrorCode :: ierr<br>
  PetscSection :: Section<br>
  PetscScalar, pointer :: FArray(:)<br>
  Vec :: gVec, lVec<br>
  PetscScalar, pointer :: lArray(:)<br>
  integer :: numFieldDof,offset,dofId<br>
  integer :: CellId,CellIdStart,CellIdEnd<br>
  PetscScalar :: zero<br>
<br>
  zero = 0.0_WP<br>
  call DMGetGlobalVector(meshDM, gVec, ierr)<br>
  call DMGetLocalVector(meshDM, lVec, ierr)<br>
  call Vecset(lVec, zero, ierr)<br>
<br>
  ! Get acess to default section<br>
  call DMGetDefaultSection(meshDM,section,ierr)<br>
  ! Get chart<br>
  call  PetscSectionGetChart(section,cellIdstart,CellIdEnd,ierr)<br>
  call vecGetArrayF90(F, FArray, ierr)<br>
  call VecGetArrayF90(lVec, lArray, ierr)<br>
  do cellId=cellIdstart,cellIdEnd-1<br>
    call PetscSectionGetFieldDof(section,CellId,fieldId,numFieldDof,ierr)<br>
    if(numFieldDof > 0) then<br>
      call PetscSectionGetFieldOffset(section,cellId,fieldId,offset,ierr)<br>
      do dofId=0,numFieldDof-1<br>
        ic = offset+dofId<br>
        lArray(ic) = FArray(ic)<br>
      end do<br>
    end if<br>
  end do<br>
  call VecRestoreArrayF90(lVec, lArray, ierr)<br>
  call vecRestoreArrayF90(F, Farray,ierr)<br>
<br>
  call DMLocalToGlobalBegin(meshDM, lVec, INSERT_VALUES, gVec, ierr)<br>
  call DMLocalToGlobalEnd(meshDM, lVec, INSERT_VALUES, gVec, ierr)<br>
<br>
  call DMRestoreLocalVector(meshDM, lVec, ierr)<br>
<br>
  call VecNorm( gVec, NORM_2, fNorm, ierr)<br>
<br>
  call DMRestoreGlobalVector(meshDM, gVec, ierr)<br>
end subroutine computeResidualNorm<br>
<div class="HOEnZb"><div class="h5"><br>
On Mon, Jan 13, 2014 at 4:15 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br>
> On Mon, Jan 13, 2014 at 4:11 PM, Dharmendar Reddy <<a href="mailto:dharmareddy84@gmail.com">dharmareddy84@gmail.com</a>><br>
> wrote:<br>
>><br>
>> Hello,<br>
>>         Is there a petsc function which can help accomplish the following:<br>
>> I am solving a coupled system of equations:<br>
>><br>
>> F1(x1,x2,x3) = 0<br>
>><br>
>> F2(x1,x2,x3) = 0<br>
>><br>
>> F3(x1,x2,x3) = 0<br>
>><br>
>><br>
>> I am using DMPlex to handle to mesh data. I would like to compute the<br>
>> norm of residual of each function.<br>
>><br>
>> I assemble the function element wise. I am looking for some interface<br>
>> like,<br>
>> ComputNorm(dm, F, fieldIndex).<br>
><br>
><br>
> Our thinking is generally that you would use<br>
><br>
>   DMCreateSubDM()<br>
><br>
> and then use VecScatterCreate() with the IS it returns. This will map F into<br>
> one of the F1/F2/F3.<br>
><br>
> Does this make sense?<br>
><br>
>    Matt<br>
><br>
>><br>
>> Thanks<br>
>> Reddy<br>
><br>
><br>
><br>
><br>
> --<br>
> What most experimenters take for granted before they begin their experiments<br>
> is infinitely more interesting than any results to which their experiments<br>
> lead.<br>
> -- Norbert Wiener<br>
<br>
<br>
<br>
</div></div><span class="HOEnZb"><font color="#888888">--<br>
-----------------------------------------------------<br>
Dharmendar Reddy Palle<br>
Graduate Student<br>
Microelectronics Research center,<br>
University of Texas at Austin,<br>
10100 Burnet Road, Bldg. 160<br>
MER 2.608F, TX 78758-4445<br>
e-mail: <a href="mailto:dharmareddy84@gmail.com">dharmareddy84@gmail.com</a><br>
Phone: <a href="tel:%2B1-512-350-9082" value="+15123509082">+1-512-350-9082</a><br>
United States of America.<br>
Homepage: <a href="https://webspace.utexas.edu/~dpr342" target="_blank">https://webspace.utexas.edu/~dpr342</a><br>
</font></span></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener
</div></div>