<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Apr 22, 2015 at 9:57 PM, Justin Chang <span dir="ltr"><<a href="mailto:jchang27@uh.edu" target="_blank">jchang27@uh.edu</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><div><div><div>Hello,<br><br></div>I am using DMPlexProjectField(...) to postprocess my final solution u into another vector sol. I have the following pointwise functions:<br><br></div><div>/* x-component of gradient */<br></div><div>void velx(const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], const PetscReal x[], PetscScalar v[])<br>{<br>  v[0] = -a[0]*mu*u_x[0];<br>}<br><br>/* y-component of gradient */<br>void vely(const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], const PetscReal x[], PetscScalar v[])<br>{<br>  v[0] = -a[0]*mu*u_x[1];<br>}<br><br>/* z-component of gradient */<br>void velz(const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], const PetscReal x[], PetscScalar v[])<br>{<br>  v[0] = -a[0]*mu*u_x[2];<br>}<br><br></div>Where a[0] is the cell-wise permeability and mu is the (inverse of) viscosity. I am solving the diffusion equation (P1 elements). At the end of my simulation I have these lines:<br><br>ierr = VecDuplicate(user.u,&user.sol);CHKERRQ(ierr);<br>ierr = PetscObjectSetName((PetscObject) user.sol, "velocity");CHKERRQ(ierr);<br>void (*vX[1])(const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscReal *, PetscScalar *) = {velx};<br>ierr = DMPlexProjectField(dm,user.u,vX,INSERT_ALL_VALUES,user.sol);CHKERRQ(ierr);<br>ierr = VecViewFromOptions(user.sol,NULL,"-vec_view_x");CHKERRQ(ierr);<br><br>void (*vY[1])(const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscReal *, PetscScalar *) = {vely};<br>ierr = DMPlexProjectField(dm,user.u,vY,INSERT_ALL_VALUES,user.sol);CHKERRQ(ierr);<br>ierr = VecViewFromOptions(user.sol,NULL,"-vec_view_y");CHKERRQ(ierr);<br><br>if (spatialDim == 3) {<br>  void (*vZ[1])(const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscScalar *, const PetscReal *, PetscScalar *) = {velz};<br>  ierr = DMPlexProjectField(dm,user.u,vZ,INSERT_ALL_VALUES,user.sol);CHKERRQ(ierr);<br>  ierr = VecViewFromOptions(user.sol,NULL,"-vec_view_z");CHKERRQ(ierr);<br>}<br>ierr = VecDestroy(&user.sol);CHKERRQ(ierr);<br><br></div>But when i look at the results, the dirichlet constrained values of u were projected into the respective locations in sol. I thought putting INSERT_ALL_VALUES would insert the pointwise functions into all locations including the constrained ones. How can I project the above pointwise functions into the constrained fields? Or what am I potentially missing here?<br></div></div></blockquote><div><br></div><div>There are no Dirichlet constrained values in the global vector user.sol. When you call VecView() it gets projected to</div><div>an enlarged vector including the BC and output. What you want is to project that vector into a space with no constraints,</div><div>so you really need another DM/Section combination.</div><div><br></div><div>I think this will work:</div><div><br></div><div>  1) Get a local vector from dm</div><div><br></div><div>  2) Use DMPlexProjectFieldLocal(INSERT_ALL_VALUES) to fill it up</div><div><br></div><div>  3) DMClone(dm) and give it a Section with no constraints</div><div><br></div><div>  4) Use DMLocalToGlobal() into the new DM</div><div><br></div><div>  5) VecView() that new global vector</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</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>Thanks,<span class="HOEnZb"><font color="#888888"><br><br clear="all"><br>-- <br><div><div dir="ltr"><div><div><div>Justin Chang<br></div>PhD Candidate, Civil Engineering - Computational Sciences<br></div>University of Houston, Department of Civil and Environmental Engineering<br></div>Houston, TX 77004<br><a href="tel:%28512%29%20963-3262" value="+15129633262" target="_blank">(512) 963-3262</a><br></div></div></font></span></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">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>