<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Mon, Dec 30, 2013 at 7:23 PM, Yaakoub El Khamra <span dir="ltr"><<a href="mailto:yaakoub@tacc.utexas.edu" target="_blank">yaakoub@tacc.utexas.edu</a>></span> wrote:<br>
<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"><div dir="ltr"><div><br></div>Greetings<div>I am try to use DMRefine to refine a dm I create with DMPlexCreateBoxMesh with tonight's checkout. A quick glance at the DMView output from the DM's shows that the refinedMesh looks a lot like the original dm. DMGetRefineLevel returns 1 though. I must be missing something elementary here. Can someone please confirm if there are any glaring mistakes I made? Also if I want to refine and distribute, should I refine on a single processor then distribute or can I refine a couple of times, distribute, and refine some more?</div>
</div></blockquote><div><br></div><div>1) Sorry this took a long time.</div><div><br></div><div>2) Refinement is somewhat complicated right now. If you use a mesh generator to make the mesh (as BoxMesh does), then</div><div>
the default refinement is to use the generator itself. However, then you must call</div><div><br></div><div> ierr = DMPlexSetRefinementLimit(dm, refinementLimit);CHKERRQ(ierr);</div><div><br></div><div> to set the maximum cell volume after refinement. If, however, the mesh is constructed by hand the default is to use uniform</div>
<div> refinement. You can also get uniform refinement using</div><div><br></div><div><div> ierr = DMPlexSetRefinementUniform(dm, refinementUniform);CHKERRQ(ierr);</div></div><div><br></div><div>3) Only uniform refinement works in parallel, so I tend to refine using the generator in serial, distribute, and then uniformly refine.</div>
<div> This is what SNES ex12 does.</div><div><br></div><div>4) Did you have another problem that I was supposed to look at on your machine? I have lost my notes.</div><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"><div dir="ltr"><div>The test code looks as follows:</div>
<div><br></div><div><div> program testRefine</div><div> implicit none</div><div>#include "finclude/petsc.h90"</div><div> DM :: dm</div><div>
DM :: refinedMesh</div><div> PetscInt :: r</div><div><br></div><div> PetscInt :: dim</div><div> PetscBool :: interpolate</div><div> PetscErrorCode :: ierr</div><div><br></div><div> call PetscInitialize(PETSC_NULL_CHARACTER, ierr)</div>
<div> CHKERRQ(ierr)</div><div> dim = 2</div><div> interpolate = PETSC_TRUE</div><div>! Create a mesh</div><div> call DMPlexCreateBoxMesh(PETSC_COMM_WORLD, dim, interpolate, dm,</div><div> $ ierr)</div>
<div> CHKERRQ(ierr)</div><div><br></div><div> call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr)</div><div> CHKERRQ(ierr)</div><div><br></div><div><br></div><div> call DMRefine(dm, PETSC_COMM_WORLD, refinedMesh, ierr)</div>
<div> CHKERRQ(ierr)</div><div><br></div><div> call DMView(refinedMesh, PETSC_VIEWER_STDOUT_WORLD, ierr)</div><div> CHKERRQ(ierr)</div><div><br></div><div> call DMGetRefineLevel(refinedMesh, r, ierr)</div>
<div> write (*,*) "The value of r is:", r</div><div><br></div><div> call DMDestroy(dm, ierr)</div><div> call DMDestroy(refinedMesh, ierr)</div><div> CHKERRQ(ierr)</div><div><br></div><div> call PetscFinalize(ierr)</div>
<div> end program testRefine</div><div><br></div><div><br clear="all"><div><br></div><div>And the output:</div><div><div>DM_0x1122bb0_0 in 2 dimensions:</div><div> 0-cells: 9</div><div> 1-cells: 16</div><div> 2-cells: 8</div>
<div>Labels:</div><div> marker: 2 strata of sizes (16, 5)</div><div> depth: 3 strata of sizes (9, 16, 8)</div><div>DM_0x1122bb0_1 in 2 dimensions:</div><div> 0-cells: 9</div><div> 1-cells: 16</div><div> 2-cells: 8</div>
<div>Labels:</div><div> marker: 2 strata of sizes (16, 1)</div><div> depth: 3 strata of sizes (9, 16, 8)</div><div> The value of r is: 1</div></div><div><br></div><div><br></div><div><br></div><div><br></div>
<div>
<br></div><div><br></div><div><br>Regards<span class=""><font color="#888888"><br>Yaakoub El Khamra<br></font></span></div>
</div>
</div></div>
</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>