1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Tests binary I/O of matrices and illustrates user-defined event logging.\n\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown #include <petscmat.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown /* Note: Most applications would not read and write the same matrix within 7c4762a1bSJed Brown the same program. This example is intended only to demonstrate 8c4762a1bSJed Brown both input and output. */ 9c4762a1bSJed Brown 10c4762a1bSJed Brown int main(int argc,char **args) 11c4762a1bSJed Brown { 12c4762a1bSJed Brown Mat C; 13c4762a1bSJed Brown PetscScalar v; 14c4762a1bSJed Brown PetscInt i,j,Ii,J,Istart,Iend,N,m = 4,n = 4; 15c4762a1bSJed Brown PetscMPIInt rank,size; 16c4762a1bSJed Brown PetscErrorCode ierr; 17c4762a1bSJed Brown PetscViewer viewer; 18c4762a1bSJed Brown #if defined(PETSC_USE_LOG) 19c4762a1bSJed Brown PetscLogEvent MATRIX_GENERATE,MATRIX_READ; 20c4762a1bSJed Brown #endif 21c4762a1bSJed Brown 22c4762a1bSJed Brown ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 23*5f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 24*5f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 25*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL)); 26*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL)); 27c4762a1bSJed Brown N = m*n; 28c4762a1bSJed Brown 29c4762a1bSJed Brown /* PART 1: Generate matrix, then write it in binary format */ 30c4762a1bSJed Brown 31*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventRegister("Generate Matrix",0,&MATRIX_GENERATE)); 32*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventBegin(MATRIX_GENERATE,0,0,0,0)); 33c4762a1bSJed Brown 34c4762a1bSJed Brown /* Generate matrix */ 35*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatCreate(PETSC_COMM_WORLD,&C)); 36*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,N,N)); 37*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetFromOptions(C)); 38*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetUp(C)); 39*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatGetOwnershipRange(C,&Istart,&Iend)); 40c4762a1bSJed Brown for (Ii=Istart; Ii<Iend; Ii++) { 41c4762a1bSJed Brown v = -1.0; i = Ii/n; j = Ii - i*n; 42*5f80ce2aSJacob Faibussowitsch if (i>0) {J = Ii - n; CHKERRQ(MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES));} 43*5f80ce2aSJacob Faibussowitsch if (i<m-1) {J = Ii + n; CHKERRQ(MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES));} 44*5f80ce2aSJacob Faibussowitsch if (j>0) {J = Ii - 1; CHKERRQ(MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES));} 45*5f80ce2aSJacob Faibussowitsch if (j<n-1) {J = Ii + 1; CHKERRQ(MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES));} 46*5f80ce2aSJacob Faibussowitsch v = 4.0; CHKERRQ(MatSetValues(C,1,&Ii,1,&Ii,&v,ADD_VALUES)); 47c4762a1bSJed Brown } 48*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY)); 49*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY)); 50*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatView(C,PETSC_VIEWER_STDOUT_WORLD)); 51c4762a1bSJed Brown 52*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"writing matrix in binary to matrix.dat ...\n")); 53*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_WRITE,&viewer)); 54*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatView(C,viewer)); 55*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 56*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&C)); 57*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventEnd(MATRIX_GENERATE,0,0,0,0)); 58c4762a1bSJed Brown 59c4762a1bSJed Brown /* PART 2: Read in matrix in binary format */ 60c4762a1bSJed Brown 61c4762a1bSJed Brown /* All processors wait until test matrix has been dumped */ 62*5f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Barrier(PETSC_COMM_WORLD)); 63c4762a1bSJed Brown 64*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventRegister("Read Matrix",0,&MATRIX_READ)); 65*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventBegin(MATRIX_READ,0,0,0,0)); 66*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"reading matrix in binary from matrix.dat ...\n")); 67*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&viewer)); 68*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatCreate(PETSC_COMM_WORLD,&C)); 69*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatLoad(C,viewer)); 70*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 71*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscLogEventEnd(MATRIX_READ,0,0,0,0)); 72*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatView(C,PETSC_VIEWER_STDOUT_WORLD)); 73c4762a1bSJed Brown 74c4762a1bSJed Brown /* Free data structures */ 75*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&C)); 76c4762a1bSJed Brown 77c4762a1bSJed Brown ierr = PetscFinalize(); 78c4762a1bSJed Brown return ierr; 79c4762a1bSJed Brown } 80c4762a1bSJed Brown 81c4762a1bSJed Brown /*TEST 82c4762a1bSJed Brown 83c4762a1bSJed Brown test: 84c4762a1bSJed Brown filter: grep -v "MPI processes" 85c4762a1bSJed Brown 86c4762a1bSJed Brown TEST*/ 87