• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    公众号
登陆 注册

C++ KSPSolve函数代码示例

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

本文整理汇总了C++中KSPSolve函数的典型用法代码示例。如果您正苦于以下问题：C++ KSPSolve函数的具体用法？C++ KSPSolve怎么用？C++ KSPSolve使用的例子？那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。



在下文中一共展示了KSPSolve函数的20个代码示例，这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞，您的评价将有助于我们的系统推荐出更棒的C++代码示例。

示例1: KSPSolve_SpecEst

static PetscErrorCode  KSPSolve_SpecEst(KSP ksp)
{
  PetscErrorCode ierr;
  KSP_SpecEst    *spec = (KSP_SpecEst*)ksp->data;

  PetscFunctionBegin;
  if (spec->current) {
    ierr = KSPSolve(spec->kspcheap,ksp->vec_rhs,ksp->vec_sol);CHKERRQ(ierr);
    ierr = KSPSpecEstPropagateUp(ksp,spec->kspcheap);CHKERRQ(ierr);
  } else {
    PetscInt  i,its,neig;
    PetscReal *real,*imag,rad = 0;
    ierr = KSPSolve(spec->kspest,ksp->vec_rhs,ksp->vec_sol);CHKERRQ(ierr);
    ierr = KSPSpecEstPropagateUp(ksp,spec->kspest);CHKERRQ(ierr);
    ierr = KSPComputeExtremeSingularValues(spec->kspest,&spec->max,&spec->min);CHKERRQ(ierr);

    ierr = KSPGetIterationNumber(spec->kspest,&its);CHKERRQ(ierr);
    ierr = PetscMalloc2(its,PetscReal,&real,its,PetscReal,&imag);CHKERRQ(ierr);
    ierr = KSPComputeEigenvalues(spec->kspest,its,real,imag,&neig);CHKERRQ(ierr);
    for (i=0; i<neig; i++) {
      /* We would really like to compute w (nominally 1/radius) to minimize |1-wB|.  Empirically it
         is better to compute rad = |1-B| than rad = |B|.  There must be a cheap way to do better. */
      rad = PetscMax(rad,PetscRealPart(PetscSqrtScalar((PetscScalar)(PetscSqr(real[i]-1.) + PetscSqr(imag[i])))));
    }
    ierr = PetscFree2(real,imag);CHKERRQ(ierr);
    spec->radius = rad;

    ierr = KSPChebyshevSetEigenvalues(spec->kspcheap,spec->max*spec->maxfactor,spec->min*spec->minfactor);CHKERRQ(ierr);
    ierr = KSPRichardsonSetScale(spec->kspcheap,spec->richfactor/spec->radius);
    ierr = PetscInfo3(ksp,"Estimated singular value min=%G max=%G, spectral radius=%G",spec->min,spec->max,spec->radius);CHKERRQ(ierr);
    spec->current = PETSC_TRUE;
  }
  PetscFunctionReturn(0);
}
开发者ID:Kun-Qu，项目名称:petsc，代码行数:34，代码来源:specest.c


示例2: PCMGKCycle_Private

PetscErrorCode PCMGKCycle_Private(PC pc,PC_MG_Levels **mglevels)
{
  PetscErrorCode ierr;
  PetscInt       i,l = mglevels[0]->levels;

  PetscFunctionBegin;
  /* restrict the RHS through all levels to coarsest. */
  for (i=l-1; i>0; i--){
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
  }
  
  /* work our way up through the levels */
  ierr = VecSet(mglevels[0]->x,0.0);CHKERRQ(ierr); 
  for (i=0; i<l-1; i++) {
    if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
    ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr);
    if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
    if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatInterpolate(mglevels[i+1]->interpolate,mglevels[i]->x,mglevels[i+1]->x);CHKERRQ(ierr);
    if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
  }
  if (mglevels[l-1]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[l-1]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
  ierr = KSPSolve(mglevels[l-1]->smoothd,mglevels[l-1]->b,mglevels[l-1]->x);CHKERRQ(ierr);
  if (mglevels[l-1]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[l-1]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}

  PetscFunctionReturn(0);
}
开发者ID:Kun-Qu，项目名称:petsc，代码行数:29，代码来源:fmg.c


示例3: G

/*
S^{-1} = ( G^T G )^{-1} G^T K G ( G^T G )^{-1}
       = A C A
S^{-T} = A^T (A C)^T
       = A^T C^T A^T, but A = G^T G which is symmetric
       = A C^T A
       = A G^T ( G^T K )^T A
       = A G^T K^T G A

*/
PetscErrorCode BSSCR_PCApplyTranspose_GtKG( PC pc, Vec x, Vec y )
{
	
	PC_GtKG ctx = (PC_GtKG)pc->data;
	
	KSP ksp;
	Mat K, G;
	Vec s,t,X;
	PetscLogDouble t0,t1;
	
	ksp = ctx->ksp;
	K = ctx->K;
	G = ctx->G;
	s = ctx->s;
	t = ctx->t;
	X = ctx->X;
	
	if (ctx->monitor_rhs_consistency) {		BSSCRBSSCR_Lp_monitor_check_rhs_consistency(ksp,x,1);		}
	PetscGetTime(&t0);
	KSPSolve( ksp, x, t ); /* t <- GtG_inv x */
	PetscGetTime(&t1);
	if (ctx->monitor_activated) {	BSSCR_PCBFBTSubKSPMonitor(ksp,1,(t1-t0));		}
	
	MatMult( G, t, s ); /* s <- G t */
	MatMultTranspose( K, s, X ); /* X <- K^T s */
	MatMultTranspose( G, X, t ); /* t <- Gt X */
	
	if (ctx->monitor_rhs_consistency) {		BSSCRBSSCR_Lp_monitor_check_rhs_consistency(ksp,t,2);		}
	PetscGetTime(&t0);
	KSPSolve( ksp, t, y ); /* y <- GtG_inv t */
	PetscGetTime(&t1);
	if (ctx->monitor_activated) {	BSSCR_PCBFBTSubKSPMonitor(ksp,2,(t1-t0));		}
	
	PetscFunctionReturn(0);
}
开发者ID:AngelValverdePerez，项目名称:underworld2，代码行数:45，代码来源:pc_GtKG.c


示例4: BSSCR_BSSCR_PCApplyTranspose_GtKG_diagonal_scaling

PetscErrorCode BSSCR_BSSCR_PCApplyTranspose_GtKG_diagonal_scaling( PC pc, Vec x, Vec y )
{
	PC_GtKG ctx = (PC_GtKG)pc->data;
	
	KSP ksp;
	Mat K, G;
	Vec s,t,X,di;
	
	
	ksp = ctx->ksp;
	K = ctx->K;
	G = ctx->G;
	di = ctx->inv_diag_M;
	s = ctx->s;
	t = ctx->t;
	X = ctx->X;
	
	if (ctx->monitor_rhs_consistency) {		BSSCRBSSCR_Lp_monitor_check_rhs_consistency(ksp,x,1);		}
	KSPSolve( ksp, x, t );       /* t <- GtG_inv x */
	if (ctx->monitor_activated) {	BSSCR_Lp_monitor(ksp,1);		}
	
	MatMult( G, t, s );          /* s <- G t */
	VecPointwiseMult( s, s,di ); /* s <- s * di */
	
	MatMultTranspose( K, s, X ); /* X <- K^T s */
	VecPointwiseMult( X, X,di ); /* X <- X * di */
	
	MatMultTranspose( G, X, t ); /* t <- Gt X */
	
	if (ctx->monitor_rhs_consistency) {		BSSCRBSSCR_Lp_monitor_check_rhs_consistency(ksp,t,2);		}
	KSPSolve( ksp, t, y );       /* y <- GtG_inv t */
	if (ctx->monitor_activated) {	BSSCR_Lp_monitor(ksp,2);		}
	
	PetscFunctionReturn(0);
}
开发者ID:AngelValverdePerez，项目名称:underworld2，代码行数:35，代码来源:pc_GtKG.c


示例5: PressurePoisson

PetscErrorCode PressurePoisson( UserContext *uc )
{
  Mat A;
  Vec b, px, py, p, u, v, ss, c;
  KSP ksp;
  PetscErrorCode ierr;
  
  
  PetscFunctionBegin;
  PetscLogEventBegin(EVENT_PressurePoisson,0,0,0,0);
  PetscLogEventRegister(&EVENT_PressurePoisson,"PressurePoisson", 0);
  
  /*  Assemble and Solve for pressure */
ierr = PetscPrintf(PETSC_COMM_WORLD, "***********************\nPRESSURE\n"); CHKERRQ(ierr);
  ierr = AssemblePressureMatrx( uc ); CHKERRQ(ierr);
  ierr = AssemblePressureRHS( uc ); CHKERRQ(ierr);
  ierr = SolvePressure( uc ); CHKERRQ(ierr);
  ierr = WriteResult( uc, uc->p, "p_pressure" ); CHKERRQ(ierr);
    
  /*  Assemble and Solve for Velocity  */
ierr = PetscPrintf(PETSC_COMM_WORLD, "***********************\nVELOCITY\n"); CHKERRQ(ierr);
  ierr = AssembleVelocityRHS( uc ); CHKERRQ(ierr);
  ierr = AssembleVelocityMatrix( uc ); CHKERRQ(ierr);
  
  ierr = KSPSetOperators(uc->ksp,uc->A, uc->A, SAME_PRECONDITIONER); CHKERRQ(ierr);
  ierr = KSPSolve(uc->ksp,uc->px,uc->u);CHKERRQ(ierr);
  ierr = KSPSolve(uc->ksp,uc->py,uc->v);CHKERRQ(ierr);
  ierr = WriteResult( uc, uc->u, "u_velocity" ); CHKERRQ(ierr);
  ierr = WriteResult( uc, uc->v, "v_velocity" ); CHKERRQ(ierr);

  ierr = ComputeShearStress(uc); CHKERRQ(ierr);
  ierr = WriteResult( uc, uc->ss,"shear_stress" ); CHKERRQ(ierr);
  ierr = ConservationTest(uc); CHKERRQ(ierr);
  
  /*  Output indexing  */
  PetscReal *idx;
  VecGetArray(uc->b,&idx);
  for( int i = 0; i < uc->numNodes; i++)
    idx[i] = i;
  VecRestoreArray(uc->b,&idx);
  WriteResult(uc, uc->b, "indexes");
  VecSet(uc->b, 0.);
  
  ierr = VecDestroy(uc->c); CHKERRQ(ierr);
  ierr = VecDestroy(uc->ss); CHKERRQ(ierr);
  ierr = MatDestroy(uc->A); CHKERRQ(ierr);
  ierr = KSPDestroy(uc->ksp); CHKERRQ(ierr);
  ierr = VecDestroy(uc->b); CHKERRQ(ierr);
  ierr = VecDestroy(uc->p); CHKERRQ(ierr);
  ierr = VecDestroy(uc->px); CHKERRQ(ierr);
  ierr = VecDestroy(uc->py); CHKERRQ(ierr);
  ierr = VecDestroy(uc->u); CHKERRQ(ierr);
  ierr = VecDestroy(uc->v); CHKERRQ(ierr);
  PetscLogEventEnd(EVENT_PressurePoisson,0,0,0,0);
  PetscFunctionReturn(0);
}
开发者ID:adrielb，项目名称:DCell，代码行数:56，代码来源:PressurePoisson.c


示例6: PCMGMCycle_Private

PetscErrorCode PCMGMCycle_Private(PC pc,PC_MG_Levels **mglevelsin,PCRichardsonConvergedReason *reason)
{
  PC_MG          *mg = (PC_MG*)pc->data;
  PC_MG_Levels   *mgc,*mglevels = *mglevelsin;
  PetscErrorCode ierr;
  PetscInt       cycles = (mglevels->level == 1) ? 1 : (PetscInt) mglevels->cycles;
  PC             subpc;
  PCFailedReason pcreason;

  PetscFunctionBegin;
  if (mglevels->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
  ierr = KSPSolve(mglevels->smoothd,mglevels->b,mglevels->x);CHKERRQ(ierr);  /* pre-smooth */
  ierr = KSPGetPC(mglevels->smoothd,&subpc);CHKERRQ(ierr);
  ierr = PCGetSetUpFailedReason(subpc,&pcreason);CHKERRQ(ierr); 
  if (pcreason) {
    pc->failedreason = PC_SUBPC_ERROR;
  }
  if (mglevels->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
  if (mglevels->level) {  /* not the coarsest grid */
    if (mglevels->eventresidual) {ierr = PetscLogEventBegin(mglevels->eventresidual,0,0,0,0);CHKERRQ(ierr);}
    ierr = (*mglevels->residual)(mglevels->A,mglevels->b,mglevels->x,mglevels->r);CHKERRQ(ierr);
    if (mglevels->eventresidual) {ierr = PetscLogEventEnd(mglevels->eventresidual,0,0,0,0);CHKERRQ(ierr);}

    /* if on finest level and have convergence criteria set */
    if (mglevels->level == mglevels->levels-1 && mg->ttol && reason) {
      PetscReal rnorm;
      ierr = VecNorm(mglevels->r,NORM_2,&rnorm);CHKERRQ(ierr);
      if (rnorm <= mg->ttol) {
        if (rnorm < mg->abstol) {
          *reason = PCRICHARDSON_CONVERGED_ATOL;
          ierr    = PetscInfo2(pc,"Linear solver has converged. Residual norm %g is less than absolute tolerance %g\n",(double)rnorm,(double)mg->abstol);CHKERRQ(ierr);
        } else {
          *reason = PCRICHARDSON_CONVERGED_RTOL;
          ierr    = PetscInfo2(pc,"Linear solver has converged. Residual norm %g is less than relative tolerance times initial residual norm %g\n",(double)rnorm,(double)mg->ttol);CHKERRQ(ierr);
        }
        PetscFunctionReturn(0);
      }
    }

    mgc = *(mglevelsin - 1);
    if (mglevels->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatRestrict(mglevels->restrct,mglevels->r,mgc->b);CHKERRQ(ierr);
    if (mglevels->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = VecSet(mgc->x,0.0);CHKERRQ(ierr);
    while (cycles--) {
      ierr = PCMGMCycle_Private(pc,mglevelsin-1,reason);CHKERRQ(ierr);
    }
    if (mglevels->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatInterpolateAdd(mglevels->interpolate,mgc->x,mglevels->x,mglevels->x);CHKERRQ(ierr);
    if (mglevels->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    if (mglevels->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
    ierr = KSPSolve(mglevels->smoothu,mglevels->b,mglevels->x);CHKERRQ(ierr);    /* post smooth */
    if (mglevels->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
  }
  PetscFunctionReturn(0);
}
开发者ID:ziolai，项目名称:petsc，代码行数:56，代码来源:mg.c


示例7: PCApply_NN

static PetscErrorCode PCApply_NN(PC pc,Vec r,Vec z)
{
  PC_IS          *pcis = (PC_IS*)(pc->data);
  PetscErrorCode ierr;
  PetscScalar    m_one = -1.0;
  Vec            w = pcis->vec1_global;

  PetscFunctionBegin;
  /*
    Dirichlet solvers.
    Solving $ B_I^{(i)}r_I^{(i)} $ at each processor.
    Storing the local results at vec2_D
  */
  ierr = VecScatterBegin(pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = VecScatterEnd  (pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = KSPSolve(pcis->ksp_D,pcis->vec1_D,pcis->vec2_D);CHKERRQ(ierr);
  
  /*
    Computing $ r_B - \sum_j \tilde R_j^T A_{BI}^{(j)} (B_I^{(j)}r_I^{(j)}) $ .
    Storing the result in the interface portion of the global vector w.
  */
  ierr = MatMult(pcis->A_BI,pcis->vec2_D,pcis->vec1_B);CHKERRQ(ierr);
  ierr = VecScale(pcis->vec1_B,m_one);CHKERRQ(ierr);
  ierr = VecCopy(r,w);CHKERRQ(ierr);
  ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,w,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  ierr = VecScatterEnd  (pcis->global_to_B,pcis->vec1_B,w,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);

  /*
    Apply the interface preconditioner
  */
  ierr = PCNNApplyInterfacePreconditioner(pc,w,z,pcis->work_N,pcis->vec1_B,pcis->vec2_B,pcis->vec3_B,pcis->vec1_D,
                                          pcis->vec3_D,pcis->vec1_N,pcis->vec2_N);CHKERRQ(ierr);

  /*
    Computing $ t_I^{(i)} = A_{IB}^{(i)} \tilde R_i z_B $
    The result is stored in vec1_D.
  */
  ierr = VecScatterBegin(pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = VecScatterEnd  (pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = MatMult(pcis->A_IB,pcis->vec1_B,pcis->vec1_D);CHKERRQ(ierr);

  /*
    Dirichlet solvers.
    Computing $ B_I^{(i)}t_I^{(i)} $ and sticking into the global vector the blocks
    $ B_I^{(i)}r_I^{(i)} - B_I^{(i)}t_I^{(i)} $.
  */
  ierr = VecScatterBegin(pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  ierr = VecScatterEnd  (pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  ierr = KSPSolve(pcis->ksp_D,pcis->vec1_D,pcis->vec2_D);CHKERRQ(ierr);
  ierr = VecScale(pcis->vec2_D,m_one);CHKERRQ(ierr);
  ierr = VecScatterBegin(pcis->global_to_D,pcis->vec2_D,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  ierr = VecScatterEnd  (pcis->global_to_D,pcis->vec2_D,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}
开发者ID:Kun-Qu，项目名称:petsc，代码行数:54，代码来源:nn.c


示例8: PCMGACycle_Private

PetscErrorCode PCMGACycle_Private(PC pc,PC_MG_Levels **mglevels)
{
  PetscErrorCode ierr;
  PetscInt       i,l = mglevels[0]->levels;

  PetscFunctionBegin;
  /* compute RHS on each level */
  for (i=l-1; i>0; i--) {
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr);
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
  }
  /* solve separately on each level */
  for (i=0; i<l; i++) {
    ierr = VecSet(mglevels[i]->x,0.0);CHKERRQ(ierr);
    if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
    ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr);
    ierr = KSPCheckSolve(mglevels[i]->smoothd,pc,mglevels[i]->x);CHKERRQ(ierr);
    if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);}
  }
  for (i=1; i<l; i++) {
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
    ierr = MatInterpolateAdd(mglevels[i]->interpolate,mglevels[i-1]->x,mglevels[i]->x,mglevels[i]->x);CHKERRQ(ierr);
    if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);}
  }
  PetscFunctionReturn(0);
}
开发者ID:firedrakeproject，项目名称:petsc，代码行数:27，代码来源:smg.c


示例9: VecScale

bool PotentialSolve::Solve(Vec delta, Vec pot, double bias) {
    PetscInt its;
    KSPConvergedReason reason;
    bool retval;

    // Delta to -Delta
    VecScale(delta, -1.0/bias);

    // Actually solve
    KSPSolve(solver, delta, pot);
    KSPGetConvergedReason(solver,&reason);
    if (reason<0) {
        PetscPrintf(PETSC_COMM_WORLD,"Diverged : %d.\n",reason);
        retval=false;
    } else {
        KSPGetIterationNumber(solver,&its);
        PetscPrintf(PETSC_COMM_WORLD,"\nConvergence in %d iterations.\n",(int)its);
        retval=true;
    }


    // Remove any values that might have crept in here
    _dg1.ZeroPad(pot);

    // Clean up
    VecScale(delta, -1.0*bias);

    return retval;
}
开发者ID:ajcuesta，项目名称:baorecon，代码行数:29，代码来源:PotentialSolve.cpp


示例10: swap

void SingleLongPipe::applyAdvection() {
  swap(rhs, concentration);
  PetscScalar* localRhs = rhs.getLocalArray();
  PetscScalar alpha = flux*dt*T0/M_PI/(L*L*L);
  const PetscScalar *localRadii = radii.getLocalArray();
  PetscInt mStart, mEnd, jStart = 0;
  MatGetOwnershipRange(M, &mStart, &mEnd);
  if (rank == 0) {
    int row = 0;
    int column = 0;
    PetscScalar val = alpha / (localRadii[0] * localRadii[0]);
    PetscScalar val1 = 1 + val;
    MatSetValues(M, 1, &row, 1, &column, &val1, INSERT_VALUES);
    localRhs[0] += val*C0;
    ++mStart;
    ++jStart;
  }
  for (int i = mStart, j = jStart; i < mEnd; ++i, ++j) {
    int row = i;
    int column = i;
    PetscScalar val = alpha / (localRadii[j] * localRadii[j]);
    PetscScalar val1 = 1 + val;
    MatSetValues(M, 1, &row, 1, &column, &val1, INSERT_VALUES);
    val1 = -val;
    --column;
    MatSetValues(M, 1, &row, 1, &column, &val1, INSERT_VALUES);
  }
  MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY);
  MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY);
  rhs.restoreVec(localRhs);
  KSPSolve(ksp, rhs.getVec(), concentration.getVec());
}
开发者ID:whyzidane，项目名称:CO2_project，代码行数:32，代码来源:single_long_pipe_class.cpp


示例11: SFieldSolveFor

double SFieldSolveFor(SField sfv, double *Y, unsigned int yCount) {
    mySField sf = static_cast<mySField>(sfv);
    assert(yCount <= sf->maxN);
    assert(Y);
    assert(sf->running);
    sf->Y = Y;
    sf->curN = yCount;

    // -------------- SOLVE
    PetscErrorCode ierr;
    PetscLogDouble tic,toc;
    PetscTime(&tic);
    int pt[sf->d];

    ierr = MatZeroEntries(sf->J); CHKERRQ(ierr);
    JacobianOnD(sf->J, sf->F, 0, pt, sf);

    ierr = MatAssemblyBegin(sf->J,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
    ierr = MatAssemblyEnd(sf->J,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);

    PetscTime(&toc);
    sf->timeAssembly += toc-tic;
    PetscTime(&tic);
    ierr = VecZeroEntries(sf->U);  CHKERRQ(ierr);
    ierr = KSPSetOperators(sf->ksp, sf->J, sf->J); CHKERRQ(ierr);
    ierr = KSPSetUp(sf->ksp);  CHKERRQ(ierr);
    ierr = KSPSolve(sf->ksp,sf->F,sf->U); CHKERRQ(ierr);
    PetscTime(&toc);
    sf->timeSolver  += toc-tic;
    return Integrate(sf->U,pt,0,sf);
}
开发者ID:StochasticNumerics，项目名称:mimclib，代码行数:31，代码来源:matern.cpp


示例12: main

int main(int argc,char **argv)
{
  PetscErrorCode ierr;
  KSP            ksp;
  DM             da,shell;
  PetscInt       levels;

  PetscInitialize(&argc,&argv,(char*)0,help);

  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
  ierr = DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,-129,1,1,0,&da);CHKERRQ(ierr);
  ierr = MyDMShellCreate(PETSC_COMM_WORLD,da,&shell);CHKERRQ(ierr);
  /* these two lines are not needed but allow PCMG to automatically know how many multigrid levels the user wants */
  ierr = DMGetRefineLevel(da,&levels);CHKERRQ(ierr);
  ierr = DMSetRefineLevel(shell,levels);CHKERRQ(ierr);

  ierr = KSPSetDM(ksp,shell);CHKERRQ(ierr);
  ierr = KSPSetComputeRHS(ksp,ComputeRHS,NULL);CHKERRQ(ierr);
  ierr = KSPSetComputeOperators(ksp,ComputeMatrix,NULL);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
  ierr = KSPSolve(ksp,NULL,NULL);CHKERRQ(ierr);

  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = DMDestroy(&shell);CHKERRQ(ierr);
  ierr = DMDestroy(&da);CHKERRQ(ierr);
  ierr = PetscFinalize();

  return 0;
}
开发者ID:000Justin000，项目名称:ATPESC，代码行数:29，代码来源:ex65.c


示例13: MatAssemblyBegin

void PETSc::Solve(void)
{
        
  //start_clock("Before Assemble matrix and vector");
	ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
  	ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
  	
  	ierr = VecAssemblyBegin(x);
  	ierr = VecAssemblyEnd(x);
  	
  	ierr = VecAssemblyBegin(b);
  	ierr = VecAssemblyEnd(b);
	//stop_clock("After Assembly matrix and vector");


        KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);
        KSPSetType(ksp,KSPBCGSL);
	KSPBCGSLSetEll(ksp,2);

	//KSPGetPC(ksp, &pc);
	//PCSetType(pc, PCJACOBI);

        KSPSetFromOptions(ksp);
        KSPSetUp(ksp);

	//start_clock("Before KSPSolve");
        KSPSolve(ksp,b,x);
	//stop_clock("After KSPSolve");
}
开发者ID:tonyguo1，项目名称:LSGFD，代码行数:29，代码来源:solver_petsc.cpp


示例14: SolvePressure

PetscErrorCode SolvePressure(UserContext *uc)
{
  PetscErrorCode ierr;
  
  PetscFunctionBegin;
  PetscLogEventBegin(EVENT_SolvePressure,0,0,0,0);

  ierr = KSPCreate(PETSC_COMM_WORLD, &uc->ksp); CHKERRQ(ierr);
  ierr = KSPSetType(uc->ksp,KSPCG); CHKERRQ(ierr);
  ierr = KSPSetOperators(uc->ksp,uc->A,uc->A,SAME_NONZERO_PATTERN); CHKERRQ(ierr);
  ierr = KSPSetFromOptions(uc->ksp); CHKERRQ(ierr);
  ierr = KSPSetTolerances(uc->ksp, 0.000001, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT); CHKERRQ(ierr);
  
  ierr = KSPSetType( uc->ksp, KSPPREONLY); CHKERRQ(ierr);
  PC pc;
  ierr = KSPGetPC(uc->ksp, &pc); CHKERRQ(ierr);
  ierr = PCSetType(pc, PCLU); CHKERRQ(ierr);
  
  ierr = KSPSolve(uc->ksp,uc->b,uc->p);CHKERRQ(ierr);
  
  KSPConvergedReason reason;
  KSPGetConvergedReason(uc->ksp,&reason);
  PetscPrintf(PETSC_COMM_WORLD,"Pressure KSPConvergedReason: %D\n", reason);
  // if( reason < 0 )    SETERRQ(PETSC_ERR_CONV_FAILED, "Exiting: Failed to converge\n");
  
  PetscLogEventEnd(EVENT_SolvePressure,0,0,0,0);
  PetscFunctionReturn(0);
}
开发者ID:adrielb，项目名称:DCell，代码行数:28，代码来源:PressurePoisson.c


示例15: main

int main(int argc, char **args) {
  PetscErrorCode ierr;
  ierr = DCellInit(); CHKERRQ(ierr);

  PetscReal dx = 1;
  iCoor size = {1625,1145,0};
//  iCoor size = {253,341,0};
  int fd;
  Coor dh = {dx,dx,dx};
  iCoor pos = {0,0,0};
  Grid chip;
  ierr = GridCreate(dh,pos,size,1,&chip); CHKERRQ(ierr);
  ierr = PetscInfo(0,"Reading image file\n"); CHKERRQ(ierr);
  ierr = PetscBinaryOpen("/scratch/n/BL Big",FILE_MODE_READ,&fd); CHKERRQ(ierr);
  ierr = PetscBinaryRead(fd,chip->v1,size.x*size.y,PETSC_DOUBLE); CHKERRQ(ierr);
  ierr = PetscBinaryClose(fd); CHKERRQ(ierr);
  ierr = PetscInfo(0,"Writing image file\n"); CHKERRQ(ierr);
  ierr = GridWrite(chip,0); CHKERRQ(ierr);

  FluidField fluid;
  ierr = FluidFieldCreate(PETSC_COMM_WORLD, &fluid);  CHKERRQ(ierr);
  ierr = FluidFieldSetDims(fluid,size); CHKERRQ(ierr);
  ierr = FluidFieldSetDx(fluid,dx); CHKERRQ(ierr);
  ierr = FluidFieldSetMask(fluid, chip); CHKERRQ(ierr);
  ierr = FluidFieldSetup(fluid); CHKERRQ(ierr);
  ierr = SetPressureBC(fluid); CHKERRQ(ierr);
  ierr = KSPSolve(fluid->ksp,fluid->rhs,fluid->vel); CHKERRQ(ierr);
  ierr = FluidFieldWrite( fluid,0); CHKERRQ(ierr);

  ierr = FluidFieldDestroy(fluid); CHKERRQ(ierr);
  ierr = GridDestroy(chip); CHKERRQ(ierr);
  ierr = DCellFinalize(); CHKERRQ(ierr);
  return 0;
}
开发者ID:adrielb，项目名称:DCell，代码行数:34，代码来源:Microfluidics.c


示例16: main

int main(int argc,char **argv)
{
  KSP            ksp;
  DM             da;
  UserContext    user;
  PetscInt       bc;
  PetscErrorCode ierr;

  PetscInitialize(&argc,&argv,(char *)0,help);
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
  ierr = DMDACreate2d(PETSC_COMM_WORLD, DMDA_BOUNDARY_NONE, DMDA_BOUNDARY_NONE,DMDA_STENCIL_STAR,-11,-11,PETSC_DECIDE,PETSC_DECIDE,1,1,PETSC_NULL,PETSC_NULL,&da);CHKERRQ(ierr);
  ierr = KSPSetDM(ksp,(DM)da);
  ierr = DMSetApplicationContext(da,&user);CHKERRQ(ierr);

  user.uu = 1.0;
  user.tt = 1.0;
  bc   = (PetscInt)NEUMANN; // Use Neumann Boundary Conditions
  user.bcType = (BCType)bc;


  ierr = KSPSetComputeRHS(ksp,ComputeRHS,&user);CHKERRQ(ierr);
  ierr = KSPSetComputeOperators(ksp,ComputeJacobian,&user);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
  ierr = KSPSolve(ksp,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);

  ierr = DMDestroy(&da);CHKERRQ(ierr);
  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = PetscFinalize();CHKERRQ(ierr);
  return 0;
}
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:30，代码来源:ex50.c


示例17: START_TEST

END_TEST

START_TEST( Generate2DLapacian_test )
{
  PetscErrorCode ierr;
  PetscInt d = 50, len = d*d;
  Mat m;
  PetscViewer view;
  Generate2DLapacian( d, d, &m);
  
  PetscViewerASCIIOpen(PETSC_COMM_SELF, "mat.dat", &view);
  MatView(m, view);
  
  KSP ksp;
  PC pc;
  ierr = KSPCreate(PETSC_COMM_SELF, &ksp); CHKERRQ(ierr);
  KSPSetType(ksp, KSPPREONLY);
  KSPGetPC(ksp, &pc);
  PCSetType(pc, PCCHOLESKY);
  KSPSetOperators(ksp, m, m, SAME_PRECONDITIONER);
  Vec v, s;
  VecCreateSeq(PETSC_COMM_SELF,len,&v);
  VecDuplicate(v, &s);
  VecSet(v, 1);
  ierr = KSPSolve(ksp, v, s); CHKERRQ(ierr);
  
}
开发者ID:adrielb，项目名称:DCell，代码行数:27，代码来源:check.c


示例18: _check_if_constant_nullsp_present

Bool _check_if_constant_nullsp_present( Stokes_SLE_UzawaSolver* self, Mat K, Mat G, Mat M, Vec t1, Vec ustar,  Vec r, Vec l, KSP ksp )
{
    PetscInt N;
    PetscScalar sum;
    PetscReal nrm;
    Bool nullsp_present;

    VecGetSize(l,&N);
    sum  = 1.0/N;
    VecSet(l,sum);

    /* [S] {l} = {r} */
    MatMult( G,l, t1 );
    KSPSolve( ksp, t1, ustar );
    MatMultTranspose( G, ustar, r );     
    if ( M ) {
	VecScale( r, -1.0 );
	MatMultAdd( M,l, r, r );
	VecScale( r, -1.0 );
    }

    VecNorm(r,NORM_2,&nrm);
    if (nrm < 1.e-7) {
	Journal_PrintfL( self->info, 1, "Constant null space detected, " ); 
	nullsp_present = True;
    }
    else {
	Journal_PrintfL( self->info, 1, "Constant null space not present, " );
	nullsp_present = False;
    }
    Journal_PrintfL( self->info, 1, "|| [S]{1} || = %G\n", nrm );


    return nullsp_present;
}
开发者ID:OlympusMonds，项目名称:EarthByte_Underworld，代码行数:35，代码来源:Stokes_SLE_UzawaSolver.c


示例19: _Stokes_SLE_UzawaSolver_GetRhs

PetscErrorCode _Stokes_SLE_UzawaSolver_GetRhs( void *stokesSLE, void *solver, Vec rhs )
{
        Stokes_SLE_UzawaSolver *self = (Stokes_SLE_UzawaSolver*)solver;
	Stokes_SLE *sle = (Stokes_SLE*)stokesSLE;
	/* stg linear algebra */
	KSP     A11_solver = self->velSolver;
        Mat         A12 = sle->gStiffMat->matrix;
        Vec         b1  = sle->fForceVec->vector;
        Vec         b2  = sle->hForceVec->vector;
        Vec         u_star = self->vStarVec;
	
	/* petsc variables */
	KSP ksp_A11;

	/* check operations will be valid */
	if (sle->dStiffMat!=NULL) {   Stg_SETERRQ( PETSC_ERR_SUP, "A21 must be NULL" ); }
	if (A11_solver==NULL){    Stg_SETERRQ( PETSC_ERR_ARG_NULL, "vel_solver is NULL" ); }
        if (A12==NULL){           Stg_SETERRQ( PETSC_ERR_ARG_NULL, "A12 is NULL" ); }
        if (b1==NULL){       Stg_SETERRQ( PETSC_ERR_ARG_NULL, "b1 is NULL" ); }
        if (b2==NULL){       Stg_SETERRQ( PETSC_ERR_ARG_NULL, "b2 is NULL" ); }
        if (u_star==NULL){   Stg_SETERRQ( PETSC_ERR_ARG_NULL, "u* is NULL" ); }

	/* Extract petsc objects */
	ksp_A11 = A11_solver;

	/* compute rhs = A12^T A11^{-1} b1 - b2 */
	KSPSolve( ksp_A11, b1, u_star );         /* u* = A11^{-1} b1 */
	MatMultTranspose( A12, u_star, rhs );    /* b2 = A12^T u* */
	VecAXPY( rhs, -1.0, b2 );  /* rhs <- rhs - b2 */
}
开发者ID:OlympusMonds，项目名称:EarthByte_Underworld，代码行数:30，代码来源:Stokes_SLE_UzawaSolver.c


示例20: PCApply_Redistribute

static PetscErrorCode PCApply_Redistribute(PC pc,Vec b,Vec x)
{
  PC_Redistribute   *red = (PC_Redistribute*)pc->data;
  PetscErrorCode    ierr;
  PetscInt          dcnt = red->dcnt,i;
  const PetscInt    *drows = red->drows;
  PetscScalar       *xwork;
  const PetscScalar *bwork,*diag = red->diag;

  PetscFunctionBegin;
  if (!red->work) {
    ierr = VecDuplicate(b,&red->work);CHKERRQ(ierr);
  }
  /* compute the rows of solution that have diagonal entries only */
  ierr = VecSet(x,0.0);CHKERRQ(ierr);         /* x = diag(A)^{-1} b */
  ierr = VecGetArray(x,&xwork);CHKERRQ(ierr);
  ierr = VecGetArrayRead(b,&bwork);CHKERRQ(ierr);
  for (i=0; i<dcnt; i++) {
    xwork[drows[i]] = diag[i]*bwork[drows[i]];
  }
  ierr = PetscLogFlops(dcnt);CHKERRQ(ierr);
  ierr = VecRestoreArray(red->work,&xwork);CHKERRQ(ierr);
  ierr = VecRestoreArrayRead(b,&bwork);CHKERRQ(ierr);
  /* update the right hand side for the reduced system with diagonal rows (and corresponding columns) removed */
  ierr = MatMult(pc->pmat,x,red->work);CHKERRQ(ierr);
  ierr = VecAYPX(red->work,-1.0,b);CHKERRQ(ierr);   /* red->work = b - A x */

  ierr = VecScatterBegin(red->scatter,red->work,red->b,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = VecScatterEnd(red->scatter,red->work,red->b,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
  ierr = KSPSolve(red->ksp,red->b,red->x);CHKERRQ(ierr);
  ierr = VecScatterBegin(red->scatter,red->x,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  ierr = VecScatterEnd(red->scatter,red->x,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:34，代码来源:redistribute.c



注：本文中的KSPSolve函数示例由纯净天空整理自Github/MSDocs等源码及文档管理平台，相关代码片段筛选自各路编程大神贡献的开源项目，源码版权归原作者所有，传播和使用请参考对应项目的License；未经允许，请勿转载。

鲜花
握手
雷人
路过
鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
More+
上一篇：
C++ KSP_MatMult函数代码示例发布时间：2022-05-30
下一篇：
C++ KSPSetUp函数代码示例发布时间：2022-05-30
热门推荐
More+
热门话题
More+
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服（服务时间 9:00～18:00）

在线QQ客服
地址：深圳市南山区西丽大学城创智工业园
电邮：jeky_zhao#qq.com
移动电话：139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap