C++ KSPSetOperators函数代码示例

OStack程序员社区-中国程序员成长平台 › 门户 › 编程› C++›C++教程

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

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

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

示例1: StokesSetupPC

PetscErrorCode StokesSetupPC(Stokes *s, KSP ksp)
{
  KSP            *subksp;
  PC             pc;
  PetscInt       n = 1;
  PetscErrorCode ierr;

  PetscFunctionBeginUser;
  ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
  ierr = PCFieldSplitSetIS(pc, "0", s->isg[0]);CHKERRQ(ierr);
  ierr = PCFieldSplitSetIS(pc, "1", s->isg[1]);CHKERRQ(ierr);
  if (s->userPC) {
    ierr = PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_USER, s->myS);CHKERRQ(ierr);
  }
  if (s->userKSP) {
    ierr = PCSetUp(pc);CHKERRQ(ierr);
    ierr = PCFieldSplitGetSubKSP(pc, &n, &subksp);CHKERRQ(ierr);
    ierr = KSPSetOperators(subksp[1], s->myS, s->myS);CHKERRQ(ierr);
    ierr = PetscFree(subksp);CHKERRQ(ierr);
  }
  PetscFunctionReturn(0);
}

开发者ID:plguhur，项目名称:petsc，代码行数:22，代码来源:ex70.c

示例2: ComputeKSPFETIDP

static PetscErrorCode ComputeKSPFETIDP(DomainData dd, KSP ksp_bddc, KSP *ksp_fetidp)
{
  PetscErrorCode ierr;
  KSP            temp_ksp;
  PC             pc,D;
  Mat            F;

  PetscFunctionBeginUser;
  ierr        = KSPGetPC(ksp_bddc,&pc);CHKERRQ(ierr);
  ierr        = PCBDDCCreateFETIDPOperators(pc,&F,&D);CHKERRQ(ierr);
  ierr        = KSPCreate(PetscObjectComm((PetscObject)F),&temp_ksp);CHKERRQ(ierr);
  ierr        = KSPSetOperators(temp_ksp,F,F);CHKERRQ(ierr);
  ierr        = KSPSetType(temp_ksp,KSPCG);CHKERRQ(ierr);
  ierr        = KSPSetPC(temp_ksp,D);CHKERRQ(ierr);
  ierr        = KSPSetComputeSingularValues(temp_ksp,PETSC_TRUE);CHKERRQ(ierr);
  ierr        = KSPSetFromOptions(temp_ksp);CHKERRQ(ierr);
  ierr        = KSPSetUp(temp_ksp);CHKERRQ(ierr);
  *ksp_fetidp = temp_ksp;
  ierr        = MatDestroy(&F);CHKERRQ(ierr);
  ierr        = PCDestroy(&D);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

开发者ID:tom-klotz，项目名称:petsc，代码行数:22，代码来源:ex59.c

示例3: Solver

HeatSolverBTCS::HeatSolverBTCS(int ny_, double dy_, int nx_, double dx_):
	Solver(ny_, dy_, nx_, dx_)
{
	// Create a sparse matrix A
	two_d_heat_BTCS(A, ny, dy, nx, dx);
	dt = 1.0;
	// set up linear solver context (ksp) and preconditioner (pc)
	KSPCreate(PETSC_COMM_WORLD,&ksp);
	KSPSetOperators(ksp,A,A);
	KSPGetPC(ksp,&pc);
	PCSetType(pc,PCLU);
	KSPSetFromOptions(ksp);

	// create rhs
	VecCreateSeq(PETSC_COMM_SELF, nx*ny, &rhs);
	VecSetFromOptions(rhs);
	VecSet(rhs,0.0);

	// prepare temp
	VecCreateSeq(PETSC_COMM_SELF, nx*ny, &temp);
	VecSetFromOptions(temp);
}

开发者ID:shigh，项目名称:pswr-charm，代码行数:22，代码来源:solver.cpp

示例4: main

int main(int argc, char **argv){
	PetscErrorCode ierr;
	Vec x,b;
	Mat A;
	KSP ksp;

	ierr=PetscInitialize(&argc,&argv,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
	PetscPrintf(PETSC_COMM_WORLD,"]> Initializing PETSc/SLEPc\n");
	
	/*Load data*/
	ierr=loadInputs(&A,&b,&x);CHKERRQ(ierr);
	PetscPrintf(PETSC_COMM_WORLD,"]> Data loaded\n");

	/*Create the KSP context and setup*/
	ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);	
	ierr = KSPSetType(ksp,KSPFGMRES);CHKERRQ(ierr);	
	ierr = KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);	
	ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);	
	ierr = KSPSetUp(ksp);CHKERRQ(ierr);	
	PetscPrintf(PETSC_COMM_WORLD,"]> Krylov Solver settings done\n");

	/*Solve the system*/
	PetscPrintf(PETSC_COMM_WORLD,"]> Krylov Solver Launching solving process\n");
	ierr = KSPSolve(ksp, b, x); CHKERRQ(ierr);
	PetscPrintf(PETSC_COMM_WORLD,"]> Krylov Solver System solved\n");

	/*Clean*/
	ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
	ierr = VecDestroy(&b);CHKERRQ(ierr);
	ierr = VecDestroy(&x);CHKERRQ(ierr);
	ierr = MatDestroy(&A);CHKERRQ(ierr);
	PetscPrintf(PETSC_COMM_WORLD,"]> Cleaned structures, finalizing\n");

	/*Finalize PETSc*/
	PetscFinalize(); 

	return 0;
}

开发者ID:Perif，项目名称:PETScHPCBench，代码行数:38，代码来源:hpc_petsc_bench.c

示例5: PressurePoissonCreate

PetscErrorCode PressurePoissonCreate(  )
{
  PetscErrorCode ierr;
  
  
  PetscFunctionBegin;
  PetscLogEventBegin(EVENT_PressurePoissonCreate,0,0,0,0);
//  PetscLogEventRegister(&EVENT_PressurePoissonCreate,"PressurePoissonCreate", 0);
  ierr = KSPCreate(PETSC_COMM_SELF, &ksp); CHKERRQ(ierr);
  KSPGetPC(ksp, &pc);
//  KSPSetType(ksp, KSPCG);
//  PCSetType(pc, PCICC);
//  PCFactorSetLevels(pc, 0);
//  KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
  KSPSetType(ksp, KSPPREONLY);
  PCSetType(pc, PCCHOLESKY);
  PCFactorSetMatOrderingType(pc, MATORDERING_ND);
  KSPSetFromOptions(ksp);
  KSPSetOperators(ksp, m, m, SAME_PRECONDITIONER);
  PCSetUp(pc);
  
  ierr = IIMCreate(&iim, 12); CHKERRQ(ierr);
//  IIMSetForceComponents(iim, ForceComponentNormalSimple, ForceComponentTangentialSimple);
  CreateGrid2D(d1, d2, &rhsC);
  CreateGrid2D(d1, d2, &p);
  CreateGrid2D(d1, d2, &px);
  CreateGrid2D(d1, d2, &py);
  CreateGrid2D(d1, d2, &u);
  CreateGrid2D(d1, d2, &v);
  ierr = CreateLevelSet2D(d1,d2,&ls); CHKERRQ(ierr);  
  CreateLevelSet2D(d1, d2, &lstemp);
  LevelSetInitializeToStar(ls,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
//  LevelSetInitializeToCircle(ls, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE);
  
  PetscLogEventEnd(EVENT_PressurePoissonCreate,0,0,0,0);
  PetscFunctionReturn(0);
}

开发者ID:adrielb，项目名称:DCell，代码行数:37，代码来源:PressurePoisson.c

示例6: MatAssemblyBegin

void PETSc::Solve_withPureNeumann(void)
{
	ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
  	ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
  	
  	ierr = VecAssemblyBegin(x);
  	ierr = VecAssemblyEnd(x);
  	
  	ierr = VecAssemblyBegin(b);
  	ierr = VecAssemblyEnd(b);
  	
	
	MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,0,PETSC_NULL,&nullsp);
        KSPSetNullSpace(ksp,nullsp);
	
        KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);
        
	//KSPSetType(ksp,KSPMINRES);
	//KSPSetType(ksp,KSPGMRES);
	//KSPSetType(ksp,KSPBCGS);
	KSPSetType(ksp,KSPBCGSL);
	KSPBCGSLSetEll(ksp,2);

	//KSPGetPC(ksp, &pc);
	//PCSetType(pc, PCASM);
	//PCSetType(pc, PCMG);
	//PCMGSetLevels(pc, 3, &PETSC_COMM_WORLD);
	//PCMGSetType(pc,PC_MG_MULTIPLICATIVE);
	//PCMGSetCycleType(pc,PC_MG_CYCLE_V);
	//
        KSPSetFromOptions(ksp);
        KSPSetUp(ksp);

	//start_clock("Before Petsc Solve in pure neumann solver");
        KSPSolve(ksp,b,x);
	//stop_clock("After Petsc Solve in pure neumann solver");
}

开发者ID:tonyguo1，项目名称:LSGFD，代码行数:37，代码来源:solver_petsc.cpp

示例7: main

int main(int argc, char** argv)
{
  Mat            Q;
  Vec            v, a, se;
  KSP            QRsolver;
  PC             pc;
  PetscErrorCode ierr;

  ierr = PetscInitialize(&argc, &argv, NULL, NULL);if (ierr) return ierr;

  ierr = VecCreate(PETSC_COMM_WORLD, &v);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD, &Q);CHKERRQ(ierr);
  ierr = MatSetType(Q, MATDENSE);CHKERRQ(ierr);
  ierr = fill(Q, v);CHKERRQ(ierr);

  ierr = MatCreateVecs(Q, &a, NULL);CHKERRQ(ierr);
  ierr = KSPCreate(PETSC_COMM_WORLD, &QRsolver);CHKERRQ(ierr);
  ierr = KSPGetPC(QRsolver, &pc);CHKERRQ(ierr);
  ierr = PCSetType(pc, PCNONE);CHKERRQ(ierr);
  ierr = KSPSetType(QRsolver, KSPLSQR);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(QRsolver);CHKERRQ(ierr);
  ierr = KSPSetOperators(QRsolver, Q, Q);CHKERRQ(ierr);
  ierr = MatViewFromOptions(Q, NULL, "-sys_view");CHKERRQ(ierr);
  ierr = VecViewFromOptions(a, NULL, "-rhs_view");CHKERRQ(ierr);
  ierr = KSPSolve(QRsolver, v, a);CHKERRQ(ierr);
  ierr = KSPLSQRGetStandardErrorVec(QRsolver, &se);CHKERRQ(ierr);
  if (se) {
    ierr = VecViewFromOptions(se, NULL, "-se_view");CHKERRQ(ierr);
  }
  ierr = KSPDestroy(&QRsolver);CHKERRQ(ierr);
  ierr = VecDestroy(&a);CHKERRQ(ierr);
  ierr = VecDestroy(&v);CHKERRQ(ierr);
  ierr = MatDestroy(&Q);CHKERRQ(ierr);

  ierr = PetscFinalize();
  return ierr;
}

开发者ID:firedrakeproject，项目名称:petsc，代码行数:37，代码来源:ex54.c

示例8: Petsc_KSPSolve

PetscErrorCode Petsc_KSPSolve(AppCtx *obj)
{
  PetscErrorCode ierr;
  KSP            ksp;
  PC             pc;

  /*create the ksp context and set the operators,that is, associate the system matrix with it*/
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
  ierr = KSPSetOperators(ksp,obj->Amat,obj->Amat);CHKERRQ(ierr);

  /*get the preconditioner context, set its type and the tolerances*/
  ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
  ierr = PCSetType(pc,PCLU);CHKERRQ(ierr);
  ierr = KSPSetTolerances(ksp,1.e-7,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr);

  /*get the command line options if there are any and set them*/
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);

  /*get the linear system (ksp) solve*/
  ierr = KSPSolve(ksp,obj->ksp_rhs,obj->ksp_sol);CHKERRQ(ierr);

  KSPDestroy(&ksp);
  return 0;
}

开发者ID:wgapl，项目名称:petsc，代码行数:24，代码来源:ex3.c

示例9: KSPCreate

void Field_solver::create_solver_and_preconditioner( KSP *ksp, PC *pc, Mat *A )
{
    PetscReal rtol = 1.e-12;
    // Default.     
    // Possible to specify from command line using '-ksp_rtol' option.
    
    PetscErrorCode ierr;
    ierr = KSPCreate( PETSC_COMM_WORLD, ksp ); CHKERRXX(ierr);
    ierr = KSPSetOperators( *ksp, *A, *A, DIFFERENT_NONZERO_PATTERN ); CHKERRXX(ierr);
    //ierr = KSPSetOperators( *ksp, *A, *A ); CHKERRXX(ierr);
    ierr = KSPGetPC( *ksp, pc ); CHKERRXX(ierr);
    ierr = PCSetType( *pc, PCGAMG ); CHKERRXX(ierr);
    ierr = KSPSetType( *ksp, KSPGMRES ); CHKERRXX(ierr);
    ierr = KSPSetTolerances( *ksp, rtol,
			     PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT); CHKERRXX(ierr);
    ierr = KSPSetFromOptions( *ksp ); CHKERRXX(ierr);    
    ierr = KSPSetInitialGuessNonzero( *ksp, PETSC_TRUE ); CHKERRXX( ierr );

    // For test purposes
    //ierr = KSPSetInitialGuessNonzero( *ksp, PETSC_FALSE ); CHKERRXX( ierr );

    ierr = KSPSetUp( *ksp ); CHKERRXX(ierr);
    return;
}

开发者ID:epicf，项目名称:ef，代码行数:24，代码来源:field_solver.cpp

示例10: main

int main(int argc,char **args)
{
  Vec            x,b,u;      /* approx solution, RHS, exact solution */
  Mat            A;            /* linear system matrix */
  KSP            ksp;         /* KSP context */
  KSP            *subksp;     /* array of local KSP contexts on this processor */
  PC             pc;           /* PC context */
  PC             subpc;        /* PC context for subdomain */
  PetscReal      norm;         /* norm of solution error */
  PetscErrorCode ierr;
  PetscInt       i,j,Ii,J,*blks,m = 4,n;
  PetscMPIInt    rank,size;
  PetscInt       its,nlocal,first,Istart,Iend;
  PetscScalar    v,one = 1.0,none = -1.0;
  PetscBool      isbjacobi;

  PetscInitialize(&argc,&args,(char*)0,help);
  ierr = PetscOptionsGetInt(NULL,"-m",&m,NULL);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
  n    = m+2;

  /* -------------------------------------------------------------------
         Compute the matrix and right-hand-side vector that define
         the linear system, Ax = b.
     ------------------------------------------------------------------- */

  /*
     Create and assemble parallel matrix
  */
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatMPIAIJSetPreallocation(A,5,NULL,5,NULL);CHKERRQ(ierr);
  ierr = MatSeqAIJSetPreallocation(A,5,NULL);CHKERRQ(ierr);
  ierr = MatGetOwnershipRange(A,&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(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);}
    if (i<m-1) {J = Ii + n; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);}
    if (j>0)   {J = Ii - 1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);}
    if (j<n-1) {J = Ii + 1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);}
    v = 4.0; ierr = MatSetValues(A,1,&Ii,1,&Ii,&v,ADD_VALUES);CHKERRQ(ierr);
  }
  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);

  /*
     Create parallel vectors
  */
  ierr = VecCreate(PETSC_COMM_WORLD,&u);CHKERRQ(ierr);
  ierr = VecSetSizes(u,PETSC_DECIDE,m*n);CHKERRQ(ierr);
  ierr = VecSetFromOptions(u);CHKERRQ(ierr);
  ierr = VecDuplicate(u,&b);CHKERRQ(ierr);
  ierr = VecDuplicate(b,&x);CHKERRQ(ierr);

  /*
     Set exact solution; then compute right-hand-side vector.
  */
  ierr = VecSet(u,one);CHKERRQ(ierr);
  ierr = MatMult(A,u,b);CHKERRQ(ierr);

  /*
     Create linear solver context
  */
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);

  /*
     Set operators. Here the matrix that defines the linear system
     also serves as the preconditioning matrix.
  */
  ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);

  /*
     Set default preconditioner for this program to be block Jacobi.
     This choice can be overridden at runtime with the option
        -pc_type <type>
  */
  ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
  ierr = PCSetType(pc,PCBJACOBI);CHKERRQ(ierr);


  /* -------------------------------------------------------------------
                   Define the problem decomposition
     ------------------------------------------------------------------- */

  /*
     Call PCBJacobiSetTotalBlocks() to set individually the size of
     each block in the preconditioner.  This could also be done with
     the runtime option
         -pc_bjacobi_blocks <blocks>
     Also, see the command PCBJacobiSetLocalBlocks() to set the
     local blocks.

      Note: The default decomposition is 1 block per processor.
  */
  ierr = PetscMalloc1(m,&blks);CHKERRQ(ierr);
  for (i=0; i<m; i++) blks[i] = n;
  ierr = PCBJacobiSetTotalBlocks(pc,m,blks);CHKERRQ(ierr);
  ierr = PetscFree(blks);CHKERRQ(ierr);
//.........这里部分代码省略.........

开发者ID:fengyuqi，项目名称:petsc，代码行数:101，代码来源:ex7.c

示例11: TaylorGalerkinStepIIMomentum


//.........这里部分代码省略.........
  PetscReal      hx, hy, area;
  const PetscInt *necon;
  PetscInt       j, k, e, ne, nc, mx, my;
  PetscErrorCode ierr;

  PetscFunctionBeginUser;
  ierr = PetscObjectGetComm((PetscObject) da, &comm);CHKERRQ(ierr);
  ierr = DMSetMatType(da,MATAIJ);CHKERRQ(ierr);
  ierr = DMCreateMatrix(da, &mat);CHKERRQ(ierr);
  ierr = MatSetOption(mat,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da, &rhs_u);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da, &rhs_v);CHKERRQ(ierr);
  ierr = KSPCreate(comm, &ksp);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);

  ierr = DMDAGetInfo(da, 0, &mx, &my, 0,0,0,0,0,0,0,0,0,0);CHKERRQ(ierr);
  hx   = 1.0 / (PetscReal)(mx-1);
  hy   = 1.0 / (PetscReal)(my-1);
  area = 0.5*hx*hy;
  ierr = VecGetArray(user->sol_n.u,       &u_n);CHKERRQ(ierr);
  ierr = VecGetArray(user->sol_n.v,       &v_n);CHKERRQ(ierr);
  ierr = VecGetArray(user->mu,            &mu_n);CHKERRQ(ierr);
  ierr = VecGetArray(user->sol_phi.u,     &u_phi);CHKERRQ(ierr);
  ierr = VecGetArray(user->sol_phi.v,     &v_phi);CHKERRQ(ierr);
  ierr = VecGetArray(user->sol_phi.rho_u, &rho_u_phi);CHKERRQ(ierr);
  ierr = VecGetArray(user->sol_phi.rho_v, &rho_v_phi);CHKERRQ(ierr);
  ierr = DMDAGetElements(da, &ne, &nc, &necon);CHKERRQ(ierr);
  for (e = 0; e < ne; e++) {
    for (j = 0; j < 3; j++) {
      idx[j]      = necon[3*e+j];
      values_u[j] = 0.0;
      values_v[j] = 0.0;
    }
    /* Get basis function deriatives (we need the orientation of the element here) */
    if (idx[1] > idx[0]) {
      psi_x[0] = -hy; psi_x[1] =  hy; psi_x[2] = 0.0;
      psi_y[0] = -hx; psi_y[1] = 0.0; psi_y[2] =  hx;
    } else {
      psi_x[0] =  hy; psi_x[1] = -hy; psi_x[2] = 0.0;
      psi_y[0] =  hx; psi_y[1] = 0.0; psi_y[2] = -hx;
    }
    /*  <\nabla\psi, F^{n+\phi}_e>: Divergence of the element-averaged convective fluxes */
    for (j = 0; j < 3; j++) {
      values_u[j] += psi_x[j]*rho_u_phi[e]*u_phi[e] + psi_y[j]*rho_u_phi[e]*v_phi[e];
      values_v[j] += psi_x[j]*rho_v_phi[e]*u_phi[e] + psi_y[j]*rho_v_phi[e]*v_phi[e];
    }
    /*  -<\nabla\psi, F^n_v>: Divergence of the viscous fluxes */
    for (j = 0; j < 3; j++) {
      /* \tau_{xx} = 2/3 \mu(T) (2 {\partial u\over\partial x} - {\partial v\over\partial y}) */
      /* \tau_{xy} =     \mu(T) (  {\partial u\over\partial y} + {\partial v\over\partial x}) */
      /* \tau_{yy} = 2/3 \mu(T) (2 {\partial v\over\partial y} - {\partial u\over\partial x}) */
      mu     = 0.0;
      tau_xx = 0.0;
      tau_xy = 0.0;
      tau_yy = 0.0;
      for (k = 0; k < 3; k++) {
        mu     += mu_n[idx[k]];
        tau_xx += 2.0*psi_x[k]*u_n[idx[k]] - psi_y[k]*v_n[idx[k]];
        tau_xy +=     psi_y[k]*u_n[idx[k]] + psi_x[k]*v_n[idx[k]];
        tau_yy += 2.0*psi_y[k]*v_n[idx[k]] - psi_x[k]*u_n[idx[k]];
      }
      mu          /= 3.0;
      tau_xx      *= (2.0/3.0)*mu;
      tau_xy      *= mu;
      tau_yy      *= (2.0/3.0)*mu;
      values_u[j] -= area*(psi_x[j]*tau_xx + psi_y[j]*tau_xy);
      values_v[j] -= area*(psi_x[j]*tau_xy + psi_y[j]*tau_yy);
    }
    /* Accumulate to global structures */
    ierr = VecSetValuesLocal(rhs_u, 3, idx, values_u, ADD_VALUES);CHKERRQ(ierr);
    ierr = VecSetValuesLocal(rhs_v, 3, idx, values_v, ADD_VALUES);CHKERRQ(ierr);
    ierr = MatSetValuesLocal(mat, 3, idx, 3, idx, identity, ADD_VALUES);CHKERRQ(ierr);
  }
  ierr = DMDARestoreElements(da, &ne, &nc, &necon);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_n.u,       &u_n);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_n.v,       &v_n);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->mu,            &mu_n);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_phi.u,     &u_phi);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_phi.v,     &v_phi);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_phi.rho_u, &rho_u_phi);CHKERRQ(ierr);
  ierr = VecRestoreArray(user->sol_phi.rho_v, &rho_v_phi);CHKERRQ(ierr);

  ierr = VecAssemblyBegin(rhs_u);CHKERRQ(ierr);
  ierr = VecAssemblyBegin(rhs_v);CHKERRQ(ierr);
  ierr = MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(rhs_u);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(rhs_v);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = VecScale(rhs_u,user->dt);CHKERRQ(ierr);
  ierr = VecScale(rhs_v,user->dt);CHKERRQ(ierr);

  ierr = KSPSetOperators(ksp, mat, mat, DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
  ierr = KSPSolve(ksp, rhs_u, user->sol_np1.rho_u);CHKERRQ(ierr);
  ierr = KSPSolve(ksp, rhs_v, user->sol_np1.rho_v);CHKERRQ(ierr);
  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = MatDestroy(&mat);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da, &rhs_u);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da, &rhs_v);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

开发者ID:hsahasra，项目名称:petsc-magma-dense-mat，代码行数:101，代码来源:ex31.c

示例12: PCSetUp_MG

PetscErrorCode PCSetUp_MG(PC pc)
{
  PC_MG          *mg        = (PC_MG*)pc->data;
  PC_MG_Levels   **mglevels = mg->levels;
  PetscErrorCode ierr;
  PetscInt       i,n = mglevels[0]->levels;
  PC             cpc;
  PetscBool      dump = PETSC_FALSE,opsset,use_amat,missinginterpolate = PETSC_FALSE;
  Mat            dA,dB;
  Vec            tvec;
  DM             *dms;
  PetscViewer    viewer = 0;

  PetscFunctionBegin;
  /* FIX: Move this to PCSetFromOptions_MG? */
  if (mg->usedmfornumberoflevels) {
    PetscInt levels;
    ierr = DMGetRefineLevel(pc->dm,&levels);CHKERRQ(ierr);
    levels++;
    if (levels > n) { /* the problem is now being solved on a finer grid */
      ierr     = PCMGSetLevels(pc,levels,NULL);CHKERRQ(ierr);
      n        = levels;
      ierr     = PCSetFromOptions(pc);CHKERRQ(ierr); /* it is bad to call this here, but otherwise will never be called for the new hierarchy */
      mglevels =  mg->levels;
    }
  }
  ierr = KSPGetPC(mglevels[0]->smoothd,&cpc);CHKERRQ(ierr);


  /* If user did not provide fine grid operators OR operator was not updated since last global KSPSetOperators() */
  /* so use those from global PC */
  /* Is this what we always want? What if user wants to keep old one? */
  ierr = KSPGetOperatorsSet(mglevels[n-1]->smoothd,NULL,&opsset);CHKERRQ(ierr);
  if (opsset) {
    Mat mmat;
    ierr = KSPGetOperators(mglevels[n-1]->smoothd,NULL,&mmat);CHKERRQ(ierr);
    if (mmat == pc->pmat) opsset = PETSC_FALSE;
  }

  if (!opsset) {
    ierr = PCGetUseAmat(pc,&use_amat);CHKERRQ(ierr);
    if(use_amat){
      ierr = PetscInfo(pc,"Using outer operators to define finest grid operator \n  because PCMGGetSmoother(pc,nlevels-1,&ksp);KSPSetOperators(ksp,...); was not called.\n");CHKERRQ(ierr);
      ierr = KSPSetOperators(mglevels[n-1]->smoothd,pc->mat,pc->pmat);CHKERRQ(ierr);
    }
    else {
      ierr = PetscInfo(pc,"Using matrix (pmat) operators to define finest grid operator \n  because PCMGGetSmoother(pc,nlevels-1,&ksp);KSPSetOperators(ksp,...); was not called.\n");CHKERRQ(ierr);
      ierr = KSPSetOperators(mglevels[n-1]->smoothd,pc->pmat,pc->pmat);CHKERRQ(ierr);
    }
  }

  for (i=n-1; i>0; i--) {
    if (!(mglevels[i]->interpolate || mglevels[i]->restrct)) {
      missinginterpolate = PETSC_TRUE;
      continue;
    }
  }
  /*
   Skipping if user has provided all interpolation/restriction needed (since DM might not be able to produce them (when coming from SNES/TS)
   Skipping for galerkin==2 (externally managed hierarchy such as ML and GAMG). Cleaner logic here would be great. Wrap ML/GAMG as DMs?
  */
  if (missinginterpolate && pc->dm && mg->galerkin != 2 && !pc->setupcalled) {
    /* construct the interpolation from the DMs */
    Mat p;
    Vec rscale;
    ierr     = PetscMalloc1(n,&dms);CHKERRQ(ierr);
    dms[n-1] = pc->dm;
    /* Separately create them so we do not get DMKSP interference between levels */
    for (i=n-2; i>-1; i--) {ierr = DMCoarsen(dms[i+1],MPI_COMM_NULL,&dms[i]);CHKERRQ(ierr);}
    for (i=n-2; i>-1; i--) {
      DMKSP     kdm;
      PetscBool dmhasrestrict;
      ierr = KSPSetDM(mglevels[i]->smoothd,dms[i]);CHKERRQ(ierr);
      if (mg->galerkin) {ierr = KSPSetDMActive(mglevels[i]->smoothd,PETSC_FALSE);CHKERRQ(ierr);}
      ierr = DMGetDMKSPWrite(dms[i],&kdm);CHKERRQ(ierr);
      /* Ugly hack so that the next KSPSetUp() will use the RHS that we set. A better fix is to change dmActive to take
       * a bitwise OR of computing the matrix, RHS, and initial iterate. */
      kdm->ops->computerhs = NULL;
      kdm->rhsctx          = NULL;
      if (!mglevels[i+1]->interpolate) {
        ierr = DMCreateInterpolation(dms[i],dms[i+1],&p,&rscale);CHKERRQ(ierr);
        ierr = PCMGSetInterpolation(pc,i+1,p);CHKERRQ(ierr);
        if (rscale) {ierr = PCMGSetRScale(pc,i+1,rscale);CHKERRQ(ierr);}
        ierr = VecDestroy(&rscale);CHKERRQ(ierr);
        ierr = MatDestroy(&p);CHKERRQ(ierr);
      }
      ierr = DMHasCreateRestriction(dms[i],&dmhasrestrict);CHKERRQ(ierr);
      if (dmhasrestrict && !mglevels[i+1]->restrct){
        ierr = DMCreateRestriction(dms[i],dms[i+1],&p);CHKERRQ(ierr);
        ierr = PCMGSetRestriction(pc,i+1,p);CHKERRQ(ierr);
        ierr = MatDestroy(&p);CHKERRQ(ierr);
      }
    }

    for (i=n-2; i>-1; i--) {ierr = DMDestroy(&dms[i]);CHKERRQ(ierr);}
    ierr = PetscFree(dms);CHKERRQ(ierr);
  }

  if (pc->dm && !pc->setupcalled) {
    /* finest smoother also gets DM but it is not active, independent of whether galerkin==2 */
//.........这里部分代码省略.........

开发者ID:ziolai，项目名称:petsc，代码行数:101，代码来源:mg.c

示例13: main

int main(int argc,char **args)
{
    Mat            C;
    PetscErrorCode ierr;
    PetscInt       N = 2,rowidx,colidx;
    Vec            u,b,r;
    KSP            ksp;
    PetscReal      norm;
    PetscMPIInt    rank,size;
    PetscScalar    v;

    PetscInitialize(&argc,&args,(char*)0,help);
    ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
    CHKERRQ(ierr);
    ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);
    CHKERRQ(ierr);

    /* create stiffness matrix C = [1 2; 2 3] */
    ierr = MatCreate(PETSC_COMM_WORLD,&C);
    CHKERRQ(ierr);
    ierr = MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,N,N);
    CHKERRQ(ierr);
    ierr = MatSetFromOptions(C);
    CHKERRQ(ierr);
    ierr = MatSetUp(C);
    CHKERRQ(ierr);
    if (!rank) {
        rowidx = 0;
        colidx = 0;
        v = 1.0;
        ierr   = MatSetValues(C,1,&rowidx,1,&colidx,&v,INSERT_VALUES);
        CHKERRQ(ierr);
        rowidx = 0;
        colidx = 1;
        v = 2.0;
        ierr   = MatSetValues(C,1,&rowidx,1,&colidx,&v,INSERT_VALUES);
        CHKERRQ(ierr);

        rowidx = 1;
        colidx = 0;
        v = 2.0;
        ierr   = MatSetValues(C,1,&rowidx,1,&colidx,&v,INSERT_VALUES);
        CHKERRQ(ierr);
        rowidx = 1;
        colidx = 1;
        v = 3.0;
        ierr   = MatSetValues(C,1,&rowidx,1,&colidx,&v,INSERT_VALUES);
        CHKERRQ(ierr);
    }
    ierr = MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
    CHKERRQ(ierr);
    ierr = MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
    CHKERRQ(ierr);

    /* create right hand side and solution */
    ierr = VecCreate(PETSC_COMM_WORLD,&u);
    CHKERRQ(ierr);
    ierr = VecSetSizes(u,PETSC_DECIDE,N);
    CHKERRQ(ierr);
    ierr = VecSetFromOptions(u);
    CHKERRQ(ierr);
    ierr = VecDuplicate(u,&b);
    CHKERRQ(ierr);
    ierr = VecDuplicate(u,&r);
    CHKERRQ(ierr);
    ierr = VecSet(u,0.0);
    CHKERRQ(ierr);
    ierr = VecSet(b,1.0);
    CHKERRQ(ierr);

    /* solve linear system C*u = b */
    ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);
    CHKERRQ(ierr);
    ierr = KSPSetOperators(ksp,C,C);
    CHKERRQ(ierr);
    ierr = KSPSetFromOptions(ksp);
    CHKERRQ(ierr);
    ierr = KSPSolve(ksp,b,u);
    CHKERRQ(ierr);

    /* check residual r = C*u - b */
    ierr = MatMult(C,u,r);
    CHKERRQ(ierr);
    ierr = VecAXPY(r,-1.0,b);
    CHKERRQ(ierr);
    ierr = VecNorm(r,NORM_2,&norm);
    CHKERRQ(ierr);
    ierr = PetscPrintf(PETSC_COMM_WORLD,"|| C*u - b|| = %g\n",(double)norm);
    CHKERRQ(ierr);

    /* solve C^T*u = b twice */
    ierr = KSPSolveTranspose(ksp,b,u);
    CHKERRQ(ierr);
    /* check residual r = C^T*u - b */
    ierr = MatMultTranspose(C,u,r);
    CHKERRQ(ierr);
    ierr = VecAXPY(r,-1.0,b);
    CHKERRQ(ierr);
    ierr = VecNorm(r,NORM_2,&norm);
    CHKERRQ(ierr);
//.........这里部分代码省略.........

开发者ID:00liujj，项目名称:petsc，代码行数:101，代码来源:ex2.c

示例14: main


//.........这里部分代码省略.........
    ierr = MatAssemblyEnd(C1,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);

    /*
       Indicate same nonzero structure of successive linear system matrices
    */
    ierr = MatSetOption(C1,MAT_NEW_NONZERO_LOCATIONS,PETSC_TRUE);CHKERRQ(ierr);

    /*
       Compute right-hand-side vector
    */
    ierr = MatMult(C1,u,b1);CHKERRQ(ierr);

    /*
       Set operators. Here the matrix that defines the linear system
       also serves as the preconditioning matrix.
        - The flag SAME_NONZERO_PATTERN indicates that the
          preconditioning matrix has identical nonzero structure
          as during the last linear solve (although the values of
          the entries have changed). Thus, we can save some
          work in setting up the preconditioner (e.g., no need to
          redo symbolic factorization for ILU/ICC preconditioners).
        - If the nonzero structure of the matrix is different during
          the second linear solve, then the flag DIFFERENT_NONZERO_PATTERN
          must be used instead.  If you are unsure whether the
          matrix structure has changed or not, use the flag
          DIFFERENT_NONZERO_PATTERN.
        - Caution:  If you specify SAME_NONZERO_PATTERN, PETSc
          believes your assertion and does not check the structure
          of the matrix.  If you erroneously claim that the structure
          is the same when it actually is not, the new preconditioner
          will not function correctly.  Thus, use this optimization
          feature with caution!
    */
    ierr = KSPSetOperators(ksp1,C1,C1,SAME_NONZERO_PATTERN);CHKERRQ(ierr);

    /*
       Use the previous solution of linear system #1 as the initial
       guess for the next solve of linear system #1.  The user MUST
       call KSPSetInitialGuessNonzero() in indicate use of an initial
       guess vector; otherwise, an initial guess of zero is used.
    */
    if (t>0) {
      ierr = KSPSetInitialGuessNonzero(ksp1,PETSC_TRUE);CHKERRQ(ierr);
    }

    /*
       Solve the first linear system.  Here we explicitly call
       KSPSetUp() for more detailed performance monitoring of
       certain preconditioners, such as ICC and ILU.  This call
       is optional, ase KSPSetUp() will automatically be called
       within KSPSolve() if it hasn't been called already.
    */
    ierr = KSPSetUp(ksp1);CHKERRQ(ierr);
    ierr = KSPSolve(ksp1,b1,x1);CHKERRQ(ierr);
    ierr = KSPGetIterationNumber(ksp1,&its);CHKERRQ(ierr);

    /*
       Check error of solution to first linear system
    */
    ierr = CheckError(u,x1,b1,its,1.e-4,CHECK_ERROR);CHKERRQ(ierr);

    /* - - - - - - - - - - - - Stage 2: - - - - - - - - - - - - - -
                 Assemble and solve second linear system
       - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

    /*

开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:67，代码来源:ex9.c

示例15: internal_solve

  long internal_solve()
  {

    /* 
       Create linear solver context
    */
    ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
    ierr = KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
    ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
    ierr = PCSetType(pc,PCJACOBI);CHKERRQ(ierr);
    ierr = KSPSetTolerances(ksp,1.e-7,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr);
    ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
 
    /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
       Solve the linear system
       - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
    ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr); 

    PetscScalar*   solution_vector;
    VecGetArray(x,&solution_vector);
    for (long i =0; i < n; ++i)
      {
	sol_vec.push_back(solution_vector[i]);
      }

    int reason=0;
    KSPGetConvergedReason(ksp,(KSPConvergedReason*)&reason);
    switch(reason)
      {
      case 2:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_RTOL .. " << std::endl;
	break;
      case 3:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_ATOL .. " << std::endl;
	break;
      case 4:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_ITS .. " << std::endl;
	break;
      case 5:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_QCG_NEG_CURVE .. " << std::endl;
	break;
      case 6:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_QCG_CONSTRAINED .. " << std::endl;
	break;
      case 7:
	std::cout << "### CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_CONVERGED_STEP_LENGTH .. " << std::endl;
	break;
      case -3:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_ITS .. " << std::endl;
	break;
      case -4:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_DTOL .. " << std::endl;
	break;
      case -5:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_BREAKDOWN  .. " << std::endl;
	break;
      case -6:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_BREAKDOWN_BICG .. " << std::endl;
	break;
      case -7:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_NONSYMMETRIC .. " << std::endl;
	break;
      case -8:
	std::cout << "### !! N O T  !! CONVERGED ### " << std::endl;
	std::cout << "### KSP ..KSP_DIVERGED_INDEFINITE_PC .. " << std::endl;
	break;

      default:
	std::cout << "### KSP .. no problem detected.. " << std::endl;
	break;

      }

    

    /* 
    ierr = KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);



       Check the error
    ierr = VecAXPY(x,neg_one,u);CHKERRQ(ierr);
    ierr  = VecNorm(x,NORM_2,&norm);CHKERRQ(ierr);
    ierr = KSPGetIterationNumber(ksp,&its);CHKERRQ(ierr);
    ierr = PetscPrintf(PETSC_COMM_WORLD,"Norm of error %A, Iterations %D\n",
		       norm,its);CHKERRQ(ierr);
    */
    return 0;
  }

开发者ID:AlexanderToifl，项目名称:viennamesh-dev，代码行数:100，代码来源:petsc_solver_interface.hpp

示例16: init_fft2d_


//.........这里部分代码省略.........
  
  if(((ly-topo_.mylc)!=0) || topo_.npzc>1) {
  printf("Error,npz should equal to 1, or %d\t%d\n",irankx,ly-topo_.mylc);
  MPI_Abort(commx,1);
  }

  minp = fftw_alloc_complex(alloc_ly);
  mout = fftw_alloc_complex(alloc_ly);
  
  if( !(freq = r3tensor(topo_.mxlc, topo_.mylc*topo_.mzlc, 2)) ) 
  	printf("Malloc error!\n");  
  if( !(data = r3tensor(topo_.mxlc, topo_.mylc*topo_.mzlc, 2)) ) 
  	printf("Malloc error!\n");
/*  if( !(dar = r3tensor(topo_.mxlc, topo_.mylc, topo_.mzlc)) ) 
  	printf("Malloc error!\n");
  if( !(dai = r3tensor(topo_.mxlc, topo_.mylc, topo_.mzlc)) ) 
  	printf("Malloc error!\n");
*/
//***********
  mplanF = fftw_mpi_plan_dft_2d(my, mz, minp, mout, commx, FFTW_FORWARD, FFTW_MEASURE);
  mplanR = fftw_mpi_plan_dft_2d(my, mz, minp, mout, commx, FFTW_BACKWARD, FFTW_MEASURE);
/*  mplanF = fftw_mpi_plan_many_dft(rnk, myz, howmany,
        FFTW_MPI_DEFAULT_BLOCK,FFTW_MPI_DEFAULT_BLOCK,
        minp, mout, commx, FFTW_FORWARD, FFTW_MEASURE);
  mplanR = fftw_mpi_plan_many_dft(rnk, myz, howmany,
        FFTW_MPI_DEFAULT_BLOCK,FFTW_MPI_DEFAULT_BLOCK,
        minp, mout, commx, FFTW_BACKWARD, FFTW_MEASURE);
*/
//***********
//***** Solver part ******

  dxs = topo_.dx0*topo_.dx0;
  dys = topo_.dy0*topo_.dy0;
  dzs = topo_.dz0*topo_.dz0;

  vm2=-1.0/12.0;
  vm1=16.0/12.0;
  v  =-30.0/12.0;
  vp1=16.0/12.0;
  vp2=-1.0/12.0;

  MatCreateMPIAIJ(commyz, PETSC_DECIDE, PETSC_DECIDE, mx, mx, 5, PETSC_NULL, 5, PETSC_NULL, &A);

  ierr = MatGetOwnershipRange(A,&Istart,&Iend);

  
  for (Ii=Istart; Ii<Iend; Ii++) {
    i = Ii; j = Ii;
    if ((i>1)&&(i<mx-2))   {

	J = Ii - 2; MatSetValues(A,1,&Ii,1,&J,&vm2,INSERT_VALUES);
	J = Ii - 1; ierr = MatSetValues(A,1,&Ii,1,&J,&vm1,INSERT_VALUES);CHKERRQ(ierr);
	J = Ii; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);
	J = Ii + 1; ierr = MatSetValues(A,1,&Ii,1,&J,&vp1,INSERT_VALUES);CHKERRQ(ierr);
	J = Ii + 2; ierr = MatSetValues(A,1,&Ii,1,&J,&vp2,INSERT_VALUES);CHKERRQ(ierr);
    }

    if (i==0) {
	J = Ii; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);}
    if (i==1)   {
	J = Ii - 1; vb = 11.0/12.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);CHKERRQ(ierr);
	J = Ii ; vb = -5.0/3.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
	J = Ii + 1; vb = 0.5; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
	J = Ii + 2; vb = 1.0/3.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
	J = Ii + 3; vb = -1.0/12.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
}

    if (i==mx-2) {
        J = Ii + 1; vb = 11.0/12.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);CHKERRQ(ierr);
        J = Ii ; vb = -5.0/3.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
        J = Ii - 1; vb = 0.5; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
        J = Ii - 2; vb = 1.0/3.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
        J = Ii - 3; vb = -1.0/12.0; ierr = MatSetValues(A,1,&Ii,1,&J,&vb,INSERT_VALUES);
}
    if (i==mx-1) {J = Ii; ierr = MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);} 

    }

  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);

  ierr = VecCreate(commyz,&br);CHKERRQ(ierr);
  ierr = VecSetSizes(br,PETSC_DECIDE,mx);CHKERRQ(ierr);
  ierr = VecSetFromOptions(br);CHKERRQ(ierr);
  ierr = VecDuplicate(br,&xr);CHKERRQ(ierr);
  ierr = VecDuplicate(br,&bi);CHKERRQ(ierr);
  ierr = VecDuplicate(br,&xi);CHKERRQ(ierr);


  ierr = KSPCreate(commyz,&ksp);CHKERRQ(ierr);

  ierr = KSPSetOperators(ksp,A,A,SAME_PRECONDITIONER);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
  ierr = KSPGetPC(ksp,&pc);
  PCSetType(pc,PCJACOBI);
  ierr = KSPSetTolerances(ksp,1.e-7,1.e-50,PETSC_DEFAULT,
                          PETSC_DEFAULT);CHKERRQ(ierr);

//******* End  *******
}

开发者ID:LeeJH，项目名称:pde，代码行数:101，代码来源:fft2d.c

示例17: main

int main(int argc,char **argv)
{
  PetscErrorCode ierr;
  KSP            ksp;
  PC             pc;
  Vec            x,b;
  DM             da;
  Mat            A;
  PetscInt       dof=1;
  PetscBool      flg;
  PetscScalar    zero=0.0;

  PetscInitialize(&argc,&argv,(char *)0,help);
  ierr = PetscOptionsGetInt(PETSC_NULL,"-dof",&dof,PETSC_NULL);CHKERRQ(ierr);

  ierr = DMDACreate(PETSC_COMM_WORLD,&da);CHKERRQ(ierr);
  ierr = DMDASetDim(da,3);CHKERRQ(ierr);
  ierr = DMDASetBoundaryType(da,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE);CHKERRQ(ierr);
  ierr = DMDASetStencilType(da,DMDA_STENCIL_STAR);CHKERRQ(ierr);
  ierr = DMDASetSizes(da,3,3,3);CHKERRQ(ierr);
  ierr = DMDASetNumProcs(da,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
  ierr = DMDASetDof(da,dof);CHKERRQ(ierr);
  ierr = DMDASetStencilWidth(da,1);CHKERRQ(ierr);
  ierr = DMDASetOwnershipRanges(da,PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
  ierr = DMSetFromOptions(da);CHKERRQ(ierr);
  ierr = DMSetUp(da);CHKERRQ(ierr);

  ierr = DMCreateGlobalVector(da,&x);CHKERRQ(ierr);
  ierr = DMCreateGlobalVector(da,&b);CHKERRQ(ierr);
  ierr = DMCreateMatrix(da,MATAIJ,&A);CHKERRQ(ierr);
  ierr = VecSet(b,zero);CHKERRQ(ierr);

  /* Test sbaij matrix */
  flg = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL,"-test_sbaij",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) {
    Mat sA;
    ierr = MatSetOption(A,MAT_SYMMETRIC,PETSC_TRUE);CHKERRQ(ierr);
    ierr = MatConvert(A,MATSBAIJ,MAT_INITIAL_MATRIX,&sA);CHKERRQ(ierr);
    ierr = MatDestroy(&A);CHKERRQ(ierr);
    A = sA;
  }

  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
  ierr = KSPSetOperators(ksp,A,A,SAME_NONZERO_PATTERN);CHKERRQ(ierr);
  ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
  ierr = PCSetDM(pc,(DM)da);CHKERRQ(ierr);

  ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr);

  /* check final residual */
  flg  = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL, "-check_final_residual", &flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg){
    Vec            b1;
    PetscReal      norm;
    ierr = KSPGetSolution(ksp,&x);CHKERRQ(ierr);
    ierr = VecDuplicate(b,&b1);CHKERRQ(ierr);
    ierr = MatMult(A,x,b1);CHKERRQ(ierr);
    ierr = VecAXPY(b1,-1.0,b);CHKERRQ(ierr);
    ierr = VecNorm(b1,NORM_2,&norm);CHKERRQ(ierr);
    ierr = PetscPrintf(PETSC_COMM_WORLD,"Final residual %g\n",norm);CHKERRQ(ierr);
    ierr = VecDestroy(&b1);CHKERRQ(ierr);
  }

  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = VecDestroy(&x);CHKERRQ(ierr);
  ierr = VecDestroy(&b);CHKERRQ(ierr);
  ierr = MatDestroy(&A);CHKERRQ(ierr);
  ierr = DMDestroy(&da);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:74，代码来源:ex35.c

示例18: SNESSolve_NEWTONTR

鲜花

握手

雷人

路过

鸡蛋

该文章已有0人参与评论

请发表评论

全部评论

专题导读

More+

10-27 六六分期app的软件客服如何联系？(六六分期

11-06 可心卡盟:win10系统火狐flash插件崩溃怎么

11-06 亲亲特价:怎么删除回收站图标

11-06 济南大学虚拟社区:鲁大师节能降温的具体办

11-06 xlueops.exe:无线网络安装向导

11-06 女斗合众国:win7系统cf与主机连接不稳定怎

11-06 0xc000022-[cf烟雾头]cf怎么调烟雾头

11-06 qizideyouhuo:应用程序无法正常启动0xc0000

11-06 ipz-185:win7系统vcf文件怎么打开

11-06 傻哥蹦迪:win10系统s4怎么打开usb调试

11-06 八神浩树gtaste:回收站清空了怎么恢复

11-06 妖尾之黑色守护:win10系统电脑没有1440x900

11-06 校园至尊魔王小说:win7系统浏览网页时字体

11-06 女斗合众国:win10系统访问共享文件夹提示请

11-06 tokyo hot n0654:恢复win7系统默认字体一招

11-06 雨酷仙境:设置win7系统转移临时文件夹腾出

11-06 阿穆纳伊之杖:win7系统开始菜单在右边还原

11-06 tunespotting:win10系统火狐flash插件总是

11-06 甘尔葛分析师：计谋网站seo关键词暴涨有什

11-06 蔡贵霖: 计谋网站seo关键词暴涨有什么秘密

11-06 博益网首页:ao3网页版进入不了解决方法

11-06 漏斗子专栏: 网站数据分析小白易懂精华篇

11-06 见证双虹怎么做:win7系统开启telnet命令的

11-06 颾狐蝶蜋:系统资源不足无法完成请求的服务

11-06 国光中学校歌:提交网站到alexa查询详细步骤

11-06 西安有情天:静态网页和动态网页的区别

11-06 红木雅尚斋:外部链接构造对网站的好处

11-06 前官礼遇：防止域名劫持–增强域安全性的10

11-06 密传二转答案: 中文分词算法有哪些

11-06 金泉家园邮编:百度快照劫持的表现及应对方

C++ KSPSetTolerances函数代码示例发布时间：2022-05-30

C++ KSPSetFromOptions函数代码示例发布时间：2022-05-30

剪的笔顺,诠释剪的笔画,认识剪的部首

1 六六分期app的软件客服如何联系？(六六分期

六六分期app的软件客服如何联系？不知道吗？加qq群【895510560】即可！标题：六六分期

阅读：19232|2023-10-27

2 可心卡盟:win10系统火狐flash插件崩溃怎么

今天小编告诉大家如何处理win10系统火狐flash插件总是崩溃的问题，可能很多用户都不知

阅读：9998|2022-11-06

3 亲亲特价:怎么删除回收站图标

今天小编告诉大家如何对win10系统删除桌面回收站图标进行设置，可能很多用户都不知道

阅读：8331|2022-11-06

4 济南大学虚拟社区:鲁大师节能降温的具体办

今天小编告诉大家如何对win10系统电脑设置节能降温的设置方法，想必大家都遇到过需要

阅读：8701|2022-11-06

5 xlueops.exe:无线网络安装向导

我们在使用xp系统的过程中,经常需要对xp系统无线网络安装向导设置进行设置，可能很多

阅读：8646|2022-11-06

6 女斗合众国:win7系统cf与主机连接不稳定怎

今天小编告诉大家如何处理win7系统玩cf老是与主机连接不稳定的问题，可能很多用户都不

阅读：9667|2022-11-06

7 0xc000022-[cf烟雾头]cf怎么调烟雾头

电脑对日常生活的重要性小编就不多说了，可是一旦碰到win7系统设置cf烟雾头的问题，很

阅读：8631|2022-11-06

8 qizideyouhuo:应用程序无法正常启动0xc0000

我们在日常使用电脑的时候，有的小伙伴们可能在打开应用的时候会遇见提示应用程序无法

阅读：8005|2022-11-06

9 ipz-185:win7系统vcf文件怎么打开

今天小编告诉大家如何对win7系统打开vcf文件进行设置，可能很多用户都不知道怎么对win

阅读：8665|2022-11-06

10 傻哥蹦迪:win10系统s4怎么打开usb调试

今天小编告诉大家如何对win10系统s4开启USB调试模式进行设置，可能很多用户都不知道怎

阅读：7540|2022-11-06

客服电话

电子邮件

C++ KSPSetOperators函数代码示例

示例1: StokesSetupPC

示例2: ComputeKSPFETIDP

示例3: Solver

示例4: main

示例5: PressurePoissonCreate

示例6: MatAssemblyBegin

示例7: main

示例8: Petsc_KSPSolve

示例9: KSPCreate

示例10: main

示例11: TaylorGalerkinStepIIMomentum

示例12: PCSetUp_MG

示例13: main

示例14: main

示例15: internal_solve

示例16: init_fft2d_

示例17: main

示例18: SNESSolve_NEWTONTR

请发表评论

全部评论

上一篇：

下一篇：

kojino/Harvard-Robust-Machine-Learning:

scotch-io/simple-laravel-crud: Code for

cescoffier/puppet-nexus: A Puppet Module

文的笔顺,诠释文的笔画,解读文的部首

Cassolotl/sentient.cloud · GitHub

剪的笔顺,诠释剪的笔画,认识剪的部首

六六分期app的软件客服如何联系？(六六分期

florent37/ViewAnimator: A fluent Android

florent37/Shrine-MaterialDesign2: implem

CVE-2020-36276

SimpleSoftwareIO/simple-sms: Send and re

关于我们

产品与服务

解决方案

139-2527-9053