<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Fri, Jan 10, 2014 at 8:38 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">Hello,<br>
        I am getting an error using VecSetValuesSectionF90. I had this<br>
working before, i do not know whats wrong now. The call is deep inside<br>
my code so i am not able to create a test case.<br></blockquote><div><br></div><div>I believe this is our bug, from the name change. If you look in</div><div><br></div><div>  include/finclude/ftn-custom/petscvec.h90:43</div>
<div><br></div><div>you see that the function name is wrong. This causes the Fortran compiler to emit wrong code and</div><div>never tell us about it. If you change the name here and recompile it should work. I will push this fix</div>
<div>tomorrow.</div><div><br></div><div>  Thanks,</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">
I am using: Petsc Development GIT revision:<br>
b45093d79c132b908a729253318f899155602205  GIT Date: 2013-11-11 1<br>
6:19:16 -0600<br>
Hope the code chunk below, can help identify the issue.<br>
<br>
 The following is the sequence of calls:<br>
<br>
    DM :: dm<br>
  Vec :: localX<br>
  InsertMode :: mode<br>
  PetscErrorCode :: ierr<br>
  !<br>
  PetscSection :: section ! defualt section<br>
  PetscScalar,allocatable :: values(:) ! size numDofPerPoint<br>
  integer :: vertexId<br>
  integer :: vertexIdStart, vertexIdEnd<br>
  integer :: vertextdim<br>
  integer :: numConstrainedDOF<br>
  integer :: rgnId,bcId,offset<br>
  integer :: IdStart,IdEnd<br>
  integer :: fc<br>
  vertextdim = 0<br>
  call DMPlexGetDepthStratum(dm,vertextdim,vertexIdStart,vertexIdEnd, ierr)<br>
  call DMGetDefaultSection(dm, section, ierr)<br>
  call PetscSectionView(section,PETSC_VIEWER_STDOUT_SELF,ierr)<br>
  allocate(values(this%dofMap%numFieldComp))<br>
  print*, 'size of values is',this%dofMap%numFieldComp<br>
  do vertexId=vertexIdStart,vertexIdEnd-1<br>
    call PetscSectionGetConstraintDof(section,vertexId,numConstrainedDOF,ierr)<br>
    print*,'vertex dof',vertexId,numConstrainedDOF<br>
    if(numConstrainedDOF > 0) then<br>
      call DMPlexGetLabelValue(dm,'Boundary',vertexId,rgnId, ierr)<br>
        values = 0.0_WP<br>
        offset = 0<br>
        do fc = 1,this%numField<br>
        ! check if Field at vertextId is contained<br>
        IdStart = offset + 1<br>
        IdEnd = offset + this%numComp(fc)<br>
        !if(fieldIdisconstrained) then<br>
          bcId = this%getBCId(fc,rgnId)<br>
          values(IdStart:IdEnd) = this%BC(bcId)%getBC()<br>
        !end if<br>
        offset = offset + this%numComp(fc)<br>
        end do<br>
        print*,'Setting bc value',values(1),'at vertex',vertexId,rgnId<br>
        call VecSetValuesSectionF90(localX, section, vertexId, values,Mode,ierr)<br>
      end if<br>
  end do<br>
  deallocate(values)<br>
<br>
 The three print statements above print the following lines.<br>
<br>
 size of values is           1<br>
 vertex dof         600           1<br>
 Setting bc value  -35.7688912272406      at vertex         600           1<br>
<br>
 So, i am passing the correct size array to VecSetValuesSectionF90.<br>
However, when i run the code through intel debugger, i get the<br>
following:<br>
<br>
 size of values is           1<br>
 vertex dof         600           1<br>
 Setting bc value  -35.7688912272406      at vertex         600           1<br>
Program received signal SIGSEGV<br>
VecSetValuesSection (v=0x2bf5620, s=0x2b41a60, point=600,<br>
values=0x447e3d3b071996fc, mode=INSERT_BC_VALUES<br>
) at /home1/00924/Reddy135/LocalApps/petsc/src/vec/vec/utils/vsection.c:185<br>
185                 if (doBC) array[i] = values[i]; /* Constrained update */<br>
(idb) print values[0]<br>
$1 = <no value><br>
<br>
<br>
The value i passed from Fortran code is not missing.<br>
<br>
The relevant call trace is:<br>
<br>
(idb) bt<br>
#0  0x00002b477a299b2d in VecSetValuesSection (v=0x2bf5620,<br>
s=0x2b41a60, point=600, values=0x447e3d3b071996fc, mode<br>
=INSERT_BC_VALUES) at<br>
/home1/00924/Reddy135/LocalApps/petsc/src/vec/vec/utils/vsection.c:185<br>
#1  0x00002b477a29aead in vecsetvaluessectionf90_ (v=0x7fff16adce30,<br>
section=0x7fff16adc790, point=0x7fff16adc6f4,<br>
ptr=0x2b12740, mode=0x9aa990, __ierr=0x7fff16adce08) at<br>
/home1/00924/Reddy135/LocalApps/petsc/src/vec/vec/utils/f90<br>
-custom/zvsectionf90.c:17<br>
<br>
<br>
Am i doing some thing wrong ?<br>
<br>
Thanks<br>
<span class="HOEnZb"><font color="#888888">Reddy<br>
<br>
--<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>