MODULE create_limit_unstruct_mod PRIVATE INTEGER, PARAMETER :: lmdep = 12 PUBLIC create_limit_unstruct CONTAINS SUBROUTINE create_limit_unstruct USE dimphy USE lmdz_xios USE ioipsl, ONLY: ioget_year_len USE time_phylmdz_mod, ONLY: annee_ref USE indice_sol_mod USE phys_state_var_mod USE lmdz_phys_para USE lmdz_abort_physic, ONLY: abort_physic USE lmdz_iniprint, ONLY: lunout, prt_level IMPLICIT NONE REAL, DIMENSION(:, :), ALLOCATABLE :: sic REAL, DIMENSION(:, :), ALLOCATABLE :: sst REAL, DIMENSION(klon, lmdep) :: rugos REAL, DIMENSION(klon, lmdep) :: albedo REAL, DIMENSION(:, :), ALLOCATABLE :: sic_mpi REAL, DIMENSION(:, :), ALLOCATABLE :: sst_mpi REAL, DIMENSION(klon_mpi, lmdep) :: rugos_mpi REAL, DIMENSION(klon_mpi, lmdep) :: albedo_mpi INTEGER :: ndays REAL :: fi_ice(klon) REAL, ALLOCATABLE :: sic_year(:, :) REAL, ALLOCATABLE :: sst_year(:, :) REAL, ALLOCATABLE :: rugos_year(:, :) REAL, ALLOCATABLE :: albedo_year(:, :) REAL, ALLOCATABLE :: pctsrf_t(:, :, :) REAL, ALLOCATABLE :: phy_bil(:, :) REAL, ALLOCATABLE :: sst_year_mpi(:, :) REAL, ALLOCATABLE :: rugos_year_mpi(:, :) REAL, ALLOCATABLE :: albedo_year_mpi(:, :) REAL, ALLOCATABLE :: pctsrf_t_mpi(:, :, :) REAL, ALLOCATABLE :: phy_bil_mpi(:, :) INTEGER :: l, k INTEGER :: nbad INTEGER :: sic_time_axis_size INTEGER :: sst_time_axis_size CHARACTER(LEN = 99) :: mess ! error message ndays = ioget_year_len(annee_ref) IF (is_omp_master) CALL xios_get_axis_attr("time_sic", n_glo = sic_time_axis_size) CALL bcast_omp(sic_time_axis_size) ALLOCATE(sic_mpi(klon_mpi, sic_time_axis_size)) ALLOCATE(sic(klon, sic_time_axis_size)) IF (is_omp_master) CALL xios_get_axis_attr("time_sst", n_glo = sst_time_axis_size) CALL bcast_omp(sst_time_axis_size) ALLOCATE(sst_mpi(klon_mpi, sst_time_axis_size)) ALLOCATE(sst(klon, sst_time_axis_size)) IF (is_omp_master) THEN CALL xios_recv_field("sic_limit", sic_mpi) CALL xios_recv_field("sst_limit", sst_mpi) CALL xios_recv_field("rugos_limit", rugos_mpi) CALL xios_recv_field("albedo_limit", albedo_mpi) ENDIF CALL scatter_omp(sic_mpi, sic) CALL scatter_omp(sst_mpi, sst) CALL scatter_omp(rugos_mpi, rugos) CALL scatter_omp(albedo_mpi, albedo) ALLOCATE(sic_year(klon, ndays)) ALLOCATE(sst_year(klon, ndays)) ALLOCATE(rugos_year(klon, ndays)) ALLOCATE(albedo_year(klon, ndays)) ALLOCATE(pctsrf_t(klon, nbsrf, ndays)) ALLOCATE(phy_bil(klon, ndays)); phy_bil = 0.0 ! sic IF (sic_time_axis_size==lmdep) THEN CALL time_interpolation(ndays, sic, 'gregorian', sic_year) ELSE IF (sic_time_axis_size==ndays) THEN sic_year = sic ELSE WRITE(mess, *) 'sic time axis is nor montly, nor daily. sic time interpolation ', & 'is requiered but is not currently managed' CALL abort_physic('create_limit_unstruct', TRIM(mess), 1) ENDIF sic_year(:, :) = sic_year(:, :) / 100. ! convert percent to fraction WHERE(sic_year(:, :)>1.0) sic_year(:, :) = 1.0 ! Some fractions have some time large negative values WHERE(sic_year(:, :)<0.0) sic_year(:, :) = 0.0 ! probably better to apply alse this filter before horizontal interpolation ! sst IF (sst_time_axis_size==lmdep) THEN CALL time_interpolation(ndays, sst, 'gregorian', sst_year) ELSE IF (sst_time_axis_size==ndays) THEN sst_year = sst ELSE WRITE(mess, *)'sic time axis is nor montly, nor daily. sic time interpolation ', & 'is requiered but is not currently managed' CALL abort_physic('create_limit_unstruct', TRIM(mess), 1) ENDIF WHERE(sst_year(:, :)<271.38) sst_year(:, :) = 271.38 ! rugos DO l = 1, lmdep WHERE(NINT(zmasq(:))/=1) rugos(:, l) = 0.001 ENDDO CALL time_interpolation(ndays, rugos, '360_day', rugos_year) ! albedo CALL time_interpolation(ndays, albedo, '360_day', albedo_year) DO k = 1, ndays fi_ice = sic_year(:, k) WHERE(fi_ice>=1.0) fi_ice = 1.0 WHERE(fi_ice=1.0 - zmasq) pctsrf_t(:, is_sic, k) = 1.0 - zmasq pctsrf_t(:, is_oce, k) = 0.0 ELSEWHERE pctsrf_t(:, is_oce, k) = 1.0 - zmasq - pctsrf_t(:, is_sic, k) WHERE(pctsrf_t(:, is_oce, k)0) WRITE(lunout, *) 'pb sous maille pour nb point = ', nbad nbad = COUNT(abs(sum(pctsrf_t(:, :, k), dim = 2) - 1.0)>EPSFRA) IF(nbad>0) WRITE(lunout, *) 'pb sous surface pour nb points = ', nbad END DO ALLOCATE(sst_year_mpi(klon_mpi, ndays)) ALLOCATE(rugos_year_mpi(klon_mpi, ndays)) ALLOCATE(albedo_year_mpi(klon_mpi, ndays)) ALLOCATE(pctsrf_t_mpi(klon_mpi, nbsrf, ndays)) ALLOCATE(phy_bil_mpi(klon_mpi, ndays)) CALL gather_omp(pctsrf_t, pctsrf_t_mpi) CALL gather_omp(sst_year, sst_year_mpi) CALL gather_omp(phy_bil, phy_bil_mpi) CALL gather_omp(albedo_year, albedo_year_mpi) CALL gather_omp(rugos_year, rugos_year_mpi) IF (is_omp_master) THEN CALL xios_send_field("foce_limout", pctsrf_t_mpi(:, is_oce, :)) CALL xios_send_field("fsic_limout", pctsrf_t_mpi(:, is_sic, :)) CALL xios_send_field("fter_limout", pctsrf_t_mpi(:, is_ter, :)) CALL xios_send_field("flic_limout", pctsrf_t_mpi(:, is_lic, :)) CALL xios_send_field("sst_limout", sst_year_mpi) CALL xios_send_field("bils_limout", phy_bil_mpi) CALL xios_send_field("alb_limout", albedo_year_mpi) CALL xios_send_field("rug_limout", rugos_year_mpi) ENDIF END SUBROUTINE create_limit_unstruct SUBROUTINE time_interpolation(ndays, field_in, calendar, field_out) USE lmdz_libmath_pch, ONLY: pchsp_95, pchfe_95 USE lmdz_arth, ONLY: arth USE dimphy, ONLY: klon USE ioipsl, ONLY: ioget_year_len USE time_phylmdz_mod, ONLY: annee_ref USE lmdz_phys_para USE lmdz_abort_physic, ONLY: abort_physic USE lmdz_iniprint, ONLY: lunout, prt_level IMPLICIT NONE INTEGER, INTENT(IN) :: ndays REAL, INTENT(IN) :: field_in(klon, lmdep) CHARACTER(LEN = *), INTENT(IN) :: calendar REAL, INTENT(OUT) :: field_out(klon, ndays) INTEGER :: ndays_in REAL :: timeyear(lmdep) REAL :: yder(lmdep) INTEGER :: ij, ierr, n_extrap LOGICAL :: skip CHARACTER (len = 50) :: modname = 'create_limit_unstruct.time_interpolation' CHARACTER (len = 80) :: abort_message IF (is_omp_master) ndays_in = year_len(annee_ref, calendar) CALL bcast_omp(ndays_in) IF (is_omp_master) timeyear = mid_months(annee_ref, calendar, lmdep) CALL bcast_omp(timeyear) n_extrap = 0 skip = .FALSE. DO ij = 1, klon yder = pchsp_95(timeyear, field_in(ij, :), ibeg = 2, iend = 2, vc_beg = 0., vc_end = 0.) CALL pchfe_95(timeyear, field_in(ij, :), yder, skip, arth(0., real(ndays_in) / ndays, ndays), field_out(ij, :), ierr) IF (ierr < 0) THEN abort_message = 'error in pchfe_95' CALL abort_physic(modname, abort_message, 1) endif n_extrap = n_extrap + ierr END DO IF (n_extrap /= 0) THEN WRITE(lunout, *) "get_2Dfield pchfe_95: n_extrap = ", n_extrap ENDIF END SUBROUTINE time_interpolation !------------------------------------------------------------------------------- FUNCTION year_len(y, cal_in) !------------------------------------------------------------------------------- USE ioipsl, ONLY: ioget_calendar, ioconf_calendar, lock_calendar, ioget_year_len IMPLICIT NONE !------------------------------------------------------------------------------- ! Arguments: INTEGER :: year_len INTEGER, INTENT(IN) :: y CHARACTER(LEN = *), INTENT(IN) :: cal_in !------------------------------------------------------------------------------- ! Local variables: CHARACTER(LEN = 20) :: cal_out ! calendar (for outputs) !------------------------------------------------------------------------------- !--- Getting the input calendar to reset at the end of the function CALL ioget_calendar(cal_out) !--- Unlocking calendar and setting it to wanted one CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) !--- Getting the number of days in this year year_len = ioget_year_len(y) !--- Back to original calendar CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) END FUNCTION year_len !------------------------------------------------------------------------------- !------------------------------------------------------------------------------- FUNCTION mid_months(y, cal_in, nm) !------------------------------------------------------------------------------- USE ioipsl, ONLY: ioget_calendar, ioconf_calendar, lock_calendar, ioget_mon_len USE lmdz_abort_physic, ONLY: abort_physic IMPLICIT NONE !------------------------------------------------------------------------------- ! Arguments: INTEGER, INTENT(IN) :: y ! year CHARACTER(LEN = *), INTENT(IN) :: cal_in ! calendar INTEGER, INTENT(IN) :: nm ! months/year number REAL, DIMENSION(nm) :: mid_months ! mid-month times !------------------------------------------------------------------------------- ! Local variables: CHARACTER(LEN = 99) :: mess ! error message CHARACTER(LEN = 20) :: cal_out ! calendar (for outputs) INTEGER, DIMENSION(nm) :: mnth ! months lengths (days) INTEGER :: m ! months counter INTEGER :: nd ! number of days INTEGER :: k !------------------------------------------------------------------------------- nd = year_len(y, cal_in) IF(nm==12) THEN !--- Getting the input calendar to reset at the end of the function CALL ioget_calendar(cal_out) !--- Unlocking calendar and setting it to wanted one CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) !--- Getting the length of each month DO m = 1, nm; mnth(m) = ioget_mon_len(y, m); END DO !--- Back to original calendar CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) ELSE IF(MODULO(nd, nm)/=0) THEN WRITE(mess, '(a,i3,a,i3,a)')'Unconsistent calendar: ', nd, ' days/year, but ', & nm, ' months/year. Months number should divide days number.' CALL abort_physic('mid_months', TRIM(mess), 1) ELSE mnth = (/(m, m = 1, nm, nd / nm)/) END IF !--- Mid-months times mid_months(1) = 0.5 * REAL(mnth(1)) DO k = 2, nm mid_months(k) = mid_months(k - 1) + 0.5 * REAL(mnth(k - 1) + mnth(k)) END DO END FUNCTION mid_months END MODULE create_limit_unstruct_mod