[petsc-users] how to change KSP of A00 inside the Schur complement?
Klaij, Christiaan
C.Klaij at marin.nl
Tue Sep 9 02:31:23 CDT 2014
In the program below, I'm using PCFieldSplitGetSubKSP to get the
sub KSP's of a Schur fieldsplit preconditioner. I'm setting
fieldsplit_0 to BICG+ILU and fieldsplit_1 to CG+ICC. Running
$ ./fieldsplittry -ksp_view
shows that this works as expected (full output below). Now, I
would like to change the KSP of A00 inside the Schur complement,
so I'm running
$ ./fieldsplittry -fieldsplit_1_inner_ksp_type preonly -fieldsplit_1_inner_pc_type jacobi -ksp_view
(full output below). To my surprise, this shows the
fieldsplit_1_inner_ KSP to be BICG+ILU while the fieldsplit_0_
KSP is changed to preonly; the fieldsplit_0_ PC is still ILU (no
Jacobi anywhere). What am I doing wrong this time?
The bottom line is: how do I set the KSP of A00 inside the Schur
complement to preonly+jacobi while keeping the other settings?
Preferably directly in the code without command-line arguments.
Chris
$ cat fieldsplittry.F90
program fieldsplittry
use petscksp
implicit none
#include <finclude/petsckspdef.h>
PetscErrorCode :: ierr
PetscInt :: size,i,j,start,end,n=4,numsplit=1
PetscScalar :: zero=0.0,one=1.0
Vec :: diag3,x,b
Mat :: A,subA(4),myS
PC :: pc,subpc(2)
KSP :: ksp,subksp(2)
IS :: isg(2)
call PetscInitialize(PETSC_NULL_CHARACTER,ierr); CHKERRQ(ierr)
call MPI_Comm_size(PETSC_COMM_WORLD,size,ierr); CHKERRQ(ierr);
! vectors
call VecCreateMPI(MPI_COMM_WORLD,3*n,PETSC_DECIDE,diag3,ierr); CHKERRQ(ierr)
call VecSet(diag3,one,ierr); CHKERRQ(ierr)
call VecCreateMPI(MPI_COMM_WORLD,4*n,PETSC_DECIDE,x,ierr); CHKERRQ(ierr)
call VecSet(x,zero,ierr); CHKERRQ(ierr)
call VecDuplicate(x,b,ierr); CHKERRQ(ierr)
call VecSet(b,one,ierr); CHKERRQ(ierr)
! matrix a00
call MatCreateAIJ(MPI_COMM_WORLD,3*n,3*n,PETSC_DECIDE,PETSC_DECIDE,1,PETSC_NULL_INTEGER,0,PETSC_NULL_INTEGER,subA(1),ierr);CHKERRQ(ierr)
call MatDiagonalSet(subA(1),diag3,INSERT_VALUES,ierr);CHKERRQ(ierr)
call MatAssemblyBegin(subA(1),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
call MatAssemblyEnd(subA(1),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
! matrix a01
call MatCreateAIJ(MPI_COMM_WORLD,3*n,n,PETSC_DECIDE,PETSC_DECIDE,1,PETSC_NULL_INTEGER,1,PETSC_NULL_INTEGER,subA(2),ierr);CHKERRQ(ierr)
call MatGetOwnershipRange(subA(2),start,end,ierr);CHKERRQ(ierr);
do i=start,end-1
j=mod(i,size*n)
call MatSetValue(subA(2),i,j,one,INSERT_VALUES,ierr);CHKERRQ(ierr)
end do
call MatAssemblyBegin(subA(2),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
call MatAssemblyEnd(subA(2),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
! matrix a10
call MatTranspose(subA(2),MAT_INITIAL_MATRIX,subA(3),ierr);CHKERRQ(ierr)
! matrix a11 (empty)
call MatCreateAIJ(MPI_COMM_WORLD,n,n,PETSC_DECIDE,PETSC_DECIDE,0,PETSC_NULL_INTEGER,0,PETSC_NULL_INTEGER,subA(4),ierr);CHKERRQ(ierr)
call MatAssemblyBegin(subA(4),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
call MatAssemblyEnd(subA(4),MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
! nested mat [a00,a01;a10,a11]
call MatCreateNest(MPI_COMM_WORLD,2,PETSC_NULL_OBJECT,2,PETSC_NULL_OBJECT,subA,A,ierr);CHKERRQ(ierr)
call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr);CHKERRQ(ierr)
call MatNestGetISs(A,isg,PETSC_NULL_OBJECT,ierr);CHKERRQ(ierr);
! KSP and PC
call KSPCreate(MPI_COMM_WORLD,ksp,ierr);CHKERRQ(ierr)
call KSPSetOperators(ksp,A,A,ierr);CHKERRQ(ierr)
call KSPSetType(ksp,KSPGCR,ierr);CHKERRQ(ierr)
call KSPGetPC(ksp,pc,ierr);CHKERRQ(ierr)
call PCSetType(pc,PCFIELDSPLIT,ierr);CHKERRQ(ierr)
call PCFieldSplitSetType(pc,PC_COMPOSITE_SCHUR,ierr);CHKERRQ(ierr)
call PCFieldSplitSetIS(pc,"0",isg(1),ierr);CHKERRQ(ierr)
call PCFieldSplitSetIS(pc,"1",isg(2),ierr);CHKERRQ(ierr)
call PCFieldSplitSetSchurFactType(pc,PC_FIELDSPLIT_SCHUR_FACT_LOWER,ierr);CHKERRQ(ierr)
call PCFieldSplitSetSchurPre(pc,PC_FIELDSPLIT_SCHUR_PRE_SELFP,PETSC_NULL_OBJECT,ierr);CHKERRQ(ierr)
call KSPSetUp(ksp,ierr);CHKERRQ(ierr);
call PCFieldSplitGetSubKSP(pc,numsplit,subksp,ierr);CHKERRQ(ierr)
call KSPSetType(subksp(1),KSPBICG,ierr);CHKERRQ(ierr)
call KSPGetPC(subksp(1),subpc(1),ierr);CHKERRQ(ierr)
call PCSetType(subpc(1),PCILU,ierr);CHKERRQ(ierr)
call KSPSetType(subksp(2),KSPCG,ierr);CHKERRQ(ierr)
call KSPGetPC(subksp(2),subpc(2),ierr);CHKERRQ(ierr)
call PCSetType(subpc(2),PCICC,ierr);CHKERRQ(ierr)
! call PetscFree(subksp);CHKERRQ(ierr);
call KSPSetFromOptions(ksp,ierr);CHKERRQ(ierr)
call KSPSolve(ksp,b,x,ierr);CHKERRQ(ierr)
call KSPGetSolution(ksp,x,ierr);CHKERRQ(ierr)
! call VecView(x,PETSC_VIEWER_STDOUT_WORLD,ierr);CHKERRQ(ierr)
call PetscFinalize(ierr)
end program fieldsplittry
$ ./fieldsplittry -ksp_view
KSP Object: 1 MPI processes
type: gcr
GCR: restart = 30
GCR: restarts performed = 1
maximum iterations=10000, 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: bicg
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_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=12, cols=12
package used to perform factorization: petsc
total: nonzeros=12, allocated nonzeros=12
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=12, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP solver for S = A11 - A10 inv(A00) A01
KSP Object: (fieldsplit_1_) 1 MPI processes
type: cg
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_1_) 1 MPI processes
type: icc
0 levels of fill
tolerance for zero pivot 2.22045e-14
using Manteuffel shift [POSITIVE_DEFINITE]
matrix ordering: natural
factor fill ratio given 1, needed 1
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqsbaij
rows=4, cols=4
package used to perform factorization: petsc
total: nonzeros=4, allocated nonzeros=4
total number of mallocs used during MatSetValues calls =0
block size is 1
linear system matrix followed by preconditioner matrix:
Mat Object: (fieldsplit_1_) 1 MPI processes
type: schurcomplement
rows=4, cols=4
Schur complement A11 - A10 inv(A00) A01
A11
Mat Object: (fieldsplit_1_) 1 MPI processes
type: seqaij
rows=4, cols=4
total: nonzeros=0, allocated nonzeros=0
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 1 nodes, limit used is 5
A10
Mat Object: 1 MPI processes
type: seqaij
rows=4, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP of A00
KSP Object: (fieldsplit_0_) 1 MPI processes
type: bicg
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_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=12, cols=12
package used to perform factorization: petsc
total: nonzeros=12, allocated nonzeros=12
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=12, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
A01
Mat Object: 1 MPI processes
type: seqaij
rows=12, cols=4
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
Mat Object: 1 MPI processes
type: seqaij
rows=4, cols=4
total: nonzeros=4, allocated nonzeros=4
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: nest
rows=16, cols=16
Matrix object:
type=nest, rows=2, cols=2
MatNest structure:
(0,0) : prefix="fieldsplit_0_", type=seqaij, rows=12, cols=12
(0,1) : type=seqaij, rows=12, cols=4
(1,0) : type=seqaij, rows=4, cols=12
(1,1) : prefix="fieldsplit_1_", type=seqaij, rows=4, cols=4
$ ./fieldsplittry -fieldsplit_1_inner_ksp_type preonly -fieldsplit_1_inner_pc_type jacobi -ksp_view
KSP Object: 1 MPI processes
type: gcr
GCR: restart = 30
GCR: restarts performed = 1
maximum iterations=10000, 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=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_) 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=12, cols=12
package used to perform factorization: petsc
total: nonzeros=12, allocated nonzeros=12
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=12, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP solver for S = A11 - A10 inv(A00) A01
KSP Object: (fieldsplit_1_) 1 MPI processes
type: cg
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_1_) 1 MPI processes
type: icc
0 levels of fill
tolerance for zero pivot 2.22045e-14
using Manteuffel shift [POSITIVE_DEFINITE]
matrix ordering: natural
factor fill ratio given 1, needed 1
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqsbaij
rows=4, cols=4
package used to perform factorization: petsc
total: nonzeros=4, allocated nonzeros=4
total number of mallocs used during MatSetValues calls =0
block size is 1
linear system matrix followed by preconditioner matrix:
Mat Object: (fieldsplit_1_) 1 MPI processes
type: schurcomplement
rows=4, cols=4
Schur complement A11 - A10 inv(A00) A01
A11
Mat Object: (fieldsplit_1_) 1 MPI processes
type: seqaij
rows=4, cols=4
total: nonzeros=0, allocated nonzeros=0
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 1 nodes, limit used is 5
A10
Mat Object: 1 MPI processes
type: seqaij
rows=4, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP of A00
KSP Object: (fieldsplit_1_inner_) 1 MPI processes
type: bicg
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_1_inner_) 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=12, cols=12
package used to perform factorization: petsc
total: nonzeros=12, allocated nonzeros=12
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=12, cols=12
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
A01
Mat Object: 1 MPI processes
type: seqaij
rows=12, cols=4
total: nonzeros=12, allocated nonzeros=12
total number of mallocs used during MatSetValues calls =0
not using I-node routines
Mat Object: 1 MPI processes
type: seqaij
rows=4, cols=4
total: nonzeros=4, allocated nonzeros=4
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: nest
rows=16, cols=16
Matrix object:
type=nest, rows=2, cols=2
MatNest structure:
(0,0) : prefix="fieldsplit_0_", type=seqaij, rows=12, cols=12
(0,1) : type=seqaij, rows=12, cols=4
(1,0) : type=seqaij, rows=4, cols=12
(1,1) : prefix="fieldsplit_1_", type=seqaij, rows=4, cols=4
$
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
More information about the petsc-users
mailing list