Actual source code: ex43.c

  1: static char help[] = "Tests VecMDot(),VecDot(),VecMTDot(), and VecTDot()\n";

  3: #include <petscvec.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Vec         *V, t;
  8:   PetscInt     i, j, reps, n = 15, k = 6;
  9:   PetscRandom  rctx;
 10:   PetscScalar *val_dot, *val_mdot, *tval_dot, *tval_mdot;

 12:   PetscFunctionBeginUser;
 13:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 14:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
 15:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-k", &k, NULL));
 16:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test with %" PetscInt_FMT " random vectors of length %" PetscInt_FMT "\n", k, n));
 17:   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
 18:   PetscCall(PetscRandomSetFromOptions(rctx));
 19: #if defined(PETSC_USE_COMPLEX)
 20:   PetscCall(PetscRandomSetInterval(rctx, -1. + 4. * PETSC_i, 1. + 5. * PETSC_i));
 21: #else
 22:   PetscCall(PetscRandomSetInterval(rctx, -1., 1.));
 23: #endif
 24:   PetscCall(VecCreate(PETSC_COMM_WORLD, &t));
 25:   PetscCall(VecSetSizes(t, n, PETSC_DECIDE));
 26:   PetscCall(VecSetFromOptions(t));
 27:   PetscCall(VecDuplicateVecs(t, k, &V));
 28:   PetscCall(VecSetRandom(t, rctx));
 29:   PetscCall(VecViewFromOptions(t, NULL, "-t_view"));
 30:   PetscCall(PetscMalloc1(k, &val_dot));
 31:   PetscCall(PetscMalloc1(k, &val_mdot));
 32:   PetscCall(PetscMalloc1(k, &tval_dot));
 33:   PetscCall(PetscMalloc1(k, &tval_mdot));
 34:   for (i = 0; i < k; i++) PetscCall(VecSetRandom(V[i], rctx));
 35:   for (reps = 0; reps < 20; reps++) {
 36:     for (i = 1; i < k; i++) {
 37:       PetscCall(VecMDot(t, i, V, val_mdot));
 38:       PetscCall(VecMTDot(t, i, V, tval_mdot));
 39:       for (j = 0; j < i; j++) {
 40:         PetscCall(VecDot(t, V[j], &val_dot[j]));
 41:         PetscCall(VecTDot(t, V[j], &tval_dot[j]));
 42:       }
 43:       /* Check result */
 44:       for (j = 0; j < i; j++) {
 45:         if (PetscAbsScalar(val_mdot[j] - val_dot[j]) / PetscAbsScalar(val_dot[j]) > 1e-5) {
 46:           PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%" PetscInt_FMT ", j=%" PetscInt_FMT ", val_mdot[j]=%g, val_dot[j]=%g\n", i, j, (double)PetscAbsScalar(val_mdot[j]), (double)PetscAbsScalar(val_dot[j])));
 47:           break;
 48:         }
 49:         if (PetscAbsScalar(tval_mdot[j] - tval_dot[j]) / PetscAbsScalar(tval_dot[j]) > 1e-5) {
 50:           PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%" PetscInt_FMT ", j=%" PetscInt_FMT ", tval_mdot[j]=%g, tval_dot[j]=%g\n", i, j, (double)PetscAbsScalar(tval_mdot[j]), (double)PetscAbsScalar(tval_dot[j])));
 51:           break;
 52:         }
 53:       }
 54:     }
 55:   }
 56:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test completed successfully!\n"));
 57:   PetscCall(PetscFree(val_dot));
 58:   PetscCall(PetscFree(val_mdot));
 59:   PetscCall(PetscFree(tval_dot));
 60:   PetscCall(PetscFree(tval_mdot));
 61:   PetscCall(VecDestroyVecs(k, &V));
 62:   PetscCall(VecDestroy(&t));
 63:   PetscCall(PetscRandomDestroy(&rctx));
 64:   PetscCall(PetscFinalize());
 65:   return 0;
 66: }

 68: /*TEST

 70:    test:

 72:    testset:
 73:       output_file: output/ex43_1.out

 75:       test:
 76:          suffix: cuda
 77:          args: -vec_type cuda -random_type curand
 78:          requires: cuda

 80:       test:
 81:          suffix: kokkos
 82:          args: -vec_type kokkos
 83:          requires: kokkos_kernels

 85:       test:
 86:          suffix: hip
 87:          args: -vec_type hip
 88:          requires: hip
 89: TEST*/