On Tue, Mar 18, 2014 at 2:28 AM, TAY wee-beng <[email protected]> wrote:
> On 18/3/2014 3:07 PM, Matthew Knepley wrote: > > On Tue, Mar 18, 2014 at 1:58 AM, TAY wee-beng <[email protected]> wrote: > >> Hi Barry, >> >> My command line is : >> >> mpiexec -np 46 ./a.out -options_left -poisson_ksp_view -poisson_ksp_type >> gmres -poisson_pc_type hypre -poisson_pc_type_hypre boomeramg > log >> > > > > ^^^^^^^^^^^^^^^^^^^^^^^^^^^ > > It should be -pisson_pc_hypre_type boomeramg, but BoomerAMG is the > default. > > So in other words, -poisson_pc_type hypre is sufficient, is that so? > > Also if only my RHS of the Poisson eqn's changes at every timestep, do I > set the ksp and pc only once at the start? E.g. : > Yes, only at the start. Matt > call KSPSetType(ksp,ksptype,ierr) > > ksptype=KSPGMRES > > call PCSetType(pc,'hypre',ierr) > > call PCHYPREGetType(pc,'boomeramg',ierr) > > or do I have to do it at every time step? > > > Matt > > My result is : >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> >> *KSP Object:(poisson_) 46 MPI processes type: gmres GMRES: >> restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization >> with no iterative refinement GMRES: happy breakdown tolerance 1e-30 >> maximum iterations=10000, initial guess is zero tolerances: >> relative=1e-05, absolute=1e-50, divergence=10000 left preconditioning >> using PRECONDITIONED norm type for convergence test PC Object:(poisson_) 46 >> MPI processes type: hypre HYPRE BoomerAMG preconditioning HYPRE >> BoomerAMG: Cycle type V HYPRE BoomerAMG: Maximum number of levels 25 >> HYPRE BoomerAMG: Maximum number of iterations PER hypre call 1 >> HYPRE BoomerAMG: Convergence tolerance PER hypre call 0 HYPRE >> BoomerAMG: Threshold for strong coupling 0.25 HYPRE BoomerAMG: >> Interpolation truncation factor 0 HYPRE BoomerAMG: Interpolation: max >> elements per row 0 HYPRE BoomerAMG: Number of levels of aggressive >> coarsening 0 HYPRE BoomerAMG: Number of paths for aggressive coarsening >> 1 HYPRE BoomerAMG: Maximum row sums 0.9 HYPRE BoomerAMG: Sweeps >> down 1 HYPRE BoomerAMG: Sweeps up 1 HYPRE >> BoomerAMG: Sweeps on coarse 1 HYPRE BoomerAMG: Relax down >> symmetric-SOR/Jacobi HYPRE BoomerAMG: Relax up >> symmetric-SOR/Jacobi HYPRE BoomerAMG: Relax on coarse >> Gaussian-elimination HYPRE BoomerAMG: Relax weight (all) 1 >> HYPRE BoomerAMG: Outer relax weight (all) 1 HYPRE BoomerAMG: Using >> CF-relaxation HYPRE BoomerAMG: Measure type local HYPRE >> BoomerAMG: Coarsen type Falgout HYPRE BoomerAMG: Interpolation >> type classical linear system matrix = precond matrix: Matrix Object: >> 46 MPI processes .... ... * >> >> >> >> >> >> >> *#PETSc Option Table entries: -options_left -poisson_ksp_type gmres >> -poisson_pc_type hypre -poisson_pc_type_hypre boomeramg #End of PETSc >> Option Table entries There is one unused database option. It is: Option >> left: name:-poisson_pc_type_hypre value: boomeramg* >> >> It seems that it is using boomeramg but why does it say "one unused >> database option"? >> >> Did I do something wrong? >> >> Also if only my RHS of the Poisson eqn's changes, do I set the ksp and pc >> once at the start? E.g. : >> >> call KSPSetType(ksp,ksptype,ierr) >> >> ksptype=KSPGMRES >> >> call PCSetType(pc,'hypre',ierr) >> >> call PCHYPREGetType(pc,'boomeramg',ierr) >> >> or do I have to do it at each time step? >> >> Thank you >> >> Yours sincerely, >> >> TAY wee-beng >> >> On 18/3/2014 2:23 AM, Barry Smith wrote: >> >> Yes. You can run with -poisson_ksp_view and -options_left to make sure the >> options you provide are actually used. >> >> Barry >> >> On Mar 17, 2014, at 9:39 AM, TAY wee-beng <[email protected]> >> <[email protected]> wrote: >> >> >> Hi, >> >> I use >> >> call KSPSetOptionsPrefix(ksp_semi_xyz,"momentum_",ierr) >> >> and >> >> call KSPSetOptionsPrefix(ksp,"poisson_",ierr) >> >> so that I can choose separate ksp/pc options for my momentum and poisson >> equations through command line e.g. >> >> -poisson_ksp_type gmres -poisson_pc_type hypre -poisson_pc_type_hypre >> boomeramg >> >> In general, I need to use boomeramg as the preconditioner and gmres as the >> solver for my poisson eqn, separate from my momentum eqn, which has its own >> default pc and ksp. Is the above the correct way? >> >> Thanks! >> >> -- >> Yours sincerely, >> >> TAY wee-beng >> >> >> >> > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
