<html xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<meta name="Generator" content="Microsoft Word 15 (filtered medium)">
<style><!--
/* Font Definitions */
@font-face
        {font-family:"Cambria Math";
        panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
        {font-family:Calibri;
        panose-1:2 15 5 2 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
        {margin:0in;
        margin-bottom:.0001pt;
        font-size:11.0pt;
        font-family:"Calibri",sans-serif;}
a:link, span.MsoHyperlink
        {mso-style-priority:99;
        color:blue;
        text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
        {mso-style-priority:99;
        color:purple;
        text-decoration:underline;}
p.msonormal0, li.msonormal0, div.msonormal0
        {mso-style-name:msonormal;
        mso-margin-top-alt:auto;
        margin-right:0in;
        mso-margin-bottom-alt:auto;
        margin-left:0in;
        font-size:11.0pt;
        font-family:"Calibri",sans-serif;}
span.EmailStyle18
        {mso-style-type:personal-reply;
        font-family:"Calibri",sans-serif;
        color:windowtext;}
.MsoChpDefault
        {mso-style-type:export-only;
        font-size:10.0pt;}
@page WordSection1
        {size:8.5in 11.0in;
        margin:1.0in 1.0in 1.0in 1.0in;}
div.WordSection1
        {page:WordSection1;}
--></style>
</head>
<body lang="EN-US" link="blue" vlink="purple">
<div class="WordSection1">
<p class="MsoNormal">What about ‘bending’ the physics a bit?<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Q_i-1/2 proportional to sqrt(x[i-1] + z[i-1] - (x[i] + z[i])),                  if             x[i-1]+z[i-1]  >  x[i]+z[i] + eps<o:p></o:p></p>
<p class="MsoNormal">Q_i-1/2 proportional to -1.0*sqrt(x[i] + z[i] - (x[i-1] + z[i-1])),         if            x[i]+z[i]  >  x[i-1]+z[i-1] + eps<o:p></o:p></p>
<p class="MsoNormal">Q_i-1/2 proportional to x[i-1] + z[i-1] - (x[i] + z[i])                               in between<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">in which, eps is a very small positive number.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">-Ling<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<div style="border:none;border-top:solid #B5C4DF 1.0pt;padding:3.0pt 0in 0in 0in">
<p class="MsoNormal"><b><span style="font-size:12.0pt;color:black">From: </span></b><span style="font-size:12.0pt;color:black">petsc-users <petsc-users-bounces@mcs.anl.gov> on behalf of Alexander B Prescott <alexprescott@email.arizona.edu><br>
<b>Date: </b>Monday, April 6, 2020 at 1:06 PM<br>
<b>To: </b>PETSc <petsc-users@mcs.anl.gov><br>
<b>Subject: </b>[petsc-users] Discontinuities in the Jacobian matrix for nonlinear problem<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Hello, <o:p></o:p></p>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">The non-linear boundary-value problem I am applying PETSc to is a relatively simple steady-state flow routing algorithm based on the continuity equation, such that Div(Q) = 0 everywhere (Q=discharge). I use a finite volume approach to calculate
 flow between nodes, with Q calculated as a piecewise smooth function of the local flow depth and the water-surface slope. In 1D, the residual is calculated as R(x_i)=Q_i-1/2 - Q_i+1/2.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">For example, Q_i-1/2 at x[i]:<o:p></o:p></p>
</div>
<div>
<blockquote style="margin-left:30.0pt;margin-right:0in">
<div>
<p class="MsoNormal">Q_i-1/2 proportional to sqrt(x[i-1] + z[i-1] - (x[i] + z[i])),              if  x[i-1]+z[i-1]  >  x[i]+z[i]<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">Q_i-1/2 proportional to -1.0*sqrt(x[i] + z[i] - (x[i-1] + z[i-1])),      if         x[i]+z[i]  >  x[i-1]+z[i-1] <o:p></o:p></p>
</div>
</blockquote>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Where z[i] is local topography and doesn't change over the iterations, and Q_i+1/2 is computed analogously. So the residual derivatives with respect to x[i-1], x[i] and x[i+1] are not continuous when the water-surface slope = 0.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Are there intelligent ways to handle this problem? My 1D trial runs naively fix any zero-valued water-surface slopes to a small non-zero positive value (e.g. 1e-12). Solver convergence has been mixed and highly dependent on the initial
 guess. So far, FAS with QN coarse solver has been the most robust.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Restricting x[i] to be non-negative is a separate issue, to which I have applied the SNES_VI solvers. They perform modestly but have been less robust.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Best,<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">Alexander<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal"><br clear="all">
<o:p></o:p></p>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<p class="MsoNormal">-- <o:p></o:p></p>
<div>
<div>
<div>
<div>
<div>
<div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif">Alexander Prescott</span><span style="font-size:9.5pt"><o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif"><a href="mailto:alexprescott@email.arizona.edu" target="_blank">alexprescott@email.arizona.edu</a></span><span style="font-size:9.5pt"><o:p></o:p></span></p>
</div>
<div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif">PhD Candidate, The University of Arizona<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif">Department of Geosciences<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif">1040 E. 4th Street<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><span style="font-size:9.5pt;font-family:"Arial",sans-serif">Tucson, AZ, 85721<o:p></o:p></span></p>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</body>
</html>