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 PetscCall(PetscInitialize(&argc,&argv,0,help)); 12 PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 13 PetscCall(PetscOptionsGetInt(NULL,NULL,"-lda",&lda,NULL)); 14 PetscCheck(lda>=n,PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"lda %" PetscInt_FMT " < n %" PetscInt_FMT,lda,n); 15 16 /* 17 * Create two identical matrices (MatDuplicate does not preserve lda) 18 */ 19 PetscCall(PetscCalloc2(lda*n,&a,lda*n,&b)); 20 for (i=0; i<n; i++) { 21 a[i+i*lda] = 1.0+2.0*PETSC_i; 22 if (i>0) a[i+(i-1)*lda] = 3.0-0.5*PETSC_i; 23 b[i+i*lda] = 1.0+2.0*PETSC_i; 24 if (i>0) b[i+(i-1)*lda] = 3.0-0.5*PETSC_i; 25 } 26 PetscCall(MatCreate(PETSC_COMM_SELF,&A)); 27 PetscCall(MatSetSizes(A,n,n,n,n)); 28 PetscCall(MatSetType(A,MATSEQDENSE)); 29 PetscCall(MatSeqDenseSetPreallocation(A,a)); 30 PetscCall(MatDenseSetLDA(A,lda)); 31 32 PetscCall(MatCreate(PETSC_COMM_SELF,&B)); 33 PetscCall(MatSetSizes(B,n,n,n,n)); 34 PetscCall(MatSetType(B,MATSEQDENSE)); 35 PetscCall(MatSeqDenseSetPreallocation(B,b)); 36 PetscCall(MatDenseSetLDA(B,lda)); 37 38 PetscCall(MatView(A,NULL)); 39 PetscCall(MatConjugate(A)); 40 PetscCall(MatView(A,NULL)); 41 PetscCall(MatRealPart(A)); 42 PetscCall(MatView(A,NULL)); 43 PetscCall(MatImaginaryPart(B)); 44 PetscCall(MatView(B,NULL)); 45 46 PetscCall(PetscFree2(a,b)); 47 PetscCall(MatDestroy(&A)); 48 PetscCall(MatDestroy(&B)); 49 PetscCall(PetscFinalize()); 50 return 0; 51 } 52 53 /*TEST 54 55 build: 56 requires: complex 57 58 test: 59 60 TEST*/ 61