<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> Feng,<div class=""><br class=""></div><div class=""> The first thing to check is that for each linear solve that involves a new operator (values in the base vector) the MFFD matrix knows it is using a new operator. </div><div class=""><br class=""></div><div class="">The easiest way is to call <font face="Calibri, Helvetica, sans-serif" class="">MatMFFDSetBase() before each solve that involves a new operator (new values in the base vector). Also be careful about </font><span style="font-family: Calibri, Helvetica, sans-serif;" class=""> </span><span style="font-family: Calibri, Helvetica, sans-serif;" class="">petsc_baserhs, when you change the base vector's values you also need to change the </span><span style="font-family: Calibri, Helvetica, sans-serif;" class="">petsc_baserhs values to the function evaluation at that point.</span></div><div class=""><span style="font-family: Calibri, Helvetica, sans-serif;" class=""><br class=""></span></div><div class=""><span style="font-family: Calibri, Helvetica, sans-serif;" class="">If that is correct I would check with a trivial function evaluator to make sure the infrastructure is all set up correctly. For examples use for the matrix free a 1 4 1 operator applied matrix free. </span></div><div class=""><br class=""></div><div class="">Barry</div><div class=""><br class=""></div><div class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Mar 11, 2021, at 7:35 AM, feng wang <<a href="mailto:snailsoar@hotmail.com" class="">snailsoar@hotmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">Dear All,</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">I am new to petsc and trying to implement a matrix-free GMRES. I have assembled an approximate Jacobian matrix just for preconditioning. After reading some previous questions on this topic, my approach is:</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">the matrix-free matrix is created as:</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""> ierr = MatCreateMFFD(*A_COMM_WORLD, iqe*blocksize, iqe*blocksize, PETSC_DETERMINE, PETSC_DETERMINE, &petsc_A_mf); CHKERRQ(ierr);<div class=""> ierr = MatMFFDSetFunction(petsc_A_mf, FormFunction_mf, this); CHKERRQ(ierr);</div><div class=""></div><div class=""><br class=""></div><div class="">KSP linear operator is set up as:</div><div class=""><br class=""></div><div class=""> ierr = KSPSetOperators(petsc_ksp, petsc_A_mf, petsc_A_pre); CHKERRQ(ierr); //petsc_A_pre is my assembled pre-conditioning matrix<br class=""></div><div class=""><br class=""></div><div class="">Before calling KSPSolve, I do:</div><div class=""><br class=""></div> ierr = MatMFFDSetBase(petsc_A_mf, petsc_csv, petsc_baserhs); CHKERRQ(ierr); //petsc_csv is the flow states, petsc_baserhs is the pre-computed right hand side</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">The call back function is defined as:</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""> PetscErrorCode cFdDomain::FormFunction_mf(void *ctx, Vec in_vec, Vec out_vec)<div class=""> {</div><div class=""> PetscErrorCode ierr;</div><div class=""> cFdDomain *user_ctx;</div><div class=""><br class=""></div><div class=""> cout << "FormFunction_mf called\n";</div><div class=""><br class=""></div><div class=""> //in_vec: flow states</div><div class=""> //out_vec: right hand side + diagonal contributions from CFL number</div><div class=""><br class=""></div><div class=""> user_ctx = (cFdDomain*)ctx;</div><div class=""><br class=""></div><div class=""> //get perturbed conservative variables from petsc</div><div class=""> user_ctx->petsc_getcsv(in_vec);</div><div class=""><br class=""></div><div class=""> //get new right side</div><div class=""> user_ctx->petsc_fd_rhs();</div><div class=""><br class=""></div><div class=""> //set new right hand side to the output vector</div><div class=""> user_ctx->petsc_setrhs(out_vec);</div><div class=""><br class=""></div><div class=""> ierr = 0;</div><div class=""> return ierr;</div><div class=""> }</div></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">The linear system I am solving is (J+D)x=RHS. J is the Jacobian matrix. D is a diagonal matrix and it is used to stabilise the solution at the start but reduced gradually when the solution moves on to recover Newton's method. I add D*x to the true right side when non-linear function is computed to work out finite difference Jacobian, so when finite difference is used, it actually computes (J+D)*dx.<span class="Apple-converted-space"> </span><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">The code runs but diverges in the end. If I don't do matrix-free and use my approximate Jacobian matrix, GMRES works. So something is wrong with my matrix-free implementation. Have I missed something in my implementation? Besides, is there a way to check if the finite difference Jacobian matrix is computed correctly in a matrix-free implementation?</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class=""><br class=""></div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">Thanks for your help in advance.</div><div style="font-style: normal; font-variant-caps: normal; font-weight: normal; letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; text-decoration: none; font-family: Calibri, Helvetica, sans-serif; font-size: 12pt;" class="">Feng <span class="Apple-converted-space"> </span></div></div></blockquote></div><br class=""></div></body></html>