From b95d84d7d992609c0f0d8a76ad2a3e8f5b5bf9a3 Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Wed, 3 Jan 2024 12:38:48 +0100
Subject: [PATCH] V.Masson 03/01/2024: optimization of NRAD_AGG options (now
 radiations is not called on already-called columns when using the
 aggregation)

---
 src/MNH/default_desfmn.f90     |   3 +-
 src/MNH/ini_modeln.f90         |   8 ++-
 src/MNH/ini_radiations.f90     |  13 ++--
 src/MNH/ini_radiations_agg.f90 | 126 ++++++++++++++++++++++++++++++++-
 src/MNH/modd_param_radn.f90    |   7 +-
 src/MNH/modd_radiationsn.f90   |   8 ++-
 src/MNH/modn_param_radn.f90    |   7 +-
 src/MNH/phys_paramn.f90        |   6 +-
 src/MNH/radiations.f90         |   1 +
 src/MNH/radiations_agg.f90     |  69 ++++++++++++++++--
 src/MNH/read_exsegn.f90        |   1 +
 11 files changed, 217 insertions(+), 32 deletions(-)

diff --git a/src/MNH/default_desfmn.f90 b/src/MNH/default_desfmn.f90
index 9ccf63878..d3a764ee9 100644
--- a/src/MNH/default_desfmn.f90
+++ b/src/MNH/default_desfmn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1994-2023 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -226,6 +226,7 @@ END MODULE MODI_DEFAULT_DESFM_n
 !  Delbeke/Vie    03/2022: KHKO option in LIMA
 !  P. Wautelet 27/04/2022: add namelist for profilers
 !  PA. Joulin     04/2023: add EOL/ADR
+!  V. Masson      01/2024: aggregation of columns for radiation
 !-------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90
index 4277d47e8..729df6d04 100644
--- a/src/MNH/ini_modeln.f90
+++ b/src/MNH/ini_modeln.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1994-2023 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -299,6 +299,7 @@ END MODULE MODI_INI_MODEL_n
 !  A. Costes      12/2021: Blaze fire model
 !  H. Toumi       09/2022: add EOL/ADR
 !  C. Barthe      03/2023: if cloud electricity is activated, both ini_micron and ini_elecn are called
+!  V. Masson      01/2024: aggregation of columns for radiation
 !---------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
@@ -1556,6 +1557,7 @@ IF (CRAD /= 'NONE') THEN
   ALLOCATE(XDTHRADSW(IIU,IJU,IKU))
   ALLOCATE(XDTHRADLW(IIU,IJU,IKU))
   ALLOCATE(XRADEFF(IIU,IJU,IKU))
+  ALLOCATE(NRAD_AGG_FLAG(IIU,IJU))
 ELSE
   ALLOCATE(XSLOPANG(0,0))
   ALLOCATE(XSLOPAZI(0,0))
@@ -1575,6 +1577,7 @@ ELSE
   ALLOCATE(XDTHRADSW(0,0,0))
   ALLOCATE(XDTHRADLW(0,0,0))
   ALLOCATE(XRADEFF(0,0,0))
+  ALLOCATE(NRAD_AGG_FLAG(0,0))
 END IF
 
 IF (CRAD == 'ECMW' .OR. CRAD == 'ECRA') THEN
@@ -2535,7 +2538,8 @@ IF (CRAD   /= 'NONE') THEN
                       XRADEFF,XSWU,XSWD,XLWU,             &
                       XLWD,XDTHRADSW,XDTHRADLW,           &
                       NRAD_AGG,NI_RAD_AGG,NJ_RAD_AGG,     &
-                      NIOR_RAD_AGG,NJOR_RAD_AGG           )
+                      NIOR_RAD_AGG,NJOR_RAD_AGG,          &
+                      NRAD_AGG_FLAG                       )
   !
   IF (GINIRAD) CALL SUNPOS_n(XZENITH,PAZIMSOL=XAZIM)
   CALL SURF_SOLAR_GEOM    (XZS, XZS_XY)
diff --git a/src/MNH/ini_radiations.f90 b/src/MNH/ini_radiations.f90
index 2b11a5890..a19a3b71f 100644
--- a/src/MNH/ini_radiations.f90
+++ b/src/MNH/ini_radiations.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 2003-2020 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2003-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -16,7 +16,8 @@ INTERFACE
          PFLALWD,PDIRSRFSWD,KCLEARCOL_TM1,                             &
          PZENITH, PAZIM, TPDTRAD_FULL,TPDTRAD_CLONLY,TPINITHALO2D_ll,  &
          PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW,              &
-         KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG      )
+         KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,     &
+         KRAD_AGG_FLAG                                                 )
 !
 USE MODD_ARGSLIST_ll, ONLY : LIST_ll
 USE MODD_IO,       ONLY : TFILEDATA
@@ -66,6 +67,7 @@ INTEGER, INTENT(OUT) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(OUT) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(OUT) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(OUT) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(OUT) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
 
 END SUBROUTINE INI_RADIATIONS
 !
@@ -82,7 +84,8 @@ END MODULE MODI_INI_RADIATIONS
          PFLALWD,PDIRSRFSWD,KCLEARCOL_TM1,                             &
          PZENITH, PAZIM, TPDTRAD_FULL,TPDTRAD_CLONLY,TPINITHALO2D_ll,  &
          PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW,              &
-         KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG      )
+         KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,     &
+         KRAD_AGG_FLAG                                                 )
 !   ####################################################################
 !
 !!****  *INI_RADIATION_TIME * - initialisation for radiation scheme in the MesoNH framework
@@ -119,6 +122,7 @@ END MODULE MODI_INI_RADIATIONS
 !  P. Wautelet 14/02/2019: remove CLUOUT/CLUOUT0 and associated variables
 !  P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
 !  P. Wautelet 20/05/2019: add name argument to ADDnFIELD_ll + new ADD4DFIELD_ll subroutine
+!       V. Masson 03/01/2024: aggregation of columns for radiation
 !-------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
@@ -187,6 +191,7 @@ INTEGER, INTENT(OUT) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(OUT) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(OUT) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(OUT) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(OUT) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
 !
 !*       0.2   declarations of local variables
 !
@@ -343,7 +348,7 @@ END IF
 !*       10.         INITIALIZE COLUMN AGGREGATION FOR RADIATION CALL
 !	            -------------------------------------------------
 
-CALL INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+CALL INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,KRAD_AGG_FLAG)
 !
 !-------------------------------------------------------------------------------
 !
diff --git a/src/MNH/ini_radiations_agg.f90 b/src/MNH/ini_radiations_agg.f90
index 0ec198677..71679e0dd 100644
--- a/src/MNH/ini_radiations_agg.f90
+++ b/src/MNH/ini_radiations_agg.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1995-2022 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2023-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -9,12 +9,13 @@
 !
 INTERFACE
 
-    SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+    SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,KRAD_AGG_FLAG)
 INTEGER, INTENT(IN)  :: KRAD_AGG      ! number of aggregated points
 INTEGER, INTENT(OUT) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(OUT) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(OUT) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(OUT) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(OUT) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
 END SUBROUTINE INI_RADIATIONS_AGG
 
 END INTERFACE
@@ -22,7 +23,7 @@ END INTERFACE
 END MODULE MODI_INI_RADIATIONS_AGG
 !
 !   ############################################################################
-    SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+    SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,KRAD_AGG_FLAG)
 !   ############################################################################
 !
 !!****  *INI_RADIATIONS_AGG * - routine to call the SW and LW radiation calculations
@@ -95,6 +96,7 @@ INTEGER, INTENT(OUT) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(OUT) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(OUT) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(OUT) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(OUT) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
 !
 !
 !*       0.2   DECLARATIONS OF LOCAL VARIABLES
@@ -119,6 +121,9 @@ INTEGER :: IJOR_ll      ! index of first point in the processor relative to the
 !
 INTEGER :: ILUOUT       ! Logical unit 
 INTEGER :: IMI
+
+LOGICAL :: LREMOVE_SOUTH, LREMOVE_NORTH, LREMOVE_EAST, LREMOVE_WEST ! flags to not keep packed column in current processor
+INTEGER :: ISOUTH, INORTH, IEAST, IWEST ! inner limits of packed colums on each side of the processor (in local coordinate)
 !-------------------------------------------------------------------------
 !-------------------------------------------------------------------------
 !-------------------------------------------------------------------------
@@ -191,6 +196,7 @@ KJOR_RAD_AGG = IJOR_RAD_AGG_ll - IJOR_ll + 1
 KIOR_RAD_AGG = KIOR_RAD_AGG + ((IIB-KIOR_RAD_AGG)/KRAD_AGG) * KRAD_AGG
 KJOR_RAD_AGG = KJOR_RAD_AGG + ((IJB-KJOR_RAD_AGG)/KRAD_AGG) * KRAD_AGG
 !
+!-------------------------------------------------------------------------------
 !
 ! Number of PACKED radiative subdomains inside current processor domain
 KI_RAD_AGG = (IIE - KIOR_RAD_AGG) / KRAD_AGG + 1
@@ -198,5 +204,119 @@ KJ_RAD_AGG = (IJE - KJOR_RAD_AGG) / KRAD_AGG + 1
 !
 !-------------------------------------------------------------------------------
 !
+! REMOVES Aggregated columns that are duplicated in several processors
+!
+! Checks if the middle of the packed column is in the current processor
+! 
+LREMOVE_SOUTH = .FALSE.
+LREMOVE_NORTH = .FALSE.
+LREMOVE_EAST  = .FALSE.
+LREMOVE_WEST  = .FALSE.
+!
+KRAD_AGG_FLAG(:,:) = 0.
+!
+! inner limits of packed colums on each side of the processor (in local coordinate)
+IWEST  = KIOR_RAD_AGG +                KRAD_AGG-1
+IEAST  = KIOR_RAD_AGG + (KI_RAD_AGG-1)*KRAD_AGG
+ISOUTH = KJOR_RAD_AGG +                KRAD_AGG-1
+INORTH = KJOR_RAD_AGG + (KJ_RAD_AGG-1)*KRAD_AGG
+!
+! Eastern side of processor (if NOT of whole domain) (checks on last X index of aggregated columns)
+IF (IEAST+KRAD_AGG/2 > IIE .AND. .NOT. (LEAST_ll() .AND. CLBCX(2)=='OPEN') ) THEN
+   KRAD_AGG_FLAG(IEAST:IIE,:) = 3  ! points located there will take values computed by processor towards east
+   LREMOVE_EAST = .TRUE.
+END IF
+! Western side of processor (if NOT of whole domain) (checks on first X index of aggregated columns)
+IF (KIOR_RAD_AGG+KRAD_AGG/2 < IIB  .AND. .NOT. (LWEST_ll() .AND. CLBCX(1)=='OPEN') ) THEN
+   KRAD_AGG_FLAG(IIB:IWEST,:) = 1  ! points located there will take values computed by processor towards west
+   LREMOVE_WEST = .TRUE.
+END IF
+! Northern side of processor (if NOT of whole domain) (checks on last Y index of aggregated columns)
+IF (INORTH+KRAD_AGG/2 > IJE .AND. .NOT. (LNORTH_ll() .AND. CLBCY(2)=='OPEN') ) THEN
+   KRAD_AGG_FLAG(:,INORTH:IJE) = 4  ! points located there will take values computed by processor towards north
+   LREMOVE_NORTH= .TRUE.
+END IF
+! Southern side of processor (if NOT of whole domain) (checks on first X index of aggregated columns)
+IF (KJOR_RAD_AGG+KRAD_AGG/2 < IJB  .AND. .NOT. (LSOUTH_ll() .AND. CLBCY(1)=='OPEN') ) THEN
+   KRAD_AGG_FLAG(:,IJB:ISOUTH) = 2  ! points located there will take values computed by processor towards south
+   LREMOVE_SOUTH= .TRUE.
+END IF
+!
+! North-Eastern corner of processor (if NOT of whole domain)
+IF (       IEAST+KRAD_AGG/2 > IIE .AND.  INORTH+KRAD_AGG/2 > IJE  ) THEN
+   IF (.NOT. (LEAST_ll()  .AND. CLBCX(2)=='OPEN') .AND. .NOT. (LNORTH_ll() .AND. CLBCY(2)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IEAST:IIE,INORTH:IJE) = 7  
+     ! points located there will take values computed by processor towards NE
+     LREMOVE_EAST  = .TRUE.
+     LREMOVE_NORTH = .TRUE.
+   ELSE IF ( (LEAST_ll()  .AND. CLBCX(2)=='OPEN') .AND. .NOT. (LNORTH_ll() .AND. CLBCY(2)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(:,INORTH:IJE) = 4  ! points located there will take values computed by processor towards north
+     LREMOVE_NORTH = .TRUE.
+   ELSE IF ( (LNORTH_ll()  .AND. CLBCY(2)=='OPEN') .AND. .NOT. (LEAST_ll()  .AND. CLBCX(2)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IEAST:IIE,:) = 3  ! points located there will take values computed by processor towards east
+     LREMOVE_EAST = .TRUE.
+   END IF
+END IF
+!
+! North-Western corner of processor (if NOT of whole domain)
+IF (       KIOR_RAD_AGG+KRAD_AGG/2 < IIB .AND.  INORTH+KRAD_AGG/2 > IJE  ) THEN
+   IF (.NOT. (LWEST_ll()  .AND. CLBCX(1)=='OPEN') .AND. .NOT. (LNORTH_ll() .AND. CLBCY(2)=='OPEN')) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,INORTH:IJE) = 8  ! points located there will take values computed by processor towards NW
+     LREMOVE_WEST  = .TRUE.
+     LREMOVE_NORTH = .TRUE.
+   ELSE IF ( (LWEST_ll()  .AND. CLBCX(1)=='OPEN') .AND. .NOT. (LNORTH_ll() .AND. CLBCY(2)=='OPEN')) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,INORTH:IJE) = 4  ! points located there will take values computed by processor towards north
+     LREMOVE_NORTH = .TRUE.
+   ELSE IF ( (LNORTH_ll()  .AND. CLBCY(2)=='OPEN') .AND. .NOT. (LWEST_ll()  .AND. CLBCX(1)=='OPEN')) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,INORTH:IJE) = 1  ! points located there will take values computed by processor towards west
+     LREMOVE_WEST = .TRUE.
+   END IF
+END IF
+!
+! South-Western corner of processor (if NOT of whole domain)
+IF (       KIOR_RAD_AGG+KRAD_AGG/2 < IIB .AND.  KJOR_RAD_AGG+KRAD_AGG/2 < IJB ) THEN
+   IF (.NOT. (LWEST_ll()  .AND. CLBCX(1)=='OPEN') .AND. .NOT. (LSOUTH_ll() .AND. CLBCY(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,IIB:ISOUTH) = 5  ! points located there will take values computed by processor towards SW
+     LREMOVE_WEST  = .TRUE.
+     LREMOVE_SOUTH = .TRUE.
+   ELSE IF ( (LWEST_ll()  .AND. CLBCX(1)=='OPEN')  .AND. .NOT. (LSOUTH_ll() .AND. CLBCY(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,IIB:ISOUTH) = 2  ! points located there will take values computed by processor towards south
+     LREMOVE_SOUTH = .TRUE.
+   ELSE IF ( (LSOUTH_ll()  .AND. CLBCY(2)=='OPEN') .AND. .NOT. (LWEST_ll()  .AND. CLBCX(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IIB:IWEST,IIB:ISOUTH) = 1  ! points located there will take values computed by processor towards west
+     LREMOVE_WEST = .TRUE.
+   END IF
+END IF
+!
+! South-Eastern corner of processor (if NOT of whole domain)
+IF (       IEAST+KRAD_AGG/2 > IIE .AND.  KJOR_RAD_AGG+KRAD_AGG/2 < IJB ) THEN
+   IF (.NOT. (LEAST_ll()  .AND. CLBCX(1)=='OPEN') .AND. .NOT. (LSOUTH_ll() .AND. CLBCY(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IEAST:IIE,IIB:ISOUTH) = 6  ! points located there will take values computed by processor towards SW
+     LREMOVE_EAST  = .TRUE.
+     LREMOVE_SOUTH = .TRUE.
+   ELSE IF ( (LEAST_ll()  .AND. CLBCX(1)=='OPEN')  .AND. .NOT. (LSOUTH_ll() .AND. CLBCY(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IEAST:IIE,IIB:ISOUTH) = 2  ! points located there will take values computed by processor towards south
+     LREMOVE_SOUTH = .TRUE.
+   ELSE IF ( (LSOUTH_ll()  .AND. CLBCY(2)=='OPEN') .AND. .NOT. (LEAST_ll()  .AND. CLBCX(1)=='OPEN') ) THEN
+     KRAD_AGG_FLAG(IEAST:IIE,IIB:ISOUTH) = 3  ! points located there will take values computed by processor towards west
+     LREMOVE_EAST = .TRUE.
+   END IF
+END IF
+!
+! removes from current processor the column that was partially (and majoritaly) in the other processor
+!
+IF (LREMOVE_EAST) KI_RAD_AGG = KI_RAD_AGG -1 
+IF (LREMOVE_WEST) THEN
+   KIOR_RAD_AGG = KIOR_RAD_AGG + KRAD_AGG
+   KI_RAD_AGG = KI_RAD_AGG -1
+END IF
+IF (LREMOVE_NORTH) KJ_RAD_AGG = KJ_RAD_AGG -1
+IF (LREMOVE_SOUTH) THEN
+   KJOR_RAD_AGG = KJOR_RAD_AGG + KRAD_AGG
+   KJ_RAD_AGG = KJ_RAD_AGG -1
+END IF
+!
+!-------------------------------------------------------------------------------
+!
 END SUBROUTINE INI_RADIATIONS_AGG
 
diff --git a/src/MNH/modd_param_radn.f90 b/src/MNH/modd_param_radn.f90
index f429c6938..613a2de1e 100644
--- a/src/MNH/modd_param_radn.f90
+++ b/src/MNH/modd_param_radn.f90
@@ -1,12 +1,8 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
 !MNH_LIC for details. version 1.
 !-----------------------------------------------------------------
-!--------------- special set of characters for RCS information
-!-----------------------------------------------------------------
-! $Source$ $Revision$
-! MASDEV4_7 modd 2006/11/23 17:28:26
 !-----------------------------------------------------------------
 !     ######################## 
       MODULE MODD_PARAM_RAD_n
@@ -43,6 +39,7 @@
 !!      F.Solmon      15/03/02  add the control parameter for aerosol and cloud radiative 
 !!                              properties. Remove the NSPOT option.
 !!      B.Aouizerats  07/11     add aerosol optical properties CAOP
+!       V. Masson 03/01/2024: aggregation of columns for radiation
 !-------------------------------------------------------------------------------
 !
 !*       0.   DECLARATIONS
diff --git a/src/MNH/modd_radiationsn.f90 b/src/MNH/modd_radiationsn.f90
index 64be9ce63..4cac01d2d 100644
--- a/src/MNH/modd_radiationsn.f90
+++ b/src/MNH/modd_radiationsn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1995-2021 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1995-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -38,7 +38,7 @@
 !!                                         multiple wavelengths for surface SW
 !!  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
 !  P. Wautelet 08/02/2019: add missing NULL association for pointers
-!!
+!       V. Masson 03/01/2024: aggregation of columns for radiation
 !-------------------------------------------------------------------------------
 !
 !*       0.   DECLARATIONS
@@ -132,6 +132,7 @@ TYPE RADIATIONS_t
   INTEGER :: NJ_RAD_AGG    ! reformatted Y array size
   INTEGER :: NIOR_RAD_AGG  ! index of first point of packed array according to current domain
   INTEGER :: NJOR_RAD_AGG  ! index of first point of packed array according to current domain
+  INTEGER, DIMENSION(:,:),     POINTER :: NRAD_AGG_FLAG=>NULL() ! status on which processor is calculated the aggregated column
 !
 END TYPE RADIATIONS_t
 
@@ -196,6 +197,7 @@ INTEGER, POINTER :: NI_RAD_AGG=>NULL()    ! reformatted X array size
 INTEGER, POINTER :: NJ_RAD_AGG=>NULL()    ! reformatted Y array size
 INTEGER, POINTER :: NIOR_RAD_AGG=>NULL()  ! index of first point of packed array according to current domain
 INTEGER, POINTER :: NJOR_RAD_AGG=>NULL()  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:),POINTER :: NRAD_AGG_FLAG=>NULL()
 
 CONTAINS
 
@@ -236,6 +238,7 @@ RADIATIONS_MODEL(KFROM)%XLWD=>XLWD
 RADIATIONS_MODEL(KFROM)%XDTHRADSW=>XDTHRADSW
 RADIATIONS_MODEL(KFROM)%XDTHRADLW=>XDTHRADLW
 RADIATIONS_MODEL(KFROM)%XRADEFF=>XRADEFF
+RADIATIONS_MODEL(KFROM)%NRAD_AGG_FLAG=>NRAD_AGG_FLAG
 !
 ! Current model is set to model KTO
 NDLON=>RADIATIONS_MODEL(KTO)%NDLON
@@ -297,6 +300,7 @@ NI_RAD_AGG=>RADIATIONS_MODEL(KTO)%NI_RAD_AGG
 NJ_RAD_AGG=>RADIATIONS_MODEL(KTO)%NJ_RAD_AGG
 NIOR_RAD_AGG=>RADIATIONS_MODEL(KTO)%NIOR_RAD_AGG
 NJOR_RAD_AGG=>RADIATIONS_MODEL(KTO)%NJOR_RAD_AGG
+NRAD_AGG_FLAG=>RADIATIONS_MODEL(KTO)%NRAD_AGG_FLAG
 
 END SUBROUTINE RADIATIONS_GOTO_MODEL
 
diff --git a/src/MNH/modn_param_radn.f90 b/src/MNH/modn_param_radn.f90
index 477e14df9..ddcc8c11a 100644
--- a/src/MNH/modn_param_radn.f90
+++ b/src/MNH/modn_param_radn.f90
@@ -1,13 +1,8 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
 !MNH_LIC for details. version 1.
 !-----------------------------------------------------------------
-!--------------- special set of characters for RCS information
-!-----------------------------------------------------------------
-! $Source$ $Revision$
-! MASDEV4_7 modn 2006/11/23 17:22:54
-!-----------------------------------------------------------------
 !     ########################  
       MODULE MODN_PARAM_RAD_n
 !     ########################
diff --git a/src/MNH/phys_paramn.f90 b/src/MNH/phys_paramn.f90
index 1c237e478..edd731e42 100644
--- a/src/MNH/phys_paramn.f90
+++ b/src/MNH/phys_paramn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1995-2023 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1995-2024 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -242,6 +242,7 @@ END MODULE MODI_PHYS_PARAM_n
 !  A. Costes      12/2021: add Blaze fire model
 !  Q. Rodier      2022   : integration with PHYEX
 !  C. Barthe      03/2023: add CELEC in call to turbulence
+!  V. Masson      01/2024: aggregation of columns for radiation
 !!-------------------------------------------------------------------------------
 !
 !*       0.     DECLARATIONS
@@ -326,6 +327,7 @@ USE MODD_TIME_n
 USE MODD_TIME, ONLY : TDTEXP  ! Ajout PP
 USE MODD_TURB_FLUX_AIRCRAFT_BALLOON, ONLY : XTHW_FLUX, XRCW_FLUX, XSVW_FLUX
 USE MODD_TURB_n
+USE MODD_VAR_ll,      ONLY: IP
 
 USE MODE_AERO_PSD
 use mode_budget,            only: Budget_store_end, Budget_store_init
@@ -777,7 +779,7 @@ CALL SUNPOS_n   ( XZENITH, ZCOSZEN, ZSINZEN, ZAZIMSOL )
       XLWD(:,:,:)=0.0
       XDTHRADSW(:,:,:)=0.0
       XDTHRADLW(:,:,:)=0.0
-      CALL RADIATIONS_AGG(NRAD_AGG,NI_RAD_AGG,NJ_RAD_AGG,NIOR_RAD_AGG,NJOR_RAD_AGG, TPFILE,      &
+      CALL RADIATIONS_AGG(NRAD_AGG,NI_RAD_AGG,NJ_RAD_AGG,NIOR_RAD_AGG,NJOR_RAD_AGG, NRAD_AGG_FLAG, TPFILE,      &
                        LCLEAR_SKY, OCLOUD_ONLY, NCLEARCOL_TM1, CEFRADL, CEFRADI, COPWSW, COPISW, &
                        COPWLW, COPILW, XFUDG,                                                    &
                        NDLON, NFLEV, NRAD_DIAG, NFLUX, NRAD, NAER, NSWB_OLD, NSWB_MNH, NLWB_MNH, &
diff --git a/src/MNH/radiations.f90 b/src/MNH/radiations.f90
index 7917ae111..246c30f19 100644
--- a/src/MNH/radiations.f90
+++ b/src/MNH/radiations.f90
@@ -123,6 +123,7 @@ CONTAINS
 !  P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg
 !  P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
 !  P. Wautelet 06/09/2022: small fix: GSURF_CLOUD was not set outside of physical domain
+!       V. Masson 03/01/2024: aggregation of columns for radiation
 !-------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
diff --git a/src/MNH/radiations_agg.f90 b/src/MNH/radiations_agg.f90
index bbc5d383f..f06055180 100644
--- a/src/MNH/radiations_agg.f90
+++ b/src/MNH/radiations_agg.f90
@@ -9,7 +9,7 @@
 !
 INTERFACE
 
-    SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG, &
+    SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,KRAD_AGG_FLAG, &
                TPFILE,OCLEAR_SKY,OCLOUD_ONLY,                                  &
                KCLEARCOL_TM1,HEFRADL,HEFRADI,HOPWSW,HOPISW,HOPWLW,HOPILW,      &
                PFUDG, KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD,    &
@@ -26,6 +26,8 @@ INTEGER, INTENT(IN) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(IN) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(IN) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(IN) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(IN) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
+
 
 TYPE(TFILEDATA),  INTENT(IN)         :: TPFILE    ! Output file
 LOGICAL, INTENT(IN)                  :: OCLOUD_ONLY! flag for the cloud column
@@ -113,7 +115,7 @@ END INTERFACE
 END MODULE MODI_RADIATIONS_AGG
 !
 !   ############################################################################
-    SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG, &
+    SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG,KRAD_AGG_FLAG,  &
                TPFILE,OCLEAR_SKY,OCLOUD_ONLY,                                  &
                KCLEARCOL_TM1,HEFRADL,HEFRADI,HOPWSW,HOPISW,HOPWLW,HOPILW,      &
                PFUDG, KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD,    &
@@ -192,6 +194,7 @@ INTEGER, INTENT(IN) :: KI_RAD_AGG    ! reformatted X array size
 INTEGER, INTENT(IN) :: KJ_RAD_AGG    ! reformatted Y array size
 INTEGER, INTENT(IN) :: KIOR_RAD_AGG  ! index of first point of packed array according to current domain
 INTEGER, INTENT(IN) :: KJOR_RAD_AGG  ! index of first point of packed array according to current domain
+INTEGER, DIMENSION(:,:), INTENT(IN) :: KRAD_AGG_FLAG  ! flag to know if aggregated column is computed in this processor or another one
 
 TYPE(TFILEDATA),  INTENT(IN)         :: TPFILE    ! Output file
 LOGICAL, INTENT(IN)                  :: OCLOUD_ONLY! flag for the cloud column
@@ -630,35 +633,87 @@ REAL, DIMENSION(:,:), INTENT(IN)  :: PFULL
 REAL, DIMENSION(:,:), INTENT(OUT) :: PPACK
 LOGICAL                           :: OEXCH
 
-!  DO JJP=1,MIN(KJ_RAD_AGG,IJMAX)
-!    DO JIP=1,MIN(KI_RAD_AGG,IIMAX)
   DO JJP=1,KJ_RAD_AGG
     DO JIP=1,KI_RAD_AGG
       IXORP = KIOR_RAD_AGG + (JIP-1) * KRAD_AGG
       IYORP = KJOR_RAD_AGG + (JJP-1) * KRAD_AGG
-      PPACK(JIP,JJP) = PFULL(MIN(IXORP + KRAD_AGG/2,IImax),MIN(IYORP + KRAD_AGG/2,Ijmax) ) 
+      PPACK(JIP,JJP) = PFULL(MIN(IXORP + KRAD_AGG/2,IIMAX),MIN(IYORP + KRAD_AGG/2,IJMAX) )
     END DO
   END DO
 END SUBROUTINE PACK_RAD_AGG_MID
-!
 !-------------------------------------------------------------------------------
 !
 SUBROUTINE UNPACK_RAD_AGG_2D(PFULL,PPACK)
 REAL, DIMENSION(:,:), INTENT(OUT)  :: PFULL
 REAL, DIMENSION(:,:), INTENT(IN)   :: PPACK
 
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll   ! list of fields to exchange
+REAL, DIMENSION(SIZE(PFULL,1),SIZE(PFULL,2)) :: ZFULL
+
+  ZFULL = 0.
+  ! unpacks columns whose center is located in the current processor
   DO JJP=1,KJ_RAD_AGG
     DO JIP=1,KI_RAD_AGG
       IXORP = KIOR_RAD_AGG + (JIP-1) * KRAD_AGG
       IYORP = KJOR_RAD_AGG + (JJP-1) * KRAD_AGG
       DO JJ=IYORP,MIN(IYORP+KRAD_AGG-1,IJMAX)
         DO JI=IXORP,MIN(IXORP+KRAD_AGG-1,IIMAX)
-          PFULL(JI,JJ) = PPACK(JIP,JJP)
+          ZFULL(JI,JJ) = PPACK(JIP,JJP)
         END DO
       END DO
     END DO
   END DO
 
+  ! stores the output field
+  PFULL = 0.
+  PFULL(IIB:IIE,IJB:IJE) = ZFULL(IIB:IIE,IJB:IJE)
+
+  ! updates boundaries outside processor (that's because some aggregated columns were computed another processors)
+  NULLIFY(TZFIELDS_ll)
+  CALL ADD2DFIELD_ll( TZFIELDS_ll, ZFULL, 'RADIATION_AGG: UNPACK_RAD_AGG_2D' )
+  CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+  CALL CLEANLIST_ll(TZFIELDS_ll)
+
+  ! Get the values of the packed column for points in the HALO if those were indeed computed by the other processors
+!  DO JJ=1,IJB-1
+!    DO JI=1,IIB-1
+!      IF (ZFULL(JI,JJ) /= 0.) PFULL(JI,JJ) = ZFULL(JI,JJ)
+!    END DO
+!  END DO
+!  DO JJ=IJE+1,IJU
+!    DO JI=1,IIB-1
+!      IF (ZFULL(JI,JJ) /= 0.) PFULL(JI,JJ) = ZFULL(JI,JJ)
+!    END DO
+!  END DO
+!  DO JJ=1,IJB-1
+!    DO JI=IIE+1,IIU
+!      IF (ZFULL(JI,JJ) /= 0.) PFULL(JI,JJ) = ZFULL(JI,JJ)
+!    END DO
+!  END DO
+!  DO JJ=IJE+1,IJU
+!    DO JI=IIE+1,IIU
+!      IF (ZFULL(JI,JJ) /= 0.) PFULL(JI,JJ) = ZFULL(JI,JJ)
+!    END DO
+!  END DO
+
+
+  ! Get the values of the packed column for points whose column center is outside the processor
+  ! In that case, the column is on one of the sides outside the processor, but within the NHALO bands
+  !
+  DO JJ=IJB,IJE
+    DO JI=IIB,IIE
+      IF (KRAD_AGG_FLAG(JI,JJ)==0) CYCLE  ! points whose aggregated columns were computed within this processor
+      IF (KRAD_AGG_FLAG(JI,JJ)==1) PFULL(JI,JJ) = ZFULL(IIB-1,JJ) ! computed column was towards West side of processor
+      IF (KRAD_AGG_FLAG(JI,JJ)==2) PFULL(JI,JJ) = ZFULL(JI,IJB-1) ! computed column was towards South side of processor
+      IF (KRAD_AGG_FLAG(JI,JJ)==3) PFULL(JI,JJ) = ZFULL(IIE+1,JJ) ! computed column was towards East side of processor
+      IF (KRAD_AGG_FLAG(JI,JJ)==4) PFULL(JI,JJ) = ZFULL(JI,IJE+1) ! computed column was towards North side of processor
+      IF (KRAD_AGG_FLAG(JI,JJ)==5) PFULL(JI,JJ) = ZFULL(IIB-1,IJB-1) ! computed column was towards South-West corner
+      IF (KRAD_AGG_FLAG(JI,JJ)==6) PFULL(JI,JJ) = ZFULL(IIE+1,IJB-1) ! computed column was towards South-East corner
+      IF (KRAD_AGG_FLAG(JI,JJ)==7) PFULL(JI,JJ) = ZFULL(IIE+1,IJE+1) ! computed column was towards North-East corner
+      IF (KRAD_AGG_FLAG(JI,JJ)==8) PFULL(JI,JJ) = ZFULL(IIB-1,IJE+1) ! computed column was towards North-West corner
+    END DO
+  END DO
+  !
 END SUBROUTINE UNPACK_RAD_AGG_2D
 !-------------------------------------------------------------------------------
 !
diff --git a/src/MNH/read_exsegn.f90 b/src/MNH/read_exsegn.f90
index 8cc0cc0c4..41a70bbf1 100644
--- a/src/MNH/read_exsegn.f90
+++ b/src/MNH/read_exsegn.f90
@@ -312,6 +312,7 @@ END MODULE MODI_READ_EXSEG_n
 !  P. Wautelet 19/08/2022: add namelist for aircrafts
 !  H. Toumi       09/2022: add EOL/ADR
 !  C. Barthe   11/07/2023: ELEC: only some combinations of microphysical options are allowed
+!  V. Masson   03/01/2024: aggregation of columns for radiation
 !------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
-- 
GitLab