[petsc-users] fieldsplit_0_ monitor in combination with selfp
Matthew Knepley
knepley at gmail.com
Fri Sep 5 10:38:37 CDT 2014
On Fri, Sep 5, 2014 at 7:31 AM, Klaij, Christiaan <C.Klaij at marin.nl> wrote:
> Thanks! I've spotted another difference: you are setting the
> fieldsplit_0_ksp_type and I'm not, just relying on the default
> instead. If I add -fieldsplit_0_ksp_type gmres then is also get
> the correct answer. Probably, you will get my problem if you
> remove -fieldsplit_velocity.
>
This is not a bug. The default solver for A00 is preonly, unless it is used
as the inner
solver as well, in which case it defaults to GMRES so as not to give an
inexact
Schur complement by default.
ex62 -run_type full -refinement_limit 0.00625 -bc_type dirichlet
-interpolate 1 -vel_petscspace_order 2 -pres_petscspace_order 1 -ksp_type
fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit
-pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full
-fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
-snes_monitor_short -ksp_monitor_short -snes_converged_reason -snes_view
-show_solution 0 -fieldsplit_pressure_inner_ksp_type gmres
-fieldsplit_pressure_inner_ksp_max_it 1 -fieldsplit_pressure_inner_pc_type
jacobi -pc_fieldsplit_schur_precondition selfp
SNES Object: 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=77
total number of function evaluations=2
SNESLineSearch Object: 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: 1 MPI processes
type: fgmres
GMRES: restart=100, 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-09, absolute=1e-50, divergence=10000
right preconditioning
has attached null space
using UNPRECONDITIONED norm type for convergence test
PC Object: 1 MPI processes
type: fieldsplit
FieldSplit with Schur preconditioner, factorization FULL
Preconditioner for the Schur complement formed from Sp, an assembled
approximation to S, which uses (the lumped) A00's diagonal's inverse
Split info:
Split number 0 Defined by IS
Split number 1 Defined by IS
KSP solver for A00 block
KSP Object: (fieldsplit_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: (fieldsplit_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=962, cols=962
package used to perform factorization: petsc
total: nonzeros=19908, allocated nonzeros=19908
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 481 nodes, limit used is 5
linear system matrix = precond matrix:
Mat Object: (fieldsplit_velocity_) 1 MPI
processes
type: seqaij
rows=962, cols=962
total: nonzeros=19908, allocated nonzeros=19908
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 481 nodes, limit used is 5
KSP solver for S = A11 - A10 inv(A00) A01
KSP Object: (fieldsplit_pressure_) 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-10, absolute=1e-50, divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: (fieldsplit_pressure_) 1 MPI processes
type: jacobi
linear system matrix followed by preconditioner matrix:
Mat Object: (fieldsplit_pressure_) 1 MPI
processes
type: schurcomplement
rows=145, cols=145
has attached null space
Schur complement A11 - A10 inv(A00) A01
A11
Mat Object: (fieldsplit_pressure_)
1 MPI processes
type: seqaij
rows=145, cols=145
total: nonzeros=945, allocated nonzeros=945
total number of mallocs used during MatSetValues calls =0
has attached null space
not using I-node routines
A10
Mat Object: 1 MPI processes
type: seqaij
rows=145, cols=962
total: nonzeros=4466, allocated nonzeros=4466
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP of A00
KSP Object: (fieldsplit_pressure_inner_)
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=1, initial guess is zero
tolerances: relative=1e-09, absolute=1e-50,
divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: (fieldsplit_pressure_inner_)
1 MPI processes
type: jacobi
linear system matrix = precond matrix:
Mat Object: (fieldsplit_velocity_)
1 MPI processes
type: seqaij
rows=962, cols=962
total: nonzeros=19908, allocated nonzeros=19908
total number of mallocs used during MatSetValues calls
=0
using I-node routines: found 481 nodes, limit used is
5
A01
Mat Object: 1 MPI processes
type: seqaij
rows=962, cols=145
total: nonzeros=4466, allocated nonzeros=4466
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 481 nodes, limit used is 5
Mat Object: 1 MPI processes
type: seqaij
rows=145, cols=145
total: nonzeros=2601, allocated nonzeros=2601
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=1107, cols=1107
total: nonzeros=29785, allocated nonzeros=29785
total number of mallocs used during MatSetValues calls =0
has attached null space
using I-node routines: found 513 nodes, limit used is 5
Matt
> mpiexec -n 2 ./ex70 -nx 4 -ny 6 \
> -ksp_type fgmres \
> -pc_type fieldsplit \
> -pc_fieldsplit_type schur \
> -pc_fieldsplit_schur_fact_type lower \
> -pc_fieldsplit_schur_precondition selfp \
> -fieldsplit_1_inner_ksp_type preonly \
> -fieldsplit_1_inner_pc_type jacobi \
> -fieldsplit_0_ksp_monitor -fieldsplit_0_ksp_max_it 1 \
> -fieldsplit_1_ksp_monitor -fieldsplit_1_ksp_max_it 1 \
> -ksp_monitor -ksp_max_it 1 \
> -fieldsplit_0_ksp_type gmres -ksp_view
>
>
> KSP Object: 2 MPI processes
> type: fgmres
> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
> Orthogonalization with no iterative refinement
> GMRES: happy breakdown tolerance 1e-30
> maximum iterations=1, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> right preconditioning
> using UNPRECONDITIONED norm type for convergence test
> PC Object: 2 MPI processes
> type: fieldsplit
> FieldSplit with Schur preconditioner, factorization LOWER
> Preconditioner for the Schur complement formed from Sp, an assembled
> approximation to S, which uses (the lumped) A00's diagonal's inverse
> Split info:
> Split number 0 Defined by IS
> Split number 1 Defined by IS
> KSP solver for A00 block
> KSP Object: (fieldsplit_0_) 2 MPI processes
> type: gmres
>
>
>
>
>
>
> MARIN news: Development of a Scaled-Down Floating Wind Turbine for
> Offshore Basin Testing
> <http://www.marin.nl/web/News/News-items/Development-of-a-ScaledDown-Floating-Wind-Turbine-for-Offshore-Basin-Testing-1.htm>
>
> This e-mail may be confidential, privileged and/or protected by copyright.
> If you are not the intended recipient, you should return it to the sender
> immediately and delete your copy from your system.
>
>
>
> ------------------------------
> *From:* Matthew Knepley <knepley at gmail.com>
> *Sent:* Friday, September 05, 2014 2:10 PM
> *To:* Klaij, Christiaan; petsc-users at mcs.anl.gov
> *Subject:* Re: [petsc-users] fieldsplit_0_ monitor in combination with
> selfp
>
> On Fri, Sep 5, 2014 at 1:34 AM, Klaij, Christiaan <C.Klaij at marin.nl>
> wrote:
>
>> Matt,
>>
>> I think the problem is somehow related to
>> -pc_fieldsplit_schur_precondition selfp. In the example below your are not
>> using that option.
>>
>
> Here is the selfp output. It retains the A00 solver.
>
> ex62 -run_type full -refinement_limit 0.00625 -bc_type dirichlet
> -interpolate 1 -vel_petscspace_order 2 -pres_petscspace_order 1 -ksp_type
> fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit
> -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full
> -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres
> -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi
> -snes_monitor_short -ksp_monitor_short -snes_converged_reason
> -ksp_converged_reason -snes_view -show_solution 0
> -fieldsplit_pressure_inner_ksp_type preonly
> -fieldsplit_pressure_inner_pc_type jacobi -pc_fieldsplit_schur_precondition
> selfp
>
> SNES Object: 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=20
> total number of function evaluations=2
> SNESLineSearch Object: 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: 1 MPI processes
> type: fgmres
> GMRES: restart=100, 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-09, absolute=1e-50, divergence=10000
> right preconditioning
> has attached null space
> using UNPRECONDITIONED norm type for convergence test
> PC Object: 1 MPI processes
> type: fieldsplit
> FieldSplit with Schur preconditioner, factorization FULL
> Preconditioner for the Schur complement formed from Sp, an assembled
> approximation to S, which uses (the lumped) A00's diagonal's inverse
> Split info:
> Split number 0 Defined by IS
> Split number 1 Defined by IS
> KSP solver for A00 block
> KSP Object: (fieldsplit_velocity_) 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: (fieldsplit_velocity_) 1 MPI processes
> type: lu
> LU: out-of-place factorization
> tolerance for zero pivot 2.22045e-14
> matrix ordering: nd
> factor fill ratio given 5, needed 3.45047
> Factored matrix follows:
> Mat Object: 1 MPI processes
> type: seqaij
> rows=962, cols=962
> package used to perform factorization: petsc
> total: nonzeros=68692, allocated nonzeros=68692
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 456 nodes, limit used is 5
> linear system matrix = precond matrix:
> Mat Object: (fieldsplit_velocity_) 1 MPI
> processes
> type: seqaij
> rows=962, cols=962
> total: nonzeros=19908, allocated nonzeros=19908
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 481 nodes, limit used is 5
> KSP solver for S = A11 - A10 inv(A00) A01
> KSP Object: (fieldsplit_pressure_) 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-10, absolute=1e-50, divergence=10000
> left preconditioning
> using PRECONDITIONED norm type for convergence test
> PC Object: (fieldsplit_pressure_) 1 MPI processes
> type: jacobi
> linear system matrix followed by preconditioner matrix:
> Mat Object: (fieldsplit_pressure_) 1 MPI
> processes
> type: schurcomplement
> rows=145, cols=145
> has attached null space
> Schur complement A11 - A10 inv(A00) A01
> A11
> Mat Object: (fieldsplit_pressure_)
> 1 MPI processes
> type: seqaij
> rows=145, cols=145
> total: nonzeros=945, allocated nonzeros=945
> total number of mallocs used during MatSetValues calls =0
> has attached null space
> not using I-node routines
> A10
> Mat Object: 1 MPI processes
> type: seqaij
> rows=145, cols=962
> total: nonzeros=4466, allocated nonzeros=4466
> total number of mallocs used during MatSetValues calls =0
> not using I-node routines
> KSP of A00
> KSP Object: (fieldsplit_pressure_inner_)
> 1 MPI processes
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-09, absolute=1e-50,
> divergence=10000
> left preconditioning
> using NONE norm type for convergence test
> PC Object: (fieldsplit_pressure_inner_)
> 1 MPI processes
> type: jacobi
> linear system matrix = precond matrix:
> Mat Object: (fieldsplit_velocity_)
> 1 MPI processes
> type: seqaij
> rows=962, cols=962
> total: nonzeros=19908, allocated nonzeros=19908
> total number of mallocs used during MatSetValues calls
> =0
> using I-node routines: found 481 nodes, limit used
> is 5
> A01
> Mat Object: 1 MPI processes
> type: seqaij
> rows=962, cols=145
> total: nonzeros=4466, allocated nonzeros=4466
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 481 nodes, limit used is 5
> Mat Object: 1 MPI processes
> type: seqaij
> rows=145, cols=145
> total: nonzeros=2601, allocated nonzeros=2601
> 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=1107, cols=1107
> total: nonzeros=29785, allocated nonzeros=29785
> total number of mallocs used during MatSetValues calls =0
> has attached null space
> using I-node routines: found 513 nodes, limit used is 5
>
> Thanks,
>
> Matt
>
>
>> Chris
>>
>> dr. ir. Christiaan Klaij
>>
>> CFD Researcher
>> Research & Development
>>
>>
>>
>> *MARIN*
>>
>>
>> 2, Haagsteeg E C.Klaij at marin.nl P.O. Box 28 T +31 317 49 39 11 6700
>> AA Wageningen F +31 317 49 32 45 T +31 317 49 33 44 The Netherlands I
>> www.marin.nl
>>
>>
>>
>> MARIN news: MARIN at SMM, Hamburg, September 9-12
>> <http://www.marin.nl/web/News/News-items/MARIN-at-SMM-Hamburg-September-912.htm>
>>
>> This e-mail may be confidential, privileged and/or protected by
>> copyright. If you are not the intended recipient, you should return it to
>> the sender immediately and delete your copy from your system.
>>
>>
>>
>> ------------------------------
>> *From:* Matthew Knepley <knepley at gmail.com>
>> *Sent:* Friday, September 05, 2014 12:36 AM
>> *To:* Klaij, Christiaan
>> *Cc:* petsc-users at mcs.anl.gov
>> *Subject:* Re: [petsc-users] fieldsplit_0_ monitor in combination with
>> selfp
>>
>> On Thu, Sep 4, 2014 at 7:26 AM, Klaij, Christiaan <C.Klaij at marin.nl>
>> wrote:
>>
>>> Sorry, here's the ksp_view. I'm expecting
>>>
>>> -fieldsplit_1_inner_ksp_type preonly
>>>
>>> to set the ksp(A00) in the Schur complement only, but it seems to set it
>>> in the inv(A00) of the diagonal as well.
>>>
>>
>> I think something is wrong in your example (we strongly advise against
>> using MatNest directly). I cannot reproduce this using SNES ex62:
>>
>> ./config/builder2.py check src/snes/examples/tutorials/ex62.c
>> --testnum=36 --args="-fieldsplit_pressure_inner_ksp_type preonly
>> -fieldsplit_pressure_inner_pc_type jacobi"
>>
>> which translates to
>>
>> ex62 -run_type full -refinement_limit 0.00625 -bc_type dirichlet
>> -interpolate 1 -vel_petscspace_order 2 -pres_petscspace_order 1 -ksp_type
>> fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit
>> -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full
>> -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres
>> -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi
>> -snes_monitor_short -ksp_monitor_short -snes_converged_reason
>> -ksp_converged_reason -snes_view -show_solution 0
>> -fieldsplit_pressure_inner_ksp_type preonly
>> -fieldsplit_pressure_inner_pc_type jacobi
>>
>> gives
>>
>> Nonlinear solve converged due to CONVERGED_FNORM_RELATIVE iterations 1
>> SNES Object: 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=20
>> total number of function evaluations=2
>> SNESLineSearch Object: 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: 1 MPI processes
>> type: fgmres
>> GMRES: restart=100, 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-09, absolute=1e-50, divergence=10000
>> right preconditioning
>> has attached null space
>> using UNPRECONDITIONED norm type for convergence test
>> PC Object: 1 MPI processes
>> type: fieldsplit
>> FieldSplit with Schur preconditioner, factorization FULL
>> 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: (fieldsplit_velocity_) 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: (fieldsplit_velocity_) 1 MPI processes
>> type: lu
>> LU: out-of-place factorization
>> tolerance for zero pivot 2.22045e-14
>> matrix ordering: nd
>> factor fill ratio given 5, needed 3.45047
>> Factored matrix follows:
>> Mat Object: 1 MPI processes
>> type: seqaij
>> rows=962, cols=962
>> package used to perform factorization: petsc
>> total: nonzeros=68692, allocated nonzeros=68692
>> total number of mallocs used during MatSetValues calls
>> =0
>> using I-node routines: found 456 nodes, limit used is
>> 5
>> linear system matrix = precond matrix:
>> Mat Object: (fieldsplit_velocity_) 1 MPI
>> processes
>> type: seqaij
>> rows=962, cols=962
>> total: nonzeros=19908, allocated nonzeros=19908
>> total number of mallocs used during MatSetValues calls =0
>> using I-node routines: found 481 nodes, limit used is 5
>> KSP solver for S = A11 - A10 inv(A00) A01
>> KSP Object: (fieldsplit_pressure_) 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-10, absolute=1e-50, divergence=10000
>> left preconditioning
>> has attached null space
>> using PRECONDITIONED norm type for convergence test
>> PC Object: (fieldsplit_pressure_) 1 MPI processes
>> type: jacobi
>> linear system matrix followed by preconditioner matrix:
>> Mat Object: (fieldsplit_pressure_) 1 MPI
>> processes
>> type: schurcomplement
>> rows=145, cols=145
>> has attached null space
>> Schur complement A11 - A10 inv(A00) A01
>> A11
>> Mat Object: (fieldsplit_pressure_)
>> 1 MPI processes
>> type: seqaij
>> rows=145, cols=145
>> total: nonzeros=945, allocated nonzeros=945
>> total number of mallocs used during MatSetValues calls
>> =0
>> has attached null space
>> not using I-node routines
>> A10
>> Mat Object: 1 MPI processes
>> type: seqaij
>> rows=145, cols=962
>> total: nonzeros=4466, allocated nonzeros=4466
>> total number of mallocs used during MatSetValues calls
>> =0
>> not using I-node routines
>> KSP of A00
>> KSP Object: (fieldsplit_pressure_inner_)
>> 1 MPI processes
>> type: preonly
>> maximum iterations=10000, initial guess is zero
>> tolerances: relative=1e-09, absolute=1e-50,
>> divergence=10000
>> left preconditioning
>> using NONE norm type for convergence test
>> PC Object: (fieldsplit_pressure_inner_)
>> 1 MPI processes
>> type: jacobi
>> linear system matrix = precond matrix:
>> Mat Object: (fieldsplit_velocity_)
>> 1 MPI processes
>> type: seqaij
>> rows=962, cols=962
>> total: nonzeros=19908, allocated nonzeros=19908
>> total number of mallocs used during MatSetValues
>> calls =0
>> using I-node routines: found 481 nodes, limit used
>> is 5
>> A01
>> Mat Object: 1 MPI processes
>> type: seqaij
>> rows=962, cols=145
>> total: nonzeros=4466, allocated nonzeros=4466
>> total number of mallocs used during MatSetValues calls
>> =0
>> using I-node routines: found 481 nodes, limit used is
>> 5
>> Mat Object: (fieldsplit_pressure_) 1 MPI
>> processes
>> type: seqaij
>> rows=145, cols=145
>> total: nonzeros=945, allocated nonzeros=945
>> total number of mallocs used during MatSetValues calls =0
>> has attached null space
>> not using I-node routines
>> linear system matrix = precond matrix:
>> Mat Object: 1 MPI processes
>> type: seqaij
>> rows=1107, cols=1107
>> total: nonzeros=29785, allocated nonzeros=29785
>> total number of mallocs used during MatSetValues calls =0
>> has attached null space
>> using I-node routines: found 513 nodes, limit used is 5
>>
>> Matt
>>
>>
>>> Chris
>>>
>>> 0 KSP Residual norm 1.229687498638e+00
>>> Residual norms for fieldsplit_1_ solve.
>>> 0 KSP Residual norm 7.185799114488e+01
>>> 1 KSP Residual norm 3.873274154012e+01
>>> 1 KSP Residual norm 1.107969383366e+00
>>> KSP Object: 1 MPI processes
>>> type: fgmres
>>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
>>> Orthogonalization with no iterative refinement
>>> GMRES: happy breakdown tolerance 1e-30
>>> maximum iterations=1, initial guess is zero
>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>> right preconditioning
>>> using UNPRECONDITIONED norm type for convergence test
>>> PC Object: 1 MPI processes
>>> type: fieldsplit
>>> FieldSplit with Schur preconditioner, factorization LOWER
>>> Preconditioner for the Schur complement formed from Sp, an assembled
>>> approximation to S, which uses (the lumped) A00's diagonal's inverse
>>> Split info:
>>> Split number 0 Defined by IS
>>> Split number 1 Defined by IS
>>> KSP solver for A00 block
>>> KSP Object: (fieldsplit_0_) 1 MPI processes
>>> type: preonly
>>> maximum iterations=1, initial guess is zero
>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>> left preconditioning
>>> using NONE norm type for convergence test
>>> PC Object: (fieldsplit_0_) 1 MPI processes
>>> type: bjacobi
>>> block Jacobi: number of blocks = 1
>>> Local solve is same for all blocks, in the following KSP and
>>> PC objects:
>>> KSP Object: (fieldsplit_0_sub_) 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: (fieldsplit_0_sub_) 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=200, allocated nonzeros=200
>>> total number of mallocs used during MatSetValues
>>> calls =0
>>> not using I-node routines
>>> linear system matrix = precond matrix:
>>> Mat Object: (fieldsplit_0_) 1 MPI
>>> processes
>>> type: seqaij
>>> rows=48, cols=48
>>> total: nonzeros=200, allocated nonzeros=240
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node routines
>>> linear system matrix = precond matrix:
>>> Mat Object: (fieldsplit_0_) 1 MPI processes
>>> type: mpiaij
>>> rows=48, cols=48
>>> total: nonzeros=200, allocated nonzeros=480
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node (on process 0) routines
>>> KSP solver for S = A11 - A10 inv(A00) A01
>>> KSP Object: (fieldsplit_1_) 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=1, initial guess is zero
>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>> left preconditioning
>>> using PRECONDITIONED norm type for convergence test
>>> PC Object: (fieldsplit_1_) 1 MPI processes
>>> type: bjacobi
>>> block Jacobi: number of blocks = 1
>>> Local solve is same for all blocks, in the following KSP and
>>> PC objects:
>>> KSP Object: (fieldsplit_1_sub_) 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: (fieldsplit_1_sub_) 1 MPI
>>> processes
>>> type: bjacobi
>>> block Jacobi: number of blocks = 1
>>> Local solve is same for all blocks, in the following KSP
>>> and PC objects:
>>> KSP Object:
>>> (fieldsplit_1_sub_sub_) 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:
>>> (fieldsplit_1_sub_sub_) 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=24, cols=24
>>> package used to perform factorization: petsc
>>> total: nonzeros=120, allocated nonzeros=120
>>> 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=24, cols=24
>>> total: nonzeros=120, allocated nonzeros=120
>>> 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: mpiaij
>>> rows=24, cols=24
>>> total: nonzeros=120, allocated nonzeros=120
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node (on process 0) routines
>>> linear system matrix followed by preconditioner matrix:
>>> Mat Object: (fieldsplit_1_) 1 MPI processes
>>> type: schurcomplement
>>> rows=24, cols=24
>>> Schur complement A11 - A10 inv(A00) A01
>>> A11
>>> Mat Object: (fieldsplit_1_) 1
>>> MPI processes
>>> type: mpiaij
>>> rows=24, cols=24
>>> total: nonzeros=0, allocated nonzeros=0
>>> total number of mallocs used during MatSetValues calls =0
>>> using I-node (on process 0) routines: found 5 nodes,
>>> limit used is 5
>>> A10
>>> Mat Object: (a10_) 1 MPI
>>> processes
>>> type: mpiaij
>>> rows=24, cols=48
>>> total: nonzeros=96, allocated nonzeros=96
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node (on process 0) routines
>>> KSP of A00
>>> KSP Object:
>>> (fieldsplit_1_inner_) 1 MPI processes
>>> type: preonly
>>> maximum iterations=1, initial guess is zero
>>> tolerances: relative=1e-05, absolute=1e-50,
>>> divergence=10000
>>> left preconditioning
>>> using NONE norm type for convergence test
>>> PC Object:
>>> (fieldsplit_1_inner_) 1 MPI processes
>>> type: jacobi
>>> linear system matrix = precond matrix:
>>> Mat Object:
>>> (fieldsplit_0_) 1 MPI processes
>>> type: mpiaij
>>> rows=48, cols=48
>>> total: nonzeros=200, allocated nonzeros=480
>>> total number of mallocs used during MatSetValues calls
>>> =0
>>> not using I-node (on process 0) routines
>>> A01
>>> Mat Object: (a01_) 1 MPI
>>> processes
>>> type: mpiaij
>>> rows=48, cols=24
>>> total: nonzeros=96, allocated nonzeros=480
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node (on process 0) routines
>>> Mat Object: 1 MPI processes
>>> type: mpiaij
>>> rows=24, cols=24
>>> total: nonzeros=120, allocated nonzeros=120
>>> total number of mallocs used during MatSetValues calls =0
>>> not using I-node (on process 0) routines
>>> linear system matrix = precond matrix:
>>> Mat Object: 1 MPI processes
>>> type: nest
>>> rows=72, cols=72
>>> Matrix object:
>>> type=nest, rows=2, cols=2
>>> MatNest structure:
>>> (0,0) : prefix="fieldsplit_0_", type=mpiaij, rows=48, cols=48
>>> (0,1) : prefix="a01_", type=mpiaij, rows=48, cols=24
>>> (1,0) : prefix="a10_", type=mpiaij, rows=24, cols=48
>>> (1,1) : prefix="fieldsplit_1_", type=mpiaij, rows=24, cols=24
>>>
>>>
>>> From: Matthew Knepley <knepley at gmail.com>
>>> Sent: Thursday, September 04, 2014 2:20 PM
>>> To: Klaij, Christiaan
>>> Cc: petsc-users at mcs.anl.gov
>>> Subject: Re: [petsc-users] fieldsplit_0_ monitor in combination with
>>> selfp
>>>
>>>
>>>
>>>
>>> On Thu, Sep 4, 2014 at 7:06 AM, Klaij, Christiaan <C.Klaij at marin.nl>
>>> wrote:
>>> I'm playing with the selfp option in fieldsplit using
>>> snes/examples/tutorials/ex70.c. For example:
>>>
>>> mpiexec -n 2 ./ex70 -nx 4 -ny 6 \
>>> -ksp_type fgmres \
>>> -pc_type fieldsplit \
>>> -pc_fieldsplit_type schur \
>>> -pc_fieldsplit_schur_fact_type lower \
>>> -pc_fieldsplit_schur_precondition selfp \
>>> -fieldsplit_1_inner_ksp_type preonly \
>>> -fieldsplit_1_inner_pc_type jacobi \
>>> -fieldsplit_0_ksp_monitor -fieldsplit_0_ksp_max_it 1 \
>>> -fieldsplit_1_ksp_monitor -fieldsplit_1_ksp_max_it 1 \
>>> -ksp_monitor -ksp_max_it 1
>>>
>>> gives the following output
>>>
>>> 0 KSP Residual norm 1.229687498638e+00
>>> Residual norms for fieldsplit_1_ solve.
>>> 0 KSP Residual norm 2.330138480101e+01
>>> 1 KSP Residual norm 1.609000846751e+01
>>> 1 KSP Residual norm 1.180287268335e+00
>>>
>>> To my suprise I don't see anything for the fieldsplit_0_ solve,
>>> why?
>>>
>>>
>>>
>>> Always run with -ksp_view for any solver question.
>>>
>>>
>>> Thanks,
>>>
>>>
>>> Matt
>>> Furthermore, if I understand correctly the above should be
>>> exactly equivalent with
>>>
>>> mpiexec -n 2 ./ex70 -nx 4 -ny 6 \
>>> -ksp_type fgmres \
>>> -pc_type fieldsplit \
>>> -pc_fieldsplit_type schur \
>>> -pc_fieldsplit_schur_fact_type lower \
>>> -user_ksp \
>>> -fieldsplit_0_ksp_monitor -fieldsplit_0_ksp_max_it 1 \
>>> -fieldsplit_1_ksp_monitor -fieldsplit_1_ksp_max_it 1 \
>>> -ksp_monitor -ksp_max_it 1
>>>
>>> 0 KSP Residual norm 1.229687498638e+00
>>> Residual norms for fieldsplit_0_ solve.
>>> 0 KSP Residual norm 5.486639587672e-01
>>> 1 KSP Residual norm 6.348354253703e-02
>>> Residual norms for fieldsplit_1_ solve.
>>> 0 KSP Residual norm 2.321938107977e+01
>>> 1 KSP Residual norm 1.605484031258e+01
>>> 1 KSP Residual norm 1.183225251166e+00
>>>
>>> because -user_ksp replaces the Schur complement by the simple
>>> approximation A11 - A10 inv(diag(A00)) A01. Beside the missing
>>> fielsplit_0_ part, the numbers are pretty close but not exactly
>>> the same. Any explanation?
>>>
>>> Chris
>>>
>>>
>>> dr. ir. Christiaan Klaij
>>> CFD Researcher
>>> Research & Development
>>> E mailto:C.Klaij at marin.nl
>>> T +31 317 49 33 44
>>>
>>>
>>> MARIN
>>> 2, Haagsteeg, P.O. Box 28, 6700 AA Wageningen, The Netherlands
>>> T +31 317 49 39 11, F +31 317 49 32 45, I www.marin.nl
>>>
>>>
>>>
>>>
>>>
>>> --
>>> 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
>>
>>
>
>
> --
> 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/20140905/1693045a/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image4a95d2.JPG
Type: image/jpeg
Size: 1622 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140905/1693045a/attachment-0002.jpe>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: imagea7316f.JPG
Type: image/jpeg
Size: 1069 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140905/1693045a/attachment-0003.jpe>
More information about the petsc-users
mailing list