本文整理汇总了C++中MatGetDiagonal函数的典型用法代码示例。如果您正苦于以下问题:C++ MatGetDiagonal函数的具体用法?C++ MatGetDiagonal怎么用?C++ MatGetDiagonal使用的例子?那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。
在下文中一共展示了MatGetDiagonal函数的20个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的C++代码示例。
示例1: BSSCR_MatStokesMVBlockDefaultBuildScaling
// updated
PetscErrorCode BSSCR_MatStokesMVBlockDefaultBuildScaling( MatStokesBlockScaling BA, Mat A, Vec b, Vec x, PetscTruth sym )
{
Mat K,G,D,C;
Vec rG;
PetscScalar rg2, rg, ra;
PetscInt N;
Vec rA, rC;
Vec L1,L2, R1,R2;
Mat S;
VecNestGetSubVec( BA->Lz, 0, &L1 );
VecNestGetSubVec( BA->Lz, 1, &L2 );
VecNestGetSubVec( BA->Rz, 0, &R1 );
VecNestGetSubVec( BA->Rz, 1, &R2 );
rA = L1;
rC = L2;
MatNestGetSubMat( A, 0,0, &K );
MatNestGetSubMat( A, 0,1, &G );
MatNestGetSubMat( A, 1,0, &D );
MatNestGetSubMat( A, 1,1, &C );
VecDuplicate( rA, &rG );
/* Get magnitude of K */
//px_MatGetAbsRowSum( K, rA );
//MatGetRowMax( K, rA, PETSC_NULL );
MatGetDiagonal( K, rA);
VecSqrt( rA );
VecReciprocal( rA );
/* VecDot( rA,rA, &ra ); */
/* VecGetSize( rA, &N ); */
/* ra = PetscSqrtScalar( ra/N ); */
/* Get magnitude of G */
//px_MatGetAbsRowSum( G, rG );
//MatGetRowMax( G, rG, PETSC_NULL );
Mat A21_cpy;
Mat Shat;
Vec diag; /* same as rA*rA */
MatGetVecs( K, &diag, PETSC_NULL );
MatGetDiagonal( K, diag );
VecReciprocal( diag );
if( sym )
#if( PETSC_VERSION_MAJOR <= 2 )
MatTranspose( G, &A21_cpy );
#else
MatTranspose( G, MAT_INITIAL_MATRIX, &A21_cpy );
#endif
else {
开发者ID:AngelValverdePerez,项目名称:underworld2,代码行数:61,代码来源:stokes_mvblock_scaling.c
示例2: main
int main(int argc,char **args)
{
Mat A;
Vec x;
PetscErrorCode ierr;
PetscViewer fd; /* viewer */
char file[PETSC_MAX_PATH_LEN]; /* input file name */
PetscReal norm;
PetscBool flg;
PetscInitialize(&argc,&args,(char*)0,help);
/* Determine file from which we read the matrix A */
ierr = PetscOptionsGetString(NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr);
if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate binary file with the -f option");
/* Load matrix A */
ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);CHKERRQ(ierr);
ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
ierr = MatLoad(A,fd);CHKERRQ(ierr);
ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr);
ierr = MatCreateVecs(A,&x,NULL);CHKERRQ(ierr);
ierr = MatGetDiagonal(A,x);CHKERRQ(ierr);
ierr = VecScale(x,-1.0);CHKERRQ(ierr);
ierr = MatDiagonalSet(A,x,ADD_VALUES);CHKERRQ(ierr);
ierr = MatGetDiagonal(A,x);CHKERRQ(ierr);
ierr = VecNorm(x,NORM_2,&norm);CHKERRQ(ierr);
ierr = PetscPrintf(PETSC_COMM_WORLD,"Norm %g\n",(double)norm);CHKERRQ(ierr);
/* Free data structures */
ierr = MatDestroy(&A);CHKERRQ(ierr);
ierr = VecDestroy(&x);CHKERRQ(ierr);
ierr = PetscFinalize();
return 0;
}
开发者ID:PeiLiu90,项目名称:petsc,代码行数:35,代码来源:ex171.c
示例3: bsscr_DGMiGtD
PetscErrorCode bsscr_DGMiGtD( Mat *_K2, Mat K, Mat G, Mat M){
Mat K2;
Vec diag;
Vec Mdiag;
Mat MinvGt;
Mat Gtrans;
PetscErrorCode ierr;
PetscFunctionBegin;
MatGetVecs( K, &diag, PETSC_NULL );
MatGetDiagonal( K, diag );
VecSqrt(diag);
//VecReciprocal(diag);/* trying something different here */
MatGetVecs( M, &Mdiag, PETSC_NULL );
MatGetDiagonal( M, Mdiag );
VecReciprocal(Mdiag);
#if( PETSC_VERSION_MAJOR <= 2 )
ierr=MatTranspose(G, &Gtrans);CHKERRQ(ierr);
#else
ierr=MatTranspose(G, MAT_INITIAL_MATRIX,&Gtrans);CHKERRQ(ierr);
#endif
ierr=MatConvert(Gtrans, MATSAME, MAT_INITIAL_MATRIX, &MinvGt);CHKERRQ(ierr);/* copy Gtrans -> MinvGt */
MatDiagonalScale(MinvGt, Mdiag, PETSC_NULL);/* Minv*Gtrans */
/* MAT_INITIAL_MATRIX -> creates K2 matrix : PETSC_DEFAULT for fill ratio: run with -info to find what it should be*/
ierr=MatMatMult( G, MinvGt, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &K2);CHKERRQ(ierr);/* K2 = G*Minv*Gtrans */
MatDiagonalScale(K2, diag, diag );/* K2 = D*K2*D = D*G*Minv*Gtrans*D */
Stg_MatDestroy(&Gtrans);
Stg_MatDestroy(&MinvGt);
Stg_VecDestroy(&diag);
Stg_VecDestroy(&Mdiag);
*_K2=K2;
PetscFunctionReturn(0);
}
开发者ID:dansand,项目名称:underworld2,代码行数:33,代码来源:createK2.c
示例4: StokesSetupApproxSchur
PetscErrorCode StokesSetupApproxSchur(Stokes *s)
{
Vec diag;
PetscErrorCode ierr;
PetscFunctionBeginUser;
/* Schur complement approximation: myS = A11 - A10 diag(A00)^(-1) A01 */
/* note: A11 is zero */
/* note: in real life this matrix would be build directly, */
/* i.e. without MatMatMult */
/* inverse of diagonal of A00 */
ierr = VecCreate(PETSC_COMM_WORLD,&diag);CHKERRQ(ierr);
ierr = VecSetSizes(diag,PETSC_DECIDE,2*s->nx*s->ny);CHKERRQ(ierr);
ierr = VecSetType(diag,VECMPI);CHKERRQ(ierr);
ierr = MatGetDiagonal(s->subA[0],diag);
ierr = VecReciprocal(diag);
/* compute: - A10 diag(A00)^(-1) A01 */
ierr = MatDiagonalScale(s->subA[1],diag,NULL); /* (*warning* overwrites subA[1]) */
ierr = MatMatMult(s->subA[2],s->subA[1],MAT_INITIAL_MATRIX,PETSC_DEFAULT,&s->myS);CHKERRQ(ierr);
ierr = MatScale(s->myS,-1.0);CHKERRQ(ierr);
/* restore A10 */
ierr = MatGetDiagonal(s->subA[0],diag);
ierr = MatDiagonalScale(s->subA[1],diag,NULL);
ierr = VecDestroy(&diag);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:rjayawar,项目名称:femus,代码行数:29,代码来源:ex2.cpp
示例5: PCSetUp_Eisenstat
static PetscErrorCode PCSetUp_Eisenstat(PC pc)
{
PetscErrorCode ierr;
PetscInt M,N,m,n;
PC_Eisenstat *eis = (PC_Eisenstat*)pc->data;
PetscFunctionBegin;
if (!pc->setupcalled) {
ierr = MatGetSize(pc->mat,&M,&N);CHKERRQ(ierr);
ierr = MatGetLocalSize(pc->mat,&m,&n);CHKERRQ(ierr);
ierr = MatCreate(((PetscObject)pc)->comm,&eis->shell);CHKERRQ(ierr);
ierr = MatSetSizes(eis->shell,m,n,M,N);CHKERRQ(ierr);
ierr = MatSetType(eis->shell,MATSHELL);CHKERRQ(ierr);
ierr = MatSetUp(eis->shell);CHKERRQ(ierr);
ierr = MatShellSetContext(eis->shell,(void*)pc);CHKERRQ(ierr);
ierr = PetscLogObjectParent(pc,eis->shell);CHKERRQ(ierr);
ierr = MatShellSetOperation(eis->shell,MATOP_MULT,(void(*)(void))PCMult_Eisenstat);CHKERRQ(ierr);
}
if (!eis->usediag) PetscFunctionReturn(0);
if (!pc->setupcalled) {
ierr = MatGetVecs(pc->pmat,&eis->diag,0);CHKERRQ(ierr);
ierr = PetscLogObjectParent(pc,eis->diag);CHKERRQ(ierr);
}
ierr = MatGetDiagonal(pc->pmat,eis->diag);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:erdc-cm,项目名称:petsc-dev,代码行数:26,代码来源:eisen.c
示例6: PCSetUp_LSC
static PetscErrorCode PCSetUp_LSC(PC pc)
{
PC_LSC *lsc = (PC_LSC*)pc->data;
Mat L,Lp,B,C;
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = PCLSCAllocate_Private(pc);CHKERRQ(ierr);
ierr = PetscObjectQuery((PetscObject)pc->mat,"LSC_L",(PetscObject*)&L);CHKERRQ(ierr);
if (!L) {ierr = PetscObjectQuery((PetscObject)pc->pmat,"LSC_L",(PetscObject*)&L);CHKERRQ(ierr);}
ierr = PetscObjectQuery((PetscObject)pc->pmat,"LSC_Lp",(PetscObject*)&Lp);CHKERRQ(ierr);
if (!Lp) {ierr = PetscObjectQuery((PetscObject)pc->mat,"LSC_Lp",(PetscObject*)&Lp);CHKERRQ(ierr);}
if (!L) {
ierr = MatSchurComplementGetSubMatrices(pc->mat,NULL,NULL,&B,&C,NULL);CHKERRQ(ierr);
if (!lsc->L) {
ierr = MatMatMult(C,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&lsc->L);CHKERRQ(ierr);
} else {
ierr = MatMatMult(C,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&lsc->L);CHKERRQ(ierr);
}
Lp = L = lsc->L;
}
if (lsc->scale) {
Mat Ap;
ierr = MatSchurComplementGetSubMatrices(pc->mat,NULL,&Ap,NULL,NULL,NULL);CHKERRQ(ierr);
ierr = MatGetDiagonal(Ap,lsc->scale);CHKERRQ(ierr); /* Should be the mass matrix, but we don't have plumbing for that yet */
ierr = VecReciprocal(lsc->scale);CHKERRQ(ierr);
}
ierr = KSPSetOperators(lsc->kspL,L,Lp);CHKERRQ(ierr);
ierr = KSPSetFromOptions(lsc->kspL);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:pombredanne,项目名称:petsc,代码行数:31,代码来源:lsc.c
示例7: homemade_assert_msg
void FETI_Operations::set_jacobi_precond_vector()
{
homemade_assert_msg(m_bC_RR_MatrixSet,"Preconditioner matrix not set yet!");
// Create and set the vector
VecCreate(m_comm.get(),&m_coupling_jacobi_precond_vec);
VecSetSizes(m_coupling_jacobi_precond_vec,m_C_RR_M_local,m_C_RR_M);
VecSetFromOptions(m_coupling_jacobi_precond_vec);
// Get the diagonal
MatGetDiagonal(m_C_RR,m_coupling_jacobi_precond_vec);
// Calculate the reciprocal
VecReciprocal(m_coupling_jacobi_precond_vec);
// Export it
write_PETSC_vector(m_coupling_jacobi_precond_vec,m_scratch_folder_path + "/precond_Jacobi_vector.petscvec",m_comm.rank(),m_comm.get());
// Print MatLab debugging output? Variable defined at "carl_headers.h"
#ifdef PRINT_MATLAB_DEBUG
write_PETSC_vector_MATLAB(m_coupling_jacobi_precond_vec,m_scratch_folder_path + "/precond_Jacobi_vector.m",m_comm.get());
#endif
// Set flag
m_bCreatedPrecondJacobiVec = true;
}
开发者ID:cottereau,项目名称:CArl,代码行数:26,代码来源:FETI_operations_setup.cpp
示例8: MatGetDiagonal
void PetscSparseStorage::getDiagonal( OoqpVector& vec_in )
{
PetscVector & vec = dynamic_cast<PetscVector &>(vec_in);
int ierr;
ierr = MatGetDiagonal(M, vec.pv); assert( ierr == 0);
}
开发者ID:aig-lchion,项目名称:OOQP,代码行数:7,代码来源:PetscSparseStorage.C
示例9: applyConstraintsToLoadIncrement
void
NRSolver :: applyConstraintsToLoadIncrement(int nite, const SparseMtrx &k, FloatArray &R,
referenceLoadInputModeType rlm, TimeStep *tStep)
{
double factor = engngModel->giveDomain(1)->giveFunction(prescribedDisplacementTF)->evaluateAtTime( tStep->giveTargetTime() );
if ( ( rlm == rlm_total ) && ( !tStep->isTheFirstStep() ) ) {
//factor -= engngModel->giveDomain(1)->giveFunction(prescribedDisplacementTF)->
// at(tStep->givePreviousStep()->giveTime()) ;
factor -= engngModel->giveDomain(1)->giveFunction(prescribedDisplacementTF)->
evaluateAtTime( tStep->giveTargetTime() - tStep->giveTimeIncrement() );
}
if ( nite == 0 ) {
#if 0
#ifdef __PETSC_MODULE
if ( solverType == ST_Petsc ) {
//Natural2LocalOrdering* n2lpm = engngModel->giveParallelContext(1)->giveN2Lmap();
//IntArray* map = n2lpm->giveN2Lmap();
for ( i = 1; i <= prescribedEqs.giveSize(); i++ ) {
eq = prescribedEqs.at(i);
R.at(eq) = prescribedDofsValues.at(i) * factor; // local eq
}
return;
}
#endif
#else
#ifdef __PETSC_MODULE
const PetscSparseMtrx *lhs = dynamic_cast< const PetscSparseMtrx * >(&k);
if ( lhs ) {
Vec diag;
PetscScalar *ptr;
lhs->createVecGlobal(& diag);
MatGetDiagonal(* ( const_cast< PetscSparseMtrx * >(lhs)->giveMtrx() ), diag);
VecGetArray(diag, & ptr);
for ( int i = 1; i <= numberOfPrescribedDofs; i++ ) {
int eq = prescribedEqs.at(i) - 1;
R.at(eq + 1) = ptr [ eq ] * prescribedDofsValues.at(i) * factor;
}
VecRestoreArray(diag, & ptr);
VecDestroy(& diag);
return;
}
#endif
#endif
for ( int i = 1; i <= numberOfPrescribedDofs; i++ ) {
int eq = prescribedEqs.at(i);
R.at(eq) = k.at(eq, eq) * prescribedDofsValues.at(i) * factor;
}
} else {
for ( int i = 1; i <= numberOfPrescribedDofs; i++ ) {
R.at( prescribedEqs.at(i) ) = 0.0;
}
}
}
开发者ID:srinath-chakravarthy,项目名称:oofem_dd,代码行数:58,代码来源:nrsolver.C
示例10: MatGetDiagonal_SMF
PetscErrorCode MatGetDiagonal_SMF(Mat mat,Vec v)
{
PetscErrorCode ierr;
MatSubMatFreeCtx ctx;
PetscFunctionBegin;
ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr);
ierr = MatGetDiagonal(ctx->A,v);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:tom-klotz,项目名称:petsc,代码行数:10,代码来源:submatfree.c
示例11: MatGetDiagonal_User
static PetscErrorCode MatGetDiagonal_User(Mat A,Vec X)
{
User user;
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = MatShellGetContext(A,&user);CHKERRQ(ierr);
ierr = MatGetDiagonal(user->B,X);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:petsc,项目名称:petsc,代码行数:10,代码来源:ex221.c
示例12: PCGAMGOptProl_Classical_Jacobi
PetscErrorCode PCGAMGOptProl_Classical_Jacobi(PC pc,const Mat A,Mat *P)
{
PetscErrorCode ierr;
PetscInt f,s,n,cf,cs,i,idx;
PetscInt *coarserows;
PetscInt ncols;
const PetscInt *pcols;
const PetscScalar *pvals;
Mat Pnew;
Vec diag;
PC_MG *mg = (PC_MG*)pc->data;
PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx;
PC_GAMG_Classical *cls = (PC_GAMG_Classical*)pc_gamg->subctx;
PetscFunctionBegin;
if (cls->nsmooths == 0) {
ierr = PCGAMGTruncateProlongator_Private(pc,P);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
ierr = MatGetOwnershipRange(*P,&s,&f);CHKERRQ(ierr);
n = f-s;
ierr = MatGetOwnershipRangeColumn(*P,&cs,&cf);CHKERRQ(ierr);
ierr = PetscMalloc(sizeof(PetscInt)*n,&coarserows);CHKERRQ(ierr);
/* identify the rows corresponding to coarse unknowns */
idx = 0;
for (i=s;i<f;i++) {
ierr = MatGetRow(*P,i,&ncols,&pcols,&pvals);CHKERRQ(ierr);
/* assume, for now, that it's a coarse unknown if it has a single unit entry */
if (ncols == 1) {
if (pvals[0] == 1.) {
coarserows[idx] = i;
idx++;
}
}
ierr = MatRestoreRow(*P,i,&ncols,&pcols,&pvals);CHKERRQ(ierr);
}
ierr = MatGetVecs(A,&diag,0);CHKERRQ(ierr);
ierr = MatGetDiagonal(A,diag);CHKERRQ(ierr);
ierr = VecReciprocal(diag);CHKERRQ(ierr);
for (i=0;i<cls->nsmooths;i++) {
ierr = MatMatMult(A,*P,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&Pnew);CHKERRQ(ierr);
ierr = MatZeroRows(Pnew,idx,coarserows,0.,NULL,NULL);CHKERRQ(ierr);
ierr = MatDiagonalScale(Pnew,diag,0);CHKERRQ(ierr);
ierr = MatAYPX(Pnew,-1.0,*P,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
ierr = MatDestroy(P);CHKERRQ(ierr);
*P = Pnew;
Pnew = NULL;
}
ierr = VecDestroy(&diag);CHKERRQ(ierr);
ierr = PetscFree(coarserows);CHKERRQ(ierr);
ierr = PCGAMGTruncateProlongator_Private(pc,P);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:ZJLi2013,项目名称:petsc,代码行数:54,代码来源:classical.c
示例13: SampleShellPCApply
/*
SampleShellPCSetUp - This routine sets up a user-defined
preconditioner context.
Input Parameters:
. pc - preconditioner object
. pmat - preconditioner matrix
. x - vector
Output Parameter:
. shell - fully set up user-defined preconditioner context
Notes:
In this example, we define the shell preconditioner to be Jacobi's
method. Thus, here we create a work vector for storing the reciprocal
of the diagonal of the preconditioner matrix; this vector is then
used within the routine SampleShellPCApply().
*/
PetscErrorCode SampleShellPCSetUp(PC pc,Mat pmat,Vec x)
{
SampleShellPC *shell;
Vec diag;
PetscErrorCode ierr;
ierr = PCShellGetContext(pc,(void**)&shell);CHKERRQ(ierr);
ierr = VecDuplicate(x,&diag);CHKERRQ(ierr);
ierr = MatGetDiagonal(pmat,diag);CHKERRQ(ierr);
ierr = VecReciprocal(diag);CHKERRQ(ierr);
shell->diag = diag;
return 0;
}
开发者ID:00liujj,项目名称:petsc,代码行数:32,代码来源:ex15.c
示例14: 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
示例15: MatGetDiagonal_IS
PetscErrorCode MatGetDiagonal_IS(Mat A, Vec v)
{
Mat_IS *is = (Mat_IS*)A->data;
PetscErrorCode ierr;
PetscFunctionBegin;
/* get diagonal of the local matrix */
ierr = MatGetDiagonal(is->A,is->y);CHKERRQ(ierr);
/* scatter diagonal back into global vector */
ierr = VecSet(v,0);CHKERRQ(ierr);
ierr = VecScatterBegin(is->rctx,is->y,v,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
ierr = VecScatterEnd(is->rctx,is->y,v,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:masa-ito,项目名称:PETScToPoisson,代码行数:15,代码来源:matis.c
示例16: applyConstraintsToStiffness
void
NRSolver :: applyConstraintsToStiffness(SparseMtrx &k)
{
if ( this->smConstraintVersion == k.giveVersion() ) {
return;
}
#ifdef __PETSC_MODULE
PetscSparseMtrx *lhs = dynamic_cast< PetscSparseMtrx * >(&k);
if ( lhs ) {
Vec diag;
PetscScalar *ptr;
int eq;
lhs->createVecGlobal(& diag);
MatGetDiagonal(* lhs->giveMtrx(), diag);
VecGetArray(diag, & ptr);
for ( int i = 1; i <= numberOfPrescribedDofs; i++ ) {
eq = prescribedEqs.at(i) - 1;
MatSetValue(* ( lhs->giveMtrx() ), eq, eq, ptr [ eq ] * 1.e6, INSERT_VALUES);
}
MatAssemblyBegin(* lhs->giveMtrx(), MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(* lhs->giveMtrx(), MAT_FINAL_ASSEMBLY);
VecRestoreArray(diag, & ptr);
VecDestroy(& diag);
if ( numberOfPrescribedDofs ) {
this->smConstraintVersion = k.giveVersion();
}
return;
}
#endif // __PETSC_MODULE
for ( int i = 1; i <= numberOfPrescribedDofs; i++ ) {
k.at( prescribedEqs.at(i), prescribedEqs.at(i) ) *= 1.e6;
}
if ( numberOfPrescribedDofs ) {
this->smConstraintVersion = k.giveVersion();
}
}
开发者ID:srinath-chakravarthy,项目名称:oofem_dd,代码行数:42,代码来源:nrsolver.C
示例17: PETSC_VERSION_LESS_THAN
void PetscMatrix<T>::get_diagonal (NumericVector<T>& dest) const
{
// Make sure the NumericVector passed in is really a PetscVector
PetscVector<T>& petsc_dest = libmesh_cast_ref<PetscVector<T>&>(dest);
// Call PETSc function.
#if PETSC_VERSION_LESS_THAN(2,3,1)
libMesh::out << "This method has been developed with PETSc 2.3.1. "
<< "No one has made it backwards compatible with older "
<< "versions of PETSc so far; however, it might work "
<< "without any change with some older version." << std::endl;
libmesh_error();
#else
// Needs a const_cast since PETSc does not work with const.
PetscErrorCode ierr =
MatGetDiagonal(const_cast<PetscMatrix<T>*>(this)->mat(),petsc_dest.vec()); LIBMESH_CHKERRABORT(ierr);
#endif
}
开发者ID:dengchangtao,项目名称:libmesh,代码行数:24,代码来源:petsc_matrix.C
示例18: VecSetSizes
int Epetra_PETScAIJMatrix::ExtractDiagonalCopy(Epetra_Vector & Diagonal) const
{
//TODO optimization: only get this diagonal once
Vec petscDiag;
double *vals=0;
int length;
int ierr=VecCreate(Comm_->Comm(),&petscDiag);CHKERRQ(ierr);
VecSetSizes(petscDiag,NumMyRows_,NumGlobalRows_);
# ifdef HAVE_MPI
ierr = VecSetType(petscDiag,VECMPI);CHKERRQ(ierr);
# else //TODO untested!!
VecSetType(petscDiag,VECSEQ);
# endif
MatGetDiagonal(Amat_, petscDiag);
VecGetArray(petscDiag,&vals);
VecGetLocalSize(petscDiag,&length);
for (int i=0; i<length; i++) Diagonal[i] = vals[i];
VecRestoreArray(petscDiag,&vals);
VecDestroy(petscDiag);
return(0);
}
开发者ID:00liujj,项目名称:trilinos,代码行数:24,代码来源:EpetraExt_PETScAIJMatrix.cpp
示例19: MatGetSchurComplement_Basic
/* Developer Notes: This should be implemented with a MatCreate_SchurComplement() as that is the standard design for new Mat classes. */
PetscErrorCode MatGetSchurComplement_Basic(Mat mat,IS isrow0,IS iscol0,IS isrow1,IS iscol1,MatReuse mreuse,Mat *newmat,MatReuse preuse,Mat *newpmat)
{
PetscErrorCode ierr;
Mat A=0,Ap=0,B=0,C=0,D=0;
PetscFunctionBegin;
PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
PetscValidHeaderSpecific(isrow0,IS_CLASSID,2);
PetscValidHeaderSpecific(iscol0,IS_CLASSID,3);
PetscValidHeaderSpecific(isrow1,IS_CLASSID,4);
PetscValidHeaderSpecific(iscol1,IS_CLASSID,5);
if (mreuse == MAT_REUSE_MATRIX) PetscValidHeaderSpecific(*newmat,MAT_CLASSID,7);
if (preuse == MAT_REUSE_MATRIX) PetscValidHeaderSpecific(*newpmat,MAT_CLASSID,9);
PetscValidType(mat,1);
if (mat->factortype) SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_ARG_WRONGSTATE,"Not for factored matrix");
if (mreuse != MAT_IGNORE_MATRIX) {
/* Use MatSchurComplement */
if (mreuse == MAT_REUSE_MATRIX) {
ierr = MatSchurComplementGetSubmatrices(*newmat,&A,&Ap,&B,&C,&D);CHKERRQ(ierr);
if (!A || !Ap || !B || !C) SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_ARG_WRONGSTATE,"Attempting to reuse matrix but Schur complement matrices unset");
if (A != Ap) SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_ARG_WRONGSTATE,"Preconditioning matrix does not match operator");
ierr = MatDestroy(&Ap);CHKERRQ(ierr); /* get rid of extra reference */
}
ierr = MatGetSubMatrix(mat,isrow0,iscol0,mreuse,&A);CHKERRQ(ierr);
ierr = MatGetSubMatrix(mat,isrow0,iscol1,mreuse,&B);CHKERRQ(ierr);
ierr = MatGetSubMatrix(mat,isrow1,iscol0,mreuse,&C);CHKERRQ(ierr);
ierr = MatGetSubMatrix(mat,isrow1,iscol1,mreuse,&D);CHKERRQ(ierr);
switch (mreuse) {
case MAT_INITIAL_MATRIX:
ierr = MatCreateSchurComplement(A,A,B,C,D,newmat);CHKERRQ(ierr);
break;
case MAT_REUSE_MATRIX:
ierr = MatSchurComplementUpdate(*newmat,A,A,B,C,D,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
break;
default:
SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_SUP,"Unrecognized value of mreuse");
}
}
if (preuse != MAT_IGNORE_MATRIX) {
/* Use the diagonal part of A to form D - C inv(diag(A)) B */
Mat Ad,AdB,S;
Vec diag;
PetscInt i,m,n,mstart,mend;
PetscScalar *x;
/* We could compose these with newpmat so that the matrices can be reused. */
if (!A) {ierr = MatGetSubMatrix(mat,isrow0,iscol0,MAT_INITIAL_MATRIX,&A);CHKERRQ(ierr);}
if (!B) {ierr = MatGetSubMatrix(mat,isrow0,iscol1,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr);}
if (!C) {ierr = MatGetSubMatrix(mat,isrow1,iscol0,MAT_INITIAL_MATRIX,&C);CHKERRQ(ierr);}
if (!D) {ierr = MatGetSubMatrix(mat,isrow1,iscol1,MAT_INITIAL_MATRIX,&D);CHKERRQ(ierr);}
ierr = MatGetVecs(A,&diag,PETSC_NULL);CHKERRQ(ierr);
ierr = MatGetDiagonal(A,diag);CHKERRQ(ierr);
ierr = VecReciprocal(diag);CHKERRQ(ierr);
ierr = MatGetLocalSize(A,&m,&n);CHKERRQ(ierr);
/* We need to compute S = D - C inv(diag(A)) B. For row-oriented formats, it is easy to scale the rows of B and
* for column-oriented formats the columns of C can be scaled. Would skip creating a silly diagonal matrix. */
ierr = MatCreate(((PetscObject)A)->comm,&Ad);CHKERRQ(ierr);
ierr = MatSetSizes(Ad,m,n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
ierr = MatSetOptionsPrefix(Ad,((PetscObject)mat)->prefix);CHKERRQ(ierr);
ierr = MatAppendOptionsPrefix(Ad,"diag_");CHKERRQ(ierr);
ierr = MatSetFromOptions(Ad);CHKERRQ(ierr);
ierr = MatSeqAIJSetPreallocation(Ad,1,PETSC_NULL);CHKERRQ(ierr);
ierr = MatMPIAIJSetPreallocation(Ad,1,PETSC_NULL,0,PETSC_NULL);CHKERRQ(ierr);
ierr = MatGetOwnershipRange(Ad,&mstart,&mend);CHKERRQ(ierr);
ierr = VecGetArray(diag,&x);CHKERRQ(ierr);
for (i=mstart; i<mend; i++) {
ierr = MatSetValue(Ad,i,i,x[i-mstart],INSERT_VALUES);CHKERRQ(ierr);
}
ierr = VecRestoreArray(diag,&x);CHKERRQ(ierr);
ierr = MatAssemblyBegin(Ad,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = MatAssemblyEnd(Ad,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
ierr = VecDestroy(&diag);CHKERRQ(ierr);
ierr = MatMatMult(Ad,B,MAT_INITIAL_MATRIX,1,&AdB);CHKERRQ(ierr);
S = (preuse == MAT_REUSE_MATRIX) ? *newpmat : (Mat)0;
ierr = MatMatMult(C,AdB,preuse,PETSC_DEFAULT,&S);CHKERRQ(ierr);
ierr = MatAYPX(S,-1,D,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
*newpmat = S;
ierr = MatDestroy(&Ad);CHKERRQ(ierr);
ierr = MatDestroy(&AdB);CHKERRQ(ierr);
}
ierr = MatDestroy(&A);CHKERRQ(ierr);
ierr = MatDestroy(&B);CHKERRQ(ierr);
ierr = MatDestroy(&C);CHKERRQ(ierr);
ierr = MatDestroy(&D);CHKERRQ(ierr);
PetscFunctionReturn(0);
}
开发者ID:erdc-cm,项目名称:petsc-dev,代码行数:90,代码来源:schurm.c
示例20: main
//.........这里部分代码省略.........
CHKERRQ(ierr);
ierr = MatNorm(sA,NORM_1,&r2);
CHKERRQ(ierr);
rnorm = PetscAbsReal(r1-r2)/r2;
if (rnorm > tol && !rank) {
ierr = PetscPrintf(PETSC_COMM_WORLD,"Error: MatNorm_1(), Anorm=%16.14e, sAnorm=%16.14e bs=%D\n",r1,r2,bs);
CHKERRQ(ierr);
}
/* Test MatGetOwnershipRange() */
ierr = MatGetOwnershipRange(sA,&rstart,&rend);
CHKERRQ(ierr);
ierr = MatGetOwnershipRange(A,&i2,&j2);
CHKERRQ(ierr);
i2 -= rstart;
j2 -= rend;
if (i2 || j2) {
ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d], Error: MaGetOwnershipRange()\n",rank);
CHKERRQ(ierr);
ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
CHKERRQ(ierr);
}
/* Test MatDiagonalScale() */
ierr = MatDiagonalScale(A,x,x);
CHKERRQ(ierr);
ierr = MatDiagonalScale(sA,x,x);
CHKERRQ(ierr);
ierr = MatMultEqual(A,sA,10,&flg);
CHKERRQ(ierr);
if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in MatDiagonalScale");
/* Test MatGetDiagonal(), MatScale() */
ierr = MatGetDiagonal(A,s1);
CHKERRQ(ierr);
ierr = MatGetDiagonal(sA,s2);
CHKERRQ(ierr);
ierr = VecNorm(s1,NORM_1,&r1);
CHKERRQ(ierr);
ierr = VecNorm(s2,NORM_1,&r2);
CHKERRQ(ierr);
r1 -= r2;
if (r1<-tol || r1>tol) {
ierr = PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d], Error: MatDiagonalScale() or MatGetDiagonal(), r1=%g \n",rank,(double)r1);
CHKERRQ(ierr);
ierr = PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
CHKERRQ(ierr);
}
ierr = MatScale(A,alpha);
CHKERRQ(ierr);
ierr = MatScale(sA,alpha);
CHKERRQ(ierr);
/* Test MatGetRowMaxAbs() */
ierr = MatGetRowMaxAbs(A,s1,NULL);
CHKERRQ(ierr);
ierr = MatGetRowMaxAbs(sA,s2,NULL);
CHKERRQ(ierr);
ierr = VecNorm(s1,NORM_1,&r1);
CHKERRQ(ierr);
ierr = VecNorm(s2,NORM_1,&r2);
CHKERRQ(ierr);
r1 -= r2;
if (r1<-tol || r1>tol) {
开发者ID:petsc,项目名称:petsc,代码行数:67,代码来源:ex75.c
注:本文中的MatGetDiagonal函数示例由纯净天空整理自Github/MSDocs等源码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。 |
请发表评论