--- trunk/phylmd/Interface_surf/calcul_fluxs.f	2018/07/12 17:53:18	278
+++ trunk/phylmd/Interface_surf/calcul_fluxs.f	2018/07/20 14:30:23	279
@@ -4,18 +4,21 @@
 
 contains
 
-  SUBROUTINE calcul_fluxs(dtime, tsurf, p1lay, cal, beta, coef1lay, ps, &
-       qsurf, radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, petAcoef, &
-       peqAcoef, petBcoef, peqBcoef, tsurf_new, evap, fluxlat, flux_t, &
-       dflux_s, dflux_l)
+  SUBROUTINE calcul_fluxs(dtime, tsurf, p1lay, cal, beta, coef1lay, ps, qsurf, &
+       radsol, dif_grnd, t1lay, q1lay, u1lay, v1lay, petAcoef, peqAcoef, &
+       petBcoef, peqBcoef, tsurf_new, evap, fluxlat, flux_t, dflux_s, dflux_l)
 
     ! Cette routine calcule les flux en h et q à l'interface et une
     ! température de surface.
 
     ! L. Fairhead, April 2000
 
-    USE fcttre, ONLY: foede, foeew
+    ! Note that, if cal = 0, beta = 1 and dif_grnd = 0, then tsurf_new
+    ! = tsurf and qsurf = qsat.
+
     use nr_util, only: assert_eq
+
+    USE fcttre, ONLY: foede, foeew
     USE suphec_m, ONLY: rcpd, rd, retv, rlstt, rlvtt, rtt
     USE yoethf_m, ONLY: r2es, r5ies, r5les, rvtmp2
 
@@ -57,7 +60,7 @@
 
     ! Local:
     integer i
-    integer knon ! nombre de points a traiter
+    integer knon ! nombre de points \`a traiter
     real, dimension(size(ps)):: mh, oh, mq, nq, oq, dq_s_dt, coef ! (knon)
     real qsat(size(ps)) ! (knon) mass fraction
     real sl(size(ps)) ! (knon) chaleur latente d'évaporation ou de sublimation
@@ -67,12 +70,12 @@
 
     !---------------------------------------------------------------------
 
-    knon = assert_eq((/size(tsurf), size(p1lay), size(cal), size(beta), &
+    knon = assert_eq([size(tsurf), size(p1lay), size(cal), size(beta), &
          size(coef1lay), size(ps), size(qsurf), size(radsol), size(dif_grnd), &
          size(t1lay), size(q1lay), size(u1lay), size(v1lay), size(petAcoef), &
          size(peqAcoef), size(petBcoef), size(peqBcoef), size(tsurf_new), &
          size(evap), size(fluxlat), size(flux_t), size(dflux_s), &
-         size(dflux_l)/), "calcul_fluxs knon")
+         size(dflux_l)], "calcul_fluxs knon")
 
     ! Traitement de l'humidité du sol
 
@@ -92,24 +95,20 @@
     ! Q
     oq = 1. - beta * coef * peqBcoef * dtime
     mq = beta * coef * (peqAcoef - qsat + dq_s_dt * tsurf) / oq
-    nq = beta * coef * (- 1. * dq_s_dt) / oq
+    nq = - beta * coef * dq_s_dt / oq
 
     ! H
-    oh = 1. - (coef * petBcoef * dtime)
+    oh = 1. - coef * petBcoef * dtime
     mh = coef * petAcoef / oh
-    dflux_s = - (coef * RCPD)/ oh
+    dflux_s = - coef * RCPD / oh
 
-    ! Tsurface
     tsurf_new = (tsurf + cal / RCPD * dtime * (radsol + mh + sl * mq) &
          + dif_grnd * t_grnd * dtime) / (1. - dtime * cal / RCPD * (dflux_s &
          + sl * nq) + dtime * dif_grnd)
-
     evap = - mq - nq * tsurf_new
     fluxlat = - evap * sl
     flux_t = mh + dflux_s * tsurf_new
     dflux_l = sl * nq
-
-    ! Nouvelle valeur de l'humidité au dessus du sol :
     qsurf = (peqAcoef - peqBcoef * evap * dtime) * (1. - beta) + beta * (qsat &
          + dq_s_dt * (tsurf_new - tsurf))