1 2 static char help[] = "Tests MatILUFactorSymbolic() on matrix with missing diagonal.\n\n"; 3 4 #include <petscmat.h> 5 #include <petscpc.h> 6 7 int main(int argc,char **args) 8 { 9 Mat C,A; 10 PetscInt i,j; 11 PetscErrorCode ierr; 12 PetscScalar v; 13 PC pc; 14 Vec xtmp; 15 16 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 17 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&C)); 18 CHKERRQ(MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,3,3)); 19 CHKERRQ(MatSetFromOptions(C)); 20 CHKERRQ(MatSetUp(C)); 21 CHKERRQ(VecCreateSeq(PETSC_COMM_WORLD,3,&xtmp)); 22 i = 0; j = 0; v = 4; 23 CHKERRQ(MatSetValues(C,1,&i,1,&j,&v,INSERT_VALUES)); 24 i = 0; j = 2; v = 1; 25 CHKERRQ(MatSetValues(C,1,&i,1,&j,&v,INSERT_VALUES)); 26 i = 1; j = 0; v = 1; 27 CHKERRQ(MatSetValues(C,1,&i,1,&j,&v,INSERT_VALUES)); 28 i = 1; j = 1; v = 4; 29 CHKERRQ(MatSetValues(C,1,&i,1,&j,&v,INSERT_VALUES)); 30 i = 2; j = 1; v = 1; 31 CHKERRQ(MatSetValues(C,1,&i,1,&j,&v,INSERT_VALUES)); 32 33 CHKERRQ(MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY)); 34 CHKERRQ(MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY)); 35 36 CHKERRQ(MatView(C,PETSC_VIEWER_STDOUT_WORLD)); 37 CHKERRQ(PCCreate(PETSC_COMM_WORLD,&pc)); 38 CHKERRQ(PCSetFromOptions(pc)); 39 CHKERRQ(PCSetOperators(pc,C,C)); 40 CHKERRQ(PCSetUp(pc)); 41 CHKERRQ(PCFactorGetMatrix(pc,&A)); 42 CHKERRQ(MatView(A,PETSC_VIEWER_STDOUT_WORLD)); 43 44 CHKERRQ(PCDestroy(&pc)); 45 CHKERRQ(VecDestroy(&xtmp)); 46 CHKERRQ(MatDestroy(&C)); 47 48 ierr = PetscFinalize(); 49 return ierr; 50 } 51