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; 171a0cae60SJed Brown PetscInt i,j,nz,m,n,rstart,rend,lm,ln,prbs,pcbs,cstart,cend,*dnz,*onz; 18c1ac3661SBarry Smith const PetscInt *cwork; 19*d3bf5c86SHong Zhang PetscBool isseqsbaij,ismpisbaij,isseqbaij,ismpibaij,isseqdense,ismpidense; 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); 291a0cae60SJed Brown ierr = MatSetBlockSizes(M,mat->rmap->bs,mat->cmap->bs);CHKERRQ(ierr); 30676c34cdSKris Buschelman ierr = MatSetType(M,newtype);CHKERRQ(ierr); 31*d3bf5c86SHong Zhang ierr = MatGetOwnershipRange(mat,&rstart,&rend);CHKERRQ(ierr); 32*d3bf5c86SHong Zhang 33251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)M,MATSEQSBAIJ,&isseqsbaij);CHKERRQ(ierr); 34251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)M,MATMPISBAIJ,&ismpisbaij);CHKERRQ(ierr); 35d1bf975fSJed Brown if (isseqsbaij || ismpisbaij) {ierr = MatSetOption(M,MAT_IGNORE_LOWER_TRIANGULAR,PETSC_TRUE);CHKERRQ(ierr);} 361a0cae60SJed Brown ierr = PetscObjectTypeCompare((PetscObject)M,MATSEQBAIJ,&isseqbaij);CHKERRQ(ierr); 371a0cae60SJed Brown ierr = PetscObjectTypeCompare((PetscObject)M,MATMPIBAIJ,&ismpibaij);CHKERRQ(ierr); 38*d3bf5c86SHong Zhang ierr = PetscObjectTypeCompare((PetscObject)M,MATSEQDENSE,&isseqdense);CHKERRQ(ierr); 39*d3bf5c86SHong Zhang ierr = PetscObjectTypeCompare((PetscObject)M,MATMPIDENSE,&ismpidense);CHKERRQ(ierr); 40273d9f13SBarry Smith 41*d3bf5c86SHong Zhang if (isseqdense) { 42*d3bf5c86SHong Zhang ierr = MatSeqDenseSetPreallocation(M,NULL);CHKERRQ(ierr); 43*d3bf5c86SHong Zhang } else if (ismpidense) { 44*d3bf5c86SHong Zhang ierr = MatMPIDenseSetPreallocation(M,NULL);CHKERRQ(ierr); 45*d3bf5c86SHong Zhang } else { 461a0cae60SJed Brown /* Preallocation block sizes. (S)BAIJ matrices will have one index per block. */ 471a0cae60SJed Brown prbs = (isseqbaij || ismpibaij || isseqsbaij || ismpisbaij) ? M->rmap->bs : 1; 481a0cae60SJed Brown pcbs = (isseqbaij || ismpibaij || isseqsbaij || ismpisbaij) ? M->cmap->bs : 1; 491a0cae60SJed Brown 501a0cae60SJed Brown ierr = PetscMalloc2(lm/prbs,PetscInt,&dnz,lm/prbs,PetscInt,&onz);CHKERRQ(ierr); 511a0cae60SJed Brown ierr = MatGetOwnershipRangeColumn(mat,&cstart,&cend);CHKERRQ(ierr); 521a0cae60SJed Brown for (i=0; i<lm; i+=prbs) { 531a0cae60SJed Brown ierr = MatGetRow(mat,rstart+i,&nz,&cwork,NULL);CHKERRQ(ierr); 541a0cae60SJed Brown dnz[i] = 0; 551a0cae60SJed Brown onz[i] = 0; 561a0cae60SJed Brown for (j=0; j<nz; j+=pcbs) { 571a0cae60SJed Brown if ((isseqsbaij || ismpisbaij) && cwork[j] < rstart+i) continue; 581a0cae60SJed Brown if (cstart <= cwork[j] && cwork[j] < cend) dnz[i]++; 591a0cae60SJed Brown else onz[i]++; 601a0cae60SJed Brown } 611a0cae60SJed Brown ierr = MatRestoreRow(mat,rstart+i,&nz,&cwork,NULL);CHKERRQ(ierr); 621a0cae60SJed Brown } 631a0cae60SJed Brown ierr = MatXAIJSetPreallocation(M,M->rmap->bs,dnz,onz,dnz,onz);CHKERRQ(ierr); 641a0cae60SJed Brown ierr = PetscFree2(dnz,onz);CHKERRQ(ierr); 65*d3bf5c86SHong Zhang } 661a0cae60SJed Brown 678b6375c0SLois Curfman McInnes for (i=rstart; i<rend; i++) { 688b6375c0SLois Curfman McInnes ierr = MatGetRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 69676c34cdSKris Buschelman ierr = MatSetValues(M,1,&i,nz,cwork,vwork,INSERT_VALUES);CHKERRQ(ierr); 708b6375c0SLois Curfman McInnes ierr = MatRestoreRow(mat,i,&nz,&cwork,&vwork);CHKERRQ(ierr); 718b6375c0SLois Curfman McInnes } 72676c34cdSKris Buschelman ierr = MatAssemblyBegin(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 73676c34cdSKris Buschelman ierr = MatAssemblyEnd(M,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 74676c34cdSKris Buschelman 75ceb03754SKris Buschelman if (reuse == MAT_REUSE_MATRIX) { 768ab5b326SKris Buschelman ierr = MatHeaderReplace(mat,M);CHKERRQ(ierr); 77c3d102feSKris Buschelman } else { 7818f87311SHong Zhang *newmat = M; 79c3d102feSKris Buschelman } 803a40ed3dSBarry Smith PetscFunctionReturn(0); 818b6375c0SLois Curfman McInnes } 82