<div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr">Thank you very much, Matt.<br></div><div dir="ltr"> I followed your advice and learned some MATMFFD related functions. Since I really wanted to use the Matrix-Free method, I learned the few examples of this method. I would like to enhance my understanding of the Matrix-Free method in PETSc according to my questions about these examples and hope to get your answer. </div><div dir="ltr">1. In the example src/snes/example/tutorials/ex22.c</div></div></div></div><blockquote style="margin:0px 0px 0px 40px;border:none;padding:0px"><div><div><div><div><div><font size="1">......</font></div></div></div></div></div><div><div><div><div><div><font size="1">107: packer->ops->creatematrix = DMCreateMatrix_MF;</font></div></div></div></div></div><div><div><div><div><div><br></div></div></div></div></div><div><div><div><div><div><font size="1">110: SNESCreate(PETSC_COMM_WORLD,&snes);</font></div></div></div></div></div><div><div><div><div><div><font size="1">111: SNESSetDM(snes,packer);</font></div></div></div></div></div><div><div><div><div><div><font size="1">112: SNESSetFunction(snes,NULL,ComputeFunction,NULL);</font></div></div></div></div></div><div><div><div><div><div><font size="1">113: SNESSetJacobian(snes,NULL, NULL,ComputeJacobian_MF,NULL);</font></div></div></div></div></div><div><div><div><div><div><font size="1"><br></font></div></div></div></div></div><div><div><div><div><div><font size="1">115: SNESSetFromOptions(snes);</font></div></div></div></div></div><div><div><div><div><div><font size="1">......</font></div></div></div></div></div><div><div><div><div><div><font size="1">285: PetscErrorCode DMCreateMatrix_MF(DM packer,Mat *A)</font></div></div></div></div></div><div><div><div><div><div><font size="1">286: {</font></div></div></div></div></div><div><div><div><div><div><font size="1">288: Vec t;</font></div></div></div></div></div><div><div><div><div><div><font size="1">289: PetscInt m;</font></div></div></div></div></div><div><div><div><div><div><font size="1"><br></font></div></div></div></div></div><div><div><div><div><div><font size="1">292: DMGetGlobalVector(packer,&t);</font></div></div></div></div></div><div><div><div><div><div><font size="1">293: VecGetLocalSize(t,&m);</font></div></div></div></div></div><div><div><div><div><div><font size="1">294: DMRestoreGlobalVector(packer,&t);</font></div></div></div></div></div><div><div><div><div><div><font size="1">295: MatCreateMFFD(PETSC_COMM_WORLD,m,m,PETSC_DETERMINE,PETSC_DETERMINE,A);</font></div></div></div></div></div><div><div><div><div><div><font size="1">296: MatSetUp(*A);</font></div></div></div></div></div><div><div><div><div><div><font size="1">297: return(0);</font></div></div></div></div></div><div><div><div><div><div><font size="1">298: }</font></div></div></div></div></div><div><div><div><div><div><font size="1"><br></font></div></div></div></div></div><div><div><div><div><div><font size="1">300: PetscErrorCode ComputeJacobian_MF(SNES snes,Vec x,Mat A,Mat B,void *ctx)</font></div></div></div></div></div><div><div><div><div><div><font size="1">301: {</font></div></div></div></div></div><div><div><div><div><div><font size="1"><br></font></div></div></div></div></div><div><div><div><div><div><font size="1">305: MatMFFDSetFunction(A,(PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction,snes);</font></div></div></div></div></div><div><div><div><div><div><font size="1">306: MatMFFDSetBase(A,x,NULL);</font></div></div></div></div></div><div><div><div><div><div><font size="1">307: return(0);</font></div></div></div></div></div><div><div><div><div><div><font size="1">308: }</font></div></div></div></div></div></blockquote>Why are the solution vector, Jacobian matrix, and precondition matrix set to NULL in SNESSetJacobian? Is it because they are automatically generated in SNESSetDM? </div><div dir="ltr">There is no precondition matrix B in ComputeJacobian_MF (Jacobian evaluate program). If I want to add my own precondition matrix, then the Assembly matrix B in this function, is that right? </div><div dir="ltr">What is the function of MatMFFDSetBase? I read the documentation but didn't understand it. What would happen if I removed it?</div><div dir="ltr"><br></div><div dir="ltr">2. In the example /src/snes/examples/tutorials/ex14.c <br></div></div></div><blockquote style="margin:0px 0px 0px 40px;border:none;padding:0px"><div><div><div><div><font size="1">114: SNESSetFunction(snes,r,FormFunction,(void*)&user);</font></div></div></div></div><div><div><div><div><font size="1"><br></font></div></div></div></div><div><div><div><div><font size="1">131: PetscOptionsGetBool(NULL,NULL,"-snes_mf",&matrix_free,NULL);</font></div></div></div></div><div><div><div><div><font size="1">132: PetscOptionsGetBool(NULL,NULL,"-fdcoloring",&coloring,NULL);</font></div></div></div></div><div><div><div><div><font size="1">133: PetscOptionsGetBool(NULL,NULL,"-fdcoloring_ds",&coloring_ds,NULL);</font></div></div></div></div><div><div><div><div><font size="1">134: PetscOptionsGetBool(NULL,NULL,"-fdcoloring_local",&local_coloring,NULL);</font></div></div></div></div><div><div><div><div><font size="1">135: if (!matrix_free) {</font></div></div></div></div><div><div><div><div><font size="1">136: DMSetMatType(user.da,MATAIJ);</font></div></div></div></div><div><div><div><div><font size="1">137: DMCreateMatrix(user.da,&J);</font></div></div></div></div><div><div><div><div><font size="1">138: if (coloring) {</font></div></div></div></div><div><div><div><div><font size="1">139: ISColoring iscoloring;</font></div></div></div></div><div><div><div><div><font size="1">140: if (!local_coloring) {</font></div></div></div></div><div><div><div><div><font size="1">141: DMCreateColoring(user.da,IS_COLORING_GLOBAL,&iscoloring);</font></div></div></div></div><div><div><div><div><font size="1">142: MatFDColoringCreate(J,iscoloring,&matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">143: MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))FormFunction,&user);</font></div></div></div></div><div><div><div><div><font size="1">144: } else {</font></div></div></div></div><div><div><div><div><font size="1">145: DMCreateColoring(user.da,IS_COLORING_LOCAL,&iscoloring);</font></div></div></div></div><div><div><div><div><font size="1">146: MatFDColoringCreate(J,iscoloring,&matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">147: MatFDColoringUseDM(J,matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">148: MatFDColoringSetFunction(matfdcoloring,(PetscErrorCode (*)(void))FormFunctionLocal,&user);</font></div></div></div></div><div><div><div><div><font size="1">149: }</font></div></div></div></div><div><div><div><div><font size="1">150: if (coloring_ds) {</font></div></div></div></div><div><div><div><div><font size="1">151: MatFDColoringSetType(matfdcoloring,MATMFFD_DS);</font></div></div></div></div><div><div><div><div><font size="1">152: }</font></div></div></div></div><div><div><div><div><font size="1">153: MatFDColoringSetFromOptions(matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">154: MatFDColoringSetUp(J,iscoloring,matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">155: SNESSetJacobian(snes,J,J,SNESComputeJacobianDefaultColor,matfdcoloring);</font></div></div></div></div><div><div><div><div><font size="1">156: ISColoringDestroy(&iscoloring);</font></div></div></div></div><div><div><div><div><font size="1">157: } else {</font></div></div></div></div><div><div><div><div><font size="1">158: SNESSetJacobian(snes,J,J,FormJacobian,&user);</font></div></div></div></div><div><div><div><div><font size="1">159: }</font></div></div></div></div><div><div><div><div><font size="1">160: }</font></div></div></div></div></blockquote>If the command "- snes_mf" is used, the function SNESSetJacobian is not called, and the Jacobian matrix type is initialized as MATMFFD. The function in SNESSetFunction is used to approximate the matrix vector product. Do I understand this correctly? </div><div dir="ltr">The description of my question may be complicated, but as a beginner, I really want to learn something about matrix-free method from the examples, and I hope to get your reply. <br></div><div dir="ltr"><br></div><div>Thanks,</div><div>Yingjie</div></div><br><div class="gmail_quote"><div dir="ltr">Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> 于2018年10月9日周二 下午11:34写道：<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Mon, Oct 8, 2018 at 11:33 PM Yingjie Wu <<a href="mailto:yjwu16@gmail.com" target="_blank">yjwu16@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr">Dear Petsc developer: </div><div dir="ltr">Hi,</div><div dir="ltr"> </div><div dir="ltr">I've been studying Petsc recently about Precontioner and Metrix-Free, and I have some questions that puzzle me. </div><div dir="ltr"> </div><div dir="ltr">1. I want to test block Jacobi preconditioner, so I chose /snes/example/tutorial/ex3.c as an example. According to the reference in the example, the input parameters are: </div><div dir="ltr">mpiexec -n 8./ex3 -nox -n 10000 -ksp_type fgmres -pc_type bjacobi -pc_bjacobi_blocks 4 -sub_ksp_type gmres -sub_ksp_max_it 3 -post_setsubksp -sub_ksp_rtol 1.e-16 </div></div></div></div></div></div></div></blockquote><div><br></div><div>You do not care about recursive blocks, so just use</div><div><br></div><div> <span style="font-family:Menlo;font-size:13px">$MPIEXEC -n 8 ./ex3 -nox -n 10000 -ksp_type fgmres -pc_type bjacobi -pc_bjacobi_blocks 4 -sub_ksp_type gmres -sub_ksp_max_it 3 -snes_view -ksp_view</span></div><div><span style="font-family:Menlo;font-size:13px"><br></span></div><div>and I get the attached output.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr">I want to export each block of KSP and PC information : </div><div dir="ltr"> -snes_view -ksp_view </div><div dir="ltr">However, the procedure is wrong and the wrong information is as follows: </div><blockquote style="margin:0px 0px 0px 40px;border:none;padding:0px"><div> <font size="1"> <span style="font-family:等线;text-align:justify">[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></font></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: Different communicators in the two objects: Argument # 1 and 2 flag 3</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: [1]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[2]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[2]PETSC ERROR: [3]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[3]PETSC ERROR: Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[3]PETSC ERROR: [6]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[6]PETSC ERROR: Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[6]PETSC ERROR: Different communicators in the two objects: Argument # 1 and 2 flag 3</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[6]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[7]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[7]PETSC ERROR: Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[7]PETSC ERROR: Different communicators in the two objects: Argument # 1 and 2 flag 3</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[7]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[7]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: ./ex3 on a arch-linux2-c-debug named yjwu-XPS-8910 by yjwu Mon Oct 8 22:35:34 2018</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-mpich --download-fblaslapack</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: #1 KSPView() line 118 in /home/yjwu/petsc-3.10.1/src/ksp/ksp/interface/itcreate.c</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: #2 PCView_BJacobi() line 232 in /home/yjwu/petsc-3.10.1/src/ksp/pc/impls/bjacobi/bjacobi.c</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[0]PETSC ERROR: #3 PCView() line 1651 in /home/yjwu/petsc-3.10.1/src/ksp/pc/interface/precon.c</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: Different communicators in the two objects: Argument # 1 and 2 flag 3</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: ./ex3 on a arch-linux2-c-debug named yjwu-XPS-8910 by yjwu Mon Oct 8 22:35:34 2018</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-mpich --download-fblaslapack</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[1]PETSC ERROR: #1 KSPView() line 118 in /home/yjwu/petsc-3.10.1/src/ksp/ksp/interface/itcreate.c</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">Arguments must have same communicators</span></div><div><span style="font-family:等线;font-size:x-small;text-align:justify">[2]PETSC ERROR: Different communicators in the two objects: Argument # 1 and 2 flag 3</span></div><div><font size="1"><span style="font-family:等线;text-align:justify">[2]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</span> </font></div></blockquote><p class="m_-5073886653601292629m_3708552985005253687gmail-15" align="justify" style="text-align:justify"><font size="1"><span style="font-family:等线"></span></font></p><p class="m_-5073886653601292629m_3708552985005253687gmail-15" align="justify" style="text-align:justify"><font size="1"><span style="font-family:等线"></span></font></p><div dir="ltr">If I want to get the information of KSP and PC in each block, what should I do? </div><div dir="ltr"><br></div><div dir="ltr">2. There is a requirement in my program to use Matrix-Free method to approximate Jacobian matrix by finite difference of residual function in solving nonlinear equations. But I'll also provide an analytic( incomplete, some terms are missing or approximate) for preconditioning. Because my problem is about a large set of equations composed of several physical fields, I want to use block Jacobian precondition for each subfield(block), and ILU sub_pc for each subfield. After reading the Users'Guide, I found that using - snex_mf_operator can do the above, so I added:</div><div dir="ltr"> - snes_mf_operator </div></div></div></div></div></div></div></blockquote><div><br></div><div>You need to be careful what matrix you are adding values to in your FormJacobian() routine. The primal matrix J is of type MFFD (finite difference)</div><div>and thus cannot accept values. You put your approximate values in the preconditioner M.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr">after the example above. However, the procedure is wrong and does not seem to support this. </div><div dir="ltr"> The wrong information is: </div></div><blockquote style="margin:0px 0px 0px 40px;border:none;padding:0px"><div><div> <font size="1">[1]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</font></div></div><div><div><font size="1">[2]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</font></div></div><div><div><font size="1">[2]PETSC ERROR: No support for this operation for this object type</font></div></div><div><div><font size="1">[2]PETSC ERROR: Mat type mffd</font></div></div><div><div><font size="1">[2]PETSC ERROR: [3]PETSC ERROR: [4]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</font></div></div><div><div><font size="1">[4]PETSC ERROR: No support for this operation for this object type</font></div></div><div><div><font size="1">[4]PETSC ERROR: Mat type mffd</font></div></div><div><div><font size="1">[4]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</font></div></div><div><div><font size="1">[4]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018 </font></div></div><div><div><font size="1">[6]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</font></div></div><div><div><font size="1">[6]PETSC ERROR: No support for this operation for this object type</font></div></div><div><div><font size="1">[6]PETSC ERROR: Mat type mffd</font></div></div><div><div><font size="1">[6]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</font></div></div><div><div><font size="1">[6]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018 </font></div></div><div><div><font size="1">[6]PETSC ERROR: ./ex3 on a arch-linux2-c-debug named yjwu-XPS-8910 by yjwu Mon Oct 8 23:01:57 2018</font></div></div><div><div><font size="1">[6]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-mpich --download-fblaslapack</font></div></div><div><div><font size="1">[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</font></div></div><div><div><font size="1">[0]PETSC ERROR: No support for this operation for this object type</font></div></div><div><div><font size="1">[0]PETSC ERROR: Mat type mffd</font></div></div><div><div><font size="1">[0]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</font></div></div><div><div><font size="1">[0]PETSC ERROR: Petsc Release Version 3.10.1, Sep, 26, 2018 </font></div></div><div><div><font size="1">[0]PETSC ERROR: ./ex3 on a arch-linux2-c-debug named yjwu-XPS-8910 by yjwu Mon Oct 8 23:01:57 2018</font></div></div><div><div><font size="1">[0]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-mpich --download-fblaslapack</font></div></div><div><div><font size="1">[0]PETSC ERROR: #1 MatSetValues() line 1330 in /home/yjwu/petsc-3.10.1/src/mat/interface/matrix.c</font></div></div><div><div><font size="1">[1]PETSC ERROR: No support for this operation for this object type</font></div></div><div><div><font size="1">[1]PETSC ERROR: Mat type mffd</font></div></div><div><div><font size="1">[1]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting. </font></div></div></blockquote><div dir="ltr"><div dir="ltr">If I want to implement matrix-free and block jacobi precondition in my program, what do I need to pay attention or add ?<br></div><div dir="ltr"><br></div><div dir="ltr">The second question is very important for my research and I hope to get your answer. <br></div><div dir="ltr"><br></div><div>Thanks,</div><div>Yingjie</div></div></div></div></div></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="m_-5073886653601292629gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>
</blockquote></div>