[petsc-users] fieldsplit_0_ monitor in combination with selfp
Matthew Knepley
knepley at gmail.com
Fri Sep 5 07:10:26 CDT 2014
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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140905/eabc0118/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/eabc0118/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/eabc0118/attachment-0003.jpe>
More information about the petsc-users
mailing list