<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Jun 14, 2017 at 11:48 PM, Adrian Croucher <span dir="ltr"><<a href="mailto:a.croucher@auckland.ac.nz" target="_blank">a.croucher@auckland.ac.nz</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">On 14/06/17 07:45, Jed Brown wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Barry Smith <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> writes:<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
On Jun 13, 2017, at 10:06 AM, Jed Brown <<a href="mailto:jed@jedbrown.org" target="_blank">jed@jedbrown.org</a>> wrote:<br>
<br>
Adrian Croucher <<a href="mailto:a.croucher@auckland.ac.nz" target="_blank">a.croucher@auckland.ac.nz</a>> writes:<br>
<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
One way might be to form the whole Jacobian but somehow use a modified<br>
KSP solve which would implement the reduction process, do a KSP solve on<br>
the reduced system of size n, and finally back-substitute to find the<br>
unknowns in the matrix rock cells.<br>
</blockquote>
You can do this with PCFieldSplit type Schur, but it's a lot heavier<br>
than you might like.<br>
</blockquote>
    Is it clear that it would produce much overhead compared to doing a custom "reduction to a smaller problem". Perhaps he should start with this and then profiling can show if there are any likely benefits to "specializing more"?<br>
</blockquote>
Yeah, that would be reasonable.  We don't have a concept of sparsity for<br>
preconditioners so don't have a clean way to produce the exact (sparse)<br>
Schur complement.  Computing this matrix using coloring should be<br>
relatively inexpensive due to the independence in each cell and its<br>
tridiagonal structure.<br>
</blockquote>
<br>
Thanks for those ideas, very helpful.<br>
<br>
If I try this approach (forming whole Jacobian matrix and using PCFieldSplit Schur), I guess I will first need to set up a modified DMPlex for the whole fracture + matrix mesh- so I can use it to create vectors and the Jacobian matrix (with the right sparsity pattern), and also to work out the coloring for finite differencing.<br>
<br>
Would that be straight-forward to do? Currently my DM just comes from DMPlexCreateFromFile(). Presumably you can use DMPlexInsertCone() or similar to add points into it?</blockquote><div><br></div><div>You can certainly modify the mesh. I need to get a better idea what kind of modification, and then</div><div>I can suggest a way to do it. What do you start with, and what exactly do you want to add?</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"><span class="HOEnZb"><font color="#888888"><br>
- Adrian<br>
<br>
-- <br>
Dr Adrian Croucher<br>
Senior Research Fellow<br>
Department of Engineering Science<br>
University of Auckland, New Zealand<br>
email: <a href="mailto:a.croucher@auckland.ac.nz" target="_blank">a.croucher@auckland.ac.nz</a><br>
tel: <a href="tel:%2B64%20%280%299%20923%204611" value="+6499234611" target="_blank">+64 (0)9 923 4611</a><br>
<br>
</font></span></blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div>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><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">http://www.caam.rice.edu/~mk51/</a><br></div></div></div>
</div></div>