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

[Kaushik Vijaykumar]

Matt,

Thanks for the comment, I see your point now. I was using the linear
elastic [K] from the previous load step as the guess for the jacobian. The
stiffness matrix was being built outside formfunction using displacements
from the previous load step. What is the best approach, Build [K] abd {F}
inside the form function and compute the residual [K]_{i}{U}_{i} - {F}_{i},
where i designates the iterated matrices and vector for a given load step ?


Thanks
Kaushik

On Wed, May 19, 2021 at 11:51 AM Matthew Knepley <knepley at gmail.com> wrote:

> On Wed, May 19, 2021 at 11:48 AM Kaushik Vijaykumar <kaushikv318 at gmail.com>
> wrote:
>
>> 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
>>
>
> 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.
>
>   Thanks,
>
>     Matt
>
>
>> 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
>>>
>>>
>>>
>>>
>>>
>>>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which their
> experiments lead.
> -- Norbert Wiener
>
> https://www.cse.buffalo.edu/~knepley/
> <http://www.cse.buffalo.edu/~knepley/>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210519/48ed3468/attachment.html>