Index: codes/icosagcm/devel/src/dynamics/caldyn_hevi.f90
===================================================================
--- codes/icosagcm/devel/src/dynamics/caldyn_hevi.f90	(revision 920)
+++ codes/icosagcm/devel/src/dynamics/caldyn_hevi.f90	(revision 921)
@@ -9,8 +9,8 @@
   USE compute_caldyn_kv_mod, ONLY : compute_caldyn_kv
   USE compute_caldyn_Coriolis_mod, ONLY : compute_caldyn_Coriolis
-  USE compute_caldyn_slow_hydro_mod, ONLY : compute_caldyn_slow_hydro
+  USE compute_caldyn_slow_hydro_mod, ONLY : compute_caldyn_slow_hydro => compute_caldyn_slow_hydro_manual
   USE compute_caldyn_slow_NH_mod, ONLY : compute_caldyn_slow_NH
   USE compute_caldyn_solver_mod, ONLY : compute_caldyn_solver
-  USE compute_caldyn_fast_mod, ONLY : compute_caldyn_fast
+  USE compute_caldyn_fast_mod, ONLY : compute_caldyn_fast => compute_caldyn_fast_manual
   USE compute_NH_geopot_mod, ONLY : compute_NH_geopot
   IMPLICIT NONE
@@ -136,5 +136,6 @@
        END IF
        u=f_u(ind)
-       CALL compute_caldyn_fast(tau,u,mass,theta,pk,geopot,du) ! computes du_fast and updates u
+       Kv=f_Kv(ind) ! buffer for partial Bernoulli function
+       CALL compute_caldyn_fast(tau,theta,geopot, pk,Kv,du,u) ! computes du_fast and updates u
     ENDDO
     
@@ -185,5 +186,5 @@
           hv=f_hv(ind)
           Kv=f_Kv(ind)
-          CALL compute_caldyn_slow_hydro(u,mass,hv, hflux,Kv,du, .TRUE.)
+          CALL compute_caldyn_slow_hydro(.TRUE., u,mass,hv,Kv, Kv, hflux,du)
        ELSE
           W = f_W(ind)
Index: codes/icosagcm/devel/src/dynamics/compute_caldyn_fast.F90
===================================================================
--- codes/icosagcm/devel/src/dynamics/compute_caldyn_fast.F90	(revision 920)
+++ codes/icosagcm/devel/src/dynamics/compute_caldyn_fast.F90	(revision 921)
@@ -1,4 +1,9 @@
 MODULE compute_caldyn_fast_mod
-  USE grid_param, ONLY : llm
+  USE prec, ONLY : rstd
+  USE grid_param
+  USE earth_const
+  USE disvert_mod
+  USE omp_para
+  USE trace
   IMPLICIT NONE
   PRIVATE
@@ -6,5 +11,5 @@
 #include "../unstructured/unstructured.h90"
 
-  PUBLIC :: compute_caldyn_fast
+  PUBLIC :: compute_caldyn_fast_unst, compute_caldyn_fast_hex, compute_caldyn_fast_manual
 
 CONTAINS
@@ -65,8 +70,6 @@
 #endif END_DYSL
 
-  SUBROUTINE compute_caldyn_fast_unst(tau, pk,berni,theta,geopot, du,u)
+  SUBROUTINE compute_caldyn_fast_unst(tau,theta,geopot, pk,berni,du,u)
     USE ISO_C_BINDING, only : C_DOUBLE, C_FLOAT
-    USE earth_const
-    USE grid_param, ONLY : nqdyn
     USE data_unstructured_mod, ONLY : enter_trace, exit_trace, &
          id_fast, primal_num, dual_num, edge_num, &
@@ -86,27 +89,38 @@
   END SUBROUTINE compute_caldyn_fast_unst
 
-  SUBROUTINE compute_caldyn_fast(tau,u,rhodz,theta,pk,geopot,du)
+  SUBROUTINE compute_caldyn_fast_hex(tau,theta,geopot, pk,berni,du,u)
     USE icosa
-    USE trace
-    USE caldyn_vars_mod
-    USE omp_para, ONLY : ll_begin, ll_end
     REAL(rstd),INTENT(IN)    :: tau                ! "solve" u-tau*du/dt = rhs
-    REAL(rstd),INTENT(INOUT) :: u(iim*3*jjm,llm)   ! OUT if tau>0
-    REAL(rstd),INTENT(IN)    :: rhodz(iim*jjm,llm)
     REAL(rstd),INTENT(IN)    :: theta(iim*jjm,llm,nqdyn)
+    REAL(rstd),INTENT(IN)    :: geopot(iim*jjm,llm+1)
     REAL(rstd),INTENT(INOUT) :: pk(iim*jjm,llm)
-    REAL(rstd),INTENT(INOUT) :: geopot(iim*jjm,llm+1)
-    REAL(rstd),INTENT(INOUT)   :: du(iim*3*jjm,llm)
-    REAL(rstd) :: berni(iim*jjm,llm)  ! Bernoulli function
+    REAL(rstd),INTENT(INOUT) :: berni(iim*jjm,llm)  ! partial Bernoulli function
+    REAL(rstd),INTENT(INOUT) :: du(iim*3*jjm,llm)
+    REAL(rstd),INTENT(INOUT) :: u(iim*3*jjm,llm)   ! INOUT if tau>0
     REAL(rstd) :: berniv(iim*jjm,llm)  ! moist Bernoulli function
 
-    INTEGER :: i,j,ij,l
+    INTEGER :: ij,l
     REAL(rstd) :: cp_ik, qv, temp, chi, log_p_preff, due, due_right, due_lup, due_ldown
 
     CALL trace_start("compute_caldyn_fast")
-
-    IF(dysl_caldyn_fast) THEN
 #include "../kernels_hex/caldyn_fast.k90"
-    ELSE
+    CALL trace_end("compute_caldyn_fast")
+  END SUBROUTINE compute_caldyn_fast_hex
+
+  SUBROUTINE compute_caldyn_fast_manual(tau,theta,geopot, pk,berni,du,u)
+    USE icosa
+    REAL(rstd),INTENT(IN)    :: tau                ! "solve" u-tau*du/dt = rhs
+    REAL(rstd),INTENT(IN)    :: theta(iim*jjm,llm,nqdyn)
+    REAL(rstd),INTENT(IN)    :: geopot(iim*jjm,llm+1)
+    REAL(rstd),INTENT(INOUT) :: pk(iim*jjm,llm)
+    REAL(rstd),INTENT(INOUT) :: berni(iim*jjm,llm)  ! partial Bernoulli function
+    REAL(rstd),INTENT(INOUT) :: du(iim*3*jjm,llm)
+    REAL(rstd),INTENT(INOUT) :: u(iim*3*jjm,llm)   ! INOUT if tau>0
+    REAL(rstd) :: berniv(iim*jjm,llm)  ! moist Bernoulli function
+
+    INTEGER :: ij,l
+    REAL(rstd) :: cp_ik, qv, temp, chi, log_p_preff, due, due_right, due_lup, due_ldown
+
+    CALL trace_start("compute_caldyn_fast")
 
     ! Compute Bernoulli term
@@ -208,8 +222,7 @@
     END DO
 
-    END IF ! dysl
     CALL trace_end("compute_caldyn_fast")
 
-  END SUBROUTINE compute_caldyn_fast
+  END SUBROUTINE compute_caldyn_fast_manual
 
 END MODULE compute_caldyn_fast_mod
Index: codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_hydro.F90
===================================================================
--- codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_hydro.F90	(revision 920)
+++ codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_hydro.F90	(revision 921)
@@ -1,4 +1,9 @@
 MODULE compute_caldyn_slow_hydro_mod
-  USE grid_param, ONLY : llm
+  USE prec, ONLY : rstd
+  USE grid_param
+  USE earth_const
+  USE disvert_mod
+  USE omp_para
+  USE trace
   IMPLICIT NONE
   PRIVATE
@@ -6,5 +11,6 @@
 #include "../unstructured/unstructured.h90"
 
-  PUBLIC :: compute_caldyn_slow_hydro
+  PUBLIC :: compute_caldyn_slow_hydro_unst, &
+       compute_caldyn_slow_hydro_hex, compute_caldyn_slow_hydro_manual
 
 CONTAINS
@@ -47,15 +53,13 @@
 #endif END_DYSL
 
-  SUBROUTINE compute_caldyn_slow_hydro_unst(rhodz,theta,u, berni,hflux,du)
+  SUBROUTINE compute_caldyn_slow_hydro_unst(zero, u,rhodz,hv,Kv, berni, hflux,du)
     USE ISO_C_BINDING, only : C_DOUBLE, C_FLOAT
-    USE grid_param, ONLY : nqdyn
     USE data_unstructured_mod, ONLY :left,right,edge_num,primal_num,dual_num,id_slow_hydro, &
         primal_deg,le_de,primal_edge,Ai,enter_trace, exit_trace
-    FIELD_MASS  :: rhodz,berni! IN
-    FIELD_THETA :: theta      ! IN
+    LOGICAL, INTENT(IN) :: zero
+    FIELD_MASS  :: rhodz, hv, Kv, berni  ! IN, IN, IN, BUF
     FIELD_U     :: u,hflux,du ! IN, OUT, OUT
     DECLARE_INDICES
     DECLARE_EDGES
-    LOGICAL, PARAMETER :: zero=.TRUE.
     NUM :: ke, uu
     START_TRACE(id_slow_hydro, 3,0,3)
@@ -64,8 +68,6 @@
   END SUBROUTINE compute_caldyn_slow_hydro_unst
 
-  SUBROUTINE compute_caldyn_slow_hydro(u,rhodz,hv, hflux,Kv,du, zero)
+  SUBROUTINE compute_caldyn_slow_hydro_hex(zero, u,rhodz,hv,Kv, berni, hflux,du)
     USE icosa
-    USE trace
-    USE omp_para, ONLY : ll_begin, ll_end
     USE caldyn_vars_mod
     LOGICAL, INTENT(IN) :: zero
@@ -83,12 +85,28 @@
 
     CALL trace_start("compute_caldyn_slow_hydro")
-
-    IF(dysl_slow_hydro) THEN
-
 #define BERNI(ij,l) berni(ij,l)
 #include "../kernels_hex/caldyn_slow_hydro.k90"
 #undef BERNI
+    CALL trace_end("compute_caldyn_slow_hydro")    
 
-     ELSE
+  END SUBROUTINE compute_caldyn_slow_hydro_hex
+
+  SUBROUTINE compute_caldyn_slow_hydro_manual(zero, u,rhodz,hv,Kv, berni, hflux,du)
+    USE icosa
+    USE caldyn_vars_mod
+    LOGICAL, INTENT(IN) :: zero
+    REAL(rstd),INTENT(IN)  :: u(3*iim*jjm,llm)    ! prognostic "velocity"
+    REAL(rstd),INTENT(IN)  :: Kv(2*iim*jjm,llm)   ! kinetic energy at vertices
+    REAL(rstd),INTENT(IN)  :: hv(2*iim*jjm,llm)   ! height/mass averaged to vertices
+    REAL(rstd),INTENT(IN)  :: rhodz(iim*jjm,llm)
+    REAL(rstd),INTENT(OUT) :: hflux(3*iim*jjm,llm) ! hflux in kg/s
+    REAL(rstd),INTENT(INOUT) :: du(3*iim*jjm,llm)
+    
+    REAL(rstd) :: berni(iim*jjm,llm)  ! Bernoulli function
+    REAL(rstd) :: berni1(iim*jjm)  ! Bernoulli function
+    REAL(rstd) :: uu_right, uu_lup, uu_ldown, ke, uu
+    INTEGER :: ij,l
+
+    CALL trace_start("compute_caldyn_slow_hydro")
 
 #define BERNI(ij) berni1(ij)
@@ -172,7 +190,6 @@
 #undef BERNI
 
-    END IF ! dysl
     CALL trace_end("compute_caldyn_slow_hydro")    
-  END SUBROUTINE compute_caldyn_slow_hydro
+  END SUBROUTINE compute_caldyn_slow_hydro_manual
 
 END MODULE compute_caldyn_slow_hydro_mod
