Actual source code: ex66.c

  1: const char help[] = "Test VecPointwiseSign()";

  3: #include <petscvec.h>

  5: int main(int argc, char **argv)
  6: {
  7:   MPI_Comm    comm;
  8:   PetscRandom rand;
  9:   PetscInt    n = 12;
 10:   Vec         x, y, z, test;
 11:   VecSignMode types[] = {VEC_SIGN_ZERO_TO_ZERO, VEC_SIGN_ZERO_TO_SIGNED_ZERO, VEC_SIGN_ZERO_TO_SIGNED_UNIT};

 13:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 14:   comm = PETSC_COMM_WORLD;
 15:   PetscCall(PetscRandomCreate(comm, &rand));
 16:   PetscCall(PetscRandomSetInterval(rand, 0.5, 1.5));
 17:   PetscCall(VecCreate(comm, &x));
 18:   PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
 19:   PetscCall(VecSetFromOptions(x));
 20:   PetscCall(VecDuplicate(x, &y));
 21:   PetscCall(VecDuplicate(x, &z));
 22:   PetscCall(VecDuplicate(x, &test));
 23:   {
 24:     PetscScalar *_x;
 25:     PetscCall(VecGetArrayWrite(x, &_x));
 26:     for (PetscInt i = 0; i < n; i++) {
 27:       PetscReal r;

 29:       PetscCall(PetscRandomGetValueReal(rand, &r));
 30:       switch (i % 4) {
 31:       case 0:
 32:         _x[i] = -r;
 33:         break;
 34:       case 1:
 35:         _x[i] = -0.0;
 36:         break;
 37:       case 2:
 38:         _x[i] = 0.0;
 39:         break;
 40:       default:
 41:         _x[i] = r;
 42:         break;
 43:       }
 44:     }
 45:     PetscCall(VecRestoreArrayWrite(x, &_x));
 46:   }
 47:   for (size_t t = 0; t < PETSC_STATIC_ARRAY_LENGTH(types); t++) {
 48:     PetscScalar *_test;
 49:     VecSignMode  sign_type = types[t];

 51:     PetscCall(VecGetArrayWrite(test, &_test));
 52:     for (PetscInt i = 0; i < n; i++) {
 53:       switch (i % 4) {
 54:       case 0:
 55:         _test[i] = -1.0;
 56:         break;
 57:       case 1:
 58:         switch (sign_type) {
 59:         case VEC_SIGN_ZERO_TO_ZERO:
 60:           _test[i] = 0.0;
 61:           break;
 62:         case VEC_SIGN_ZERO_TO_SIGNED_ZERO:
 63:           _test[i] = -0.0;
 64:           break;
 65:         case VEC_SIGN_ZERO_TO_SIGNED_UNIT:
 66:           _test[i] = -1.0;
 67:           break;
 68:         }
 69:         break;
 70:       case 2:
 71:         switch (sign_type) {
 72:         case VEC_SIGN_ZERO_TO_ZERO:
 73:           _test[i] = 0.0;
 74:           break;
 75:         case VEC_SIGN_ZERO_TO_SIGNED_ZERO:
 76:           _test[i] = 0.0;
 77:           break;
 78:         case VEC_SIGN_ZERO_TO_SIGNED_UNIT:
 79:           _test[i] = 1.0;
 80:           break;
 81:         }
 82:         break;
 83:       default:
 84:         _test[i] = 1.0;
 85:         break;
 86:       }
 87:     }
 88:     PetscCall(VecRestoreArrayWrite(test, &_test));
 89:     for (PetscInt j = 0; j < 2; j++) {
 90:       Vec                out = j ? z : y;
 91:       const PetscScalar *_out;
 92:       const PetscScalar *_test;

 94:       PetscCall(VecCopy(x, y));
 95:       PetscCall(VecPointwiseSign(out, y, sign_type));
 96:       PetscCall(VecGetArrayRead(out, &_out));
 97:       PetscCall(VecGetArrayRead(test, &_test));
 98:       for (PetscInt i = 0; i < n; i++) {
 99:         PetscScalar _o = _out[i];
100:         PetscScalar _t = _test[i];

102:         PetscCheck(PetscImaginaryPart(_o) == 0.0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Nonzero imaginary part");
103:         PetscCheck(PetscAbsScalar(_o) == PetscAbsScalar(_t), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Unexpected magnitude");
104:         PetscCheck(PetscCopysignReal(1.0, PetscRealPart(_o)) == PetscCopysignReal(1.0, PetscRealPart(_t)), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Unexpected sign");
105:       }
106:       PetscCall(VecRestoreArrayRead(test, &_test));
107:       PetscCall(VecRestoreArrayRead(out, &_out));
108:     }
109:   }
110:   PetscCall(VecDestroy(&test));
111:   PetscCall(VecDestroy(&z));
112:   PetscCall(VecDestroy(&y));
113:   PetscCall(VecDestroy(&x));
114:   PetscCall(PetscRandomDestroy(&rand));
115:   PetscCall(PetscFinalize());
116:   return 0;
117: }

119: /*TEST

121:   test:
122:     nsize: {{1 2}}
123:     suffix: 0
124:     output_file: output/empty.out

126:   test:
127:     requires: cuda
128:     nsize: {{1 2}}
129:     suffix: cuda
130:     output_file: output/empty.out
131:     args: -vec_type cuda

133:   test:
134:     requires: hip
135:     nsize: {{1 2}}
136:     suffix: hip
137:     output_file: output/empty.out
138:     args: -vec_type hip

140: TEST*/