1 static char help[] = "Tests MatSetValuesBlockedStencil() in 3d.\n\n"; 2 3 #include <petscmat.h> 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 7 int main(int argc, char **argv) 8 { 9 PetscInt M = 3, N = 4, P = 2, s = 1, w = 2, i, m = PETSC_DECIDE, n = PETSC_DECIDE, p = PETSC_DECIDE; 10 DM da; 11 Mat mat; 12 DMDAStencilType stencil_type = DMDA_STENCIL_BOX; 13 PetscBool flg = PETSC_FALSE; 14 MatStencil idx[2], idy[2]; 15 PetscScalar *values; 16 17 PetscFunctionBeginUser; 18 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 19 PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL)); 20 PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL)); 21 PetscCall(PetscOptionsGetInt(NULL, NULL, "-P", &P, NULL)); 22 PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL)); 23 PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL)); 24 PetscCall(PetscOptionsGetInt(NULL, NULL, "-p", &p, NULL)); 25 PetscCall(PetscOptionsGetInt(NULL, NULL, "-s", &s, NULL)); 26 PetscCall(PetscOptionsGetInt(NULL, NULL, "-w", &w, NULL)); 27 PetscCall(PetscOptionsGetBool(NULL, NULL, "-star", &flg, NULL)); 28 if (flg) stencil_type = DMDA_STENCIL_STAR; 29 30 /* Create distributed array and get vectors */ 31 PetscCall(DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, stencil_type, M, N, P, m, n, p, w, s, 0, 0, 0, &da)); 32 PetscCall(DMSetFromOptions(da)); 33 PetscCall(DMSetUp(da)); 34 PetscCall(DMCreateMatrix(da, &mat)); 35 PetscCall(MatSetFromOptions(mat)); 36 37 idx[0].i = 1; 38 idx[0].j = 1; 39 idx[0].k = 0; 40 idx[1].i = 2; 41 idx[1].j = 1; 42 idx[1].k = 0; 43 idy[0].i = 1; 44 idy[0].j = 2; 45 idy[0].k = 0; 46 idy[1].i = 2; 47 idy[1].j = 2; 48 idy[1].k = 0; 49 PetscCall(PetscMalloc1(2 * 2 * w * w, &values)); 50 for (i = 0; i < 2 * 2 * w * w; i++) values[i] = i; 51 PetscCall(MatSetValuesBlockedStencil(mat, 2, idx, 2, idy, values, INSERT_VALUES)); 52 PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY)); 53 PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY)); 54 55 flg = PETSC_FALSE; 56 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_matdiagonalscalelocal", &flg, NULL)); 57 if (flg) { 58 Vec vec; 59 PetscInt size; 60 PetscMPIInt rank; 61 62 PetscCall(DMGetLocalVector(da, &vec)); 63 PetscCall(VecGetLocalSize(vec, &size)); 64 PetscCall(PetscFree(values)); 65 PetscCall(VecGetArrayWrite(vec, &values)); 66 PetscCallMPI(MPI_Comm_rank(MPI_COMM_WORLD, &rank)); 67 for (i = 0; i < size; i++) values[i] = (PetscScalar)(rank + 1); 68 PetscCall(VecRestoreArrayWrite(vec, &values)); 69 PetscCall(MatDiagonalScaleLocal(mat, vec)); 70 PetscCall(DMRestoreLocalVector(da, &vec)); 71 } 72 73 /* Free memory */ 74 PetscCall(PetscFree(values)); 75 PetscCall(MatDestroy(&mat)); 76 PetscCall(DMDestroy(&da)); 77 PetscCall(PetscFinalize()); 78 return 0; 79 } 80 81 /*TEST 82 83 test: 84 suffix: baij 85 nsize: 2 86 args: -mat_type baij 87 output_file: output/empty.out 88 89 test: 90 suffix: aij 91 nsize: 2 92 args: -mat_type aij 93 output_file: output/empty.out 94 95 test: 96 suffix: baij-diagonalscalelocal 97 nsize: 2 98 args: -mat_type baij -test_matdiagonalscalelocal 99 output_file: output/empty.out 100 101 test: 102 suffix: aij-diagonalscalelocal 103 nsize: 2 104 args: -mat_type aij -test_matdiagonalscalelocal 105 output_file: output/empty.out 106 107 TEST*/ 108