[petsc-users] Jacobian matrix construction when dof > 1

Rahul Samala srahul_05 at yahoo.co.in
Thu Jan 25 08:38:02 CST 2018


Hello PetSc users,
I am a beginner to PetSc and I am trying to write Jacobian matrix for my problem involving 2 nonlinear equations in 2 unknown variables, solving in a 1D domain.Most of the petsc examples involve one variable and I am finding it difficult to understand how the Jacobian matrix is built when there are more than one variables.
For my problem the DOF is 2 (since 2 unknown variables). I am using finite difference and my stencil is 3 (i-1,i,i+1).Jacobian matrix will have 2x2 submatrix elements and with stencil 3, we have to set 4x3=12 values for each grid point.
Let the unknown variables are "u,v" and (F1_u)_i-1 represent "partial derivative of residual of equation 1 at grid point 'i' with respect to variable 'u' at i-1".
Is the following pseudo-code the correct way of filling jacobian matrix?

 FormJacobianLocal(DMDALocalInfo *info,Field *sol,Mat J,Mat Jpre, AppCtx *user){  .......
    MatStencil     col[2][3],row;    PetscReal       v[2][3];
.........for(i=xs;i<xs+xm;i++){    row.i = i; row.c=0;    { // Internal points
        col[0][0].i = i-1;     col[0][0].c = 0;     v[0][0] = jacobian value =  (F1_u)_i-1
        col[0][1].i = i;        col[0][1].c = 0;     v[0][1] = jacobian value = (F1_u)_i
        col[0][2].i = i+1;    col[0][2].c = 0;    v[0][1] = jacobian value =  (F1_u)_i+1
               col[1][0].i = i-1;     col[1][0].c = 1;     v[1][0] = jacobian value =  (F1_v)_i-1        col[1][1].i = i;        col[1][1].c = 1;     v[1][1] = jacobian value = (F1_v)_i
        col[1][2].i = i+1;    col[1][2].c = 1;    v[1][1] = jacobian value =  (F1_v)_i+1
    }  MatSetValuesStencil(Jpre,1,&row,6,col,v,INSERT_VALUES);row.i = i; row.c=1;    { // Internal points
        col[0][0].i = i-1;     col[0][0].c = 0;     v[0][0] = jacobian value =  (F2_u)_i-1
        col[0][1].i = i;        col[0][1].c = 0;     v[0][1] = jacobian value = (F2_u)_i
        col[0][2].i = i+1;    col[0][2].c = 0;    v[0][1] = jacobian value =  (F2_u)_i+1
               col[1][0].i = i-1;     col[1][0].c = 1;     v[1][0] = jacobian value =  (F2_v)_i-1        col[1][1].i = i;        col[1][1].c = 1;     v[1][1] = jacobian value = (F2_v)_i
        col[1][2].i = i+1;    col[1][2].c = 1;    v[1][1] = jacobian value =  (F2_v)_i+1
    }  MatSetValuesStencil(Jpre,1,&row,6,col,v,INSERT_VALUES);
}.......
}
Thank you,Rahul.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20180125/af2633b0/attachment.html>


More information about the petsc-users mailing list