Actual source code: convert.c

 2:  #include src/mat/matimpl.h

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

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

 21:   MatGetSize(mat,&m,&n);
 22:   MatGetLocalSize(mat,&lm,&ln);

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

 26:   MatCreate(mat->comm,lm,ln,m,n,&M);
 27:   MatSetType(M,newtype);

 29:   MatGetOwnershipRange(mat,&rstart,&rend);
 30:   for (i=rstart; i<rend; i++) {
 31:     MatGetRow(mat,i,&nz,&cwork,&vwork);
 32:     MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);
 33:     MatRestoreRow(mat,i,&nz,&cwork,&vwork);
 34:   }
 35:   MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);
 36:   MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);

 38:   /* Fake support for "inplace" convert. */
 39:   if (*newmat == mat) {
 40:     MatDestroy(mat);
 41:   }
 42:   *newmat = M;
 43:   return(0);
 44: }