Actual source code: ex132.c


  2: static char help[] = "Test MatAXPY()\n\n";

  4: #include <petscmat.h>

  6: int main(int argc, char **args)
  7: {
  8:   Mat         C, C1, C2, CU;
  9:   PetscScalar v;
 10:   PetscInt    Ii, J, Istart, Iend;
 11:   PetscInt    i, j, m = 3, n;
 12:   PetscMPIInt size;
 13:   PetscBool   mat_nonsymmetric = PETSC_FALSE, flg;
 14:   MatInfo     info;

 16:   PetscFunctionBeginUser;
 17:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 18:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
 19:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 20:   n = 2 * size;

 22:   /* Set flag if we are doing a nonsymmetric problem; the default is symmetric. */
 23:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-mat_nonsym", &mat_nonsymmetric, NULL));

 25:   PetscCall(MatCreate(PETSC_COMM_WORLD, &C));
 26:   PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m * n, m * n));
 27:   PetscCall(MatSetFromOptions(C));
 28:   PetscCall(MatSeqAIJSetPreallocation(C, 5, NULL));
 29:   PetscCall(MatMPIAIJSetPreallocation(C, 5, NULL, 5, NULL));

 31:   PetscCall(MatGetOwnershipRange(C, &Istart, &Iend));
 32:   for (Ii = Istart; Ii < Iend; Ii++) {
 33:     v = -1.0;
 34:     i = Ii / n;
 35:     j = Ii - i * n;
 36:     if (i > 0) {
 37:       J = Ii - n;
 38:       PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, ADD_VALUES));
 39:     }
 40:     if (i < m - 1) {
 41:       J = Ii + n;
 42:       PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, ADD_VALUES));
 43:     }
 44:     if (j > 0) {
 45:       J = Ii - 1;
 46:       PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, ADD_VALUES));
 47:     }
 48:     if (j < n - 1) {
 49:       J = Ii + 1;
 50:       PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, ADD_VALUES));
 51:     }
 52:     v = 4.0;
 53:     PetscCall(MatSetValues(C, 1, &Ii, 1, &Ii, &v, ADD_VALUES));
 54:   }

 56:   /* Make the matrix nonsymmetric if desired */
 57:   if (mat_nonsymmetric) {
 58:     for (Ii = Istart; Ii < Iend; Ii++) {
 59:       v = -1.5;
 60:       i = Ii / n;
 61:       if (i > 1) {
 62:         J = Ii - n - 1;
 63:         PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, ADD_VALUES));
 64:       }
 65:     }
 66:   } else {
 67:     PetscCall(MatSetOption(C, MAT_SYMMETRIC, PETSC_TRUE));
 68:     PetscCall(MatSetOption(C, MAT_SYMMETRY_ETERNAL, PETSC_TRUE));
 69:   }
 70:   PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
 71:   PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));
 72:   PetscCall(PetscObjectSetName((PetscObject)C, "C"));
 73:   PetscCall(MatViewFromOptions(C, NULL, "-view"));

 75:   /* C1 = 2.0*C1 + C, C1 is anti-diagonal and has different non-zeros than C */
 76:   PetscCall(MatCreate(PETSC_COMM_WORLD, &C1));
 77:   PetscCall(MatSetSizes(C1, PETSC_DECIDE, PETSC_DECIDE, m * n, m * n));
 78:   PetscCall(MatSetFromOptions(C1));
 79:   PetscCall(MatSeqAIJSetPreallocation(C1, 1, NULL));
 80:   PetscCall(MatMPIAIJSetPreallocation(C1, 1, NULL, 1, NULL));
 81:   for (Ii = Istart; Ii < Iend; Ii++) {
 82:     v = 1.0;
 83:     i = m * n - Ii - 1;
 84:     j = Ii;
 85:     PetscCall(MatSetValues(C1, 1, &i, 1, &j, &v, ADD_VALUES));
 86:   }
 87:   PetscCall(MatAssemblyBegin(C1, MAT_FINAL_ASSEMBLY));
 88:   PetscCall(MatAssemblyEnd(C1, MAT_FINAL_ASSEMBLY));
 89:   PetscCall(PetscObjectSetName((PetscObject)C1, "C1"));
 90:   PetscCall(MatViewFromOptions(C1, NULL, "-view"));
 91:   PetscCall(MatDuplicate(C1, MAT_COPY_VALUES, &CU));

 93:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " MatAXPY(C1,2.0,C,DIFFERENT_NONZERO_PATTERN)...\n"));
 94:   PetscCall(MatAXPY(C1, 2.0, C, DIFFERENT_NONZERO_PATTERN));
 95:   PetscCall(MatAXPY(CU, 2.0, C, UNKNOWN_NONZERO_PATTERN));
 96:   PetscCall(MatGetInfo(C1, MAT_GLOBAL_SUM, &info));
 97:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " C1: nz_allocated = %g; nz_used = %g; nz_unneeded = %g\n", info.nz_allocated, info.nz_used, info.nz_unneeded));
 98:   PetscCall(MatViewFromOptions(C1, NULL, "-view"));
 99:   PetscCall(MatMultEqual(CU, C1, 10, &flg));
100:   if (!flg) {
101:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error UNKNOWN_NONZERO_PATTERN (supposedly DIFFERENT_NONZERO_PATTERN)\n"));
102:     PetscCall(MatViewFromOptions(CU, NULL, "-view"));
103:   }
104:   PetscCall(MatDestroy(&CU));

106:   /* Secondly, compute C1 = 2.0*C2 + C1, C2 has non-zero pattern of C */
107:   PetscCall(MatDuplicate(C, MAT_DO_NOT_COPY_VALUES, &C2));
108:   PetscCall(MatDuplicate(C1, MAT_COPY_VALUES, &CU));

110:   for (Ii = Istart; Ii < Iend; Ii++) {
111:     v = 1.0;
112:     PetscCall(MatSetValues(C2, 1, &Ii, 1, &Ii, &v, ADD_VALUES));
113:   }
114:   PetscCall(MatAssemblyBegin(C2, MAT_FINAL_ASSEMBLY));
115:   PetscCall(MatAssemblyEnd(C2, MAT_FINAL_ASSEMBLY));
116:   PetscCall(PetscObjectSetName((PetscObject)C2, "C2"));
117:   PetscCall(MatViewFromOptions(C2, NULL, "-view"));
118:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " MatAXPY(C1,2.0,C2,SUBSET_NONZERO_PATTERN)...\n"));
119:   PetscCall(MatAXPY(C1, 2.0, C2, SUBSET_NONZERO_PATTERN));
120:   PetscCall(MatAXPY(CU, 2.0, C2, UNKNOWN_NONZERO_PATTERN));
121:   PetscCall(MatGetInfo(C1, MAT_GLOBAL_SUM, &info));
122:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " C1: nz_allocated = %g; nz_used = %g; nz_unneeded = %g\n", info.nz_allocated, info.nz_used, info.nz_unneeded));
123:   PetscCall(MatViewFromOptions(C1, NULL, "-view"));
124:   PetscCall(MatMultEqual(CU, C1, 10, &flg));
125:   if (!flg) {
126:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error UNKNOWN_NONZERO_PATTERN (supposedly SUBSET_NONZERO_PATTERN)\n"));
127:     PetscCall(MatViewFromOptions(CU, NULL, "-view"));
128:   }
129:   PetscCall(MatDestroy(&CU));

131:   /* Test SAME_NONZERO_PATTERN computing C2 = C2 + 2.0 * C */
132:   PetscCall(MatDuplicate(C2, MAT_COPY_VALUES, &CU));
133:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " MatAXPY(C2,2.0,C,SAME_NONZERO_PATTERN)...\n"));
134:   PetscCall(MatAXPY(C2, 2.0, C, SAME_NONZERO_PATTERN));
135:   PetscCall(MatAXPY(CU, 2.0, C, UNKNOWN_NONZERO_PATTERN));
136:   PetscCall(MatGetInfo(C2, MAT_GLOBAL_SUM, &info));
137:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, " C2: nz_allocated = %g; nz_used = %g; nz_unneeded = %g\n", info.nz_allocated, info.nz_used, info.nz_unneeded));
138:   PetscCall(MatViewFromOptions(C2, NULL, "-view"));
139:   PetscCall(MatMultEqual(CU, C2, 10, &flg));
140:   if (!flg) {
141:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error UNKNOWN_NONZERO_PATTERN (supposedly SUBSET_NONZERO_PATTERN)\n"));
142:     PetscCall(MatViewFromOptions(CU, NULL, "-view"));
143:   }
144:   PetscCall(MatDestroy(&CU));

146:   PetscCall(MatDestroy(&C1));
147:   PetscCall(MatDestroy(&C2));
148:   PetscCall(MatDestroy(&C));

150:   PetscCall(PetscFinalize());
151:   return 0;
152: }

154: /*TEST

156:    test:
157:      suffix: 1
158:      filter: grep -v " type:" | grep -v "Mat Object"
159:      args: -view
160:      diff_args: -j

162:    test:
163:      output_file: output/ex132_1.out
164:      requires: cuda
165:      suffix: 1_cuda
166:      filter: grep -v " type:" | grep -v "Mat Object"
167:      args: -view -mat_type aijcusparse
168:      diff_args: -j

170:    test:
171:      output_file: output/ex132_1.out
172:      requires: kokkos_kernels
173:      suffix: 1_kokkos
174:      filter: grep -v " type:" | grep -v "Mat Object"
175:      args: -view -mat_type aijkokkos
176:      diff_args: -j

178:    test:
179:      suffix: 2
180:      filter: grep -v " type:" | grep -v "Mat Object"
181:      args: -view -mat_nonsym
182:      diff_args: -j

184:    test:
185:      output_file: output/ex132_2.out
186:      requires: cuda
187:      suffix: 2_cuda
188:      filter: grep -v " type:" | grep -v "Mat Object"
189:      args: -view -mat_type aijcusparse -mat_nonsym
190:      diff_args: -j

192:    test:
193:      output_file: output/ex132_2.out
194:      requires: kokkos_kernels
195:      suffix: 2_kokkos
196:      filter: grep -v " type:" | grep -v "Mat Object"
197:      args: -view -mat_type aijkokkos -mat_nonsym
198:      diff_args: -j

200:    test:
201:      nsize: 2
202:      suffix: 1_par
203:      filter: grep -v " type:" | grep -v "Mat Object"
204:      args: -view
205:      diff_args: -j

207:    test:
208:      nsize: 2
209:      output_file: output/ex132_1_par.out
210:      requires: cuda
211:      suffix: 1_par_cuda
212:      filter: grep -v " type:" | grep -v "Mat Object"
213:      args: -view -mat_type aijcusparse
214:      diff_args: -j

216:    test:
217:      nsize: 2
218:      output_file: output/ex132_1_par.out
219:      requires: kokkos_kernels
220:      suffix: 1_par_kokkos
221:      filter: grep -v " type:" | grep -v "Mat Object"
222:      args: -view -mat_type aijkokkos
223:      diff_args: -j

225:    test:
226:      nsize: 2
227:      suffix: 2_par
228:      filter: grep -v " type:" | grep -v "Mat Object"
229:      args: -view -mat_nonsym
230:      diff_args: -j

232:    test:
233:      nsize: 2
234:      output_file: output/ex132_2_par.out
235:      requires: cuda
236:      suffix: 2_par_cuda
237:      filter: grep -v " type:" | grep -v "Mat Object"
238:      args: -view -mat_type aijcusparse -mat_nonsym
239:      diff_args: -j

241:    testset:
242:      nsize: 2
243:      output_file: output/ex132_2_par.out
244:      requires: kokkos_kernels
245:      filter: grep -v " type:" | grep -v "Mat Object"
246:      args: -view -mat_type aijkokkos -mat_nonsym
247:      diff_args: -j
248:      test:
249:        suffix: 2_par_kokkos_no_gpu_aware
250:        args: -use_gpu_aware_mpi 0
251:      test:
252:        requires: defined(HAVE_MPI_GPU_AWARE)
253:        suffix: 2_par_kokkos_gpu_aware
254:        args: -use_gpu_aware_mpi 1

256: TEST*/