Actual source code: ex256.c

  1: static char help[] = "Test some operations of SeqDense matrices with an LDA larger than M.\n";

  3: #include <petscmat.h>

  5: int main(int argc, char **argv)
  6: {
  7:   Mat          A, B;
  8:   PetscScalar *a, *b;
  9:   PetscInt     n = 4, lda = 5, i;

 12:   PetscInitialize(&argc, &argv, 0, help);
 13:   PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
 14:   PetscOptionsGetInt(NULL, NULL, "-lda", &lda, NULL);

 17:   /*
 18:    * Create two identical matrices (MatDuplicate does not preserve lda)
 19:    */
 20:   PetscCalloc2(lda * n, &a, lda * n, &b);
 21:   for (i = 0; i < n; i++) {
 22:     a[i + i * lda] = 1.0 + 2.0 * PETSC_i;
 23:     if (i > 0) a[i + (i - 1) * lda] = 3.0 - 0.5 * PETSC_i;
 24:     b[i + i * lda] = 1.0 + 2.0 * PETSC_i;
 25:     if (i > 0) b[i + (i - 1) * lda] = 3.0 - 0.5 * PETSC_i;
 26:   }
 27:   MatCreate(PETSC_COMM_SELF, &A);
 28:   MatSetSizes(A, n, n, n, n);
 29:   MatSetType(A, MATSEQDENSE);
 30:   MatSeqDenseSetPreallocation(A, a);
 31:   MatDenseSetLDA(A, lda);

 33:   MatCreate(PETSC_COMM_SELF, &B);
 34:   MatSetSizes(B, n, n, n, n);
 35:   MatSetType(B, MATSEQDENSE);
 36:   MatSeqDenseSetPreallocation(B, b);
 37:   MatDenseSetLDA(B, lda);

 39:   MatView(A, NULL);
 40:   MatConjugate(A);
 41:   MatView(A, NULL);
 42:   MatRealPart(A);
 43:   MatView(A, NULL);
 44:   MatImaginaryPart(B);
 45:   MatView(B, NULL);

 47:   PetscFree2(a, b);
 48:   MatDestroy(&A);
 49:   MatDestroy(&B);
 50:   PetscFinalize();
 51:   return 0;
 52: }

 54: /*TEST

 56:    build:
 57:      requires: complex

 59:    test:

 61: TEST*/