!$Author: wagner $
!$Date: 2006-08-11 18:44:39 $
!$Revision: 1.1 $
!$Source: /weru/cvs/weps/weps_ref/src/ref_pm10.f95,v $

!     +++ PURPOSE +++
! Calculate the horizontal PM10 discharge based upon
! the physical processes involved.  Each of the primary sources
! and sinks/traps are represented as individual contributing processes

      REAL FUNCTION dq10_dx (Gssen,Gssan, Gssbk, SF10en, SF10an, SF10bk) ! net PM10 discharge flux

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: Gssen  ! vert. susp. emission flux from loose aggs
      REAL, INTENT (IN) :: Gssan  ! vert. susp. flux from abrasion of surf. clods and crust
      REAL, INTENT (IN) :: Gssbk  ! vert. susp. flux from breakage of salt/creep aggs
      REAL, INTENT (IN) :: SF10en ! soil fraction of PM10 in suspension-size surface soil
      REAL, INTENT (IN) :: SF10an ! soil fraction of PM10 in susp-size aggs. from abraded clods & crust
      REAL, INTENT (IN) :: SF10bk ! soil fraction of PM10 in susp-size aggs. from salt/creep-size aggs


      ! dq10_dx - net PM10 discharge flux (kg/m2/s)
      ! Gssen  - vertical susp. emission flux from loose aggs (kg/m^2/s)
      ! Gssan  - vertical susp. emission flux from abrasion of surface clods and crust (kg/m^2/s)
      ! Gssbk  - vertical susp. emission flux from breakage of salt/creep aggs (kg/m^2/s)
      ! Gssdp  - vert. susp. flux (deposition) above a non-eroding surface (kg/m^2/s)
      ! SF10en - soil fraction of PM10 in suspension-size surface soil
      ! SF10an - soil fraction of PM10 in susp-size aggs. from abraded clods & crust
      ! SF10bk - soil fraction of PM10 in susp-size aggs. broken from salt/creep-size aggs


      dq10_dx = SF10en*Gssen + SF10an*Gssan + SF10bk*Gssbk

      END FUNCTION dq10_dx

!-----------------------------------------------------------------------

!     +++ PURPOSE +++
! Calculate the horizontal suspension discharge based upon
! the physical processes involved.  Each of the primary sources
! and sinks/traps are represented as individual contributing processes
! (Same as "dqss_dx" but with the source and sink terms re-arranged)

      REAL FUNCTION dq10_dx_1 (H, K, q) ! net PM10 discharge flux

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: H      ! intercept coefficient
      REAL, INTENT (IN) :: K      ! linear coefficient
      REAL, INTENT (IN) :: q      ! horizontal suspension discharge


      ! dq10_dx_1 - net PM10 discharge flux (kg/m2/s)
      ! H      - intercept coefficient (kg/m^2/s)
      ! K      - linear coefficient (unitless)
      ! q      - horizontal suspension discharge (kg/m/s)

      dq10_dx_1 = H + K*q

      END FUNCTION dq10_dx_1

!     +++ PURPOSE +++
! Calculate the "intercept" component of the vertical flux in the
! control volume for the PM10 component

      REAL FUNCTION H (SF10en, SFssen, Cen, qen) ! "intercept" component vert. flux eqn.

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: SF10en ! soil fraction of PM10 in suspension-size surface soil
      REAL, INTENT (IN) :: SFssen ! mass fract. of susp-size among loose aggs
      REAL, INTENT (IN) :: Cen    ! coefficient of emission
      REAL, INTENT (IN) :: qen    ! transport capacity

      ! H      - "intercept" component of vertical flux (kg/m^2/s)
      ! SF10en - soil fraction of PM10 in suspension-size surface soil
      ! SFssen - mass fraction of susp-size (<0.10 mm) among loose aggs (<2.0 mm dia.)
      ! Cen    - coefficient of emission (1/m)
      ! qen    - transport capacity (kg/m/s)

      H = (SF10en * SFssen * Cen) * qen

      END FUNCTION H


!     +++ PURPOSE +++
! Calculate the "linear" component of the vertical flux in the
! control volume for the PM10 component

      REAL FUNCTION K (SF10en, Cm, SFssen, Cen, SF10an, SFssan, Fancr,  &
     & Cancr, Fancl, Cancl, SF10bk, Cbk)! "linear" component of vert. flux eqn.

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: SF10en ! soil fraction of PM10 in suspension-size surface soil
      REAL, INTENT (IN) :: Cm     ! coefficient of mixing
      REAL, INTENT (IN) :: SFssen ! mass fract. of susp-size among loose aggs
      REAL, INTENT (IN) :: Cen    ! coefficient of emission
      REAL, INTENT (IN) :: SF10an ! soil fraction of PM10 in susp-size aggs. from abraded clods & crust
      REAL, INTENT (IN) :: SFssan ! mass fract. of susp-size from abrasion of clods and crust
      REAL, INTENT (IN) :: Fancr  ! mass fraction of saltation impacting crust
      REAL, INTENT (IN) :: Cancr  ! coefficient of abrasion on crust
      REAL, INTENT (IN) :: Fancl  ! mass fraction of saltation impacting clods
      REAL, INTENT (IN) :: Cancl  ! coefficient of abrasion on clods
      REAL, INTENT (IN) :: SF10bk ! soil fraction of PM10 in susp-size aggs. from salt/creep-size aggs
      REAL, INTENT (IN) :: Cbk    ! coefficient of breakage

      ! K      - "linear" component of vertical flux (unitless)
      ! SF10en - soil fraction of PM10 in suspension-size surface soil
      ! Cm     - coefficient of mixing (1/m)
      ! SFssen - mass fraction of susp-size (<0.10 mm) among loose aggs (<2.0 mm dia.)
      ! Cen    - coefficient of emission (1/m)
      ! SF10an - soil fraction of PM10 in susp-size aggs. from abraded clods & crust
      ! SFssan - mass fraction of susp-size (<0.10 mm) from abrasion of clods and crust
      ! Fancr  - mass fraction of saltation impacting crust
      ! Cancr  - coefficient of abrasion on crust (1/m)
      ! Fancl  - mass fraction of saltation impacting clods
      ! Cancl  - coefficient of abrasion on clods (1/m)
      ! SF10bk - soil fraction of PM10 in susp-size aggs. broken from salt/creep-size aggs
      ! Cbk    - coefficient of breakage (1/m)

      K = SF10en*(Cm - SFssen*Cen) +                                    &
     &    SF10an*SFssan*((Fancr*Cancr) + (Fancl*Cancl)) + SF10bk*Cbk

      END FUNCTION K


!-----------------------------------------------------------------------

!     +++ PURPOSE +++
! Calculate net soil loss from the
! control volume for the PM10 component

      REAL FUNCTION q10o (A, B, C, H, K, lx, q10i) ! net suspension discharge

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: A      ! intercept coefficient (salt/creep)
      REAL, INTENT (IN) :: B      ! linear coefficient (salt/creep)
      REAL, INTENT (IN) :: C      ! quadratic coefficient (salt/creep)
      REAL, INTENT (IN) :: H      ! intercept coefficient
      REAL, INTENT (IN) :: K      ! linear coefficient
      REAL, INTENT (IN) :: lx     ! length of control volume
      REAL, INTENT (IN) :: q10i   ! incoming PM10

      ! qo     - net suspension discharge (kg/m)
      ! A      - intercept coefficient (kg/m^2/s)
      ! B      - linear coefficient (unitless)
      ! C      - quadratic coefficient (s/m^2/kg)
      ! H      - intercept coefficient (kg/m^2/s)
      ! K      - linear coefficient (unitless)
      ! lx     - length of control volume (m)
      ! q10i   - incoming horizontal suspension (kg/m)

      REAL :: S                   ! temporary variables

      S = (4.0*A*C + B*B)**0.5
  
      ! Would like to know what "physical" conditions can cause p

      q10o = q10i + (1.0/(C*C)) * ( (K*S + K*B + H*H*C*C)*(lx) +        &
     &       K*K * (log(exp(S*(-lx))*(B+S) - B + S) ) )

      END FUNCTION q10o


!     +++ PURPOSE +++
! Calculate net soil loss from the
! control volume for the saltation/creep component

      REAL FUNCTION qsso_1 (A, B, C, F, G, lx, qi, qssi) ! net suspension discharge

      IMPLICIT NONE

!     +++ ARGUMENT DECLARATIONS +++
      REAL, INTENT (IN) :: A      ! intercept coefficient (salt/creep)
      REAL, INTENT (IN) :: B      ! linear coefficient (salt/creep)
      REAL, INTENT (IN) :: C      ! quadratic coefficient (salt/creep)
      REAL, INTENT (IN) :: F      ! intercept coefficient
      REAL, INTENT (IN) :: G      ! linear coefficient
      REAL, INTENT (IN) :: lx     ! length of control volume
      REAL, INTENT (IN) :: qi     ! incoming saltation/cree
      REAL, INTENT (IN) :: qssi   ! incoming horizontal suspension

      ! qo     - net suspension discharge (kg/m)
      ! A      - intercept coefficient (kg/m^2/s)
      ! B      - linear coefficient (unitless)
      ! C      - quadratic coefficient (s/m^2/kg)
      ! G      - intercept coefficient (kg/m^2/s)
      ! F      - linear coefficient (unitless)
      ! lx     - length of control volume (m)
      ! qi     - incoming horizontal saltation/creep (kg/m)
      ! qssi   - incoming horizontal suspension (kg/m)

      REAL :: S, p, t2                   ! temporary variables

      S = (4.0*A*C + B*B)**0.5
      p = (-(C+C)*qi + B) / S

      ! Would like to know what "physical" conditions can cause (S*lx)
      ! to go out of range - LEW
      IF ((S*lx) .gt. 40.0) THEN
          S = 40.0/lx
      ENDIF

      ! Would like to know what "physical" conditions can cause p
      ! to go out of range (-1 < p < 1) - LEW
      IF (p .gt. 1.0) THEN
         t2 = log(2.0)               ! Why the arbritrary values?
      ELSE IF (p .le. -1.0) THEN
         t2 = log(exp(S*lx)*2.0)     ! Why the arbritrary values?
      ELSE
         t2 = log(exp(S*lx)*(1.0-p) + p + 1.0)
      END IF 

      qsso_1 = qssi + (1.0/(C*C)) * ( (-G*S + G*B + F*F*C*C)*(lx) +      &
     &       G*G * (-log(2.0) + t2))

      !qsso_1 = (0.5 * qssi) * (1.0 + exp(-Cdp * lx))

      END FUNCTION qsso_1