[petsc-users] Using SNES to solve Non linear elasticity, passing matrix to form function

[Kaushik Vijaykumar]

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 <bsmith at petsc.dev> 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 <kaushikv318 at gmail.com>
> 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
>
>
>
>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210519/3da1d8af/attachment-0001.html>