Actual source code: convert.c

petsc-3.8.0 2017-09-26
Report Typos and Errors

  2:  #include <petsc/private/matimpl.h>

  4: /*
  5:   MatConvert_Basic - Converts from any input format to another format. For
  6:   parallel formats, the new matrix distribution is determined by PETSc.

  8:   Does not do preallocation so in general will be slow
  9:  */
 10: PETSC_INTERN PetscErrorCode MatConvert_Basic(Mat mat, MatType newtype,MatReuse reuse,Mat *newmat)
 11: {
 12:   Mat               M;
 13:   const PetscScalar *vwork;
 14:   PetscErrorCode    ierr;
 15:   PetscInt          nz,i,m,n,rstart,rend,lm,ln;
 16:   const PetscInt    *cwork;
 17:   PetscBool         isSBAIJ;

 20:   PetscObjectTypeCompare((PetscObject)mat,MATSEQSBAIJ,&isSBAIJ);
 21:   if (!isSBAIJ) {
 22:     PetscObjectTypeCompare((PetscObject)mat,MATMPISBAIJ,&isSBAIJ);
 23:   }
 24:   if (isSBAIJ) SETERRQ(PetscObjectComm((PetscObject)mat),PETSC_ERR_SUP,"Cannot convert from SBAIJ matrix since cannot obtain entire rows of matrix");


 27:   MatGetSize(mat,&m,&n);
 28:   MatGetLocalSize(mat,&lm,&ln);

 30:   if (ln == n) ln = PETSC_DECIDE; /* try to preserve column ownership */

 32:   MatCreate(PetscObjectComm((PetscObject)mat),&M);
 33:   MatSetSizes(M,lm,ln,m,n);
 34:   MatSetBlockSizesFromMats(M,mat,mat);
 35:   MatSetType(M,newtype);

 37:   MatSeqDenseSetPreallocation(M,NULL);
 38:   MatMPIDenseSetPreallocation(M,NULL);
 39:   MatSetUp(M);
 40:   MatSetOption(M,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
 41:   MatSetOption(M,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);

 43:     PetscObjectTypeCompare((PetscObject)M,MATSEQSBAIJ,&isSBAIJ);
 44:   if (!isSBAIJ) {
 45:     PetscObjectTypeCompare((PetscObject)M,MATMPISBAIJ,&isSBAIJ);
 46:   }
 47:   if (isSBAIJ) {
 48:     MatSetOption(M,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);
 49:   }

 51:   MatGetOwnershipRange(mat,&rstart,&rend);
 52:   for (i=rstart; i<rend; i++) {
 53:     MatGetRow(mat,i,&nz,&cwork,&vwork);
 54:     MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);
 55:     MatRestoreRow(mat,i,&nz,&cwork,&vwork);
 56:   }
 57:   MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);
 58:   MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);

 60:   if (reuse == MAT_INPLACE_MATRIX) {
 61:     MatHeaderReplace(mat,&M);
 62:   } else {
 63:     *newmat = M;
 64:   }
 65:   return(0);
 66: }