calc_kr_mvg Function

pure elemental function calc_kr_mvg(theta_s, theta_r, alpha1, n1, m1, l, h_crit, h) result(kr)

Arguments

Type IntentOptional Attributes Name
real(kind=real64), intent(in) :: theta_s
real(kind=real64), intent(in) :: theta_r
real(kind=real64), intent(in) :: alpha1
real(kind=real64), intent(in) :: n1
real(kind=real64), intent(in) :: m1
real(kind=real64), intent(in) :: l
real(kind=real64), intent(in) :: h_crit
real(kind=real64), intent(in) :: h

Return Value real(kind=real64)


Called by

proc~~calc_kr_mvg~~CalledByGraph proc~calc_kr_mvg calc_kr_mvg proc~calc_kr_base_mvg calc_kr_base_mvg proc~calc_kr_base_mvg->proc~calc_kr_mvg interface~calc_kr_base_mvg type_hcf_base_mvg%calc_kr_base_mvg interface~calc_kr_base_mvg->proc~calc_kr_base_mvg

Source Code

    pure elemental function calc_kr_mvg(theta_s, theta_r, alpha1, n1, m1, l, h_crit, h) result(kr)
        implicit none
        real(real64), intent(in) :: theta_s
        real(real64), intent(in) :: theta_r
        real(real64), intent(in) :: alpha1
        real(real64), intent(in) :: n1
        real(real64), intent(in) :: m1
        real(real64), intent(in) :: l
        real(real64), intent(in) :: h_crit
        real(real64), intent(in) :: h
        real(real64) :: kr
        real(real64) :: s_w, theta_m

        theta_m = theta_r + (theta_s - theta_r) * (1.0d0 + (-alpha1 * h_crit)**n1)**(-m1)

        if (h < h_crit) then
            s_w = (theta_s - theta_r) / (theta_m - theta_r) * (1.0d0 + abs(alpha1 * h)**n1)**(-m1)
            kr = s_w**l * ((1.0d0 - (1.0d0 - s_w**(1.0d0 / m1))**m1) / (1.0d0 - (1.0d0 - 1.0d0**(1.0d0 / m1))**m1))**2.0d0
        else
            kr = 1.0d0
        end if

    end function calc_kr_mvg