Actual source code: ex64.c

  1: static char help[] = "Tests VecLog().\n\n";

  3: #include <petscvec.h>

  5: static PetscErrorCode CheckLog(Vec v, PetscInt n, PetscScalar *arr, PetscScalar value)
  6: {
  7:   const PetscReal    rtol = 1e-10, atol = PETSC_SMALL;
  8:   const PetscScalar *varr;

 10:   PetscFunctionBegin;
 11:   PetscCall(VecSet(v, value));
 12:   PetscCall(VecViewFromOptions(v, NULL, "-vec_view"));
 13:   PetscCall(VecLog(v));
 14:   PetscCall(VecViewFromOptions(v, NULL, "-vec_view"));

 16:   for (PetscInt i = 0; i < n; ++i) arr[i] = PetscLogScalar(value);
 17:   PetscCall(VecGetArrayRead(v, &varr));
 18:   for (PetscInt i = 0; i < n; ++i) {
 19:     const PetscScalar lhs = varr[i];
 20:     const PetscScalar rhs = arr[i];

 22:     if (!PetscIsCloseAtTolScalar(lhs, rhs, rtol, atol)) {
 23:       const PetscReal lhs_r = PetscRealPart(lhs);
 24:       const PetscReal lhs_i = PetscImaginaryPart(lhs);
 25:       const PetscReal rhs_r = PetscRealPart(rhs);
 26:       const PetscReal rhs_i = PetscImaginaryPart(rhs);

 28:       PetscCheck(PetscIsCloseAtTol(lhs_r, rhs_r, rtol, atol), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Real component actual[%" PetscInt_FMT "] %g != expected[%" PetscInt_FMT "] %g", i, (double)lhs_r, i, (double)rhs_r);
 29:       PetscCheck(PetscIsCloseAtTol(lhs_i, rhs_i, rtol, atol), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Imaginary component actual[%" PetscInt_FMT "] %g != expected[%" PetscInt_FMT "] %g", i, (double)lhs_i, i, (double)rhs_i);
 30:     }
 31:   }
 32:   PetscCall(VecRestoreArrayRead(v, &varr));
 33:   PetscFunctionReturn(PETSC_SUCCESS);
 34: }

 36: int main(int argc, char **argv)
 37: {
 38:   Vec          v;
 39:   PetscInt     n;
 40:   PetscScalar *arr;

 42:   PetscFunctionBeginUser;
 43:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));

 45:   PetscCall(VecCreate(PETSC_COMM_WORLD, &v));
 46:   PetscCall(VecSetSizes(v, 10, PETSC_DECIDE));
 47:   PetscCall(VecSetFromOptions(v));

 49:   PetscCall(VecGetLocalSize(v, &n));
 50:   PetscCall(PetscMalloc1(n, &arr));

 52:   PetscCall(CheckLog(v, n, arr, 1.0));
 53:   PetscCall(CheckLog(v, n, arr, PetscExpScalar(1.0)));
 54:   PetscCall(CheckLog(v, n, arr, 10.0));

 56:   PetscCall(PetscFree(arr));
 57:   PetscCall(VecDestroy(&v));

 59:   PetscCall(PetscFinalize());
 60:   return 0;
 61: }

 63: /*TEST

 65:   testset:
 66:     output_file: ./output/empty.out
 67:     nsize: {{1 2}}
 68:     test:
 69:       suffix: standard
 70:       args: -vec_type standard
 71:     test:
 72:       suffix: viennacl
 73:       requires: viennacl
 74:       args: -vec_type viennacl
 75:     test:
 76:       suffix: cuda
 77:       requires: cuda
 78:       args: -vec_type cuda
 79:     test:
 80:       suffix: hip
 81:       requires: hip
 82:       args: -vec_type hip
 83:     test:
 84:       suffix: kokkos
 85:       requires: kokkos kokkos_kernels
 86:       args: -vec_type kokkos

 88: TEST*/