Actual source code: ex53.c

  1: static const char help[] = "Tests VecShift()\n\n";

  3: #include <petscvec.h>

  5: static PetscErrorCode CheckVecShift(Vec v, PetscInt n, PetscScalar *array_copy, PetscScalar shift)
  6: {
  7:   const PetscScalar *array;

  9:   PetscFunctionBegin;
 10:   for (PetscInt i = 0; i < n; ++i) array_copy[i] += shift;
 11:   PetscCall(VecShift(v, shift));
 12:   PetscCall(VecGetArrayRead(v, &array));
 13:   for (PetscInt i = 0; i < n; ++i) {
 14:     const PetscScalar actual = array[i], expected = array_copy[i];

 16:     PetscCheck(PetscIsCloseAtTolScalar(actual, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecShift() returned array[%" PetscInt_FMT "] %g + %gi != expected_array[%" PetscInt_FMT "] %g + %gi", i, (double)PetscRealPart(actual), (double)PetscImaginaryPart(actual), i, (double)PetscRealPart(expected), (double)PetscImaginaryPart(expected));
 17:   }
 18:   PetscCall(VecRestoreArrayRead(v, &array));
 19:   PetscFunctionReturn(PETSC_SUCCESS);
 20: }

 22: int main(int argc, char **argv)
 23: {
 24:   Vec                x;
 25:   PetscInt           n;
 26:   const PetscScalar *array;
 27:   PetscScalar       *array_copy;
 28:   PetscReal          norm_before, norm_after;
 29:   PetscBool          available;

 31:   PetscFunctionBeginUser;
 32:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));

 34:   PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
 35:   PetscCall(VecSetSizes(x, PETSC_DECIDE, 10));
 36:   PetscCall(VecSetFromOptions(x));

 38:   PetscCall(VecZeroEntries(x));

 40:   // get a copy of the vectors array, anything we do to the vector via VecShift we will also do
 41:   // to the copy, and hence they should always match
 42:   PetscCall(VecGetLocalSize(x, &n));
 43:   PetscCall(PetscMalloc1(n, &array_copy));
 44:   PetscCall(VecGetArrayRead(x, &array));
 45:   PetscCall(PetscArraycpy(array_copy, array, n));
 46:   PetscCall(VecRestoreArrayRead(x, &array));

 48:   PetscCall(CheckVecShift(x, n, array_copy, 0.0));
 49:   PetscCall(CheckVecShift(x, n, array_copy, 1.0));
 50:   PetscCall(CheckVecShift(x, n, array_copy, -1.0));
 51:   PetscCall(CheckVecShift(x, n, array_copy, 15.0));

 53:   PetscCall(VecNorm(x, NORM_2, &norm_before));
 54:   PetscCall(VecNormAvailable(x, NORM_2, &available, &norm_after));
 55:   PetscCheck(available, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecNormAvailable() returned FALSE right after calling VecNorm()");
 56:   // a shift of zero should not invalidate norms
 57:   PetscCall(CheckVecShift(x, n, array_copy, 0.0));
 58:   PetscCall(VecNormAvailable(x, NORM_2, &available, &norm_after));
 59:   PetscCheck(available, PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecNormAvailable() returned FALSE after calling VecShift() with a shift of 0.0!");
 60:   // these can be compared with equality as the number should not change *at all*
 61:   PetscCheck(norm_before == norm_after, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Norms differ before and after calling VecShift() with shift of 0.0: before %g after %g", (double)norm_before, (double)norm_after);

 63:   PetscCall(PetscFree(array_copy));
 64:   PetscCall(VecDestroy(&x));
 65:   PetscCall(PetscFinalize());
 66:   return 0;
 67: }

 69: /*TEST

 71:   testset:
 72:     output_file: ./output/empty.out
 73:     nsize: {{1 2}}
 74:     test:
 75:       suffix: standard
 76:     test:
 77:       requires: defined(PETSC_USE_SHARED_MEMORY)
 78:       args: -vec_type shared
 79:       suffix: shared
 80:     test:
 81:       requires: viennacl
 82:       args: -vec_type viennacl
 83:       suffix: viennacl
 84:     test:
 85:       requires: kokkos_kernels
 86:       args: -vec_type kokkos
 87:       suffix: kokkos
 88:     test:
 89:       requires: cuda
 90:       args: -vec_type cuda
 91:       suffix: cuda
 92:     test:
 93:       requires: hip
 94:       args: -vec_type hip
 95:       suffix: hip

 97: TEST*/