<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=us-ascii">
<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:Aptos;
panose-1:2 11 0 4 2 2 2 2 2 4;}
@font-face
{font-family:Consolas;
panose-1:2 11 6 9 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
{margin:0cm;
font-size:11.0pt;
font-family:"Aptos",sans-serif;
mso-ligatures:standardcontextual;
mso-fareast-language:EN-US;}
span.EmailStyle17
{mso-style-type:personal-compose;
font-family:"Aptos",sans-serif;
color:windowtext;}
.MsoChpDefault
{mso-style-type:export-only;
font-size:11.0pt;
mso-fareast-language:EN-US;}
@page WordSection1
{size:612.0pt 792.0pt;
margin:72.0pt 72.0pt 72.0pt 72.0pt;}
div.WordSection1
{page:WordSection1;}
--></style>
</head>
<body lang="en-CH" link="#467886" vlink="#96607D" style="word-wrap:break-word">
<div class="WordSection1">
<p class="MsoNormal"><span lang="EN-US">Hello,<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"><o:p> </o:p></span></p>
<p class="MsoNormal"><span lang="EN-US">I have this simple example on Firedrake to illustrate my point. I am solving for a two-component poisson equation (uncoupled). Only the first component has a non-zero residual.<br>
<br>
```<br>
</span>import firedrake as fd<o:p></o:p></p>
<p class="MsoNormal">from firedrake import inner, grad, dx, sin, pi<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">N = 10<o:p></o:p></p>
<p class="MsoNormal">mesh = fd.UnitSquareMesh(N, N)<o:p></o:p></p>
<p class="MsoNormal">V = fd.FunctionSpace(mesh, "CG", 1)<o:p></o:p></p>
<p class="MsoNormal">W = V * V<o:p></o:p></p>
<p class="MsoNormal">u = fd.Function(W)<o:p></o:p></p>
<p class="MsoNormal">v = fd.TestFunction(W)<o:p></o:p></p>
<p class="MsoNormal">a = inner(grad(u[0]), grad(v[0])) * dx + inner(grad(u[1]), grad(v[1])) * dx<o:p></o:p></p>
<p class="MsoNormal">x = fd.SpatialCoordinate(mesh)<o:p></o:p></p>
<p class="MsoNormal">f = fd.Function(V)<o:p></o:p></p>
<p class="MsoNormal">f.interpolate(fd.Constant(1e4) * sin(x[0] * pi) * sin(2 * x[1] * pi))<o:p></o:p></p>
<p class="MsoNormal">L = f * v[0] * dx<o:p></o:p></p>
<p class="MsoNormal">F = a - L<o:p></o:p></p>
<p class="MsoNormal">bcs = [fd.DirichletBC(W.sub(0), fd.Constant(2.0), (1,))]<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">def snes_firedrake_residual(F, u, bcs):<o:p></o:p></p>
<p class="MsoNormal"> for bcs_ in bcs:<o:p></o:p></p>
<p class="MsoNormal"> bcs_.apply(u)<o:p></o:p></p>
<p class="MsoNormal"> residual = fd.assemble(F, bcs=bcs, zero_bc_nodes=True)<o:p></o:p></p>
<p class="MsoNormal"> with residual.dat.vec_ro as r:<o:p></o:p></p>
<p class="MsoNormal"> print("Initial residual:", r.norm())<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">snes_firedrake_residual(F, u, bcs)<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">problem = fd.NonlinearVariationalProblem(F, u, bcs=bcs)<o:p></o:p></p>
<p class="MsoNormal">solver_mumps_assembled = {<o:p></o:p></p>
<p class="MsoNormal"> "ksp_type": "preonly",<o:p></o:p></p>
<p class="MsoNormal"> "ksp_monitor": None,<o:p></o:p></p>
<p class="MsoNormal"> "pc_type": "python",<o:p></o:p></p>
<p class="MsoNormal"> "pc_python_type": "firedrake.AssembledPC",<o:p></o:p></p>
<p class="MsoNormal"> "assembled_pc_type": "lu",<o:p></o:p></p>
<p class="MsoNormal"> "assembled_pc_factor_mat_solver_type": "mumps",<o:p></o:p></p>
<p class="MsoNormal"> "assembled_mat_mumps_icntl_14": 200,<o:p></o:p></p>
<p class="MsoNormal"> "assembled_mat_mumps_icntl_24": 1,<o:p></o:p></p>
<p class="MsoNormal">}<o:p></o:p></p>
<p class="MsoNormal">solver_fieldsplit = {<o:p></o:p></p>
<p class="MsoNormal"> "mat_type": "matfree",<o:p></o:p></p>
<p class="MsoNormal"> "snes_type": "newtonls",<o:p></o:p></p>
<p class="MsoNormal"> "ksp_type": "fgmres",<o:p></o:p></p>
<p class="MsoNormal"> "ksp_rtol": 1e-1,<o:p></o:p></p>
<p class="MsoNormal"> "ksp_monitor": None,<o:p></o:p></p>
<p class="MsoNormal"> "pc_type": "fieldsplit",<o:p></o:p></p>
<p class="MsoNormal"> "pc_fieldsplit_type": "additive",<o:p></o:p></p>
<p class="MsoNormal"> "fieldsplit_0": solver_mumps_assembled,<o:p></o:p></p>
<p class="MsoNormal"> "fieldsplit_1": solver_mumps_assembled,<o:p></o:p></p>
<p class="MsoNormal">}<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">solver = fd.NonlinearVariationalSolver(problem, solver_parameters=solver_fieldsplit)<o:p></o:p></p>
<p class="MsoNormal">solver.solve()<o:p></o:p></p>
<p class="MsoNormal"><span lang="EN-US">```<br>
<br>
The PETSc output is as follow<br>
<br>
```<br>
Initial residual: 462.13689530272404<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> 0 SNES Function norm 4.621368953027e+02<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> Residual norms for firedrake_0_ solve.<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> 0 KSP Residual norm 4.621368953027e+02<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> Residual norms for firedrake_0_fieldsplit_0_ solve.<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> </span><span lang="IT-CH">0 KSP Residual norm 1.000000000000e+00<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="IT-CH"> 1 KSP Residual norm 3.501082228626e-15<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="IT-CH"> </span><span lang="EN-US">Residual norms for firedrake_0_fieldsplit_1_ solve.<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> </span><span lang="IT-CH">0 KSP Residual norm 0.000000000000e+00<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="IT-CH"> 1 KSP Residual norm 0.000000000000e+00<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="IT-CH"> </span><span lang="EN-US">1 KSP Residual norm 1.612167203819e-12<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"> 1 SNES Function norm 1.599350481360e-12<br>
```<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US"><o:p> </o:p></span></p>
<p class="MsoNormal">Using the fieldsplit additive preconditioner, the problem converges in a single KSP iteration, as expected. However, I do not understand
<span lang="EN-US">why the</span> residual of fieldsplit_0 (1e+0) does not coincide with the outer residual (462.13689530272404)<span lang="EN-US">. It should be the case given that only
</span>fieldsplit_0 has<span lang="EN-US"> a non-zero residual contribution</span>. The fact that it is just 1 is suspicious. Is there something about how the fieldsplit works that I am missing?<span lang="EN-US"><br>
<br>
Thanks,<o:p></o:p></span></p>
<p class="MsoNormal"><span lang="EN-US">Miguel<o:p></o:p></span></p>
</div>
<p style="FONT-SIZE: 10pt; FONT-FAMILY: Arial">
</p>
<table style="WIDTH: 409px" cellspacing="0" cellpadding="0" border="0">
<tbody>
<tr>
<td style="FONT-SIZE: 10pt; FONT-FAMILY: Arial; WIDTH: 260px; PADDING-BOTTOM: 20px" valign="top" class=""><br></td><td style="FONT-SIZE: 10pt; FONT-FAMILY: Arial; WIDTH: 149px; PADDING-BOTTOM: 20px" valign="top" class=""><br></td></tr><tr><td valign="top" colspan="1" style="width: 260px; font-family: Arial; color: rgb(138, 138, 140); font-size: 8pt; line-height: 14pt; padding: 0px 0px 15px;" class=""><span data-codetwo-visible="Name"><b style="COLOR: #1b1464">MIGUEL ANGEL SALAZAR DE TROYA<br></b></span><span style="FONT-SIZE: 8pt; COLOR: #8a8a8c"><span data-codetwo-visible="Job title">Head of Software Engineering</span></span><br><span data-codetwo-visible="Email"><a style="FONT-SIZE: 8pt; TEXT-DECORATION: underline; COLOR: #8a8a8c" href="mailto:miguel.salazar@corintis.com"><em>miguel.salazar@corintis.com</em></a><br></span><span data-codetwo-visible="Phone">Corintis SA<br></span><span data-codetwo-visible="Company address">EPFL Innovation Park Building C<br>1015 Lausanne<br></span></td><td valign="top" colspan="1" style="width: 149px; font-family: Arial; line-height: 14pt; padding: 0px 0px 15px;" class=""><img src="cid:2024-08-1609_02_11-re_mailmigrationfromgoogletooffice365-sebastien.gobel@corintis.com-co_4456a1bf-dd3e-46b7-85f9-21f278e66a79.png" border="0" id="0.k3axcrzm7nf" alt="www.corintis.com" style="width: 247px; height: 97px;" width="247" height="97"></td>
</tr>
<tr data-codetwo-visible="Disclaimer">
<td style="FONT-SIZE: 7pt; FONT-FAMILY: Arial; WIDTH: 409px; COLOR: #ffffff; PADDING-BOTTOM: 3px; PADDING-TOP: 3px; PADDING-LEFT: 10px; PADDING-RIGHT: 10px; BACKGROUND-COLOR: #1b1464" valign="top" colspan="2">Here at Corintis we care for your privacy. That is why we have taken appropriate
measures to ensure that the data you have provided to us is always secure.</td>
</tr>
</tbody>
</table>
</body>
</html>