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