<div dir="ltr">On Tue, Apr 2, 2013 at 4:18 PM, Dharmendar Reddy <span dir="ltr"><<a href="mailto:dharmareddy84@gmail.com" target="_blank">dharmareddy84@gmail.com</a>></span> wrote:<br><div class="gmail_extra"><div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div><div>Hello,<br></div> I should be able to generate petscFEM struct with my code, i have the required data, I just need to prepare it in the format used by DMPLexComputeResidualFEM. I will give it a shot. <br>
</div></div></blockquote><div><br></div><div style>Cool. If its not general enough, you can always write the loop yourself, or we can look into expanding it.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div dir="ltr"><div></div>What is the pStart and pEnd in the call to DMLabelCreateIndex referring to ? Also,my code complains if i use DMLabel :: labelName<br></div></blockquote><div><br></div><div style>It builds an index, so it needs bounds on the possible points you will be mapping.</div>
<div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">This name does not have a type, and must have an explicit type.</div></blockquote><div><br></div>
<div style>Ah, I did not put DMLabel into Fortran. Goes on the list.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>For now i can do without DMLabelCreateIndex as fortran interface to DMLabel doesn't seem to exist.</div>
</div></blockquote><div><br></div><div style>Yes, its only an optimization.</div><div style><br></div><div style> Thanks,</div><div style><br></div><div style> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div class="HOEnZb"><div class="h5"><div class="gmail_extra"><div class="gmail_quote">On Mon, Apr 1, 2013 at 11:02 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@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 dir="ltr"><div>On Tue, Apr 2, 2013 at 1:47 PM, Dharmendar Reddy <span dir="ltr"><<a href="mailto:dharmareddy84@gmail.com" target="_blank">dharmareddy84@gmail.com</a>></span> wrote:<br>
</div><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"><div dir="ltr"><div><div>Hello,<br></div><div>
I looked at the source code of SNESSetFunction after i sent the email. From what i understood, SNESSetFunction calls DMSNESSetFunction internally. I was like, why not just use SNESSetFunction with out passing DM. But i get it now, not all problems have to use DM. If i am using DM i should rather use DMSNESSetFunction. A better way is DMSNESSetFunctionLocal right ? DMSNesSetFunction i called inside DMSNESSetFunctionLocal.<br>
<br> <br></div></div><div>I did look at DMPlexComputeResidualFEM but that requires userContext to have PetscFEM object. Is PetscFEM object accisable from Fortran or is it genrated using<br></div></div></blockquote>
<div><br></div><div>
The object is just a struct, and has only information you must have for FEM, namely</div><div>a quadrature rule, and tabulation of the basis functions are quadrature points.</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"><pre width="80"><font color="#B22222">bin/pythonscripts/PetscGenerateFEMQuadrature.py</font></pre></div></blockquote><div>This script just generates a header with the struct information. You can make this header</div>
<div>yourself and obviate the need for the script.</div><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><div><div><div><div><div><div><div class="gmail_extra">I structured my code such that, user provides the Linear and bilinear form subroutines:<br><br></div><div class="gmail_extra">
subroutine biLinearFrom(u,v,cell,Kelm,eqnCtx)<br></div><div class="gmail_extra">type(FunctionSpace) :: u,v<br></div><div class="gmail_extra">type(FEMCell) :: cell<br></div><div class="gmail_extra">PetscScalar :: Kelm<br>
</div>
<div class="gmail_extra">class(*),pointer :: eqnCtx<br></div><div class="gmail_extra">end subroutine<br></div><div class="gmail_extra"><br></div><div class="gmail_extra">Function space object provide the routines for evaluation of basis functions and gradients, given cell coordinates, and cell degrees of freedom.</div>
</div></div></div></div></div></div></div></div></blockquote><div><br></div></div><div>Can you just use these to create the PetscFEM struct?</div><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><div><div><div><div><div><div><div class="gmail_extra"> <br>
So i should be able to get the cell geometry information and cell degree of freedom data via DMPlexVecGet/SetClosure right ?<br></div></div></div></div></div></div></div></div></div></blockquote><div><br></div></div><div>
Yes, although tell me why DMPlexComputeCellGeometry() does not work, and if not you can look at the code</div><div>to see me getting the coordinate information.</div><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><div><div><div><div><div><div><div class="gmail_extra"></div><div class="gmail_extra">At this point i was able to create a DMPlex object from a mesh created by gmsh. I used DMPLexCreateLabel to add the physical region labels i created in GMSh . So no i have a mapping from given a region label, list of cells in that region. <br>
<br></div><div class="gmail_extra">If i am doing a physical region wise assemble of the problem, i could use the DMPlexgetStratumIS for each regionlabel. <br></div></div></div></div></div></div></div></div></div></blockquote>
<div><br></div></div><div>Yes.</div><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><div><div><div><div><div><div><div class="gmail_extra"></div><div class="gmail_extra">since i am doing a cell wise assemble, i would be interested in a mapping from cell id to regionId (i can have a local array of labels indexed by region Id) so that i can pass the material information to FEMCell. Right now i have a integer array which holds the map from cellId to regionId.<br>
<br></div><div class="gmail_extra">I was wondering if i could add a integer section to DM object ? <br></div></div></div></div></div></div></div></div></div></blockquote><div><br></div></div><div>You use DMLabelGetValue(), and to trade memory for speed call DMLabelCreateIndex() first if you want.</div>
<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><div><div><div><div><div>
<div>
<div class="gmail_extra"></div><div class="gmail_extra">I also need to have some other auxiliary scalar data associated with the DM. You suggested in one of the earlier email that i created another DM using DMplexClone ? If i clone a object is it creating another copy of the mesh in memory or is it just storing pointers ? <br>
</div></div></div></div></div></div></div></div></div></blockquote><div><br></div></div><div>Just a pointer.</div><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><div><div><div><div><div><div><div class="gmail_extra">
</div><div class="gmail_extra">I was thinking i could create a section using DMPlexCreateSection and pass that section to DMPlexVecGetClosure when i have to access auxilary data. <br></div></div></div></div></div></div></div>
</div></div></blockquote><div><br></div></div><div>This is exactly what I was suggesting, but you have the section be the "default section" of the cloned DM. That way</div><div>it interacts with all PETSc functions nicely and you only have the extra memory for a DM object which is tiny.</div>
<div><br></div><div> Matt</div><div><div><div><br></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><div><div><div><div><div><div><div class="gmail_extra"></div><div class="gmail_extra">Basically the problem is to solve Poisson equation with current continuity equation for a semiconductor. I solve the problem by solving each equation until the solutions of both are consistent. For Poisson problem, the degree freedom is potential and takes charge density as input. For continuation equation degree of freedom is charge density and take potential as input. <br>
<br></div><div class="gmail_extra">Do you have alternative approaches to what i am thinking ? <br></div><div><div><div class="gmail_extra"><br></div><div class="gmail_extra"><div class="gmail_quote">On Mon, Apr 1, 2013 at 8:24 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>On Tue, Apr 2, 2013 at 10:16 AM, Dharmendar Reddy <span dir="ltr"><<a href="mailto:dharmareddy84@gmail.com" target="_blank">dharmareddy84@gmail.com</a>></span> wrote:<br>
</div><div class="gmail_extra"><div class="gmail_quote"><div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><div><div>Hello,<br></div> I am confused about the usage of SNESSetFunction, DMSNESSetFunction, DMSNESSetFunctionLocal and the corresponding Jacobin routines.<br>
<br></div>At present my code is set up as follows. <br>
<br></div>I initiate SNES with getFunction and getJacobain routines required as per SNESSet<Function,Jacobian>. I pass a user context which has the mesh and field layout information. I assemble the residual function and Jacobin element wise using the information from user context. My code at this point is serial, so i run the element loop from 1 to numberofTotalElements.<br>
<br></div><div>Now i want to switch to using DMPlex. How should i change the interface to my code. <br></div></div></blockquote><div><br></div></div><div>You can always look at SNES ex62 <a href="http://www.mcs.anl.gov/petsc/petsc-dev/src/snes/examples/tutorials/ex62.c.html" target="_blank">http://www.mcs.anl.gov/petsc/petsc-dev/src/snes/examples/tutorials/ex62.c.html</a></div>
<div>
<div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex">
<div dir="ltr"><div>If i set the dm for snes using SNESSetDM, i can see that the DM can be accessed via SNESGetDM inside getFunction and getJacobian, then i get the elementStartID and elementEndID to run the loop for assmbler. But now i will have to assemble to local vector into global vector right ? <br>
</div></div></blockquote><div><br></div></div><div>No, here is the sequence:</div><div><br></div><div> 1) SNESSetDM(snes, dm)</div><div><br></div><div> 2) DMSNESSetFunctionLocal(dm, userResidual, userCtx)</div>
<div><br></div><div>where we have</div><div><br></div><div> userResidual(DM dm, Vec X, Vec F, void *userCtx)</div><div><br></div><div>Both X and F are local vectors which I normally interact with using DMPlexVecGet/SetClosure(). If you</div>
<div>are using FEM, you can DMPlexComputeResidualFEM() here and use DMPlexSetFEMIntegration() to</div><div>input point-wise physics functions as is done in ex62.</div><div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex">
<div dir="ltr"><div></div><div>Looks like the DMSNESSetFunctionLocal will assemble the local <br></div><div>vector into global vector. The function provided to this should just evaluate the local vector. Did i understand this right ?<br clear="all">
</div><div><div><div><div><div><div><br></div><div>I am confused about passing the DM in DMSNESSet<Function,Jacobian>. When the DM can be accessed via SNESGetDM, why do we pass it again explicitly ?</div></div></div>
</div></div></div></div></blockquote><div><br></div></div><div>You are not passing the SNES, so where would it come from in this call?</div></div></div></div></blockquote><div><br><br> <br> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex">
<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><div><br></div><div> Matt</div><div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex">
<div dir="ltr"><div><div><div><div><div><div><span><font color="#888888"><br>
</font></span></div><span><font color="#888888"><div>-- <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" target="_blank">dharmareddy84@gmail.com</a><br>Phone: <a href="tel:%2B1-512-350-9082" value="+15123509082" target="_blank">+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>
</div></font></span></div></div></div></div></div></div>
</blockquote></div></div><span><font color="#888888"><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
</font></span></div></div>
</blockquote></div><br><br clear="all"><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" target="_blank">dharmareddy84@gmail.com</a><br>Phone: <a href="tel:%2B1-512-350-9082" value="+15123509082" target="_blank">+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>
</div></div></div></div></div></div></div></div></div></div></div>
</blockquote></div></div></div><div><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></div></div>
</blockquote></div><br><br clear="all"><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" target="_blank">dharmareddy84@gmail.com</a><br>Phone: <a href="tel:%2B1-512-350-9082" value="+15123509082" target="_blank">+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>
</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>