diff --git a/src/MNH/default_desfmn.f90 b/src/MNH/default_desfmn.f90 index 47d3c118d920b33fd47f69cb36c6e608175299e8..07359fc5a4638350daa6ee14fac793f8521c511a 100644 --- a/src/MNH/default_desfmn.f90 +++ b/src/MNH/default_desfmn.f90 @@ -209,6 +209,7 @@ END MODULE MODI_DEFAULT_DESFM_n !! 11/2019 C.Lac correction in the drag formula and application to building in addition to tree ! P. Wautelet 17/04/2020: move budgets switch values into modd_budget ! P. Wautelet 30/06/2020: add NNETURSV, NNEADVSV and NNECONSV variables +! PA. Joulin 10/2020: add default variables for wind tubrines (EOL) ! F. Auguste, T. Nagel 02/2021: add IBM defaults parameters ! T. Nagel 02/2021: add turbulence recycling defaults parameters ! P-A Joulin 21/05/2021: add Wind turbines @@ -224,6 +225,7 @@ END MODULE MODI_DEFAULT_DESFM_n ! C. Barthe 03/2022: add CIBU and RDSF options in LIMA ! Delbeke/Vie 03/2022: KHKO option in LIMA ! P. Wautelet 27/04/2022: add namelist for profilers +! PA. Joulin 04/2023: add EOL/ADR !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -274,9 +276,9 @@ USE MODD_MEAN_FIELD USE MODD_DRAGTREE_n USE MODD_DRAGBLDG_n USE MODD_COUPLING_LEVELS_n -USE MODD_EOL_MAIN -USE MODD_EOL_ADNR -USE MODD_EOL_ALM +USE MODD_EOL_MAIN, ONLY: LMAIN_EOL,CMETH_EOL,CSMEAR,NMODEL_EOL +USE MODD_EOL_ADR, ONLY: NNB_RADELT, NNB_AZIELT +USE MODD_EOL_ALM, ONLY: NNB_BLAELT, LTIMESPLIT USE MODD_EOL_SHARED_IO USE MODD_ALLPROFILER_n USE MODD_ALLSTATION_n @@ -567,14 +569,20 @@ CBLADE_CSVDATA = 'data_blade.csv' CAIRFOIL_CSVDATA = 'data_airfoil.csv' ! CINTERP = 'CLS' +LTIPLOSSG = .TRUE. +LTECOUTPTS = .FALSE. +LCSVOUTFRM = .FALSE. +! +! 10d.iii) MODD_EOL_ADR +! +NNB_RADELT = 18 +NNB_AZIELT = 56 ! -! 10d.iii) MODD_EOL_ALM +! 10d.iv) MODD_EOL_ALM ! NNB_BLAELT = 42 LTIMESPLIT = .FALSE. -LTIPLOSSG = .TRUE. -LTECOUTPTS = .FALSE. -! + !------------------------------------------------------------------------------ !* 10.e SET DEFAULT VALUES FOR MODD_ALLPROFILER_n : ! ---------------------------------- diff --git a/src/MNH/eol_adr.f90 b/src/MNH/eol_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..7aa786a3761e8975f33c5baeb68c48338b6e420e --- /dev/null +++ b/src/MNH/eol_adr.f90 @@ -0,0 +1,474 @@ +!MNH_LIC Copyright 2017-2021 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. +!----------------------------------------------------------------- +! ####################### + MODULE MODI_EOL_ADR +! ####################### +! +INTERFACE +! +SUBROUTINE EOL_ADR(KTCOUNT, PTSTEP, & + PDXX, PDYY, PDZZ, & + PRHO_M, & + PUT_M, PVT_M, PWT_M, & + PFX_RG, PFY_RG, PFZ_RG ) + +! +INTEGER, INTENT(IN) :: KTCOUNT ! iteration count +REAL, INTENT(IN) :: PTSTEP ! timestep except +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHO_M ! dry Density +REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT_M,PVT_M,PWT_M ! Wind speed at mass point +! +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFX_RG ! Aerodynamic force .. +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFY_RG ! .. cartesian mesh .. +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFZ_RG ! .. global frame) +! +! +END SUBROUTINE EOL_ADR +! +END INTERFACE +! +END MODULE MODI_EOL_ADR +! +! ################################################################### + SUBROUTINE EOL_ADR(KTCOUNT, PTSTEP, & + PDXX, PDYY, PDZZ, & + PRHO_M, & + PUT_M, PVT_M, PWT_M, & + PFX_RG, PFY_RG, PFZ_RG ) +! ################################################################### +! +!!**** *MODI_EOL_ADR* - +!! +!!** PURPOSE +!! ------- +!! It is possible to include wind turbines parameterization in Meso-NH, +!! and several methods are available. One of the models is the Actuator +!! Disc with Rotation (ADR). It allows to compute aerodynamic forces +!! according the wind speed and the caracteristics of the wind turbine. +!! +!!** METHOD +!! ------ +!! The ADR consists in modeling the rotor by a disc, drawn by the blades +!! (Mikkelsen 2003). The disc is discretized in a cyldrical frame. Each +!! element apply an aerodynamic force into the flow. Each element carries +!! a two-dimensional (2D) airfoil, and its characteristics, as lift and +!! drag coefficients. Knowing these coefficients, and the angle of attack, +!! the lift and drag forces can be evaluated. +!! +!!** REFERENCE +!! --------- +!! PA. Joulin PhD Thesis. 2020. +!! H. Toumi Master Thesis. 2022. +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +!! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!!--------------------------------------------------------------- +! +! +!* 0. DECLARATIONS +! ------------ +! +! To work with wind turbines +USE MODD_EOL_ADR +USE MODD_EOL_KINE_ADR +! +USE MODD_EOL_SHARED_IO, ONLY: CINTERP, LTECOUTPTS, LTIPLOSSG +USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT +USE MODD_EOL_SHARED_IO, ONLY: XELT_RAD +USE MODD_EOL_SHARED_IO, ONLY: XAOA_GLB, XFLIFT_GLB, XFDRAG_GLB, XFAERO_RG_GLB +! +USE MODI_EOL_MATHS +USE MODI_EOL_READER, ONLY: GET_AIRFOIL_ID_ADR +USE MODI_EOL_PRINTER, ONLY: PRINT_TSPLIT +USE MODI_EOL_DEBUGGER +! Math +USE MODD_CST, ONLY: XPI +! To know the grid +USE MODD_GRID_n, ONLY: XXHAT,XYHAT,XZS,XZZ +USE MODE_ll, ONLY: GET_INDICE_ll +USE MODD_PARAMETERS, ONLY: JPVEXT +! MPI stuffs +USE MODD_VAR_ll, ONLY: NMNH_COMM_WORLD +USE MODD_PRECISION, ONLY: MNHREAL_MPI +USE MODD_MPIF, ONLY: MPI_SUM +USE MODE_SUM_ll, ONLY: MIN_ll +USE MODD_VAR_ll, ONLY: IP +! +! +IMPLICIT NONE +! +!* 0.1 Declarations of dummy arguments : +! +INTEGER, INTENT(IN) :: KTCOUNT ! iteration count +REAL, INTENT(IN) :: PTSTEP ! timestep except +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHO_M ! dry Density +REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT_M,PVT_M,PWT_M ! Wind speed at mass point +! +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFX_RG ! Aerodynamic force .. +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFY_RG ! .. cartesian mesh .. +REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PFZ_RG ! .. global frame) +! +!* 0.2 Declarations of local variables : +! +! Indicies Compteurs +INTEGER :: IIB,IJB,IKB ! Begin of a CPU domain +INTEGER :: IIE,IJE,IKE ! End of a CPU domain +INTEGER :: IKU ! Vertical size of the domain +INTEGER :: JI, JJ, JK ! Loop index +INTEGER :: JROT, JAZI, JRAD ! Rotor, azimutal, and radial element indicies +! +! Aggregated variables +REAL, DIMENSION(TFARM%NNB_TURBINES,TBLADE%NNB_AZIELT,TBLADE%NNB_RADELT) :: ZAOA_ELT ! Angle of attack of an element, hub frame [rad] +REAL, DIMENSION(TFARM%NNB_TURBINES,TBLADE%NNB_AZIELT,TBLADE%NNB_RADELT) :: ZFLIFT_ELT ! Aerodynamic lift force, parallel to Urel [N] +REAL, DIMENSION(TFARM%NNB_TURBINES,TBLADE%NNB_AZIELT,TBLADE%NNB_RADELT) :: ZFDRAG_ELT ! Aerodynamic drag force, perpendicular to Urel [N] +REAL, DIMENSION(TFARM%NNB_TURBINES,TBLADE%NNB_AZIELT,TBLADE%NNB_RADELT,3) :: ZFAERO_RA_ELT ! Aerodynamic force (lift+drag) in RA [N] +REAL, DIMENSION(TFARM%NNB_TURBINES,TBLADE%NNB_AZIELT,TBLADE%NNB_RADELT,3) :: ZFAERO_RG_ELT ! Aerodynamic force (lift+drag) in RG [N] + +! +! -- Wind -- +REAL :: ZRHO_I ! Interpolated density [kg/m3] +REAL :: ZUT_I ! Interpolated wind speed U (RG) [m/s] +REAL :: ZVT_I ! Interpolated wind speed V (RG) [m/s] +REAL :: ZWT_I ! Interpolated wind speed W (RG) [m/s] +REAL, DIMENSION(3) :: ZWIND_VEL_RG ! Wind velocity in RG frame [m/s] +REAL, DIMENSION(3) :: ZWIND_VEL_RA ! Wind velocity in RE frame [m/s] +REAL, DIMENSION(3) :: ZWINDREL_VEL_RA ! Relative wind velocity in RE frame [m/s] +REAL :: ZWINDREL_VEL ! Norm of the relative wind velocity [m/s] +REAL, DIMENSION(SIZE(PUT_M,1),SIZE(PUT_M,2),SIZE(PUT_M,3)) :: ZZH ! True heigth to interpolate 8NB +! +! -- Wind turbine -- +INTEGER :: INB_WT, INB_B ! Total numbers of turbines and blades +INTEGER :: INB_AELT, INB_RELT ! Total numbers of azimutal and radial elt +REAL :: ZRAD ! Blade radius [m] +INTEGER :: IAID ! Airfoil index [-] +! +! -- Aero -- +REAL :: ZAOA ! Attack angle of an element [rad] +REAL :: ZCDRAG ! Drag coefficient of an element [] +REAL :: ZCLIFT ! Lift coefficient of an element [] +REAL :: ZFDRAG ! Drag force of an element, parallel to Urel [N] +REAL :: ZFLIFT ! Lift force of an element, perpendicular to Urel [N] +REAL, DIMENSION(3) :: ZFAERO_RA ! Aerodynamic force (lift+drag) in RA [N] +REAL, DIMENSION(3) :: ZFAERO_RG ! Aerodynamic force (lift+drag) in RG [N] +! Tip loss +REAL :: ZFTIPL ! tip loss function +REAL :: ZPHI ! angle twist+pitch+aa +! +! Thrust, Torque and Power +REAL, DIMENSION(3) :: ZFAERO_RH ! Aerodynamic force (lift+drag) in RH [N] (thrust/torque) +REAL, DIMENSION(3) :: ZDIST_HBELT_RH ! Distance between blade element and hub, in RH [m] +REAL, DIMENSION(3) :: ZDIST_HBELT_RG ! Distances between blade element and hub, in RG [m] +REAL, DIMENSION(3) :: Z3D_TORQT ! Full torque force (3D) of the wind turbine [N] +! +! -- Numerical -- +INTEGER :: IINFO ! code info return +! +! +!* 0.3 Implicit arguments +! +! A. From MODD_EOL_ADR +!TYPE(FARM) :: TFARM +!TYPE(TURBINE) :: TTURBINE +!TYPE(BLADE) :: TBLADE +!TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE :: TAIRFOIL +! +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_AZI ! Elements azimut [rad] +!REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RA_GLB ! Aero. force (lift+drag) RA [N] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFAERO_BLEQ_RA_GLB ! Blade Eq Aero. force (lift+drag) in RA [N] +!REAL, DIMENSION(:,:), ALLOCATABLE :: XAOA_BLEQ_GLB ! Blade Eq. AoA + +! B. From MODD_EOL_SHARED_IO: +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Blade elements radius [m] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] +!REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RE_GLB ! Aerodyn. force (lift+drag) in RE [N] +!REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] +! +! for namelist NAM_EOL_ADR +!CHARACTER(LEN=3) :: CINTERP ! Interpolation method for wind speed +!LOGICAL :: LTIPLOSSG ! Flag to apply Glauert's tip loss correction +!LOGICAL :: LTECOUTPTS ! Flag to get Tecplot file output of element points +!LOGICAL :: LCSVOUTFRM ! Flag to get CSV files output of frames +! +! for output +!REAL, DIMENSION(:), ALLOCATABLE :: XTHRUT ! Thrust [N] +!REAL, DIMENSION(:), ALLOCATABLE :: XTORQT ! Torque [Nm] +!REAL, DIMENSION(:), ALLOCATABLE :: XPOWT ! Power [W] +! +! +!------------------------------------------------------------------------------- +! +! +!* 1. INITIALIZATIONS +! --------------- +! +!* 1.1 Subdomain (CPU) indices +! +CALL GET_INDICE_ll(IIB,IJB,IIE,IJE) ! Get begin and end domain index (CPU) +IKU = SIZE(PUT_M,3) ! Top of the domain end index +IKB=1+JPVEXT ! Vertical begin index +IKE=IKU-JPVEXT ! Vertical end index +! +!* 1.2 Some usefull integers +! +INB_WT = TFARM%NNB_TURBINES +INB_B = TTURBINE%NNB_BLADES +INB_AELT = TBLADE%NNB_AZIELT +INB_RELT = TBLADE%NNB_RADELT +! +!* 1.3 Vertical coordinate in case of interpolation +! +IF (CINTERP=='8NB') THEN + DO JK=1,IKU-1 + ZZH(:,:,JK) = (0.5*(XZZ(:,:,JK)+XZZ(:,:,JK+1))-XZS(:,:)) + END DO + ZZH(:,:,IKU) = 2*ZZH(:,:,IKU-1) - ZZH(:,:,IKU-2) +END IF +! +!* 1.4 Set to zeros at each MNH time steps +! +ZAOA_ELT(:,:,:) = 0. +ZFLIFT_ELT(:,:,:) = 0. +ZFDRAG_ELT(:,:,:) = 0. +ZFAERO_RA_ELT(:,:,:,:) = 0. +ZFAERO_RG_ELT(:,:,:,:) = 0. +! Global variables (seen by all CPU) +XAOA_GLB(:,:,:) = 0. +XFLIFT_GLB(:,:,:) = 0. +XFDRAG_GLB(:,:,:) = 0. +XFAERO_RA_GLB(:,:,:,:) = 0. +XFAERO_RG_GLB(:,:,:,:) = 0. +! +XFAERO_BLEQ_RA_GLB(:,:,:) = 0. +XAOA_BLEQ_GLB(:,:) = 0. +! +XTHRUT(:) = 0. +XTORQT(:) = 0. +! +!------------------------------------------------------------------------------ +! +!* 2. KINEMATICS COMPUTATIONS +! ----------------------- +! + CALL EOL_KINE_ADR(KTCOUNT, PTSTEP) +! +! +!------------------------------------------------------------------------------- +! +!* 4. COMPUTES AERODYNAMIC FORCES THAT ACTS ON THE BLADES DUE TO THE WIND +! -------------------------------------------------------------- +! +!* 4.1 Finding the position of wind turbines +! +! Loop over domain + DO JK=IKB,IKE + DO JJ=IJB,IJE + DO JI=IIB,IIE + ! Loop over wind turbines + DO JROT=1, INB_WT + DO JAZI=1, INB_AELT + DO JRAD=1, INB_RELT + + ! Position test + IF (XPOS_ELT_RG(JROT,JAZI,JRAD,1) >= XXHAT(JI) .AND. & + XPOS_ELT_RG(JROT,JAZI,JRAD,1) < XXHAT(JI) + PDXX(JI,JJ,JK)) THEN +! + IF (XPOS_ELT_RG(JROT,JAZI,JRAD,2) >= XYHAT(JJ) .AND. & + XPOS_ELT_RG(JROT,JAZI,JRAD,2) < XYHAT(JJ) + PDYY(JI,JJ,JK)) THEN +! + IF (XPOS_ELT_RG(JROT,JAZI,JRAD,3) >= XZZ(JI,JJ,JK) .AND. & + XPOS_ELT_RG(JROT,JAZI,JRAD,3) < XZZ(JI,JJ,JK) + PDZZ(JI,JJ,JK)) THEN +! +!* 4.2 Extracting the wind +! + SELECT CASE(CINTERP) + CASE('CLS') + ZUT_I = PUT_M(JI,JJ,JK) + ZVT_I = PVT_M(JI,JJ,JK) + ZWT_I = PWT_M(JI,JJ,JK) + ZRHO_I = PRHO_M(JI,JJ,JK) + CASE('8NB') + ZUT_I = INTERP_LIN8NB(XPOS_ELT_RG(JROT,JAZI,JRAD,:),& + JI,JJ,JK,PUT_M,ZZH) + ZVT_I = INTERP_LIN8NB(XPOS_ELT_RG(JROT,JAZI,JRAD,:),& + JI,JJ,JK,PVT_M,ZZH) + ZWT_I = INTERP_LIN8NB(XPOS_ELT_RG(JROT,JAZI,JRAD,:),& + JI,JJ,JK,PWT_M,ZZH) + ZRHO_I = INTERP_LIN8NB(XPOS_ELT_RG(JROT,JAZI,JRAD,:),& + JI,JJ,JK,PRHO_M,ZZH) + END SELECT + ZWIND_VEL_RG(1) = ZUT_I + ZWIND_VEL_RG(2) = ZVT_I + ZWIND_VEL_RG(3) = ZWT_I +! +!* 4.3 Calculating the wind in RA frame +! + ZWIND_VEL_RA(:) = MATMUL(XMAT_RA_RG(JROT,JAZI,JRAD,:,:), ZWIND_VEL_RG(:)) +! +!* 4.4 Calculating the relative wind speed in RE frame + norm +! + ZWINDREL_VEL_RA(:) = ZWIND_VEL_RA(:) - XTVEL_ELT_RA(JROT,JAZI,JRAD,:) + ZWINDREL_VEL = NORM(ZWINDREL_VEL_RA) +! +!* 4.5 Calculating the angle of attack +! + ZAOA = ATAN2(-ZWINDREL_VEL_RA(3), ZWINDREL_VEL_RA(2)) +! +!* 4.6 Getting aerodynamic coefficients from tabulated data +! + ZRAD = XELT_RAD(JROT,JAZI,JRAD) ! Radius of the element + IAID = GET_AIRFOIL_ID_ADR(TTURBINE,TBLADE,TAIRFOIL,ZRAD) ! ID of the airfoil + ZCLIFT = INTERP_SPLCUB(ZAOA*180/XPI, & + TAIRFOIL(IAID)%XAA,& + TAIRFOIL(IAID)%XCL) + ZCDRAG = INTERP_SPLCUB(ZAOA*180/XPI, & + TAIRFOIL(IAID)%XAA,& + TAIRFOIL(IAID)%XCD) +! +!* 4.7 Tip loss correction (Glauert) +! + IF (LTIPLOSSG) THEN + ZPHI = + ZAOA & + - TFARM%XBLA_PITCH(JROT) & + - XTWIST_ELT(JROT,JRAD) + IF (ZPHI > 0.0) THEN + ZFTIPL = (2.0/XPI)*ACOS(MIN( & + 1.0, EXP(-(TTURBINE%NNB_BLADES/2.0) & + *(TTURBINE%XR_MAX-ZRAD)/(ZRAD*SIN(ZPHI))))) + ELSE + ZFTIPL = 1.0 + END IF + ZCLIFT = ZFTIPL*ZCLIFT + ZCDRAG = ZFTIPL*ZCDRAG + END IF +! +!* 4.8 Computing aerodynamic forces in relative frame +! that act on blades (wind->blade) + ZFLIFT = 0.5*ZRHO_I*XSURF_APP_ELT(JROT,JRAD)*ZCLIFT*ZWINDREL_VEL**2 + ZFDRAG = 0.5*ZRHO_I*XSURF_APP_ELT(JROT,JRAD)*ZCDRAG*ZWINDREL_VEL**2 +! +!* 4.9 Evaluating the aerodynamiques forces in RE frame +! that act on blades (wind->blade) + ZFAERO_RA(1) = .0 + ZFAERO_RA(2) = COS(ZAOA)*ZFDRAG - SIN(ZAOA)*ZFLIFT + ZFAERO_RA(3) = - SIN(ZAOA)*ZFDRAG - COS(ZAOA)*ZFLIFT +! +!* 4.10 Evaluating the aerodynamiques forces in RG frame +! that act on blades (wind->blade) + ZFAERO_RG(:) = MATMUL(XMAT_RG_RA(JROT,JAZI,JRAD,:,:), ZFAERO_RA(:)) +! +!* 4.11 Adding it to the cell of Meso-NH + PFX_RG(JI,JJ,JK) = PFX_RG(JI,JJ,JK) + ZFAERO_RG(1) + PFY_RG(JI,JJ,JK) = PFY_RG(JI,JJ,JK) + ZFAERO_RG(2) + PFZ_RG(JI,JJ,JK) = PFZ_RG(JI,JJ,JK) + ZFAERO_RG(3) +! +!* 4.12 Storing variables (because they were local only) + ZAOA_ELT(JROT,JAZI,JRAD) = ZAOA + ZFLIFT_ELT(JROT,JAZI,JRAD) = ZFLIFT + ZFDRAG_ELT(JROT,JAZI,JRAD) = ZFDRAG + ZFAERO_RA_ELT(JROT,JAZI,JRAD,:) = ZFAERO_RA(:) + ZFAERO_RG_ELT(JROT,JAZI,JRAD,:) = ZFAERO_RG(:) + + ! End of position tests + END IF + END IF + END IF + ! End of wind turbine loops + END DO + END DO + END DO + ! End of domain loops + END DO + END DO + END DO +! +!* 4.13 Top and bottom conditions +PFX_RG(:,:,IKB-1) = PFX_RG(:,:,IKB) +PFX_RG(:,:,IKE+1) = PFX_RG(:,:,IKE) +! +PFY_RG(:,:,IKB-1) = PFY_RG(:,:,IKB) +PFY_RG(:,:,IKE+1) = PFY_RG(:,:,IKE) +! +PFZ_RG(:,:,IKB-1) = PFZ_RG(:,:,IKB) +PFZ_RG(:,:,IKE+1) = PFZ_RG(:,:,IKE) +! +! +!------------------------------------------------------------------------------- +! +!* 5. SHARING THE DATAS OVER THE CPUS +! ------------------------------- +CALL MPI_ALLREDUCE(ZAOA_ELT, XAOA_GLB, SIZE(XAOA_GLB), & + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) +CALL MPI_ALLREDUCE(ZFLIFT_ELT, XFLIFT_GLB, SIZE(XFLIFT_GLB), & + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) +CALL MPI_ALLREDUCE(ZFDRAG_ELT, XFDRAG_GLB, SIZE(XFDRAG_GLB), & + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) +CALL MPI_ALLREDUCE(ZFAERO_RA_ELT, XFAERO_RA_GLB, SIZE(XFAERO_RA_GLB),& + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) +CALL MPI_ALLREDUCE(ZFAERO_RG_ELT, XFAERO_RG_GLB, SIZE(XFAERO_RG_GLB),& + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) +! +! +!------------------------------------------------------------------------------- +! +!* 6. COMPUTING THRUST, TORQUE AND POWER +! --------------------------------- +! +IF(IP == 1) THEN + DO JROT=1,TFARM%NNB_TURBINES + + DO JRAD=1, TBLADE%NNB_RADELT + DO JAZI=1, TBLADE%NNB_AZIELT +! +!* 6.1 Preliminaries +! Mean AOA on one blade + XAOA_BLEQ_GLB(JROT,JRAD) = XAOA_BLEQ_GLB(JROT,JRAD) + XAOA_GLB(JROT,JAZI,JRAD) +! Aerodynamic load (wind->blade) in RA + XFAERO_BLEQ_RA_GLB(JROT,JRAD,:) = XFAERO_BLEQ_RA_GLB(JROT,JRAD,:) & + + XFAERO_RA_GLB(JROT,JAZI,JRAD,:)/INB_B + END DO + END DO + + + DO JAZI=1, TBLADE%NNB_AZIELT + DO JRAD=1, TBLADE%NNB_RADELT +! Aerodynamic load (wind->blade) in RH + ZFAERO_RH(:) = MATMUL(XMAT_RH_RG(JROT,:,:), & + XFAERO_RG_GLB(JROT,JAZI,JRAD,:)) +! Distance between element and hub in RG + ZDIST_HBELT_RG(:) = XPOS_ELT_RG(JROT,JAZI,JRAD,:) - XPOS_HUB_RG(JROT,:) +! Distance between element and hub in RH + ZDIST_HBELT_RH(:) = MATMUL(XMAT_RH_RG(JROT,:,:),ZDIST_HBELT_RG(:)) +! +!* 6.2 Thrust (wind->rotor): in RH + XTHRUT(JROT) = XTHRUT(JROT) + ZFAERO_RH(3) ! Only Z component +!* 6.3 Torque (wind->rotor) in RH + Z3D_TORQT = CROSS(ZDIST_HBELT_RH(:),ZFAERO_RH(:)) + XTORQT(JROT) = XTORQT(JROT) + Z3D_TORQT(3) ! Only Z component + END DO + END DO +! +!* 6.4 Power (wind->rotor) + XPOWT(JROT) = XTORQT(JROT) * TFARM%XOMEGA(JROT) + END DO +END IF +! +! +END SUBROUTINE EOL_ADR diff --git a/src/MNH/eol_alm.f90 b/src/MNH/eol_alm.f90 index 509a014b15177f20aba4a338929e8ee2ae4b4ab8..648ae6575115848cd28188bfd0ff58c70f1f0855 100644 --- a/src/MNH/eol_alm.f90 +++ b/src/MNH/eol_alm.f90 @@ -74,7 +74,8 @@ END MODULE MODI_EOL_ALM !! ------------- !! Original 24/01/17 !! Modification 14/10/20 (PA. Joulin) Updated for a main version -! P. Wautelet 23/07/2021: replace non-standard FLOAT function by REAL function +!! P. Wautelet 23/07/2021: replace non-standard FLOAT function by REAL function +!! H. Toumi 09/2022: uses specific functions for ALM !!--------------------------------------------------------------- ! ! @@ -85,11 +86,13 @@ END MODULE MODI_EOL_ALM USE MODD_EOL_ALM USE MODD_EOL_KINE_ALM ! -USE MODD_EOL_SHARED_IO, ONLY: CINTERP +USE MODD_EOL_SHARED_IO, ONLY: CINTERP, LTECOUTPTS,LTIPLOSSG USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT +USE MODD_EOL_SHARED_IO, ONLY: XELT_RAD +USE MODD_EOL_SHARED_IO, ONLY: XAOA_GLB, XFLIFT_GLB, XFDRAG_GLB, XFAERO_RG_GLB ! USE MODI_EOL_MATHS -USE MODI_EOL_READER, ONLY: GET_AIRFOIL_ID +USE MODI_EOL_READER, ONLY: GET_AIRFOIL_ID_ALM USE MODI_EOL_PRINTER, ONLY: PRINT_TSPLIT USE MODI_EOL_ERROR, ONLY: EOL_WTCFL_ERROR ! Math @@ -367,7 +370,7 @@ DO KTSUBCOUNT=1,INBSUBCOUNT !* 4.6 Getting aerodynamic coefficients from tabulated data ! ZRAD = XELT_RAD(JROT,JBLA,JBELT) ! Radius of the element - IAID = GET_AIRFOIL_ID(TTURBINE,TBLADE,TAIRFOIL,ZRAD) ! ID of the airfoil + IAID = GET_AIRFOIL_ID_ALM(TTURBINE,TBLADE,TAIRFOIL,ZRAD) ! ID of the airfoil ZCLIFT = INTERP_SPLCUB(ZAOA*180/XPI, & TAIRFOIL(IAID)%XAA,& TAIRFOIL(IAID)%XCL) diff --git a/src/MNH/eol_error.f90 b/src/MNH/eol_error.f90 index c6c221f9f4a5ba6818a6a8242084a253b3ba518d..5568718960ebe667e8c0b151cf8ed3fae18c5aec 100644 --- a/src/MNH/eol_error.f90 +++ b/src/MNH/eol_error.f90 @@ -23,22 +23,39 @@ SUBROUTINE EOL_CSVEMPTY_ERROR(HFILE,KNBLINE) END SUBROUTINE EOL_CSVEMPTY_ERROR ! ! -! *** -! ALM -! *** +! ********* +! ALM & ADR +! ********* ! SUBROUTINE EOL_AIRFOILNOTFOUND_ERROR(HFILE,HVAR) CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file read CHARACTER(LEN=*), INTENT(IN) :: HVAR ! missing data END SUBROUTINE EOL_AIRFOILNOTFOUND_ERROR ! +SUBROUTINE EOL_BLADEDATA_ERROR(PDELTARAD) + REAL, INTENT(IN) :: PDELTARAD ! Span lenght of an element +END SUBROUTINE EOL_BLADEDATA_ERROR +! +SUBROUTINE EOL_DR_ERROR(KNB_RELT_MIN) + INTEGER, INTENT(IN) :: KNB_RELT_MIN ! minimum number of rad elt for ADR & ALM +END SUBROUTINE EOL_DR_ERROR +! +! *** +! ALM +! *** +! SUBROUTINE EOL_WTCFL_ERROR(PMAXTSTEP) REAL, INTENT(IN) :: PMAXTSTEP ! maximum acceptable time-step END SUBROUTINE EOL_WTCFL_ERROR ! -SUBROUTINE EOL_BLADEDATA_ERROR(PDELTARAD) - REAL, INTENT(IN) :: PDELTARAD ! Span lenght of an element -END SUBROUTINE EOL_BLADEDATA_ERROR +! *** +! ADR +! *** +! +SUBROUTINE EOL_DA_ERROR(KNB_AELT_MIN) + INTEGER, INTENT(IN) :: KNB_AELT_MIN ! minimum number of rad elt for ADR +END SUBROUTINE EOL_DA_ERROR +! ! END INTERFACE ! @@ -148,7 +165,7 @@ USE MODE_IO_FILE, ONLY: IO_File_close USE MODE_MSG USE MODD_EOL_SHARED_IO, ONLY: CBLADE_CSVDATA ! -REAL, INTENT(IN) :: PDELTARAD ! hals section width +REAL, INTENT(IN) :: PDELTARAD ! aero section width ! CHARACTER(LEN=4) :: YDELTARAD ! @@ -164,3 +181,47 @@ CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'EOL_BLADEDATA_ERROR' ) END SUBROUTINE EOL_BLADEDATA_ERROR !######################################################### ! +!######################################################### +SUBROUTINE EOL_DR_ERROR(KNB_RELT_MIN) +! +USE MODD_LUNIT_n, ONLY: TLUOUT +USE MODE_IO_FILE, ONLY: IO_File_close +USE MODE_MSG +! +INTEGER, INTENT(IN) :: KNB_RELT_MIN ! minimum number of radial element +! +CHARACTER(LEN=4) :: YNB_RELT_MIN +! +WRITE( YNB_RELT_MIN, '( I3 )' ) KNB_RELT_MIN +! +CMNHMSG(1) = 'EOL Initialization error: error while meshing blades.' +CMNHMSG(2) = 'Considering the grid, the number of discretized radial' +CMNHMSG(3) = 'element should be at least : ' // TRIM(YNB_RELT_MIN) +CMNHMSG(4) = 'Please, modify NNB_RADELT (NAM_ADR) or NNB_BLAELT (NAM_ALM)' +CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'EOL_MESH_ERROR' ) +! +END SUBROUTINE EOL_DR_ERROR +! +!######################################################### +! +!######################################################### +SUBROUTINE EOL_DA_ERROR(KNB_AELT_MIN) +! +USE MODD_LUNIT_n, ONLY: TLUOUT +USE MODE_IO_FILE, ONLY: IO_File_close +USE MODE_MSG +! +INTEGER, INTENT(IN) :: KNB_AELT_MIN ! minimum number of radial element +! +CHARACTER(LEN=4) :: YNB_AELT_MIN +! +WRITE( YNB_AELT_MIN, '( I3 )' ) KNB_AELT_MIN +! +CMNHMSG(1) = 'EOL Initialization error: error while meshing blades.' +CMNHMSG(2) = 'Considering the grid, the number of discretized azimutal' +CMNHMSG(3) = 'element should be at least : ' // TRIM(YNB_AELT_MIN) +CMNHMSG(4) = 'Please, modify NNB_AZIELT in NAM_ADR' +CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'EOL_MESH_ERROR' ) +! +END SUBROUTINE EOL_DA_ERROR +! diff --git a/src/MNH/eol_kine_adr.f90 b/src/MNH/eol_kine_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..af5ccb1151bade03fd4ffce7d35d784606a05ec4 --- /dev/null +++ b/src/MNH/eol_kine_adr.f90 @@ -0,0 +1,343 @@ +!MNH_LIC Copyright 2017-2021 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. +!----------------------------------------------------------------- +! +!####################### +MODULE MODI_EOL_KINE_ADR +! +INTERFACE !---------------------------------------- +! +SUBROUTINE EOL_KINE_ADR(KTCOUNT,PTSTEP) +! + INTEGER, INTENT(IN) :: KTCOUNT ! iteration count + REAL, INTENT(IN) :: PTSTEP ! timestep +! +END SUBROUTINE EOL_KINE_ADR +! +END INTERFACE !------------------------------------ +! +END MODULE MODI_EOL_KINE_ADR +!####################### +! +! +! +! +!################################################################### +SUBROUTINE EOL_KINE_ADR(KTCOUNT,PTSTEP) +! +!!**** *EOL_KINEMATICS * - +!! +!! PURPOSE +!! ------- +!! Compute positions, oritentations and velocities of all the +!! elements of the wind turbine +!! +!!** METHOD +!! ------ +!! +!! REFERENCE +!! --------- +!! +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +!! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!!--------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ------------ +! +!* 0.1 Modules +USE MODD_EOL_KINE_ADR +USE MODD_EOL_ADR +USE MODI_EOL_MATHS +USE MODI_EOL_PRINTER +USE MODD_TIME_n, ONLY : TDTCUR +USE MODD_CST, ONLY : XPI +USE MODD_EOL_SHARED_IO, ONLY: LTECOUTPTS,LCSVOUTFRM +! +IMPLICIT NONE +! +!* 0.2 Declarations of dummy arguments : +! +INTEGER, INTENT(IN) :: KTCOUNT ! iteration count +REAL, INTENT(IN) :: PTSTEP ! timestep +! +!* 0.3 Local variables +REAL, DIMENSION(3,3) :: ZORI_TMAT_X, ZORI_TMAT_Y, ZORI_TMAT_Z +REAL, DIMENSION(3,3) :: ZORI_MAT_X, ZORI_MAT_Y, ZORI_MAT_Z +REAL, DIMENSION(3) :: ZADD_TO_POS +! +REAL, DIMENSION(3) :: ZDIST_TOWO_TELT_RG ! Distance between tower elmt and tower base +REAL, DIMENSION(3) :: ZDIST_TOWO_NELT_RG ! Distance between nacelle and base of tower +REAL, DIMENSION(3) :: ZDIST_NAC_HUB_RG ! Distance between hub and base of nacelle +REAL, DIMENSION(3) :: ZDIST_HUB_CYL_RG ! Distance between cylindrical base and base of hub +REAL, DIMENSION(3) :: ZDIST_CYL_ELT_RG ! Distance between cylindrical base and elements +! +REAL :: ZTIME ! TIME +INTEGER :: JROT, JTELT, JNELT, JAZI ,JRAD ! Loop control +INTEGER :: INB_WT,INB_B,INB_AELT, INB_RELT ! Total numbers +INTEGER :: INB_TELT, INB_NELT ! Total numbers +INTEGER :: ITECOUT ! Unit number for Tecplot file +INTEGER :: ICSVOUT ! Unit number for TOW.csv file +! +! +!--------------------------------------------------------------- +! +!* 1. PRELIMINARIES +! ------------- +! +!* 1.1 Some useful integers +INB_WT = TFARM%NNB_TURBINES +INB_B = TTURBINE%NNB_BLADES +INB_TELT = 2 +INB_NELT = 2 +INB_AELT = TBLADE%NNB_AZIELT +INB_RELT = TBLADE%NNB_RADELT + +! +!* 1.2 Sub-time computation +ZTIME = TDTCUR%xtime +! +!* 1.3 Tecplotfile : opening + headers +IF (LTECOUTPTS) THEN + CALL OPEN_TECOUT_ADR(ITECOUT, KTCOUNT) +END IF +IF (LCSVOUTFRM) THEN + CALL OPEN_3DCSV_ADR(ICSVOUT, KTCOUNT) +END IF +! +! +!--------------------------------------------------------------- +! +!* 2. COMPUTATIONS +! ------------ +! +DO JROT=1, INB_WT +! +! ---- TOWER ---- +! +!* T.0 Update origin positions in RG (Base position + floating) + XPOS_TOWO_RG(JROT,:) = XPOS_REF(:) + XPOSINI_TOWO_RG(JROT,:) & + + XTVEL_TOWO_RG(JROT,:)*ZTIME +! +!* T.1 Update orientation + CALL GET_ORI_MAT_X(XANGINI_TOW_RG(JROT,1) + XRVEL_RT_RG(JROT,1)*ZTIME, ZORI_MAT_X) + CALL GET_ORI_MAT_Y(XANGINI_TOW_RG(JROT,2) + XRVEL_RT_RG(JROT,2)*ZTIME, ZORI_MAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_TOW_RG(JROT,3) + XRVEL_RT_RG(JROT,3)*ZTIME, ZORI_MAT_Z) +! Compute orientation matrix + XMAT_RG_RT(JROT,:,:) = MATMUL(ZORI_MAT_X, MATMUL(ZORI_MAT_Y, ZORI_MAT_Z)) +! +!* T.2 Update positions in RG + DO JTELT=1, INB_TELT + XPOS_TELT_RG(JROT,JTELT,:) = XPOS_TOWO_RG(JROT,:) & + + MATMUL(XMAT_RG_RT(JROT,:,:),XPOS_TELT_RT(JROT,JTELT,:)) + END DO +! +!* T.3 Update structural velocities + DO JTELT=1, INB_TELT + ! Rotation of tower already in RG + ! Translation of elements + ZDIST_TOWO_TELT_RG(:) = XPOS_TELT_RG(JROT,JTELT,:) - XPOS_TOWO_RG(JROT,:) + XTVEL_TELT_RG(JROT,JTELT,:) = XTVEL_TOWO_RG(JROT,:) & + + CROSS(XRVEL_RT_RG(JROT,:),ZDIST_TOWO_TELT_RG(:)) + ENDDO +! +!* T.4 Print in tecplot file + IF (LTECOUTPTS) THEN + DO JTELT=1, INB_TELT + CALL PRINT_TECOUT(ITECOUT, XPOS_TELT_RG(JROT,JTELT,:)) + END DO + END IF +! +!* T.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'TOW',XMAT_RG_RT(JROT,:,:),XPOS_TOWO_RG(JROT,:)) + END IF +! +! +! ---- NACELLE ---- +! +!* N.0 Update origin positions in RG + XPOS_NACO_RG(JROT,:) = XPOS_TELT_RG(JROT,INB_TELT,:) & + + MATMUL(XMAT_RG_RT(JROT,:,:),XPOSINI_NACO_RT(JROT,:)) +! +!* N.1 Update orientation + CALL GET_ORI_MAT_X(XANGINI_NAC_RT(JROT,1) + XRVEL_RN_RT(JROT,1)*ZTIME, ZORI_MAT_X) + CALL GET_ORI_MAT_Y(XANGINI_NAC_RT(JROT,2) + XRVEL_RN_RT(JROT,2)*ZTIME, ZORI_MAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_NAC_RT(JROT,3) + XRVEL_RN_RT(JROT,3)*ZTIME, ZORI_MAT_Z) +! +! Orientation matrix + XMAT_RT_RN(JROT,:,:) = MATMUL(ZORI_MAT_X, MATMUL(ZORI_MAT_Y, ZORI_MAT_Z)) + XMAT_RG_RN(JROT,:,:) = MATMUL(XMAT_RG_RT(JROT,:,:), XMAT_RT_RN(JROT,:,:)) +! +!* N.2 Update positions in RG + DO JNELT=1, INB_NELT + XPOS_NELT_RG(JROT,JNELT,:) = XPOS_NACO_RG(JROT,:) & + + MATMUL(XMAT_RG_RN(JROT,:,:),XPOS_NELT_RN(JROT,JNELT,:)) + END DO +! +!* N.3 Update structural velocities + ! Rotation of nacelle in RG + XRVEL_RN_RG(JROT,:) = MATMUL(XMAT_RG_RT(JROT,:,:),XRVEL_RN_RT(JROT,:)) & + + XRVEL_RT_RG(JROT,:) + DO JNELT=1, INB_NELT + ! Translation of elements in RG + ZDIST_TOWO_NELT_RG(:) = XPOS_NELT_RG(JROT,JNELT,:) - XPOS_TOWO_RG(JROT,:) + XTVEL_NELT_RG(JROT,JNELT,:) = XTVEL_TOWO_RG(JROT,:) & + + CROSS(XRVEL_RN_RG(JROT,:),ZDIST_TOWO_NELT_RG(:)) + END DO +! +!* N.4 Print in tecplot file + IF (LTECOUTPTS) THEN + DO JNELT=1, INB_NELT + CALL PRINT_TECOUT(ITECOUT, XPOS_NELT_RG(JROT,JNELT,:)) + END DO + END IF +! +!* N.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'NAC',XMAT_RG_RN(JROT,:,:),XPOS_NACO_RG(JROT,:)) + END IF +! +! ---- HUB ---- +! +!* H.1 Update positions + XPOS_HUB_RG(JROT,:) = XPOS_NELT_RG(JROT,INB_NELT,:) & + + MATMUL(XMAT_RG_RN(JROT,:,:),XPOSINI_HUB_RN(JROT,:)) +! +!* H.2 Update orientation + CALL GET_ORI_MAT_X(XANGINI_HUB_RN(JROT,1) + XRVEL_RH_RN(JROT,1)*ZTIME, ZORI_MAT_X) + CALL GET_ORI_MAT_Y(XANGINI_HUB_RN(JROT,2) + XRVEL_RH_RN(JROT,2)*ZTIME, ZORI_MAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_HUB_RN(JROT,3) + XRVEL_RH_RN(JROT,3)*ZTIME, ZORI_MAT_Z) +! Orientation matrix + XMAT_RN_RH(JROT,:,:) = MATMUL(ZORI_MAT_X, MATMUL(ZORI_MAT_Y, ZORI_MAT_Z)) + XMAT_RG_RH(JROT,:,:) = MATMUL(XMAT_RG_RN(JROT,:,:), XMAT_RN_RH(JROT,:,:)) + XMAT_RH_RG(JROT,:,:) = TRANSPOSE(XMAT_RG_RH(JROT,:,:)) +! +!* H.3 Update structural velocities +! Rotation of hub in RG + XRVEL_RH_RG(JROT,:) = MATMUL(XMAT_RG_RH(JROT,:,:),XRVEL_RH_RN(JROT,:)) & + + XRVEL_RN_RG(JROT,:) +! Translation of hub in RG + ZDIST_NAC_HUB_RG(:) = XPOS_HUB_RG(JROT,:) - XPOS_NELT_RG(JROT,INB_NELT,:) + XTVEL_HUB_RG(JROT,:) = XTVEL_NELT_RG(JROT,INB_NELT,:) + CROSS(XRVEL_RH_RG(JROT,:),ZDIST_NAC_HUB_RG(:)) +! +!* H.4 Print in tecplot file + IF (LTECOUTPTS) THEN + CALL PRINT_TECOUT(ITECOUT, XPOS_HUB_RG(JROT,:)) + END IF +! +!* H.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'HUB',XMAT_RG_RH(JROT,:,:),XPOS_HUB_RG(JROT,:)) + END IF +! +! ---- CYLINDRICAL ---- +!* C.1 Update positions + XPOS_CYL_RG(JROT,:) = XPOS_HUB_RG(JROT,:) & + + MATMUL(XMAT_RG_RH(JROT,:,:),XPOSINI_CYL_RH(JROT,:)) +! +!* C.2 Update orientation + CALL GET_ORI_MAT_X(XANGINI_CYL_RH(JROT,1) + XRVEL_RC_RH(JROT,1)*ZTIME, ZORI_MAT_X) + CALL GET_ORI_MAT_Y(XANGINI_CYL_RH(JROT,2) + XRVEL_RC_RH(JROT,2)*ZTIME, ZORI_MAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_CYL_RH(JROT,3) + XRVEL_RC_RH(JROT,3)*ZTIME, ZORI_MAT_Z) +! Orientation matrix + XMAT_RH_RC(JROT,:,:) = MATMUL(ZORI_MAT_X, MATMUL(ZORI_MAT_Y, ZORI_MAT_Z)) + XMAT_RG_RC(JROT,:,:) = MATMUL(XMAT_RG_RH(JROT,:,:), XMAT_RH_RC(JROT,:,:)) +! +!* C.3 Update structural velocities +! Rotation of cylindrical in RG + XRVEL_RC_RG(JROT,:) = MATMUL(XMAT_RG_RC(JROT,:,:),XRVEL_RC_RH(JROT,:)) & + + XRVEL_RH_RG(JROT,:) +! Translation of hub in RG + ZDIST_HUB_CYL_RG(:) = XPOS_CYL_RG(JROT,:) - XPOS_HUB_RG(JROT,:) + XTVEL_CYL_RG(JROT,:) = XTVEL_HUB_RG(JROT,:) + CROSS(XRVEL_RC_RG(JROT,:),ZDIST_HUB_CYL_RG(:)) +! +!* C.4 Print in tecplot file + IF (LTECOUTPTS) THEN + CALL PRINT_TECOUT(ITECOUT, XPOS_CYL_RG(JROT,:)) + END IF +! +!* C.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'CYL',XMAT_RG_RC(JROT,:,:),XPOS_CYL_RG(JROT,:)) + END IF +! +! ---- ANNULAR ELEMENTS ---- +!* A.0 Positioning sections (cuts) in RG + DO JAZI=1, INB_AELT + DO JRAD=1, INB_RELT+1 + XPOS_SEC_RG(JROT,JAZI,JRAD,:) = XPOS_CYL_RG(JROT,:) & + + MATMUL(XMAT_RG_RC(JROT,:,:),GET_VEC_CART(XPOS_SEC_RC(JROT,JAZI,JRAD,:))) + END DO + + DO JRAD=1, INB_RELT +!* A.1 Positioning sections centers (application points) in RG + XPOS_ELT_RG(JROT,JAZI,JRAD,:) = XPOS_CYL_RG(JROT,:) & + + MATMUL(XMAT_RG_RC(JROT,:,:),GET_VEC_CART(XPOS_ELT_RC(JROT,JAZI,JRAD,:))) +!* A.2 Update orientation +! A.2.1 Orientation in elements transition ref. frame + CALL GET_ORI_MAT_X(XANGINI_ELT_RC(JROT,JAZI,JRAD,1) & + + XRVEL_RA_RC(JROT,JAZI,JRAD,1)*ZTIME, ZORI_TMAT_X) + CALL GET_ORI_MAT_Y(XANGINI_ELT_RC(JROT,JAZI,JRAD,2) & + + XRVEL_RA_RC(JROT,JAZI,JRAD,2)*ZTIME, ZORI_TMAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_ELT_RC(JROT,JAZI,JRAD,3) & + + XRVEL_RA_RC(JROT,JAZI,JRAD,3)*ZTIME, ZORI_TMAT_Z) +! Orientation matrix + XMAT_RC_RTA(JROT,JAZI,JRAD,:,:) = MATMUL(ZORI_TMAT_X, MATMUL(ZORI_TMAT_Y, ZORI_TMAT_Z)) + XMAT_RG_RTA(JROT,JAZI,JRAD,:,:) = MATMUL(XMAT_RG_RC(JROT,:,:), XMAT_RC_RTA(JROT,JAZI,JRAD,:,:)) + XMAT_RTA_RG(JROT,JAZI,JRAD,:,:) = TRANSPOSE(XMAT_RG_RTA(JROT,JAZI,JRAD,:,:)) +! +! A.2.2 Orientation in annular elements ref. frame + + CALL GET_ORI_MAT_X(XANGINI_ELT_RA(JROT,JAZI,JRAD,1), ZORI_MAT_X) + CALL GET_ORI_MAT_Y(XANGINI_ELT_RA(JROT,JAZI,JRAD,2), ZORI_MAT_Y) + CALL GET_ORI_MAT_Z(XANGINI_ELT_RA(JROT,JAZI,JRAD,3), ZORI_MAT_Z) + + XMAT_RTA_RA(JROT,JAZI,JRAD,:,:) = MATMUL(ZORI_MAT_X, MATMUL(ZORI_MAT_Y, ZORI_MAT_Z)) + XMAT_RG_RA(JROT,JAZI,JRAD,:,:) = MATMUL(XMAT_RG_RTA(JROT,JAZI,JRAD,:,:), XMAT_RTA_RA(JROT,JAZI,JRAD,:,:)) + XMAT_RA_RG(JROT,JAZI,JRAD,:,:) = TRANSPOSE(XMAT_RG_RA(JROT,JAZI,JRAD,:,:)) +! + +!* A.3 Update structural velocities +! Rotation of elements in RG + XRVEL_RA_RG(JROT,JAZI,JRAD,:) = XRVEL_RC_RG(JROT,:) & + + MATMUL(XMAT_RG_RA(JROT,JAZI,JRAD,:,:),& + XRVEL_RA_RC(JROT,JAZI,JRAD,:)) +! Translation of elements in RG + ZDIST_CYL_ELT_RG(:) = XPOS_ELT_RG(JROT,JAZI,JRAD,:) - XPOS_CYL_RG(JROT,:) + XTVEL_ELT_RG(JROT,JAZI,JRAD,:) = XTVEL_CYL_RG(JROT,:) & + + CROSS(XRVEL_RA_RG(JROT,JAZI,JRAD,:),ZDIST_CYL_ELT_RG(:)) + XTVEL_ELT_RA(JROT,JAZI,JRAD,:) = MATMUL(XMAT_RA_RG(JROT,JAZI,JRAD,:,:),& + XTVEL_ELT_RG(JROT,JAZI,JRAD,:)) +! +!* A.4 Print in tecplot file + IF (LTECOUTPTS) THEN + CALL PRINT_TECOUT(ITECOUT, XPOS_ELT_RG(JROT,JAZI,JRAD,:)) + END IF +! + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'ELT',XMAT_RG_RA(JROT,JAZI,JRAD,:,:),XPOS_ELT_RG(JROT,JAZI,JRAD,:)) + END IF + END DO ! Radial loop + END DO ! Azimutal loop +END DO ! Rotor loop +! +! Closing tec file +IF (LTECOUTPTS) THEN + CLOSE(ITECOUT) +END IF +IF (LCSVOUTFRM) THEN + CLOSE(ICSVOUT) +END IF +! +END SUBROUTINE EOL_KINE_ADR diff --git a/src/MNH/eol_kine_alm.f90 b/src/MNH/eol_kine_alm.f90 index 45b4cf6795942491c1fac17b449af168bf242e6f..502e92b2213b6f56ddc2a699f001cc75a81734c0 100644 --- a/src/MNH/eol_kine_alm.f90 +++ b/src/MNH/eol_kine_alm.f90 @@ -52,6 +52,7 @@ SUBROUTINE EOL_KINE_ALM(KTCOUNT,KTSUBCOUNT,PTSUBSTEP,PTSTEP) !! Original 04/2017 !! Modification 10/11/20 (PA. Joulin) Updated for a main version ! P. Wautelet 19/07/2021: replace double precision by real to allow MNH_REAL=4 compilation +! PA. Joulin 04/2023: add csv output for 3d frames !! !!--------------------------------------------------------------- ! @@ -61,7 +62,9 @@ SUBROUTINE EOL_KINE_ALM(KTCOUNT,KTSUBCOUNT,PTSUBSTEP,PTSTEP) !* 0.1 Modules USE MODD_EOL_KINE_ALM USE MODD_EOL_ALM +USE MODI_EOL_PRINTER USE MODI_EOL_MATHS +USE MODD_EOL_SHARED_IO, ONLY: LTECOUTPTS,LCSVOUTFRM USE MODD_TIME_n, ONLY : TDTCUR USE MODD_CST, ONLY : XPI ! @@ -94,6 +97,7 @@ INTEGER :: JROT, JBLA, JTELT, JNELT, JBELT ! Loop co INTEGER :: INB_WT, INB_B, INB_BELT ! Total numbers INTEGER :: INB_TELT, INB_NELT ! Total numbers INTEGER :: ITECOUT ! Unit number for Tecplot file +INTEGER :: ICSVOUT ! Unit number for TOW.csv file ! ! !--------------------------------------------------------------- @@ -113,7 +117,10 @@ ZTIME = TDTCUR%xtime+(KTSUBCOUNT)*PTSUBSTEP ! !* 1.3 Tecplotfile : opening + headers IF (LTECOUTPTS) THEN - CALL OPEN_TECOUT(ITECOUT, KTCOUNT, KTSUBCOUNT) + CALL OPEN_TECOUT_ALM(ITECOUT, KTCOUNT, KTSUBCOUNT) +END IF +IF (LCSVOUTFRM) THEN + CALL OPEN_3DCSV_ALM(ICSVOUT, KTCOUNT, KTSUBCOUNT) END IF ! ! @@ -159,6 +166,11 @@ DO JROT=1, INB_WT END DO END IF ! +!* T.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'TOW',XMAT_RG_RT(JROT,:,:), XPOS_TOWO_RG(JROT,:)) + END IF +! ! ! ---- NACELLE ---- ! @@ -199,6 +211,11 @@ DO JROT=1, INB_WT END DO END IF ! +!* N.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'NAC',XMAT_RG_RN(JROT,:,:), XPOS_NACO_RG(JROT,:)) + END IF +! ! ! ---- HUB ---- ! @@ -228,6 +245,11 @@ DO JROT=1, INB_WT CALL PRINT_TECOUT(ITECOUT, XPOS_HUB_RG(JROT,:)) END IF ! +!* H.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'HUB',XMAT_RG_RH(JROT,:,:), XPOS_HUB_RG(JROT,:)) + END IF +! ! ! ---- BLADES ---- ! @@ -258,6 +280,10 @@ DO JROT=1, INB_WT CALL PRINT_TECOUT(ITECOUT, XPOS_BLA_RG(JROT,JBLA,:)) END IF ! +!* B.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'BLA',XMAT_RG_RB(JROT,JBLA,:,:), XPOS_BLA_RG(JROT,JBLA,:)) + END IF ! ! ---- ELEMENTS ---- !* E.0 Positioning sections (cuts) in RG @@ -300,6 +326,11 @@ DO JROT=1, INB_WT CALL PRINT_TECOUT(ITECOUT, XPOS_ELT_RG(JROT,JBLA,JBELT,:)) END IF ! +!* E.5 Print in 3D csv file + IF (LCSVOUTFRM) THEN + CALL PRINT_3DFRM(ICSVOUT,'ELT',XMAT_RG_RE(JROT,JBLA,JBELT,:,:), XPOS_ELT_RG(JROT,JBLA,JBELT,:)) + END IF +! ! ---- Leading Edge and Trailing Edge ---- ! It is just to have a tecplot more beautiful ! For the moment, they are useless for computations @@ -333,5 +364,8 @@ END DO ! Rotor loop IF (LTECOUTPTS) THEN CLOSE(ITECOUT) END IF +IF (LCSVOUTFRM) THEN + CLOSE(ICSVOUT) +END IF ! END SUBROUTINE EOL_KINE_ALM diff --git a/src/MNH/eol_main.f90 b/src/MNH/eol_main.f90 index 6b0633c1dec3746119312794abaa4b75c4a16fed..a921c8ce65091174cba9c21903685d46db4481f1 100644 --- a/src/MNH/eol_main.f90 +++ b/src/MNH/eol_main.f90 @@ -59,6 +59,7 @@ END MODULE MODI_EOL_MAIN !! MODIFICATIONS !! ------------- !! 21/10/20 Original +!! 09/22 H. Toumi : add ADR model !! !!--------------------------------------------------------------- ! @@ -69,6 +70,7 @@ END MODULE MODI_EOL_MAIN ! To work with wind turbines USE MODD_EOL_MAIN USE MODI_EOL_ADNR +USE MODI_EOL_ADR USE MODI_EOL_ALM USE MODI_EOL_SMEAR ! To play with MPI @@ -205,6 +207,12 @@ SELECT CASE(CMETH_EOL) ZUT_M, ZVT_M, ZWT_M, & XFX_RG, XFY_RG, XFZ_RG ) ! + CASE('ADR') ! Actuator Disc Rotating + CALL EOL_ADR(KTCOUNT, PTSTEP, & + PDXX, PDYY, PDZZ, & + ZRHO_M, & + ZUT_M, ZVT_M, ZWT_M, & + XFX_RG, XFY_RG, XFZ_RG ) END SELECT ! !* 3.2 Sharing 3D field diff --git a/src/MNH/eol_maths.f90 b/src/MNH/eol_maths.f90 index 81e2a0a885d7012caff260547722e3829d0e7aee..e38dd234efb9f0b57461ed0a5b1fdb3da85a78d9 100644 --- a/src/MNH/eol_maths.f90 +++ b/src/MNH/eol_maths.f90 @@ -3,8 +3,24 @@ !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt !MNH_LIC for details. version 1. !----------------------------------------------------------------- -! Modifications: -! P. Wautelet 19/07/2021: replace double precision by real to allow MNH_REAL=4 compilation +!! +!!**** *MODI_EOL_MATHS* - +!! +!! PURPOSE +!! ------- +!! Contains all mathematical functions for EOL computations. +!! +!! AUTHOR +!! ------ +!! PA. Joulin *CNRM & IFPEN* +!! +!! +!! MODIFICATIONS +!! ------------- +!! Original 24/01/17 +!! P. Wautelet 19/07/2021: replace double precision by real to allow MNH_REAL=4 compilation +!! H. Toumi 04/23 : add ADR functions +!! !----------------------------------------------------------------- ! ####################### MODULE MODI_EOL_MATHS @@ -37,6 +53,16 @@ SUBROUTINE GET_ORI_MAT_Z(PTHETA, PORI_MAT_Z) REAL, DIMENSION(3,3), INTENT(OUT) :: PORI_MAT_Z ! Matrix END SUBROUTINE GET_ORI_MAT_Z ! +FUNCTION GET_VEC_CYL(VEC_CART) + REAL, DIMENSION(3), INTENT(IN) :: VEC_CART ! cartesian vector + REAL, DIMENSION(3) :: GET_VEC_CYL ! cylindrical vector +END FUNCTION GET_VEC_CYL +! +FUNCTION GET_VEC_CART(VEC_CYL) + REAL, DIMENSION(3), INTENT(IN) :: VEC_CYL ! cartesian vector + REAL, DIMENSION(3) :: GET_VEC_CART ! cylindrical vector +END FUNCTION GET_VEC_CART +! FUNCTION INTERP_SPLCUB(PAV, PX, PY) REAL :: INTERP_SPLCUB ! interface REAL, INTENT(IN) :: PAV ! Abscissa where spline is to be evaluate @@ -157,6 +183,34 @@ END SUBROUTINE GET_ORI_MAT_Z !######################################################### ! !######################################################### +FUNCTION GET_VEC_CYL(VEC_CART) +! Obtain cylindrical coordinates from a cartesian vector +! + REAL, DIMENSION(3), INTENT(IN) :: VEC_CART ! cartesian vector + REAL, DIMENSION(3) :: GET_VEC_CYL ! cylindrical vector +! + GET_VEC_CYL(1) = SQRT(VEC_CART(1)**2+VEC_CART(2)**2) + GET_VEC_CYL(2) = ATAN2(VEC_CART(2),VEC_CART(1)) + GET_VEC_CYL(3) = VEC_CART(3) +! +END FUNCTION GET_VEC_CYL +!######################################################### +! +!######################################################### +FUNCTION GET_VEC_CART(VEC_CYL) +! Obtain cylindrical coordinates from a cartesian vector +! + REAL, DIMENSION(3), INTENT(IN) :: VEC_CYL ! cartesian vector + REAL, DIMENSION(3) :: GET_VEC_CART ! cylindrical vector +! + GET_VEC_CART(1) = VEC_CYL(1)*COS(VEC_CYL(2)) + GET_VEC_CART(2) = VEC_CYL(1)*SIN(VEC_CYL(2)) + GET_VEC_CART(3) = VEC_CYL(3) +! +END FUNCTION GET_VEC_CART +!######################################################### +! +!######################################################### FUNCTION INTERP_SPLCUB(PAV, PX, PY) ! adapted from https://ww2.odu.edu/~agodunov/computing/programs/book2/Ch01/spline.f90 ! diff --git a/src/MNH/eol_printer.f90 b/src/MNH/eol_printer.f90 index b4fe7281673df7ba6b9000ed4353c7430db48423..9bd815e4826690ecdaa394da8dab7aaddfcb3c09 100644 --- a/src/MNH/eol_printer.f90 +++ b/src/MNH/eol_printer.f90 @@ -26,6 +26,40 @@ SUBROUTINE PRINT_DATA_TURBINE_ADNR(KFILE,TPTURBINE) END SUBROUTINE PRINT_DATA_TURBINE_ADNR ! ! *** +! ADR +! *** +! +SUBROUTINE PRINT_DATA_FARM_ADR(KFILE,TPFARM) + USE MODD_EOL_ADR, ONLY: FARM + INTEGER, INTENT(IN) :: KFILE ! output file + TYPE(FARM), INTENT(IN) :: TPFARM ! stored farm data +END SUBROUTINE PRINT_DATA_FARM_ADR +! +SUBROUTINE PRINT_DATA_TURBINE_ADR(KFILE,TPTURBINE) + USE MODD_EOL_ADR, ONLY : TURBINE + INTEGER, INTENT(IN) :: KFILE ! output file + TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data +END SUBROUTINE PRINT_DATA_TURBINE_ADR +! +SUBROUTINE PRINT_DATA_BLADE_ADR(KFILE,TPBLADE) + USE MODD_EOL_ADR, ONLY : BLADE + INTEGER, INTENT(IN) :: KFILE ! output file + TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data +END SUBROUTINE PRINT_DATA_BLADE_ADR +! +SUBROUTINE PRINT_DATA_AIRFOIL_ADR(KFILE,TPAIRFOIL) + USE MODD_EOL_ADR, ONLY : AIRFOIL + INTEGER, INTENT(IN) :: KFILE ! output file + TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored airfoil data +END SUBROUTINE PRINT_DATA_AIRFOIL_ADR +! +SUBROUTINE OPEN_TECOUT_ADR(KFILE, KTCOUNT) + INTEGER, INTENT(IN) :: KFILE ! File index + INTEGER, INTENT(IN) :: KTCOUNT ! Time step index +END SUBROUTINE OPEN_TECOUT_ADR +! +! +! *** ! ALM ! *** ! @@ -53,11 +87,11 @@ SUBROUTINE PRINT_DATA_AIRFOIL_ALM(KFILE,TPAIRFOIL) TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored airfoil data END SUBROUTINE PRINT_DATA_AIRFOIL_ALM ! -SUBROUTINE OPEN_TECOUT(KFILE, KTCOUNT, KTSUBCOUNT) +SUBROUTINE OPEN_TECOUT_ALM(KFILE, KTCOUNT, KTSUBCOUNT) INTEGER, INTENT(IN) :: KFILE ! File index INTEGER, INTENT(IN) :: KTCOUNT ! Time step index INTEGER, INTENT(IN) :: KTSUBCOUNT ! Subtime step index -END SUBROUTINE OPEN_TECOUT +END SUBROUTINE OPEN_TECOUT_ALM ! SUBROUTINE PRINT_TECOUT(KFILE,PVAR) INTEGER, INTENT(IN) :: KFILE ! File index @@ -69,6 +103,25 @@ SUBROUTINE PRINT_TSPLIT(KNBSUBCOUNT,PTSUBSTEP) REAL, INTENT(IN) :: PTSUBSTEP ! sub timestep END SUBROUTINE PRINT_TSPLIT ! +! 3D FRAMES +! +SUBROUTINE OPEN_3DCSV_ADR(KFILE, KTCOUNT) + INTEGER, INTENT(IN) :: KFILE ! File index + INTEGER, INTENT(IN) :: KTCOUNT ! Time step index +END SUBROUTINE OPEN_3DCSV_ADR +! +SUBROUTINE OPEN_3DCSV_ALM(KFILE, KTCOUNT, KTSUBCOUNT) + INTEGER, INTENT(IN) :: KFILE ! File index + INTEGER, INTENT(IN) :: KTCOUNT ! Time step index + INTEGER, INTENT(IN) :: KTSUBCOUNT ! Subtime step index +END SUBROUTINE OPEN_3DCSV_ALM +! +SUBROUTINE PRINT_3DFRM(KFILE,HPART,PMAT,PPOS) + INTEGER, INTENT(IN) :: KFILE ! File index + CHARACTER(LEN=3), INTENT(IN) :: HPART ! Part of the wind turbine + REAL, DIMENSION(3,3), INTENT(IN) :: PMAT ! Frame to plot + REAL, DIMENSION(3), INTENT(IN) :: PPOS ! Frame origin +END SUBROUTINE PRINT_3DFRM ! END INTERFACE ! @@ -88,6 +141,7 @@ END MODULE MODI_EOL_PRINTER !! MODIFICATIONS !! ------------- !! Original 26/10/2020 +!! H. Toumi 09/22 add ADR functions !! !!--------------------------------------------------------------- ! @@ -151,6 +205,135 @@ END IF END SUBROUTINE PRINT_DATA_TURBINE_ADNR !######################################################### ! +! ADR +! +!######################################################### +SUBROUTINE PRINT_DATA_FARM_ADR(KFILE,TPFARM) +! +USE MODD_EOL_ADR, ONLY : FARM +USE MODD_EOL_SHARED_IO, ONLY : CFARM_CSVDATA +USE MODD_VAR_ll, ONLY : IP ! only master cpu +! +IMPLICIT NONE +! +INTEGER, INTENT(IN) :: KFILE ! File index +TYPE(FARM), INTENT(IN) :: TPFARM ! dummy stored farm data +! +INTEGER :: JROT ! Loop index +! +IF (IP==1) THEN + WRITE(KFILE,*) '' + WRITE(KFILE,*) '======================== WIND TURBINE DATA ========================' + WRITE(KFILE,*) '' + WRITE(KFILE,*) '---- Farm ----' + WRITE(KFILE,*) 'Data from file : ', TRIM(CFARM_CSVDATA) + WRITE(KFILE,*) 'Number of turbines : ', TPFARM%NNB_TURBINES + WRITE(KFILE,*) 'Tower base positions (X,Y) [m] : ' + DO JROT=1, TPFARM%NNB_TURBINES + WRITE(KFILE, '(1X,A,I3,A,F10.1,A,F10.1,A)') 'n.', JROT,& + ': (', TPFARM%XPOS_X(JROT),',',TPFARM%XPOS_Y(JROT),')' + END DO + WRITE(KFILE,*) 'Working state (rad/s,rad,rad) : ' + DO JROT=1, TPFARM%NNB_TURBINES + WRITE(KFILE, '(1X,A,I3,A,F10.5,A,F10.5,A,F10.5)') 'n.', JROT,& + ': Omega = ', TPFARM%XOMEGA(JROT), & + ' ; Yaw = ', TPFARM%XNAC_YAW(JROT),& + ' ; Pitch = ', TPFARM%XBLA_PITCH(JROT) + END DO + WRITE(KFILE,*) '' +END IF +! +END SUBROUTINE PRINT_DATA_FARM_ADR +!######################################################### +! +!######################################################### +SUBROUTINE PRINT_DATA_TURBINE_ADR(KFILE,TPTURBINE) +! +USE MODD_EOL_ADR, ONLY : TURBINE +USE MODD_EOL_SHARED_IO, ONLY : CTURBINE_CSVDATA +USE MODD_VAR_ll, ONLY : IP ! only master cpu +! +IMPLICIT NONE +! +INTEGER, INTENT(IN) :: KFILE ! File index +TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! dummy stored turbine data +! +IF (IP==1) THEN + WRITE(KFILE,*) '---- Turbine ----' + WRITE(KFILE,* ) 'Data from file : ', TRIM(CTURBINE_CSVDATA) + WRITE(KFILE,'(1X,A,A10)' ) 'Wind turbine : ', TPTURBINE%CNAME + WRITE(KFILE,'(1X,A,I10)' ) 'Number of blades : ', TPTURBINE%NNB_BLADES + WRITE(KFILE,'(1X,A,F10.1)') 'Hub height [m] : ', TPTURBINE%XH_HEIGHT + WRITE(KFILE,'(1X,A,F10.3)') 'Blade min radius [m] : ', TPTURBINE%XR_MIN + WRITE(KFILE,'(1X,A,F10.3)') 'Blade max radius [m] : ', TPTURBINE%XR_MAX + WRITE(KFILE,'(1X,A,F10.3)') 'Nacelle tilt [rad] : ', TPTURBINE%XNAC_TILT + WRITE(KFILE,'(1X,A,F10.3)') 'Hub deport [m] : ', TPTURBINE%XH_DEPORT + WRITE(KFILE,*) '' +END IF +! +END SUBROUTINE PRINT_DATA_TURBINE_ADR +!######################################################### +! +!######################################################### +SUBROUTINE PRINT_DATA_BLADE_ADR(KFILE,TPBLADE) +! +USE MODD_EOL_ADR, ONLY : BLADE +USE MODD_EOL_SHARED_IO, ONLY : CBLADE_CSVDATA +USE MODD_VAR_ll, ONLY : IP ! only master cpu +! +IMPLICIT NONE +! +INTEGER, INTENT(IN) :: KFILE ! File index +TYPE(BLADE), INTENT(IN) :: TPBLADE ! dummy stored blade data +! +IF (IP==1) THEN + WRITE(KFILE,*) '---- Blade ----' + WRITE(KFILE,* ) 'Data from file : ', TRIM(CBLADE_CSVDATA) + WRITE(KFILE,'(1X,A,I10)' ) 'Nb of data (from data file) : ', TPBLADE%NNB_BLADAT + WRITE(KFILE,'(1X,A,F10.1)') 'First node radius [m] : ', TPBLADE%XRAD(1) + WRITE(KFILE,'(1X,A,F10.1)') 'Last node radius [m] : ', TPBLADE%XRAD(TPBLADE%NNB_BLADAT) + WRITE(KFILE,'(1X,A,F10.1)') 'Chord max. [m] : ', MAXVAL(TPBLADE%XCHORD(:)) + WRITE(KFILE,'(1X,A,I10)' ) 'Nb of radial element (from nam) : ', TPBLADE%NNB_RADELT + WRITE(KFILE,'(1X,A,I10)' ) 'Nb of azimutal element (from nam) : ', TPBLADE%NNB_AZIELT + WRITE(KFILE,*) '' +END IF +! +END SUBROUTINE PRINT_DATA_BLADE_ADR +!######################################################### +! +!######################################################### +SUBROUTINE PRINT_DATA_AIRFOIL_ADR(KFILE,TPAIRFOIL) +! +USE MODD_EOL_ADR, ONLY : AIRFOIL +USE MODD_EOL_SHARED_IO, ONLY : CAIRFOIL_CSVDATA +USE MODD_VAR_ll, ONLY : IP ! only master cpu +! +IMPLICIT NONE +! +INTEGER, INTENT(IN) :: KFILE ! File index +TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! dummy stored airfoil data +! +INTEGER :: JA +! +IF (IP==1) THEN + WRITE(KFILE,*) '---- Airfoils ----' + WRITE(KFILE,* ) 'Data from file : ', TRIM(CAIRFOIL_CSVDATA) + WRITE(KFILE,'(1X,A,I10)' ) 'Nb of airfoils (from data file) : ', SIZE(TPAIRFOIL) + WRITE(KFILE,'(1X,A)' ) 'Different airfoils : ' + DO JA=1,SIZE(TPAIRFOIL) + WRITE(KFILE,'(1X,A,I3,A,A)') 'Airfoil n.', JA,& + ': ', TPAIRFOIL(JA)%CNAME + END DO + WRITE(KFILE,*) '' + WRITE(KFILE,*) '===================================================================' + WRITE(KFILE,*) '' +END IF +! +END SUBROUTINE PRINT_DATA_AIRFOIL_ADR +!######################################################### +! +! ALM +! !######################################################### SUBROUTINE PRINT_DATA_FARM_ALM(KFILE,TPFARM) ! @@ -274,10 +457,47 @@ END IF ! END SUBROUTINE PRINT_DATA_AIRFOIL_ALM !######################################################### +!! +!######################################################### +SUBROUTINE OPEN_TECOUT_ADR(KFILE, KTCOUNT) +! +USE MODD_EOL_ADR, ONLY:TFARM,TTURBINE,TBLADE +! +IMPLICIT NONE +! +INTEGER, INTENT(OUT) :: KFILE ! File index +INTEGER, INTENT(IN) :: KTCOUNT ! Time step index +! +INTEGER :: INB_WT, INB_RADELT, INB_AZIELT ! Total numbers of wind turbines, radial elt, and azimutal elt +INTEGER :: INB_TELT, INB_NELT ! Total numbers of tower elt, and nacelle elt +INTEGER :: ITOTELT ! Total number of points +! +CHARACTER(LEN=1024) :: HFILE ! File name +! +INB_WT = TFARM%NNB_TURBINES +INB_RADELT = TBLADE%NNB_RADELT +INB_AZIELT = TBLADE%NNB_AZIELT +! Hard coded variables, but they will be useful in next updates +INB_TELT = 2 +INB_NELT = 2 +! +ITOTELT = INB_WT*(INB_TELT+INB_NELT+INB_AZIELT*INB_RADELT) +! +! File name and opening +WRITE(HFILE, "(A18,I4.4,A3)") "Tecplot2.0_Output_", KTCOUNT,".tp" +OPEN( NEWUNIT=KFILE, file=HFILE, form="FORMATTED") +! +! Tecplot Header +WRITE(KFILE,*) 'TITLE="Wind Turbines Points"' +WRITE(KFILE,*) 'VARIABLES="X" "Y" "Z"' +WRITE(KFILE,*) 'ZONE I=',ITOTELT,' J=3 K=1 DATAPACKING=POINT' +! +END SUBROUTINE OPEN_TECOUT_ADR +!######################################################### ! ! !######################################################### -SUBROUTINE OPEN_TECOUT(KFILE, KTCOUNT, KTSUBCOUNT) +SUBROUTINE OPEN_TECOUT_ALM(KFILE, KTCOUNT, KTSUBCOUNT) ! USE MODD_EOL_ALM, ONLY:TFARM,TTURBINE,TBLADE ! @@ -311,7 +531,7 @@ WRITE(KFILE,*) 'TITLE="Wind Turbines Points"' WRITE(KFILE,*) 'VARIABLES="X" "Y" "Z"' WRITE(KFILE,*) 'ZONE I=',ITOTELT,' J=3 K=1 DATAPACKING=POINT' ! -END SUBROUTINE OPEN_TECOUT +END SUBROUTINE OPEN_TECOUT_ALM !######################################################### ! !######################################################### @@ -354,4 +574,76 @@ WRITE(KFILE,'(A,A,A,A,A,F10.8,A)') & 'Please, turn on the time-splitting method (LTIMESPLIT=.TRUE.), ',& 'or decrease XTSTEP to a value lower than ', PMAXTSTEP, ' sec.' END SUBROUTINE PRINT_ERROR_WTCFL +!######################################################### +! +!######################################################### +SUBROUTINE OPEN_3DCSV_ALM(KFILE, KTCOUNT, KTSUBCOUNT) +! +IMPLICIT NONE +! +INTEGER, INTENT(OUT) :: KFILE ! File index +INTEGER, INTENT(IN) :: KTCOUNT ! Time step index +INTEGER, INTENT(IN) :: KTSUBCOUNT ! Subtime step index +! +! +CHARACTER(LEN=1024) :: HFILE ! File name +! +! File name and opening +WRITE(HFILE, "(A6,I4.4,I2.2,A4)") "CSV3D_", KTCOUNT, KTSUBCOUNT,".csv" +OPEN( NEWUNIT=KFILE, file=HFILE, form="FORMATTED") +! +! CSV Header +WRITE(KFILE,'(A50)') 'Frame;& + R11;R12;R13;& + R21;R22;R23;& + R31;R32;R33;& + P1;P2;P3' +! +END SUBROUTINE OPEN_3DCSV_ALM +!######################################################### +! +!######################################################### +SUBROUTINE OPEN_3DCSV_ADR(KFILE, KTCOUNT) +! +IMPLICIT NONE +! +INTEGER, INTENT(OUT) :: KFILE ! File index +INTEGER, INTENT(IN) :: KTCOUNT ! Time step index +! +! +CHARACTER(LEN=1024) :: HFILE ! File name +! +! File name and opening +WRITE(HFILE, "(A6,I4.4,A4)") "CSV3D_", KTCOUNT,".csv" +OPEN( NEWUNIT=KFILE, file=HFILE, form="FORMATTED") +! +! CSV Header +WRITE(KFILE,'(A50)') 'Frame;& + R11;R12;R13;& + R21;R22;R23;& + R31;R32;R33;& + P1;P2;P3' +! +END SUBROUTINE OPEN_3DCSV_ADR +!######################################################### +! +!######################################################### +SUBROUTINE PRINT_3DFRM(KFILE,HPART,PMAT,PPOS) +IMPLICIT NONE +INTEGER, INTENT(IN) :: KFILE ! File index +CHARACTER(LEN=3), INTENT(IN) :: HPART ! Part of the wind turbine +REAL, DIMENSION(3,3), INTENT(IN) :: PMAT ! Frame to plot +REAL, DIMENSION(3), INTENT(IN) :: PPOS ! Point to plot +! +! It writes frames informations +WRITE(KFILE, '(A3,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7,A1,E14.7)')& +!'(A3,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F,A1,F)')& + HPART, ';', & + PMAT(1,1),';',PMAT(1,2),';',PMAT(1,3),';', & + PMAT(2,1),';',PMAT(2,2),';',PMAT(2,3),';', & + PMAT(3,1),';',PMAT(3,2),';',PMAT(3,3),';', & + PPOS(1),';',PPOS(2),';',PPOS(3) +! +END SUBROUTINE PRINT_3DFRM +!######################################################### ! diff --git a/src/MNH/eol_reader.f90 b/src/MNH/eol_reader.f90 index 8d5ecc6cfea05ea462138d0f8472d23083e0d7ce..c827b82c6c5cb1106d8acc1e23fc0a2be73cd911 100644 --- a/src/MNH/eol_reader.f90 +++ b/src/MNH/eol_reader.f90 @@ -23,6 +23,45 @@ SUBROUTINE READ_CSVDATA_TURBINE_ADNR(HFILE,TPTURBINE) TYPE(TURBINE), INTENT(OUT) :: TPTURBINE ! stored turbine data END SUBROUTINE READ_CSVDATA_TURBINE_ADNR ! +! ADR +! +SUBROUTINE READ_CSVDATA_FARM_ADR(HFILE,TPFARM) + USE MODD_EOL_ADR, ONLY: FARM + CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read + TYPE(FARM), INTENT(OUT) :: TPFARM ! stored farm data +END SUBROUTINE READ_CSVDATA_FARM_ADR +! +SUBROUTINE READ_CSVDATA_TURBINE_ADR(HFILE,TPTURBINE) + USE MODD_EOL_ADR, ONLY : TURBINE + CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read + TYPE(TURBINE), INTENT(OUT) :: TPTURBINE ! stored turbine data +END SUBROUTINE READ_CSVDATA_TURBINE_ADR + +SUBROUTINE READ_CSVDATA_BLADE_ADR(HFILE,TPTURBINE,TPBLADE) + USE MODD_EOL_ADR, ONLY : TURBINE, BLADE + CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read + TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data + TYPE(BLADE), INTENT(OUT) :: TPBLADE ! stored blade data +END SUBROUTINE READ_CSVDATA_BLADE_ADR +! +SUBROUTINE READ_CSVDATA_AIRFOIL_ADR(HFILE,TPBLADE,TPAIRFOIL) + USE MODD_EOL_ADR, ONLY : BLADE, AIRFOIL + CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read + TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data (to select airfoils) + TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE, INTENT(OUT) :: TPAIRFOIL ! stored airfoil data +END SUBROUTINE READ_CSVDATA_AIRFOIL_ADR +! +FUNCTION GET_AIRFOIL_ID_ADR(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) + USE MODD_EOL_ADR, ONLY : TURBINE, BLADE, AIRFOIL + IMPLICIT NONE + INTEGER :: GET_AIRFOIL_ID_ADR + TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data + TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data + TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored airfoil data + REAL, INTENT(IN) :: PRADIUS ! Radius position studied +END FUNCTION GET_AIRFOIL_ID_ADR +! + ! ALM ! SUBROUTINE READ_CSVDATA_FARM_ALM(HFILE,TPFARM) @@ -57,16 +96,17 @@ SUBROUTINE HOW_MANY_LINES_OF(KLUNAM,HFILE,HNAME,KLINE) CHARACTER(LEN=*), INTENT(IN) :: HNAME ! turbine's name INTEGER, INTENT(OUT) :: KLINE END SUBROUTINE HOW_MANY_LINES_OF +!! ! -FUNCTION GET_AIRFOIL_ID(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) +FUNCTION GET_AIRFOIL_ID_ALM(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) USE MODD_EOL_ALM, ONLY : TURBINE, BLADE, AIRFOIL IMPLICIT NONE - INTEGER :: GET_AIRFOIL_ID + INTEGER :: GET_AIRFOIL_ID_ALM TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored airfoil data REAL, INTENT(IN) :: PRADIUS ! Radius position studied -END FUNCTION GET_AIRFOIL_ID +END FUNCTION GET_AIRFOIL_ID_ALM ! END INTERFACE ! @@ -87,6 +127,7 @@ END MODULE MODI_EOL_READER !! ------------- !! Original 05/2018 !! Modification 21/10/20 (PA. Joulin) Updated for a main version +!! H. Toumi 09/22 add ADR functions !! !!--------------------------------------------------------------- ! @@ -223,9 +264,9 @@ END SUBROUTINE READ_CSVDATA_TURBINE_ADNR !######################################################### ! !######################################################### -SUBROUTINE READ_CSVDATA_FARM_ALM(HFILE,TPFARM) +SUBROUTINE READ_CSVDATA_FARM_ADR(HFILE,TPFARM) ! -USE MODD_EOL_ALM, ONLY: FARM +USE MODD_EOL_ADR, ONLY: FARM USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR ! IMPLICIT NONE @@ -296,13 +337,13 @@ ELSE RETURN END IF ! -END SUBROUTINE READ_CSVDATA_FARM_ALM +END SUBROUTINE READ_CSVDATA_FARM_ADR !######################################################### ! !######################################################### -SUBROUTINE READ_CSVDATA_TURBINE_ALM(HFILE,TPTURBINE) +SUBROUTINE READ_CSVDATA_TURBINE_ADR(HFILE,TPTURBINE) ! -USE MODD_EOL_ALM, ONLY: TURBINE +USE MODD_EOL_ADR, ONLY: TURBINE USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR ! IMPLICIT NONE @@ -367,13 +408,13 @@ ELSE CLOSE(ILU) END IF ! -END SUBROUTINE READ_CSVDATA_TURBINE_ALM +END SUBROUTINE READ_CSVDATA_TURBINE_ADR !######################################################### ! !######################################################### -SUBROUTINE READ_CSVDATA_BLADE_ALM(HFILE,TPTURBINE,TPBLADE) +SUBROUTINE READ_CSVDATA_BLADE_ADR(HFILE,TPTURBINE,TPBLADE) ! -USE MODD_EOL_ALM, ONLY: TURBINE, BLADE, NNB_BLAELT +USE MODD_EOL_ADR, ONLY: TURBINE, BLADE, NNB_RADELT, NNB_AZIELT USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR USE MODI_EOL_ERROR, ONLY: EOL_BLADEDATA_ERROR ! @@ -416,7 +457,8 @@ END DO IF (INBLINE < 2) THEN CALL EOL_CSVEMPTY_ERROR(HFILE,INBLINE) ELSE - TPBLADE%NNB_BLAELT = NNB_BLAELT + TPBLADE%NNB_RADELT = NNB_RADELT + TPBLADE%NNB_AZIELT = NNB_AZIELT ! Saving number of data TPBLADE%NNB_BLADAT = INBLINE - 1 ALLOCATE(TPBLADE%XRAD(TPBLADE%NNB_BLADAT)) @@ -448,13 +490,13 @@ ELSE RETURN END IF ! -END SUBROUTINE READ_CSVDATA_BLADE_ALM +END SUBROUTINE READ_CSVDATA_BLADE_ADR !######################################################### ! !######################################################### -SUBROUTINE READ_CSVDATA_AIRFOIL_ALM(HFILE,TPBLADE,TPAIRFOIL) +SUBROUTINE READ_CSVDATA_AIRFOIL_ADR(HFILE,TPBLADE,TPAIRFOIL) ! -USE MODD_EOL_ALM, ONLY: BLADE, AIRFOIL +USE MODD_EOL_ADR, ONLY: BLADE, AIRFOIL USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR USE MODI_EOL_ERROR, ONLY: EOL_AIRFOILNOTFOUND_ERROR ! @@ -559,9 +601,10 @@ CLOSE(ILU) IF (INBLINE == 0) THEN CALL EOL_AIRFOILNOTFOUND_ERROR(HFILE,YAIRFOIL(IA)) END IF -END SUBROUTINE READ_CSVDATA_AIRFOIL_ALM +END SUBROUTINE READ_CSVDATA_AIRFOIL_ADR !######################################################### ! +! !######################################################### SUBROUTINE HOW_MANY_LINES_OF(KLUNAM,HFILE,HNAME,KLINE) ! @@ -601,9 +644,76 @@ END DO END IF END SUBROUTINE HOW_MANY_LINES_OF !######################################################### +!! +!######################################################### +FUNCTION GET_AIRFOIL_ID_ADR(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) +! Allows to link an airfoil from a TPBLADE, at a specific radius (PRADIUS) +! to the airfoils characteristics (TPAIRFOIL). +! The result is an integer (IAID) that should be used like that : +! IAID = GET_GET_AIRFOIL_ID(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) +! TPAIRFOIL(IAID)%MEMBER_OF_TPAIRFOIL +! +USE MODD_EOL_ADR, ONLY : TURBINE, BLADE, AIRFOIL +! +USE MODE_MSG +! +IMPLICIT NONE +! +TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data +TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data +TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored arifoil data +REAL, INTENT(IN) :: PRADIUS ! Radius position studied +INTEGER :: GET_AIRFOIL_ID_ADR +! +CHARACTER(LEN=:), ALLOCATABLE :: YMSG +CHARACTER(LEN=10) :: YRADIUS, YRMIN, YRMAX +INTEGER :: INB_BDATA ! Total number of blade data +INTEGER :: JBDATA ! Index over blade's data +INTEGER :: JA ! Index over diffetents airfoils +REAL, DIMENSION(SIZE(TPBLADE%XRAD)) :: ZDELTARAD ! 2*ZDELTARAD = section lenght +! +! Checking data +IF ((PRADIUS < TPTURBINE%XR_MIN) .OR. (PRADIUS > TPTURBINE%XR_MAX)) THEN + WRITE( YRADIUS, '( F10.2 )' ) PRADIUS + WRITE( YRMIN, '( F10.2 )' ) TPTURBINE%XR_MIN + WRITE( YRMAX, '( F10.2 )' ) TPTURBINE%XR_MAX + YMSG = 'The studied radius R=' // TRIM( YRADIUS ) // ' is out of blade range : [' // TRIM( YRMIN ) // ';' // TRIM( YRMAX ) // ']' + CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'GET_AIRFOIL_ID_ADR', YMSG ) +END IF +! +! Preliminaires +INB_BDATA = SIZE(TPBLADE%XRAD) +! +! Computes half length of sections +ZDELTARAD(1) = TPBLADE%XRAD(1) - TPTURBINE%XR_MIN +DO JBDATA=2,INB_BDATA-1 + ZDELTARAD(JBDATA) = TPBLADE%XRAD(JBDATA) & + - TPBLADE%XRAD(JBDATA-1) & + - ZDELTARAD(JBDATA-1) +END DO +ZDELTARAD(INB_BDATA) = TPTURBINE%XR_MAX - TPBLADE%XRAD(INB_BDATA) +! +! Looking for the section at r=PRADIUS +DO JBDATA=1,INB_BDATA + IF ((PRADIUS >= TPBLADE%XRAD(JBDATA)-ZDELTARAD(JBDATA)) .AND. & + (PRADIUS < TPBLADE%XRAD(JBDATA)+ZDELTARAD(JBDATA))) THEN +! Looking for the ID of the airfoil of this section + DO JA=1,SIZE(TPAIRFOIL) + IF (TRIM(TPBLADE%CAIRFOIL(JBDATA)) == TRIM(TPAIRFOIL(JA)%CNAME)) THEN + GET_AIRFOIL_ID_ADR = JA + EXIT + END IF + END DO +! + EXIT + END IF +END DO +! +END FUNCTION GET_AIRFOIL_ID_ADR +!######################################################### ! !######################################################### -FUNCTION GET_AIRFOIL_ID(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) +FUNCTION GET_AIRFOIL_ID_ALM(TPTURBINE,TPBLADE,TPAIRFOIL,PRADIUS) ! Allows to link an airfoil from a TPBLADE, at a specific radius (PRADIUS) ! to the airfoils characteristics (TPAIRFOIL). ! The result is an integer (IAID) that should be used like that : @@ -620,7 +730,7 @@ TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data TYPE(AIRFOIL), DIMENSION(:), INTENT(IN) :: TPAIRFOIL ! stored arifoil data REAL, INTENT(IN) :: PRADIUS ! Radius position studied -INTEGER :: GET_AIRFOIL_ID +INTEGER :: GET_AIRFOIL_ID_ALM ! CHARACTER(LEN=:), ALLOCATABLE :: YMSG CHARACTER(LEN=10) :: YRADIUS, YRMIN, YRMAX @@ -635,7 +745,7 @@ IF ((PRADIUS < TPTURBINE%XR_MIN) .OR. (PRADIUS > TPTURBINE%XR_MAX)) THEN WRITE( YRMIN, '( F10.2 )' ) TPTURBINE%XR_MIN WRITE( YRMAX, '( F10.2 )' ) TPTURBINE%XR_MAX YMSG = 'The studied radius R=' // TRIM( YRADIUS ) // ' is out of blade range : [' // TRIM( YRMIN ) // ';' // TRIM( YRMAX ) // ']' - CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'GET_AIRFOIL_ID', YMSG ) + CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'GET_AIRFOIL_ID_ALM', YMSG ) END IF ! ! Preliminaires @@ -657,7 +767,7 @@ DO JBDATA=1,INB_BDATA ! Looking for the ID of the airfoil of this section DO JA=1,SIZE(TPAIRFOIL) IF (TRIM(TPBLADE%CAIRFOIL(JBDATA)) == TRIM(TPAIRFOIL(JA)%CNAME)) THEN - GET_AIRFOIL_ID = JA + GET_AIRFOIL_ID_ALM = JA EXIT END IF END DO @@ -666,5 +776,346 @@ DO JBDATA=1,INB_BDATA END IF END DO ! -END FUNCTION GET_AIRFOIL_ID +END FUNCTION GET_AIRFOIL_ID_ALM +!######################################################### +! +!######################################################### +SUBROUTINE READ_CSVDATA_FARM_ALM(HFILE,TPFARM) +! +USE MODD_EOL_ALM, ONLY: FARM +USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR +! +IMPLICIT NONE +! +CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read +TYPE(FARM), INTENT(OUT) :: TPFARM ! dummy stored data blade +! +LOGICAL :: GEXIST ! Existence of file +! +INTEGER :: ILU ! logical unit of the file +INTEGER :: INBLINE ! Nb of line in csv file +! +CHARACTER(LEN=400) :: YSTRING +! +! Read data +REAL :: ZPOS_X +REAL :: ZPOS_Y +REAL :: ZOMEGA +REAL :: ZYAW +REAL :: ZPITCH +! +! Checking file existence +INQUIRE(FILE=HFILE, EXIST=GEXIST) +IF (.NOT.GEXIST) THEN + CALL EOL_CSVNOTFOUND_ERROR(HFILE) +END IF +! +! Opening the file +OPEN(NEWUNIT=ILU,FILE=HFILE, FORM='formatted', STATUS='OLD') +! Counting number of line +REWIND(ILU) +INBLINE=0 +DO + READ(ILU,END=101,FMT='(A400)') YSTRING + IF (LEN_TRIM(YSTRING) > 0) THEN + INBLINE = INBLINE + 1 + END IF +END DO +! +101 CONTINUE +IF (INBLINE < 2) THEN + CALL EOL_CSVEMPTY_ERROR(HFILE,INBLINE) +ELSE + ! Saving number of wind turbine + TPFARM%NNB_TURBINES = INBLINE - 1 + ! Allocations + ALLOCATE(TPFARM%XPOS_X (TPFARM%NNB_TURBINES)) + ALLOCATE(TPFARM%XPOS_Y (TPFARM%NNB_TURBINES)) + ALLOCATE(TPFARM%XOMEGA (TPFARM%NNB_TURBINES)) + ALLOCATE(TPFARM%XNAC_YAW (TPFARM%NNB_TURBINES)) + ALLOCATE(TPFARM%XBLA_PITCH(TPFARM%NNB_TURBINES)) + ! + ! New read + REWIND(ILU) + READ(ILU,FMT='(A400)') YSTRING ! Header reading + ! + ! Saving data + DO INBLINE=1, TPFARM%NNB_TURBINES + READ(ILU,FMT='(A400)') YSTRING + READ(YSTRING,*) ZPOS_X, ZPOS_Y, ZOMEGA, ZYAW, ZPITCH + TPFARM%XPOS_X(INBLINE) = ZPOS_X + TPFARM%XPOS_Y(INBLINE) = ZPOS_Y + TPFARM%XOMEGA(INBLINE) = ZOMEGA + TPFARM%XNAC_YAW(INBLINE) = ZYAW + TPFARM%XBLA_PITCH(INBLINE) = ZPITCH + END DO + CLOSE(ILU) + RETURN +END IF +! +END SUBROUTINE READ_CSVDATA_FARM_ALM !######################################################### +! +!######################################################### +SUBROUTINE READ_CSVDATA_TURBINE_ALM(HFILE,TPTURBINE) +! +USE MODD_EOL_ALM, ONLY: TURBINE +USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR +! +IMPLICIT NONE +! +CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read +TYPE(TURBINE), INTENT(OUT) :: TPTURBINE ! dummy stored data turbine +! +LOGICAL :: GEXIST ! Existence of file +! +INTEGER :: ILU ! logical unit of the file +INTEGER :: INBLINE ! Nb of line in csv file +! +CHARACTER(LEN=400) :: YSTRING +! +CHARACTER(LEN=80) :: YWT_NAME +INTEGER :: INB_BLADE +REAL :: ZH_HEIGHT +REAL :: ZR_MIN +REAL :: ZR_MAX +REAL :: ZNAC_TILT +REAL :: ZHUB_DEPORT +! +! Checking file existence +INQUIRE(FILE=HFILE, EXIST=GEXIST) +IF (.NOT.GEXIST) THEN + CALL EOL_CSVNOTFOUND_ERROR(HFILE) +END IF +! +! Opening +OPEN(NEWUNIT=ILU,FILE=HFILE, FORM='formatted', STATUS='OLD') +! +! Counting number of line +REWIND(ILU) +INBLINE=0 +DO + READ(ILU,END=101,FMT='(A400)') YSTRING + IF (LEN_TRIM(YSTRING) > 0) THEN + INBLINE = INBLINE + 1 + END IF +END DO +! +101 CONTINUE +IF (INBLINE /= 2) THEN + CALL EOL_CSVEMPTY_ERROR(HFILE,INBLINE) +ELSE + REWIND(ILU) + READ(ILU,FMT='(A400)') YSTRING ! Header reading + READ(ILU,FMT='(A400)') YSTRING ! Reading next line + ! Read data + READ(YSTRING,*) YWT_NAME, INB_BLADE, ZH_HEIGHT,& ! reading data + ZR_MIN, ZR_MAX, ZNAC_TILT, & + ZHUB_DEPORT + TPTURBINE%CNAME = YWT_NAME ! Saving them + TPTURBINE%NNB_BLADES = INB_BLADE + TPTURBINE%XH_HEIGHT = ZH_HEIGHT + TPTURBINE%XR_MIN = ZR_MIN + TPTURBINE%XR_MAX = ZR_MAX + TPTURBINE%XNAC_TILT = ZNAC_TILT + TPTURBINE%XH_DEPORT = ZHUB_DEPORT + REWIND(ILU) + RETURN + CLOSE(ILU) +END IF +! +END SUBROUTINE READ_CSVDATA_TURBINE_ALM +!######################################################### +! +!######################################################### +SUBROUTINE READ_CSVDATA_BLADE_ALM(HFILE,TPTURBINE,TPBLADE) +! +USE MODD_EOL_ALM, ONLY: TURBINE, BLADE, NNB_BLAELT +USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR +USE MODI_EOL_ERROR, ONLY: EOL_BLADEDATA_ERROR +! +CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read +TYPE(TURBINE), INTENT(IN) :: TPTURBINE ! stored turbine data +TYPE(BLADE), INTENT(OUT) :: TPBLADE ! dummy stored data blade +! +LOGICAL :: GEXIST ! Existence of file +! +INTEGER :: ILU ! logical unit of the file +INTEGER :: INBLINE ! Nb of line in csv file +INTEGER :: INBDATA ! Nb of data (line/section) of blade +! +CHARACTER(LEN=400) :: YSTRING +! +REAL :: ZCENTER ! Center pos. of elmt [m] +REAL :: ZCHORD ! Blade chord of elmt [m] +REAL :: ZTWIST ! Twist of elmt [rad] +CHARACTER(LEN=20) :: YAIRFOIL ! Airfoil name [-] +! +! Checking file existence +INQUIRE(FILE=HFILE, EXIST=GEXIST) +IF (.NOT.GEXIST) THEN + CALL EOL_CSVNOTFOUND_ERROR(HFILE) +END IF +! +! Ouverture +OPEN(NEWUNIT=ILU,FILE=HFILE, FORM='formatted', STATUS='OLD') +! Counting number of line +REWIND(ILU) +INBLINE=0 +DO + READ(ILU,END=101,FMT='(A400)') YSTRING + IF (LEN_TRIM(YSTRING) > 0) THEN + INBLINE = INBLINE + 1 + END IF +END DO +! +101 CONTINUE +IF (INBLINE < 2) THEN + CALL EOL_CSVEMPTY_ERROR(HFILE,INBLINE) +ELSE + TPBLADE%NNB_BLAELT = NNB_BLAELT + ! Saving number of data + TPBLADE%NNB_BLADAT = INBLINE - 1 + ALLOCATE(TPBLADE%XRAD(TPBLADE%NNB_BLADAT)) + ALLOCATE(TPBLADE%XCHORD(TPBLADE%NNB_BLADAT)) + ALLOCATE(TPBLADE%XTWIST(TPBLADE%NNB_BLADAT)) + ALLOCATE(TPBLADE%CAIRFOIL(TPBLADE%NNB_BLADAT)) + ! + ! New read + REWIND(ILU) + READ(ILU,FMT='(A400)') YSTRING ! Header reading + ! + ! Saving data + DO INBLINE=1, TPBLADE%NNB_BLADAT + READ(ILU,FMT='(A400)') YSTRING + READ(YSTRING,*) ZCENTER, ZCHORD, ZTWIST, YAIRFOIL ! Reading data + IF ((ZCENTER<=0.0) .OR. (ZCENTER>= 1.0)) THEN + ! Checking data + CALL EOL_BLADEDATA_ERROR(ZCENTER) + ELSE + ! Storing them + TPBLADE%XRAD(INBLINE) = ZCENTER*(TPTURBINE%XR_MAX-TPTURBINE%XR_MIN) & ! Data in % + + TPTURBINE%XR_MIN + TPBLADE%XCHORD(INBLINE) = ZCHORD + TPBLADE%XTWIST(INBLINE) = ZTWIST + TPBLADE%CAIRFOIL(INBLINE) = YAIRFOIL + END IF + END DO + CLOSE(ILU) + RETURN +END IF +! +END SUBROUTINE READ_CSVDATA_BLADE_ALM +!######################################################### +! +!######################################################### +SUBROUTINE READ_CSVDATA_AIRFOIL_ALM(HFILE,TPBLADE,TPAIRFOIL) +! +USE MODD_EOL_ALM, ONLY: BLADE, AIRFOIL +USE MODI_EOL_ERROR, ONLY: EOL_CSVNOTFOUND_ERROR, EOL_CSVEMPTY_ERROR +USE MODI_EOL_ERROR, ONLY: EOL_AIRFOILNOTFOUND_ERROR +! +CHARACTER(LEN=*), INTENT(IN) :: HFILE ! file to read +TYPE(BLADE), INTENT(IN) :: TPBLADE ! stored blade data (to select airfoils) +TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE, INTENT(OUT) :: TPAIRFOIL ! dummy stored data blade +! +LOGICAL :: GEXIST ! Existence of file +! +INTEGER :: ILU ! logical unit of the file +INTEGER :: INBDATA ! Nb of data (line/section) per airfoil +INTEGER :: INBLINE ! Nb of line in csv file +LOGICAL :: GAIRFLAG ! Flag for airfoil counting +! +INTEGER :: JI, JJ, IA ! loop control +INTEGER :: INBAIRFOIL ! Nb of differents airfoils on one blade +! +CHARACTER(LEN=400) :: YSTRING +! +CHARACTER(LEN=15), DIMENSION(:), ALLOCATABLE :: YAIRFOIL +! +CHARACTER(LEN=15) :: YREAD_NAME +REAL :: ZAA ! Attack Angle [rad] +REAL :: ZRE ! Reynolds Number [-] +REAL :: ZCL ! Lift Coef [-] +REAL :: ZCD ! Drag Coef [-] +REAL :: ZCM ! Moment Coef [-] +! +! Checking file existence +INQUIRE(FILE=HFILE, EXIST=GEXIST) +IF (.NOT.GEXIST) THEN + CALL EOL_CSVNOTFOUND_ERROR(HFILE) +END IF +! +! Ouverture +OPEN(NEWUNIT=ILU,FILE=HFILE, FORM='formatted', STATUS='OLD') +! +! 1. Counting number of differents airfoils along the blade and selection : +! +! Allcation of local airfoil array +ALLOCATE(YAIRFOIL(SIZE(TPBLADE%CAIRFOIL))) +! +YAIRFOIL(:) = '' +INBAIRFOIL = 1 +YAIRFOIL(1) = TRIM(TPBLADE%CAIRFOIL(1)) +! +DO JI=1, SIZE(TPBLADE%CAIRFOIL) + GAIRFLAG = .FALSE. + DO JJ=1, INBAIRFOIL + IF (TRIM(TPBLADE%CAIRFOIL(JI)) == TRIM(YAIRFOIL(JJ))) THEN + GAIRFLAG = .TRUE. + END IF + END DO + IF (GAIRFLAG .EQV. .FALSE.) THEN + INBAIRFOIL = INBAIRFOIL + 1 + YAIRFOIL(INBAIRFOIL) = TRIM(TPBLADE%CAIRFOIL(JI)) + END IF +END DO +ALLOCATE(TPAIRFOIL(INBAIRFOIL)) +! +! 2. Reading and storing data : +! +DO IA = 1, INBAIRFOIL + ! Array allocation + CALL HOW_MANY_LINES_OF(ILU,HFILE,YAIRFOIL(IA),INBDATA) + ALLOCATE(TPAIRFOIL(IA)%XAA(INBDATA)) + ALLOCATE(TPAIRFOIL(IA)%XRE(INBDATA)) + ALLOCATE(TPAIRFOIL(IA)%XCL(INBDATA)) + ALLOCATE(TPAIRFOIL(IA)%XCD(INBDATA)) + ALLOCATE(TPAIRFOIL(IA)%XCM(INBDATA)) + ! + REWIND(ILU) + INBLINE = 0 + DO + READ(ILU,END=101,FMT='(A400)') YSTRING ! Header + !* reads the string + IF (LEN_TRIM(YSTRING)>0) THEN + READ(YSTRING,FMT=*) YREAD_NAME + IF (TRIM(YREAD_NAME)==TRIM(YAIRFOIL(IA))) THEN ! Read data + INBLINE = INBLINE + 1 + READ(YSTRING,*) YREAD_NAME, ZAA, ZRE, & + ZCL, ZCD, ZCM + ! Storing data + TPAIRFOIL(IA)%CNAME = YREAD_NAME + TPAIRFOIL(IA)%XAA(INBLINE) = ZAA + TPAIRFOIL(IA)%XRE(INBLINE) = ZRE + TPAIRFOIL(IA)%XCL(INBLINE) = ZCL + TPAIRFOIL(IA)%XCD(INBLINE) = ZCD + TPAIRFOIL(IA)%XCM(INBLINE) = ZCM + ELSE ! The name doesnt appear during a new read.. + IF (INBLINE > 0) THEN ! .. but it has already been found, .. + REWIND(ILU) ! .. so it is the end of the data .. + EXIT ! .. and we can exit :) + END IF + END IF + END IF + END DO +END DO +! +CLOSE(ILU) +101 CONTINUE + IF (INBLINE == 0) THEN + CALL EOL_AIRFOILNOTFOUND_ERROR(HFILE,YAIRFOIL(IA)) + END IF +END SUBROUTINE READ_CSVDATA_AIRFOIL_ALM +!######################################################### +! diff --git a/src/MNH/ini_eol_adr.f90 b/src/MNH/ini_eol_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..2ecf8c777af6fa68543d6e8dce608a4cc55995c8 --- /dev/null +++ b/src/MNH/ini_eol_adr.f90 @@ -0,0 +1,500 @@ +!MNH_LIC Copyright 2018-2021 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. +!----------------------------------------------------------------- +! ########################## + MODULE MODI_INI_EOL_ADR +! ########################## +! +INTERFACE +! +SUBROUTINE INI_EOL_ADR(PDXX,PDYY,PDZZ) +! +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size +! +END SUBROUTINE INI_EOL_ADR +! +END INTERFACE +! +END MODULE MODI_INI_EOL_ADR +! +! ############################################################ + SUBROUTINE INI_EOL_ADR(PDXX,PDYY,PDZZ) +! ############################################################ +! +!!**** *INI_EOL_ADR* - routine to initialize the Actuator Disc +!! with Rotation to simulate wind turbines +!! +!! PURPOSE +!! ------- +!! Routine to initialized the ADR (wind turbine model) variables +!! +!!** METHOD +!! ------ +!! +!! EXTERNAL +!! -------- +!! +!! IMPLICIT ARGUMENTS +!! ------------------ +!! +!! *MODD_EOL_SHARED_IO: +!! *Namelist NAM_EOL_ADNR (INPUT) : +!! CHARACTER(LEN=100) :: CFARM_CSVDATA ! File to read, with farm data +!! CHARACTER(LEN=100) :: CTURBINE_CSVDATA ! File to read, turbine data +!! CHARACTER(LEN=100) :: CBLADE_CSVDATA ! Blade file to read +!! CHARACTER(LEN=100) :: CAIRFOIL_CSVDATA ! Airfoil file to read +!! *for ouputs : +!! REAL, DIMENSION(:), ALLOCATABLE :: XTHRUT ! Thrust [N] +!! REAL, DIMENSION(:), ALLOCATABLE :: XTORQT ! Torque [Nm] +!! REAL, DIMENSION(:), ALLOCATABLE :: XPOWT ! Power [W] +!! REAL, DIMENSION(:), ALLOCATABLE :: XTHRU_SUM ! Sum of thrust (N) +!! REAL, DIMENSION(:), ALLOCATABLE :: XTORQ_SUM ! Sum of torque (Nm) +!! REAL, DIMENSION(:), ALLOCATABLE :: XPOW_SUM ! Sum of power (W) +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Elements radius [m] +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] +!! REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] +!! +! +!! *MODD_EOL_ADR (OUTPUT): +!! TYPE(FARM) :: TFARM +!! TYPE(TURBINE) :: TTURBINE +!! TYPE(BLADE) :: TBLADE +!! TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE :: TAIRFOIL +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_AZI ! Elements azimut [rad] +!! REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RA_GLB ! Aerodyn. force (lift+drag) in RA [N] +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFAERO_BLEQ_RA ! Blade Eq Aero. force in RA [N] +!! REAL, DIMENSION(:,:), ALLOCATABLE :: XAOA_BLEQ_GLB ! Blade Eq. AoA +!! REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFAERO_RA_SUM ! Sum of aero. force in RA [N] +!! INTEGER :: NNB_RADELT ! Number of radial elements +!! INTEGER :: NNB_AZIELT ! Number of azimutal elements +!! +!! *MODD_EOL_KINEMATICS +!! Positions +!! Orientations +!! Velocities +!! +!! REFERENCE +!! --------- +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +!! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!------------------------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ------------ +! +!* 0.1 Modules +! +USE MODD_EOL_ADR +USE MODI_EOL_ERROR, ONLY: EOL_DR_ERROR, EOL_DA_ERROR +USE MODD_EOL_SHARED_IO, ONLY: CFARM_CSVDATA, CTURBINE_CSVDATA +USE MODD_EOL_SHARED_IO, ONLY: CBLADE_CSVDATA, CAIRFOIL_CSVDATA +USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT +USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM +USE MODD_EOL_SHARED_IO, ONLY: XELT_RAD +USE MODD_EOL_SHARED_IO, ONLY: XAOA_GLB, XFLIFT_GLB, XFDRAG_GLB, XFAERO_RG_GLB +USE MODD_EOL_SHARED_IO, ONLY: XAOA_SUM +USE MODI_EOL_READER, ONLY: READ_CSVDATA_FARM_ADR +USE MODI_EOL_READER, ONLY: READ_CSVDATA_TURBINE_ADR +USE MODI_EOL_READER, ONLY: READ_CSVDATA_BLADE_ADR +USE MODI_EOL_READER, ONLY: READ_CSVDATA_AIRFOIL_ADR +USE MODI_EOL_PRINTER, ONLY: PRINT_DATA_FARM_ADR +USE MODI_EOL_PRINTER, ONLY: PRINT_DATA_TURBINE_ADR +USE MODI_EOL_PRINTER, ONLY: PRINT_DATA_BLADE_ADR +USE MODI_EOL_PRINTER, ONLY: PRINT_DATA_AIRFOIL_ADR +USE MODD_EOL_KINE_ADR +USE MODI_EOL_MATHS +! To print in output listing +USE MODD_LUNIT_n, ONLY: TLUOUT +! Constant +USE MODD_CST, ONLY: XPI +! To know the grid +USE MODD_GRID_n, ONLY: XXHAT,XYHAT,XZS +USE MODE_ll, ONLY: GET_INDICE_ll +USE MODD_PARAMETERS, ONLY: JPVEXT +! MPI stuffs +USE MODD_VAR_ll, ONLY: NMNH_COMM_WORLD +USE MODD_PRECISION, ONLY: MNHREAL_MPI +USE MODD_MPIF, ONLY: MPI_SUM +USE MODE_SUM_ll, ONLY: MIN_ll +! +!* 0.2 Variables +! +IMPLICIT NONE +! Interface +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size +! +! Some loop controlers +! .. for wind turbines +INTEGER :: JROT ! Rotor index +INTEGER :: JTELT ! Tower element index +INTEGER :: JNELT ! Nacelle element index +INTEGER :: JRAD ! Radius index +INTEGER :: JAZI ! Azimut index +! .. for domain +INTEGER :: IIB,IJB,IKB ! Begin of a CPU domain +INTEGER :: IIE,IJE ! End of a CPU domain +INTEGER :: JI,JJ ! Domain index +! Some variables to be coder-friendly +INTEGER :: INB_WT, INB_B ! Total numbers of wind turbines, blade +INTEGER :: INB_RELT, INB_AELT ! Total numbers of radial elt and azimut elt +INTEGER :: INB_TELT, INB_NELT ! Total numbers of tower elt and nacelle elt +REAL :: ZRAD ! Radius along the blade +REAL :: ZDELTA_AZI, ZDELTA_RAD ! Azimuthal and radial step +REAL :: ZMIN_SIZE ! Minimum size of a mesh side +REAL :: ZDELTA_RAD_MAX, ZDELTA_AZI_MAX ! Max size of ADR element +INTEGER :: INB_RELT_MIN, INB_AELT_MIN ! Minimum number of ADR mesh elements +! Tower base folowing the terrain +REAL,DIMENSION(:,:),ALLOCATABLE :: ZPOSINI_TOWO_RG ! Initial tower origin position +INTEGER :: IINFO ! code info return +! +!------------------------------------------------------------------- +! +!* 1. READING AND ALLOCATING DATA +! --------------------------- +! Reading in csv files +! Allocation of TFARM, TTURBINE, TBLADE and TAIRFOILS inside the function +! +!* 1.1 Wind farm data +! +CALL READ_CSVDATA_FARM_ADR(TRIM(CFARM_CSVDATA),TFARM) +! +!* 1.2 Wind turbine data +! +CALL READ_CSVDATA_TURBINE_ADR(TRIM(CTURBINE_CSVDATA),TTURBINE) +! +!* 1.3 Blade data +! +CALL READ_CSVDATA_BLADE_ADR(TRIM(CBLADE_CSVDATA),TTURBINE,TBLADE) +! +!* 1.4 Airfoil data +! +CALL READ_CSVDATA_AIRFOIL_ADR(TRIM(CAIRFOIL_CSVDATA),TBLADE,TAIRFOIL) +! +! +!------------------------------------------------------------------- +! +!* 2. PRINTING DATA +! ------------- +! +!* 2.1 Wind farm data +! +CALL PRINT_DATA_FARM_ADR(TLUOUT%NLU,TFARM) +! +!* 2.2 Wind turbine data +! +CALL PRINT_DATA_TURBINE_ADR(TLUOUT%NLU,TTURBINE) +! +!* 2.3 Blade data +! +CALL PRINT_DATA_BLADE_ADR(TLUOUT%NLU,TBLADE) +! +!* 2.4 Airfoil data +! +CALL PRINT_DATA_AIRFOIL_ADR(TLUOUT%NLU,TAIRFOIL) +! +! +!------------------------------------------------------------------- +! +!* 3. ALLOCATING VARIABLES +! -------------------- +! +!* 3.0 Preliminaries +! +!* 3.0.a) Shortening the names of the most frequently used var +INB_WT = TFARM%NNB_TURBINES +INB_B = TTURBINE%NNB_BLADES +INB_AELT = TBLADE%NNB_AZIELT +INB_RELT = TBLADE%NNB_RADELT +INB_TELT = 2 ! Hard coded variable, usefull for next updates +INB_NELT = 2 ! Hard coded variable, usefull for next updates +! +!* 3.0.b) Mesh variable +ZDELTA_AZI = 2d0*XPI/INB_AELT ! Rotor disc azimutal devision +ZDELTA_RAD = (TTURBINE%XR_MAX - TTURBINE%XR_MIN)/INB_RELT ! Rotor disc radialal devision +! +!* 3.0.c) Checking mesh +! Finding the minimum size of a cell, in the whole domain. +! (Later, the research could be done only in wind turbine area) +ZMIN_SIZE = MIN(MIN_ll(PDXX(:,:,:),IINFO),& + MIN_ll(PDYY(:,:,:),IINFO),& + MIN_ll(PDZZ(:,:,:),IINFO)) +! 3.0.c)i) Check the number of radial elements: +ZDELTA_RAD_MAX = ZMIN_SIZE +IF (ZDELTA_RAD > ZDELTA_RAD_MAX) THEN + INB_RELT_MIN = INT(CEILING((TTURBINE%XR_MAX - TTURBINE%XR_MIN)/ZDELTA_RAD_MAX)) ! User proper value + CALL EOL_DR_ERROR(INB_RELT_MIN) +END IF +! 3.0.c)ii) Check the number of azimutal elements: +ZDELTA_AZI_MAX = ATAN2(ZMIN_SIZE,TTURBINE%XR_MAX) +IF (ZDELTA_AZI > ZDELTA_AZI_MAX) THEN + INB_AELT_MIN = INT(CEILING(2d0*XPI/ZDELTA_AZI_MAX)) ! User proper value + CALL EOL_DA_ERROR(INB_AELT_MIN) +END IF +! +!* 3.1 MODD_EOL_ADR variables +! at t +ALLOCATE(XELT_RAD (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XELT_AZI (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XAOA_GLB (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XAOA_BLEQ_GLB (INB_WT, INB_RELT )) +ALLOCATE(XFDRAG_GLB (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XFLIFT_GLB (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XFAERO_RA_GLB (INB_WT,INB_AELT,INB_RELT,3)) +ALLOCATE(XFAERO_BLEQ_RA_GLB (INB_WT, INB_RELT,3)) +ALLOCATE(XFAERO_RG_GLB (INB_WT,INB_AELT,INB_RELT,3)) +ALLOCATE(XTHRUT (INB_WT )) +ALLOCATE(XTORQT (INB_WT )) +ALLOCATE(XPOWT (INB_WT )) +! for mean values +ALLOCATE(XAOA_SUM (INB_WT,INB_AELT,INB_RELT )) +ALLOCATE(XFAERO_RA_SUM (INB_WT,INB_AELT,INB_RELT,3)) +ALLOCATE(XTHRU_SUM (INB_WT )) +ALLOCATE(XTORQ_SUM (INB_WT )) +ALLOCATE(XPOW_SUM (INB_WT )) +ALLOCATE(XAOA_BLEQ_SUM (INB_WT, INB_RELT )) +ALLOCATE(XFAERO_BLEQ_RA_SUM (INB_WT, INB_RELT,3)) + +! +!* 3.2 MODD_EOL_KINE_ADR variables +! +! - ORIENTATION MATRIX - +ALLOCATE(XMAT_RG_RT(INB_WT,3,3)) +ALLOCATE(XMAT_RG_RN(INB_WT,3,3)) +ALLOCATE(XMAT_RT_RN(INB_WT,3,3)) +ALLOCATE(XMAT_RG_RH(INB_WT,3,3)) +ALLOCATE(XMAT_RH_RG(INB_WT,3,3)) +ALLOCATE(XMAT_RN_RH(INB_WT,3,3)) +ALLOCATE(XMAT_RG_RC(INB_WT,3,3)) +ALLOCATE(XMAT_RH_RC(INB_WT,3,3)) +ALLOCATE(XMAT_RG_RTA(INB_WT,INB_AELT,INB_RELT,3,3)) +ALLOCATE(XMAT_RTA_RG(INB_WT,INB_AELT,INB_RELT,3,3)) +ALLOCATE(XMAT_RC_RTA(INB_WT,INB_AELT,INB_RELT,3,3)) +ALLOCATE(XMAT_RG_RA(INB_WT,INB_AELT,INB_RELT,3,3)) +ALLOCATE(XMAT_RA_RG(INB_WT,INB_AELT,INB_RELT,3,3)) +ALLOCATE(XMAT_RTA_RA(INB_WT,INB_AELT,INB_RELT,3,3)) +! +! - POSITIONS & ORIENTATIONS - +! Tower +ALLOCATE(ZPOSINI_TOWO_RG (INB_WT,3) ) ! Initial tower origin pos. in RG +ALLOCATE(XPOSINI_TOWO_RG (INB_WT,3) ) ! Initial tower origin pos. in RG +ALLOCATE(XPOS_TOWO_RG (INB_WT,3) ) ! Current tower origin pos. in RG +ALLOCATE(XPOS_TELT_RG (INB_WT,INB_TELT,3) ) ! Current tower element pos. in RG +ALLOCATE(XPOS_TELT_RT (INB_WT,INB_TELT,3) ) ! Current tower element pos. in RT +ALLOCATE(XANGINI_TOW_RG (INB_WT,3) ) ! Initial tower ori. in RG +! Nacelle +ALLOCATE(XPOSINI_NACO_RT (INB_WT,3) ) ! Initial nacelle origin pos. in RT +ALLOCATE(XPOS_NACO_RG (INB_WT,3) ) ! Current nacelle origin pos. in RG +ALLOCATE(XPOS_NELT_RG (INB_WT,INB_NELT,3) ) ! Current nacelle element pos. in RG +ALLOCATE(XPOS_NELT_RN (INB_WT,INB_NELT,3) ) ! Current nacelle element pos. in RN +ALLOCATE(XANGINI_NAC_RT (INB_WT,3) ) ! Initial nacelle ori. in RT +! Hub +ALLOCATE(XPOSINI_HUB_RN (INB_WT,3) ) ! Initial hub pos. in RN +ALLOCATE(XPOS_HUB_RG (INB_WT,3) ) ! Current hub pos. in RG +ALLOCATE(XANGINI_HUB_RN (INB_WT,3) ) ! Initial hub ori. in RN +! Cylindrical +ALLOCATE(XPOSINI_CYL_RH (INB_WT,3) ) ! Initial cylindrical center pos. in RH +ALLOCATE(XPOS_CYL_RG (INB_WT,3) ) ! Current cylindrical center pos. in RG +ALLOCATE(XANGINI_CYL_RH (INB_WT,3) ) ! Initial cylindrical center ori. RH +! Annular Element +ALLOCATE(XPOS_ELT_RC (INB_WT,INB_AELT,INB_RELT,3) ) ! Element pos. in RC +ALLOCATE(XPOS_ELT_RG (INB_WT,INB_AELT,INB_RELT,3) ) ! Element pos. in RG +ALLOCATE(XPOS_SEC_RC (INB_WT,INB_AELT,INB_RELT+1,3)) ! Section pos. in RC +ALLOCATE(XPOS_SEC_RG (INB_WT,INB_AELT,INB_RELT+1,3)) ! Section pos. in RG +ALLOCATE(XANGINI_ELT_RC (INB_WT,INB_AELT,INB_RELT,3)) ! Initial element ori. in RC +ALLOCATE(XANGINI_ELT_RA (INB_WT,INB_AELT,INB_RELT,3)) ! Initial element ori. in transition ref. +ALLOCATE(XTWIST_ELT (INB_WT,INB_RELT) ) ! Element twist in RC +ALLOCATE(XCHORD_ELT (INB_WT,INB_RELT) ) ! Element chord lenght +ALLOCATE(XSURF_ELT (INB_WT,INB_RELT) ) ! Element lift surface +ALLOCATE(XSURF_APP_ELT (INB_WT,INB_RELT) ) ! Element lift surface +! +! - STRUCTURAL VELOCITIES - +! Tower +ALLOCATE(XTVEL_TOWO_RG (INB_WT,3) ) ! Tower base trans. vel. in RG +ALLOCATE(XTVEL_TELT_RG (INB_WT,INB_TELT,3) ) ! Tower element trans. vel. in RG +ALLOCATE(XRVEL_RT_RG (INB_WT,3) ) ! RT/RG rot. vel. +! Nacelle +ALLOCATE(XTVEL_NACO_RT (INB_WT,3) ) ! Nacelle base trans. vel. in RT +ALLOCATE(XTVEL_NELT_RG (INB_WT,INB_NELT,3) ) ! Nacelle element trans. vel. in RG +ALLOCATE(XRVEL_RN_RT (INB_WT,3) ) ! RN/RT rot. vel. +ALLOCATE(XRVEL_RN_RG (INB_WT,3) ) ! RN/RG rot. vel. +! Hub +ALLOCATE(XTVEL_HUB_RN (INB_WT,3) ) ! Hub base trans. vel. in RN +ALLOCATE(XTVEL_HUB_RG (INB_WT,3) ) ! Hub base trans. vel. in RG +ALLOCATE(XRVEL_RH_RN (INB_WT,3) ) ! RH/RT rot. vel. +ALLOCATE(XRVEL_RH_RG (INB_WT,3) ) ! RH/RG rot. vel. +! Cylindrical +ALLOCATE(XTVEL_CYL_RH (INB_WT,3) ) ! Cylindrical base trans. vel. in RH +ALLOCATE(XTVEL_CYL_RG (INB_WT,3) ) ! Cylindrical base trans. vel. in RG +ALLOCATE(XRVEL_RC_RH (INB_WT,3) ) ! RC/RH rot. vel. +ALLOCATE(XRVEL_RC_RG (INB_WT,3) ) ! RC/RG rot. vel. +! Annular Elements +ALLOCATE(XTVEL_ELT_RC (INB_WT,INB_AELT,INB_RELT,3) ) ! Element base trans. vel. in RC +ALLOCATE(XTVEL_ELT_RG (INB_WT,INB_AELT,INB_RELT,3) ) ! Element base trans. vel. in RG +ALLOCATE(XTVEL_ELT_RA (INB_WT,INB_AELT,INB_RELT,3) ) ! Element base trans. vel. in RA +ALLOCATE(XRVEL_RA_RC (INB_WT,INB_AELT,INB_RELT,3) ) ! RA/RB rot. vel. +ALLOCATE(XRVEL_RA_RG (INB_WT,INB_AELT,INB_RELT,3) ) ! RA/RG rot. vel. +! +!------------------------------------------------------------------- +! +!* 4. FIRST BUILDING OF WIND TURBINES +! ------------------------------- +! +!* 4.1 Preliminaries + + +CALL GET_INDICE_ll(IIB,IJB,IIE,IJE) ! Get begin and end domain index (CPU) +IKB=1+JPVEXT ! Vertical begin index +! +XPOS_REF(:) = 0d0 ! Global Origin +! +! +!* 4.2 Tower +DO JROT=1, INB_WT +! Velocities (Rotor, (XYZ)) + XTVEL_TOWO_RG(JROT,:) = 0d0 + XRVEL_RT_RG(JROT,:) = 0d0 +! Positions + ! Init + ZPOSINI_TOWO_RG(JROT,1) = 0 + ZPOSINI_TOWO_RG(JROT,2) = 0 + ZPOSINI_TOWO_RG(JROT,3) = 0 + ! Finding the elevation at the WT position + DO JJ=IJB,IJE + DO JI=IIB,IIE + IF (TFARM%XPOS_X(JROT) >= XXHAT(JI) .AND. & + TFARM%XPOS_X(JROT) < XXHAT(JI) + PDXX(JI,JJ,IKB)) THEN + IF (TFARM%XPOS_Y(JROT) >= XYHAT(JJ) .AND. & + TFARM%XPOS_Y(JROT) < XYHAT(JJ) + PDYY(JI,JJ,IKB)) THEN + ! Horizontal position + ZPOSINI_TOWO_RG(JROT,1) = TFARM%XPOS_X(JROT) ! Tower base position + ZPOSINI_TOWO_RG(JROT,2) = TFARM%XPOS_Y(JROT) + ! Finding the elevation at the WT position + ZPOSINI_TOWO_RG(JROT,3) = XZS(JI,JJ) + END IF + END IF + END DO + END DO + ! Sharing information + CALL MPI_ALLREDUCE(ZPOSINI_TOWO_RG, XPOSINI_TOWO_RG, SIZE(XPOSINI_TOWO_RG),& + MNHREAL_MPI,MPI_SUM,NMNH_COMM_WORLD,IINFO) + ! Tower elements + DO JTELT=1, INB_TELT + XPOS_TELT_RT(JROT,JTELT,1) = 0d0 + XPOS_TELT_RT(JROT,JTELT,2) = 0d0 + XPOS_TELT_RT(JROT,JTELT,3) = (JTELT-1)*(TTURBINE%XH_HEIGHT)/(INB_TELT-1) + END DO + ! Angles + XANGINI_TOW_RG(JROT,1) = 0d0 + XANGINI_TOW_RG(JROT,2) = 0d0 + XANGINI_TOW_RG(JROT,3) = 0d0 +! +! +!* 4.3 Nacelle +! Velocities + XTVEL_NACO_RT(JROT,:) = 0d0 + XRVEL_RN_RT(JROT,:) = 0d0 +! Positions (Rotor, (XYZ)) ! From last point of tower + ! Origin + XPOSINI_NACO_RT(JROT,:) = 0d0 ! Distance between nacelle base and tower top + ! Elements + DO JNELT=1, INB_NELT + XPOS_NELT_RN(JROT,JNELT,1) = (JNELT-1)*TTURBINE%XH_DEPORT/(INB_NELT-1) + XPOS_NELT_RN(JROT,JNELT,2) = 0d0 + XPOS_NELT_RN(JROT,JNELT,3) = 0d0 + END DO + ! Angles + XANGINI_NAC_RT(JROT,1) = 0d0 + XANGINI_NAC_RT(JROT,2) = 0d0 + XANGINI_NAC_RT(JROT,3) = XPI + TFARM%XNAC_YAW(JROT) +! +! +!* 4.4 Hub +! Velocities + XTVEL_HUB_RN(JROT,:) = 0d0 + XRVEL_RH_RN(JROT,1) = 0d0 + XRVEL_RH_RN(JROT,2) = 0d0 + XRVEL_RH_RN(JROT,3) = TFARM%XOMEGA(JROT) +! Position (Rotor, (XYZ)) ! From nacelle last point + XPOSINI_HUB_RN(JROT,1) = 0d0 + XPOSINI_HUB_RN(JROT,2) = 0d0 + XPOSINI_HUB_RN(JROT,3) = 0d0 + XANGINI_HUB_RN(JROT,1) = 0d0 + XANGINI_HUB_RN(JROT,2) = XPI/2d0 - TTURBINE%XNAC_TILT + XANGINI_HUB_RN(JROT,3) = 0d0 +! +! +!* 4.5 Cylindrical +! Velocities + XRVEL_RC_RH(JROT,:) = 0d0 + XTVEL_CYL_RH(JROT,:) = 0d0 +! Position + XPOSINI_CYL_RH(JROT,:) = 0d0 + XANGINI_CYL_RH(JROT,1) = 0d0 + XANGINI_CYL_RH(JROT,2) = 0d0 + XANGINI_CYL_RH(JROT,3) = XPI +! +! +!* 4.6 Annular Elements + + DO JAZI=1, INB_AELT + + + DO JRAD=1, INB_RELT+1 + +! - Positioning of sections (cuts) + + XPOS_SEC_RC(JROT,JAZI,JRAD,1) = TTURBINE%XR_MIN + (JRAD-1) * ZDELTA_RAD + XPOS_SEC_RC(JROT,JAZI,JRAD,2) = (JAZI-1) * ZDELTA_AZI + XPOS_SEC_RC(JROT,JAZI,JRAD,3) = 0d0 + END DO + + DO JRAD=1, INB_RELT +! - Positioning of centers (points of application) + XPOS_ELT_RC(JROT,JAZI,JRAD,1) = XPOS_SEC_RC(JROT,JAZI,JRAD,1) + ZDELTA_RAD/2d0 + XPOS_ELT_RC(JROT,JAZI,JRAD,2) = XPOS_SEC_RC(JROT,JAZI,JRAD,2) + ZDELTA_AZI/2d0 + XPOS_ELT_RC(JROT,JAZI,JRAD,3) = 0d0 + XELT_RAD(JROT,JAZI,JRAD) = XPOS_ELT_RC(JROT,JAZI,JRAD,1) + XELT_AZI(JROT,JAZI,JRAD) = XPOS_ELT_RC(JROT,JAZI,JRAD,2) + + +! - Calculating chord and twist + ZRAD = XELT_RAD(JROT,JAZI,JRAD) + XCHORD_ELT(JROT,JRAD) = INTERP_SPLCUB(ZRAD, TBLADE%XRAD, TBLADE%XCHORD) + XTWIST_ELT(JROT,JRAD) = INTERP_SPLCUB(ZRAD, TBLADE%XRAD, TBLADE%XTWIST) +! - Calculating lifting surface + XSURF_ELT(JROT,JRAD) = XCHORD_ELT(JROT,JRAD) & + * (XPOS_SEC_RC(JROT,JAZI,JRAD+1,1) & + - XPOS_SEC_RC(JROT,JAZI,JRAD,1)) + XSURF_APP_ELT(JROT,JRAD) = INB_B * XSURF_ELT(JROT,JRAD) & + * ZDELTA_AZI/(2d0*XPI) + +! - Velocities + XRVEL_RA_RC(JROT,JAZI,JRAD,:) = 0d0 + XTVEL_ELT_RC(JROT,JAZI,JRAD,:) = 0d0 + +! - Orientation +! - Orientation in the transition reference frame + XANGINI_ELT_RC(JROT,JAZI,JRAD,1) = 0d0 + XANGINI_ELT_RC(JROT,JAZI,JRAD,2) = 0d0 + XANGINI_ELT_RC(JROT,JAZI,JRAD,3) = 0d0 + (JAZI-1) * ZDELTA_AZI + ZDELTA_AZI/2d0 +! - Orientation in the annular elements reference frame + XANGINI_ELT_RA(JROT,JAZI,JRAD,1) = 0d0 + TFARM%XBLA_PITCH(JROT) + XTWIST_ELT(JROT,JRAD) + XANGINI_ELT_RA(JROT,JAZI,JRAD,2) = 0d0 + XANGINI_ELT_RA(JROT,JAZI,JRAD,3) = 0d0 + END DO ! Loop radial + END DO ! Loop azimutal +END DO ! Loop turbine +! +END SUBROUTINE INI_EOL_ADR diff --git a/src/MNH/ini_eol_alm.f90 b/src/MNH/ini_eol_alm.f90 index 41938d789ae7c876168bdf7c62004aa320eb2615..301f7d35d7ee1fcd1f1efc6f0cf36813e7e5a78a 100644 --- a/src/MNH/ini_eol_alm.f90 +++ b/src/MNH/ini_eol_alm.f90 @@ -9,9 +9,9 @@ ! INTERFACE ! -SUBROUTINE INI_EOL_ALM(PDXX,PDYY) +SUBROUTINE INI_EOL_ALM(PDXX,PDYY,PDZZ) ! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY ! mesh size +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size ! END SUBROUTINE INI_EOL_ALM ! @@ -20,7 +20,7 @@ END INTERFACE END MODULE MODI_INI_EOL_ALM ! ! ############################################################ - SUBROUTINE INI_EOL_ALM(PDXX,PDYY) + SUBROUTINE INI_EOL_ALM(PDXX,PDYY,PDZZ) ! ############################################################ ! !!**** *INI_EOL_ALM* - routine to initialize the Actuator Line Model @@ -86,6 +86,7 @@ END MODULE MODI_INI_EOL_ALM !! ------------- !! Original 31/05/18 !! Modification 10/11/20 (PA. Joulin) Updated for a main version +!! Modification 04/23 (H. Toumi) Modified shared_io var !! !------------------------------------------------------------------------------- ! @@ -95,12 +96,16 @@ END MODULE MODI_INI_EOL_ALM !* 0.1 Modules ! USE MODD_EOL_ALM +USE MODI_EOL_ERROR, ONLY: EOL_DR_ERROR USE MODD_EOL_SHARED_IO, ONLY: CFARM_CSVDATA USE MODD_EOL_SHARED_IO, ONLY: CTURBINE_CSVDATA USE MODD_EOL_SHARED_IO, ONLY: CBLADE_CSVDATA USE MODD_EOL_SHARED_IO, ONLY: CAIRFOIL_CSVDATA USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM +USE MODD_EOL_SHARED_IO, ONLY: XELT_RAD +USE MODD_EOL_SHARED_IO, ONLY: XAOA_GLB, XFLIFT_GLB, XFDRAG_GLB, XFAERO_RG_GLB +USE MODD_EOL_SHARED_IO, ONLY: XAOA_SUM USE MODI_EOL_READER, ONLY: READ_CSVDATA_FARM_ALM USE MODI_EOL_READER, ONLY: READ_CSVDATA_TURBINE_ALM USE MODI_EOL_READER, ONLY: READ_CSVDATA_BLADE_ALM @@ -123,12 +128,13 @@ USE MODD_PARAMETERS, ONLY: JPVEXT USE MODD_VAR_ll, ONLY: NMNH_COMM_WORLD USE MODD_PRECISION, ONLY: MNHREAL_MPI USE MODD_MPIF, ONLY: MPI_SUM +USE MODE_SUM_ll, ONLY: MIN_ll ! !* 0.2 Variables ! IMPLICIT NONE ! Interface -REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY ! mesh size +REAL, DIMENSION(:,:,:), INTENT(IN) :: PDXX,PDYY,PDZZ ! mesh size ! ! Some loop controlers ! .. for wind turbines @@ -145,6 +151,11 @@ INTEGER :: JI,JJ ! Domain index INTEGER :: INB_WT, INB_B, INB_BELT ! Total numbers of wind turbines, blades, and blade elt INTEGER :: INB_TELT, INB_NELT ! Total numbers of tower elt, and nacelle elt REAL :: ZRAD ! Radius along the blade +REAL :: ZDELTA_RAD ! Radial size of blade element +REAL :: ZMIN_SIZE ! Minimum size of a mesh side +REAL :: ZDELTA_RAD_MAX ! Max size of ALM rad element +INTEGER :: INB_BELT_MIN ! Minimum number of ALM mesh elements +! ! Tower base folowing the terrain REAL,DIMENSION(:,:),ALLOCATABLE :: ZPOSINI_TOWO_RG ! Initial tower origin position INTEGER :: IINFO ! code info return @@ -201,12 +212,29 @@ CALL PRINT_DATA_AIRFOIL_ALM(TLUOUT%NLU,TAIRFOIL) ! -------------------- ! !* 3.0 Preliminaries +! +!* 3.0.a) Shortening the names of the most frequently used var INB_WT = TFARM%NNB_TURBINES INB_B = TTURBINE%NNB_BLADES -INB_BELT = TBLADE%NNB_BLAELT -! Hard coded variables, but they will be usefull in next updates -INB_TELT = 2 -INB_NELT = 2 +INB_BELT = TBLADE%NNB_BLAELT +INB_TELT = 2 ! Hard coded variable, usefull for next updates +INB_NELT = 2 ! Hard coded variable, usefull for next updates +! +!* 3.0.b) Mesh variable +ZDELTA_RAD = (TTURBINE%XR_MAX-TTURBINE%XR_MIN)/INB_BELT ! Radial element size +! +!* 3.0.c) Checking mesh +! Finding the minimum size of a cell, in the whole domain. +! Later, the research could be done only in wind turbine area. +ZMIN_SIZE = MIN(MIN_ll(PDXX(:,:,:),IINFO),& + MIN_ll(PDYY(:,:,:),IINFO),& + MIN_ll(PDZZ(:,:,:),IINFO)) +! 3.0.c)i) Check the number of radial elements: +ZDELTA_RAD_MAX = ZMIN_SIZE +IF (ZDELTA_RAD > ZDELTA_RAD_MAX) THEN + INB_BELT_MIN = INT(CEILING((TTURBINE%XR_MAX-TTURBINE%XR_MIN)/ZDELTA_RAD_MAX)) ! User proper value + CALL EOL_DR_ERROR(INB_BELT_MIN) +END IF ! !* 3.1 MODD_EOL_ALM variables ! at t diff --git a/src/MNH/ini_mean_field.f90 b/src/MNH/ini_mean_field.f90 index ed14412d1763587b08f421da91559639d87b75e0..146cc12c89b60381cd95060059cf8e4158bc29e7 100644 --- a/src/MNH/ini_mean_field.f90 +++ b/src/MNH/ini_mean_field.f90 @@ -52,6 +52,7 @@ END MODULE MODI_INI_MEAN_FIELD !! 05/2021 (PA.Joulin) add wind turbine variables !! 12/2021 (R. Schoetter) adds humidity and other mean diagnostics !! 11/2022 (E. Jezequel) add covariances +!! 09/2022 (H.Toumi) add EOL/ADR variables !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -62,8 +63,9 @@ USE MODD_MEAN_FIELD_n USE MODD_MEAN_FIELD USE MODD_PARAM_n USE MODD_EOL_MAIN, ONLY: LMAIN_EOL, CMETH_EOL, NMODEL_EOL -USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM -USE MODD_EOL_ALM +USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM, XAOA_SUM +USE MODD_EOL_ADR, ONLY: XFAERO_RA_SUM, XAOA_BLEQ_SUM, XFAERO_BLEQ_RA_SUM +USE MODD_EOL_ALM, ONLY: XFAERO_RE_SUM USE MODE_MODELN_HANDLER ! IMPLICIT NONE @@ -113,6 +115,14 @@ IF (LMAIN_EOL .AND. IMI==NMODEL_EOL) THEN SELECT CASE(CMETH_EOL) CASE('ADNR') ! Actuator Disc Non-Rotating XTHRU_SUM = 0.0 + CASE('ADR') ! Actuator Disc Method + XAOA_SUM = 0.0 + XFAERO_RA_SUM = 0.0 + XTHRU_SUM = 0.0 + XTORQ_SUM = 0.0 + XPOW_SUM = 0.0 + XAOA_BLEQ_SUM = 0.0 + XFAERO_BLEQ_RA_SUM = 0.0 CASE('ALM') ! Actuator Line Method XAOA_SUM = 0.0 XFAERO_RE_SUM = 0.0 diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90 index 4d22e05cdc68a9d90f1d8039051bcdf706448c09..ecaae1f6ef7aec7571272cdbd978fa948f962b43 100644 --- a/src/MNH/ini_modeln.f90 +++ b/src/MNH/ini_modeln.f90 @@ -290,12 +290,14 @@ END MODULE MODI_INI_MODEL_n ! P. Wautelet 13/09/2019: budget: simplify and modernize date/time management ! C. Lac 11/2019: correction in the drag formula and application to building in addition to tree ! S. Riette 04/2020: XHL* fields +! PA. Joulin 10/2020: add wind tubrines (EOL) ! F. Auguste 02/2021: add IBM ! T.Nigel 02/2021: add turbulence recycling ! J.L.Redelsperger 06/2011: OCEAN case ! R. Schoetter 12/2021: multi-level coupling between MesoNH and SURFEX ! R. Schoetter 12/2021: adds humidity and other mean diagnostics ! 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 !--------------------------------------------------------------------------------- ! @@ -446,6 +448,7 @@ USE MODI_INI_DRAG USE MODI_INI_DYNAMICS USE MODI_INI_ELEC_n USE MODI_INI_EOL_ADNR +USE MODI_INI_EOL_ADR USE MODI_INI_EOL_ALM USE MODI_INI_LES_N USE MODI_INI_LG @@ -2909,8 +2912,10 @@ IF (LMAIN_EOL .AND. KMI == NMODEL_EOL) THEN SELECT CASE(CMETH_EOL) CASE('ADNR') CALL INI_EOL_ADNR + CASE('ADR') + CALL INI_EOL_ADR(XDXX,XDYY,XDZZ) CASE('ALM') - CALL INI_EOL_ALM(XDXX,XDYY) + CALL INI_EOL_ALM(XDXX,XDYY,XDZZ) END SELECT END IF ! diff --git a/src/MNH/mean_field.f90 b/src/MNH/mean_field.f90 index e19ea9ed0f1ed6a6680652738eb721813259eabb..573046c23188f3d2f7edd9804864543ceebb933b 100644 --- a/src/MNH/mean_field.f90 +++ b/src/MNH/mean_field.f90 @@ -46,11 +46,12 @@ END MODULE MODI_MEAN_FIELD !! !! MODIFICATIONS !! ------------- -!! Original 07/2009 -!! (C.Lac) 09/2016 Max values -!! (PA.Joulin) 12/2020 Wind turbine variables +!! Original 07/2009 +!! (C.Lac) 09/2016 Max values +!! (PA.Joulin) 12/2020 Wind turbine variables !! (R. Schoetter) 12/2021 adds humidity and other mean diagnostics -!! (E. Jezequel) 11/2022 Welford algorithm and covariances +!! (E. Jezequel) 11/2022 Welford algorithm and covariances +!! (H. Toumi) 09/2022: add ADR !!--------------------------------------------------------------- ! ! @@ -67,11 +68,15 @@ USE MODD_PASPOL USE MODE_THERMO USE MODI_SHUMAN ! +! * EOL USE MODD_EOL_MAIN, ONLY: LMAIN_EOL, CMETH_EOL, NMODEL_EOL -USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT -USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM -USE MODD_EOL_ALM -USE MODD_EOL_ADNR +USE MODD_EOL_SHARED_IO, ONLY: XTHRUT, XTORQT, XPOWT, XAOA_GLB +USE MODD_EOL_SHARED_IO, ONLY: XTHRU_SUM, XTORQ_SUM, XPOW_SUM, XAOA_SUM +USE MODD_EOL_ADR, ONLY: XFAERO_RA_GLB, XFAERO_RA_SUM +USE MODD_EOL_ADR, ONLY: XFAERO_BLEQ_RA_GLB, XFAERO_BLEQ_RA_SUM +USE MODD_EOL_ADR, ONLY: XAOA_BLEQ_GLB, XAOA_BLEQ_SUM +USE MODD_EOL_ALM, ONLY: XFAERO_RE_SUM, XFAERO_RE_GLB +! USE MODE_MODELN_HANDLER USE MODI_UPDATE_WELFORD ! @@ -199,6 +204,14 @@ END IF SELECT CASE(CMETH_EOL) CASE('ADNR') ! Actuator Disc Non-Rotating XTHRU_SUM = XTHRUT + XTHRU_SUM + CASE('ADR') ! Actuator Disc with Rotation + XAOA_SUM = XAOA_GLB + XAOA_SUM + XFAERO_RA_SUM = XFAERO_RA_GLB + XFAERO_RA_SUM + XTHRU_SUM = XTHRUT + XTHRU_SUM + XTORQ_SUM = XTORQT + XTORQ_SUM + XPOW_SUM = XPOWT + XPOW_SUM + XAOA_BLEQ_SUM = XAOA_BLEQ_GLB + XAOA_BLEQ_SUM + XFAERO_BLEQ_RA_SUM = XFAERO_BLEQ_RA_GLB + XFAERO_BLEQ_RA_SUM CASE('ALM') ! Actuator Line Method XAOA_SUM = XAOA_GLB + XAOA_SUM XFAERO_RE_SUM = XFAERO_RE_GLB + XFAERO_RE_SUM diff --git a/src/MNH/modd_eol_adr.f90 b/src/MNH/modd_eol_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..9dcfb46b37f282d6764b7d46f8842c908dd632d0 --- /dev/null +++ b/src/MNH/modd_eol_adr.f90 @@ -0,0 +1,129 @@ +!MNH_LIC Copyright 2017-2021 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. +!----------------------------------------------------------------- +!! +!! ##################### + MODULE MODD_EOL_ADR +!! ##################### +!! +!!*** *MODD_EOL_ADR* +!! +!! PURPOSE +!! ------- +!! It is possible to include wind turbines parameterization in Meso-NH, +!! and several models are available. One of the models is the Actuator +!! Disc with Rotation (ADR). MODD_EOL_ADR contains all the declarations +!! for the ADR model. +!! +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +!! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!----------------------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ----------------- +! +IMPLICIT NONE +! +! ------ TYPES ------ +! +TYPE FARM + INTEGER :: NNB_TURBINES ! Number of wind turbines [-] + REAL, DIMENSION(:), ALLOCATABLE :: XPOS_X ! Tower base position, X coord [m] + REAL, DIMENSION(:), ALLOCATABLE :: XPOS_Y ! Tower base position, Y coord [m] + REAL, DIMENSION(:), ALLOCATABLE :: XOMEGA ! Rotor rotation speed [rad/s] + REAL, DIMENSION(:), ALLOCATABLE :: XNAC_YAW ! Nacelle yaw angle [rad] + REAL, DIMENSION(:), ALLOCATABLE :: XBLA_PITCH ! Blade picth angle [rad] +END TYPE FARM +! +TYPE TURBINE + CHARACTER(LEN=10) :: CNAME ! Wind turbine name [-] + INTEGER :: NNB_BLADES ! Number of blades [-] + REAL :: XH_HEIGHT ! Hub height [m] + REAL :: XH_DEPORT ! Hub deport [m] + REAL :: XNAC_TILT ! Tilt of the nacelle [m] + REAL :: XR_MIN ! Minimum blade radius [m] + REAL :: XR_MAX ! Maximum blade radius [m] +END TYPE TURBINE +! +TYPE BLADE + INTEGER :: NNB_RADELT ! Number of radial elements + INTEGER :: NNB_AZIELT ! Number of azimutal elements + INTEGER :: NNB_BLADAT ! Number of blade data + REAL, DIMENSION(:), ALLOCATABLE :: XRAD ! Data node radius [m] + REAL, DIMENSION(:), ALLOCATABLE :: XCHORD ! Element chord [m] + REAL, DIMENSION(:), ALLOCATABLE :: XTWIST ! Element twist [rad] + CHARACTER(LEN=20), DIMENSION(:), ALLOCATABLE :: CAIRFOIL ! Arifoil name [-] +END TYPE BLADE +! +TYPE AIRFOIL + CHARACTER(LEN=15) :: CNAME ! Airfoil name [-] + REAL, DIMENSION(:), ALLOCATABLE :: XAA ! Attack Angle [rad] + REAL, DIMENSION(:), ALLOCATABLE :: XRE ! Reynolds Number [-] + REAL, DIMENSION(:), ALLOCATABLE :: XCL ! Lift coefficient [-] + REAL, DIMENSION(:), ALLOCATABLE :: XCD ! Drag coefficient [-] + REAL, DIMENSION(:), ALLOCATABLE :: XCM ! Moment coefficient [-] +END TYPE AIRFOIL +! +! ------ VARIABLES ------ +! --- Farm data --- +TYPE(FARM) :: TFARM +TYPE(TURBINE) :: TTURBINE +TYPE(BLADE) :: TBLADE +TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE :: TAIRFOIL +! +! +! --- Global variables (Code & CPU) --- +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_AZI ! Elements azimut [rad] +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RA_GLB ! Aerodyn. force (lift+drag) in RA [N], global +! +! Blade equivalent values +REAL, DIMENSION(:,:), ALLOCATABLE :: XAOA_BLEQ_GLB ! Blade Eq. AoA +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFAERO_BLEQ_RA_GLB ! Blade Eq Aerodyn. force (lift+drag) in RA [N] +! +! Mean values +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RA_SUM ! Sum of aerodyn. force (lift+drag) in RA [N] +REAL, DIMENSION(:,:), ALLOCATABLE :: XAOA_BLEQ_SUM ! Sum of Blade Eq. AoA +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFAERO_BLEQ_RA_SUM ! Sum of Blade Eq Aero. force (lift+drag) in RA [N] + +! +! Implicit from MODD_EOL_SHARED_IO : +! --- Thruts torque and power --- +!REAL, DIMENSION(:), ALLOCATABLE :: XTHRUT ! Thrust [N] +!REAL, DIMENSION(:), ALLOCATABLE :: XTORQT ! Torque [Nm] +!REAL, DIMENSION(:), ALLOCATABLE :: XPOWT ! Power [W] +!REAL, DIMENSION(:), ALLOCATABLE :: XTHRU_SUM ! Sum of thrust (N) +!REAL, DIMENSION(:), ALLOCATABLE :: XTORQ_SUM ! Sum of torque (Nm) +!REAL, DIMENSION(:), ALLOCATABLE :: XPOW_SUM ! Sum of power (W) +! +! --- Global variables (Code & CPU) --- +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Elements radius [m] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] +!REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] +! +! Mean values +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_SUM ! Sum of angle of attack [rad] +! +! --- Namelist NAM_EOL_ADR --- +! Implicit from MODD_EOL_SHARED_IO : +!CHARACTER(LEN=100) :: CFARM_CSVDATA ! Farm file to read +!CHARACTER(LEN=100) :: CTURBINE_CSVDATA ! Turbine file to read +!CHARACTER(LEN=100) :: CBLADE_CSVDATA ! Blade file to read +!CHARACTER(LEN=100) :: CAIRFOIL_CSVDATA ! Airfoil file to read +!CHARACTER(LEN=3) :: CINTERP ! Interpolation method for wind speed +! +INTEGER :: NNB_RADELT ! Number of radial elements +INTEGER :: NNB_AZIELT ! Number of azimutal elements +! +END MODULE MODD_EOL_ADR diff --git a/src/MNH/modd_eol_alm.f90 b/src/MNH/modd_eol_alm.f90 index c04f0cb538ec0d35b16b0cb2ceb5a005401dfa81..5118bc1e3059de87e9459a9628c62372cebe00ec 100644 --- a/src/MNH/modd_eol_alm.f90 +++ b/src/MNH/modd_eol_alm.f90 @@ -26,6 +26,7 @@ !! ------------- !! Original 04/01/17 !! Modification 14/10/20 (PA. Joulin) Updated for a main version +!! Modification 04/23 (PA. Joulin) variables moved in modd_eol_shared_io !! !----------------------------------------------------------------------------- ! @@ -81,21 +82,23 @@ TYPE(BLADE) :: TBLADE TYPE(AIRFOIL), DIMENSION(:), ALLOCATABLE :: TAIRFOIL ! ! --- Global variables (Code & CPU) --- -REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Blade elements radius [m] -REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] -REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] -REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RE_GLB ! Aerodyn. force (lift+drag) in RE [N] -REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] ! ! Mean values -REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_SUM ! Sum of angle of attack [rad] REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RE_SUM ! Sum of aerodyn. force (lift+drag) in RE [N] ! ! Implicit from MODD_EOL_SHARED_IO : +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Blade elements radius [m] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] +!REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] +! !REAL, DIMENSION(:), ALLOCATABLE :: XTHRUT ! Thrust [N] !REAL, DIMENSION(:), ALLOCATABLE :: XTORQT ! Torque [Nm] !REAL, DIMENSION(:), ALLOCATABLE :: XPOWT ! Power [W] +! +!REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_SUM ! Sum of angle of attack [rad] !REAL, DIMENSION(:), ALLOCATABLE :: XTHRU_SUM ! Sum of thrust (N) !REAL, DIMENSION(:), ALLOCATABLE :: XTORQ_SUM ! Sum of torque (Nm) !REAL, DIMENSION(:), ALLOCATABLE :: XPOW_SUM ! Sum of power (W) @@ -108,11 +111,12 @@ REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RE_SUM ! Sum of aerodyn. force !CHARACTER(LEN=100) :: CBLADE_CSVDATA ! Blade file to read !CHARACTER(LEN=100) :: CAIRFOIL_CSVDATA ! Airfoil file to read !CHARACTER(LEN=3) :: CINTERP ! Interpolation method for wind speed +!LOGICAL :: LTIPLOSSG ! Flag to apply Glauert's tip loss correction +!LOGICAL :: LTECOUTPTS ! Flag to get Tecplot file output of element points +!LOGICAL :: LCSVOUTFRM ! Flag to get CSV files output of frames ! INTEGER :: NNB_BLAELT ! Number of blade elements ! LOGICAL :: LTIMESPLIT ! Flag to apply Time splitting method -LOGICAL :: LTIPLOSSG ! Flag to apply Glauert's tip loss correction -LOGICAL :: LTECOUTPTS ! Flag to get Tecplot file output of element points ! END MODULE MODD_EOL_ALM diff --git a/src/MNH/modd_eol_kine_adr.f90 b/src/MNH/modd_eol_kine_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..e10f891b82ce50b46ef5aacf4a44a7aeedfba161 --- /dev/null +++ b/src/MNH/modd_eol_kine_adr.f90 @@ -0,0 +1,112 @@ +!MNH_LIC Copyright 2018-2021 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. +!----------------------------------------------------------------- +!! +!! ##################### + MODULE MODD_EOL_KINE_ADR +!! ##################### +!! +!!*** *MODD_EOL_KINE_ADR* +!! +!! PURPOSE +!! ------- +!! Declaration to take into account wind turbine motion +!! +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +!! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!----------------------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ----------------- +! +! +! - Matrix to move from one fram to an other - +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XMAT_RG_RT ! RG = Repere GLOBAL, RT = Repere TOWER +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XMAT_RG_RN, XMAT_RT_RN ! RN = Repere NACELLE +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XMAT_RG_RH, XMAT_RH_RG, XMAT_RN_RH ! RH = Repere HUB +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XMAT_RG_RC, XMAT_RH_RC ! RC = Repere CYLINDRICAL +REAL, DIMENSION(:,:,:,:,:), ALLOCATABLE :: XMAT_RG_RTA, XMAT_RTA_RG, XMAT_RC_RTA ! RTA = Repere TRANSITION ANNULAR +REAL, DIMENSION(:,:,:,:,:), ALLOCATABLE :: XMAT_RG_RA, XMAT_RA_RG, XMAT_RTA_RA ! RE = Repere ANNULAR + +! - POSITIONS & ORIENTATIONS - +REAL, DIMENSION(3) :: XPOS_REF ! Reference position + +! Tower +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOSINI_TOWO_RG ! Initial tower origin position, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOS_TOWO_RG ! Current tower origin real position, in global reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XPOS_TELT_RG ! Current tower element position, in global reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XPOS_TELT_RT ! Current tower element position, in tower frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XANGINI_TOW_RG ! Initial tower orientation in global ref frame + +! Nacelle +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOSINI_NACO_RT ! Initial nacelle position, in tower reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOS_NACO_RG ! Initial nacelle position, in global reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XPOS_NELT_RG ! Initial nacelle position, in global reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XPOS_NELT_RN ! Initial nacelle position, in nacelle reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XANGINI_NAC_RT ! Initial nacelle orientation, in tower reference frame + +! Hub +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOSINI_HUB_RN ! Initial hub position, in nacelle reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOS_HUB_RG ! Current hub position, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XANGINI_HUB_RN ! Initial hub orientation, in nacelle reference frame + +! Cylindrical +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOSINI_CYL_RH ! Initial blade root position, in hub reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XPOS_CYL_RG ! Current blade root position, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XANGINI_CYL_RH ! Initial blade orientation, in hub reference frame + +! Annular Element +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XPOS_ELT_RC ! Element position, in cylindrical reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XPOS_ELT_RG ! Element position, in global reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XPOS_SEC_RC ! Section position, in cylindrical reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XPOS_SEC_RG ! Section position, in global reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XANGINI_ELT_RC ! Initial element orientation in cyl reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XANGINI_ELT_RA ! Initial element orientation in annular elt transition ref frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XTWIST_ELT ! Element twist, interpolated from data +REAL, DIMENSION(:,:), ALLOCATABLE :: XCHORD_ELT ! Element chord lenght, interpolated from data +REAL, DIMENSION(:,:), ALLOCATABLE :: XSURF_ELT ! Element lift surface +REAL, DIMENSION(:,:), ALLOCATABLE :: XSURF_APP_ELT ! Element lift surface +! +! +! - STRUCTURAL VELOCITIES - +! Tower +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_TOWO_RG ! Tower base translation velocity, in global reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XTVEL_TELT_RG ! Tower element velocity, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RT_RG ! RT/RG rotational velocity + +! Nacelle +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_NACO_RT ! Nacelle base translation velocity, in tower reference frame +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XTVEL_NELT_RG ! Nacelle element translation velocity, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RN_RT ! RN/RT rotational velocity +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RN_RG ! RN/RG rotational velocity + +! Hub +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_HUB_RN ! Hub base translation velocity, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_HUB_RG ! Hub base translation velocity, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RH_RN ! RH/RN rotational velocity +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RH_RG ! RH/RG rotational velocity + +! Cylindrical +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_CYL_RH ! Cyllindrical base translation velocity, in hub reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XTVEL_CYL_RG ! Cylindrical base translation velocity, in global reference frame +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RC_RH ! RC/RH rotational velocity +REAL, DIMENSION(:,:), ALLOCATABLE :: XRVEL_RC_RG ! RC/RG rotational velocity + +! Annular Elements +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XTVEL_ELT_RC ! Element base translation velocity, in cylindrical reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XTVEL_ELT_RG ! Element base translation velocity, in global reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XTVEL_ELT_RA ! Element base translation velocity, in annular element reference frame +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XRVEL_RA_RC ! RA/RC rotational velocity +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XRVEL_RA_RG ! RA/RG rotational velocity + +END MODULE MODD_EOL_KINE_ADR diff --git a/src/MNH/modd_eol_shared_io.f90 b/src/MNH/modd_eol_shared_io.f90 index 6132fb28f19890e2f30b56539b55735b70f26be6..a80dd1fa8e84de29de2dacd84e4df54410755b2f 100644 --- a/src/MNH/modd_eol_shared_io.f90 +++ b/src/MNH/modd_eol_shared_io.f90 @@ -26,6 +26,7 @@ !! MODIFICATIONS !! ------------- !! Original 17/11/20 +!! Modification 09/22 (H. Toumi) ADR integration : new shared var !! !----------------------------------------------------------------------------- ! @@ -42,17 +43,30 @@ CHARACTER(LEN=100) :: CBLADE_CSVDATA ! Blade file to read CHARACTER(LEN=100) :: CAIRFOIL_CSVDATA ! Airfoil file to read ! * flags CHARACTER(LEN=3) :: CINTERP ! Interpolation method for wind speed +LOGICAL :: LTIPLOSSG ! Flag to apply Glauert's tip loss correction +LOGICAL :: LTECOUTPTS ! Flag to get Tecplot file output of element points +LOGICAL :: LCSVOUTFRM ! Flag to get CSV files output of frames ! ! !* 2. OUTPUTS VAR ! ----------------- ! +! --- Global variables of aerodynamic models (Code & CPU) --- +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XELT_RAD ! Elements radius [m] +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_GLB ! Angle of attack of an element [rad] +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFLIFT_GLB ! Lift force, parallel to Urel [N] +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XFDRAG_GLB ! Drag force, perpendicular to Urel [N] +REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: XFAERO_RG_GLB ! Aerodyn. force (lift+drag) in RG [N] +! ! --- Thruts torque and power --- REAL, DIMENSION(:), ALLOCATABLE :: XTHRUT ! Thrust [N] REAL, DIMENSION(:), ALLOCATABLE :: XTORQT ! Torque [Nm] REAL, DIMENSION(:), ALLOCATABLE :: XPOWT ! Power [W] +! +! Mean values REAL, DIMENSION(:), ALLOCATABLE :: XTHRU_SUM ! Sum of thrust (N) REAL, DIMENSION(:), ALLOCATABLE :: XTORQ_SUM ! Sum of torque (Nm) REAL, DIMENSION(:), ALLOCATABLE :: XPOW_SUM ! Sum of power (W) +REAL, DIMENSION(:,:,:), ALLOCATABLE :: XAOA_SUM ! Sum of angle of attack [rad] ! END MODULE MODD_EOL_SHARED_IO diff --git a/src/MNH/modn_eol_adr.f90 b/src/MNH/modn_eol_adr.f90 new file mode 100644 index 0000000000000000000000000000000000000000..db55c0cfcec70db6f3bd1b060fc5cc2a5df1636b --- /dev/null +++ b/src/MNH/modn_eol_adr.f90 @@ -0,0 +1,49 @@ +!MNH_LIC Copyright 2017-2021 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. +!----------------------------------------------------------------- +!! +!! ##################### + MODULE MODN_EOL_ADR +!! ##################### +!! +!!*** *MODN_EOL_ADR* +!! +!! PURPOSE +!! ------- +!! NAM_EOL activate the parameterization of wind turbines, and several +!! models are available. One of the models is the Actuator Disc with +!! Rotation (ADR). +!! The aim of NAM_EOL_ADR is to specify ADR parameters. +!! +!!** AUTHOR +!! ------ +!! H. Toumi *IFPEN* +! +!! MODIFICATIONS +!! ------------- +!! Original 09/22 +!! Modification 05/04/23 (PA. Joulin) Updated for a main version +!! +!! IMPLICIT ARGUMENTS +!! ------------------ +USE MODD_EOL_ADR , ONLY : NNB_AZIELT, NNB_RADELT +USE MODD_EOL_SHARED_IO, ONLY : CFARM_CSVDATA,CTURBINE_CSVDATA +USE MODD_EOL_SHARED_IO, ONLY : CBLADE_CSVDATA, CAIRFOIL_CSVDATA +USE MODD_EOL_SHARED_IO, ONLY : CINTERP, LTIPLOSSG +USE MODD_EOL_SHARED_IO, ONLY : LTECOUTPTS, LCSVOUTFRM +!! +!----------------------------------------------------------------------------- +! +!* 0. DECLARATIONS +! ----------------- +IMPLICIT NONE +SAVE +NAMELIST /NAM_EOL_ADR/ & + CFARM_CSVDATA, CTURBINE_CSVDATA, CBLADE_CSVDATA, CAIRFOIL_CSVDATA, & + NNB_AZIELT, NNB_RADELT, & + CINTERP, LTIPLOSSG, & + LTECOUTPTS,LCSVOUTFRM +! +END MODULE MODN_EOL_ADR diff --git a/src/MNH/modn_eol_alm.f90 b/src/MNH/modn_eol_alm.f90 index 54f7a389bda529ce0ea80fac0fd6ff62124d4b18..473d357e1486e1acb07bc77fa645cd39a37722fa 100644 --- a/src/MNH/modn_eol_alm.f90 +++ b/src/MNH/modn_eol_alm.f90 @@ -25,6 +25,7 @@ !! ------------- !! Original 24/01/17 !! Modification 14/10/20 (PA. Joulin) Updated for a main version +!! Modification 04/23 (PA. Joulin) Add csv outputs of frames !! !! IMPLICIT ARGUMENTS !! ------------------ @@ -41,6 +42,6 @@ NAMELIST /NAM_EOL_ALM/ & CFARM_CSVDATA, CTURBINE_CSVDATA, CBLADE_CSVDATA, CAIRFOIL_CSVDATA, & NNB_BLAELT, & CINTERP, LTIMESPLIT, LTIPLOSSG, & - LTECOUTPTS + LTECOUTPTS,LCSVOUTFRM ! END MODULE MODN_EOL_ALM diff --git a/src/MNH/read_exsegn.f90 b/src/MNH/read_exsegn.f90 index a5c44aee26580ff1bcd278192a9f83db8f8c5991..cab38bf983a4f1c15ac01c52f5685dfdfb2b158b 100644 --- a/src/MNH/read_exsegn.f90 +++ b/src/MNH/read_exsegn.f90 @@ -295,6 +295,7 @@ END MODULE MODI_READ_EXSEG_n ! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg ! C. Lac 11/2019: correction in the drag formula and application to building in addition to tree ! Q. Rodier 03/2020: add abort if use of any LHORELAX and cyclic conditions +! PA. Joulin 10/2020: add namelists for wind tubrines (EOL) ! F.Auguste 02/2021: add IBM ! T.Nagel 02/2021: add turbulence recycling ! E.Jezequel 02/2021: add stations read from CSV file @@ -309,6 +310,7 @@ END MODULE MODI_READ_EXSEG_n ! P. Wautelet 24/06/2022: remove check on CSTORAGE_TYPE for restart of ForeFire variables ! P. Wautelet 13/07/2022: add namelist for flyers and balloons ! 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 !------------------------------------------------------------------------------ ! @@ -383,6 +385,7 @@ USE MODN_DYN_n ! to avoid the duplication of this routine for each model. USE MODN_ELEC USE MODN_EOL USE MODN_EOL_ADNR +USE MODN_EOL_ADR USE MODN_EOL_ALM USE MODN_FIRE_n USE MODN_FLYERS @@ -570,6 +573,8 @@ CALL POSNAM( TPEXSEGFILE, 'NAM_EOL', GFOUND ) IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_EOL) CALL POSNAM( TPEXSEGFILE, 'NAM_EOL_ADNR', GFOUND ) IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_EOL_ADNR) +CALL POSNAM( TPEXSEGFILE, 'NAM_EOL_ADR', GFOUND ) +IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_EOL_ADR) CALL POSNAM( TPEXSEGFILE, 'NAM_EOL_ALM', GFOUND ) IF (GFOUND) READ(UNIT=ILUSEG,NML=NAM_EOL_ALM) CALL POSNAM( TPEXSEGFILE, 'NAM_PROFILERN', GFOUND ) diff --git a/src/MNH/write_lfin.f90 b/src/MNH/write_lfin.f90 index 9e99477d3cb000f97fb64c031bf0ec03982a0032..e663656a9bf48a9989a6a450ba908151d4d67755 100644 --- a/src/MNH/write_lfin.f90 +++ b/src/MNH/write_lfin.f90 @@ -182,6 +182,8 @@ END MODULE MODI_WRITE_LFIFM_n ! A. Costes 12/2021: add Blaze fire model ! P. Wautelet 04/02/2022: use TSVLIST to manage metadata of scalar variables ! E. Jezequel 11/2022: add covariances from MEAN fields +! H. Toumi 09/2022: add ADR +! PA. Joulin 04/2023: update EOL metadata management !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS @@ -214,6 +216,7 @@ USE MODD_ELEC_DESCR, ONLY: LLNOX_EXPLICIT USE MODD_ELEC_FLASH USE MODD_ELEC_n USE MODD_EOL_ADNR +USE MODD_EOL_ADR USE MODD_EOL_ALM USE MODD_EOL_MAIN USE MODD_EOL_SHARED_IO @@ -2067,7 +2070,7 @@ END IF ! IF (LMAIN_EOL .AND. IMI == NMODEL_EOL) THEN TZFIELD = TFIELDMETADATA( & - CMNHNAME = 'generic for wind turbine variables', & !Temporary name to ease identification + CMNHNAME = 'generic for wind turbine var', & !Temporary name to ease identification CSTDNAME = '', & CUNITS = 'N', & CDIR = 'XY', & @@ -2149,18 +2152,20 @@ SELECT CASE(CMETH_EOL) CALL IO_Field_write(TPFILE,TZFIELD,XTHRU_SUM/MEAN_COUNT) ! END IF -! iii) Actuator Line Model ! - CASE('ALM') ! Actuator Line Method +! iii) Actuator Disc with Rotation Model +! + CASE('ADR') ! Actuator Disc with Rotation ! +! * 1D Variables (rotor id) TZFIELD = TFIELDMETADATA( & - CMNHNAME = 'generic for ALM variables', & !Temporary name to ease identification + CMNHNAME = '1D ADR var: (rot)', & !Temporary name to ease identification CSTDNAME = '', & CDIR = '--', & NGRID = 1, & NTYPE = TYPEREAL, & NDIMS = 1, & - LTIMEDEP = .TRUE. ) + LTIMEDEP = .TRUE. ) ! TZFIELD%CMNHNAME = 'THRUT' TZFIELD%CLONGNAME = 'THRUSTT_EOL' @@ -2180,14 +2185,222 @@ SELECT CASE(CMETH_EOL) TZFIELD%CCOMMENT = 'RID instantaneous power (W) of wind turbines' CALL IO_Field_write(TPFILE,TZFIELD,XPOWT) ! +! +! * 3D Variables (rotor id, azimuthal id, radial id) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = '3D ADR var: (rot,azi,rad)', & + CSTDNAME = '', & + CDIR = '--', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 3, & + LTIMEDEP = .TRUE. ) +! + TZFIELD%CMNHNAME = 'ELT_RAD' + TZFIELD%CLONGNAME = 'ELT_RAD' + TZFIELD%CUNITS = 'm' + TZFIELD%CCOMMENT = 'RID_AZI_RAD radius (m) of wind turbine blade elements' + CALL IO_Field_write(TPFILE,TZFIELD,XELT_RAD) +! + TZFIELD%CMNHNAME = 'ELT_AZI' + TZFIELD%CLONGNAME = 'ELT_AZI' + TZFIELD%CUNITS = 'rad' + TZFIELD%CCOMMENT = 'RID_AZI_RAD Azimutal angle (rad) of wind turbine blade elements' + CALL IO_Field_write(TPFILE,TZFIELD,XELT_AZI) +! + TZFIELD%CMNHNAME = 'AOA' + TZFIELD%CLONGNAME = 'ANGLE OF ATTACK' + TZFIELD%CUNITS = 'rad' + TZFIELD%CCOMMENT = 'RID_AZI_RAD instantaneous angle of attack (rad)' + CALL IO_Field_write(TPFILE,TZFIELD,XAOA_GLB) +! + TZFIELD%CMNHNAME = 'FLIFT' + TZFIELD%CLONGNAME = 'LIFT FORCE' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID_AZI_RAD instantaneous lift (N) in relative frame' + CALL IO_Field_write(TPFILE,TZFIELD,XFLIFT_GLB) +! + TZFIELD%CMNHNAME = 'FDRAG' + TZFIELD%CLONGNAME = 'DRAG FORCE' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID_AZI_RAD instantaneous drag (N) in relative frame' + CALL IO_Field_write(TPFILE,TZFIELD,XFDRAG_GLB) +! +! * 4D Variables (rotor id, azimuthal id, radial id, xyz) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = '4D ADR var: (rot,azi,rad,xyz)', & + CSTDNAME = '', & + CDIR = '--', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 4, & + LTIMEDEP = .TRUE. ) +! + TZFIELD%CMNHNAME = 'FAERO_RA' + TZFIELD%CLONGNAME = 'AERODYNAMIC FORCE RA' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID_AZI_RAD_XYZ instantaneous forces (N) in RA' + CALL IO_Field_write(TPFILE,TZFIELD,XFAERO_RA_GLB) +! + TZFIELD%CMNHNAME = 'FAERO_RG' + TZFIELD%CLONGNAME = 'AERODYNAMIC FORCE RG' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID_AZI_RAD_XYZ instantaneous forces (N) in RG' + CALL IO_Field_write(TPFILE,TZFIELD,XFAERO_RG_GLB) +! +! * Blade Equivalent Variables (rotor id, radial id, xyz) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'AOA_BLEQ', & + CSTDNAME = '', & + CLONGNAME = 'BLADE_EQ_ANGLE_OF_ATTACK', & + CUNITS = 'rad', & + CDIR = '--', & + CCOMMENT = 'RID_RAD blade eq. angle of attack (rad)', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 2, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XAOA_BLEQ_GLB) + + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'FAERO_BLEQ_RA', & + CSTDNAME = '', & + CLONGNAME = 'BLADE_EQ_AERO_FORCE_RA', & + CUNITS = 'N', & + CDIR = '--', & + CCOMMENT = 'RID_RAD_XYZ blade eq. forces (N) in RA', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 3, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XFAERO_BLEQ_RA_GLB) +! + IF (MEAN_COUNT /= 0) THEN +! +! * 1D Variables (rotor id) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = '1D ADR mean var: (rot)', & !Temporary name to ease identification + CSTDNAME = '', & + CDIR = '--', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 1, & + LTIMEDEP = .TRUE. ) +! + TZFIELD%CMNHNAME = 'THRUMME' + TZFIELD%CLONGNAME = 'MEAN_THRUST_EOL' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID mean thrust of the wind turbines (N)' + CALL IO_Field_write(TPFILE,TZFIELD,XTHRU_SUM/MEAN_COUNT) +! + TZFIELD%CMNHNAME = 'TORQMME' + TZFIELD%CLONGNAME = 'MEAN_TORQUE_EOL' + TZFIELD%CUNITS = 'Nm' + TZFIELD%CCOMMENT = 'RID mean torque of the wind turbines (Nm)' + CALL IO_Field_write(TPFILE,TZFIELD,XTORQ_SUM/MEAN_COUNT) +! + TZFIELD%CMNHNAME = 'POWMME' + TZFIELD%CLONGNAME = 'MEAN_POWER_EOL' + TZFIELD%CUNITS = 'W' + TZFIELD%CCOMMENT = 'RID mean power of the wind turbines (W)' + CALL IO_Field_write(TPFILE,TZFIELD,XPOW_SUM/MEAN_COUNT) +! + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'AOAMME', & + CSTDNAME = '', & + CLONGNAME = 'MEAN_ANGLE_OF_ATTACK', & + CUNITS = 'rad', & + CDIR = '--', & + CCOMMENT = 'RID_AZI_RAD mean angle of attack (rad)', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 3, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XAOA_SUM/MEAN_COUNT) +! + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'FAEROMME_RA', & + CSTDNAME = '', & + CLONGNAME = 'MEAN_AERODYNAMIC_FORCE_RA', & + CUNITS = 'N', & + CDIR = '--', & + CCOMMENT = 'RID_AZI_RAD_XYZ mean forces (N) in RA', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 4, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XFAERO_RA_SUM/MEAN_COUNT) +! +! * Blade Equivalent Variables (rotor id, radial id, xyz) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'AOAMME_BLEQ', & + CSTDNAME = '', & + CLONGNAME = 'MEAN_BLADE_EQ_AOA', & + CUNITS = 'rad', & + CDIR = '--', & + CCOMMENT = 'RID_RAD mean blade eq. angle of attack (rad)', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 2, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XAOA_BLEQ_SUM/MEAN_COUNT) + + TZFIELD = TFIELDMETADATA( & + CMNHNAME = 'FAEROMME_BLEQ_RA', & + CSTDNAME = '', & + CLONGNAME = 'MEAN_BLADE_EQ_AERO_F_RA', & + CUNITS = 'N', & + CDIR = '--', & + CCOMMENT = 'RID_RAD_XYZ mean blade eq. forces (N) in RA', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 3, & + LTIMEDEP = .TRUE. ) + CALL IO_Field_write(TPFILE,TZFIELD,XFAERO_BLEQ_RA_SUM/MEAN_COUNT) +! + END IF +! +! iv) Actuator Line Model +! + CASE('ALM') ! Actuator Line Method +! +! * 1D Variables (rotor id) TZFIELD = TFIELDMETADATA( & - CMNHNAME = 'generic for ALM variables', & !Temporary name to ease identification + CMNHNAME = '1D ALM var: (rot)', & !Temporary name to ease identification CSTDNAME = '', & CDIR = '--', & NGRID = 1, & NTYPE = TYPEREAL, & - NDIMS = 3, & + NDIMS = 1, & LTIMEDEP = .TRUE. ) +! + TZFIELD%CMNHNAME = 'THRUT' + TZFIELD%CLONGNAME = 'THRUSTT_EOL' + TZFIELD%CUNITS = 'N' + TZFIELD%CCOMMENT = 'RID instantaneous thrust (N) of wind turbines' + CALL IO_Field_write(TPFILE,TZFIELD,XTHRUT) +! + TZFIELD%CMNHNAME = 'TORQT' + TZFIELD%CLONGNAME = 'TORQUET_EOL' + TZFIELD%CUNITS = 'Nm' + TZFIELD%CCOMMENT = 'RID instantaneous torque (Nm) of wind turbines' + CALL IO_Field_write(TPFILE,TZFIELD,XTORQT) +! + TZFIELD%CMNHNAME = 'POWT' + TZFIELD%CLONGNAME = 'POWERT_EOL' + TZFIELD%CUNITS = 'W' + TZFIELD%CCOMMENT = 'RID instantaneous power (W) of wind turbines' + CALL IO_Field_write(TPFILE,TZFIELD,XPOWT) +! +! * 3D Variables (rotor id, blade id, radial id) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = '3D ALM var: (rot,blade,rad)', & + CSTDNAME = '', & + CDIR = '--', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 3, & + LTIMEDEP = .TRUE. ) ! TZFIELD%CMNHNAME = 'ELT_RAD' TZFIELD%CLONGNAME = 'ELT_RAD' @@ -2213,14 +2426,15 @@ SELECT CASE(CMETH_EOL) TZFIELD%CCOMMENT = 'RID_BID_EID instantaneous drag (N) in relative frame' CALL IO_Field_write(TPFILE,TZFIELD,XFDRAG_GLB) ! - TZFIELD = TFIELDMETADATA( & - CMNHNAME = 'generic for ALM variables', & !Temporary name to ease identification - CSTDNAME = '', & - CDIR = '--', & - NGRID = 1, & - NTYPE = TYPEREAL, & - NDIMS = 4, & - LTIMEDEP = .TRUE. ) +! * 4D Variables (rotor id, azimuthal id, radial id, xyz) + TZFIELD = TFIELDMETADATA( & + CMNHNAME = '4D ALM var: (rot,blade,rad,xyz)', & + CSTDNAME = '', & + CDIR = '--', & + NGRID = 1, & + NTYPE = TYPEREAL, & + NDIMS = 4, & + LTIMEDEP = .TRUE. ) ! TZFIELD%CMNHNAME = 'FAERO_RE' TZFIELD%CLONGNAME = 'AERODYNAMIC FORCE RE' @@ -2236,8 +2450,9 @@ SELECT CASE(CMETH_EOL) ! IF (MEAN_COUNT /= 0) THEN ! +! * 1D Variables (rotor id) TZFIELD = TFIELDMETADATA( & - CMNHNAME = 'generic for ALM mean variables', & !Temporary name to ease identification + CMNHNAME = '1D ALM mean var: (rot)', & !Temporary name to ease identification CSTDNAME = '', & CDIR = '--', & NGRID = 1, & diff --git a/src/configure b/src/configure index c031c125297153eb6cab33a6597aa3b99c89790b..d911ff218651940a977c9e3f96f83a30aaf3379b 100755 --- a/src/configure +++ b/src/configure @@ -257,6 +257,30 @@ export ARMCI_SHR_BUF_METHOD=COPY ;; esac ;; +'Linux topaze'*) # Topaze TGCC + export MNH_ARCH=`echo $ARCH | grep LX` + export ARCH=${MNH_ARCH:-LXifort} + export VER_MPI=${VER_MPI:-MPIAUTO} + export OPTLEVEL=${OPTLEVEL:-O3} + export MVWORK=${MVWORK:-NO} + export VER_CDF=${VER_CDF:-CDFAUTO} + export MNHENV=${MNHENV:-" +module purge +module load cmake/3.20.3 +module load inteloneapi/21.4.0 +module load mpi/openmpi/4.1.4 + +export SLURM_CPU_BIND=none +# Set some openmpi variable for pb with nb of cores >> 1024 +export OMPI_MCA_coll_hcoll_enable=0 +export HCOLL_ENABLE_MCAST_ALL=0 +export OMPI_MCA_coll_tuned_barrier_algorithm=2 +# For GA version set GA/ARMCI variables +export ARMCI_VERBOSE=1 +export ARMCI_STRIDED_METHOD=IOV ARMCI_IOV_METHOD=BATCHED +export ARMCI_SHR_BUF_METHOD=COPY +"} +;; 'Linux belenos'*|'Linux taranis'*) export ARCH=${ARCH:-LXifort} export VER_MPI=${VER_MPI:-MPIAUTO} @@ -507,7 +531,30 @@ module load intel/20.0.015 intelmpi/20.0.015 export SLURM_CPU_BIND=none export I_MPI_PIN_PROCESSOR_LIST=all:map=spread "} - ;; + ;; +'Linux ener'*) # Cluster IFPEN + export ARCH=${ARCH:-LXifort} + export VER_MPI=${VER_MPI:-MPIAUTO} + export OPTLEVEL=${OPTLEVEL:-O2} + export VER_CDF=${VER_CDF:-CDFAUTO} + export MNHENV=${MNHENV:-" +module purge +source /soft/irsrvsoft1/expl/eb/r17/$(echo el_`lsb_release -r | awk -F ":" ' {print $2}' | awk -F "." '{print $1}' | tr -d '\t'`-`uname -i`)/envs/mesonh.sh +module load intel-compilers/2021.4.0 +export OMPI_FC=ifort +"} + ;; + +'Linux irlin'*|'Linux islin'*) # PC IFPEN + export ARCH=${ARCH:-LXgfortran} + export VER_MPI=${VER_MPI:-MPIAUTO} + export OPTLEVEL=${OPTLEVEL:-O2} + export VER_CDF=${VER_CDF:-CDFAUTO} + export MNHENV=${MNHENV:-" +module purge +source /soft/irsrvsoft1/expl/eb/r17/$(echo el_`lsb_release -r | awk -F ":" ' {print $2}' | awk -F "." '{print $1}' | tr -d '\t'`-`uname -i`)/envs/mesonh.sh +"} + ;; Linux*) export ARCH=${ARCH:-LXgfortran} diff --git a/src/job_make_mesonh_dev_pc_ifpen b/src/job_make_mesonh_dev_pc_ifpen new file mode 100755 index 0000000000000000000000000000000000000000..d735db715dd27d80153515264580b0432e309900 --- /dev/null +++ b/src/job_make_mesonh_dev_pc_ifpen @@ -0,0 +1,25 @@ +#!/bin/bash +#MNH_LIC Copyright 1994-2019 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. + +pwd + +set -x + +export ARCH=LXgfortran +export VERSION_XYZ=MNH-V5-6-1 +export VER_MPI=MPIAUTO +export OPTLEVEL=DEBUG + +# Test de l'existence du profile +if ! [ -e ../conf/profile_mesonh-${ARCH}-R8I4-${VERSION_XYZ}-${VER_MPI}-${OPTLEVEL} ]; then + ./configure +fi + +# Chargement du profil (qui contient le source vers l'env ifpen-mesonh) +. ../conf/profile_mesonh-${ARCH}-R8I4-${VERSION_XYZ}-${VER_MPI}-${OPTLEVEL} + +time gmake -j 10 +time gmake -j 1 installmaster diff --git a/src/job_make_mesonh_ener440 b/src/job_make_mesonh_ener440 new file mode 100755 index 0000000000000000000000000000000000000000..1b2675aecaafe054fdd769bf86ad3681dcb7f9ad --- /dev/null +++ b/src/job_make_mesonh_ener440 @@ -0,0 +1,23 @@ +#!/bin/bash +#MNH_LIC Copyright 1994-2019 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. + +#SBATCH -N 1 +#SBATCH -n 10 +#SBATCH -o MasterI.eo%j # +#SBATCH -e MasterI.eo%j # +#SBATCH -p debug +#SBATCH -J "MesoNH-compile" +#SBATCH --wckey MPX13003 # Numero du projet +#SBATCH -t 02:00:00 + +# On va lancer la compilation dans le répertoire de lancement du job +cd ${SLURM_SUBMIT_DIR} + +# Chargement du profil (qui contient le source vers l'env ifpen-mesonh) +. ../conf/profile_mesonh-LXifort-R8I4-MNH-V5-6-1-MPIAUTO-O2 + +time gmake -j 10 +time gmake -j 1 installmaster diff --git a/src/job_make_mesonh_pc_ifpen b/src/job_make_mesonh_pc_ifpen new file mode 100755 index 0000000000000000000000000000000000000000..559c64dd8550e7fb40e6054635c7ae1cea02d90d --- /dev/null +++ b/src/job_make_mesonh_pc_ifpen @@ -0,0 +1,14 @@ +#!/bin/bash +#MNH_LIC Copyright 1994-2019 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. + +# On va lancer la compilation dans le répertoire de lancement du job +pwd + +# Chargement du profil (qui contient le source vers l'env ifpen-mesonh) +. ../conf/profile_mesonh-LXgfortran-R8I4-MNH-V5-6-1-MPIAUTO-O2 + +time gmake -j 10 +time gmake -j 1 installmaster