MODULE compute_caldyn_mod
  USE prec, ONLY : rstd
  IMPLICIT NONE
  SAVE
  
  ! fake array dimensions, for interfaces
  INTEGER, PARAMETER :: iim_jjm_i=1, iim_jjm_u=1, iim_jjm_v =1, llm_=1, llm1=1, nqdyn_=1

  INTERFACE
     
    SUBROUTINE comp_pvort_only(u,rhodz,qu,qv, hv)
      IMPORT
      REAL(rstd),INTENT(IN)  :: u(iim_jjm_i, llm_)
      REAL(rstd),INTENT(INOUT) :: rhodz(iim_jjm_i, llm_)
      REAL(rstd),INTENT(OUT) :: qu(iim_jjm_u, llm_)
      REAL(rstd),INTENT(OUT) :: qv(iim_jjm_v, llm_)
      REAL(rstd),INTENT(OUT) :: hv(iim_jjm_v, llm_)
    END SUBROUTINE comp_pvort_only

    SUBROUTINE comp_theta(mass_col,theta_rhodz, rhodz,theta)
      IMPORT
      REAL(rstd),INTENT(IN)    :: mass_col(iim_jjm_i)
      REAL(rstd),INTENT(IN)    :: theta_rhodz(iim_jjm_i, llm_, nqdyn_)
      REAL(rstd),INTENT(INOUT) :: rhodz(iim_jjm_i, llm_)
      REAL(rstd),INTENT(OUT)   :: theta(iim_jjm_i, llm_, nqdyn_)
    END SUBROUTINE comp_theta

    SUBROUTINE comp_geopot(rhodz,theta, ps,pk,geopot) 
      IMPORT
      REAL(rstd),INTENT(IN)    :: rhodz(iim_jjm_i, llm_)         ! rho*dz = mass per unit surface in each full model level
      REAL(rstd),INTENT(IN)    :: theta(iim_jjm_i, llm_, nqdyn_) ! active scalars : theta/entropy, moisture, ...
      REAL(rstd),INTENT(INOUT) :: ps(iim_jjm_i)                  ! surface pressure
      REAL(rstd),INTENT(OUT)   :: pk(iim_jjm_i, llm_)            ! Exner function (compressible) /Lagrange multiplier (Boussinesq)
      REAL(rstd),INTENT(INOUT) :: geopot(iim_jjm_i, llm1)        ! geopotential
    END SUBROUTINE comp_geopot

    SUBROUTINE comp_caldyn_fast(tau,theta,geopot, pk,berni,du,u)
      IMPORT
      REAL(rstd),INTENT(IN)    :: tau                ! "solve" u-tau*du/dt = rhs
      REAL(rstd),INTENT(IN)    :: theta(iim_jjm_i, llm_, nqdyn_)
      REAL(rstd),INTENT(IN)    :: geopot(iim_jjm_i, llm1)
      REAL(rstd),INTENT(INOUT) :: pk(iim_jjm_i, llm_)
      REAL(rstd),INTENT(INOUT) :: berni(iim_jjm_i, llm_)  ! partial Bernoulli function
      REAL(rstd),INTENT(INOUT) :: du(iim_jjm_u, llm_)
      REAL(rstd),INTENT(INOUT) :: u(iim_jjm_u, llm_)   ! INOUT if tau>0
    END SUBROUTINE comp_caldyn_fast

    SUBROUTINE comp_caldyn_slow_hydro(zero, u,rhodz,hv,Kv, berni, hflux,du)
      IMPORT
      LOGICAL, INTENT(IN) :: zero
      REAL(rstd),INTENT(IN)  :: u(iim_jjm_u, llm_)    ! prognostic "velocity"
      REAL(rstd),INTENT(IN)  :: rhodz(iim_jjm_i, llm_)
      REAL(rstd),INTENT(IN)  :: hv(iim_jjm_v, llm_)   ! height/mass averaged to vertices
      REAL(rstd),INTENT(IN)  :: Kv(iim_jjm_v, llm_)   ! kinetic energy at vertices
      REAL(rstd), INTENT(OUT) :: berni(iim_jjm_i, llm_)  ! Bernoulli function
      REAL(rstd),INTENT(OUT) :: hflux(iim_jjm_u, llm_) ! hflux in kg/s 
      REAL(rstd),INTENT(INOUT) :: du(iim_jjm_u, llm_)
    END SUBROUTINE comp_caldyn_slow_hydro
    
    SUBROUTINE comp_caldyn_coriolis(hflux,theta,qu, Ftheta, convm,dtheta_rhodz,du)
      IMPORT
      REAL(rstd),INTENT(IN)    :: hflux(iim_jjm_u, llm_)  ! hflux in kg/s
      REAL(rstd),INTENT(IN)    :: theta(iim_jjm_i, llm_, nqdyn_) ! active scalars
      REAL(rstd),INTENT(IN)    :: qu(iim_jjm_u, llm_)
      REAL(rstd), INTENT(OUT)  :: Ftheta(iim_jjm_u, llm_)  ! potential temperature flux
      REAL(rstd),INTENT(OUT)   :: convm(iim_jjm_i, llm_)  ! mass flux convergence
      REAL(rstd),INTENT(OUT)   :: dtheta_rhodz(iim_jjm_i, llm_, nqdyn_)
      REAL(rstd),INTENT(INOUT) :: du(iim_jjm_u, llm_)
    END SUBROUTINE comp_caldyn_coriolis

  END INTERFACE

  PROCEDURE(comp_pvort_only),        POINTER :: compute_pvort_only        => NULL()
  PROCEDURE(comp_theta),             POINTER :: compute_theta             => NULL()
  PROCEDURE(comp_geopot),            POINTER :: compute_geopot            => NULL()
  PROCEDURE(comp_caldyn_fast),       POINTER :: compute_caldyn_fast       => NULL()
  PROCEDURE(comp_caldyn_slow_hydro), POINTER :: compute_caldyn_slow_hydro => NULL()
  PROCEDURE(comp_caldyn_coriolis),   POINTER :: compute_caldyn_coriolis   => NULL()

END MODULE compute_caldyn_mod