1 static char help[] = "Test some operations of SeqDense matrices with an LDA larger than M.\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc, char **argv) 6 { 7 Mat A, B; 8 PetscScalar *a, *b; 9 PetscInt n = 4, lda = 5, i; 10 11 PetscFunctionBeginUser; 12 PetscCall(PetscInitialize(&argc, &argv, 0, help)); 13 PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL)); 14 PetscCall(PetscOptionsGetInt(NULL, NULL, "-lda", &lda, NULL)); 15 PetscCheck(lda >= n, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "lda %" PetscInt_FMT " < n %" PetscInt_FMT, lda, n); 16 17 /* 18 * Create two identical matrices (MatDuplicate does not preserve lda) 19 */ 20 PetscCall(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 PetscCall(MatCreate(PETSC_COMM_SELF, &A)); 28 PetscCall(MatSetSizes(A, n, n, n, n)); 29 PetscCall(MatSetType(A, MATSEQDENSE)); 30 PetscCall(MatSeqDenseSetPreallocation(A, a)); 31 PetscCall(MatDenseSetLDA(A, lda)); 32 33 PetscCall(MatCreate(PETSC_COMM_SELF, &B)); 34 PetscCall(MatSetSizes(B, n, n, n, n)); 35 PetscCall(MatSetType(B, MATSEQDENSE)); 36 PetscCall(MatSeqDenseSetPreallocation(B, b)); 37 PetscCall(MatDenseSetLDA(B, lda)); 38 39 PetscCall(MatView(A, NULL)); 40 PetscCall(MatConjugate(A)); 41 PetscCall(MatView(A, NULL)); 42 PetscCall(MatRealPart(A)); 43 PetscCall(MatView(A, NULL)); 44 PetscCall(MatImaginaryPart(B)); 45 PetscCall(MatView(B, NULL)); 46 47 PetscCall(PetscFree2(a, b)); 48 PetscCall(MatDestroy(&A)); 49 PetscCall(MatDestroy(&B)); 50 PetscCall(PetscFinalize()); 51 return 0; 52 } 53 54 /*TEST 55 56 build: 57 requires: complex 58 59 test: 60 61 TEST*/ 62