Barry, 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. 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. 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(). Please let me know, if this unclear and I can clarify with more details. Thanks Kaushik On Wed, May 19, 2021 at 11:25 AM Barry Smith <[email protected]> wrote:
*ierr = MatSeqAIJGetArray(ptr->K,&Kc);CHKERRQ(ierr); *
*ierr = MatGetLocalSize(ptr->K,&m,&n); CHKERRQ(ierr);*
*for(i=0;i<m;i++) { for(j=0;j<n;j++) { v = *(Kc + (i*n+j)); ierr = MatSetValues(Kloc,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr); }* * } *
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.
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 https://www.mcs.anl.gov/petsc/documentation/faq.html#functionjacobian nor do you need to pass into FormFunction in the context the F or Jacobian matrix.
Maybe if you explained your use case a bit more we could make suggestions on how to accomplish it.
Barry
On May 19, 2021, at 10:00 AM, Kaushik Vijaykumar <[email protected]> wrote:
Hi everyone,
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 " *formfunction*" to calculate the residual and "*formjacobian*" to calculate the tangent stiffness matrix. To do so, I have defined a struct:
*extern PetscErrorCode FormJacobian(SNES,Vec,Mat,Mat,void*);extern PetscErrorCode FormFunction(SNES,Vec,Vec,void*);*
// define a struct type to pass K and f as "user context"
*typedef struct { Mat K; Vec f; } K_and_f;*
In the main program, the struct is declared *int main()* *{* *K_and_f main_K_and_f;* // declare the struct
// SNES - Populate K and f into the struct *main_K_and_f.K = K;* // K matrix *main_K_and_f.f = f; *// f vector .... .... *}* In form function
*PetscErrorCode FormFunction(SNES snes, Vec x, Vec F, void* ctx) { PetscErrorCode ierr; PetscReal *ax,*c; PetscReal *aF; PetscScalar *Kc,v; PetscInt nlocal,m,n,i,j,index[100000]; * // Create local Vec and Mat * Mat Kloc;* * Vec Floc, Uloc, KUloc, resloc;*
* K_and_f* ptr = (K_and_f*) ctx; *// cast the pointer to void into pointer to struct // Get local F array, FLOC * ierr = VecGhostGetLocalForm(ptr->f,&Floc);CHKERRQ(ierr);*
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.
The following trial gives me incorrect K values in Kloc: *ierr = MatSeqAIJGetArray(ptr->K,&Kc);CHKERRQ(ierr); *
*ierr = MatGetLocalSize(ptr->K,&m,&n); CHKERRQ(ierr);*
*for(i=0;i<m;i++) { for(j=0;j<n;j++) { v = *(Kc + (i*n+j)); ierr = MatSetValues(Kloc,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr); }* * } *
When I compare K and Kloc, they are not identical
Please let me know if there is an equivalent function like *VecGhostGetLocalForm()* for Matrices. If not, is there a better way to do this.
Any guidance/help is greatly appreciated.
Thanks Kaushik