<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Feb 22, 2018 at 6:10 PM, Danyang Su <span dir="ltr"><<a href="mailto:danyang.su@gmail.com" target="_blank">danyang.su@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<p>Hi Matt,<br>
</p>
<p>Just to let you know that after updating to PETSc 3.8.3 version,
the DMPlexCreateSection in my code now works.</p>
<p>One more question, what is the PETSC_NULL_XXX for IS pointer, as
shown below, in C code, it just pass NULL, but in fortran, what is
the name of null object for pBcCompIS and pBcPointIS. <br></p></div></blockquote><div>Fortran does not "use pointers", so we need to pass a real object that we then convert to a NULL pointer for C.</div><div>PETSC_NULL_IS is a real IS object in Fortran that then gets converted to NULL before calling the C function.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div text="#000000" bgcolor="#FFFFFF"><p>
</p><span class="">
call
DMPlexCreateSection(dmda_flow%<wbr>da,dmda_flow%dim, &<br>
<wbr>
numFields,pNumComp,pNumDof, <wbr> &<br></span>
<wbr>
numBC,pBcField, <wbr> &<br>
<wbr>
pBcCompIS,pBcPointIS, <wbr> &<br>
<wbr> PETSC_NULL_IS,section,ierr)<br>
CHKERRQ(ierr)<br>
<br>
Thanks,<br>
<br>
Danyang<div><div class="h5"><br>
<br>
<div class="m_-4502761336514138491moz-cite-prefix">On 18-02-21 09:22 AM, Danyang Su wrote:<br>
</div>
<blockquote type="cite">
<p>Hi Matt,</p>
<p>To test the Segmentation Violation problem in my code, I
modified the example ex1f90.F to reproduce the problem I have in
my own code. <br>
</p>
<p>If use DMPlexCreateBoxMesh to generate the mesh, the code works
fine. However, if I use DMPlexCreateGmshFromFile, using the same
mesh exported from "DMPlexCreateBoxMesh", it gives Segmentation
Violation error.</p>
<p>Did I miss something in the input mesh file? My first guess is
the label "marker" used in the code, but I couldn't find any
place to set this label.</p>
<p>Would you please let me know how to solve this problem. My code
is done in a similar way as ex1f90, it reads mesh from external
file or creates from cell list, distributes the mesh (these
already work), and then creates sections and sets ndof to the
nodes.</p>
<p>Thanks,</p>
<p>Danyang<br>
</p>
<br>
<div class="m_-4502761336514138491moz-cite-prefix">On 18-02-20 10:07 AM, Danyang Su
wrote:<br>
</div>
<blockquote type="cite">
On 18-02-20 09:52 AM, Matthew Knepley wrote:<br>
<blockquote type="cite">
<div dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">On Tue, Feb 20, 2018 at 12:30 PM,
Danyang Su <span dir="ltr"><<a href="mailto:danyang.su@gmail.com" target="_blank">danyang.su@gmail.com</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div bgcolor="#FFFFFF">
<p>Hi All,</p>
<p>I tried to compile the DMPlexCreateSection code
but got error information as shown below.<br>
</p>
<p>Error: Symbol 'petsc_null_is' at (1) has no
IMPLICIT type</p>
<p>I tried to use PETSC_NULL_OBJECT instead of
PETSC_NULL_IS, then the code can be compiled but
run into Segmentation Violation error in
DMPlexCreateSection.<br>
</p>
</div>
</blockquote>
<div>From the webpage</div>
<div><br>
</div>
<div> <a href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/DMPLEX/DMPlexCreateSection.html" target="_blank">http://www.mcs.anl.gov/petsc/<wbr>petsc-current/docs/<wbr>manualpages/DMPLEX/<wbr>DMPlexCreateSection.html</a> <br>
</div>
<div><br>
</div>
<div>The F90 version is DMPlexCreateSectionF90. Doing
this with F77 arrays would have been too painful.</div>
</div>
</div>
</div>
</blockquote>
Hi Matt,<br>
<br>
Sorry, I still cannot compile the code if use
DMPlexCreateSectionF90 instead of DMPlexCreateSection. Would you
please tell me in more details? <br>
<br>
undefined reference to `dmplexcreatesectionf90_'<br>
<br>
then I #include <petsc/finclude/petscdmplex.<wbr>h90>, but this
throws more error during compilation.<br>
<br>
<br>
Included at
/home/dsu/Soft/PETSc/petsc-3.<wbr>7.5/include/petsc/finclude/<wbr>petscdmplex.h90:6:<br>
Included at ../../solver/solver_ddmethod.<wbr>F90:62:<br>
<br>
PETSCSECTION_HIDE section<br>
1<br>
Error: Unclassifiable statement at (1)<br>
/home/dsu/Soft/PETSc/petsc-3.<wbr>7.5/include/petsc/finclude/<wbr>ftn-custom/petscdmplex.h90:<wbr>167.10:<br>
Included at
/home/dsu/Soft/PETSc/petsc-3.<wbr>7.5/include/petsc/finclude/<wbr>petscdmplex.h90:6:<br>
Included at ../../solver/solver_ddmethod.<wbr>F90:62:<br>
<br>
PETSCSECTION_HIDE section<br>
1<br>
Error: Unclassifiable statement at (1)<br>
/home/dsu/Soft/PETSc/petsc-3.<wbr>7.5/include/petsc/finclude/<wbr>ftn-custom/petscdmplex.h90:<wbr>179.10:<br>
Included at
/home/dsu/Soft/PETSc/petsc-3.<wbr>7.5/include/petsc/finclude/<wbr>petscdmplex.h90:6:<br>
Included at ../../solver/solver_ddmethod.<wbr>F90:62:<br>
<br>
<blockquote type="cite">
<div dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">
<div><br>
</div>
<div> Thanks,</div>
<div><br>
</div>
<div> Matt</div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div bgcolor="#FFFFFF">
<p> </p>
<p>dmda_flow%da is distributed dm object that works
fine. </p>
<p>The fortran example I follow is <a class="m_-4502761336514138491gmail-m_8962138310016109093moz-txt-link-freetext" href="http://www.mcs.anl.gov/petsc/petsc-dev/src/dm/impls/plex/examples/tutorials/ex1f90.F90" target="_blank">http://www.mcs.anl.gov/petsc/p<wbr>etsc-dev/src/dm/impls/plex/exa<wbr>mples/tutorials/ex1f90.F90</a>.
<br>
</p>
<p>What parameters should I use if passing null to
bcField, bcComps, bcPoints and perm.<br>
</p>
<pre style="color:rgb(0,0,0);font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;word-spacing:0px;text-decoration-style:initial;text-decoration-color:initial"><a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscErrorCode.html#PetscErrorCode" target="_blank">PetscErrorCode</a> <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/DMPLEX/DMPlexCreateSection.html#DMPlexCreateSection" target="_blank">DMPlexCreateSection</a>(<a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/DM/DM.html#DM" target="_blank">DM</a> dm, <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> dim, <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> numFields,const <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> numComp[],const <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> numDof[], <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> numBC,const <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/Vec/PetscInt.html#PetscInt" target="_blank">PetscInt</a> bcField[],
const <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/IS/IS.html#IS" target="_blank">IS</a> bcComps[], const <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/IS/IS.html#IS" target="_blank">IS</a> bcPoints[], <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/IS/IS.html#IS" target="_blank">IS</a> perm, <a href="https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/IS/PetscSection.html#PetscSection" target="_blank">PetscSection</a> *section)
</pre>
<p>#include <petsc/finclude/petscis.h><br>
#include <petsc/finclude/petscis.h90><br>
#include <petsc/finclude/petscdmplex.h></p>
<p>...<br>
</p>
<p>#ifdef USG<br>
numFields = 1<br>
numComp(1) = 1<br>
pNumComp => numComp<br>
<br>
do i = 1, numFields*(dmda_flow%dim+1)<br>
numDof(i) = 0<br>
end do<br>
numDof(0*(dmda_flow%dim+1)+1) =
dmda_flow%dof<br>
pNumDof => numDof<br>
<br>
numBC = 0<br>
<br>
call DMPlexCreateSection(dmda_flow%<wbr>da,dmda_flow%dim,
&<br>
<wbr>
numFields,pNumComp,pNumDof, <wbr>
&<br>
<wbr>
numBC,PETSC_NULL_INTEGER, <wbr>
&<br>
<wbr>
PETSC_NULL_IS,PETSC_NULL_IS, <wbr>
& !Error here<br>
<wbr>
PETSC_NULL_IS,section,ierr)<br>
CHKERRQ(ierr)<br>
<br>
call PetscSectionSetFieldName(secti<wbr>on,0,'flow',ierr)<br>
CHKERRQ(ierr)<br>
<br>
call DMSetDefaultSection(dmda_flow%<wbr>da,section,ierr)<br>
CHKERRQ(ierr)<br>
<br>
call PetscSectionDestroy(section,ie<wbr>rr)<br>
CHKERRQ(ierr)<br>
#endif</p>
<p>Thanks,</p>
<p>Danyang<br>
</p>
</div>
</blockquote>
</div>
<br>
<br clear="all">
<div><br>
</div>
-- <br>
<div class="m_-4502761336514138491gmail_signature">
<div dir="ltr">
<div>
<div dir="ltr">
<div>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><br>
</div>
<div><a href="http://www.caam.rice.edu/%7Emk51/" target="_blank">https://www.cse.buffalo.edu/~<wbr>knepley/</a><br>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</blockquote>
<br>
</blockquote>
<br>
</blockquote>
<br>
</div></div></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div>
</div></div>