<html><head><meta http-equiv="Content-Type" content="text/html; charset=us-ascii"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><div class=""><br class=""></div><div class=""> Everybody is correct.</div><div class=""><br class=""></div> Using the extremely sparse, nearly diagonal MPI matrix will have a bit more overhead than doing the rotations on the elements due to the needed MPI check and general use of a sparse-sparse-sparse matrix vector product. But if a code does not have support for the local transformations the extra cost of the the sparse matrix operation will be trivial relative to the rest of computations and so is fine to use. <div class=""><br class=""></div><div class=""> Barry</div><div class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On May 31, 2021, at 1:00 PM, Sanjay Govindjee <<a href="mailto:s_g@berkeley.edu" class="">s_g@berkeley.edu</a>> wrote:</div><br class="Apple-interchange-newline"><div class="">
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8" class="">
<div class="">
In our FEA code we perform the rotations at the local level, before
assembly so that it is easy to apply the boundary conditions, then
unrotate locally after solution to get the usual Cartesian
components. Somehow this seems more efficient than doing this
globally, but perhaps I am missing something.<br class="">
-sanjay<br class="">
<pre class="moz-signature" cols="72">
</pre>
<div class="moz-cite-prefix">On 5/31/21 9:33 AM, Matthew Knepley
wrote:<br class="">
</div>
<blockquote type="cite" cite="mid:CAMYG4GmjmbTg8_6TU9FGEfkwooRrfeJtiwtxaYwg-GDAQF2MAQ@mail.gmail.com" class="">
<meta http-equiv="content-type" content="text/html; charset=UTF-8" class="">
<div dir="ltr" class="">
<div dir="ltr" class="">On Mon, May 31, 2021 at 11:12 AM Stefano Zampini
<<a href="mailto:stefano.zampini@gmail.com" moz-do-not-send="true" class="">stefano.zampini@gmail.com</a>>
wrote:<br class="">
</div>
<div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left:1px solid
rgb(204,204,204);padding-left:1ex">
<div dir="ltr" class="">Mike
<div class=""><br class="">
</div>
<div class="">as long as P is a sparse matrix with compatible rows
and cols (i.e. rows(P)= cols(A) = rows (A)) , MatPtAP
will compute the result. </div>
</div>
</blockquote>
<div class=""><br class="">
</div>
<div class="">Stefano and Mark are correct. This will work.</div>
<div class=""><br class="">
</div>
<div class="">I implemented the same thing in my code in a different
way. I put this transformation into the mapping between
local and global vector spaces. The global degrees of</div>
<div class="">freedom are the ones you want for boundary conditions
(normal and tangential to the boundary), and I eliminate the
ones that are constrained. The local degrees of</div>
<div class="">freedom are the normal Caresian ones, and these are used
for assembly. The map is used when I
execute DMGlobalToLocal() and DMLocalToGlobal(). There is an</div>
<div class="">example of me doing this in SNES ex71, Poiseuille flow in
a tilted channel.</div>
<div class=""><br class="">
</div>
<div class=""> Thanks,</div>
<div class=""><br class="">
</div>
<div class=""> Matt</div>
<div class=""> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left:1px solid
rgb(204,204,204);padding-left:1ex">
<div class="gmail_quote">
<div dir="ltr" class="gmail_attr">Il giorno lun 31 mag
2021 alle ore 16:52 Mark Adams <<a href="mailto:mfadams@lbl.gov" target="_blank" moz-do-not-send="true" class="">mfadams@lbl.gov</a>> ha
scritto:<br class="">
</div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px
0.8ex;border-left:1px solid
rgb(204,204,204);padding-left:1ex">
<div dir="ltr" class="">
<div dir="ltr" class=""><br class="">
</div>
<br class="">
<div class="gmail_quote">
<div dir="ltr" class="gmail_attr">On Mon, May 31,
2021 at 9:20 AM Michael Wick <<a href="mailto:michael.wick.1980@gmail.com" target="_blank" moz-do-not-send="true" class="">michael.wick.1980@gmail.com</a>>
wrote:<br class="">
</div>
<blockquote class="gmail_quote" style="margin:0px
0px 0px 0.8ex;border-left:1px solid
rgb(204,204,204);padding-left:1ex">
<div dir="ltr" class="">Hi PETSc team:
<div class=""><br class="">
</div>
<div class="">I am considering implementing a skew roller
boundary condition for my elasticity problem.
The method is based on this journal paper: <a href="http://inside.mines.edu/~vgriffit/pubs/All_J_Pubs/18.pdf" target="_blank" moz-do-not-send="true" class="">http://inside.mines.edu/~vgriffit/pubs/All_J_Pubs/18.pdf</a></div>
<div class=""><br class="">
</div>
<div class="">Or you may find the method in the attached
Bathe's slides, pages 9 -10.</div>
<div class=""><br class="">
</div>
<div class="">Roughly speaking, a (very) sparse matrix T
will be created which takes the shape [ I, O;
O, R], where R is a 3x3 rotation matrix. And
the original linear problem K U = F will be
modified into (T^t K T) (T^t U) = T^t F. In
doing so, one can enforce a roller boundary
condition on a slanted surface.</div>
<div class=""><br class="">
</div>
<div class="">I think it can be an easy option if I can
generate the T matrix and do two matrix
multiplications to get T^t K T. I noticed that
there is a MatPtAP function. Yet, after
reading a previous discussion, it seems that
this function is not designed for this
purposes (<a href="https://lists.mcs.anl.gov/pipermail/petsc-users/2018-June/035477.html" target="_blank" moz-do-not-send="true" class="">https://lists.mcs.anl.gov/pipermail/petsc-users/2018-June/035477.html</a>).</div>
</div>
</blockquote>
<div class=""><br class="">
</div>
<div class="">Yes, and no. It is motivated and optimized for
a Galerkin coarse grid operator for AMG solvers,
but it is a projection and it should be fine. If
not, we will fix it.</div>
<div class=""><br class="">
</div>
<div class="">We try to test our methods of "empty" operators
, but I don't know if MatPtAP has ever been tested
for super sparse P. Give it a shot and see what
happens.</div>
<div class=""><br class="">
</div>
<div class="">Mark</div>
<div class=""> </div>
<blockquote class="gmail_quote" style="margin:0px
0px 0px 0.8ex;border-left:1px solid
rgb(204,204,204);padding-left:1ex">
<div dir="ltr" class="">
<div class=""><br class="">
</div>
<div class="">I assume I can only call MatMatMult &
MatTransposeMatMult to do this job, correct?
Is there any existingly PETSc function to do
T^t K T in one call?</div>
<div class=""><br class="">
</div>
<div class="">Thanks,</div>
<div class=""><br class="">
</div>
<div class="">Mike</div>
<div class=""><br class="">
</div>
</div>
</blockquote>
</div>
</div>
</blockquote>
</div>
<br clear="all" class="">
<div class=""><br class="">
</div>
-- <br class="">
<div dir="ltr" class="">Stefano</div>
</blockquote>
</div>
<br clear="all" class="">
<div class=""><br class="">
</div>
-- <br class="">
<div dir="ltr" class="gmail_signature">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">
<div dir="ltr" class="">
<div class="">What most experimenters take for granted before
they begin their experiments is infinitely more
interesting than any results to which their
experiments lead.<br class="">
-- Norbert Wiener</div>
<div class=""><br class="">
</div>
<div class=""><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank" moz-do-not-send="true" class="">https://www.cse.buffalo.edu/~knepley/</a><br class="">
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</blockquote>
<br class="">
</div>
</div></blockquote></div><br class=""></div></body></html>