Hi All,
I have many linear equations with the same matrix structure (same
non-zero entries) that are derived from a flow problem at different time
steps. I feel puzzled that the results are a little different when the
solver run repeatedly and one by one. Say, I have three equations, I can
get the following results if running three equations together
Equation 1: Iterations 1 norm 0.9457E-02 Result error PETSc vs
Solver2, max 0.4362E-02 min -0.2277E-04 norm 0.9458E-02
Equation 2: Iterations 2 norm 0.2994E-05 Result error PETSc vs
Solver2, max 0.1381E-05 min -0.7209E-08 norm 0.2994E-05
Equation 3: Iterations 2 norm 0.3919E-04 Result error PETSc vs
Solver2, max 0.9435E-07 min -0.1808E-04 norm 0.3919E-04
But if I solve only one equation every time, then restart the program to
run another one, the results are like this:
Equation 1: Iterations 1 norm 0.9457E-02 Result error PETSc vs
Solver2, max 0.4362E-02 min -0.2277E-04 norm 0.9458E-02
Equation 2: Iterations 1 norm 0.7949E-05 Result error PETSc vs
Solver2, max 0.3501E-05 min -0.8377E-06 norm 0.7949E-05
Equation 3: Iterations 1 norm 0.1980E-04 Result error PETSc vs
Solver2, max 0.4168E-08 min -0.9085E-05 norm 0.1980E-04
Note: Solver2 is the original sequential solver used in this flow model.
Though there are no big difference in the solution for the above
equations, I want to know why?
For another large linear equations with more than 400,000 unknowns and
10,000,000 non-zero entries, if the equations are solved repeatedly,
they need a lot of iterations or fail, but if the equations are solved
one by one, it only needs 1 to 2 iterations.
How does this difference come from?
The sample codes are attached bellow.
Thanks and regards,
Danyang
!***************************************************************************!
!Create matrix, rhs and solver
call MatCreateAIJ(Petsc_Comm_World, Petsc_Decide, Petsc_Decide, nb, nb,
nd_nzrow, &
Petsc_Null_Integer, nd_nzrow, Petsc_Null_Integer, a,
ierr)
call MatSetOption(a,Mat_New_Nonzero_Allocation_Err,Petsc_False,ierr)
call VecCreateMPI(Petsc_Comm_World, Petsc_Decide, nb, b, ierr)
call VecDuplicate(b, x, ierr)
call VecDuplicate(x, u, ierr)
call KSPCreate(Petsc_Comm_World,ksp,ierr)
call KSPSetTolerances(ksp,tol, &
PETSC_DEFAULT_DOUBLE_PRECISION, &
PETSC_DEFAULT_DOUBLE_PRECISION, &
100,ierr)
call KSPSetFromOptions(ksp,ierr)
!Do time loop
do i = 1, nTimeStep
call MatGetOwnershipRange(a,istart,iend,ierr)
do i = istart, iend - 1
ii = ia_in(i+1)
jj = ia_in(i+2)
call MatSetValues(a, ione, i, jj-ii, ja_in(ii:jj-1)-1,
a_in(ii:jj-1), Insert_Values, ierr)
end do
call MatAssemblyBegin(a, Mat_Final_Assembly, ierr)
call MatAssemblyEnd(a, Mat_Final_Assembly, ierr)
call VecGetOwnershipRange(b,istart,iend,ierr)
call VecSetValues(b, iend-istart, ix(istart+1:iend),
b_in(istart+1:iend), Insert_Values, ierr)
call VecAssemblyBegin(b,ierr)
call VecAssemblyEnd(b,ierr)
if(i == 1) then
call MatConvert(a,MATSAME,MAT_INITIAL_MATRIX,a2,ierr)
end if
!call KSPSetOperators(ksp,a,a2,SAME_PRECONDITIONER,ierr)
call KSPSetOperators(ksp,a,a2,SAME_NONZERO_PATTERN,ierr) !These
three patterns make no difference in current codes
!call KSPSetOperators(ksp,a,a2,DIFFERENT_NONZERO_PATTERN,ierr)
call KSPSolve(ksp,b,x,ierr)
call KSPGetResidualNorm(ksp,norm,ierr)
call KSPGetIterationNumber(ksp,its,ierr)
end do
!Destroy objects
!...
!***************************************************************************!
