1 2 #include <petsc/private/matimpl.h> 3 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; 15 16 PetscFunctionBegin; 17 if (!mat->ops->getrow) { /* missing get row, use matvecs */ 18 CHKERRQ(MatConvert_Shell(mat,newtype,reuse,newmat)); 19 PetscFunctionReturn(0); 20 } 21 CHKERRQ(PetscObjectTypeCompare((PetscObject)mat,MATSEQSBAIJ,&isSBAIJ)); 22 if (!isSBAIJ) { 23 CHKERRQ(PetscObjectTypeCompare((PetscObject)mat,MATMPISBAIJ,&isSBAIJ)); 24 } 25 PetscCheck(!isSBAIJ,PetscObjectComm((PetscObject)mat),PETSC_ERR_SUP,"Cannot convert from SBAIJ matrix since cannot obtain entire rows of matrix"); 26 27 if (reuse == MAT_REUSE_MATRIX) { 28 M = *newmat; 29 } else { 30 PetscInt m,n,lm,ln; 31 CHKERRQ(MatGetSize(mat,&m,&n)); 32 CHKERRQ(MatGetLocalSize(mat,&lm,&ln)); 33 CHKERRQ(MatCreate(PetscObjectComm((PetscObject)mat),&M)); 34 CHKERRQ(MatSetSizes(M,lm,ln,m,n)); 35 CHKERRQ(MatSetBlockSizesFromMats(M,mat,mat)); 36 CHKERRQ(MatSetType(M,newtype)); 37 CHKERRQ(MatSetUp(M)); 38 39 CHKERRQ(MatSetOption(M,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE)); 40 CHKERRQ(MatSetOption(M,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE)); 41 CHKERRQ(PetscObjectTypeCompare((PetscObject)M,MATSEQSBAIJ,&isSBAIJ)); 42 if (!isSBAIJ) { 43 CHKERRQ(PetscObjectTypeCompare((PetscObject)M,MATMPISBAIJ,&isSBAIJ)); 44 } 45 if (isSBAIJ) { 46 CHKERRQ(MatSetOption(M,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE)); 47 } 48 } 49 50 CHKERRQ(MatGetOwnershipRange(mat,&rstart,&rend)); 51 for (i=rstart; i<rend; i++) { 52 CHKERRQ(MatGetRow(mat,i,&nz,&cwork,&vwork)); 53 CHKERRQ(MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES)); 54 CHKERRQ(MatRestoreRow(mat,i,&nz,&cwork,&vwork)); 55 } 56 CHKERRQ(MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY)); 57 CHKERRQ(MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY)); 58 59 if (reuse == MAT_INPLACE_MATRIX) { 60 CHKERRQ(MatHeaderReplace(mat,&M)); 61 } else { 62 *newmat = M; 63 } 64 PetscFunctionReturn(0); 65 } 66