Actual source code: convert.c


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

  4: /*
  5:   MatConvert_Basic - Converts from any input format to another format.
  6:   Does not do preallocation so in general will be slow
  7:  */
  8: PETSC_INTERN PetscErrorCode MatConvert_Basic(Mat mat,MatType newtype,MatReuse reuse,Mat *newmat)
  9: {
 10:   Mat               M;
 11:   const PetscScalar *vwork;
 12:   PetscInt          i,rstart,rend,nz;
 13:   const PetscInt    *cwork;
 14:   PetscBool         isSBAIJ;

 16:   if (!mat->ops->getrow) { /* missing get row, use matvecs */
 17:     MatConvert_Shell(mat,newtype,reuse,newmat);
 18:     return 0;
 19:   }
 20:   PetscObjectTypeCompare((PetscObject)mat,MATSEQSBAIJ,&isSBAIJ);
 21:   if (!isSBAIJ) {
 22:     PetscObjectTypeCompare((PetscObject)mat,MATMPISBAIJ,&isSBAIJ);
 23:   }

 26:   if (reuse == MAT_REUSE_MATRIX) {
 27:     M = *newmat;
 28:   } else {
 29:     PetscInt m,n,lm,ln;
 30:     MatGetSize(mat,&m,&n);
 31:     MatGetLocalSize(mat,&lm,&ln);
 32:     MatCreate(PetscObjectComm((PetscObject)mat),&M);
 33:     MatSetSizes(M,lm,ln,m,n);
 34:     MatSetBlockSizesFromMats(M,mat,mat);
 35:     MatSetType(M,newtype);
 36:     MatSetUp(M);

 38:     MatSetOption(M,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
 39:     MatSetOption(M,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
 40:     PetscObjectTypeCompare((PetscObject)M,MATSEQSBAIJ,&isSBAIJ);
 41:     if (!isSBAIJ) {
 42:       PetscObjectTypeCompare((PetscObject)M,MATMPISBAIJ,&isSBAIJ);
 43:     }
 44:     if (isSBAIJ) {
 45:       MatSetOption(M,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);
 46:     }
 47:   }

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

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