Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

fix(par): fix backtracking and residual norm calculation for certain cases in parallel #1780

Merged
merged 2 commits into from
May 8, 2024
Merged

fix(par): fix backtracking and residual norm calculation for certain cases in parallel #1780

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
8beb2ef
- update petsc matrix before multiplication
mjr-deltares May 8, 2024
d41bdb1
- fprettify
mjr-deltares May 8, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
90 changes: 58 additions & 32 deletions src/Solution/NumericalSolution.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -168,6 +168,8 @@ module NumericalSolutionModule
procedure :: sln_calc_ptc
procedure :: sln_underrelax
procedure :: sln_backtracking_xupdate
procedure :: get_backtracking_flag
procedure :: apply_backtracking

! private
procedure, private :: sln_connect
Expand All @@ -179,7 +181,7 @@ module NumericalSolutionModule
procedure, private :: sln_calcdx
procedure, private :: sln_calc_residual
procedure, private :: sln_l2norm
procedure, private :: sln_outer_check
procedure, private :: sln_get_dxmax
procedure, private :: sln_get_loc
procedure, private :: sln_get_nodeu
procedure, private :: allocate_scalars
Expand Down Expand Up @@ -1617,7 +1619,7 @@ subroutine solve(this, kiter)
!-------------------------------------------------------
!
! -- check convergence of solution
call this%sln_outer_check(this%hncg(kiter), this%lrch(1, kiter))
call this%sln_get_dxmax(this%hncg(kiter), this%lrch(1, kiter))
this%icnvg = 0
if (this%sln_has_converged(this%hncg(kiter))) then
this%icnvg = 1
Expand Down Expand Up @@ -1752,7 +1754,7 @@ subroutine solve(this, kiter)
this%icnvg = 1
!
! -- reset outer dependent-variable change and location for output
call this%sln_outer_check(this%hncg(kiter), this%lrch(1, kiter))
call this%sln_get_dxmax(this%hncg(kiter), this%lrch(1, kiter))
!
! -- write revised dependent-variable change data after
! newton under-relaxation
Expand Down Expand Up @@ -2790,41 +2792,65 @@ end subroutine sln_backtracking
!! update exceeds the dependent variable closure criteria.
!!
!<
subroutine sln_backtracking_xupdate(this, btflag)
subroutine sln_backtracking_xupdate(this, bt_flag)
! -- dummy variables
class(NumericalSolutionType), intent(inout) :: this !< NumericalSolutionType instance
integer(I4B), intent(inout) :: btflag !< backtracking flag (1) backtracking performed (0) backtracking not performed
! -- local variables
integer(I4B), intent(inout) :: bt_flag !< backtracking flag (1) backtracking performed (0) backtracking not performed

bt_flag = this%get_backtracking_flag()

! perform backtracking if ...
if (bt_flag > 0) then
call this%apply_backtracking()
end if

end subroutine sln_backtracking_xupdate

!> @brief Check if backtracking should be applied for this solution,
!< returns 1: yes, 0: no
function get_backtracking_flag(this) result(bt_flag)
class(NumericalSolutionType) :: this !< NumericalSolutionType instance
integer(I4B) :: bt_flag !< backtracking flag (1) backtracking performed (0) backtracking not performed
! local
integer(I4B) :: n
real(DP) :: dx
real(DP) :: dx_abs
real(DP) :: dx_abs_max

! default is off
bt_flag = 0

! find max. change
dx_abs_max = 0.0
do n = 1, this%neq
if (this%active(n) < 1) cycle
dx = this%x(n) - this%xtemp(n)
dx_abs = abs(dx)
if (dx_abs > dx_abs_max) dx_abs_max = dx_abs
end do

! if backtracking, set flag
if (this%breduc * dx_abs_max >= this%dvclose) then
bt_flag = 1
end if

end function get_backtracking_flag

!> @brief Update x with backtracking
!<
subroutine apply_backtracking(this)
class(NumericalSolutionType) :: this !< NumericalSolutionType instance
! local
integer(I4B) :: n
real(DP) :: delx
real(DP) :: absdelx
real(DP) :: chmax
!
! -- initialize dummy variables
btflag = 0
!
! -- no backtracking if maximum change is less than closure so return
chmax = 0.0

do n = 1, this%neq
if (this%active(n) < 1) cycle
delx = this%breduc * (this%x(n) - this%xtemp(n))
absdelx = abs(delx)
if (absdelx > chmax) chmax = absdelx
this%x(n) = this%xtemp(n) + delx
end do
!
! -- perform backtracking if free of constraints and set counter and flag
if (chmax >= this%dvclose) then
btflag = 1
do n = 1, this%neq
if (this%active(n) < 1) cycle
delx = this%breduc * (this%x(n) - this%xtemp(n))
this%x(n) = this%xtemp(n) + delx
end do
end if
!
! -- return
return
end subroutine sln_backtracking_xupdate

end subroutine

!> @ brief Calculate the solution L-2 norm for all
!! active cells using
Expand Down Expand Up @@ -3116,7 +3142,7 @@ end subroutine sln_underrelax
!! Picard iteration.
!!
!<
subroutine sln_outer_check(this, hncg, lrch)
subroutine sln_get_dxmax(this, hncg, lrch)
! -- dummy variables
class(NumericalSolutionType), intent(inout) :: this !< NumericalSolutionType instance
real(DP), intent(inout) :: hncg !< maximum dependent-variable change
Expand Down Expand Up @@ -3150,7 +3176,7 @@ subroutine sln_outer_check(this, hncg, lrch)
!
! -- return
return
end subroutine sln_outer_check
end subroutine sln_get_dxmax

function sln_has_converged(this, max_dvc) result(has_converged)
class(NumericalSolutionType) :: this !< NumericalSolutionType instance
Expand Down
16 changes: 11 additions & 5 deletions src/Solution/ParallelSolution.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -183,24 +183,30 @@ end subroutine par_underrelax

!> @brief synchronize backtracking flag over processes
!<
subroutine par_backtracking_xupdate(this, btflag)
subroutine par_backtracking_xupdate(this, bt_flag)
! -- dummy variables
class(ParallelSolutionType), intent(inout) :: this !< ParallelSolutionType instance
integer(I4B), intent(inout) :: btflag !< global backtracking flag (1) backtracking performed (0) backtracking not performed
integer(I4B), intent(inout) :: bt_flag !< global backtracking flag (1) backtracking performed (0) backtracking not performed
! -- local variables
integer(I4B) :: btflag_local
type(MpiWorldType), pointer :: mpi_world
integer :: ierr

mpi_world => get_mpi_world()

btflag_local = 0
call this%NumericalSolutionType%sln_backtracking_xupdate(btflag_local)
! get local bt flag
btflag_local = this%NumericalSolutionType%get_backtracking_flag()

call MPI_Allreduce(btflag_local, btflag, 1, MPI_INTEGER, &
! reduce into global decision (if any, then all)
call MPI_Allreduce(btflag_local, bt_flag, 1, MPI_INTEGER, &
MPI_MAX, mpi_world%comm, ierr)
call CHECK_MPI(ierr)

! perform backtracking if ...
if (bt_flag > 0) then
call this%NumericalSolutionType%apply_backtracking()
end if

end subroutine par_backtracking_xupdate

end module ParallelSolutionModule
7 changes: 5 additions & 2 deletions src/Utilities/Matrix/PetscMatrix.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,8 @@ module PetscMatrixModule
private

type, public, extends(MatrixBaseType) :: PetscMatrixType
Mat :: mat
! offset in the global matrix
Mat :: mat !< the PETSc matrix object, NOTE: update() should be called before using this,
!! in case the matrix CSR array has changed!!!
integer(I4B) :: nrow !< number of rows in this portion of the global matrix
integer(I4B) :: ncol !< number of columns in the matrix
integer(I4B) :: nnz !< number of nonzeros in the matrix
Expand Down Expand Up @@ -446,6 +446,9 @@ subroutine pm_multiply(this, vec_x, vec_y)
y => vec_y
end select

! copy data into petsc object
call this%update()
! and multiply
call MatMult(this%mat, x%vec_impl, y%vec_impl, ierr)
CHKERRQ(ierr)

Expand Down
Loading