1: static const char help[] = "Test VecMAXPY() coefficient lifetime with short local vectors.\n";

  3: #include <petscvec.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Vec          x, y, ref, ys[9];
  8:   PetscScalar *alpha;
  9:   struct {
 10:     PetscReal   padding;
 11:     PetscScalar values[9];
 12:   }                 *host;
 13:   const PetscScalar *a, *r;
 14:   const PetscInt     counts[] = {1, 2, 3, 4, 7, 8, 9};
 15:   PetscInt           n        = 1;
 16:   PetscMPIInt        rank;
 17:   PetscBool          empty = PETSC_FALSE;

 19:   PetscFunctionBeginUser;
 20:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 21:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
 22:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-empty_rank", &empty, NULL));
 23:   if (empty && rank) n = 0;
 24:   PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
 25:   PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
 26:   PetscCall(VecSetFromOptions(x));
 27:   PetscCall(VecDuplicate(x, &y));
 28:   PetscCall(VecSet(y, 2));
 29:   PetscCall(VecCreate(PETSC_COMM_WORLD, &ref));
 30:   PetscCall(VecSetSizes(ref, n, PETSC_DECIDE));
 31:   PetscCall(VecSetType(ref, VECSTANDARD));
 32:   // A scalar array inside a C struct need not satisfy CUDA/HIP complex alignment.
 33:   PetscCall(PetscMalloc1(1, &host));
 34:   alpha = host->values;
 35:   for (PetscInt j = 0; j < 9; ++j) ys[j] = y;
 36:   for (size_t k = 0; k < PETSC_STATIC_ARRAY_LENGTH(counts); ++k) {
 37:     PetscScalar expected = 0;
 38:     PetscInt    nv       = counts[k];

 40:     PetscCall(VecSet(x, 0));
 41:     for (PetscInt i = 0; i < 4; ++i) {
 42:       for (PetscInt j = 0; j < nv; ++j) {
 43:         alpha[j] = i + j + 1;
 44: #if PetscDefined(USE_COMPLEX)
 45:         alpha[j] += PETSC_i * (j + 1);
 46: #endif
 47:         expected += 2 * alpha[j];
 48:       }
 49:       PetscCall(VecMAXPY(x, nv, alpha, ys));
 50:     }
 51:     for (PetscInt j = 0; j < nv; ++j) alpha[j] = -100;
 52:     PetscCall(VecSet(ref, expected));
 53:     PetscCall(VecGetArrayRead(x, &a));
 54:     PetscCall(VecGetArrayRead(ref, &r));
 55:     for (PetscInt i = 0; i < n; ++i) PetscCheck(PetscAbsScalar(a[i] - r[i]) <= 10 * PETSC_MACHINE_EPSILON * PetscMax(PetscAbsScalar(r[i]), 1), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incorrect VecMAXPY() result for %" PetscInt_FMT " coefficients", nv);
 56:     PetscCall(VecRestoreArrayRead(ref, &r));
 57:     PetscCall(VecRestoreArrayRead(x, &a));
 58:   }
 59:   PetscCall(PetscFree(host));
 60:   PetscCall(VecDestroy(&ref));
 61:   PetscCall(VecDestroy(&y));
 62:   PetscCall(VecDestroy(&x));
 63:   PetscCall(PetscFinalize());
 64:   return 0;
 65: }

 67: /*TEST

 69:   testset:
 70:     nsize: {{1 2}}
 71:     output_file: output/empty.out
 72:     args: -empty_rank {{0 1}}
 73:     test:
 74:       suffix: cuda
 75:       requires: cuda
 76:       args: -vec_type cuda
 77:     test:
 78:       suffix: hip
 79:       requires: hip
 80:       args: -vec_type hip

 82: TEST*/