tem_reduction_transient_update Subroutine

public subroutine tem_reduction_transient_update(me, res)

Apply time reduction on given res

Arguments

Type IntentOptional Attributes Name
type(tem_reduction_transient_type), intent(inout) :: me

Time reduction data to update

real(kind=rk), intent(in) :: res(:)

Current values of the variable to reduce.


Calls

proc~~tem_reduction_transient_update~~CallsGraph proc~tem_reduction_transient_update tem_reduction_transient_update proc~tem_reduction_transient_reset tem_reduction_transient_reset proc~tem_reduction_transient_update->proc~tem_reduction_transient_reset

Called by

proc~~tem_reduction_transient_update~~CalledByGraph proc~tem_reduction_transient_update tem_reduction_transient_update proc~tem_opvar_reduction_transient_update tem_opVar_reduction_transient_update proc~tem_opvar_reduction_transient_update->proc~tem_reduction_transient_update

Contents


Source Code

  subroutine tem_reduction_transient_update(me, res)
    ! -------------------------------------------------------------------- !
    !> Time reduction data to update
    type(tem_reduction_transient_type), intent(inout) :: me

    !> Current values of the variable to reduce.
    real(kind=rk), intent(in) :: res(:)
    ! -------------------------------------------------------------------- !
    integer :: idx
    ! -------------------------------------------------------------------- !
    me%nTimes = me%nTimes + 1

    select case( trim(me%config%reduceType) )
      case('min')
        do idx=1,me%nEntries
          me%val(idx, me%curr) = min( me%val(idx, me%curr),  res(idx) )
        end do

      case('max')
        do idx=1,me%nEntries
          me%val(idx, me%curr) = max( me%val(idx, me%curr),  res(idx) )
        end do

      case('sum', 'average')
        me%val(:, me%curr) = me%val(:, me%curr) + res(:)
    end select

    ! Check whether interval is completed
    if (me%nTimes == me%config%nRecord) then
      select case( trim(me%config%reduceType) )
        case('average')
          me%val(:, me%curr) = me%val(:, me%curr) / me%nTimes
      end select

      ! swap curr and last
      me%curr = mod(me%curr,2)+1
      me%last = mod(me%last,2)+1

      ! Reset current val reduction when nTimes has reached nRecord
      call tem_reduction_transient_reset(me)
    end if

  end subroutine tem_reduction_transient_update