Ignore:
Timestamp:
Nov 17, 2025, 3:50:06 PM (6 weeks ago)
Author:
yann meurdesoif
Message:

Separate pbl_surface into 3 subroutines for GPU port

  • pbl_surface_uncompress_pre : prepare computation for sub subsurface before compressing
  • pbl_surface_subsrf : each sub-surface is called one after other (horizontal = knon)
  • pbl_surface_uncompress_post : sub-surface are uncompressed, computation is done on whole domain (horizontal = klon)

pbl_surface_main becomes the driver, calling pbl_surface_uncompress_pre, and then looping under sub-surface (and calling pbl_surface_subsrf) and then calling pbl_surface_uncompress_post.

YM

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/PBLSURF_GPUPORT/libf/phylmd/ocean_forced_mod.F90

    r5662 r5868  
    5959!****************************************************************************************
    6060    INTEGER, INTENT(IN)                      :: itime, jour, knon
    61     INTEGER, DIMENSION(klon), INTENT(IN)     :: knindex
     61    INTEGER, DIMENSION(knon), INTENT(IN)     :: knindex
    6262    REAL, INTENT(IN)                         :: dtime
    63     REAL, DIMENSION(klon), INTENT(IN)        :: p1lay
    64     REAL, DIMENSION(klon), INTENT(IN)        :: cdragh, cdragq, cdragm
    65     REAL, DIMENSION(klon), INTENT(IN)        :: precip_rain, precip_snow
    66     REAL, DIMENSION(klon), INTENT(IN)        :: temp_air, spechum
    67     REAL, DIMENSION(klon), INTENT(IN)        :: AcoefH, AcoefQ, BcoefH, BcoefQ
    68     REAL, DIMENSION(klon), INTENT(IN)        :: AcoefU, AcoefV, BcoefU, BcoefV
    69     REAL, DIMENSION(klon), INTENT(IN)        :: ps
    70     REAL, DIMENSION(klon), INTENT(IN)        :: u1, v1, gustiness
    71     REAL, DIMENSION(klon), INTENT(IN)        :: tsurf_in
    72     real, intent(in):: rhoa(:) ! (knon) density of moist air  (kg / m3)
    73 !GG
    74      REAL, DIMENSION(klon), INTENT(IN)        :: dthetadz300
     63    REAL, DIMENSION(knon), INTENT(IN)        :: p1lay
     64    REAL, DIMENSION(knon), INTENT(IN)        :: cdragh, cdragq, cdragm
     65    REAL, DIMENSION(knon), INTENT(IN)        :: precip_rain, precip_snow
     66    REAL, DIMENSION(knon), INTENT(IN)        :: temp_air, spechum
     67    REAL, DIMENSION(knon), INTENT(IN)        :: AcoefH, AcoefQ, BcoefH, BcoefQ
     68    REAL, DIMENSION(knon), INTENT(IN)        :: AcoefU, AcoefV, BcoefU, BcoefV
     69    REAL, DIMENSION(knon), INTENT(IN)        :: ps
     70    REAL, DIMENSION(knon), INTENT(IN)        :: u1, v1, gustiness
     71    REAL, DIMENSION(knon), INTENT(IN)        :: tsurf_in
     72    real, intent(in):: rhoa(knon) ! (knon) density of moist air  (kg / m3)
     73!GG
     74     REAL, DIMENSION(knon), INTENT(IN)        :: dthetadz300
    7575     REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf
    7676!
    7777
    7878#ifdef ISO
    79     REAL, DIMENSION(ntiso,klon), INTENT(IN)  :: xtprecip_rain, xtprecip_snow
    80     REAL, DIMENSION(ntiso,klon), INTENT(IN)  :: xtspechum
     79    REAL, DIMENSION(ntiso,knon), INTENT(IN)  :: xtprecip_rain, xtprecip_snow
     80    REAL, DIMENSION(ntiso,knon), INTENT(IN)  :: xtspechum
    8181    REAL, DIMENSION(klon),       INTENT(IN)  :: rlat
    8282#endif
    … …  
    8484! In/Output arguments
    8585!****************************************************************************************
    86     REAL, DIMENSION(klon), INTENT(INOUT)     :: radsol
    87     REAL, DIMENSION(klon), INTENT(INOUT)     :: snow
    88     REAL, DIMENSION(klon), INTENT(INOUT)     :: agesno !? put to 0 in ocean
     86    REAL, DIMENSION(knon), INTENT(INOUT)     :: radsol
     87    REAL, DIMENSION(knon), INTENT(INOUT)     :: snow
     88    REAL, DIMENSION(knon), INTENT(INOUT)     :: agesno !? put to 0 in ocean
    8989#ifdef ISO     
    90     REAL, DIMENSION(niso,klon), INTENT(IN)   :: xtsnow
    91     REAL, DIMENSION(niso,klon), INTENT(INOUT):: Roce
     90    REAL, DIMENSION(niso,knon), INTENT(IN)   :: xtsnow
     91    REAL, DIMENSION(niso,knon), INTENT(INOUT):: Roce
    9292#endif
    9393
    9494! Output arguments
    9595!****************************************************************************************
    96     REAL, DIMENSION(klon), INTENT(OUT)       :: qsurf
    97     REAL, DIMENSION(klon), INTENT(OUT)       :: evap, fluxsens, fluxlat
    98     REAL, DIMENSION(klon), INTENT(OUT)       :: flux_u1, flux_v1
    99     REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    100     REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l
     96    REAL, DIMENSION(knon), INTENT(OUT)       :: qsurf
     97    REAL, DIMENSION(knon), INTENT(OUT)       :: evap, fluxsens, fluxlat
     98    REAL, DIMENSION(knon), INTENT(OUT)       :: flux_u1, flux_v1
     99    REAL, DIMENSION(knon), INTENT(OUT)       :: tsurf_new
     100    REAL, DIMENSION(knon), INTENT(OUT)       :: dflux_s, dflux_l
    101101    REAL, INTENT(out):: sens_prec_liq(:) ! (knon)
    102102!GG
    103      REAL, DIMENSION(klon), INTENT(OUT)       :: Ampl
     103     REAL, DIMENSION(knon), INTENT(OUT)       :: Ampl
    104104!
    105105
    106106#ifdef ISO     
    107     REAL, DIMENSION(ntiso,klon), INTENT(OUT) :: xtevap ! isotopes in evaporation flux
    108     REAL, DIMENSION(klon),       INTENT(OUT) :: h1 ! just a diagnostic, not useful for the simulation
     107    REAL, DIMENSION(ntiso,knon), INTENT(OUT) :: xtevap ! isotopes in evaporation flux
     108    REAL, DIMENSION(knon),       INTENT(OUT) :: h1 ! just a diagnostic, not useful for the simulation
    109109#endif
    110110
    … …  
    112112!****************************************************************************************
    113113    INTEGER                     :: i, j
    114     REAL, DIMENSION(klon)       :: cal, beta, dif_grnd
    115     REAL, DIMENSION(klon)       :: alb_neig, tsurf_lim, zx_sl
    116     REAL, DIMENSION(klon)       :: u0, v0
    117     REAL, DIMENSION(klon)       :: u1_lay, v1_lay
     114    REAL, DIMENSION(knon)       :: cal, beta, dif_grnd
     115    REAL, DIMENSION(knon)       :: alb_neig, tsurf_lim, zx_sl
     116    REAL, DIMENSION(knon)       :: u0, v0
     117    REAL, DIMENSION(knon)       :: u1_lay, v1_lay
    118118    LOGICAL                     :: check=.FALSE.
    119119    REAL, DIMENSION(knon)       :: sens_prec_sol
    120     REAL, DIMENSION(klon)       :: lat_prec_liq, lat_prec_sol   
     120    REAL, DIMENSION(knon)       :: lat_prec_liq, lat_prec_sol   
    121121! GG
    122     REAL, DIMENSION(klon)       :: l_CBL, sicfra
     122    REAL, DIMENSION(knon)       :: l_CBL, sicfra
    123123!
    124124#ifdef ISO   
    … …  
    198198
    199199! Calcul de tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l and qsurf
    200     CALL calcul_fluxs(knon, is_oce, dtime, &
    201          merge(tsurf_in, tsurf_lim, activate_ocean_skin == 2), p1lay, cal, &
     200    IF (activate_ocean_skin == 2) THEN
     201      CALL calcul_fluxs(knon, is_oce, dtime, &
     202         tsurf_in, p1lay, cal, &
    202203         beta, cdragh, cdragq, ps, &
    203204         precip_rain, precip_snow, snow, qsurf,  &
    … …  
    206207         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
    207208         sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol, rhoa)
    208     if (activate_ocean_skin == 2) tsurf_new = tsurf_lim
     209         tsurf_new = tsurf_lim
     210    ELSE
     211      CALL calcul_fluxs(knon, is_oce, dtime, &
     212         tsurf_lim, p1lay, cal, &
     213         beta, cdragh, cdragq, ps, &
     214         precip_rain, precip_snow, snow, qsurf,  &
     215         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, gustiness, &
     216         f_qsat_oce,AcoefH, AcoefQ, BcoefH, BcoefQ, &
     217         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l, &
     218         sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol, rhoa)
     219    ENDIF
    209220
    210221    do j = 1, knon
    … …  
    217228
    218229!GG
     230
    219231    if (iflag_leads == 1) then
    220       l_CBL = -52381.*dthetadz300 + 3008.1
    221       Ampl = 6.012e-08*l_CBL**2 - 4.036e-04*l_CBL + 1.4979
    222       WHERE(Ampl(:)>1.2) Ampl(:)=1.2
    223       sicfra(:)=pctsrf(:,is_sic)/(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter))
    224       WHERE(pctsrf(:,is_sic)+pctsrf(:,is_oce)<EPSFRA) sicfra(:)=0.
    225       WHERE(sicfra<0.7) Ampl(:)=1.
    226       WHERE((sicfra>0.7).and.(sicfra<0.9)) Ampl=((sicfra-0.7)/0.2)*Ampl+((0.9-sicfra)/0.2)
    227       fluxsens=Ampl*fluxsens
    228       dflux_s=Ampl*dflux_s
     232!ym hyper faux, melange de champ compresse et non compresse, je rétabli...
     233!ym      l_CBL = -52381.*dthetadz300 + 3008.1
     234!ym      Ampl = 6.012e-08*l_CBL**2 - 4.036e-04*l_CBL + 1.4979
     235!ym      WHERE(Ampl(:)>1.2) Ampl(:)=1.2
     236!ym      sicfra(:)=pctsrf(:,is_sic)/(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter))
     237!ym      WHERE(pctsrf(:,is_sic)+pctsrf(:,is_oce)<EPSFRA) sicfra(:)=0.
     238!ym      WHERE(sicfra<0.7) Ampl(:)=1.
     239!ym      WHERE((sicfra>0.7).and.(sicfra<0.9)) Ampl=((sicfra-0.7)/0.2)*Ampl+((0.9-sicfra)/0.2)
     240!ym      fluxsens=Ampl*fluxsens
     241!ym      dflux_s=Ampl*dflux_s
     242     
     243      do j = 1, knon
     244        i = knindex(j)
     245        l_CBL(j) = -52381.*dthetadz300(j) + 3008.1
     246        Ampl(j) = 6.012e-08*l_CBL(j)**2 - 4.036e-04*l_CBL(j) + 1.4979
     247        IF (Ampl(j)>1.2) Ampl(j)=1.2
     248        sicfra(j)=pctsrf(i,is_sic)/(1.-pctsrf(i,is_lic)-pctsrf(i,is_ter))
     249        IF (pctsrf(i,is_sic)+pctsrf(i,is_oce)<EPSFRA) sicfra(j)=0.
     250        IF (sicfra(j)<0.7) Ampl(j)=1.
     251        IF (sicfra(j)>0.7 .and. sicfra(j)<0.9) Ampl(j)=((sicfra(j)-0.7)/0.2)*Ampl(j)+((0.9-sicfra(j))/0.2)
     252        fluxsens(j)=Ampl(j)*fluxsens(j)
     253        dflux_s(j)=Ampl(j)*dflux_s(j)
     254      enddo
    229255    endif
    230256
    … …  
    238264
    239265#ifdef ISO     
    240     CALL calcul_iso_surf_oce_vectall(klon, knon,t_coup, &
     266    CALL calcul_iso_surf_oce_vectall(knon, knon,t_coup, &
    241267     &    ps,tsurf_new,spechum,u1_lay, v1_lay, xtspechum, &
    242268     &    evap, Roce,xtevap,h1 &
    … …  
    324350!****************************************************************************************
    325351    INTEGER, INTENT(IN)                  :: itime, jour, knon
    326     INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
     352    INTEGER, DIMENSION(knon), INTENT(IN) :: knindex
    327353    REAL, INTENT(IN)                     :: dtime
    328     REAL, DIMENSION(klon), INTENT(IN)    :: tsurf_in
    329     REAL, DIMENSION(klon), INTENT(IN)    :: p1lay
    330     REAL, DIMENSION(klon), INTENT(IN)    :: cdragh, cdragm
    331     REAL, DIMENSION(klon), INTENT(IN)    :: precip_rain, precip_snow
    332     REAL, DIMENSION(klon), INTENT(IN)    :: temp_air, spechum
    333     REAL, DIMENSION(klon), INTENT(IN)    :: AcoefH, AcoefQ, BcoefH, BcoefQ
    334     REAL, DIMENSION(klon), INTENT(IN)    :: AcoefU, AcoefV, BcoefU, BcoefV
    335     REAL, DIMENSION(klon), INTENT(IN)    :: ps
    336     REAL, DIMENSION(klon), INTENT(IN)    :: u1, v1, gustiness
    337     real, intent(in):: rhoa(:) ! (knon) density of moist air  (kg / m3)
    338 !GG
    339     REAL, DIMENSION(klon), INTENT(IN)    :: swnet
     354    REAL, DIMENSION(knon), INTENT(IN)    :: tsurf_in
     355    REAL, DIMENSION(knon), INTENT(IN)    :: p1lay
     356    REAL, DIMENSION(knon), INTENT(IN)    :: cdragh, cdragm
     357    REAL, DIMENSION(knon), INTENT(IN)    :: precip_rain, precip_snow
     358    REAL, DIMENSION(knon), INTENT(IN)    :: temp_air, spechum
     359    REAL, DIMENSION(knon), INTENT(IN)    :: AcoefH, AcoefQ, BcoefH, BcoefQ
     360    REAL, DIMENSION(knon), INTENT(IN)    :: AcoefU, AcoefV, BcoefU, BcoefV
     361    REAL, DIMENSION(knon), INTENT(IN)    :: ps
     362    REAL, DIMENSION(knon), INTENT(IN)    :: u1, v1, gustiness
     363    real, intent(in):: rhoa(knon) ! (knon) density of moist air  (kg / m3)
     364!GG
     365    REAL, DIMENSION(knon), INTENT(IN)    :: swnet
    340366    REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf
    341367!GG
    342368#ifdef ISO
    343     REAL, DIMENSION(ntiso,klon), INTENT(IN) :: xtprecip_rain, xtprecip_snow
    344     REAL, DIMENSION(ntiso,klon), INTENT(IN) :: xtspechum
    345     REAL, DIMENSION(niso,klon),  INTENT(IN) :: Roce
    346     REAL, DIMENSION(niso,klon),  INTENT(IN) :: Rland_ice
     369    REAL, DIMENSION(ntiso,knon), INTENT(IN) :: xtprecip_rain, xtprecip_snow
     370    REAL, DIMENSION(ntiso,knon), INTENT(IN) :: xtspechum
     371    REAL, DIMENSION(niso,knon),  INTENT(IN) :: Roce
     372    REAL, DIMENSION(niso,knon),  INTENT(IN) :: Rland_ice
    347373#endif
    348374
    349375! In/Output arguments
    350376!****************************************************************************************
    351     REAL, DIMENSION(klon), INTENT(INOUT)          :: radsol
    352     REAL, DIMENSION(klon), INTENT(INOUT)          :: snow, qsol
    353     REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
    354     REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
    355 !GG
    356     REAL, DIMENSION(klon), INTENT(INOUT)          :: hice
    357     REAL, DIMENSION(klon), INTENT(INOUT)          :: tice
    358     REAL, DIMENSION(klon), INTENT(INOUT)          :: bilg_cumul
    359     REAL, DIMENSION(klon), INTENT(INOUT)          :: fcds
    360     REAL, DIMENSION(klon), INTENT(INOUT)          :: fcdi
    361     REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_basal_growth
    362     REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_basal_melt
    363     REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_top_melt
    364     REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_snow2sic
    365     REAL, DIMENSION(klon), INTENT(INOUT)          :: dtice_melt
    366     REAL, DIMENSION(klon), INTENT(INOUT)          :: dtice_snow2sic
     377    REAL, DIMENSION(knon), INTENT(INOUT)          :: radsol
     378    REAL, DIMENSION(knon), INTENT(INOUT)          :: snow, qsol
     379    REAL, DIMENSION(knon), INTENT(INOUT)          :: agesno
     380    REAL, DIMENSION(knon, nsoilmx), INTENT(INOUT) :: tsoil
     381!GG
     382    REAL, DIMENSION(klon), INTENT(INOUT)          :: hice !ym WARNING uncompressed
     383    REAL, DIMENSION(klon), INTENT(INOUT)          :: tice !ym WARNING uncompressed
     384    REAL, DIMENSION(klon), INTENT(INOUT)          :: bilg_cumul !ym WARNING uncompressed
     385    REAL, DIMENSION(klon), INTENT(INOUT)          :: fcds !ym WARNING uncompressed
     386    REAL, DIMENSION(klon), INTENT(INOUT)          :: fcdi !ym WARNING uncompressed
     387    REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_basal_growth !ym WARNING uncompressed
     388    REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_basal_melt !ym WARNING uncompressed
     389    REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_top_melt !ym WARNING uncompressed
     390    REAL, DIMENSION(klon), INTENT(INOUT)          :: dh_snow2sic !ym WARNING uncompressed
     391    REAL, DIMENSION(klon), INTENT(INOUT)          :: dtice_melt !ym WARNING uncompressed
     392    REAL, DIMENSION(klon), INTENT(INOUT)          :: dtice_snow2sic !ym WARNING uncompressed
    367393!GG
    368394#ifdef ISO     
    369     REAL, DIMENSION(niso,klon), INTENT(INOUT)     :: xtsnow
    370     REAL, DIMENSION(niso,klon), INTENT(IN)        :: xtsol
     395    REAL, DIMENSION(niso,knon), INTENT(INOUT)     :: xtsnow
     396    REAL, DIMENSION(niso,knon), INTENT(IN)        :: xtsol
    371397#endif
    372398
    373399! Output arguments
    374400!****************************************************************************************
    375     REAL, DIMENSION(klon), INTENT(OUT)            :: qsurf
    376     REAL, DIMENSION(klon), INTENT(OUT)            :: alb1_new  ! new albedo in visible SW interval
    377     REAL, DIMENSION(klon), INTENT(OUT)            :: alb2_new  ! new albedo in near IR interval
    378     REAL, DIMENSION(klon), INTENT(OUT)            :: evap, fluxsens, fluxlat
    379     REAL, DIMENSION(klon), INTENT(OUT)            :: flux_u1, flux_v1
    380     REAL, DIMENSION(klon), INTENT(OUT)            :: tsurf_new
    381     REAL, DIMENSION(klon), INTENT(OUT)            :: dflux_s, dflux_l     
     401    REAL, DIMENSION(knon), INTENT(OUT)            :: qsurf
     402    REAL, DIMENSION(knon), INTENT(OUT)            :: alb1_new  ! new albedo in visible SW interval
     403    REAL, DIMENSION(knon), INTENT(OUT)            :: alb2_new  ! new albedo in near IR interval
     404    REAL, DIMENSION(knon), INTENT(OUT)            :: evap, fluxsens, fluxlat
     405    REAL, DIMENSION(knon), INTENT(OUT)            :: flux_u1, flux_v1
     406    REAL, DIMENSION(knon), INTENT(OUT)            :: tsurf_new
     407    REAL, DIMENSION(knon), INTENT(OUT)            :: dflux_s, dflux_l     
    382408#ifdef ISO     
    383     REAL, DIMENSION(ntiso,klon), INTENT(OUT)      :: xtevap
     409    REAL, DIMENSION(ntiso,knon), INTENT(OUT)      :: xtevap
    384410#endif     
    385411
    … …  
    390416    REAL                        :: zfra
    391417    REAL, PARAMETER             :: t_grnd=271.35
    392     REAL, DIMENSION(klon)       :: cal, beta, dif_grnd, capsol, icesub
    393     REAL, DIMENSION(klon)       :: alb_neig, tsurf_tmp
    394     REAL, DIMENSION(klon)       :: soilcap, soilflux
    395     REAL, DIMENSION(klon)       :: u0, v0
    396     REAL, DIMENSION(klon)       :: u1_lay, v1_lay
     418    REAL, DIMENSION(knon)       :: cal, beta, dif_grnd, capsol, icesub
     419    REAL, DIMENSION(knon)       :: alb_neig, tsurf_tmp
     420    REAL, DIMENSION(knon)       :: soilcap, soilflux
     421    REAL, DIMENSION(knon)       :: u0, v0
     422    REAL, DIMENSION(knon)       :: u1_lay, v1_lay
    397423    REAL, DIMENSION(knon)       :: sens_prec_liq, sens_prec_sol
    398     REAL, DIMENSION(klon)       :: lat_prec_liq, lat_prec_sol   
     424    REAL, DIMENSION(knon)       :: lat_prec_liq, lat_prec_sol   
     425    REAL, DIMENSION(knon)       :: yhice ! compressed version of hice
    399426!GG
    400427    INTEGER                     :: ki
    401428    INTEGER                     :: cpl_pas
    402     REAL, DIMENSION(klon)       :: bilg, fsic, f_bot
     429    REAL, DIMENSION(knon)       :: bilg
     430    REAL, DIMENSION(knon)       :: fsic
     431    REAL, DIMENSION(knon)       :: f_bot
    403432    REAL, PARAMETER             :: latent_ice = 334.0e3
    404433    REAL, PARAMETER             :: rau_ice = 917.0
    … …  
    471500#ifdef ISO
    472501    REAL, PARAMETER :: t_coup = 273.15
    473     REAL, DIMENSION(klon) :: fq_fonte_diag
    474     REAL, DIMENSION(klon) :: fqfonte_diag
    475     REAL, DIMENSION(klon) :: snow_evap_diag
    476     REAL, DIMENSION(klon) :: fqcalving_diag
    477     REAL, DIMENSION(klon) :: run_off_lic_diag
     502    REAL, DIMENSION(knon) :: fq_fonte_diag
     503    REAL, DIMENSION(knon) :: fqfonte_diag
     504    REAL, DIMENSION(knon) :: snow_evap_diag
     505    REAL, DIMENSION(knon) :: fqcalving_diag
     506    REAL, DIMENSION(knon) :: run_off_lic_diag
    478507    REAL :: coeff_rel_diag
    479508    REAL :: max_eau_sol_diag 
    480     REAL, DIMENSION(klon) :: runoff_diag   
     509    REAL, DIMENSION(knon) :: runoff_diag   
    481510    INTEGER IXT
    482     REAL, DIMENSION(niso,klon) :: xtsnow_prec, xtsol_prec
    483     REAL, DIMENSION(klon) :: snow_prec, qsol_prec 
     511    REAL, DIMENSION(niso,knon) :: xtsnow_prec, xtsol_prec
     512    REAL, DIMENSION(knon) :: snow_prec, qsol_prec 
    484513#endif
    485514
    … …  
    536565    dif_grnd(:)=0.
    537566    beta(:) = 1.
    538     fsic(:) = pctsrf(:,is_sic)
    539567    cpl_pas =  NINT(86400./dtime * 1.0) ! une fois par jour
    540568
    … …  
    545573    DO i=1,knon
    546574        ki=knindex(i)
     575        fsic(i) = pctsrf(ki,is_sic)
    547576        IF (snow(i).GT.snow_min) THEN
    548577            !  1 / snow-layer heat capacity
    … …  
    566595            ! GG no conductive flux in this case?
    567596        END IF
    568         bilg(ki)=f_swpen-f_cond_s
     597        bilg(i)=f_swpen-f_cond_s
    569598    END DO
    570599
    … …  
    652681!        write(*,*) 'klon,knon=',klon,knon
    653682!#endif
    654     CALL calcul_iso_surf_sic_vectall(klon,knon, &
     683    CALL calcul_iso_surf_sic_vectall(knon,knon, &
    655684     &          evap,snow_evap_diag,Tsurf_new,Roce,snow, &
    656685     &          fq_fonte_diag,fqfonte_diag,dtime,t_coup, &
    … …  
    676705! Calculation of albedo at snow (alb_neig) and update the age of snow (agesno)
    677706!
    678     CALL albsno(klon, knon, dtime, agesno(:), alb_neig(:), precip_snow(:)) 
     707    CALL albsno(knon, knon, dtime, agesno(:), alb_neig(:), precip_snow(:)) 
    679708
    680709    WHERE (snow(1:knon) .LT. 0.0001) agesno(1:knon) = 0.
    … …  
    719748           fcdi(ki) = f_cond_i - fsea
    720749           fcds(i) = f_cond_s
    721            !bilg(ki) = bilg(ki)+f_cond_i
     750           !bilg(i) = bilg(i)+f_cond_i
    722751        END DO
    723752
    … …  
    728757    IF (iflag_seaice==1) THEN
    729758!   Read from limit
    730     CALL limit_read_hice(knon,knindex,hice)
     759!ym   totally wrong since "hice" is an uncompressed field (klon) and limit_read_hice return a compressed field (knon)
     760!ym    CALL limit_read_hice(knon,knindex,hice)
     761      CALL limit_read_hice(knon, knindex, yhice)
     762      DO i=1,knon
     763        ki=knindex(i)
     764        hice(ki) = yhice(i)
     765      ENDDO
    731766    ENDIF
    732767!   Formula Krinner et al. 1997 : h = (0.2 + 3.8(f_min**2))(1 + 2(f- f_min))
    … …  
    736771    DO i=1,knon
    737772        ki=knindex(i)
     773!ym WARNING : fsic uninitilazed, take initialisation similar to  iflag_seaice==2
     774        fsic(i) = pctsrf(ki,is_sic)       
    738775        IF (precip_snow(i) > 0.) THEN
    739             snow(i) = snow(i)+precip_snow(i)*dtime*(1.-snow_wfact*(1.-fsic(ki)))
     776            snow(i) = snow(i)+precip_snow(i)*dtime*(1.-snow_wfact*(1.-fsic(i)))
    740777        END IF
    741778! snow and ice sublimation
    … …  
    774811                  / (h_ice_max-h_ice_thick)))
    775812           ! quantity of new ice (formed at mean ice temp)
    776            e_melt= -f_bot(i) * dtime * fsic(ki) &
     813           e_melt= -f_bot(i) * dtime * fsic(i) &
    777814                   / (ice_lat+ice_cap/2.*(t_freeze-tice(ki)))
    778815           ! first increase height to h_thick
    779            dhsic=MAX(0.,MIN(h_ice_thick-hice(ki),e_melt/(fsic(ki)*ice_den)))
     816           dhsic=MAX(0.,MIN(h_ice_thick-hice(ki),e_melt/(fsic(i)*ice_den)))
    780817           hice_i=hice(ki)
    781818           hice(ki)=dhsic+hice(ki)
    782            e_melt=e_melt-fsic(ki)*dhsic
     819           e_melt=e_melt-fsic(i)*dhsic
    783820           IF (e_melt.GT.0.) THEN
    784821           ! frac_mf fraction used for lateral increase
    785            dfsic=MIN(amax-fsic(ki),e_melt*frac_mf/ (hice(ki)*ice_den) )
     822           dfsic=MIN(amax-fsic(i),e_melt*frac_mf/ (hice(ki)*ice_den) )
    786823           ! No lateral growth -> forced ocean
    787824           !fsic(ki)=fsic(ki)+dfsic
    788825           e_melt=e_melt-dfsic*hice(ki)*ice_den
    789826           ! rest used to increase height
    790            hice(ki)=MIN(h_ice_max,hice(ki)+e_melt/( fsic(ki) * ice_den ) )
     827           hice(ki)=MIN(h_ice_max,hice(ki)+e_melt/( fsic(i) * ice_den ) )
    791828           END IF
    792829           dh_basal_growth(ki)=(hice(ki)-hice_i)/dtime
    … …  
    799836           ! bring ice to freezing and melt from below
    800837           ! quantity of melted ice
    801            e_melt= f_bot(i) * dtime * fsic(ki) &
     838           e_melt= f_bot(i) * dtime * fsic(i) &
    802839                   / (ice_lat+ice_cap/2.*(tice(ki)-t_freeze))
    803840           ! first decrease height to h_thick
    804841           hice_i=hice(ki)
    805            dhsic=MAX(0.,MIN(hice(ki)-h_ice_thick,e_melt/(fsic(ki)*ice_den)))
     842           dhsic=MAX(0.,MIN(hice(ki)-h_ice_thick,e_melt/(fsic(i)*ice_den)))
    806843           hice(ki)=hice(ki)-dhsic
    807            e_melt=e_melt-fsic(ki)*dhsic*ice_den
     844           e_melt=e_melt-fsic(i)*dhsic*ice_den
    808845
    809846           IF (e_melt.GT.0) THEN
    810847           ! frac_mf fraction used for height decrease
    811            dhsic=MAX(0.,MIN(hice(ki)-h_ice_min,e_melt/ice_den*frac_mf/fsic(ki)))
     848           dhsic=MAX(0.,MIN(hice(ki)-h_ice_min,e_melt/ice_den*frac_mf/fsic(i)))
    812849           hice(ki)=hice(ki)-dhsic
    813            e_melt=e_melt-fsic(ki)*dhsic*ice_den
     850           e_melt=e_melt-fsic(i)*dhsic*ice_den
    814851           ! rest used to decrease fraction (up to 0!)
    815            dfsic=MIN(fsic(ki),e_melt/(hice(ki)*ice_den))
     852           dfsic=MIN(fsic(i),e_melt/(hice(ki)*ice_den))
    816853           ! Remaining heat not used if everything melted
    817854           e_melt=e_melt-dfsic*hice(ki)*ice_den
    … …  
    836873            IF (e_melt.GT.0) THEN
    837874              ! lateral melt if ice too thin
    838               dfsic=MAX(fsic(ki)-ice_frac_min,e_melt/(h_ice_min*ice_den)*fsic(ki))
     875              dfsic=MAX(fsic(i)-ice_frac_min,e_melt/(h_ice_min*ice_den)*fsic(i))
    839876              ! if all melted do nothing with remaining heat
    840               e_melt=MAX(0.,e_melt*fsic(ki)-dfsic*h_ice_min*ice_den)
     877              e_melt=MAX(0.,e_melt*fsic(i)-dfsic*h_ice_min*ice_den)
    841878              ! slab_bilg(ki) = slab_bilg(ki) + e_melt*ice_lat/dtime
    842879            END IF
    … …  
    881918!***********************************************o*******************************
    882919    !cumul fluxes.
    883     bilg_cumul(:)=bilg_cumul(:)+bilg(:)/float(cpl_pas)
    884     IF (MOD(itime,cpl_pas).EQ.0) THEN ! time to update tslab
    885         bilg_cumul(:)=0.
    886     END IF ! coupling time
     920    DO j=i,knon
     921      ki=knindex(i)
     922      bilg_cumul(ki)=bilg_cumul(ki) + bilg(j)/float(cpl_pas)
     923    ENDDO
     924
     925!YM Pb for GPU port, done on an compressed field inside a compressed Kernel
     926!YM ==> move outside of subroutine
     927!YM CALL ocean_forced_ice_reset_bilg_cumul()
     928!YM    IF (MOD(itime,cpl_pas).EQ.0) THEN ! time to update tslab
     929!YM        bilg_cumul(1:klon)=0.
     930!YM    END IF ! coupling time
    887931
    888932!   write(*,*) 'hice 3',hice(1:100)
    889933!   write(*,*) 'tice 3',tice(1:100)
    890934    !tests ice fraction
    891     WHERE (fsic.LT.ice_frac_min)
    892         tice=t_melt
    893         hice=h_ice_min
    894     END WHERE
     935!YM Pb for GPU port : done on uncompressed field inside a compressed kernel
     936!YM update fsic only on compressed index
     937!YM    WHERE (fsic.LT.ice_frac_min)
     938!YM        tice=t_melt
     939!YM        hice=h_ice_min
     940!YM    END WHERE
     941
     942    DO j=i,knon
     943      ki=knindex(i)
     944      IF (fsic(i).LT.ice_frac_min ) THEN
     945        tice(ki)=t_melt
     946        hice(ki)=h_ice_min
     947      ENDIF
     948    ENDDO
    895949
    896950    !write(*,*) 'hice 4',hice(1:100)
    … …  
    10031057
    10041058      CASE(3)
    1005       CALL albsno(klon, knon, dtime, agesno(:), alb_neig(:), precip_snow(:))
     1059      CALL albsno(knon, knon, dtime, agesno(:), alb_neig(:), precip_snow(:))
    10061060      WHERE (snow(1:knon) .LT. 0.0001) agesno(1:knon) = 0.
    10071061      alb1_new(:) = 0.0
    … …  
    10241078  END SUBROUTINE ocean_forced_ice
    10251079
     1080  SUBROUTINE ocean_forced_ice_reset_bilg_cumul(itime, dtime, bilg_cumul)
     1081    USE dimphy, ONLY : klon
     1082    USE surface_data, ONLY : iflag_seaice, type_ocean, version_ocean
     1083    IMPLICIT NONE
     1084    INTEGER, INTENT(IN) :: itime
     1085    REAL, INTENT(IN)    :: dtime
     1086    REAL, DIMENSION(klon), INTENT(INOUT) :: bilg_cumul
     1087    INTEGER :: cpl_pas
     1088   
     1089    IF (iflag_seaice/=0) THEN
     1090      cpl_pas =  NINT(86400./dtime * 1.0) ! une fois par jour
     1091      IF (MOD(itime,cpl_pas).EQ.0) THEN ! time to update tslab
     1092        IF ( .NOT. (type_ocean == 'couple' .OR. (type_ocean == 'slab'.AND.version_ocean=='sicINT'))) THEN
     1093          bilg_cumul(1:klon)=0.
     1094        ENDIF
     1095      END IF ! coupling time
     1096    ENDIF
     1097
     1098  END SUBROUTINE ocean_forced_ice_reset_bilg_cumul
    10261099!************************************************************************
    10271100! 1D case
Note: See TracChangeset for help on using the changeset viewer.