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:   PetscFunctionBegin;
 17:   if (!mat->ops->getrow) { /* missing get row, use matvecs */
 18:     PetscCall(MatConvert_Shell(mat, newtype, reuse, newmat));
 19:     PetscFunctionReturn(PETSC_SUCCESS);
 20:   }
 21:   PetscCall(PetscObjectTypeCompare((PetscObject)mat, MATSEQSBAIJ, &isSBAIJ));
 22:   if (!isSBAIJ) PetscCall(PetscObjectTypeCompare((PetscObject)mat, MATMPISBAIJ, &isSBAIJ));
 23:   PetscCheck(!isSBAIJ, PetscObjectComm((PetscObject)mat), PETSC_ERR_SUP, "Cannot convert from SBAIJ matrix since cannot obtain entire rows of matrix");

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

 37:     PetscCall(MatSetOption(M, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE));
 38:     PetscCall(MatSetOption(M, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
 39:     PetscCall(PetscObjectTypeCompare((PetscObject)M, MATSEQSBAIJ, &isSBAIJ));
 40:     if (!isSBAIJ) PetscCall(PetscObjectTypeCompare((PetscObject)M, MATMPISBAIJ, &isSBAIJ));
 41:     if (isSBAIJ) PetscCall(MatSetOption(M, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE));
 42:   }

 44:   PetscCall(MatGetOwnershipRange(mat, &rstart, &rend));
 45:   for (i = rstart; i < rend; i++) {
 46:     PetscCall(MatGetRow(mat, i, &nz, &cwork, &vwork));
 47:     PetscCall(MatSetValues(M, 1, &i, nz, cwork, vwork, INSERT_VALUES));
 48:     PetscCall(MatRestoreRow(mat, i, &nz, &cwork, &vwork));
 49:   }
 50:   PetscCall(MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY));
 51:   PetscCall(MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY));

 53:   if (reuse == MAT_INPLACE_MATRIX) {
 54:     PetscCall(MatHeaderReplace(mat, &M));
 55:   } else {
 56:     *newmat = M;
 57:   }
 58:   PetscFunctionReturn(PETSC_SUCCESS);
 59: }