1 #include <petsc/private/matimpl.h> 2 3 /* 4 MatConvert_Basic - Converts from any input format to another format. 5 Does not do preallocation so in general will be slow 6 */ 7 PetscErrorCode MatConvert_Basic(Mat mat, MatType newtype, MatReuse reuse, Mat *newmat) 8 { 9 Mat M; 10 const PetscScalar *vwork; 11 PetscInt i, rstart, rend, nz; 12 const PetscInt *cwork; 13 PetscBool isSBAIJ; 14 15 PetscFunctionBegin; 16 if (!mat->ops->getrow) { /* missing get row, use matvecs */ 17 PetscCall(MatConvert_Shell(mat, newtype, reuse, newmat)); 18 PetscFunctionReturn(PETSC_SUCCESS); 19 } 20 PetscCall(PetscObjectTypeCompare((PetscObject)mat, MATSEQSBAIJ, &isSBAIJ)); 21 if (!isSBAIJ) PetscCall(PetscObjectTypeCompare((PetscObject)mat, MATMPISBAIJ, &isSBAIJ)); 22 PetscCheck(!isSBAIJ, PetscObjectComm((PetscObject)mat), PETSC_ERR_SUP, "Cannot convert from SBAIJ matrix since cannot obtain entire rows of matrix"); 23 24 if (reuse == MAT_REUSE_MATRIX) { 25 M = *newmat; 26 } else { 27 PetscInt m, n, lm, ln; 28 PetscCall(MatGetSize(mat, &m, &n)); 29 PetscCall(MatGetLocalSize(mat, &lm, &ln)); 30 PetscCall(MatCreate(PetscObjectComm((PetscObject)mat), &M)); 31 PetscCall(MatSetSizes(M, lm, ln, m, n)); 32 PetscCall(MatSetBlockSizesFromMats(M, mat, mat)); 33 PetscCall(MatSetType(M, newtype)); 34 PetscCall(MatSetUp(M)); 35 36 PetscCall(MatSetOption(M, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE)); 37 PetscCall(MatSetOption(M, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE)); 38 PetscCall(MatSetOption(M, MAT_NO_OFF_PROC_ENTRIES, PETSC_TRUE)); 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 } 43 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)); 52 if (reuse != MAT_REUSE_MATRIX) PetscCall(MatSetOption(M, MAT_NO_OFF_PROC_ENTRIES, PETSC_FALSE)); 53 54 if (reuse == MAT_INPLACE_MATRIX) { 55 PetscCall(MatHeaderReplace(mat, &M)); 56 } else { 57 *newmat = M; 58 } 59 PetscFunctionReturn(PETSC_SUCCESS); 60 } 61