本文整理汇总了C++中MatDuplicate函数的典型用法代码示例。如果您正苦于以下问题:C++ MatDuplicate函数的具体用法?C++ MatDuplicate怎么用?C++ MatDuplicate使用的例子?那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。
在下文中一共展示了MatDuplicate函数的20个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的C++代码示例。
示例1: main
int main(int argc,char **args)
{
Mat A1,A2,A3,A4,nest;
Mat mata[4];
Mat aij;
MPI_Comm comm;
PetscInt m,n,istart,iend,ii,i,J,j;
PetscScalar v;
PetscMPIInt size;
PetscErrorCode ierr;
ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
comm = PETSC_COMM_WORLD;
ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr);
/*
Assemble the matrix for the five point stencil, YET AGAIN
*/
ierr = MatCreate(comm,&A1);CHKERRQ(ierr);
m=2,n=2;
ierr = MatSetSizes(A1,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr);
ierr = MatSetFromOptions(A1);CHKERRQ(ierr);
ierr = MatSetUp(A1);CHKERRQ(ierr);
ierr = MatGetOwnershipRange(A1,&istart,&iend);CHKERRQ(ierr);
for (ii=istart; ii<iend; ii++) {
v = -1.0; i = ii/n; j = ii - i*n;
if (i>0) {J = ii - n; ierr = MatSetValues(A1,1,&ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);}
if (i<m-1) {J = ii + n; ierr = MatSetValues(A1,1,&ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);}
if (j>0) {J = ii - 1; ierr = MatSetValues(A1,1,&ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);}
if (j<n-1) {J = ii + 1; ierr = MatSetValues(A1,1,&ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);}
v = 4.0; ierr = MatSetValues(A1,1,&ii,1,&ii,&v,INSERT_VALUES);CHKERRQ(ierr);
}
ierr = MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatView(A1,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
ierr = MatDuplicate(A1,MAT_COPY_VALUES,&A2);CHKERRQ(ierr);
ierr = MatDuplicate(A1,MAT_COPY_VALUES,&A3);CHKERRQ(ierr);
ierr = MatDuplicate(A1,MAT_COPY_VALUES,&A4);CHKERRQ(ierr);
/*create a nest matrix */
ierr = MatCreate(comm,&nest);CHKERRQ(ierr);
ierr = MatSetType(nest,MATNEST);CHKERRQ(ierr);
mata[0]=A1,mata[1]=A2,mata[2]=A3,mata[3]=A4;
ierr = MatNestSetSubMats(nest,2,NULL,2,NULL,mata);CHKERRQ(ierr);
ierr = MatSetUp(nest);CHKERRQ(ierr);
ierr = MatConvert(nest,MATAIJ,MAT_INITIAL_MATRIX,&aij);CHKERRQ(ierr);
ierr = MatView(aij,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
ierr = MatDestroy(&nest);CHKERRQ(ierr);
ierr = MatDestroy(&aij);CHKERRQ(ierr);
ierr = MatDestroy(&A1);CHKERRQ(ierr);
ierr = MatDestroy(&A2);CHKERRQ(ierr);
ierr = MatDestroy(&A3);CHKERRQ(ierr);
ierr = MatDestroy(&A4);CHKERRQ(ierr);
ierr = PetscFinalize();
return ierr;
}
开发者ID:firedrakeproject,项目名称:petsc,代码行数:54,代码来源:ex195.c
示例2: MatConvert_MPISBAIJ_MPISBSTRM
PETSC_EXTERN PetscErrorCode MatConvert_MPISBAIJ_MPISBSTRM(Mat A,MatType type,MatReuse reuse,Mat *newmat)
{
PetscErrorCode ierr;
Mat B = *newmat;
Mat_SeqSBSTRM *sbstrm;
PetscFunctionBegin;
if (reuse == MAT_INITIAL_MATRIX) {
ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
}
/* printf(" --- in MatConvert_MPISBAIJ_MPISBSTRM -- 1 \n"); */
ierr = PetscNewLog(B, Mat_SeqSBSTRM,&sbstrm);CHKERRQ(ierr);
B->spptr = (void*)sbstrm;
/* Set function pointers for methods that we inherit from AIJ but override.
B->ops->duplicate = MatDuplicate_SBSTRM;
B->ops->mult = MatMult_SBSTRM;
B->ops->destroy = MatDestroy_MPISBSTRM;
*/
B->ops->assemblyend = MatAssemblyEnd_MPISBSTRM;
/* If A has already been assembled, compute the permutation. */
if (A->assembled) {
ierr = MPISBSTRM_create_sbstrm(B);CHKERRQ(ierr);
}
ierr = PetscObjectChangeTypeName((PetscObject) B, MATMPISBSTRM);CHKERRQ(ierr);
ierr = PetscObjectComposeFunction((PetscObject)B,"MatMPISBAIJSetPreallocation_C",MatMPISBAIJSetPreallocation_MPISBSTRM);CHKERRQ(ierr);
*newmat = B;
PetscFunctionReturn(0);
}
开发者ID:hansec,项目名称:petsc,代码行数:32,代码来源:mpisbstream.c
示例3: MatConvert_MPIAIJ_MPIAIJCRL
PETSC_EXTERN PetscErrorCode MatConvert_MPIAIJ_MPIAIJCRL(Mat A,MatType type,MatReuse reuse,Mat *newmat)
{
PetscErrorCode ierr;
Mat B = *newmat;
Mat_AIJCRL *aijcrl;
PetscFunctionBegin;
if (reuse == MAT_INITIAL_MATRIX) {
ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
}
ierr = PetscNewLog(B,Mat_AIJCRL,&aijcrl);CHKERRQ(ierr);
B->spptr = (void*) aijcrl;
/* Set function pointers for methods that we inherit from AIJ but override. */
B->ops->duplicate = MatDuplicate_AIJCRL;
B->ops->assemblyend = MatAssemblyEnd_MPIAIJCRL;
B->ops->destroy = MatDestroy_MPIAIJCRL;
B->ops->mult = MatMult_AIJCRL;
/* If A has already been assembled, compute the permutation. */
if (A->assembled) {
ierr = MatMPIAIJCRL_create_aijcrl(B);CHKERRQ(ierr);
}
ierr = PetscObjectChangeTypeName((PetscObject)B,MATMPIAIJCRL);CHKERRQ(ierr);
*newmat = B;
PetscFunctionReturn(0);
}
开发者ID:feelpp,项目名称:debian-petsc,代码行数:28,代码来源:mcrl.c
示例4: SNESMonitorJacUpdateSpectrum
PetscErrorCode SNESMonitorJacUpdateSpectrum(SNES snes,PetscInt it,PetscReal fnorm,void *ctx)
{
#if defined(PETSC_MISSING_LAPACK_GEEV)
SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_SUP,"GEEV - Lapack routine is unavailable\nNot able to provide eigen values.");
#elif defined(PETSC_HAVE_ESSL)
SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_SUP,"GEEV - No support for ESSL Lapack Routines");
#else
Vec X;
Mat J,dJ,dJdense;
PetscErrorCode ierr;
PetscErrorCode (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
PetscInt n,i;
PetscBLASInt nb,lwork;
PetscReal *eigr,*eigi;
MatStructure flg = DIFFERENT_NONZERO_PATTERN;
PetscScalar *work;
PetscScalar *a;
PetscFunctionBegin;
if (it == 0) PetscFunctionReturn(0);
/* create the difference between the current update and the current jacobian */
ierr = SNESGetSolution(snes,&X);CHKERRQ(ierr);
ierr = SNESGetJacobian(snes,&J,NULL,&func,NULL);CHKERRQ(ierr);
ierr = MatDuplicate(J,MAT_COPY_VALUES,&dJ);CHKERRQ(ierr);
ierr = SNESComputeJacobian(snes,X,&dJ,&dJ,&flg);CHKERRQ(ierr);
ierr = MatAXPY(dJ,-1.0,J,SAME_NONZERO_PATTERN);CHKERRQ(ierr);
/* compute the spectrum directly */
ierr = MatConvert(dJ,MATSEQDENSE,MAT_INITIAL_MATRIX,&dJdense);CHKERRQ(ierr);
ierr = MatGetSize(dJ,&n,NULL);CHKERRQ(ierr);
ierr = PetscBLASIntCast(n,&nb);CHKERRQ(ierr);
lwork = 3*nb;
ierr = PetscMalloc(n*sizeof(PetscReal),&eigr);CHKERRQ(ierr);
ierr = PetscMalloc(n*sizeof(PetscReal),&eigi);CHKERRQ(ierr);
ierr = PetscMalloc(lwork*sizeof(PetscScalar),&work);CHKERRQ(ierr);
ierr = MatDenseGetArray(dJdense,&a);CHKERRQ(ierr);
#if !defined(PETSC_USE_COMPLEX)
{
PetscBLASInt lierr;
ierr = PetscFPTrapPush(PETSC_FP_TRAP_OFF);CHKERRQ(ierr);
PetscStackCall("LAPACKgeev",LAPACKgeev_("N","N",&nb,a,&nb,eigr,eigi,NULL,&nb,NULL,&nb,work,&lwork,&lierr));
if (lierr) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"geev() error %d",lierr);
ierr = PetscFPTrapPop();CHKERRQ(ierr);
}
#else
SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not coded for complex");
#endif
PetscPrintf(PetscObjectComm((PetscObject)snes),"Eigenvalues of J_%d - J_%d:\n",it,it-1);CHKERRQ(ierr);
for (i=0;i<n;i++) {
PetscPrintf(PetscObjectComm((PetscObject)snes),"%5d: %20.5g + %20.5gi\n",i,eigr[i],eigi[i]);CHKERRQ(ierr);
}
ierr = MatDenseRestoreArray(dJdense,&a);CHKERRQ(ierr);
ierr = MatDestroy(&dJ);CHKERRQ(ierr);
ierr = MatDestroy(&dJdense);CHKERRQ(ierr);
ierr = PetscFree(eigr);CHKERRQ(ierr);
ierr = PetscFree(eigi);CHKERRQ(ierr);
ierr = PetscFree(work);CHKERRQ(ierr);
PetscFunctionReturn(0);
#endif
}
开发者ID:hansec,项目名称:petsc,代码行数:60,代码来源:snesut.c
示例5: MatConvert_SeqSBAIJ_SeqSBSTRM
PETSC_EXTERN PetscErrorCode MatConvert_SeqSBAIJ_SeqSBSTRM(Mat A,MatType type,MatReuse reuse,Mat *newmat)
{
PetscErrorCode ierr;
Mat B = *newmat;
Mat_SeqSBSTRM *sbstrm;
/* PetscInt bs = A->rmap->bs; */
PetscFunctionBegin;
if (reuse == MAT_INITIAL_MATRIX) {
ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
}
ierr = PetscNewLog(B,&sbstrm);CHKERRQ(ierr);
B->spptr = (void*) sbstrm;
/* Set function pointers for methods that we inherit from BAIJ but override. */
B->ops->duplicate = MatDuplicate_SeqSBSTRM;
B->ops->assemblyend = MatAssemblyEnd_SeqSBSTRM;
B->ops->destroy = MatDestroy_SeqSBSTRM;
/*B->ops->choleskyfactorsymbolic = MatCholeskyFactorSymbolic_sbstrm;
B->ops->choleskyfactornumeric = MatCholeskyFactorNumeric_sbstrm; */
/* If A has already been assembled, compute the permutation. */
if (A->assembled) {
ierr = SeqSBSTRM_create_sbstrm(B);CHKERRQ(ierr);
}
ierr = PetscObjectChangeTypeName((PetscObject)B,MATSEQSBSTRM);CHKERRQ(ierr);
*newmat = B;
PetscFunctionReturn(0);
}
开发者ID:ZJLi2013,项目名称:petsc,代码行数:32,代码来源:sbstream.c
示例6: MatCreateSchurComplement
/*@
MatSchurComplementComputeExplicitOperator - Compute the Schur complement matrix explicitly
Collective on Mat
Input Parameter:
. M - the matrix obtained with MatCreateSchurComplement()
Output Parameter:
. S - the Schur complement matrix
Note: This can be expensive, so it is mainly for testing
Level: advanced
.seealso: MatCreateSchurComplement(), MatSchurComplementUpdate()
@*/
PetscErrorCode MatSchurComplementComputeExplicitOperator(Mat M, Mat *S)
{
Mat B, C, D;
KSP ksp;
PC pc;
PetscBool isLU, isILU;
PetscReal fill = 2.0;
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = MatSchurComplementGetSubMatrices(M, NULL, NULL, &B, &C, &D);CHKERRQ(ierr);
ierr = MatSchurComplementGetKSP(M, &ksp);CHKERRQ(ierr);
ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
ierr = PetscObjectTypeCompare((PetscObject) pc, PCLU, &isLU);CHKERRQ(ierr);
ierr = PetscObjectTypeCompare((PetscObject) pc, PCILU, &isILU);CHKERRQ(ierr);
if (isLU || isILU) {
Mat fact, Bd, AinvB, AinvBd;
PetscReal eps = 1.0e-10;
/* This can be sped up for banded LU */
ierr = KSPSetUp(ksp);CHKERRQ(ierr);
ierr = PCFactorGetMatrix(pc, &fact);CHKERRQ(ierr);
ierr = MatConvert(B, MATDENSE, MAT_INITIAL_MATRIX, &Bd);CHKERRQ(ierr);
ierr = MatDuplicate(Bd, MAT_DO_NOT_COPY_VALUES, &AinvBd);CHKERRQ(ierr);
ierr = MatMatSolve(fact, Bd, AinvBd);CHKERRQ(ierr);
ierr = MatDestroy(&Bd);CHKERRQ(ierr);
ierr = MatChop(AinvBd, eps);CHKERRQ(ierr);
ierr = MatConvert(AinvBd, MATAIJ, MAT_INITIAL_MATRIX, &AinvB);CHKERRQ(ierr);
ierr = MatDestroy(&AinvBd);CHKERRQ(ierr);
ierr = MatMatMult(C, AinvB, MAT_INITIAL_MATRIX, fill, S);CHKERRQ(ierr);
ierr = MatDestroy(&AinvB);CHKERRQ(ierr);
} else {
Mat Ainvd, Ainv;
ierr = PCComputeExplicitOperator(pc, &Ainvd);CHKERRQ(ierr);
ierr = MatConvert(Ainvd, MATAIJ, MAT_INITIAL_MATRIX, &Ainv);CHKERRQ(ierr);
ierr = MatDestroy(&Ainvd);CHKERRQ(ierr);
#if 0
/* Symmetric version */
ierr = MatPtAP(Ainv, B, MAT_INITIAL_MATRIX, fill, S);CHKERRQ(ierr);
#else
/* Nonsymmetric version */
ierr = MatMatMatMult(C, Ainv, B, MAT_INITIAL_MATRIX, fill, S);CHKERRQ(ierr);
#endif
ierr = MatDestroy(&Ainv);CHKERRQ(ierr);
}
ierr = PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr);
ierr = MatView(*S, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
if (D) {
MatInfo info;
ierr = MatGetInfo(D, MAT_GLOBAL_SUM, &info);CHKERRQ(ierr);
if (info.nz_used) SETERRQ(PetscObjectComm((PetscObject) M), PETSC_ERR_SUP, "Not yet implemented");
}
PetscFunctionReturn(0);
}
开发者ID:fengyuqi,项目名称:petsc,代码行数:76,代码来源:schurm.c
示例7: MatDuplicate
PetscMatrix<Scalar>* PetscMatrix<Scalar>::duplicate() const
{
PetscMatrix<Scalar>*ptscmatrix = new PetscMatrix<Scalar>();
MatDuplicate(matrix, MAT_COPY_VALUES, &(ptscmatrix->matrix));
ptscmatrix->size = this->size;
ptscmatrix->nnz = nnz;
return ptscmatrix;
};
开发者ID:KirillPlasma,项目名称:hermes,代码行数:8,代码来源:petsc_solver.cpp
示例8: DMCreateMatrix_Shell
static PetscErrorCode DMCreateMatrix_Shell(DM dm,Mat *J)
{
PetscErrorCode ierr;
DM_Shell *shell = (DM_Shell*)dm->data;
Mat A;
PetscFunctionBegin;
PetscValidHeaderSpecific(dm,DM_CLASSID,1);
PetscValidPointer(J,3);
if (!shell->A) {
if (shell->Xglobal) {
PetscInt m,M;
ierr = PetscInfo(dm,"Naively creating matrix using global vector distribution without preallocation\n");
CHKERRQ(ierr);
ierr = VecGetSize(shell->Xglobal,&M);
CHKERRQ(ierr);
ierr = VecGetLocalSize(shell->Xglobal,&m);
CHKERRQ(ierr);
ierr = MatCreate(PetscObjectComm((PetscObject)dm),&shell->A);
CHKERRQ(ierr);
ierr = MatSetSizes(shell->A,m,m,M,M);
CHKERRQ(ierr);
ierr = MatSetType(shell->A,dm->mattype);
CHKERRQ(ierr);
ierr = MatSetUp(shell->A);
CHKERRQ(ierr);
} else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_USER,"Must call DMShellSetMatrix(), DMShellSetCreateMatrix(), or provide a vector");
}
A = shell->A;
/* the check below is tacky and incomplete */
if (dm->mattype) {
PetscBool flg,aij,seqaij,mpiaij;
ierr = PetscObjectTypeCompare((PetscObject)A,dm->mattype,&flg);
CHKERRQ(ierr);
ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&seqaij);
CHKERRQ(ierr);
ierr = PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&mpiaij);
CHKERRQ(ierr);
ierr = PetscStrcmp(dm->mattype,MATAIJ,&aij);
CHKERRQ(ierr);
if (!flg) {
if (!(aij && (seqaij || mpiaij))) SETERRQ2(PetscObjectComm((PetscObject)dm),PETSC_ERR_ARG_NOTSAMETYPE,"Requested matrix of type %s, but only %s available",dm->mattype,((PetscObject)A)->type_name);
}
}
if (((PetscObject)A)->refct < 2) { /* We have an exclusive reference so we can give it out */
ierr = PetscObjectReference((PetscObject)A);
CHKERRQ(ierr);
ierr = MatZeroEntries(A);
CHKERRQ(ierr);
*J = A;
} else { /* Need to create a copy, could use MAT_SHARE_NONZERO_PATTERN in most cases */
ierr = MatDuplicate(A,MAT_DO_NOT_COPY_VALUES,J);
CHKERRQ(ierr);
ierr = MatZeroEntries(*J);
CHKERRQ(ierr);
}
PetscFunctionReturn(0);
}
开发者ID:petsc,项目名称:petsc,代码行数:58,代码来源:dmshell.c
示例9: NEPSolve_Interpol
PetscErrorCode NEPSolve_Interpol(NEP nep)
{
PetscErrorCode ierr;
NEP_INTERPOL *ctx = (NEP_INTERPOL*)nep->data;
Mat *A; /*T=nep->function,Tp=nep->jacobian;*/
PetscScalar *x,*fx,t;
PetscReal *cs,a,b,s;
PetscInt i,j,k,deg=ctx->deg;
PetscFunctionBegin;
ierr = PetscMalloc4(deg+1,&A,(deg+1)*(deg+1),&cs,deg+1,&x,(deg+1)*nep->nt,&fx);CHKERRQ(ierr);
ierr = RGIntervalGetEndpoints(nep->rg,&a,&b,NULL,NULL);CHKERRQ(ierr);
ierr = ChebyshevNodes(deg,a,b,x,cs);CHKERRQ(ierr);
for (j=0;j<nep->nt;j++) {
for (i=0;i<=deg;i++) {
ierr = FNEvaluateFunction(nep->f[j],x[i],&fx[i+j*(deg+1)]);CHKERRQ(ierr);
}
}
/* Polynomial coefficients */
for (k=0;k<=deg;k++) {
ierr = MatDuplicate(nep->A[0],MAT_COPY_VALUES,&A[k]);CHKERRQ(ierr);
t = 0.0;
for (i=0;i<deg+1;i++) t += fx[i]*cs[i*(deg+1)+k];
t *= 2.0/(deg+1);
if (k==0) t /= 2.0;
ierr = MatScale(A[k],t);CHKERRQ(ierr);
for (j=1;j<nep->nt;j++) {
t = 0.0;
for (i=0;i<deg+1;i++) t += fx[i+j*(deg+1)]*cs[i*(deg+1)+k];
t *= 2.0/(deg+1);
if (k==0) t /= 2.0;
ierr = MatAXPY(A[k],t,nep->A[j],SUBSET_NONZERO_PATTERN);CHKERRQ(ierr);
}
}
ierr = PEPSetOperators(ctx->pep,deg+1,A);CHKERRQ(ierr);
for (k=0;k<=deg;k++) {
ierr = MatDestroy(&A[k]);CHKERRQ(ierr);
}
ierr = PetscFree4(A,cs,x,fx);CHKERRQ(ierr);
/* Solve polynomial eigenproblem */
ierr = PEPSolve(ctx->pep);CHKERRQ(ierr);
ierr = PEPGetConverged(ctx->pep,&nep->nconv);CHKERRQ(ierr);
ierr = PEPGetIterationNumber(ctx->pep,&nep->its);CHKERRQ(ierr);
ierr = PEPGetConvergedReason(ctx->pep,(PEPConvergedReason*)&nep->reason);CHKERRQ(ierr);
s = 2.0/(b-a);
for (i=0;i<nep->nconv;i++) {
ierr = PEPGetEigenpair(ctx->pep,i,&nep->eigr[i],&nep->eigi[i],NULL,NULL);CHKERRQ(ierr);
nep->eigr[i] /= s;
nep->eigr[i] += (a+b)/2.0;
nep->eigi[i] /= s;
}
nep->state = NEP_STATE_EIGENVECTORS;
PetscFunctionReturn(0);
}
开发者ID:OpenCMISS-Dependencies,项目名称:slepc,代码行数:57,代码来源:interpol.c
示例10: MatDuplicate_Transpose
PetscErrorCode MatDuplicate_Transpose(Mat N, MatDuplicateOption op, Mat* m)
{
Mat_Transpose *Na = (Mat_Transpose*)N->data;
PetscErrorCode ierr;
PetscFunctionBegin;
if (op == MAT_COPY_VALUES) {
ierr = MatTranspose(Na->A,MAT_INITIAL_MATRIX,m);CHKERRQ(ierr);
} else if (op == MAT_DO_NOT_COPY_VALUES) {
ierr = MatDuplicate(Na->A,MAT_DO_NOT_COPY_VALUES,m);CHKERRQ(ierr);
ierr = MatTranspose(*m,MAT_INPLACE_MATRIX,m);CHKERRQ(ierr);
} else SETERRQ(PetscObjectComm((PetscObject)N),PETSC_ERR_SUP,"MAT_SHARE_NONZERO_PATTERN not supported for this matrix type");
PetscFunctionReturn(0);
}
开发者ID:petsc,项目名称:petsc,代码行数:14,代码来源:transm.c
示例11: MatCopy_User
static PetscErrorCode MatCopy_User(Mat A,Mat B,MatStructure str)
{
User userA,userB;
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(A,&userA);CHKERRQ(ierr);
if (userA) {
ierr = PetscNew(&userB);CHKERRQ(ierr);
ierr = MatDuplicate(userA->A,MAT_COPY_VALUES,&userB->A);CHKERRQ(ierr);
ierr = MatShellSetContext(B, userB);CHKERRQ(ierr);
}
PetscFunctionReturn(0);
}
开发者ID:firedrakeproject,项目名称:petsc,代码行数:14,代码来源:ex205.c
示例12: matrix
/*@C
TaoMatGetSubMat - Gets a submatrix using the IS
Input Parameters:
+ M - the full matrix (n x n)
. is - the index set for the submatrix (both row and column index sets need to be the same)
. v1 - work vector of dimension n, needed for TAO_SUBSET_MASK option
- subset_type - the method TAO is using for subsetting (TAO_SUBSET_SUBVEC, TAO_SUBSET_MASK,
TAO_SUBSET_MATRIXFREE)
Output Parameters:
. Msub - the submatrix
@*/
PetscErrorCode TaoMatGetSubMat(Mat M, IS is, Vec v1, TaoSubsetType subset_type, Mat *Msub)
{
PetscErrorCode ierr;
IS iscomp;
PetscBool flg = PETSC_FALSE;
PetscFunctionBegin;
PetscValidHeaderSpecific(M,MAT_CLASSID,1);
PetscValidHeaderSpecific(is,IS_CLASSID,2);
ierr = MatDestroy(Msub);CHKERRQ(ierr);
switch (subset_type) {
case TAO_SUBSET_SUBVEC:
ierr = MatGetSubMatrix(M, is, is, MAT_INITIAL_MATRIX, Msub);CHKERRQ(ierr);
break;
case TAO_SUBSET_MASK:
/* Get Reduced Hessian
Msub[i,j] = M[i,j] if i,j in Free_Local or i==j
Msub[i,j] = 0 if i!=j and i or j not in Free_Local
*/
ierr = PetscOptionsBegin(PetscObjectComm((PetscObject)M),NULL,NULL,NULL);CHKERRQ(ierr);
ierr = PetscOptionsBool("-different_submatrix","use separate hessian matrix when computing submatrices","TaoSubsetType",flg,&flg,NULL);CHKERRQ(ierr);
ierr = PetscOptionsEnd();CHKERRQ(ierr);
if (flg) {
ierr = MatDuplicate(M, MAT_COPY_VALUES, Msub);CHKERRQ(ierr);
} else {
/* Act on hessian directly (default) */
ierr = PetscObjectReference((PetscObject)M);CHKERRQ(ierr);
*Msub = M;
}
/* Save the diagonal to temporary vector */
ierr = MatGetDiagonal(*Msub,v1);CHKERRQ(ierr);
/* Zero out rows and columns */
ierr = ISComplementVec(is,v1,&iscomp);CHKERRQ(ierr);
/* Use v1 instead of 0 here because of PETSc bug */
ierr = MatZeroRowsColumnsIS(*Msub,iscomp,1.0,v1,v1);CHKERRQ(ierr);
ierr = ISDestroy(&iscomp);CHKERRQ(ierr);
break;
case TAO_SUBSET_MATRIXFREE:
ierr = ISComplementVec(is,v1,&iscomp);CHKERRQ(ierr);
ierr = MatCreateSubMatrixFree(M,iscomp,iscomp,Msub);CHKERRQ(ierr);
ierr = ISDestroy(&iscomp);CHKERRQ(ierr);
break;
}
PetscFunctionReturn(0);
}
开发者ID:tom-klotz,项目名称:petsc,代码行数:62,代码来源:isutil.c
示例13: main
int main(int argc,char **args)
{
PetscErrorCode ierr;
Mat A,B,C;
PetscBool different=PETSC_FALSE,skip=PETSC_FALSE;
PetscInt m0,m1,n=128,i;
PetscInitialize(&argc,&args,(char*)0,help);
ierr = PetscOptionsGetBool(NULL,"-different",&different,NULL);CHKERRQ(ierr);
ierr = PetscOptionsGetBool(NULL,"-skip",&skip,NULL);CHKERRQ(ierr);
/*
Create matrices
A = tridiag(1,-2,1) and B = diag(7);
*/
ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
ierr = MatCreate(PETSC_COMM_WORLD,&B);CHKERRQ(ierr);
ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
ierr = MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,n,n);CHKERRQ(ierr);
ierr = MatSetFromOptions(A);CHKERRQ(ierr);
ierr = MatSetFromOptions(B);CHKERRQ(ierr);
ierr = MatSetUp(A);CHKERRQ(ierr);
ierr = MatSetUp(B);CHKERRQ(ierr);
ierr = MatGetOwnershipRange(A,&m0,&m1);CHKERRQ(ierr);
for (i=m0;i<m1;i++) {
if (i>0) { ierr = MatSetValue(A,i,i-1,-1.0,INSERT_VALUES);CHKERRQ(ierr); }
if (i<n-1) { ierr = MatSetValue(A,i,i+1,-1.0,INSERT_VALUES);CHKERRQ(ierr); }
ierr = MatSetValue(A,i,i,2.0,INSERT_VALUES);CHKERRQ(ierr);
ierr = MatSetValue(B,i,i,7.0,INSERT_VALUES);CHKERRQ(ierr);
}
ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatDuplicate(A,MAT_COPY_VALUES,&C);CHKERRQ(ierr);
/* Add B */
ierr = MatAXPY(C,1.0,B,(different)?DIFFERENT_NONZERO_PATTERN:SUBSET_NONZERO_PATTERN);CHKERRQ(ierr);
/* Add A */
if (!skip) { ierr = MatAXPY(C,1.0,A,SUBSET_NONZERO_PATTERN);CHKERRQ(ierr); }
/*
Free memory
*/
ierr = MatDestroy(&A);CHKERRQ(ierr);
ierr = MatDestroy(&B);CHKERRQ(ierr);
ierr = MatDestroy(&C);CHKERRQ(ierr);
ierr = PetscFinalize();
return 0;
}
开发者ID:fengyuqi,项目名称:petsc,代码行数:49,代码来源:ex172.c
示例14: TestMatZeroRows_with_no_allocation
PetscErrorCode TestMatZeroRows_with_no_allocation(Mat A,IS is,PetscScalar diag)
{
Mat B;
PetscErrorCode ierr;
/* Now copy A into B, and test it with MatZeroRows() */
ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
/* Set this flag after assembly. This way, it affects only MatZeroRows() */
ierr = MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr);
ierr = MatZeroRowsIS(B,is,diag,0,0);CHKERRQ(ierr);
ierr = MatView(B,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
ierr = MatDestroy(&B);CHKERRQ(ierr);
return 0;
}
开发者ID:tom-klotz,项目名称:petsc,代码行数:15,代码来源:ex12.c
示例15: MatConvert_MPIAIJ_MPIAIJMKL
PETSC_INTERN PetscErrorCode MatConvert_MPIAIJ_MPIAIJMKL(Mat A,MatType type,MatReuse reuse,Mat *newmat)
{
PetscErrorCode ierr;
Mat B = *newmat;
PetscFunctionBegin;
if (reuse == MAT_INITIAL_MATRIX) {
ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
}
ierr = PetscObjectChangeTypeName((PetscObject) B, MATMPIAIJMKL);CHKERRQ(ierr);
ierr = PetscObjectComposeFunction((PetscObject)B,"MatMPIAIJSetPreallocation_C",MatMPIAIJSetPreallocation_MPIAIJMKL);CHKERRQ(ierr);
*newmat = B;
PetscFunctionReturn(0);
}
开发者ID:firedrakeproject,项目名称:petsc,代码行数:15,代码来源:mpiaijmkl.c
示例16: VecCreateMPI
void StaggeredStokesPETScLevelSolver::initializeSolverStateSpecialized(
const SAMRAIVectorReal<NDIM, double>& x,
const SAMRAIVectorReal<NDIM, double>& /*b*/)
{
// Allocate DOF index data.
Pointer<PatchLevel<NDIM> > level = d_hierarchy->getPatchLevel(d_level_num);
if (!level->checkAllocated(d_u_dof_index_idx)) level->allocatePatchData(d_u_dof_index_idx);
if (!level->checkAllocated(d_p_dof_index_idx)) level->allocatePatchData(d_p_dof_index_idx);
// Setup PETSc objects.
int ierr;
StaggeredStokesPETScVecUtilities::constructPatchLevelDOFIndices(
d_num_dofs_per_proc, d_u_dof_index_idx, d_p_dof_index_idx, level);
const int mpi_rank = SAMRAI_MPI::getRank();
ierr = VecCreateMPI(
PETSC_COMM_WORLD, d_num_dofs_per_proc[mpi_rank], PETSC_DETERMINE, &d_petsc_x);
IBTK_CHKERRQ(ierr);
ierr = VecCreateMPI(
PETSC_COMM_WORLD, d_num_dofs_per_proc[mpi_rank], PETSC_DETERMINE, &d_petsc_b);
IBTK_CHKERRQ(ierr);
StaggeredStokesPETScMatUtilities::constructPatchLevelMACStokesOp(d_petsc_mat,
d_U_problem_coefs,
d_U_bc_coefs,
d_new_time,
d_num_dofs_per_proc,
d_u_dof_index_idx,
d_p_dof_index_idx,
level);
ierr = MatDuplicate(d_petsc_mat, MAT_COPY_VALUES, &d_petsc_pc);
IBTK_CHKERRQ(ierr);
HierarchyDataOpsManager<NDIM>* hier_ops_manager =
HierarchyDataOpsManager<NDIM>::getManager();
Pointer<HierarchyDataOpsInteger<NDIM> > hier_p_dof_index_ops =
hier_ops_manager->getOperationsInteger(d_p_dof_index_var, d_hierarchy, true);
hier_p_dof_index_ops->resetLevels(d_level_num, d_level_num);
const int min_p_idx = hier_p_dof_index_ops->min(
d_p_dof_index_idx); // NOTE: HierarchyDataOpsInteger::max() is broken
ierr = MatZeroRowsColumns(d_petsc_pc, 1, &min_p_idx, 1.0, NULL, NULL);
IBTK_CHKERRQ(ierr);
d_petsc_ksp_ops_flag = SAME_PRECONDITIONER;
const int u_idx = x.getComponentDescriptorIndex(0);
const int p_idx = x.getComponentDescriptorIndex(1);
d_data_synch_sched =
StaggeredStokesPETScVecUtilities::constructDataSynchSchedule(u_idx, p_idx, level);
d_ghost_fill_sched =
StaggeredStokesPETScVecUtilities::constructGhostFillSchedule(u_idx, p_idx, level);
return;
} // initializeSolverStateSpecialized
开发者ID:BijanZarif,项目名称:IBAMR,代码行数:48,代码来源:StaggeredStokesPETScLevelSolver.cpp
示例17: RawGraph
RawGraph(Mat new_global_mat) {
global_mat = new_global_mat;
//For this code to work, we're going to need the matrix structure for the sequential portion of this matrix.
//Therefore, extract a sequential AIJ matrix.
MatType type;
MatGetType(global_mat, &type);
if (!strcmp(type,MATSEQAIJ)) {
MatDuplicate(global_mat, MAT_DO_NOT_COPY_VALUES, &local_mat);
} else {
MatGetLocalMat(global_mat, MAT_INITIAL_MATRIX, &local_mat);
}
MatZeroEntries(local_mat);
MatGetOwnershipRange(global_mat, &row_begin, &row_end);
//MatView(local_mat, PETSC_VIEWER_DRAW_SELF);
seq_raw.create(local_mat);
}
开发者ID:rblake,项目名称:petsc_amg,代码行数:17,代码来源:mglib.c
示例18: main
int main(int argc,char **args)
{
Mat A,Ar,C;
PetscViewer fd; /* viewer */
char file[PETSC_MAX_PATH_LEN]; /* input file name */
PetscErrorCode ierr;
PetscInt ns=2,np;
PetscSubcomm subc;
PetscBool flg;
PetscInitialize(&argc,&args,(char*)0,help);
/*
Determine files from which we read the two linear systems
(matrix and right-hand-side vector).
*/
ierr = PetscOptionsGetString(NULL,"-f0",file,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr);
if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate binary file with the -f0 option");
ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);CHKERRQ(ierr);
ierr = MPI_Comm_size(PETSC_COMM_WORLD,&np);CHKERRQ(ierr);
ierr = PetscPrintf(PETSC_COMM_WORLD,"Reading matrix with %d processors\n",np);CHKERRQ(ierr);
ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
ierr = MatLoad(A,fd);CHKERRQ(ierr);
ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr);
/*
Determines amount of subcomunicators
*/
ierr = PetscOptionsGetInt(NULL,"-nsub",&ns,NULL);CHKERRQ(ierr);
ierr = PetscPrintf(PETSC_COMM_WORLD,"Splitting in %d subcommunicators\n",ns);CHKERRQ(ierr);
ierr = PetscSubcommCreate(PetscObjectComm((PetscObject)A),&subc);CHKERRQ(ierr);
ierr = PetscSubcommSetNumber(subc,ns);CHKERRQ(ierr);
ierr = PetscSubcommSetType(subc,PETSC_SUBCOMM_CONTIGUOUS);CHKERRQ(ierr);
ierr = PetscSubcommSetFromOptions(subc);CHKERRQ(ierr);
ierr = MatCreateRedundantMatrix(A,0,PetscSubcommChild(subc),MAT_INITIAL_MATRIX,&Ar);CHKERRQ(ierr);
ierr = PetscPrintf(PETSC_COMM_WORLD,"Copying matrix\n",ns);CHKERRQ(ierr);
ierr = MatDuplicate(Ar,MAT_COPY_VALUES,&C);CHKERRQ(ierr);
ierr = PetscSubcommDestroy(&subc);CHKERRQ(ierr);
/*
Free memory
*/
ierr = MatDestroy(&A);CHKERRQ(ierr);
ierr = MatDestroy(&Ar);CHKERRQ(ierr);
ierr = MatDestroy(&C);CHKERRQ(ierr);
ierr = PetscFinalize();
return 0;
}
开发者ID:fengyuqi,项目名称:petsc,代码行数:46,代码来源:ex169.c
示例19: interpPeriodicMatrix
PetscErrorCode interpPeriodicMatrix(PetscScalar tc, Mat *A, PetscScalar cyclePeriod,
PetscInt numPeriods, PetscScalar *tdp,
PeriodicMat *user, char *filename)
{
/* Function to interpolate a matrix that is periodic in time with period cyclePeriod. */
/* tc is the current time and numPeriods is the number of instances per period */
/* at which data are available (to be read from files). */
/* IMPORTANT: Matrix *A MUST have been created and preallocated before */
/* calling this routine. All other matrices will be created using *A as a template. */
#include <math.h>
PetscScalar t,t1;
PetscInt im,it0,it1;
/* static PetscInt iCurrTimeReadLast=-1; */
PetscErrorCode ierr;
PetscScalar alpha[2];
char tmpFile[PETSC_MAX_PATH_LEN];
if (user->firstTime) {
if (numPeriods>MAX_MATRIX_PERIODS) {
SETERRQ(1,"Number of allowable matrices in PeriodicMat struct exceeded by requested number ! Increase MAX_MATRIX_PERIODS.");
}
user->numPeriods = numPeriods;
for (im=0; im<numPeriods; im++) {
ierr = MatDuplicate(*A,MAT_DO_NOT_COPY_VALUES,&user->Ap[im]);CHKERRQ(ierr);
strcpy(tmpFile,"");
sprintf(tmpFile,"%s%02d",filename,im);
ierr = PetscPrintf(PETSC_COMM_WORLD,"Reading matrix from file %s\n", tmpFile);CHKERRQ(ierr);
ierr = MatLoadIntoMatrix3(tmpFile,user->Ap[im]);
}
user->firstTime = PETSC_FALSE;
}
t=tc; /* current time */
if (t<0.) t=cyclePeriod+t;
t1=t-cyclePeriod*floor(t/cyclePeriod);
ierr=calcPeriodicInterpFactor(numPeriods,t1,tdp,&it0,&it1,&alpha[0],&alpha[1]); CHKERRQ(ierr);
/* ierr = PetscPrintf(PETSC_COMM_WORLD,"tc=%lf,t1=%lf,it0=%d,it1=%d,a1=%17.16lf,a2=%17.16lf\n",tc,t1,it0,it1,alpha[0],alpha[1]);CHKERRQ(ierr); */
/* interpolate to current time */
ierr = MatAXPBYmy(alpha[0],alpha[1],user->Ap[it0],user->Ap[it1],A);CHKERRQ(ierr);
return 0;
}
开发者ID:samarkhatiwala,项目名称:tmm,代码行数:45,代码来源:tmm_forcing_utils.c
示例20: MatDuplicate_IS
PetscErrorCode MatDuplicate_IS(Mat mat,MatDuplicateOption op,Mat *newmat)
{
PetscErrorCode ierr;
Mat_IS *matis = (Mat_IS*)(mat->data);
PetscInt bs,m,n,M,N;
Mat B,localmat;
PetscFunctionBegin;
ierr = MatGetBlockSize(mat,&bs);CHKERRQ(ierr);
ierr = MatGetSize(mat,&M,&N);CHKERRQ(ierr);
ierr = MatGetLocalSize(mat,&m,&n);CHKERRQ(ierr);
ierr = MatCreateIS(PetscObjectComm((PetscObject)mat),bs,m,n,M,N,matis->mapping,&B);CHKERRQ(ierr);
ierr = MatDuplicate(matis->A,op,&localmat);CHKERRQ(ierr);
ierr = MatISSetLocalMat(B,localmat);CHKERRQ(ierr);
ierr = MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
*newmat = B;
PetscFunctionReturn(0);
}
开发者ID:PeiLiu90,项目名称:petsc,代码行数:19,代码来源:matis.c
注:本文中的MatDuplicate函数示例由纯净天空整理自Github/MSDocs等源码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。 |
请发表评论