<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> <br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Sep 29, 2021, at 5:37 PM, Ramakrishnan Thirumalaisamy <<a href="mailto:rthirumalaisam1857@sdsu.edu" class="">rthirumalaisam1857@sdsu.edu</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class="">Hi all,<div class=""><br class=""></div><div class="">I am trying to solve the Helmholtz equation for temperature T:</div><div class=""><br class=""></div><div class="">(C I + Div D grad) T = f</div><div class=""><br class=""></div><div class="">in IBAMR, in which C is the spatially varying diagonal entries, and D is the spatially varying diffusion coefficient. I use a matrix-free solver with matrix-based PETSc preconditioner. For the matrix-free solver, I use gmres solver and for the matrix based preconditioner, I use Richardson ksp + Jacobi as a preconditioner. As the simulation progresses, the iterations start to increase. To understand the cause, I set D to be zero, which results in a diagonal system:</div><div class=""> </div><div class="">C T = f.</div><div class=""><br class=""></div><div class="">This should result in convergence within a single iteration, but I get convergence in 3 iterations.</div><div class=""><br class=""></div><div class=""><div style="margin: 0px; font-stretch: normal; font-size: 15px; line-height: normal; font-family: Monaco; color: rgb(242, 242, 242); background-color: rgba(0, 0, 0, 0.85);" class=""><span style="font-variant-ligatures:no-common-ligatures" class="">Residual norms for temperature_ solve.</span></div><div style="margin: 0px; font-stretch: normal; font-size: 15px; line-height: normal; font-family: Monaco; color: rgb(242, 242, 242); background-color: rgba(0, 0, 0, 0.85);" class=""><span style="font-variant-ligatures:no-common-ligatures" class=""><span class="gmail-Apple-converted-space"> </span>0 KSP preconditioned resid norm 4.590811647875e-02 true resid norm 2.406067589273e+09 ||r(i)||/||b|| 4.455533946945e-05</span></div><div style="margin: 0px; font-stretch: normal; font-size: 15px; line-height: normal; font-family: Monaco; color: rgb(242, 242, 242); background-color: rgba(0, 0, 0, 0.85);" class=""><span style="font-variant-ligatures:no-common-ligatures" class=""><span class="gmail-Apple-converted-space"> </span>1 KSP preconditioned resid norm 2.347767895880e-06 true resid norm 1.210763896685e+05 ||r(i)||/||b|| 2.242081505717e-09</span></div><div style="margin: 0px; font-stretch: normal; font-size: 15px; line-height: normal; font-family: Monaco; color: rgb(242, 242, 242); background-color: rgba(0, 0, 0, 0.85);" class=""><span style="font-variant-ligatures:no-common-ligatures" class=""><span class="gmail-Apple-converted-space"> </span>2 KSP preconditioned resid norm 1.245406571896e-10 true resid norm 6.328828824310e+00 ||r(i)||/||b|| 1.171966730978e-13</span></div><div style="margin: 0px; font-stretch: normal; font-size: 15px; line-height: normal; font-family: Monaco; color: rgb(242, 242, 242); background-color: rgba(0, 0, 0, 0.85);" class=""><span style="font-variant-ligatures:no-common-ligatures" class="">Linear temperature_ solve converged due to CONVERGED_RTOL iterations 2</span></div></div><div class=""><br class=""></div></div></div></blockquote><div><br class=""></div> What is the result of -ksp_view on the solve? </div><div><br class=""></div><div> The way you describe your implementation it does not sound like standard PETSc practice. </div><div><br class=""></div><div>With PETSc using a matrix-free operation mA and a matrix from which KSP will build the preconditioner A one uses KSPSetOperator(ksp,mA,A); and then just selects the preconditioner with -pc_type xxx For example to use Jacobi preconditioning one uses -pc_type jacobi (note that this only uses the diagonal of A, the rest of A is never used).</div><div><br class=""></div><div>If you wish to precondition mA by fully solving with the matrix A one can use -ksp_monitor_true_residual -pc_type ksp -ksp_ksp_type yyy -ksp_pc_type xxx -ksp_ksp_monitor_true_residual with, for example, yyy of richardson and xxx of jacobi</div><div><br class=""></div><div> Barry</div><div><br class=""></div><div><br class=""></div><div><br class=""></div><div><br class=""><blockquote type="cite" class=""><div class=""><div dir="ltr" class=""><div class="">To verify that I am indeed solving a diagonal system I printed the PETSc matrix from the preconditioner and viewed it in Matlab. It indeed shows it to be a diagonal system. Attached is the plot of the spy command on the printed matrix. The matrix in binary form is also attached. </div><div class=""><br class=""></div><div class="">My understanding is that because the C coefficient is varying in 4 orders of magnitude, i.e., Max(C)/Min(C) ~ 10^4, the matrix is poorly scaled. When I rescale my matrix by 1/C then the system converges in 1 iteration as expected. Is my understanding correct, and that scaling 1/C should be done even for a diagonal system?</div><div class=""><br class=""></div><div class="">When D is non-zero, then scaling by 1/C seems to be very inconvenient as D is stored as side-centered data for the matrix free solver. </div><div class=""><br class=""></div><div class="">In the case that I do not scale my equations by 1/C, is there some solver setting that improves the convergence rate? (With D as non-zero, I have also tried gmres as the ksp solver in the matrix-based preconditioner to get better performance, but it didn't matter much.)</div><div class=""><br class=""></div><div class=""><br class=""></div><div class="">Thanks,</div><div class="">Ramakrishnan Thirumalaisamy</div><div class="">San Diego State University.</div></div>
<span id="cid:f_ku60x9u71"><Temperature_fill.pdf></span><span id="cid:f_ku60x9tq0"><matrix_temperature></span></div></blockquote></div><br class=""></body></html>