1*ebf8cefbSJunchao Zhang static char help[] = "Tests MatConvert from AIJ to MATIS with a block size greater than 1.\n"; 2*ebf8cefbSJunchao Zhang 3*ebf8cefbSJunchao Zhang #include <petscmat.h> 4*ebf8cefbSJunchao Zhang int main(int argc,char **args) 5*ebf8cefbSJunchao Zhang { 6*ebf8cefbSJunchao Zhang Mat A,B; 7*ebf8cefbSJunchao Zhang char file[PETSC_MAX_PATH_LEN]; 8*ebf8cefbSJunchao Zhang PetscViewer fd; 9*ebf8cefbSJunchao Zhang PetscBool flg,equal; 10*ebf8cefbSJunchao Zhang 11*ebf8cefbSJunchao Zhang PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 12*ebf8cefbSJunchao Zhang 13*ebf8cefbSJunchao Zhang /* Load an AIJ matrix */ 14*ebf8cefbSJunchao Zhang PetscCall(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg)); 15*ebf8cefbSJunchao Zhang PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option"); 16*ebf8cefbSJunchao Zhang PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd)); 17*ebf8cefbSJunchao Zhang PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 18*ebf8cefbSJunchao Zhang PetscCall(MatSetFromOptions(A)); 19*ebf8cefbSJunchao Zhang PetscCall(MatLoad(A,fd)); 20*ebf8cefbSJunchao Zhang 21*ebf8cefbSJunchao Zhang /* Convert it to MATIS */ 22*ebf8cefbSJunchao Zhang PetscCall(MatConvert(A,MATIS,MAT_INITIAL_MATRIX,&B)); 23*ebf8cefbSJunchao Zhang 24*ebf8cefbSJunchao Zhang /* Check they are equal */ 25*ebf8cefbSJunchao Zhang PetscCall(MatEqual(A,B,&equal)); 26*ebf8cefbSJunchao Zhang PetscCheck(equal,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"A and B are not equal"); 27*ebf8cefbSJunchao Zhang 28*ebf8cefbSJunchao Zhang PetscCall(MatDestroy(&A)); 29*ebf8cefbSJunchao Zhang PetscCall(MatDestroy(&B)); 30*ebf8cefbSJunchao Zhang PetscCall(PetscViewerDestroy(&fd)); 31*ebf8cefbSJunchao Zhang PetscCall(PetscFinalize()); 32*ebf8cefbSJunchao Zhang } 33*ebf8cefbSJunchao Zhang 34*ebf8cefbSJunchao Zhang /*TEST 35*ebf8cefbSJunchao Zhang test: 36*ebf8cefbSJunchao Zhang requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 37*ebf8cefbSJunchao Zhang args: -mat_type aij -matload_block_size {{1 2}} -f ${DATAFILESPATH}/matrices/smallbs2 38*ebf8cefbSJunchao Zhang output_file: output/empty.out 39*ebf8cefbSJunchao Zhang 40*ebf8cefbSJunchao Zhang TEST*/ 41*ebf8cefbSJunchao Zhang 42