<html><head><meta http-equiv="content-type" content="text/html; charset=utf-8"></head><body style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;">Dear all,<br><br>I tried to use MatZeroRowsColumns to eliminate Dirichlet boundary nodes. However, it cannot eliminate correctly in parallel.<br>Please see the attached code which uses DMDA to create the matrix.<br>When I used one process, it works as expected.<br>For two processes, the domain is split in the x direction. But the 10th row, 20th column is not eliminated as observed when using one process. The results for two processes are also attached.<br>I have input the same rows to be eliminated for both processes.<br>Thank you for any help.<br><br>#include <stdio.h><br>#include <petscdmda.h><br><br>int main(int argc, char **argv)<br>{<br> PetscInt M = 5, N = 3, m = PETSC_DECIDE, n = PETSC_DECIDE, ncomp = 2;<br> PetscInt i, j;<br> DMDALocalInfo daInfo;<br> DM da;<br> Mat A;<br> Vec x, b;<br> MatStencil row, col[5];<br> PetscScalar v[5];<br> PetscInt n_dirichlet_rows = 0, dirichlet_rows[2*(M+N)];<br><br> PetscFunctionBeginUser;<br> PetscCall(PetscInitialize(&argc, &argv, (char *)0, NULL));<br> PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_GHOSTED, DM_BOUNDARY_GHOSTED, DMDA_STENCIL_BOX, M, N, m, n, ncomp, 1, NULL, NULL, &da));<br> PetscCall(DMSetFromOptions(da));<br> PetscCall(DMSetUp(da));<br><br> PetscCall(DMView(da, PETSC_VIEWER_STDOUT_WORLD));<br><br> PetscCall(DMDAGetLocalInfo(da, &daInfo));<br> PetscCall(DMSetMatrixPreallocateOnly(da, PETSC_TRUE));<br><br> PetscCall(DMCreateMatrix(da, &A));<br> PetscCall(MatCreateVecs(A, &x, &b));<br> PetscCall(MatZeroEntries(A));<br> PetscCall(VecZeroEntries(x));<br> PetscCall(VecZeroEntries(b));<br><br> for (j = daInfo.ys; j < daInfo.ys + daInfo.ym; ++j) {<br> for (i = daInfo.xs; i < daInfo.xs + daInfo.xm; ++i) {<br> row.j = j;<br> row.i = i;<br> row.c = 0;<br><br> col[0].j = j;<br> col[0].i = i;<br> col[0].c = 0;<br> v[0] = row.i + col[0].i + 2;<br> col[1].j = j;<br> col[1].i = i - 1;<br> col[1].c = 0;<br> v[1] = row.i + col[1].i + 2;<br> col[2].j = j;<br> col[2].i = i + 1;<br> col[2].c = 0;<br> v[2] = row.i + col[2].i + 2;<br> col[3].j = j - 1;<br> col[3].i = i;<br> col[3].c = 0;<br> v[3] = row.i + col[2].i + 2;<br> col[4].j = j + 1;<br> col[4].i = i;<br> col[4].c = 0;<br> v[4] = row.i + col[2].i + 2;<br><br> PetscCall(MatSetValuesStencil(A, 1, &row, 5, col, v, ADD_VALUES));<br><br> row.j = j;<br> row.i = i;<br> row.c = 1;<br> col[0].j = j;<br> col[0].i = i;<br> col[0].c = 1;<br> v[1] = row.j + col[1].j + 2;<br> col[1].j = j - 1;<br> col[1].i = i;<br> col[1].c = 1;<br> v[1] = row.j + col[1].j + 2;<br> col[2].j = j + 1;<br> col[2].i = i;<br> col[2].c = 1;<br> v[2] = row.j + col[2].j + 2;<br> col[3].j = j;<br> col[3].i = i - 1;<br> col[3].c = 1;<br> v[3] = row.j + col[1].j + 2;<br> col[4].j = j;<br> col[4].i = i + 1;<br> col[4].c = 1;<br> v[4] = row.j + col[2].j + 2;<br><br> PetscCall(MatSetValuesStencil(A, 1, &row, 5, col, v, ADD_VALUES));<br><br> }<br> }<br> PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));<br> PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));<br> MatView(A, 0);<br><br> for (j = 0; j < daInfo.my; ++j) {<br> dirichlet_rows[n_dirichlet_rows++] = j * <a href="https://urldefense.us/v3/__http://dainfo.mx/__;!!G_uCfscf7eWS!faAU5SYXVdbeL4y_c26LP7pbZnC1p3U_6-f4qNMpQMw2dXi4plc1QMmEh9_6S1xCYurJg4RIuXOrnKz85IkV6m3YcLROpRRf$">daInfo.mx</a> * ncomp;<br> dirichlet_rows[n_dirichlet_rows++] = (j+1) * <a href="https://urldefense.us/v3/__http://dainfo.mx/__;!!G_uCfscf7eWS!faAU5SYXVdbeL4y_c26LP7pbZnC1p3U_6-f4qNMpQMw2dXi4plc1QMmEh9_6S1xCYurJg4RIuXOrnKz85IkV6m3YcLROpRRf$">daInfo.mx</a> * ncomp - ncomp;<br> }<br> PetscCall(PetscPrintf(PETSC_COMM_SELF, "n_dirichlet_rows: %d\n", n_dirichlet_rows));<br> for (j = 0; j < n_dirichlet_rows; ++j) {<br> PetscCall(PetscPrintf(PETSC_COMM_SELF, "%d, ", dirichlet_rows[j]));<br> }<br> PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));<br> PetscCall(MatZeroRowsColumns(A, n_dirichlet_rows, dirichlet_rows, 1, NULL, NULL));<br> MatView(A, 0);<br><br> PetscCall(VecDestroy(&x));<br> PetscCall(VecDestroy(&b));<br> PetscCall(DMDestroy(&da));<br> PetscCall(PetscFinalize());<br> return 0;<br>}<br><br>———————————————————<br><br>n_dirichlet_rows: 6<br>0, 8, 10, 18, 20, 28,<br>Mat Object: 2 MPI processes<br> type: mpiaij<br> row 0: (0, 1.) (2, 0.) (10, 0.)<br> row 1: (1, 2.) (3, 3.) (11, 3.)<br> row 2: (0, 0.) (2, 4.) (4, 5.) (12, 0.)<br> row 3: (1, 1.) (3, 4.) (5, 3.) (13, 3.)<br> row 4: (2, 5.) (4, 6.) (6, 0.) (14, 0.)<br> row 5: (3, 1.) (5, 6.) (7, 3.) (15, 3.)<br> row 6: (4, 0.) (6, 1.) (8, 0.) (16, 0.)<br> row 7: (5, 1.) (7, 8.) (9, 3.) (17, 3.)<br> row 8: (6, 0.) (8, 1.) (18, 0.)<br> row 9: (7, 1.) (9, 10.) (19, 3.)<br> row 10: (0, 0.) (10, 2.) (12, 0.) (20, 3.)<br> row 11: (1, 3.) (11, 2.) (13, 5.) (21, 5.)<br> row 12: (2, 0.) (10, 0.) (12, 1.) (14, 0.) (22, 0.)<br> row 13: (3, 3.) (11, 3.) (13, 4.) (15, 5.) (23, 5.)<br> row 14: (4, 0.) (12, 0.) (14, 1.) (16, 0.) (24, 0.)<br> row 15: (5, 3.) (13, 3.) (15, 6.) (17, 5.) (25, 5.)<br> row 16: (6, 0.) (14, 0.) (16, 8.) (18, 9.) (26, 9.)<br> row 17: (7, 3.) (15, 3.) (17, 8.) (19, 5.) (27, 5.)<br> row 18: (8, 0.) (16, 9.) (18, 10.) (28, 0.)<br> row 19: (9, 3.) (17, 3.) (19, 10.) (29, 5.)<br> row 20: (10, 3.) (20, 2.) (22, 3.)<br> row 21: (11, 5.) (21, 2.) (23, 7.)<br> row 22: (12, 0.) (20, 3.) (22, 4.) (24, 5.)<br> row 23: (13, 5.) (21, 5.) (23, 4.) (25, 7.)<br> row 24: (14, 0.) (22, 5.) (24, 6.) (26, 7.)<br> row 25: (15, 5.) (23, 5.) (25, 6.) (27, 7.)<br> row 26: (16, 9.) (24, 7.) (26, 8.) (28, 0.)<br> row 27: (17, 5.) (25, 5.) (27, 8.) (29, 7.)<br> row 28: (18, 0.) (26, 0.) (28, 1.)<br> row 29: (19, 5.) (27, 5.) (29, 10.)<br><br>Junming</body></html>