Actual source code: ex197.c

  1: static char help[] = "Test MatMultHermitianTranspose() and MatMultHermitianTransposeAdd().\n\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **args)
  6: {
  7:   Mat         A, B, C;
  8:   Vec         x, y, ys;
  9:   PetscInt    i, j;
 10:   PetscScalar v;
 11:   PetscBool   flg;

 14:   PetscInitialize(&argc, &args, (char *)0, help);
 15:   MatCreate(PETSC_COMM_WORLD, &A);
 16:   MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 2, 2);
 17:   MatSetType(A, MATAIJ);
 18:   MatSetFromOptions(A);
 19:   MatSetUp(A);

 21:   i = 0;
 22:   j = 0;
 23:   v = 2.0;
 24:   MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES);
 25:   i = 0;
 26:   j = 1;
 27:   v = 3.0 + 4.0 * PETSC_i;
 28:   MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES);
 29:   i = 1;
 30:   j = 0;
 31:   v = 5.0 + 6.0 * PETSC_i;
 32:   MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES);
 33:   i = 1;
 34:   j = 1;
 35:   v = 7.0 + 8.0 * PETSC_i;
 36:   MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES);
 37:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 38:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);

 40:   /* Create vectors */
 41:   VecCreate(PETSC_COMM_WORLD, &y);
 42:   VecSetSizes(y, PETSC_DECIDE, 2);
 43:   VecSetFromOptions(y);
 44:   VecDuplicate(y, &ys);
 45:   VecDuplicate(y, &x);

 47:   i = 0;
 48:   v = 10.0 + 11.0 * PETSC_i;
 49:   VecSetValues(x, 1, &i, &v, INSERT_VALUES);
 50:   i = 1;
 51:   v = 100.0 + 120.0 * PETSC_i;
 52:   VecSetValues(x, 1, &i, &v, INSERT_VALUES);
 53:   VecAssemblyBegin(x);
 54:   VecAssemblyEnd(x);

 56:   MatMultHermitianTranspose(A, x, y);
 57:   VecView(y, PETSC_VIEWER_STDOUT_WORLD);
 58:   MatMultHermitianTransposeAdd(A, x, y, ys);
 59:   VecView(ys, PETSC_VIEWER_STDOUT_WORLD);

 61:   MatHermitianTranspose(A, MAT_INITIAL_MATRIX, &B);
 62:   MatCreateHermitianTranspose(A, &C);
 63:   MatMultHermitianTransposeEqual(B, C, 4, &flg);
 65:   MatMultHermitianTransposeAddEqual(B, C, 4, &flg);
 67:   MatDestroy(&C);
 68:   MatDestroy(&B);

 70:   MatDestroy(&A);

 72:   VecDestroy(&x);
 73:   VecDestroy(&y);
 74:   VecDestroy(&ys);
 75:   PetscFinalize();
 76:   return 0;
 77: }

 79: /*TEST

 81:    build:
 82:       requires: complex
 83:    test:

 85:    test:
 86:       suffix: 2
 87:       nsize: 2

 89: TEST*/