Skip to content

Commit

Permalink
Mat & Vec: make block sizes lazy
Browse files Browse the repository at this point in the history 
In the past, MatSetUp eagerly set the block sizes (to 1) if not already
set.  This prevents use of the "convenience" constructors because they
call MatSetUp before returning.  We could add bs as an explicit argument
to those constructors, but block size is an extraneous concept for most
users.  The default value of the block size is -1, so now we use
PetscAbs(map->bs) anywhere it is not known in advance that bs has been
explicitly set.  This allows block sizes to be lazy, and they can be set
arbitrarily late, so long as the new size is compatible with the
existing distribution.  For (S)BAIJ matrices, the block size is
structural so it gets set explicitly and can thus be referenced directly
rather than via MatGetBlockSize or wrapped in PetscAbs.
  • Loading branch information
@jedbrown
jedbrown committed Feb 15, 2014
1 parent 69bbac9 commit 33d5767
Show file tree
Hide file tree
Showing 38 changed files with 174 additions and 156 deletions.
5 changes: 4 additions & 1 deletion include/petscis.h
View file
Expand Up @@ -238,7 +238,10 @@ struct _n_PetscLayout{
PetscInt n,N; /* local, global vector size */
PetscInt rstart,rend; /* local start, local end + 1 */
PetscInt *range; /* the offset of each processor */
PetscInt bs; /* number of elements in each block (generally for multi-component problems) Do NOT multiply above numbers by bs */
PetscInt bs; /* number of elements in each block (generally for multi-component
* problems). Defaults to -1 and can be arbitrarily lazy so always use
* PetscAbs(map->bs) when accessing directly and expecting result to be
* positive. Do NOT multiply above numbers by bs */
PetscInt refcnt; /* MPI Vecs obtained with VecDuplicate() and from MatGetVecs() reuse map of input object */
ISLocalToGlobalMapping mapping; /* mapping used in Vec/MatSetValuesLocal() */
ISLocalToGlobalMapping bmapping; /* mapping used in Vec/MatSetValuesBlockedLocal() */
Expand Down
1 change: 1 addition & 0 deletions include/petscmat.h
View file
Expand Up @@ -362,6 +362,7 @@ PETSC_EXTERN PetscErrorCode MatGetBlockSize(Mat,PetscInt *);
PETSC_EXTERN PetscErrorCode MatSetBlockSize(Mat,PetscInt);
PETSC_EXTERN PetscErrorCode MatGetBlockSizes(Mat,PetscInt *,PetscInt *);
PETSC_EXTERN PetscErrorCode MatSetBlockSizes(Mat,PetscInt,PetscInt);
PETSC_EXTERN PetscErrorCode MatSetBlockSizesFromMats(Mat,Mat,Mat);
PETSC_EXTERN PetscErrorCode MatSetNThreads(Mat,PetscInt);
PETSC_EXTERN PetscErrorCode MatGetNThreads(Mat,PetscInt*);

Expand Down
6 changes: 3 additions & 3 deletions src/ksp/pc/impls/pbjacobi/pbjacobi.c
View file
Expand Up @@ -213,9 +213,9 @@ static PetscErrorCode PCSetUp_PBJacobi(PC pc)
PetscFunctionBegin;
if (A->rmap->n != A->cmap->n) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Supported only for square matrices and square storage");

ierr = MatInvertBlockDiagonal(A,&jac->diag);CHKERRQ(ierr);
jac->bs = A->rmap->bs;
jac->mbs = A->rmap->n/A->rmap->bs;
ierr = MatInvertBlockDiagonal(A,&jac->diag);CHKERRQ(ierr);
ierr = MatGetBlockSize(A,&jac->bs);CHKERRQ(ierr);
jac->mbs = A->rmap->n/jac->bs;
switch (jac->bs) {
case 1:
pc->ops->apply = PCApply_PBJacobi_1;
Expand Down
2 changes: 1 addition & 1 deletion src/mat/impls/aij/mpi/mmaij.c
View file
Expand Up @@ -174,7 +174,7 @@ PetscErrorCode MatDisAssemble_MPIAIJ(Mat A)
}
ierr = MatCreate(PETSC_COMM_SELF,&Bnew);CHKERRQ(ierr);
ierr = MatSetSizes(Bnew,m,n,m,n);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Bnew,A->rmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(Bnew,A,A);CHKERRQ(ierr);
ierr = MatSetType(Bnew,((PetscObject)B)->type_name);CHKERRQ(ierr);
ierr = MatSeqAIJSetPreallocation(Bnew,0,nz);CHKERRQ(ierr);

Expand Down
4 changes: 1 addition & 3 deletions src/mat/impls/aij/mpi/mpb_aij.c
View file
Expand Up @@ -20,11 +20,9 @@ PetscErrorCode MatGetMultiProcBlock_MPIAIJ(Mat mat, MPI_Comm subComm, MatReuse
ierr = MatCreate(subComm,subMat);CHKERRQ(ierr);
ierr = MatSetType(*subMat,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatSetSizes(*subMat,mat->rmap->n,mat->cmap->n,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(*subMat,mat->rmap->bs,mat->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(*subMat,mat,mat);CHKERRQ(ierr);

/* need to setup rmap and cmap before Preallocation */
ierr = PetscLayoutSetBlockSize((*subMat)->rmap,mat->rmap->bs);CHKERRQ(ierr);
ierr = PetscLayoutSetBlockSize((*subMat)->cmap,mat->cmap->bs);CHKERRQ(ierr);
ierr = PetscLayoutSetUp((*subMat)->rmap);CHKERRQ(ierr);
ierr = PetscLayoutSetUp((*subMat)->cmap);CHKERRQ(ierr);
}
Expand Down
27 changes: 11 additions & 16 deletions src/mat/impls/aij/mpi/mpiaij.c
View file
Expand Up @@ -202,10 +202,7 @@ PETSC_EXTERN PetscErrorCode MatDistribute_MPIAIJ(MPI_Comm comm,Mat gmat,PetscInt
if (reuse == MAT_INITIAL_MATRIX) {
ierr = MatCreate(comm,&mat);CHKERRQ(ierr);
ierr = MatSetSizes(mat,m,m,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
if (!rank) {
bses[0] = gmat->rmap->bs;
bses[1] = gmat->cmap->bs;
}
ierr = MatGetBlockSizes(gmat,&bses[0],&bses[1]);CHKERRQ(ierr);
ierr = MPI_Bcast(bses,2,MPIU_INT,0,comm);CHKERRQ(ierr);
ierr = MatSetBlockSizes(mat,bses[0],bses[1]);CHKERRQ(ierr);
ierr = MatSetType(mat,MATAIJ);CHKERRQ(ierr);
Expand Down Expand Up @@ -1265,7 +1262,7 @@ PetscErrorCode MatView_MPIAIJ_Binary(Mat mat,PetscViewer viewer)
ierr = PetscFree(column_values);CHKERRQ(ierr);

ierr = PetscViewerBinaryGetInfoPointer(viewer,&file);CHKERRQ(ierr);
if (file) fprintf(file,"-matload_block_size %d\n",(int)mat->rmap->bs);
if (file) fprintf(file,"-matload_block_size %d\n",(int)PetscAbs(mat->rmap->bs));
PetscFunctionReturn(0);
}

Expand Down Expand Up @@ -1943,7 +1940,7 @@ PetscErrorCode MatTranspose_MPIAIJ(Mat A,MatReuse reuse,Mat *matout)

ierr = MatCreate(PetscObjectComm((PetscObject)A),&B);CHKERRQ(ierr);
ierr = MatSetSizes(B,A->cmap->n,A->rmap->n,N,M);CHKERRQ(ierr);
ierr = MatSetBlockSizes(B,A->cmap->bs,A->rmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizes(B,PetscAbs(A->cmap->bs),PetscAbs(A->rmap->bs));CHKERRQ(ierr);
ierr = MatSetType(B,((PetscObject)A)->type_name);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(B,0,d_nnz,0,o_nnz);CHKERRQ(ierr);
ierr = PetscFree4(d_nnz,o_nnz,g_nnz,oloc);CHKERRQ(ierr);
Expand Down Expand Up @@ -2160,7 +2157,7 @@ PetscErrorCode MatAXPY_MPIAIJ(Mat Y,PetscScalar a,Mat X,MatStructure str)
ierr = MatCreate(PetscObjectComm((PetscObject)Y),&B);CHKERRQ(ierr);
ierr = PetscObjectSetName((PetscObject)B,((PetscObject)Y)->name);CHKERRQ(ierr);
ierr = MatSetSizes(B,Y->rmap->n,Y->cmap->n,Y->rmap->N,Y->cmap->N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(B,Y->rmap->bs,Y->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(B,Y,Y);CHKERRQ(ierr);
ierr = MatSetType(B,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatAXPYGetPreallocation_SeqAIJ(yy->A,xx->A,nnz_d);CHKERRQ(ierr);
ierr = MatAXPYGetPreallocation_MPIAIJ(yy->B,yy->garray,xx->B,xx->garray,nnz_o);CHKERRQ(ierr);
Expand Down Expand Up @@ -2266,7 +2263,7 @@ PetscErrorCode MatILUFactorSymbolic_MPIAIJ(Mat fact,Mat A, IS isrow, IS iscol, c
ierr = MatGetLocalSize(A, &m, &n);CHKERRQ(ierr);
ierr = MatGetSize(A, &M, &N);CHKERRQ(ierr);
ierr = MatSetSizes(fact, m, n, M, N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(fact,A->rmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(fact,A,A);CHKERRQ(ierr);
ierr = MatSetType(fact, ((PetscObject)A)->type_name);CHKERRQ(ierr);
ierr = MatAssemblyBegin(fact, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(fact, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
Expand Down Expand Up @@ -2666,7 +2663,7 @@ PetscErrorCode MatGetRedundantMatrix_MPIAIJ_interlaced(Mat mat,PetscInt nsubcomm
} else { /* non-square matrix */
ierr = MatSetSizes(C,mloc_sub,PETSC_DECIDE,PETSC_DECIDE,mat->cmap->N);CHKERRQ(ierr);
}
ierr = MatSetBlockSizes(C,mat->rmap->bs,mat->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(C,mat,mat);CHKERRQ(ierr);
ierr = MatSetFromOptions(C);CHKERRQ(ierr);
ierr = MatSeqAIJSetPreallocation(C,rownz_max,NULL);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(C,rownz_max,NULL,rownz_max,NULL);CHKERRQ(ierr);
Expand Down Expand Up @@ -3242,12 +3239,12 @@ PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJ(Mat B,PetscInt d_nz,const Petsc
/* Explicitly create 2 MATSEQAIJ matrices. */
ierr = MatCreate(PETSC_COMM_SELF,&b->A);CHKERRQ(ierr);
ierr = MatSetSizes(b->A,B->rmap->n,B->cmap->n,B->rmap->n,B->cmap->n);CHKERRQ(ierr);
ierr = MatSetBlockSizes(b->A,B->rmap->bs,B->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(b->A,B,B);CHKERRQ(ierr);
ierr = MatSetType(b->A,MATSEQAIJ);CHKERRQ(ierr);
ierr = PetscLogObjectParent((PetscObject)B,(PetscObject)b->A);CHKERRQ(ierr);
ierr = MatCreate(PETSC_COMM_SELF,&b->B);CHKERRQ(ierr);
ierr = MatSetSizes(b->B,B->rmap->n,B->cmap->N,B->rmap->n,B->cmap->N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(b->B,B->rmap->bs,B->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(b->B,B,B);CHKERRQ(ierr);
ierr = MatSetType(b->B,MATSEQAIJ);CHKERRQ(ierr);
ierr = PetscLogObjectParent((PetscObject)B,(PetscObject)b->B);CHKERRQ(ierr);
}
Expand All @@ -3270,14 +3267,12 @@ PetscErrorCode MatDuplicate_MPIAIJ(Mat matin,MatDuplicateOption cpvalues,Mat *ne
*newmat = 0;
ierr = MatCreate(PetscObjectComm((PetscObject)matin),&mat);CHKERRQ(ierr);
ierr = MatSetSizes(mat,matin->rmap->n,matin->cmap->n,matin->rmap->N,matin->cmap->N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(mat,matin->rmap->bs,matin->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(mat,matin,matin);CHKERRQ(ierr);
ierr = MatSetType(mat,((PetscObject)matin)->type_name);CHKERRQ(ierr);
ierr = PetscMemcpy(mat->ops,matin->ops,sizeof(struct _MatOps));CHKERRQ(ierr);
a = (Mat_MPIAIJ*)mat->data;

mat->factortype = matin->factortype;
mat->rmap->bs = matin->rmap->bs;
mat->cmap->bs = matin->cmap->bs;
mat->assembled = PETSC_TRUE;
mat->insertmode = NOT_SET_VALUES;
mat->preallocated = PETSC_TRUE;
Expand Down Expand Up @@ -4440,7 +4435,7 @@ PetscErrorCode MatFileSplit(Mat A,char *outfile)
/* Should this be the type of the diagonal block of A? */
ierr = MatCreate(PETSC_COMM_SELF,&B);CHKERRQ(ierr);
ierr = MatSetSizes(B,m,N,m,N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(B,A->rmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(B,A,A);CHKERRQ(ierr);
ierr = MatSetType(B,MATSEQAIJ);CHKERRQ(ierr);
ierr = MatSeqAIJSetPreallocation(B,0,NULL);CHKERRQ(ierr);
ierr = MatGetOwnershipRange(A,&rstart,0);CHKERRQ(ierr);
Expand Down Expand Up @@ -5486,7 +5481,7 @@ PetscErrorCode MatMatMultSymbolic_MPIDense_MPIAIJ(Mat A,Mat B,PetscReal fill,Mat
if (A->cmap->n != B->rmap->n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"A->cmap->n %d != B->rmap->n %d\n",A->cmap->n,B->rmap->n);
ierr = MatCreate(PetscObjectComm((PetscObject)A),&Cmat);CHKERRQ(ierr);
ierr = MatSetSizes(Cmat,m,n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmat,A->rmap->bs,B->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(Cmat,A,B);CHKERRQ(ierr);
ierr = MatSetType(Cmat,MATMPIDENSE);CHKERRQ(ierr);
ierr = MatMPIDenseSetPreallocation(Cmat,NULL);CHKERRQ(ierr);
ierr = MatAssemblyBegin(Cmat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
Expand Down
6 changes: 4 additions & 2 deletions src/mat/impls/aij/mpi/mpicusp/mpiaijcusp.cu
View file
Expand Up @@ -52,19 +52,21 @@ PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJCUSP(Mat B,PetscInt d_nz,const P
PetscErrorCode MatGetVecs_MPIAIJCUSP(Mat mat,Vec *right,Vec *left)
{
PetscErrorCode ierr;
PetscInt rbs,cbs;

PetscFunctionBegin;
ierr = MatGetBlockSizes(mat,&rbs,&cbs);CHKERRQ(ierr);
if (right) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),right);CHKERRQ(ierr);
ierr = VecSetSizes(*right,mat->cmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,cbs);CHKERRQ(ierr);
ierr = VecSetType(*right,VECCUSP);CHKERRQ(ierr);
ierr = VecSetLayout(*right,mat->cmap);CHKERRQ(ierr);
}
if (left) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),left);CHKERRQ(ierr);
ierr = VecSetSizes(*left,mat->rmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,rbs);CHKERRQ(ierr);
ierr = VecSetType(*left,VECCUSP);CHKERRQ(ierr);
ierr = VecSetLayout(*left,mat->rmap);CHKERRQ(ierr);
}
Expand Down
6 changes: 4 additions & 2 deletions src/mat/impls/aij/mpi/mpicusparse/mpiaijcusparse.cu
View file
Expand Up @@ -61,19 +61,21 @@ PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJCUSPARSE(Mat B,PetscInt d_nz,con
PetscErrorCode MatGetVecs_MPIAIJCUSPARSE(Mat mat,Vec *right,Vec *left)
{
PetscErrorCode ierr;
PetscInt rbs,cbs;

PetscFunctionBegin;
ierr = MatGetBlockSizes(mat,&rbs,&cbs);CHKERRQ(ierr);
if (right) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),right);CHKERRQ(ierr);
ierr = VecSetSizes(*right,mat->cmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,cbs);CHKERRQ(ierr);
ierr = VecSetType(*right,VECCUSP);CHKERRQ(ierr);
ierr = VecSetLayout(*right,mat->cmap);CHKERRQ(ierr);
}
if (left) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),left);CHKERRQ(ierr);
ierr = VecSetSizes(*left,mat->rmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,rbs);CHKERRQ(ierr);
ierr = VecSetType(*left,VECCUSP);CHKERRQ(ierr);
ierr = VecSetLayout(*left,mat->rmap);CHKERRQ(ierr);

Expand Down
10 changes: 5 additions & 5 deletions src/mat/impls/aij/mpi/mpimatmatmult.c
View file
Expand Up @@ -303,7 +303,7 @@ PetscErrorCode MatMatMultSymbolic_MPIAIJ_MPIAIJ_nonscalable(Mat A,Mat P,PetscRea
/*----------------------------------------------------*/
ierr = MatCreate(comm,&Cmpi);CHKERRQ(ierr);
ierr = MatSetSizes(Cmpi,am,pn,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,A->rmap->bs,P->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(Cmpi,A,P);CHKERRQ(ierr);

ierr = MatSetType(Cmpi,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Cmpi,0,dnz,0,onz);CHKERRQ(ierr);
Expand Down Expand Up @@ -401,7 +401,7 @@ PetscErrorCode MatMatMultSymbolic_MPIAIJ_MPIDense(Mat A,Mat B,PetscReal fill,Mat
PetscFunctionBegin;
ierr = MatCreate(PetscObjectComm((PetscObject)B),C);CHKERRQ(ierr);
ierr = MatSetSizes(*C,m,n,A->rmap->N,B->cmap->N);CHKERRQ(ierr);
ierr = MatSetBlockSizes(*C,A->rmap->bs,B->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(*C,A,B);CHKERRQ(ierr);
ierr = MatSetType(*C,MATMPIDENSE);CHKERRQ(ierr);
ierr = MatMPIDenseSetPreallocation(*C,NULL);CHKERRQ(ierr);
ierr = MatAssemblyBegin(*C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
Expand Down Expand Up @@ -752,7 +752,7 @@ PetscErrorCode MatMatMultSymbolic_MPIAIJ_MPIAIJ(Mat A,Mat P,PetscReal fill,Mat *
/*----------------------------------------------------*/
ierr = MatCreate(comm,&Cmpi);CHKERRQ(ierr);
ierr = MatSetSizes(Cmpi,am,pn,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,A->rmap->bs,P->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(Cmpi,A,P);CHKERRQ(ierr);
ierr = MatSetType(Cmpi,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Cmpi,0,dnz,0,onz);CHKERRQ(ierr);
ierr = MatPreallocateFinalize(dnz,onz);CHKERRQ(ierr);
Expand Down Expand Up @@ -1328,7 +1328,7 @@ PetscErrorCode MatTransposeMatMultSymbolic_MPIAIJ_MPIAIJ_nonscalable(Mat P,Mat A

ierr = MatCreate(comm,&Cmpi);CHKERRQ(ierr);
ierr = MatSetSizes(Cmpi,pn,A->cmap->n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,P->cmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,PetscAbs(P->cmap->bs),PetscAbs(A->cmap->bs));CHKERRQ(ierr);
ierr = MatSetType(Cmpi,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Cmpi,0,dnz,0,onz);CHKERRQ(ierr);
ierr = MatPreallocateFinalize(dnz,onz);CHKERRQ(ierr);
Expand Down Expand Up @@ -1834,7 +1834,7 @@ PetscErrorCode MatTransposeMatMultSymbolic_MPIAIJ_MPIAIJ(Mat P,Mat A,PetscReal f

ierr = MatCreate(comm,&Cmpi);CHKERRQ(ierr);
ierr = MatSetSizes(Cmpi,pn,A->cmap->n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,P->cmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,PetscAbs(P->cmap->bs),PetscAbs(A->cmap->bs));CHKERRQ(ierr);
ierr = MatSetType(Cmpi,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Cmpi,0,dnz,0,onz);CHKERRQ(ierr);
ierr = MatPreallocateFinalize(dnz,onz);CHKERRQ(ierr);
Expand Down
6 changes: 3 additions & 3 deletions src/mat/impls/aij/mpi/mpiov.c
View file
Expand Up @@ -566,7 +566,7 @@ PetscErrorCode MatGetSubMatrix_MPIAIJ_All(Mat A,MatGetSubMatrixOption flag,MatRe
*/
ierr = MatCreate(PETSC_COMM_SELF,&B);CHKERRQ(ierr);
ierr = MatSetSizes(B,A->rmap->N,A->cmap->N,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(B,A->rmap->bs,A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(B,A,A);CHKERRQ(ierr);
ierr = MatSetType(B,((PetscObject)a->A)->type_name);CHKERRQ(ierr);
ierr = MatSeqAIJSetPreallocation(B,0,lens);CHKERRQ(ierr);
ierr = PetscMalloc(sizeof(Mat),Bin);CHKERRQ(ierr);
Expand Down Expand Up @@ -1542,11 +1542,11 @@ PetscErrorCode MatCreateMPIAIJFromSeqMatrices_Private(MPI_Comm comm, Mat A, Mat
if (!seqaij) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Diagonal matrix is of wrong type");
ierr = PetscObjectTypeCompare((PetscObject)B, MATSEQAIJ, &seqaij);CHKERRQ(ierr);
if (!seqaij) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Off-diagonal matrix is of wrong type");
if (A->rmap->n != B->rmap->n || A->rmap->bs != B->rmap->bs || A->cmap->bs != B->cmap->bs) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incompatible component matrices of an MPIAIJ matrix");
if (PetscAbs(A->rmap->n) != PetscAbs(B->rmap->n) || PetscAbs(A->rmap->bs) != PetscAbs(B->rmap->bs) || PetscAbs(A->cmap->bs) != PetscAbs(B->cmap->bs)) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incompatible component matrices of an MPIAIJ matrix");

ierr = MatCreate(comm, C);CHKERRQ(ierr);
ierr = MatSetSizes(*C,A->rmap->n, A->cmap->n, PETSC_DECIDE, PETSC_DECIDE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(*C,A->rmap->bs, A->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizesFromMats(*C,A,A);CHKERRQ(ierr);
ierr = PetscLayoutSetUp((*C)->rmap);CHKERRQ(ierr);
ierr = PetscLayoutSetUp((*C)->cmap);CHKERRQ(ierr);
if ((*C)->cmap->N != A->cmap->n + B->cmap->n) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Incompatible component matrices of an MPIAIJ matrix");
Expand Down
2 changes: 1 addition & 1 deletion src/mat/impls/aij/mpi/mpiptap.c
View file
Expand Up @@ -448,7 +448,7 @@ PetscErrorCode MatPtAPSymbolic_MPIAIJ_MPIAIJ(Mat A,Mat P,PetscReal fill,Mat *C)
/*--------------------------------------*/
ierr = MatCreate(comm,&Cmpi);CHKERRQ(ierr);
ierr = MatSetSizes(Cmpi,pn,pn,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,P->cmap->bs,P->cmap->bs);CHKERRQ(ierr);
ierr = MatSetBlockSizes(Cmpi,PetscAbs(P->cmap->bs),PetscAbs(P->cmap->bs));CHKERRQ(ierr);
ierr = MatSetType(Cmpi,MATMPIAIJ);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Cmpi,0,dnz,0,onz);CHKERRQ(ierr);
ierr = MatPreallocateFinalize(dnz,onz);CHKERRQ(ierr);
Expand Down
6 changes: 4 additions & 2 deletions src/mat/impls/aij/mpi/mpiviennacl/mpiaijviennacl.cxx
View file
Expand Up @@ -34,19 +34,21 @@ PetscErrorCode MatMPIAIJSetPreallocation_MPIAIJViennaCL(Mat B,PetscInt d_nz,con
PetscErrorCode MatGetVecs_MPIAIJViennaCL(Mat mat,Vec *right,Vec *left)
{
PetscErrorCode ierr;
PetscInt rbs,cbs;

PetscFunctionBegin;
ierr = MatGetBlockSizes(mat,&rbs,&cbs);CHKERRQ(ierr);
if (right) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),right);CHKERRQ(ierr);
ierr = VecSetSizes(*right,mat->cmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*right,cbs);CHKERRQ(ierr);
ierr = VecSetType(*right,VECVIENNACL);CHKERRQ(ierr);
ierr = VecSetLayout(*right,mat->cmap);CHKERRQ(ierr);
}
if (left) {
ierr = VecCreate(PetscObjectComm((PetscObject)mat),left);CHKERRQ(ierr);
ierr = VecSetSizes(*left,mat->rmap->n,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,mat->rmap->bs);CHKERRQ(ierr);
ierr = VecSetBlockSize(*left,rbs);CHKERRQ(ierr);
ierr = VecSetType(*left,VECVIENNACL);CHKERRQ(ierr);
ierr = VecSetLayout(*left,mat->rmap);CHKERRQ(ierr);
}
Expand Down
11 changes: 7 additions & 4 deletions src/mat/impls/aij/mpi/mumps/mumps.c
View file
Expand Up @@ -139,12 +139,14 @@ PetscErrorCode MatConvertToTriples_seqaij_seqaij(Mat A,int shift,MatReuse reuse,
PetscErrorCode MatConvertToTriples_seqbaij_seqaij(Mat A,int shift,MatReuse reuse,int *nnz,int **r, int **c, PetscScalar **v)
{
Mat_SeqBAIJ *aa=(Mat_SeqBAIJ*)A->data;
const PetscInt *ai,*aj,*ajj,bs=A->rmap->bs,bs2=aa->bs2,M=A->rmap->N/bs;
PetscInt nz,idx=0,rnz,i,j,k,m;
const PetscInt *ai,*aj,*ajj,bs2 = aa->bs2;
PetscInt bs,M,nz,idx=0,rnz,i,j,k,m;
PetscErrorCode ierr;
PetscInt *row,*col;

PetscFunctionBegin;
ierr = MatGetBlockSize(A,&bs);CHKERRQ(ierr);
M = A->rmap->N/bs;
*v = aa->a;
if (reuse == MAT_INITIAL_MATRIX) {
ai = aa->i; aj = aa->j;
Expand Down Expand Up @@ -381,14 +383,15 @@ PetscErrorCode MatConvertToTriples_mpibaij_mpiaij(Mat A,int shift,MatReuse reuse
Mat_SeqBAIJ *bb = (Mat_SeqBAIJ*)(mat->B)->data;
const PetscInt *ai = aa->i, *bi = bb->i, *aj = aa->j, *bj = bb->j,*ajj, *bjj;
const PetscInt *garray = mat->garray,mbs=mat->mbs,rstart=A->rmap->rstart;
const PetscInt bs = A->rmap->bs,bs2=mat->bs2;
const PetscInt bs2=mat->bs2;
PetscErrorCode ierr;
PetscInt nz,i,j,k,n,jj,irow,countA,countB,idx;
PetscInt bs,nz,i,j,k,n,jj,irow,countA,countB,idx;
PetscInt *row,*col;
const PetscScalar *av=aa->a, *bv=bb->a,*v1,*v2;
PetscScalar *val;

PetscFunctionBegin;
ierr = MatGetBlockSize(A,&bs);CHKERRQ(ierr);
if (reuse == MAT_INITIAL_MATRIX) {
nz = bs2*(aa->nz + bb->nz);
*nnz = nz;
Expand Down

0 comments on commit 33d5767

Please sign in to comment.