<div dir="ltr"><div dir="ltr">On Wed, May 19, 2021 at 11:48 AM Kaushik Vijaykumar <<a href="mailto:kaushikv318@gmail.com">kaushikv318@gmail.com</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Barry,<div><br></div><div>Thanks for your response. I am trying to copy K into Kloc, the purpose of doing that is to calculate the residual in formfunction - Residual = [Kloc]{X} - {Floc}, where Kloc and Floc are the copies of assembled global K matrix and F vector, and X is the current iterated solution. This is the reason that I need to access the global stiffness matrix and force vector in formfunction.</div><div>I agree that formjacobian only needs access to global stiffness matrix K to form the jacobian and does not need to access the force vector.</div><div><br></div><div>Therefore, to be able to pass both K and F to formfunction, I defined a struct that contains a vector F and matrix K and populated them in the main().</div><div><br></div><div>Please let me know, if this unclear and I can clarify with more details</div></div></blockquote><div><br></div><div>This still does not make sense. If you have a nonlinear problem, how can you know the Jacobian up front? It changes with the input point.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Thanks</div><div>Kaushik</div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, May 19, 2021 at 11:25 AM Barry Smith <<a href="mailto:bsmith@petsc.dev" target="_blank">bsmith@petsc.dev</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div><blockquote type="cite"><div dir="ltr"><div><b>ierr = MatSeqAIJGetArray(ptr->K,&Kc);CHKERRQ(ierr); </b></div><div><b>ierr = MatGetLocalSize(ptr->K,&m,&n); CHKERRQ(ierr);<br></b></div><b>for(i=0;i<m;i++)<br> { <br> for(j=0;j<n;j++)<br> {<br> v = *(Kc + (i*n+j));<br> ierr = MatSetValues(Kloc,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr);<br> }<br></b><div><b> } </b><br><div> </div></div></div></blockquote><div><br></div>I don't understand the purpose of this code fragment. Are you wanting to copy your K into Kloc? What is the purpose of Kloc? With standard usage of SNES one only needs to provide PETSc a global matrix, one does not need to work with "ghosted" matrices at all. <div><br></div><div>Also normally one fills the f vector in FormFunction and the Jacobian in FormJacobian and there is no need to access the Jacobians at all in SNES <a href="https://www.mcs.anl.gov/petsc/documentation/faq.html#functionjacobian" target="_blank">https://www.mcs.anl.gov/petsc/documentation/faq.html#functionjacobian</a> nor do you need to pass into FormFunction in the context the F or Jacobian matrix. </div><div><br></div><div>Maybe if you explained your use case a bit more we could make suggestions on how to accomplish it.</div><div><br></div><div>Barry</div><div><br></div><div><div><br><div><br><blockquote type="cite"><div>On May 19, 2021, at 10:00 AM, Kaushik Vijaykumar <<a href="mailto:kaushikv318@gmail.com" target="_blank">kaushikv318@gmail.com</a>> wrote:</div><br><div><div dir="ltr">Hi everyone,<div><br></div><div>We are in the process of integrating nonlinear solution to our FE code through SNES. The first steps that I understood that need to be done is to be able to pass the assembled stiffness matrix K and force vector F to the "<b>formfunction</b>" to calculate the residual and "<b>formjacobian</b>" to calculate the tangent stiffness matrix. To do so, I have defined a struct:</div><div><br></div><div><b>extern PetscErrorCode FormJacobian(SNES,Vec,Mat,Mat,void*);<br>extern PetscErrorCode FormFunction(SNES,Vec,Vec,void*);</b><br><br>// define a struct type to pass K and f as "user context"<br><b>typedef struct<br> {<br> Mat K;<br> Vec f;<br> } K_and_f;</b><br></div><div><br></div><div>In the main program, the struct is declared</div><div><b>int main()</b></div><div><b>{</b></div><div><b>K_and_f main_K_and_f;</b> // declare the struct</div><div><br></div>// SNES - Populate K and f into the struct<br><b>main_K_and_f.K = K;</b> // K matrix<br><div><b>main_K_and_f.f = f; </b>// f vector <br></div><div>....</div><div>....</div><div><b>}</b></div><div>In form function</div><div><br></div><div><b>PetscErrorCode FormFunction(SNES snes, Vec x, Vec F, void* ctx) {<br> PetscErrorCode ierr;<br> PetscReal *ax,*c;<br> PetscReal *aF;<br> PetscScalar *Kc,v; <br> PetscInt nlocal,m,n,i,j,index[100000];<br><br> </b> // Create local Vec and Mat<br><b> Mat Kloc;</b><br><b> Vec Floc, Uloc, KUloc, resloc;</b><br><b> </b><br><b> K_and_f* ptr = (K_and_f*) ctx; </b>// cast the pointer to void into pointer to struct<br></div><div><b> </b>// Get local F array, FLOC</div><b> ierr = VecGhostGetLocalForm(ptr->f,&Floc);CHKERRQ(ierr);</b><div><br></div><div><br></div><div>I am able to get the f array from the main program using "
VecGhostGetLocalForm" and the vec is correct. However, I am having trouble with reading the K matrix in formfunction. </div><div><br></div><div>The following trial gives me incorrect K values in Kloc:</div><div><b>ierr = MatSeqAIJGetArray(ptr->K,&Kc);CHKERRQ(ierr); </b></div><div><b>ierr = MatGetLocalSize(ptr->K,&m,&n); CHKERRQ(ierr);<br></b></div><b>for(i=0;i<m;i++)<br> { <br> for(j=0;j<n;j++)<br> {<br> v = *(Kc + (i*n+j));<br> ierr = MatSetValues(Kloc,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr);<br> }<br></b><div><b> } </b><br><div> </div><div> When I compare K and Kloc, they are not identical </div><div><br></div><div><br></div><div>Please let me know if there is an equivalent function like <b>VecGhostGetLocalForm()</b> for Matrices. If not, is there a better way to do this.</div><div><br></div><div>Any guidance/help is greatly appreciated.</div><div><br></div><div>Thanks</div><div>Kaushik</div><div><br></div><div><br></div><div><br></div><div><br></div></div></div>
</div></blockquote></div><br></div></div></div></blockquote></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="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>