[petsc-users] Help needed on DMDA dm_boundary_none and matsetvaluesstencil

[Karl Yang]

Hi, all

I've encountered some issues with DM_BOUDNARY_NONE and MatSetValuesStencil. I had a code with DM_BOUNDARY_PERIODIC which was working fine. And I simply change the boundary, and find it not working any more. I wonder is there any difference in terms of indexing and stencil for DM_BOUNDARY_NONE.
The following is a simplified code to demonstrate what I was doing. It is basically assembling for finite elements.
But the MatSetValuesStencil seems not adding values into the matrix as I expected and some entries disappeared after MatSetValuesStencil for the second time. I've attached the output at the two different matview location. Some entries in the matrix disappeared after the second values add. And the order of matrix is wired to me, the matrix output is not in a ascending order. Appreciate if anyone would help.

/////////////////////////demo code ///////////////////////////////////
DM dm;
Mat A;
MatStencil s_u[4];

DMDACreate2d(PETSC_COMM_SELF, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 5, 5, PETSC_DECIDE, PETSC_DECIDE, 3, 1, NULL, NULL, &dm);
DMSetMatType(dm, MATAIJ);
DMSetFromOptions(dm);
DMSetUp(dm);
DMDASetUniformCoordinates(dm, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0);

DMSetMatrixPreallocateOnly(dm, PETSC_TRUE);
DMCreateMatrix(dm, &A);

s_u[0].i = 0; s_u[0].j = 0; s_u[0].c = 1;
s_u[1].i = 0; s_u[1].j = 0+1; s_u[1].c = 1;
s_u[2].i = 0+1; s_u[2].j = 0+1; s_u[2].c = 1;
s_u[3].i = 0+1; s_u[3].j = 0; s_u[3].c = 1;

double Ke[16];
for (int n=0;n<16;++n){Ke[n]=1;};
MatSetValuesStencil(A,4,s_u,4,s_u,Ke,ADD_VALUES);

// MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
// MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
// MatView(A, PETSC_VIEWER_STDOUT_WORLD); //first matview

s_u[0].i = 1; s_u[0].j = 0; s_u[0].c = 1;
s_u[1].i = 1; s_u[1].j = 0+1; s_u[1].c = 1;
s_u[2].i = 1+1; s_u[2].j = 0+1; s_u[2].c = 1;
s_u[3].i = 1+1; s_u[3].j = 0; s_u[3].c = 1;

MatSetValuesStencil(A,4,s_u,4,s_u,Ke,ADD_VALUES);
MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
MatView(A, PETSC_VIEWER_STDOUT_WORLD); //second matview

////////////////first matview ////////////////////////
row 0:
row 1: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 2:
row 3:
row 4: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 5:
row 6:
row 7:
row 8:
row 9:
row 10:
row 11:
row 12:
row 13:
row 14:
row 15:
row 16: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 17:
row 18:
row 19: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 20:
row 21:
row 22:
row 23:
row 24:
row 25:
row 26:
row 27:
row 28:
row 29:
row 30:
row 31:
row 32:
row 33:
row 34:
row 35:
row 36:
row 37:
row 38:
row 39:
row 40:
row 41:
row 42:
row 43:
row 44:
row 45:
row 46:
row 47:
row 48:
row 49:
row 50:
row 51:
row 52:
row 53:
row 54:
row 55:
row 56:
row 57:
row 58:
row 59:
row 60:
row 61:
row 62:
row 63:
row 64:
row 65:
row 66:
row 67:
row 68:
row 69:
row 70:
row 71:
row 72:
row 73:
row 74:

///////////////second matview/////////////////////
row 0:
row 1: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 2:
row 3:
row 4: (4, 1.) (19, 1.) (22, 1.) (7, 1.)
row 5:
row 6:
row 7: (4, 1.) (19, 1.) (22, 1.) (7, 1.)
row 8:
row 9:
row 10:
row 11:
row 12:
row 13:
row 14:
row 15:
row 16: (1, 1.) (16, 1.) (19, 1.) (4, 1.)
row 17:
row 18:
row 19: (4, 1.) (19, 1.) (22, 1.) (7, 1.)
row 20:
row 21:
row 22: (4, 1.) (19, 1.) (22, 1.) (7, 1.)
row 23:
row 24:
row 25:
row 26:
row 27:
row 28:
row 29:
row 30:
row 31:
row 32:
row 33:
row 34:
row 35:
row 36:
row 37:
row 38:
row 39:
row 40:
row 41:
row 42:
row 43:
row 44:
row 45:
row 46:
row 47:
row 48:
row 49:
row 50:
row 51:
row 52:
row 53:
row 54:
row 55:
row 56:
row 57:
row 58:
row 59:
row 60:
row 61:
row 62:
row 63:
row 64:
row 65:
row 66:
row 67:
row 68:
row 69:
row 70:
row 71:
row 72:
row 73:
row 74:

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210209/13345f63/attachment.html>