<div>A bit more information:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Looking at the solution of the system of equations, knowing the exact RHS vector:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;">[-0.00197, 0.00203, -0.00197]</span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">SNES gives the solution:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 3.316262462871E-03, -4.189244774965E-03, 2.468317119413E-03</span><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">which is indeed (closer to) the solution for the jacobian shown in MatView (the "lower precision" one). The "exact" solution would be:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---<br><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">3.282607248309093b-3, -4.241572005990172b-3,, 2.425030835360701b-3</span></div>---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">which is already different form the second decimal digit.</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">***</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Requesting a low precision to have SNES do a few iterations, shows that the built jacobian seems to deviate more and more form the computed one after some iterations, e.g.:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Computed (printed from in-code):</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.536488494840 0.066474939660 0.198813816284</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.066474939660 0.529482312614 -0.002835445071</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.198813816284 -0.002835445071 0.530589055263</span></div><span></span></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">MatView:</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span>row 0: </span>(0, 0.538491) (1, 0.0662019) (2, 0.198808) </span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span>row 1: </span>(0, 0.0662019) (1, 0.526981) (2, -0.00257554) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"><span>row 2: </span>(0, 0.198808) (1, -0.00257554) (2, 0.529883) </span></div><span></span>---</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Only one digit is correct/the same in entries [1,2], [2,1] and [2,2].</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Thank you.</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);">Noam</div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px; color: rgb(0, 0, 0); background-color: rgb(255, 255, 255);"><br></div><div class="protonmail_quote">
On Sunday, March 15th, 2026 at 3:29 PM, Noam T. <dontbugthedevs@proton.me> wrote:<br>
<blockquote class="protonmail_quote" type="cite">
<div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Hello,</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Looking at the assembly of the Jacobian using MatView, I noticed values were somewhat different from the computed ones in the jacobian function (DMSNESSetJacobian). For the test, I used a single Q1 element, with 4 quadrature points, so the exact Jacobian matrix can be computed analytically. </span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">In the jacobian function, after looping through all quadrature points, printing the matrix to stdout shows:</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">----</span><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.531360000000 0.066670000000 0.197333333327</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.066670000000 0.535360000000 -0.003333333327</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> 0.197333333327 -0.003333333327 0.527413333333</span></div><span></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">----</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">(showing only the 3 free DoFs). The 2x2 submatrix in [0:1] has exactly four/five nonzero digits, as shown above. The rest of the elements have a periodic 3 ending (a couple of digits at the end are off, but that's fine).</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">This same matrix is the one added to the global jacobian during the quadrature loop, with "</span><span style="font-family: Times New Roman, serif;">DMPlexMatSetClosure". After looping, calling "MatAssemblyBegin/End", then MatView: </span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><span style="font-family: Times New Roman, serif;">----</span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><span style="font-family: Menlo, Consolas, Courier New, monospace;">Mat Object: 1 MPI process</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> type: seqaij</span></div></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 0: (0, 0.531387) (1, 0.0666733) (2, 0.197338) </span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 1: (0, 0.0666733) (1, 0.535387) (2, -0.0033385) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 2: (0, 0.197338) (1, -0.0033385) (2, 0.52744) </span></div><span></span><span style="font-family: Times New Roman, serif;">----</span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Values are close, but definitely not the same as computed. Are these values from "MatView" the ones actually stored in the matrix? It seems there is some precision loss in the process.</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Interestingly, computing the jacobian through finite differences, shows the same result for the "hand-coded" jacobian, whereas the finite differences ones is the "exact" one (same as first one above):</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">----</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><span style="font-family: Menlo, Consolas, Courier New, monospace;"> ||J - Jfd||_F/||J||_F = 4.90736e-05, ||J - Jfd||_F = 4.74265e-05</span><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> Hand-coded Jacobian ----------</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">Mat Object: 1 MPI process</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> type: seqaij</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 0: (0, 0.531387) (1, 0.0666733) (2, 0.197338) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 1: (0, 0.0666733) (1, 0.535387) (2, -0.0033385) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 2: (0, 0.197338) (1, -0.0033385) (2, 0.52744) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> Finite difference Jacobian ----------</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">Mat Object: 1 MPI process</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> type: seqaij</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 0: (0, 0.53136) (1, 0.06667) (2, 0.197333) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 1: (0, 0.06667) (1, 0.53536) (2, -0.00333333) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 2: (0, 0.197333) (1, -0.00333333) (2, 0.527413) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> Hand-coded minus finite-difference Jacobian with tolerance 1e-05 ----------</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">Mat Object: 1 MPI process</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;"> type: seqaij</span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 0: (0, 2.66554e-05) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 1: (1, 2.66554e-05) </span></div><div><span style="font-family: Menlo, Consolas, Courier New, monospace;">row 2: (2, 2.6558e-05) </span></div><span></span></span><span style="font-family: Times New Roman, serif;"> ----</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">The main code is in Fortran, using double precision variables, in case it is relevant.</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Thank you.</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span><br></span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><span style="font-family: Times New Roman, serif;">Noam</span></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div><div style="font-family: Arial, sans-serif; font-size: 14px;"><br></div>
</blockquote><br>
</div>