173f4d377SMatthew Knepley /*$Id: convert.c,v 1.76 2001/08/07 03:03:20 balay Exp $*/ 256fe5c5cSLois Curfman McInnes 394a9d846SBarry Smith #include "src/mat/matimpl.h" 48b6375c0SLois Curfman McInnes 54a2ae208SSatish Balay #undef __FUNCT__ 64a2ae208SSatish Balay #define __FUNCT__ "MatConvert_Basic" 78b6375c0SLois Curfman McInnes /* 88b6375c0SLois Curfman McInnes MatConvert_Basic - Converts from any input format to another format. For 98b6375c0SLois Curfman McInnes parallel formats, the new matrix distribution is determined by PETSc. 10273d9f13SBarry Smith 11273d9f13SBarry Smith Does not do preallocation so in general will be slow 128b6375c0SLois Curfman McInnes */ 13676c34cdSKris Buschelman int MatConvert_Basic(Mat mat,MatType newtype,Mat *newmat) { 14676c34cdSKris Buschelman Mat M; 1587828ca2SBarry Smith PetscScalar *vwork; 16435da068SBarry Smith int ierr,i,nz,m,n,*cwork,rstart,rend,lm,ln; 1725cdf11fSBarry Smith 183a40ed3dSBarry Smith PetscFunctionBegin; 198b6375c0SLois Curfman McInnes ierr = MatGetSize(mat,&m,&n);CHKERRQ(ierr); 20435da068SBarry Smith ierr = MatGetLocalSize(mat,&lm,&ln);CHKERRQ(ierr); 21273d9f13SBarry Smith 22435da068SBarry Smith if (ln == n) ln = PETSC_DECIDE; /* try to preserve column ownership */ 23435da068SBarry Smith 24676c34cdSKris Buschelman ierr = MatCreate(mat->comm,lm,ln,m,n,&M);CHKERRQ(ierr); 25676c34cdSKris Buschelman ierr = MatSetType(M,newtype);CHKERRQ(ierr); 26273d9f13SBarry Smith 27f1af5d2fSBarry Smith ierr = MatGetOwnershipRange(mat,&rstart,&rend);CHKERRQ(ierr); 288b6375c0SLois Curfman McInnes for (i=rstart; i<rend; i++) { 298b6375c0SLois Curfman McInnes ierr = MatGetRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 30676c34cdSKris Buschelman ierr = MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);CHKERRQ(ierr); 318b6375c0SLois Curfman McInnes ierr = MatRestoreRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 328b6375c0SLois Curfman McInnes } 33676c34cdSKris Buschelman ierr = MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 34676c34cdSKris Buschelman ierr = MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 35676c34cdSKris Buschelman 36676c34cdSKris Buschelman /* Fake support for "inplace" convert. */ 37*efe6450aSKris Buschelman if (*newmat == mat) { 38*efe6450aSKris Buschelman ierr = MatDestroy(mat);CHKERRQ(ierr); 39676c34cdSKris Buschelman } 40676c34cdSKris Buschelman *newmat = M; 413a40ed3dSBarry Smith PetscFunctionReturn(0); 428b6375c0SLois Curfman McInnes } 43