From f286215a8984c13c68af4b372604acec623c7d5e Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Tue, 7 May 2024 16:43:14 +0200
Subject: [PATCH] Jean-Luc R. 07/05/2024: minor update of atm-ocean LCOUPLES
 option (reading forcings, init and variables names)

---
 src/MNH/forcing.f90     | 11 ++++---
 src/MNH/ini_modeln.f90  | 37 +++++++++++++-----------
 src/MNH/modd_oceanh.f90 |  3 +-
 src/MNH/nhoa_coupln.f90 | 63 ++++++++++++++++++++---------------------
 src/MNH/phys_paramn.f90 | 37 ++++++++++++++----------
 src/MNH/read_field.f90  |  2 +-
 src/MNH/set_rsou.f90    | 51 ++++++++++++++++++++-------------
 7 files changed, 112 insertions(+), 92 deletions(-)

diff --git a/src/MNH/forcing.f90 b/src/MNH/forcing.f90
index 14a65ad60..522ac749e 100644
--- a/src/MNH/forcing.f90
+++ b/src/MNH/forcing.f90
@@ -163,6 +163,7 @@ use modd_budget,     only: lbudget_u,  lbudget_v,  lbudget_w,  lbudget_th, lbudg
 USE MODD_CONF
 USE MODD_CST
 USE MODD_DYN
+USE MODD_DYN_n, ONLY : LOCEAN
 USE MODD_FRC
 USE MODD_LUNIT
 USE MODD_PARAMETERS
@@ -737,11 +738,13 @@ END IF
 !
 !*       4.2.1    integration of the tendency forcing for uv
 !
-IF ( LTEND_UV_FRC ) THEN
- PRUS(:,:,:) = PRUS(:,:,:) +  MXM(PRHODJ) * ZTENDUF(:,:,:)
- PRVS(:,:,:) = PRVS(:,:,:) +  MYM(PRHODJ) * ZTENDVF(:,:,:)
+IF (.NOT. LOCEAN) THEN
+ IF ( LTEND_UV_FRC ) THEN
+  PRUS(:,:,:) = PRUS(:,:,:) +  MXM(PRHODJ) * ZTENDUF(:,:,:)
+  PRVS(:,:,:) = PRVS(:,:,:) +  MYM(PRHODJ) * ZTENDVF(:,:,:)
+ END IF
 END IF
-!
+ !
 !*       4.3    integration of the thermal and geostrophic wind
 !
 IF( LCORIO ) THEN
diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90
index fffa848fc..aad37b9c5 100644
--- a/src/MNH/ini_modeln.f90
+++ b/src/MNH/ini_modeln.f90
@@ -604,18 +604,19 @@ ILUOUT = TLUOUT%NLU
 !             --------------
 !*       2.1  Read number of forcing fields
 !
-IF (LFORCING) THEN ! Retrieve the number of time-dependent forcings.
+IF (.NOT. LOCEAN) THEN
+ IF (LFORCING) THEN ! Retrieve the number of time-dependent forcings.
   CALL IO_Field_read(TPINIFILE,'FRC',NFRC,IRESP)
   IF ( (IRESP /= 0) .OR. (NFRC <=0) ) THEN
     WRITE(ILUOUT,'(A/A)') &
      "INI_MODEL_n ERROR: you want to read forcing variables from FMfile", &
      "                   but no fields have been found by IO_Field_read"
-!callabortstop
+ !callabortstop
     CALL PRINT_MSG(NVERB_FATAL,'GEN','INI_MODEL_n','')
   END IF
-END IF
-!
-! Modif PP for time evolving adv forcing
+ END IF
+ !
+ ! Modif PP for time evolving adv forcing
   IF ( L2D_ADV_FRC ) THEN ! Retrieve the number of time-dependent forcings.
     WRITE(ILUOUT,FMT=*) "INI_MODEL_n ENTER ADV_FORCING"
     CALL IO_Field_read(TPINIFILE,'NADVFRC1',NADVFRC,IRESP)
@@ -627,9 +628,9 @@ END IF
       CALL PRINT_MSG(NVERB_FATAL,'GEN','INI_MODEL_n','')
     END IF
     WRITE(ILUOUT,*) 'NADVFRC = ', NADVFRC
-END IF
-!
-IF ( L2D_REL_FRC ) THEN ! Retrieve the number of time-dependent forcings.
+ END IF
+ !
+ IF ( L2D_REL_FRC ) THEN ! Retrieve the number of time-dependent forcings.
     WRITE(ILUOUT,FMT=*) "INI_MODEL_n ENTER REL_FORCING"
     CALL IO_Field_read(TPINIFILE,'NRELFRC1',NRELFRC,IRESP)
     IF ( (IRESP /= 0) .OR. (NRELFRC <=0) ) THEN
@@ -640,6 +641,7 @@ IF ( L2D_REL_FRC ) THEN ! Retrieve the number of time-dependent forcings.
       CALL PRINT_MSG(NVERB_FATAL,'GEN','INI_MODEL_n','')
     END IF
     WRITE(ILUOUT,*) 'NRELFRC = ', NRELFRC
+ END IF
 END IF
 !*       2.2  Checks the position of vertical absorbing layer
 !
@@ -2978,16 +2980,17 @@ END IF
 !
 !*     33.  Auto-coupling Atmos-Ocean LES NH
 !
-IF (LCOUPLES) THEN
- ALLOCATE(XSSUFL_C(IIU,IJU,1)); XSSUFL_C=0.0
- ALLOCATE(XSSVFL_C(IIU,IJU,1)); XSSVFL_C=0.0
- ALLOCATE(XSSTFL_C(IIU,IJU,1)); XSSTFL_C=0.0
- ALLOCATE(XSSRFL_C(IIU,IJU,1)); XSSRFL_C=0.
+! Atmos Flux at interface
+IF (LCOUPLES.AND.(.NOT.LOCEAN)) THEN
+ ALLOCATE(XSSUFL(IIU,IJU)); XSSUFL=0.0
+ ALLOCATE(XSSVFL(IIU,IJU)); XSSVFL=0.0
+ ALLOCATE(XSSTFL(IIU,IJU)); XSSTFL=0.0
+ ALLOCATE(XSSRFL(IIU,IJU)); XSSRFL=0.0
 ELSE
- ALLOCATE(XSSUFL_C(0,0,0))
- ALLOCATE(XSSVFL_C(0,0,0))
- ALLOCATE(XSSTFL_C(0,0,0))
- ALLOCATE(XSSRFL_C(0,0,0))
+ ALLOCATE(XSSUFL(0,0))
+ ALLOCATE(XSSVFL(0,0))
+ ALLOCATE(XSSTFL(0,0))
+ ALLOCATE(XSSRFL(0,0))
 END IF
 !
 END SUBROUTINE INI_MODEL_n
diff --git a/src/MNH/modd_oceanh.f90 b/src/MNH/modd_oceanh.f90
index e9936173e..0716105b6 100644
--- a/src/MNH/modd_oceanh.f90
+++ b/src/MNH/modd_oceanh.f90
@@ -40,8 +40,7 @@ IMPLICIT NONE
 INTEGER,          SAVE                  :: NFRCLT     ! number of sea surface forcings PLUS 1
 INTEGER,          SAVE                  :: NINFRT     ! Interval in second between forcings
 TYPE (DATE_TIME), SAVE, DIMENSION(:), ALLOCATABLE :: TFRCLT ! date/time of sea surface forcings
-REAL, SAVE, DIMENSION(:,:), ALLOCATABLE ::  XSSUFL,XSSVFL,XSSTFL,XSSOLA ! Time evol Flux U V T Solar_Rad at sea surface
-REAL, SAVE, DIMENSION(:,:), ALLOCATABLE :: XSSUFL_XY,XSSVFL_XY,XSSTFL_XY! XY flux shape
+REAL, SAVE, DIMENSION(:,:), ALLOCATABLE ::  XSSUFL,XSSVFL,XSSTFL,XSSRFL,XSSOLA ! Time evol Flux U,V,T,SW at sea surface
 REAL, SAVE, DIMENSION(:), ALLOCATABLE :: XSSUFL_T,XSSVFL_T,XSSTFL_T,XSSOLA_T ! given time forcing fluxes
 !
 END MODULE MODD_OCEANH
diff --git a/src/MNH/nhoa_coupln.f90 b/src/MNH/nhoa_coupln.f90
index 072dfbeb3..b156d2a4c 100644
--- a/src/MNH/nhoa_coupln.f90
+++ b/src/MNH/nhoa_coupln.f90
@@ -60,13 +60,15 @@ USE MODE_MODELN_HANDLER
 USE MODD_PARAMETERS
 USE MODD_NESTING
 USE MODD_CST
+USE MODD_REF, ONLY: XRHODREFZ
 USE MODD_REF_n        ! modules relative to the outer model $n
 USE MODD_FIELD_n
+USE MODD_METRICS_n
 USE MODD_CONF
 USE MODD_PARAM_n
 USE MODD_TURB_n
 USE MODD_DYN_n, ONLY : LOCEAN
-USE MODD_REF, ONLY: LCOUPLES
+USE MODD_OCEANH 
 !
 IMPLICIT NONE
 !
@@ -80,7 +82,7 @@ INTEGER,          INTENT(IN)    :: KKU      !
 !
 !*       0.2   declarations of local variables
 !
-INTEGER                :: IIB,IIE,IJB,IJE,IIU,IJU
+INTEGER                :: IIU,IJU
 INTEGER :: IKE
 INTEGER                :: ILBX,ILBY,ILBX2,ILBY2
 INTEGER,          INTENT(IN)    :: KTCOUNT  !  Temporal loop COUNTer
@@ -93,7 +95,8 @@ INTEGER :: IINFO_ll, IDIMX, IDIMY
 REAL, DIMENSION(:,:), ALLOCATABLE :: ZCOUPUA,ZCOUPVA,ZCOUPTA
 REAL, DIMENSION(:,:), ALLOCATABLE :: ZCOUPUO,ZCOUPVO,ZCOUPTO
 !surf flux  local work space
-REAL, DIMENSION(:,:), ALLOCATABLE :: ZCOUPTFL,ZCOUPUFL,ZCOUPVFL
+REAL, DIMENSION(:,:), ALLOCATABLE :: ZCOUPTFL,ZCOUPUFL,ZCOUPVFL,ZCOUPQFL
+REAL                              :: ZCDN ! CONSTANT DRAG COEFF
 CHARACTER(LEN=4)                    :: ZINIT_TYPE
 !
 !---Coupled OA MesoNH----------------------------------------------------------------------------
@@ -101,11 +104,7 @@ CHARACTER(LEN=4)                    :: ZINIT_TYPE
 !            --------------
 IIU = SIZE(XRHODJ,1)
 IJU = SIZE(XRHODJ,2)
-
-! allocate flux local array
-ALLOCATE(ZCOUPTFL(IIU,IJU))
-ALLOCATE(ZCOUPUFL(IIU,IJU))
-ALLOCATE(ZCOUPVFL(IIU,IJU))
+ZCDN = 1.2E-3
 ! allocate sfc variable local array
 ALLOCATE(ZCOUPUA(IIU,IJU))
 ALLOCATE(ZCOUPVA(IIU,IJU))
@@ -113,7 +112,11 @@ ALLOCATE(ZCOUPTA(IIU,IJU))
 ALLOCATE(ZCOUPUO(IIU,IJU))
 ALLOCATE(ZCOUPVO(IIU,IJU))
 ALLOCATE(ZCOUPTO(IIU,IJU))
-! values in ocean sfc
+ALLOCATE(ZCOUPTFL(IIU,IJU))
+ALLOCATE(ZCOUPUFL(IIU,IJU))
+ALLOCATE(ZCOUPVFL(IIU,IJU))
+ALLOCATE(ZCOUPQFL(IIU,IJU))
+! loading values in ocean sfc
 IKE=KKU-JPVEXT
 ZCOUPUO(:,:)= XUT(:,:,IKE)
 ZCOUPVO(:,:)= XVT(:,:,IKE)
@@ -121,37 +124,31 @@ ZCOUPTO(:,:)= XTHT(:,:,IKE)
 !
 ! we are going to the atmos model i.e. Model 1
 CALL GOTO_MODEL(KDAD)
-IIB=1
-IIE=IIU
-IJB=1
-IJE=IJU
-!
 ! compute gradient between ocean & atmosphere
 ZCOUPUA(:,:)= XUT(:,:,2)-ZCOUPUO(:,:)
 ZCOUPVA(:,:)= XVT(:,:,2)-ZCOUPVO(:,:)
 ZCOUPTA(:,:)= XTHT(:,:,2)-ZCOUPTO(:,:)
-!
-! sfc flux computation  * RHO AIR !!!!
-! flux vu atmosp
-!
-ZCOUPTFL(:,:) = -1.2*1.E-3* SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPTA(:,:)
-ZCOUPUFL(:,:) = -1.2*1.E-3* SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPUA(:,:)
-ZCOUPVFL(:,:) = -1.2*1.E-3* SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPVA(:,:)
-!
-XSSUFL_C(:,:,1)= ZCOUPUFL(:,:)
-XSSVFL_C(:,:,1)= ZCOUPVFL(:,:)
-XSSTFL_C(:,:,1)= ZCOUPTFL(:,:)
-!
-!
-! We are going back in the ocean model  
-! same sign & unit at the top of ocean model and the bottom of atmospheric model
+!!!!!!!!!!!!
+!  Simple neutral bulk flux with Z0 fixed
+ZCOUPTFL(:,:) = -ZCDN*SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPTA(:,:)&
+                *XCPD*XRHODREFZ(2)        ! in W/M2 (as WASP)              
+ZCOUPUFL(:,:) = -ZCDN*XRHODREFZ(2)* SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPUA(:,:)
+ZCOUPVFL(:,:) = -ZCDN*XRHODREFZ(2)* SQRT(ZCOUPUA(:,:)**2 +ZCOUPVA(:,:)**2) * ZCOUPVA(:,:)
+!  EVAP-PR TO BE ADDED HERE 
+ZCOUPQFL(:,:)= 0.
+! keep in W/M2, N/m2, ..
+XSSUFL= ZCOUPUFL
+XSSVFL= ZCOUPVFL
+XSSTFL= ZCOUPTFL
+XSSRFL= ZCOUPQFL
+!
+!  back in the ocean model  
+!oce flux are directly computed in phys_param routine
+! same sign & unit fluxes at top of ocean model & bottom of atmospheric model
 !  rho_atmos * (w'u')_atmos = rho_ocean * (u'w')_ocean
-!  rho_atmos *Cp_atmos* (u'w')_atmos = rho_ocean * CP_ocean * (u'w')_ocean
+!  rho_atmos *Cp_atmos* (T'w')_atmos = rho_ocean * CP_ocean * (T'w')_ocean
 !
 CALL GOTO_MODEL(KMI)
-XSSUFL_C(:,:,1)= ZCOUPUFL(:,:)/XRH00OCEAN
-XSSVFL_C(:,:,1)= ZCOUPVFL(:,:)/XRH00OCEAN
-XSSTFL_C(:,:,1)= ZCOUPTFL(:,:)*1004./(3900.*XRH00OCEAN)
 DEALLOCATE(ZCOUPUA,ZCOUPVA,ZCOUPUO,ZCOUPVO,ZCOUPTA,ZCOUPTO)
 DEALLOCATE(ZCOUPTFL)
 !
diff --git a/src/MNH/phys_paramn.f90 b/src/MNH/phys_paramn.f90
index edd731e42..ec037233c 100644
--- a/src/MNH/phys_paramn.f90
+++ b/src/MNH/phys_paramn.f90
@@ -317,7 +317,7 @@ USE MODD_PRECIP_n
 use modd_precision,        only: MNHTIME
 USE MODD_RADIATIONS_n
 USE MODD_RAIN_ICE_DESCR_n, ONLY: XRTMIN
-USE MODD_REF,              ONLY: LCOUPLES
+USE MODD_REF,              ONLY: LCOUPLES,XRHODREFZ
 USE MODD_REF_n
 USE MODD_SALT
 USE MODD_SHADOWS_n
@@ -865,24 +865,32 @@ END IF
 ! Sfc turbulent fluxes & Radiative tendency due to SW penetrating ocean
 ! 
 IF (LCOUPLES) THEN
-ZSFU(:,:)= XSSUFL_C(:,:,1)
-ZSFV(:,:)= XSSVFL_C(:,:,1)
-ZSFTH(:,:)= XSSTFL_C(:,:,1)
-ZSFRV(:,:)=XSSRFL_C(:,:,1)
-ELSE 
-IF (LOCEAN) THEN
-!
+! Flux computed in NHOA/WASP and are given in W/M2
+! convert in SI as in MESO-NH & for Ocean & Atmos model
+ IF (LOCEAN) THEN
+  ZSFU=  XSSUFL/CST%XRH00OCEAN
+  ZSFV=  XSSVFL/CST%XRH00OCEAN
+  ZSFTH= XSSTFL/(CST%XCL*CST%XRH00OCEAN)
+!to do: to add factor s/(1-s)
+  ZSFRV=-XSSRFL/CST%XRH00OCEAN
+ ELSE
+  ZSFU= XSSUFL/XRHODREFZ(2)
+  ZSFV= XSSVFL/XRHODREFZ(2)
+  ZSFTH=XSSTFL/(XCPD*XRHODREFZ(2))
+  ZSFRV=XSSRFL/(CST%XLVTT*XRHODREFZ(2))
+ ENDIF 
+ELSE ! no coupled 
+ IF (LOCEAN) THEN
   ALLOCATE( ZIZOCE(IKU)); ZIZOCE(:)=0. 
   ALLOCATE( ZPROSOL1(IKU))
   ALLOCATE( ZPROSOL2(IKU))
-  ALLOCATE(XSSOLA(IIU,IJU))
   ! Time interpolation
   JSW     = INT(TDTCUR%xtime/REAL(NINFRT))
   ZSWA    = TDTCUR%xtime/REAL(NINFRT)-REAL(JSW)
   ZSFRV = 0.
   ZSFTH  = (XSSTFL_T(JSW+1)*(1.-ZSWA)+XSSTFL_T(JSW+2)*ZSWA) 
-  ZSFU = (XSSUFL_T(JSW+1)*(1.-ZSWA)+XSSUFL_T(JSW+2)*ZSWA)
-  ZSFV = (XSSVFL_T(JSW+1)*(1.-ZSWA)+XSSVFL_T(JSW+2)*ZSWA)
+  ZSFU   = (XSSUFL_T(JSW+1)*(1.-ZSWA)+XSSUFL_T(JSW+2)*ZSWA)
+  ZSFV   = (XSSVFL_T(JSW+1)*(1.-ZSWA)+XSSVFL_T(JSW+2)*ZSWA)
 !
   ZIZOCE(IKU)   = XSSOLA_T(JSW+1)*(1.-ZSWA)+XSSOLA_T(JSW+2)*ZSWA
   ZPROSOL1(IKU) = CST%XROC*ZIZOCE(IKU)
@@ -896,12 +904,11 @@ IF (LOCEAN) THEN
     XRTHS(:,:,JKM) = XRTHS(:,:,JKM) + XRHODJ(:,:,JKM)*ZIZOCE(JKM)
   END DO
   if ( TBUCONF%LBUDGET_th ) call Budget_store_end ( TBUDGETS(NBUDGET_TH), 'OCEAN', xrths(:, :, :) )
-  DEALLOCATE (XSSOLA)
   DEALLOCATE( ZIZOCE) 
   DEALLOCATE (ZPROSOL1)
   DEALLOCATE (ZPROSOL2)
-END IF! LOCEAN NO LCOUPLES
-END IF!NO LCOUPLES
+END IF! TEST LOCEAN NO LCOUPLES
+END IF! TEST LCOUPLES
 !
 !
 CALL SECOND_MNH2(ZTIME2)
@@ -1324,7 +1331,7 @@ ELSE ! case no SURFEX (CSURF logical)
   ZCD_ROOF   = 0.
   ZSFRV_WALL = 0.
   ZSFRV_ROOF = 0.
-  IF (.NOT.LOCEAN) THEN
+  IF (.NOT.(LOCEAN.OR.LCOUPLES))THEN
     ZSFTH    = 0.
     ZSFRV    = 0.
     ZSFSV    = 0.
diff --git a/src/MNH/read_field.f90 b/src/MNH/read_field.f90
index 877437197..412e19bf3 100644
--- a/src/MNH/read_field.f90
+++ b/src/MNH/read_field.f90
@@ -1321,7 +1321,7 @@ IF (LOCEAN .AND. (.NOT.LCOUPLES) .AND. (KOCEMI==1)) THEN
     CLONGNAME  = 'SSOLA',                             &
     CUNITS     = 'kg m3 K m s-1',                     &
     CDIR       = '--',                                &
-    CCOMMENT   = 'sfc solar flux to force ocean LES', &
+    CCOMMENT   = 'sfc solar flux at sfc to force ocean LES', &
     NGRID      = 0,                                   &
     NTYPE      = TYPEREAL,                            &
     NDIMS      = 1,                                   &
diff --git a/src/MNH/set_rsou.f90 b/src/MNH/set_rsou.f90
index 452625183..e9a4dc78f 100644
--- a/src/MNH/set_rsou.f90
+++ b/src/MNH/set_rsou.f90
@@ -272,6 +272,7 @@ USE MODD_NETCDF
 USE MODD_OCEANH
 USE MODD_PARAMETERS, ONLY: JPHEXT
 USE MODD_TYPE_DATE
+USE MODD_TIME_n,     ONLY : TDTCUR
 !
 USE MODE_ll
 USE MODE_MSG
@@ -400,7 +401,8 @@ REAL, DIMENSION(:),     ALLOCATABLE :: ZOC_DEPTH
 REAL, DIMENSION(:),     ALLOCATABLE :: ZOC_LE,ZOC_H
 REAL, DIMENSION(:),     ALLOCATABLE :: ZOC_SW_DOWN,ZOC_SW_UP,ZOC_LW_DOWN,ZOC_LW_UP
 REAL, DIMENSION(:),     ALLOCATABLE :: ZOC_TAUX,ZOC_TAUY
-
+!
+REAL  :: ZJZTIME ! TIME(HOUR) READ in PRE_IDEA1.NAM
 !--------------------------------------------------------------------------------
 !
 !*	 1.     PROLOGUE : INITIALIZE SOME CONSTANTS, RETRIEVE LOGICAL
@@ -565,9 +567,10 @@ SELECT CASE(YKIND)
       TFRCLT(JKT)= ZFRCLT(JKT)
       XSSUFL_T(JKT)=ZSSUFL_T(JKT)/XRH00OCEAN
       XSSVFL_T(JKT)=ZSSVFL_T(JKT)/XRH00OCEAN
-      ! working in SI
-      XSSTFL_T(JKT)=ZSSTFL_T(JKT) /(3900.*XRH00OCEAN)
-      XSSOLA_T(JKT)=ZSSOLA_T(JKT) /(3900.*XRH00OCEAN)
+                 
+    ! working in SI
+      XSSTFL_T(JKT)=ZSSTFL_T(JKT) /(XCL*XRH00OCEAN)
+      XSSOLA_T(JKT)=ZSSOLA_T(JKT) /(XCL*XRH00OCEAN)
     END DO   
     DEALLOCATE(ZFRCLT)
     DEALLOCATE(ZSSUFL_T)
@@ -579,11 +582,12 @@ SELECT CASE(YKIND)
 ! 2.0.2  Ocean standard initialize from netcdf files
 !        U,V,T,S at Z levels + Forcings at model TOP (sea surface) 
 !--------------------------------------------------------------------------------   
-!
+!  WRITTEN TO READ NETCDF file for a DWL case of DYNAMO-CINDY
   CASE ('STANDOCE')
-!   
     XP00=XP00OCEAN
-    READ(ILUPRE,*) ZPTOP           ! P_atmosphere at sfc =P top domain
+! HOUR OF FIRST FORCING on 14/11 TO USE AS READ in PRE_IDEA1.NAM
+    ZJZTIME= INT(TDTCUR%xtime)
+    READ(ILUPRE,*) ZPTOP           ! P_atmosphere at sfc =P top ocean domain
     READ(ILUPRE,*) YINFILE, YINFISF
     WRITE(ILUOUT,FMT=*) 'Netcdf files to read:', YINFILE, YINFISF
     ! Open file containing initial profiles
@@ -670,6 +674,8 @@ SELECT CASE(YKIND)
       ! ZHEIGHT used only in set_ rsou, defined as such ZHEIGHT(ILEVELM)=H_model
       ! Z oriented in same time to have a model domain axis going
       ! from 0m (ocean bottom/model bottom) towards H (ocean sfc/model top)
+      WRITE(ILUOUT,FMT=*) 'End gridmodel comput in trans domain: JKM  U V ZHEIGHTU', &
+      JKM,ZU(JKM),ZV(JKM),ZHEIGHTU(JKM)
     END DO
 !
     DEALLOCATE(ZOC_TEMPERATURE)
@@ -728,30 +734,36 @@ SELECT CASE(YKIND)
       WRITE(ILUOUT,FMT=*) JKM, ZOC_LE(JKM), ZOC_H(JKM),ZOC_SW_DOWN(JKM),ZOC_SW_UP(JKM),&
                           ZOC_LW_DOWN(JKM),ZOC_LW_UP(JKM),ZOC_TAUX(JKM),ZOC_TAUY(JKM)   
     ENDDO
-    ! IFRCLT FORCINGS at sea surface
-    IFRCLT=idimlen
+    ! IFRCLT FORCINGS at sea surface are givean each 10 min (6/hour)
+! start at ZJZTIME read in PRE_IDEA1.nam; skip forcings given from 13/11 22LT
+ZJZTIME=(ZJZTIME+2)*6
+   IFRCLT=idimlen-ZJZTIME
+WRITE(ILUOUT,FMT=*) 'FORCINGS GIVEN FROM 13/11/11 22LT;  FORCING Number USED ZJZTIME=', ZJZTIME
+WRITE(ILUOUT,FMT=*) idimlen,'NB of given forcing', 'nb of first used forcing= ',IFRCLT
     ALLOCATE(ZFRCLT(IFRCLT)) 
     ALLOCATE(ZSSUFL_T(IFRCLT)); ZSSUFL_T = 0.0
     ALLOCATE(ZSSVFL_T(IFRCLT)); ZSSVFL_T = 0.0
     ALLOCATE(ZSSTFL_T(IFRCLT)); ZSSTFL_T = 0.0
     ALLOCATE(ZSSOLA_T(IFRCLT)); ZSSOLA_T = 0.0
-    DO JKT=1,IFRCLT
+    DO JKM=1,IFRCLT
+    JKT=JKM +ZJZTIME
       ! Initial file for CINDY-DYNAMO: all fluxes correspond to the absolute value (>0)
       ! modele ocean: axe z dirigé du bas vers la sfc de l'océan
       ! => flux dirigé vers le haut (positif ocean vers l'atmopshere i.e. bas vers le haut)
-      ZSSOLA_T(JKT)=ZOC_SW_DOWN(JKT)-ZOC_SW_UP(JKT)
-      ZSSTFL_T(JKT)=(ZOC_LW_DOWN(JKT)-ZOC_LW_UP(JKT)-ZOC_LE(JKT)-ZOC_H(JKT))
+      ZSSOLA_T(JKM)=ZOC_SW_DOWN(JKT)-ZOC_SW_UP(JKT)
+      ZSSTFL_T(JKM)=-(ZOC_LW_DOWN(JKT)-ZOC_LW_UP(JKT)-ZOC_LE(JKT)-ZOC_H(JKT))
       ! assume that Tau given on file is along Ox
       ! rho_air UW_air = rho_ocean UW_ocean= N/m2
       ! uw_ocean
-      ZSSUFL_T(JKT)=ZOC_TAUX(JKT)
-      ZSSVFL_T(JKT)=ZOC_TAUY(JKT)
-      WRITE(ILUOUT,FMT=*) 'Forcing Nb Sol NSol UW_oc VW',&
-                          JKT,ZSSOLA_T(JKT),ZSSTFL_T(JKT),ZSSUFL_T(JKT),ZSSVFL_T(JKT) 
+      ZSSUFL_T(JKM)=ZOC_TAUX(JKT)
+      ZSSVFL_T(JKM)=ZOC_TAUY(JKT) 
+     WRITE(ILUOUT,FMT=*) 'Forcing Nb Sol NSol UW_oc VW',&
+                          JKM,ZSSOLA_T(JKM),ZSSTFL_T(JKM),ZSSUFL_T(JKM),ZSSVFL_T(JKM) 
     ENDDO
     ! Allocate and Writing the corresponding variables in module MODD_OCEAN_FRC
     NFRCLT=IFRCLT
     ! value to read later on file ? 
+! interval in s between 2 given forcings
     NINFRT=600
     ALLOCATE(TFRCLT(NFRCLT))
     ALLOCATE(XSSUFL_T(NFRCLT));XSSUFL_T(:)=0.
@@ -759,14 +771,13 @@ SELECT CASE(YKIND)
     ALLOCATE(XSSTFL_T(NFRCLT));XSSTFL_T(:)=0.
     ALLOCATE(XSSOLA_T(NFRCLT));XSSOLA_T(:)=0.
     ! on passe en unités SI, signe, etc pour le modele ocean
-    !  W/m2 => SI :  /(CP_mer * rho_mer)
-    ! a revoir dans tt le code pour mettre de svaleurs plus exactes
+    !  W/m2 => SI : m/S *K  /(CP_mer * rho_mer)
     DO JKT=1,NFRCLT
       TFRCLT(JKT)= ZFRCLT(JKT)
       XSSUFL_T(JKT)=ZSSUFL_T(JKT)/XRH00OCEAN
       XSSVFL_T(JKT)=ZSSVFL_T(JKT)/XRH00OCEAN
-      XSSTFL_T(JKT)=ZSSTFL_T(JKT) /(3900.*XRH00OCEAN)
-      XSSOLA_T(JKT)=ZSSOLA_T(JKT) /(3900.*XRH00OCEAN)
+      XSSTFL_T(JKT)=ZSSTFL_T(JKT) /(XCL  *XRH00OCEAN)
+      XSSOLA_T(JKT)=ZSSOLA_T(JKT) /(XCL  *XRH00OCEAN)
     END DO   
     DEALLOCATE(ZFRCLT)
     DEALLOCATE(ZSSUFL_T)
-- 
GitLab