1*01a81e61SBarry Smith #define PETSCKSP_DLL 2*01a81e61SBarry Smith 3*01a81e61SBarry Smith #include "petscmat.h" 4*01a81e61SBarry Smith 5*01a81e61SBarry Smith /* 6*01a81e61SBarry Smith MatDumpSPAI - Dumps a PETSc matrix to a file in an ASCII format 7*01a81e61SBarry Smith suitable for the SPAI code of Stephen Barnard to solve. This routine 8*01a81e61SBarry Smith is simply here to allow testing of matrices directly with the SPAI 9*01a81e61SBarry Smith code, rather then through the PETSc interface. 10*01a81e61SBarry Smith 11*01a81e61SBarry Smith */ 12*01a81e61SBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT MatDumpSPAI(Mat A,FILE *file) 13*01a81e61SBarry Smith { 14*01a81e61SBarry Smith const PetscScalar *vals; 15*01a81e61SBarry Smith PetscErrorCode ierr; 16*01a81e61SBarry Smith int i,j,n,size,nz; 17*01a81e61SBarry Smith const int *cols; 18*01a81e61SBarry Smith MPI_Comm comm; 19*01a81e61SBarry Smith 20*01a81e61SBarry Smith PetscObjectGetComm((PetscObject)A,&comm); 21*01a81e61SBarry Smith 22*01a81e61SBarry Smith MPI_Comm_size(comm,&size); 23*01a81e61SBarry Smith if (size > 1) SETERRQ(PETSC_ERR_SUP,"Only single processor dumps"); 24*01a81e61SBarry Smith 25*01a81e61SBarry Smith ierr = MatGetSize(A,&n,&n);CHKERRQ(ierr); 26*01a81e61SBarry Smith 27*01a81e61SBarry Smith /* print the matrix */ 28*01a81e61SBarry Smith fprintf(file,"%d\n",n); 29*01a81e61SBarry Smith for (i=0; i<n; i++) { 30*01a81e61SBarry Smith ierr = MatGetRow(A,i,&nz,&cols,&vals);CHKERRQ(ierr); 31*01a81e61SBarry Smith for (j=0; j<nz; j++) { 32*01a81e61SBarry Smith fprintf(file,"%d %d %16.14e\n",i+1,cols[j]+1,vals[j]); 33*01a81e61SBarry Smith } 34*01a81e61SBarry Smith ierr = MatRestoreRow(A,i,&nz,&cols,&vals);CHKERRQ(ierr); 35*01a81e61SBarry Smith } 36*01a81e61SBarry Smith 37*01a81e61SBarry Smith PetscFunctionReturn(0); 38*01a81e61SBarry Smith } 39*01a81e61SBarry Smith 40*01a81e61SBarry Smith PetscErrorCode PETSCKSP_DLLEXPORT VecDumpSPAI(Vec b,FILE *file) 41*01a81e61SBarry Smith { 42*01a81e61SBarry Smith PetscErrorCode ierr; 43*01a81e61SBarry Smith int n,i; 44*01a81e61SBarry Smith PetscScalar *array; 45*01a81e61SBarry Smith 46*01a81e61SBarry Smith ierr = VecGetSize(b,&n);CHKERRQ(ierr); 47*01a81e61SBarry Smith ierr = VecGetArray(b,&array);CHKERRQ(ierr); 48*01a81e61SBarry Smith 49*01a81e61SBarry Smith fprintf(file,"%d\n",n); 50*01a81e61SBarry Smith for (i=0; i<n; i++) { 51*01a81e61SBarry Smith fprintf(file,"%d %16.14e\n",i+1,array[i]); 52*01a81e61SBarry Smith } 53*01a81e61SBarry Smith 54*01a81e61SBarry Smith PetscFunctionReturn(0); 55*01a81e61SBarry Smith } 56