<html>
<head>
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
<div class="moz-cite-prefix">On 08/28/2013 10:28 AM, Olivier
Bonnefon wrote:<br>
</div>
<blockquote cite="mid:521DB49C.6010101@avignon.inra.fr" type="cite">
<meta content="text/html; charset=ISO-8859-1"
http-equiv="Content-Type">
<div class="moz-cite-prefix">On 08/23/2013 04:42 PM, Matthew
Knepley wrote:<br>
</div>
<blockquote
cite="mid:CAMYG4G=Ey+YDUGBO5ZcCK3+=pUrdwY72YEP34aOAUw3ZQkvspA@mail.gmail.com"
type="cite">
<div dir="ltr">On Fri, Aug 23, 2013 at 9:35 AM, Olivier Bonnefon
<span dir="ltr"><<a moz-do-not-send="true"
href="mailto:olivier.bonnefon@avignon.inra.fr"
target="_blank">olivier.bonnefon@avignon.inra.fr</a>></span>
wrote:<br>
<div class="gmail_extra">
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Hello,<br>
<br>
Thanks for your answers, I'm now able to import and
distribute a mesh:<br>
</blockquote>
<div><br>
</div>
<div>You might simplify this to</div>
<div><br>
</div>
<div>if (rank) {obNbCells = 0; obNbVertex = 0;}</div>
<div>ierr =
DMPlexCreateFromCellList(comm,dim,obNbCells,obNbVertex,3,0,obCells,2,obVertex,dm);CHKERRQ(ierr);<br>
</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">if
(!rank){<br>
ierr =
DMPlexCreateFromCellList(comm,dim,obNbCells,obNbVertex,3,0,obCells,2,obVertex,dm);CHKERRQ(ierr);<br>
for (i=0;i<obNbBound;i++){<br>
ierr =DMPlexSetLabelValue(*dm, "marker",
obBoundary[i]+obNbCells, 1);CHKERRQ(ierr);<br>
}<br>
}else {<br>
ierr =
DMPlexCreateFromCellList(comm,dim,0,0,3,0,obCells,2,obVertex,dm);CHKERRQ(ierr);<br>
}<br>
<br>
ierr = DMPlexDistribute(*dm, partitioner, 0,
&distributedMesh);CHKERRQ(ierr);<br>
if (distributedMesh) {<br>
ierr = DMDestroy(dm);CHKERRQ(ierr);<br>
*dm = distributedMesh;<br>
}<br>
<br>
Is it possible to known the resulting partition ? ie,
What is the mapping between the initial cell number and
the local cell (used in DMPlexComputeResidualFEM)?<br>
I need this to write an efficient implementation of the
FEM struct functions f0 and g0, space depending.<br>
</blockquote>
<div><br>
</div>
<div>Yes, but I really do not think you want to do things
that way. I am assuming you want different material
models or something</div>
<div>in different places. The way I envision that is using
a DMLabel to mark up parts of the domain. All labels are
automatically</div>
<div>distributed with the mesh. Is that what you want?</div>
</div>
</div>
</div>
</blockquote>
Hello,<br>
<br>
It is exactly what I need: I'm mobilized a landscape, and the
parameters of the model depend of the type of crop. Therefore, I
have created a label for each type of crop and I have labeled each
triangle with the corresponding label:<br>
<br>
for (i=0;i<obNbCells;i++){<br>
if (labelCells[i]==1){<br>
ierr =DMPlexSetLabelValue(*dm, "marker1", i,
1);CHKERRQ(ierr);<br>
}else{<br>
ierr =DMPlexSetLabelValue(*dm, "marker2", i,
1);CHKERRQ(ierr);<br>
}<br>
}<br>
<br>
So, I'm able to mark the triangles, but I'm not able to get this
label in the plugin "fem.f0Funcs" and "fem.g0Funcs": These plugins
are called by looping on the triangles in the function
"FEMIntegrateResidualBatch", but the dm is not available, so I
can't use the functions DMPlexGetLabel, DMLabelGetStratumSize and
DMLabelGetStratumIS. What is the good way to get the labels in the
user plugins of the fem struct ?<br>
<br>
<br>
Thanks a lot for your help.<br>
<br>
Olivier B<br>
</blockquote>
Hello,<br>
<br>
This is the solution I implemented to get the label level in the
plugins "fem.f0Funcs" and "fem.g0Funcs":<br>
<br>
I need the DM and the index element, so i do:<br>
1) I add some static variables:<br>
static DM * spDM[128];<br>
static int scurElem[128];<br>
2) I overload the DMPlexComputeJacobianFEM with :<br>
PetscErrorCode MY_DMPlexComputeJacobianFEM(DM dm, Vec X, Mat Jac,
Mat JacP, MatStructure *str,void *user)<br>
{<br>
<br>
PetscMPIInt rank;<br>
PetscErrorCode ierr;<br>
<br>
PetscFunctionBeginUser;<br>
ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank);CHKERRQ(ierr);<br>
spDM[rank]=&dm;<br>
PetscFunctionReturn(DMPlexComputeJacobianFEM(dm, X, Jac,JacP,
str,user));<br>
<br>
}<br>
3) overload FEMIntegrateResidualBatch adding code:<br>
.<br>
.<br>
for (e = 0; e < Ne; ++e) {<br>
scurElem[rank]=e;//added ligne<br>
.<br>
.<br>
<br>
So that, I can get the label level using DMPlexHasLabel and
DMLabelGetValue<br>
<br>
I'm sure this solution is awful, and works only in this version, but
i didn't find a better way to get the label in the plugins fem
struc. Do you know the correct way to do that ??<br>
<br>
Thanks,<br>
Olivier B<br>
<blockquote cite="mid:521DB49C.6010101@avignon.inra.fr" type="cite">
<blockquote
cite="mid:CAMYG4G=Ey+YDUGBO5ZcCK3+=pUrdwY72YEP34aOAUw3ZQkvspA@mail.gmail.com"
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>
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Regards,<br>
<br>
Olivier B<span class=""><font color="#888888"><br>
<br>
-- <br>
Olivier Bonnefon<br>
INRA PACA-Avignon, Unité BioSP<br>
Tel: <a moz-do-not-send="true"
href="tel:%2B33%20%280%294%2032%2072%2021%2058"
value="+33432722158" target="_blank">+33 (0)4 32
72 21 58</a><br>
<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>
</blockquote>
<br>
<br>
<pre class="moz-signature" cols="72">--
Olivier Bonnefon
INRA PACA-Avignon, Unité BioSP
Tel: +33 (0)4 32 72 21 58</pre>
</blockquote>
<br>
<br>
<pre class="moz-signature" cols="72">--
Olivier Bonnefon
INRA PACA-Avignon, Unité BioSP
Tel: +33 (0)4 32 72 21 58</pre>
</body>
</html>