156fe5c5cSLois Curfman McInnes 2b45d2f2cSJed Brown #include <petsc-private/matimpl.h> 38b6375c0SLois Curfman McInnes 44a2ae208SSatish Balay #undef __FUNCT__ 54a2ae208SSatish Balay #define __FUNCT__ "MatConvert_Basic" 68b6375c0SLois Curfman McInnes /* 78b6375c0SLois Curfman McInnes MatConvert_Basic - Converts from any input format to another format. For 88b6375c0SLois Curfman McInnes parallel formats, the new matrix distribution is determined by PETSc. 9273d9f13SBarry Smith 10273d9f13SBarry Smith Does not do preallocation so in general will be slow 118b6375c0SLois Curfman McInnes */ 1219fd82e9SBarry Smith PetscErrorCode MatConvert_Basic(Mat mat, MatType newtype,MatReuse reuse,Mat *newmat) 13b3cc6726SBarry Smith { 14676c34cdSKris Buschelman Mat M; 15b3cc6726SBarry Smith const PetscScalar *vwork; 16dfbe8321SBarry Smith PetscErrorCode ierr; 17*1a0cae60SJed Brown PetscInt i,j,nz,m,n,rstart,rend,lm,ln,prbs,pcbs,cstart,cend,*dnz,*onz; 18c1ac3661SBarry Smith const PetscInt *cwork; 19*1a0cae60SJed Brown PetscBool isseqsbaij,ismpisbaij,isseqbaij,ismpibaij; 2025cdf11fSBarry Smith 213a40ed3dSBarry Smith PetscFunctionBegin; 228b6375c0SLois Curfman McInnes ierr = MatGetSize(mat,&m,&n);CHKERRQ(ierr); 23435da068SBarry Smith ierr = MatGetLocalSize(mat,&lm,&ln);CHKERRQ(ierr); 24273d9f13SBarry Smith 25435da068SBarry Smith if (ln == n) ln = PETSC_DECIDE; /* try to preserve column ownership */ 26435da068SBarry Smith 27ce94432eSBarry Smith ierr = MatCreate(PetscObjectComm((PetscObject)mat),&M);CHKERRQ(ierr); 28f69a0ea3SMatthew Knepley ierr = MatSetSizes(M,lm,ln,m,n);CHKERRQ(ierr); 29*1a0cae60SJed Brown ierr = MatSetBlockSizes(M,mat->rmap->bs,mat->cmap->bs);CHKERRQ(ierr); 30676c34cdSKris Buschelman ierr = MatSetType(M,newtype);CHKERRQ(ierr); 31251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)M,MATSEQSBAIJ,&isseqsbaij);CHKERRQ(ierr); 32251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)M,MATMPISBAIJ,&ismpisbaij);CHKERRQ(ierr); 33d1bf975fSJed Brown if (isseqsbaij || ismpisbaij) {ierr = MatSetOption(M,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);CHKERRQ(ierr);} 34*1a0cae60SJed Brown ierr = PetscObjectTypeCompare((PetscObject)M,MATSEQBAIJ,&isseqbaij);CHKERRQ(ierr); 35*1a0cae60SJed Brown ierr = PetscObjectTypeCompare((PetscObject)M,MATMPIBAIJ,&ismpibaij);CHKERRQ(ierr); 36273d9f13SBarry Smith 37*1a0cae60SJed Brown /* Preallocation block sizes. (S)BAIJ matrices will have one index per block. */ 38*1a0cae60SJed Brown prbs = (isseqbaij || ismpibaij || isseqsbaij || ismpisbaij) ? M->rmap->bs : 1; 39*1a0cae60SJed Brown pcbs = (isseqbaij || ismpibaij || isseqsbaij || ismpisbaij) ? M->cmap->bs : 1; 40*1a0cae60SJed Brown 41*1a0cae60SJed Brown ierr = PetscMalloc2(lm/prbs,PetscInt,&dnz,lm/prbs,PetscInt,&onz);CHKERRQ(ierr); 42f1af5d2fSBarry Smith ierr = MatGetOwnershipRange(mat,&rstart,&rend);CHKERRQ(ierr); 43*1a0cae60SJed Brown ierr = MatGetOwnershipRangeColumn(mat,&cstart,&cend);CHKERRQ(ierr); 44*1a0cae60SJed Brown for (i=0; i<lm; i+=prbs) { 45*1a0cae60SJed Brown ierr = MatGetRow(mat,rstart+i,&nz,&cwork,NULL);CHKERRQ(ierr); 46*1a0cae60SJed Brown dnz[i] = 0; 47*1a0cae60SJed Brown onz[i] = 0; 48*1a0cae60SJed Brown for (j=0; j<nz; j+=pcbs) { 49*1a0cae60SJed Brown if ((isseqsbaij || ismpisbaij) && cwork[j] < rstart+i) continue; 50*1a0cae60SJed Brown if (cstart <= cwork[j] && cwork[j] < cend) dnz[i]++; 51*1a0cae60SJed Brown else onz[i]++; 52*1a0cae60SJed Brown } 53*1a0cae60SJed Brown ierr = MatRestoreRow(mat,rstart+i,&nz,&cwork,NULL);CHKERRQ(ierr); 54*1a0cae60SJed Brown } 55*1a0cae60SJed Brown ierr = MatXAIJSetPreallocation(M,M->rmap->bs,dnz,onz,dnz,onz);CHKERRQ(ierr); 56*1a0cae60SJed Brown ierr = PetscFree2(dnz,onz);CHKERRQ(ierr); 57*1a0cae60SJed Brown 588b6375c0SLois Curfman McInnes for (i=rstart; i<rend; i++) { 598b6375c0SLois Curfman McInnes ierr = MatGetRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 60676c34cdSKris Buschelman ierr = MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);CHKERRQ(ierr); 618b6375c0SLois Curfman McInnes ierr = MatRestoreRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 628b6375c0SLois Curfman McInnes } 63676c34cdSKris Buschelman ierr = MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 64676c34cdSKris Buschelman ierr = MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 65676c34cdSKris Buschelman 66ceb03754SKris Buschelman if (reuse == MAT_REUSE_MATRIX) { 678ab5b326SKris Buschelman ierr = MatHeaderReplace(mat,M);CHKERRQ(ierr); 68c3d102feSKris Buschelman } else { 6918f87311SHong Zhang *newmat = M; 70c3d102feSKris Buschelman } 713a40ed3dSBarry Smith PetscFunctionReturn(0); 728b6375c0SLois Curfman McInnes } 73