convergence_criterion.F90

Go to the documentation of this file.

!! Copyright (C)  2020 N. Tancogne-Dejean
!!
!! This Source Code Form is subject to the terms of the Mozilla Public
!! License, v. 2.0. If a copy of the MPL was not distributed with this
!! file, You can obtain one at https://mozilla.org/MPL/2.0/.
!!
 
#include "global.h"
 
module convergence_criterion_oct_m
  use debug_oct_m
  use global_oct_m
  use linked_list_oct_m
  use profiling_oct_m
  use unit_oct_m
 
  implicit none
 
  private
  public ::                    &
    convergence_criterion_t,    &
    convergence_criterion_end,  &
    criterion_list_t,           &
    criterion_iterator_t
 
 
  type, abstract :: convergence_criterion_t
    private
    float,   public :: tol_abs 
    float,   public :: tol_rel 
    character(len=:), allocatable, public :: label
    type(unit_t), pointer, public :: unit => null()
 
    float, public   :: val_abs    
    float, public   :: val_rel
    float, pointer  :: value_diff
    float, pointer  :: norm
 
  contains
    procedure :: is_converged => convergence_criterion_is_converged
    procedure :: set_pointers => criterion_set_quantity_pointers
    procedure :: write_info => criterion_write_info
  end type convergence_criterion_t
 
  type, extends(linked_list_t) :: criterion_list_t
    private
  contains
    procedure :: add => criterion_list_add_node
  end type criterion_list_t
 
  type, extends(linked_list_iterator_t) :: criterion_iterator_t
    private
  contains
    procedure :: get_next => criterion_iterator_get_next
  end type criterion_iterator_t
 
 
contains
 
  ! ---------------------------------------------------------
  subroutine criterion_write_info(this, iunit)
    class(convergence_criterion_t), intent(inout) :: this
    integer,                        intent(in)    :: iunit
 
    push_sub(criterion_write_info)
 
    if (associated(this%unit)) then
      write(iunit, '(6x,a,a,a,es15.8,a,es15.8,4a)') 'abs_', this%label, ' = ', &
        units_from_atomic(this%unit, this%val_abs), &
        ' (', units_from_atomic(this%unit, this%tol_abs), ')', &
        ' [',  trim(units_abbrev(this%unit)), ']'
    else
      write(iunit, '(6x,a,a,a,es15.8,a,es15.8,a)') 'abs_', this%label, ' = ', this%val_abs, ' (', this%tol_abs, ')'
    end if
    write(iunit, '(6x,a,a,a,es15.8,a,es15.8,a)') 'rel_', this%label, ' = ', this%val_rel, ' (', this%tol_rel, ')'
 
    pop_sub(criterion_write_info)
  end subroutine criterion_write_info
 
  ! ---------------------------------------------------------
  subroutine convergence_criterion_is_converged(this, is_converged)
    class(convergence_criterion_t),  intent(inout) :: this
    logical,                        intent(out)   :: is_converged
 
    push_sub(convergence_criterion_is_converged)
 
    this%val_abs = abs(this%value_diff)
 
    assert(associated(this%norm))
    if (abs(this%norm) < m_tiny) then
      this%val_rel = m_huge
    else
      this%val_rel = this%val_abs / abs(this%norm)
    end if
 
    if (this%tol_abs <= m_zero) then
      is_converged = .true.
    else
      is_converged = this%val_abs < this%tol_abs
    end if
 
    if (this%tol_rel > m_zero) then
      is_converged = is_converged .and. this%val_rel < this%tol_rel
    end if
 
 
    pop_sub(convergence_criterion_is_converged)
  end subroutine convergence_criterion_is_converged
 
  ! ---------------------------------------------------------
  subroutine criterion_set_quantity_pointers(this, value_diff, value_norm)
    class(convergence_criterion_t),  intent(inout) :: this
    float, target,                  intent(in)    :: value_diff
    float, target,                  intent(in)    :: value_norm
 
    push_sub(criterion_set_quantity_pointers)
 
    this%value_diff => value_diff
    this%norm => value_norm
 
    pop_sub(criterion_set_quantity_pointers)
  end subroutine criterion_set_quantity_pointers
 
  ! ---------------------------------------------------------
  subroutine convergence_criterion_end(this)
    class(convergence_criterion_t),  intent(inout) :: this
 
    push_sub(convergence_criterion_end)
 
    nullify(this%value_diff)
    nullify(this%norm)
    nullify(this%unit)
    if (allocated(this%label)) then
      ! No safe_deallocate here, because it was not allocated with safe_allocate.
      deallocate(this%label)
    end if
 
    pop_sub(convergence_criterion_end)
  end subroutine convergence_criterion_end
 
  ! ---------------------------------------------------------
  subroutine criterion_list_add_node(this, criterion)
    class(criterion_list_t)                 :: this
    class(convergence_criterion_t), target :: criterion
 
    push_sub(criterion_list_add_node)
 
    call this%add_ptr(criterion)
 
    pop_sub(criterion_list_add_node)
  end subroutine criterion_list_add_node
 
  ! ---------------------------------------------------------
  function criterion_iterator_get_next(this) result(criterion)
    class(criterion_iterator_t),      intent(inout) :: this
    class(convergence_criterion_t),  pointer       :: criterion
 
    push_sub(criterion_iterator_get_next)
 
    select type (ptr => this%get_next_ptr())
    class is (convergence_criterion_t)
      criterion => ptr
    class default
      assert(.false.)
    end select
 
    pop_sub(criterion_iterator_get_next)
  end function criterion_iterator_get_next
 
end module convergence_criterion_oct_m
 
!! Local Variables:
!! mode: f90
!! coding: utf-8
!! End: