[petsc-users] fieldsplit doesn't pass prefix to inner ksp
anton
popov at uni-mainz.de
Fri Sep 26 10:17:34 CDT 2014
This is what I get for 3.5.0:
anton at anton ~/LIB/petsc-3.5.0-deb/src/snes/examples/tutorials $ ./ex19
-bf_pc_type fieldsplit -bf_snes_view
lid velocity = 0.0625, prandtl # = 1, grashof # = 1
Number of SNES iterations = 2
WARNING! There are options you set that were not used!
WARNING! could be spelling mistake, etc!
Option left: name:-bf_pc_type value: fieldsplit
Option left: name:-bf_snes_view (no value)
It seems like it was already corrected between 3.5.0 & 3.5.2
On 09/26/2014 04:26 PM, Matthew Knepley wrote:
> Here is the result for 3.5.2, which looks right to me:
>
> (v3.5.2)
> *:/PETSc3/petsc/release-petsc-3.5.1/src/snes/examples/tutorials$
> ./ex19 -bf_pc_type fieldsplit -bf_snes_view -bf_pc_fieldsplit_type
> schur -bf_pc_fieldsplit_0_fields 0,1,2 -bf_pc_fieldsplit_1_fields 3
> -bf_pc_fieldsplit_schur_factorization_type upper
> lid velocity = 0.0625, prandtl # = 1, grashof # = 1
> SNES Object:(bf_) 1 MPI processes
> type: newtonls
> maximum iterations=50, maximum function evaluations=10000
> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
> total number of linear solver iterations=4
> total number of function evaluations=3
> SNESLineSearch Object: (bf_) 1 MPI processes
> type: bt
> interpolation: cubic
> alpha=1.000000e-04
> maxstep=1.000000e+08, minlambda=1.000000e-12
> tolerances: relative=1.000000e-08, absolute=1.000000e-15,
> lambda=1.000000e-08
> maximum iterations=40
> KSP Object: (bf_) 1 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: (bf_) 1 MPI processes
> type: fieldsplit
> FieldSplit with Schur preconditioner, factorization UPPER
> Preconditioner for the Schur complement formed from A11
> Split info:
> Split number 0 Defined by IS
> Split number 1 Defined by IS
> KSP solver for A00 block
> KSP Object: (bf_fieldsplit_0_) 1 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: (bf_fieldsplit_0_) 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=48, cols=48
> package used to perform factorization: petsc
> total: nonzeros=576, allocated nonzeros=576
> total number of mallocs used during MatSetValues
> calls =0
> using I-node routines: found 16 nodes, limit used is 5
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_0_) 1 MPI processes
> type: seqaij
> rows=48, cols=48
> total: nonzeros=576, allocated nonzeros=576
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 16 nodes, limit used is 5
> KSP solver for S = A11 - A10 inv(A00) A01
> KSP Object: (bf_fieldsplit_temperature_) 1 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: (bf_fieldsplit_temperature_) 1 MPI
> processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=16
> package used to perform factorization: petsc
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> linear system matrix followed by preconditioner matrix:
> Mat Object: (bf_fieldsplit_temperature_) 1 MPI
> processes
> type: schurcomplement
> rows=16, cols=16
> Schur complement A11 - A10 inv(A00) A01
> A11
> Mat Object: (bf_fieldsplit_temperature_)
> 1 MPI processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> A10
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=48
> total: nonzeros=192, allocated nonzeros=192
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> KSP of A00
> KSP Object: (bf_fieldsplit_0_) 1 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: (bf_fieldsplit_0_) 1 MPI
> processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero
> pivot [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=48, cols=48
> package used to perform factorization: petsc
> total: nonzeros=576, allocated nonzeros=576
> total number of mallocs used during
> MatSetValues calls =0
> using I-node routines: found 16 nodes,
> limit used is 5
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_0_) 1
> MPI processes
> type: seqaij
> rows=48, cols=48
> total: nonzeros=576, allocated nonzeros=576
> total number of mallocs used during MatSetValues
> calls =0
> using I-node routines: found 16 nodes, limit
> used is 5
> A01
> Mat Object: 1 MPI processes
> type: seqaij
> rows=48, cols=16
> total: nonzeros=192, allocated nonzeros=192
> total number of mallocs used during MatSetValues
> calls =0
> using I-node routines: found 16 nodes, limit used is 5
> Mat Object: (bf_fieldsplit_temperature_) 1 MPI
> processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=64, cols=64, bs=4
> total: nonzeros=1024, allocated nonzeros=1024
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 16 nodes, limit used is 5
> Number of SNES iterations = 2
>
> Thanks,
>
> Matt
>
>
> On Fri, Sep 26, 2014 at 8:52 AM, Matthew Knepley <knepley at gmail.com
> <mailto:knepley at gmail.com>> wrote:
>
> On Fri, Sep 26, 2014 at 7:29 AM, anton <popov at uni-mainz.de
> <mailto:popov at uni-mainz.de>> wrote:
>
> Create preconditioner:
>
> PCCreate(PETSC_COMM_WORLD, &pc);
> PCSetOptionsPrefix(pc, "bf_");
> PCSetFromOptions(pc);
>
> Define fieldsplit options:
>
> -bf_pc_type fieldsplit
> -bf_pc_fieldsplit_type SCHUR
> -bf_pc_fieldsplit_schur_factorization_type UPPER
>
> Works OK.
>
> Set options for the first field solver:
>
> -bf_fieldsplit_0_ksp_type preonly
> -bf_fieldsplit_0_pc_type lu
>
> Doesn't work (ignored), because "bf_" prefix isn't pass to
> inner solver ksp (checked in the debugger).
>
> Indeed, the following works:
>
> -fieldsplit_0_ksp_type preonly
> -fieldsplit_0_pc_type lu
>
> Observed with 3.5 but not with 3.4
>
>
> I just tried this with master on SNES ex19, and got the correct
> result:
>
> knepley/feature-parallel-partition
> *$:/PETSc3/petsc/petsc-dev/src/snes/examples/tutorials$ ./ex19
> -bf_pc_type fieldsplit -bf_snes_view
> ./ex19 -bf_pc_type fieldsplit -bf_snes_view
> lid velocity = 0.0625, prandtl # = 1, grashof # = 1
> SNES Object:(bf_) 1 MPI processes
> type: newtonls
> maximum iterations=50, maximum function evaluations=10000
> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
> total number of linear solver iterations=13
> total number of function evaluations=3
> SNESLineSearch Object: (bf_) 1 MPI processes
> type: bt
> interpolation: cubic
> alpha=1.000000e-04
> maxstep=1.000000e+08, minlambda=1.000000e-12
> tolerances: relative=1.000000e-08, absolute=1.000000e-15,
> lambda=1.000000e-08
> maximum iterations=40
> KSP Object: (bf_) 1 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: (bf_) 1 MPI processes
> type: fieldsplit
> FieldSplit with MULTIPLICATIVE composition: total splits = 4
> Solver info for each split is in the following KSP objects:
> Split number 0 Defined by IS
> KSP Object: (bf_fieldsplit_x_velocity_) 1 MPI processes
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> using NONE norm type for convergence test
> PC Object: (bf_fieldsplit_x_velocity_) 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=16
> package used to perform factorization: petsc
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_x_velocity_) 1 MPI
> processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> Split number 1 Defined by IS
> KSP Object: (bf_fieldsplit_y_velocity_) 1 MPI processes
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> using NONE norm type for convergence test
> PC Object: (bf_fieldsplit_y_velocity_) 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=16
> package used to perform factorization: petsc
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_y_velocity_) 1 MPI
> processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> Split number 2 Defined by IS
> KSP Object: (bf_fieldsplit_Omega_) 1 MPI processes
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> using NONE norm type for convergence test
> PC Object: (bf_fieldsplit_Omega_) 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=16
> package used to perform factorization: petsc
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_Omega_) 1 MPI processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> Split number 3 Defined by IS
> KSP Object: (bf_fieldsplit_temperature_) 1 MPI processes
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> using NONE norm type for convergence test
> PC Object: (bf_fieldsplit_temperature_) 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot
> [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=16, cols=16
> package used to perform factorization: petsc
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues
> calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: (bf_fieldsplit_temperature_) 1 MPI
> processes
> type: seqaij
> rows=16, cols=16
> total: nonzeros=64, allocated nonzeros=64
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> linear system matrix = precond matrix:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=64, cols=64, bs=4
> total: nonzeros=1024, allocated nonzeros=1024
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 16 nodes, limit used is 5
> Number of SNES iterations = 2
>
> I will try with 3.5.2.
>
> Thanks,
>
> Matt
>
> Thanks.
> Anton
>
>
>
>
> --
> 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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140926/29f10a8c/attachment-0001.html>
More information about the petsc-users
mailing list