Changeset 996


Ignore:
Timestamp:
Sep 9, 2008, 3:22:23 PM (16 years ago)
Author:
lsce
Message:
  • Modifications liées au calcul des nouveau sous-fractions
  • Nettoyage de ocean slab : il reste uniquement la version avec glace de mer forcé
  • Nouveaux variables pour distiguer la version et type d'ocean : type_ocean=force/slab/couple, version_ocean=opa8/nemo pour couplé ou version_ocean=sicOBS pour slab

JG

Location:
LMDZ4/trunk/libf/phylmd
Files:
4 added
1 deleted
27 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/trunk/libf/phylmd/calcul_divers.h

    r776 r996  
    1212c surface terre
    1313       DO i=1, klon
    14          IF(pctsrf_new(i,is_ter).GT.0.) THEN
    15             paire_ter(i)=airephy(i)*pctsrf_new(i,is_ter)
     14         IF(pctsrf(i,is_ter).GT.0.) THEN
     15            paire_ter(i)=airephy(i)*pctsrf(i,is_ter)
    1616         ENDIF
    1717       ENDDO

  • LMDZ4/trunk/libf/phylmd/clesphys.h

    r900 r996  
    4747       REAL freq_ISCCP, ecrit_ISCCP
    4848       INTEGER :: ip_ebil_phy
    49        LOGICAL ok_slab_sicOBS
    5049
    5150       COMMON/clesphys/cycle_diurne, soil_model, new_oliq,              &
    … …  
    6160     &     , ecrit_mth, ecrit_tra, ecrit_reg                            &
    6261     &     , freq_ISCCP, ecrit_ISCCP, ip_ebil_phy                       &
    63      &     , ok_slab_sicOBS, ok_lic_melt, cvl_corr                      &
     62     &     , ok_lic_melt, cvl_corr                                      &
    6463     &     , qsol0
    6564     

  • LMDZ4/trunk/libf/phylmd/climb_hq_mod.F90

    r793 r996  
     1!
     2! $Header$
     3!
    14MODULE climb_hq_mod
    25!
    … …  
    274277    REAL, DIMENSION(klon,klev)               :: h_new, q_new
    275278    REAL, DIMENSION(klon)                    :: psref         
    276     INTEGER                                  :: k, i
     279    INTEGER                                  :: k, i, ierr
    277280
    278281!****************************************************************************************
    … …  
    350353!
    351354!****************************************************************************************
    352     DEALLOCATE(Ccoef_Q, Dcoef_Q, Ccoef_H, Dcoef_H)   
    353     DEALLOCATE(gamaq, gamah)
    354     DEALLOCATE(Kcoefhq)
    355 
     355    DEALLOCATE(Ccoef_Q, Dcoef_Q, Ccoef_H, Dcoef_H,stat=ierr)   
     356    IF ( ierr /= 0 )  PRINT*,' pb in dealllocate Ccoef_Q, Dcoef_Q, Ccoef_H, Dcoef_H, ierr=', ierr
     357    DEALLOCATE(gamaq, gamah,stat=ierr)
     358    IF ( ierr /= 0 )  PRINT*,' pb in dealllocate gamaq, gamah, ierr=', ierr
     359    DEALLOCATE(Kcoefhq,stat=ierr)
     360    IF ( ierr /= 0 )  PRINT*,' pb in dealllocate Kcoefhq, ierr=', ierr
    356361   
    357362  END SUBROUTINE climb_hq_up

  • LMDZ4/trunk/libf/phylmd/conf_phys.F90

    r973 r996  
    55!
    66
    7   subroutine conf_phys(ocean, ok_veget, ok_journe, ok_mensuel, ok_instan, ok_hf, &
     7  subroutine conf_phys(ok_journe, ok_mensuel, ok_instan, ok_hf, &
    88 &                     solarlong0,seuil_inversion, &
    99 &                     fact_cldcon, facttemps,ok_newmicro,iflag_radia,&
    … …  
    1616
    1717   use IOIPSL
    18 !!!!   USE surface_data,     ONLY : ocean, ok_veget
     18   USE surface_data
    1919
    2020   implicit none
    … …  
    5151! Sortie:
    5252  character (len = 6)  :: ocean
    53   logical              :: ok_veget, ok_newmicro
     53  logical              :: ok_newmicro
    5454  integer              :: iflag_radia
    5555  logical              :: ok_journe, ok_mensuel, ok_instan, ok_hf
    … …  
    116116  REAL,SAVE :: ecrit_hf_omp, ecrit_day_omp, ecrit_mth_omp, ecrit_reg_omp
    117117  REAL,SAVE :: ecrit_tra_omp
    118   LOGICAL, SAVE :: ok_slab_sicOBS_omp
    119118  REAL,SAVE :: cvl_corr_omp
    120119  LOGICAL,SAVE :: ok_lic_melt_omp
    … …  
    997996  call getin('ecrit_reg',ecrit_reg_omp)
    998997!
    999 !
    1000 !
    1001 !Config Key  = ok_slab_sicOBS
    1002 !Config Desc =
    1003 !Config Def  = .true.
    1004 !Config Help = Pour faire tourner le slab avec fraction
    1005 !              de glace de mer Observee
    1006 !
    1007   ok_slab_sicOBS_omp = .true.
    1008   call getin('ok_slab_sicOBS', ok_slab_sicOBS_omp)
    1009998!
    1010999!
    … …  
    11981187    ecrit_tra = ecrit_tra_omp
    11991188    ecrit_reg = ecrit_reg_omp
    1200     ok_slab_sicOBS = ok_slab_sicOBS_omp
    12011189    cvl_corr = cvl_corr_omp
    12021190    ok_lic_melt = ok_lic_melt_omp
    … …  
    12131201!$OMP MASTER
    12141202
     1203! Attribution of new parmeters according to parameters in .def
     1204    IF (ocean=='couple' .OR. ocean=='opa8') THEN
     1205       type_ocean='couple'
     1206       version_ocean='opa8'
     1207    ELSE IF (ocean=='nemo') THEN
     1208       type_ocean='couple'
     1209       version_ocean='nemo'
     1210    ELSE IF (ocean=='force') THEN
     1211       type_ocean='force'
     1212       version_ocean='xxxxxx'
     1213    ELSE IF (ocean=='slab') THEN
     1214       type_ocean='slab'
     1215       version_ocean='sicOBS'
     1216    ELSE
     1217       WRITE(numout,*)' ERROR ocean not valid : ocean= ', ocean
     1218       CALL abort_gcm('conf_phys','ocean not valid',1)
     1219    END IF
     1220
    12151221  write(numout,*)' ##############################################'
    12161222  write(numout,*)' Configuration des parametres de la physique: '
    12171223  write(numout,*)' Config ocean = ', ocean
     1224  write(numout,*)' Type ocean = ', type_ocean
     1225  write(numout,*)' Version ocean = ', version_ocean
    12181226  write(numout,*)' Config veget = ', ok_veget
    12191227  write(numout,*)' Sortie journaliere = ', ok_journe

  • LMDZ4/trunk/libf/phylmd/cpl_mod.F90

    r987 r996  
    3232  PRIVATE
    3333
    34   ! All subroutine are public except cpl_send_all and cpl_receive_all
    35   PUBLIC :: cpl_init, cpl_receive_ocean_fields, cpl_receive_seaice_fields, &
     34  ! All subroutine are public except cpl_send_all
     35  PUBLIC :: cpl_init, cpl_receive_frac, cpl_receive_ocean_fields, cpl_receive_seaice_fields, &
    3636       cpl_send_ocean_fields, cpl_send_seaice_fields, cpl_send_land_fields, &
    3737       cpl_send_landice_fields, gath2cpl
    … …  
    6868  !$OMP THREADPRIVATE(read_alb_sic)
    6969 
    70 ! fraction for different surface, saved during whole coupling period
    71   REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: pctsrf_sav   
    72   !$OMP THREADPRIVATE(pctsrf_sav)
    7370  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE  :: unity
    7471  !$OMP THREADPRIVATE(unity)
    … …  
    8885  !$OMP THREADPRIVATE(cpl_windsp2D)
    8986 
    90 ! variable for OPENMP parallelisation
    91 
     87! variables for OPENMP parallelisation
    9288  INTEGER,ALLOCATABLE,DIMENSION(:),SAVE :: knon_omp
    9389  REAL,ALLOCATABLE,DIMENSION(:,:),SAVE ::  buffer_omp
    … …  
    146142    ALLOCATE(unity(klon), stat = error)
    147143    sum_error = sum_error + error
    148     ALLOCATE(pctsrf_sav(klon,nbsrf), stat = error)
    149     sum_error = sum_error + error
    150144    ALLOCATE(cpl_sols(klon,2), stat = error)
    151145    sum_error = sum_error + error
    … …  
    196190       unity(ig) = ig
    197191    ENDDO
    198     pctsrf_sav = 0.
    199 
    200     cpl_sols = 0.   ; cpl_nsol = 0.  ; cpl_rain = 0.   ; cpl_snow = 0.
    201     cpl_evap = 0.   ; cpl_tsol = 0.  ; cpl_fder = 0.   ; cpl_albe = 0.
    202     cpl_taux = 0.   ; cpl_tauy = 0.  ; cpl_rriv2D = 0. ; cpl_rcoa2D = 0.
    203     cpl_rlic2D = 0. ; cpl_windsp = 0.
     192
     193!    cpl_sols = 0.   ; cpl_nsol = 0.  ; cpl_rain = 0.   ; cpl_snow = 0.
     194!    cpl_evap = 0.   ; cpl_tsol = 0.  ; cpl_fder = 0.   ; cpl_albe = 0.
     195!    cpl_taux = 0.   ; cpl_tauy = 0.  ; cpl_rriv2D = 0. ; cpl_rcoa2D = 0.
     196!    cpl_rlic2D = 0. ; cpl_windsp = 0.
    204197
    205198!*************************************************************************************
    … …  
    262255
    263256!$OMP MASTER
    264   ALLOCATE(knon_omp(0:omp_size-1))
    265   ALLOCATE(buffer_omp(klon_mpi,0:omp_size-1))       
     257    ALLOCATE(knon_omp(0:omp_size-1))
     258    ALLOCATE(buffer_omp(klon_mpi,0:omp_size-1))       
    266259!$OMP END MASTER
    267260!$OMP BARRIER
    … …  
    272265!*************************************************************************************
    273266!
    274 
    275   SUBROUTINE cpl_receive_all(itime, dtime, pctsrf)
    276 ! This subroutine reads from coupler for both ocean and seaice
     267 
     268  SUBROUTINE cpl_receive_frac(itime, dtime, pctsrf, is_modified)
     269! This subroutine receives from coupler for both ocean and seaice
    277270! 4 fields : read_sst, read_sic, read_sit and read_alb_sic.
     271! The new sea-ice-land-landice fraction is returned. The others fields
     272! are stored in this module.
     273    USE surface_data
    278274
    279275    INCLUDE "indicesol.h"
    … …  
    283279    INCLUDE "dimensions.h"
    284280
    285 ! Input arguments
     281! Arguments
    286282!************************************************************************************
    287     INTEGER, INTENT(IN)                     :: itime
    288     REAL, INTENT(IN)                        :: dtime
    289     REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf
     283    INTEGER, INTENT(IN)                        :: itime
     284    REAL, INTENT(IN)                           :: dtime
     285    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: pctsrf
     286    LOGICAL, INTENT(OUT)                       :: is_modified
    290287
    291288! Local variables
    292289!************************************************************************************
    293     INTEGER                                 :: j, ig, il_time_secs
     290    INTEGER                                 :: j, i, time_sec
    294291    INTEGER                                 :: itau_w
    295292    INTEGER, DIMENSION(iim*(jjm+1))         :: ndexcs
    296     CHARACTER(len = 20)                     :: modname = 'cpl_receive_all'
     293    CHARACTER(len = 20)                     :: modname = 'cpl_receive_frac'
    297294    CHARACTER(len = 80)                     :: abort_message
    298295    REAL, DIMENSION(klon)                   :: read_sic1D
    299     REAL, DIMENSION(iim,jj_nb,jpfldo2a)  :: tab_read_flds
    300     REAL, DIMENSION(iim,jj_nb)           :: read_sic
     296    REAL, DIMENSION(iim,jj_nb,jpfldo2a)     :: tab_read_flds
     297    REAL, DIMENSION(klon,nbsrf)             :: pctsrf_old
    301298
    302299!*************************************************************************************
    … …  
    306303!*************************************************************************************
    307304
     305    is_modified=.FALSE.
     306
     307! Check if right moment to recevie from coupler
     308    IF (MOD(itime, nexca) == 1) THEN
     309       is_modified=.TRUE.
     310 
     311       time_sec=(itime-1)*dtime
    308312#ifdef CPP_COUPLE
    309     il_time_secs=(itime-1)*dtime
    310313!$OMP MASTER
    311     CALL fromcpl(il_time_secs, tab_read_flds)
     314       CALL fromcpl(time_sec, tab_read_flds)
    312315!$OMP END MASTER
    313316#endif
    314317   
    315318! NetCDF output of received fields
    316     IF (is_sequential) THEN
    317        ndexcs(:) = 0
    318        itau_w = itau_phy + itime
    319        DO ig = 1, jpfldo2a
    320           CALL histwrite(nidcs,cl_read(ig),itau_w,tab_read_flds(:,:,ig),iim*(jjm+1),ndexcs)
    321        END DO
    322     ENDIF
    323 
     319       IF (is_sequential) THEN
     320          ndexcs(:) = 0
     321          itau_w = itau_phy + itime
     322          DO i = 1, jpfldo2a
     323             CALL histwrite(nidcs,cl_read(i),itau_w,tab_read_flds(:,:,i),iim*(jjm+1),ndexcs)
     324          END DO
     325       ENDIF
     326
     327! Save each field in a 2D array.
    324328!$OMP MASTER
    325 
    326 ! Save each field in a 2D array.
    327 
    328     IF (OPA_version=='OPA9') THEN
    329       read_sst(:,:)     = tab_read_flds(:,:,1)  ! Sea surface temperature
    330       read_sic(:,:)     = tab_read_flds(:,:,2)  ! Sea ice concentration
    331       read_sit(:,:)     = tab_read_flds(:,:,3)  ! Sea ice temperature
    332       read_alb_sic(:,:) = tab_read_flds(:,:,4)  ! Albedo at sea ice
    333     ELSE IF (OPA_version=='OPA8') THEN
    334       read_sst(:,:)     = tab_read_flds(:,:,1)  ! Sea surface temperature
    335       read_sic(:,:)     = tab_read_flds(:,:,2)  ! Sea ice concentration
    336       read_alb_sic(:,:) = tab_read_flds(:,:,3)  ! Albedo at sea ice
    337       read_sit(:,:)     = tab_read_flds(:,:,4)  ! Sea ice temperature
    338     ELSE
    339       STOP 'Bad OPA version for coupled model'
    340     ENDIF
    341 
    342 !*************************************************************************************
    343 ! Temperature and albedo are weighted with the fraction of sea-ice(read-sic)
    344 !
    345 !*************************************************************************************
    346     DO j = 1, jj_nb
    347        DO ig = 1, iim
    348           IF (ABS(1. - read_sic(ig,j)) < 0.00001) THEN
    349              read_sst(ig,j) = RTT - 1.8
    350              read_sit(ig,j) = read_sit(ig,j) / read_sic(ig,j)
    351              read_alb_sic(ig,j) = read_alb_sic(ig,j) / read_sic(ig,j)
    352           ELSE IF (ABS(read_sic(ig,j)) < 0.00001) THEN
    353              read_sst(ig,j) = read_sst(ig,j) / (1. - read_sic(ig,j))
    354              read_sit(ig,j) = read_sst(ig,j)
    355              read_alb_sic(ig,j) =  0.6
    356           ELSE
    357              read_sst(ig,j) = read_sst(ig,j) / (1. - read_sic(ig,j))
    358              read_sit(ig,j) = read_sit(ig,j) / read_sic(ig,j)
    359              read_alb_sic(ig,j) = read_alb_sic(ig,j) / read_sic(ig,j)
     329       IF (version_ocean=='nemo') THEN
     330          read_sst(:,:)     = tab_read_flds(:,:,1)  ! Sea surface temperature
     331          read_sic(:,:)     = tab_read_flds(:,:,2)  ! Sea ice concentration
     332          read_sit(:,:)     = tab_read_flds(:,:,3)  ! Sea ice temperature
     333          read_alb_sic(:,:) = tab_read_flds(:,:,4)  ! Albedo at sea ice
     334       ELSE IF (version_ocean=='opa8') THEN
     335          read_sst(:,:)     = tab_read_flds(:,:,1)  ! Sea surface temperature (multiplicated by fraction)
     336          read_sic(:,:)     = tab_read_flds(:,:,2)  ! Sea ice concentration
     337          read_alb_sic(:,:) = tab_read_flds(:,:,3)  ! Albedo at sea ice (multiplicated by fraction)
     338          read_sit(:,:)     = tab_read_flds(:,:,4)  ! Sea ice temperature (multiplicated by fraction)
     339       END IF
     340!$OMP END MASTER
     341
     342!*************************************************************************************
     343!  Transform seaice fraction (read_sic : ocean-seaice mask) into global
     344!  fraction (pctsrf : ocean-seaice-land-landice mask)
     345!
     346!*************************************************************************************
     347       CALL cpl2gath(read_sic, read_sic1D, klon, unity)
     348
     349       pctsrf_old(:,:) = pctsrf(:,:)
     350       DO i = 1, klon
     351          ! treatment only of points with ocean and/or seaice
     352          IF (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic) > 0.) THEN
     353             pctsrf(i,is_sic) = (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic)) &
     354                  * read_sic1D(i)
     355             pctsrf(i,is_oce) = (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic)) &
     356                  - pctsrf(i,is_sic)
    360357          ENDIF
    361358       ENDDO
    362     ENDDO
    363 !$OMP END MASTER
    364 
    365 !*************************************************************************************
    366 !  Transform seaice fraction, read_sic into pctsrf_sav
    367 !
    368 !*************************************************************************************
    369     CALL cpl2gath(read_sic, read_sic1D, klon, unity)
    370 
    371     DO ig = 1, klon
    372        ! treatment only of ocean and/or seaice points
    373        IF (pctsrf(ig,is_oce) > epsfra .OR. &
    374             pctsrf(ig,is_sic) > epsfra) THEN
    375           pctsrf_sav(ig,is_sic) = (pctsrf(ig,is_oce) + pctsrf(ig,is_sic)) &
    376                * read_sic1D(ig)
    377           pctsrf_sav(ig,is_oce) = (pctsrf(ig,is_oce) + pctsrf(ig,is_sic)) &
    378                - pctsrf_sav(ig,is_sic)
    379        ENDIF
    380     ENDDO
    381 
    382 !*************************************************************************************
    383 ! To avoid round up problems
    384 !
    385 !*************************************************************************************
    386     WHERE (ABS(pctsrf_sav(:,is_sic)) .LE. 2.*EPSILON(pctsrf_sav(1,is_sic)))
    387        pctsrf_sav(:,is_sic) = 0.
    388        pctsrf_sav(:,is_oce) = pctsrf(:,is_oce) + pctsrf(:,is_sic)
    389     ENDWHERE
    390     WHERE (ABS(pctsrf_sav(:,is_oce)) .LE. 2.*EPSILON(pctsrf_sav(1,is_oce)))
    391        pctsrf_sav(:,is_sic) = pctsrf(:,is_oce) + pctsrf(:,is_sic)
    392        pctsrf_sav(:,is_oce) = 0.
    393     ENDWHERE
    394     IF (MINVAL(pctsrf_sav(:,is_oce)) < 0.) THEN
    395        WRITE(*,*)'Pb fraction ocean inferieure a 0'
    396        WRITE(*,*)'au point ',MINLOC(pctsrf_sav(:,is_oce))
    397        WRITE(*,*)'valeur = ',MINVAL(pctsrf_sav(:,is_oce))
    398        abort_message = 'voir ci-dessus'
    399        CALL abort_gcm(modname,abort_message,1)
    400     ENDIF
    401     IF (MINVAL(pctsrf_sav(:,is_sic)) < 0.) THEN
    402        WRITE(*,*)'Pb fraction glace inferieure a 0'
    403        WRITE(*,*)'au point ',MINLOC(pctsrf_sav(:,is_sic))
    404        WRITE(*,*)'valeur = ',MINVAL(pctsrf_sav(:,is_sic))
    405        abort_message = 'voir ci-dessus'
    406        CALL abort_gcm(modname,abort_message,1)
    407     ENDIF
    408        
    409   END SUBROUTINE cpl_receive_all
    410 !
    411 !*************************************************************************************
    412 !
    413   SUBROUTINE cpl_receive_ocean_fields(itime, dtime, knon, knindex, pctsrf, &
    414        tsurf_new, pctsrf_oce)
    415 !
    416 ! This routine reads, if first time step in a coupling period, all fields reveived from
    417 ! coupler for all types of surfaces. It returns the fields for the ocean surface which
    418 ! are the sea surface temperature and the fraction of ocean.
    419 ! The fields are transformed into 1D arrays with valid points :
    420 ! tsurf_new(1:knon), pctsrf(1:klon).
     359
     360    END IF ! if time to receive
     361
     362  END SUBROUTINE cpl_receive_frac
     363
     364!
     365!*************************************************************************************
     366!
     367
     368  SUBROUTINE cpl_receive_ocean_fields(knon, knindex, tsurf_new)
     369!
     370! This routine returns the field for the ocean that has been read from the coupler
     371! (done earlier with cpl_receive_frac). The field is the temperature.
     372! The temperature is transformed into 1D array with valid points from index 1 to knon.
    421373!
    422374    INCLUDE "indicesol.h"
    … …  
    424376! Input arguments
    425377!*************************************************************************************
    426     INTEGER, INTENT(IN)                     :: itime
    427     REAL, INTENT(IN)                        :: dtime
    428378    INTEGER, INTENT(IN)                     :: knon
    429379    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
    430     REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf
    431380
    432381! Output arguments
    433382!*************************************************************************************
    434383    REAL, DIMENSION(klon), INTENT(OUT)      :: tsurf_new
    435     REAL, DIMENSION(klon), INTENT(OUT)      :: pctsrf_oce
    436 
    437 !*************************************************************************************
    438 ! If first time step in a coupling period receive all fields for all types
    439 ! of surfaces from coupler : read_sst, read_sit, read_alb_sic and pctsrf_sav.
    440 !
    441 !*************************************************************************************
    442 
    443     IF (MOD(itime, nexca) == 1) CALL cpl_receive_all(itime, dtime, pctsrf)
    444 
     384
     385! Local variables
     386!*************************************************************************************
     387    INTEGER               :: i
     388    REAL, DIMENSION(klon) :: sic_new
    445389
    446390!*************************************************************************************
    … …  
    449393!*************************************************************************************
    450394    CALL cpl2gath(read_sst, tsurf_new, knon, knindex)
    451     pctsrf_oce(:) = pctsrf_sav(:,is_oce)
    452 
     395    CALL cpl2gath(read_sic, sic_new, knon, knindex)
     396
     397!*************************************************************************************
     398! The fields received from the coupler have to be weighted with the fraction of ocean
     399! in relation to the total sea-ice+ocean
     400!
     401!*************************************************************************************
     402    DO i=1, knon
     403       tsurf_new(i) = tsurf_new(i)/(1. - sic_new(i))
     404    END DO
    453405
    454406  END SUBROUTINE cpl_receive_ocean_fields
    455 !
    456 !*************************************************************************************
    457 !
     407
     408!
     409!*************************************************************************************
     410!
     411
    458412  SUBROUTINE cpl_receive_seaice_fields(knon, knindex, &
    459        tsurf_new, alb_new, pctsrf_sic)
     413       tsurf_new, alb_new)
    460414!
    461415! This routine returns the fields for the seaice that have been read from the coupler
    462 ! (done earlier with cpl_receive_ocean_fields). These fields are the temperature and
     416! (done earlier with cpl_receive_frac). These fields are the temperature and
    463417! albedo at sea ice surface and fraction of sea ice.
    464 ! The fields are transformed into  1D arrays with valid points :
    465 ! tsurf_new(1:knon), alb_new(1:knon), pctsrf(1:klon).
    466 !
    467     INCLUDE "indicesol.h"
     418! The fields are transformed into 1D arrays with valid points from index 1 to knon.
     419!
    468420
    469421! Input arguments
    … …  
    476428    REAL, DIMENSION(klon), INTENT(OUT)      :: tsurf_new
    477429    REAL, DIMENSION(klon), INTENT(OUT)      :: alb_new
    478     REAL, DIMENSION(klon), INTENT(OUT)      :: pctsrf_sic
    479 
     430
     431! Local variables
     432!*************************************************************************************
     433    INTEGER               :: i
     434    REAL, DIMENSION(klon) :: sic_new
    480435
    481436!*************************************************************************************
    … …  
    485440    CALL cpl2gath(read_sit, tsurf_new, knon, knindex)
    486441    CALL cpl2gath(read_alb_sic, alb_new, knon, knindex)
    487     pctsrf_sic(:) = pctsrf_sav(:,is_sic)
     442    CALL cpl2gath(read_sic, sic_new, knon, knindex)
     443
     444!*************************************************************************************
     445! The fields received from the coupler have to be weighted with the sea-ice
     446! concentration (in relation to the total sea-ice + ocean).
     447!
     448!*************************************************************************************
     449    DO i= 1, knon
     450       tsurf_new(i) = tsurf_new(i) / sic_new(i)
     451       alb_new(i)   = alb_new(i)   / sic_new(i)
     452    END DO
    488453
    489454  END SUBROUTINE cpl_receive_seaice_fields
    … …  
    535500!*************************************************************************************
    536501    IF (MOD(itime, nexca) == 1) THEN
    537        cpl_sols(:,cpl_index) = 0.0
    538        cpl_nsol(:,cpl_index) = 0.0
    539        cpl_rain(:,cpl_index) = 0.0
    540        cpl_snow(:,cpl_index) = 0.0
    541        cpl_evap(:,cpl_index) = 0.0
    542        cpl_tsol(:,cpl_index) = 0.0
    543        cpl_fder(:,cpl_index) = 0.0
    544        cpl_albe(:,cpl_index) = 0.0
    545        cpl_taux(:,cpl_index) = 0.0
    546        cpl_tauy(:,cpl_index) = 0.0
    547        cpl_windsp(:,cpl_index) = 0.0
     502       cpl_sols(1:knon,cpl_index) = 0.0
     503       cpl_nsol(1:knon,cpl_index) = 0.0
     504       cpl_rain(1:knon,cpl_index) = 0.0
     505       cpl_snow(1:knon,cpl_index) = 0.0
     506       cpl_evap(1:knon,cpl_index) = 0.0
     507       cpl_tsol(1:knon,cpl_index) = 0.0
     508       cpl_fder(1:knon,cpl_index) = 0.0
     509       cpl_albe(1:knon,cpl_index) = 0.0
     510       cpl_taux(1:knon,cpl_index) = 0.0
     511       cpl_tauy(1:knon,cpl_index) = 0.0
     512       cpl_windsp(1:knon,cpl_index) = 0.0
    548513    ENDIF
    549514       
    … …  
    709674!*************************************************************************************
    710675    IF (MOD(itime, nexca) == 1) THEN
    711        cpl_sols(:,cpl_index) = 0.0
    712        cpl_nsol(:,cpl_index) = 0.0
    713        cpl_rain(:,cpl_index) = 0.0
    714        cpl_snow(:,cpl_index) = 0.0
    715        cpl_evap(:,cpl_index) = 0.0
    716        cpl_tsol(:,cpl_index) = 0.0
    717        cpl_fder(:,cpl_index) = 0.0
    718        cpl_albe(:,cpl_index) = 0.0
    719        cpl_taux(:,cpl_index) = 0.0
    720        cpl_tauy(:,cpl_index) = 0.0
     676       cpl_sols(1:knon,cpl_index) = 0.0
     677       cpl_nsol(1:knon,cpl_index) = 0.0
     678       cpl_rain(1:knon,cpl_index) = 0.0
     679       cpl_snow(1:knon,cpl_index) = 0.0
     680       cpl_evap(1:knon,cpl_index) = 0.0
     681       cpl_tsol(1:knon,cpl_index) = 0.0
     682       cpl_fder(1:knon,cpl_index) = 0.0
     683       cpl_albe(1:knon,cpl_index) = 0.0
     684       cpl_taux(1:knon,cpl_index) = 0.0
     685       cpl_tauy(1:knon,cpl_index) = 0.0
    721686    ENDIF
    722687       
    … …  
    893858!
    894859    INCLUDE "dimensions.h"
    895    
     860
    896861! Input varibales
    897862!*************************************************************************************
    … …  
    944909! all calculations at the different surfaces have to be done before.
    945910!   
     911    USE surface_data
    946912! Some includes
    947913!*************************************************************************************
    … …  
    962928    INTEGER                                              :: error, sum_error, j
    963929    INTEGER                                              :: itau_w
    964     INTEGER                                              :: il_time_secs
     930    INTEGER                                              :: time_sec
    965931    INTEGER, DIMENSION(iim*(jjm+1))                      :: ndexct
    966932    REAL                                                 :: Up, Down
    … …  
    994960!
    995961!*************************************************************************************
    996 !! AC >>
    997 
    998962!$OMP MASTER
    999     IF (OPA_version=='OPA9') THEN
    1000       tab_flds(:,:,7)  = cpl_windsp2D(:,:)
    1001       tab_flds(:,:,14) = cpl_sols2D(:,:,2)
    1002       tab_flds(:,:,12) = cpl_sols2D(:,:,1)
    1003       tab_flds(:,:,15) = cpl_nsol2D(:,:,2)
    1004       tab_flds(:,:,13) = cpl_nsol2D(:,:,1)
    1005       tab_flds(:,:,16) = cpl_fder2D(:,:,2)
    1006       tab_flds(:,:,11) = cpl_evap2D(:,:,2)
    1007       tab_flds(:,:,18) = cpl_rriv2D(:,:)
    1008       tab_flds(:,:,19) = cpl_rcoa2D(:,:)
    1009     ELSE IF (OPA_version=='OPA8') THEN
    1010       tab_flds(:,:,7)  = cpl_windsp2D(:,:)
    1011       tab_flds(:,:,8)  = cpl_sols2D(:,:,2)
    1012       tab_flds(:,:,9)  = cpl_sols2D(:,:,1)
    1013       tab_flds(:,:,10) = cpl_nsol2D(:,:,2)
    1014       tab_flds(:,:,11) = cpl_nsol2D(:,:,1)
    1015       tab_flds(:,:,12) = cpl_fder2D(:,:,2)
    1016       tab_flds(:,:,13) = cpl_evap2D(:,:,2)
    1017       tab_flds(:,:,14) = cpl_evap2D(:,:,1)
    1018       tab_flds(:,:,17) = cpl_rcoa2D(:,:)
    1019       tab_flds(:,:,18) = cpl_rriv2D(:,:)
    1020     ELSE
    1021       STOP 'Bad OPA version for coupled model'
    1022     ENDIF
    1023 
     963    IF (version_ocean=='nemo') THEN
     964       tab_flds(:,:,7)  = cpl_windsp2D(:,:)
     965       tab_flds(:,:,14) = cpl_sols2D(:,:,2)
     966       tab_flds(:,:,12) = cpl_sols2D(:,:,1)
     967       tab_flds(:,:,15) = cpl_nsol2D(:,:,2)
     968       tab_flds(:,:,13) = cpl_nsol2D(:,:,1)
     969       tab_flds(:,:,16) = cpl_fder2D(:,:,2)
     970       tab_flds(:,:,11) = cpl_evap2D(:,:,2)
     971       tab_flds(:,:,18) = cpl_rriv2D(:,:)
     972       tab_flds(:,:,19) = cpl_rcoa2D(:,:)
     973    ELSE IF (version_ocean=='opa8') THEN
     974       tab_flds(:,:,7)  = cpl_windsp2D(:,:)
     975       tab_flds(:,:,8)  = cpl_sols2D(:,:,2)
     976       tab_flds(:,:,9)  = cpl_sols2D(:,:,1)
     977       tab_flds(:,:,10) = cpl_nsol2D(:,:,2)
     978       tab_flds(:,:,11) = cpl_nsol2D(:,:,1)
     979       tab_flds(:,:,12) = cpl_fder2D(:,:,2)
     980       tab_flds(:,:,13) = cpl_evap2D(:,:,2)
     981       tab_flds(:,:,14) = cpl_evap2D(:,:,1)
     982       tab_flds(:,:,17) = cpl_rcoa2D(:,:)
     983       tab_flds(:,:,18) = cpl_rriv2D(:,:)
     984    END IF
     985   
    1024986!*************************************************************************************
    1025987! Transform the fraction of sub-surfaces from 1D to 2D array
    … …  
    1028990    pctsrf2D(:,:,:) = 0.
    1029991!$OMP END MASTER
    1030 
    1031992    CALL gath2cpl(pctsrf(:,is_oce), pctsrf2D(:,:,is_oce), klon, unity)
    1032993    CALL gath2cpl(pctsrf(:,is_sic), pctsrf2D(:,:,is_sic), klon, unity)
    … …  
    10401001    IF (is_omp_root) THEN
    10411002
    1042       DO j = 1, jj_nb
    1043          tmp_calv(:,j) = DOT_PRODUCT (cpl_rlic2D(1:iim,j), &
    1044               pctsrf2D(1:iim,j,is_lic)) / REAL(iim)
    1045       ENDDO
    1046    
    1047    
    1048       IF (is_parallel) THEN
     1003       DO j = 1, jj_nb
     1004          tmp_calv(:,j) = DOT_PRODUCT (cpl_rlic2D(1:iim,j), &
     1005               pctsrf2D(1:iim,j,is_lic)) / REAL(iim)
     1006       ENDDO
     1007       
     1008   
     1009       IF (is_parallel) THEN
    10491010         IF (.NOT. is_north_pole) THEN
    10501011#ifdef CPP_PARA
    … …  
    10531014#endif
    10541015         ENDIF
    1055        
     1016          
    10561017         IF (.NOT. is_south_pole) THEN
    10571018#ifdef CPP_PARA
    … …  
    10601021#endif
    10611022         ENDIF
    1062        
     1023        
    10631024         IF (.NOT. is_north_pole .AND. ii_begin /=1) THEN
    10641025            Up=Up+tmp_calv(iim,1)
    10651026            tmp_calv(:,1)=Up
    10661027         ENDIF
    1067        
     1028        
    10681029         IF (.NOT. is_south_pole .AND. ii_end /= iim) THEN
    10691030            Down=Down+tmp_calv(1,jj_nb)
    … …  
    10711032         ENDIF
    10721033      ENDIF
    1073        
    1074       IF (OPA_version=='OPA9') THEN
    1075         tab_flds(:,:,17) = tmp_calv(:,:)
    1076       ELSE IF (OPA_version=='OPA8') THEN
    1077         tab_flds(:,:,17) = tmp_calv(:,:)
    1078       ELSE
    1079         STOP 'Bad OPA version for coupled model'
    1080       ENDIF
    1081 
     1034     
     1035      IF (version_ocean=='nemo') THEN
     1036         tab_flds(:,:,17) = tmp_calv(:,:)
     1037      ELSE IF (version_ocean=='opa8') THEN
     1038         tab_flds(:,:,19) = tmp_calv(:,:)
     1039      END IF
    10821040
    10831041!*************************************************************************************
    … …  
    10891047!
    10901048!*************************************************************************************   
    1091       ! fraction oce+seaice
    1092       deno =  pctsrf2D(:,:,is_oce) + pctsrf2D(:,:,is_sic)
    1093 
    1094       IF (OPA_version=='OPA9') THEN
    1095 
    1096         tab_flds(:,:,10) = 0.0
    1097         tmp_taux(:,:)    = 0.0
    1098         tmp_tauy(:,:)    = 0.0
    1099         ! fraction oce+seaice
    1100         ! For all valid grid cells containing some fraction of ocean or sea-ice
    1101         WHERE ( deno(:,:) /= 0 )
    1102            tab_flds(:,:,10) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1103                               cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1104 
    1105            tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1106                       cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1107            tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1108                       cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1109         ENDWHERE
    1110         tab_flds(:,:,8) = (cpl_evap2D(:,:,1) - ( cpl_rain2D(:,:,1) + cpl_snow2D(:,:,1)))
    1111         tab_flds(:,:,9) = (cpl_evap2D(:,:,2) - ( cpl_rain2D(:,:,2) + cpl_snow2D(:,:,2)))
    1112 
    1113       ELSE IF (OPA_version=='OPA8') THEN
    1114 
    1115         tab_flds(:,:,15) = 0.0
    1116         tab_flds(:,:,16) = 0.0
    1117         tmp_taux(:,:)    = 0.0
    1118         tmp_tauy(:,:)    = 0.0
    1119         ! For all valid grid cells containing some fraction of ocean or sea-ice
    1120         WHERE ( deno(:,:) /= 0 )
    1121            tab_flds(:,:,15) = cpl_rain2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1122                               cpl_rain2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1123            tab_flds(:,:,16) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1124                               cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1125        
    1126            tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1127                       cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1128            tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
    1129                       cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
    1130         ENDWHERE
    1131 
    1132       ELSE
    1133         STOP 'Bad OPA version for coupled model'
    1134       ENDIF
    1135    
    1136     ENDIF ! is_omp_root 
    1137 
    1138 
    1139 ! AC <<
     1049       ! fraction oce+seaice
     1050       deno =  pctsrf2D(:,:,is_oce) + pctsrf2D(:,:,is_sic)
     1051
     1052       IF (version_ocean=='nemo') THEN
     1053          tab_flds(:,:,10) = 0.0
     1054          tmp_taux(:,:)    = 0.0
     1055          tmp_tauy(:,:)    = 0.0
     1056          ! For all valid grid cells containing some fraction of ocean or sea-ice
     1057          WHERE ( deno(:,:) /= 0 )
     1058             tab_flds(:,:,10) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1059                  cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1060             
     1061             tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1062                  cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1063             tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1064                  cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1065          ENDWHERE
     1066          tab_flds(:,:,8) = (cpl_evap2D(:,:,1) - ( cpl_rain2D(:,:,1) + cpl_snow2D(:,:,1)))
     1067          tab_flds(:,:,9) = (cpl_evap2D(:,:,2) - ( cpl_rain2D(:,:,2) + cpl_snow2D(:,:,2)))
     1068         
     1069       ELSE IF (version_ocean=='opa8') THEN
     1070          tab_flds(:,:,15) = 0.0
     1071          tab_flds(:,:,16) = 0.0
     1072          tmp_taux(:,:)    = 0.0
     1073          tmp_tauy(:,:)    = 0.0
     1074          ! For all valid grid cells containing some fraction of ocean or sea-ice
     1075          WHERE ( deno(:,:) /= 0 )
     1076             tab_flds(:,:,15) = cpl_rain2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1077                  cpl_rain2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1078             tab_flds(:,:,16) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1079                  cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1080             
     1081             tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1082                  cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1083             tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
     1084                  cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
     1085          ENDWHERE
     1086       END IF
     1087
     1088    ENDIF ! is_omp_root
     1089 
    11401090!*************************************************************************************
    11411091! Transform the wind components from local atmospheric 2D coordinates to geocentric
    … …  
    11451095
    11461096! Transform the longitudes and latitudes on 2D arrays
    1147    
    11481097    CALL gather_omp(rlon,rlon_mpi)
    11491098    CALL gather_omp(rlat,rlat_mpi)
    … …  
    12111160!
    12121161!*************************************************************************************
     1162    time_sec=(itime-1)*dtime
    12131163#ifdef CPP_COUPLE
    1214     il_time_secs=(itime-1)*dtime
    12151164!$OMP MASTER
    1216     CALL intocpl(il_time_secs, lafin, tab_flds(:,:,:))
     1165    CALL intocpl(time_sec, lafin, tab_flds(:,:,:))
    12171166!$OMP END MASTER
    12181167#endif
    … …  
    12391188!
    12401189  SUBROUTINE cpl2gath(champ_in, champ_out, knon, knindex)
    1241   USE mod_phys_lmdz_para
     1190    USE mod_phys_lmdz_para
    12421191! Cette routine ecrit un champ 'gathered' sur la grille 2D pour le passer
    12431192! au coupleur.
    … …  
    12691218!*************************************************************************************
    12701219!
    1271    
    1272 
    12731220! Transform from 2 dimensions (iim,jj_nb) to 1 dimension (klon)
    12741221!$OMP MASTER
    … …  
    12831230       champ_out(i) = temp_omp(ig)
    12841231    ENDDO
    1285  
    12861232   
    12871233  END SUBROUTINE cpl2gath

  • LMDZ4/trunk/libf/phylmd/ini_paramLMDZ_phy.h

    r879 r996  
    3131     .                "once", zstophy,zout)
    3232c
    33        CALL histdef(nid_ctesGCM, "ok_slab_sicOBS",
    34      ."ok_slab_sicOBS= 1: glace de mer observee, =0: gl.mer calculee
    35      . par le slab", "-",
    36      .                iim,jjmp1,nhori, 1,1,1, -99, 32,
    37      .                "once", zstophy,zout)
    38 c
    3933       CALL histdef(nid_ctesGCM, "type_run",
    4034     .        "Type run: 1= CLIM ou ENSP, 2= AMIP ou CFMI",

  • LMDZ4/trunk/libf/phylmd/oasis.F90

    r987 r996  
    4444  !$OMP THREADPRIVATE(out_var_id)
    4545
    46   CHARACTER(LEN=*),PARAMETER :: OPA_version='OPA9'
    4746
    4847#ifdef CPP_COUPLE
    … …  
    5857!     LF 09/2003
    5958!
     59    USE surface_data, ONLY : version_ocean
    6060    INCLUDE "dimensions.h"
    6161
    … …  
    130130!     Define symbolic name for fields exchanged from atmos to coupler,
    131131!         must be the same as (1) of the field  definition in namcouple:
    132     IF (OPA_version=='OPA9') THEN
    133       cl_writ(1)='COTAUXXU'
    134       cl_writ(2)='COTAUYYU'
    135       cl_writ(3)='COTAUZZU'
    136       cl_writ(4)='COTAUXXV'
    137       cl_writ(5)='COTAUYYV'
    138       cl_writ(6)='COTAUZZV'
    139       cl_writ(7)='COWINDSP'
    140       cl_writ(8)='COPEFWAT'
    141       cl_writ(9)='COPEFICE'
     132
     133    cl_writ(1)='COTAUXXU'
     134    cl_writ(2)='COTAUYYU'
     135    cl_writ(3)='COTAUZZU'
     136    cl_writ(4)='COTAUXXV'
     137    cl_writ(5)='COTAUYYV'
     138    cl_writ(6)='COTAUZZV'
     139    cl_writ(7)='COWINDSP'
     140   
     141    IF (version_ocean=='nemo') THEN
     142      cl_writ(8) ='COPEFWAT'
     143      cl_writ(9) ='COPEFICE'
    142144      cl_writ(10)='COTOSPSU'
    143145      cl_writ(11)='COICEVAP'
    … …  
    150152      cl_writ(18)='CRWOCERD'
    151153      cl_writ(19)='CRWOCECD'
    152     ELSE IF (OPA_version=='OPA8') THEN
    153       cl_writ(1)='COTAUXXU'
    154       cl_writ(2)='COTAUYYU'
    155       cl_writ(3)='COTAUZZU'
    156       cl_writ(4)='COTAUXXV'
    157       cl_writ(5)='COTAUYYV'
    158       cl_writ(6)='COTAUZZV'
    159       cl_writ(7)='COWINDSP'
    160       cl_writ(8)='COSHFICE'
    161       cl_writ(9)='COSHFOCE'
     154    ELSE IF (version_ocean=='opa8') THEN
     155      cl_writ(8) ='COSHFICE'
     156      cl_writ(9) ='COSHFOCE'
    162157      cl_writ(10)='CONSFICE'
    163158      cl_writ(11)='CONSFOCE'
    … …  
    170165      cl_writ(18)='CORIVFLU'
    171166      cl_writ(19)='COCALVIN'
    172     ELSE
    173       STOP 'Bad OPA version for coupled model'
    174167    ENDIF
    175168
    … …  
    178171!         must be the same as (2) of the field  definition in namcouple:
    179172!
    180     IF (OPA_version=='OPA9') THEN
    181       cl_read(1)='SISUTESW'
    182       cl_read(2)='SIICECOV'
    183       cl_read(4)='SIICEALW'
    184       cl_read(3)='SIICTEMW'
    185     ELSE IF (OPA_version=='OPA8') THEN
    186       cl_read(1)='SISUTESW'
    187       cl_read(2)='SIICECOV'
    188       cl_read(3)='SIICEALW'
    189       cl_read(4)='SIICTEMW'
    190     ELSE
    191       STOP 'Bad OPA version for coupled model'
    192     ENDIF
    193    
     173    IF (version_ocean=='nemo') THEN
     174       cl_read(1)='SISUTESW'
     175       cl_read(2)='SIICECOV'
     176       cl_read(4)='SIICEALW'
     177       cl_read(3)='SIICTEMW'
     178    ELSE IF (version_ocean=='opa8') THEN
     179       cl_read(1)='SISUTESW'
     180       cl_read(2)='SIICECOV'
     181       cl_read(3)='SIICEALW'
     182       cl_read(4)='SIICTEMW'
     183    END IF
     184
    194185    il_var_nodims(1) = 2
    195186    il_var_nodims(2) = 1

  • LMDZ4/trunk/libf/phylmd/ocean_cpl_mod.F90

    r888 r996  
    1616  PRIVATE
    1717
    18   PUBLIC :: ocean_cpl_init, ocean_cpl_get_vars, ocean_cpl_noice, ocean_cpl_ice
    19 
    20   REAL, ALLOCATABLE, DIMENSION(:), SAVE       :: tmp_flux_o
    21   !$OMP THREADPRIVATE(tmp_flux_o)
    22   REAL, ALLOCATABLE, DIMENSION(:), SAVE       :: tmp_flux_g
    23   !$OMP THREADPRIVATE(tmp_flux_g)
     18  PUBLIC :: ocean_cpl_init, ocean_cpl_noice, ocean_cpl_ice
    2419
    2520!****************************************************************************************
    … …  
    4338    CHARACTER (len = 80) :: abort_message
    4439    CHARACTER (len = 20) :: modname = 'ocean_cpl_init'
    45 
    46 
    47     ALLOCATE(tmp_flux_o(klon), stat = error)
    48     IF (error /= 0) THEN
    49        abort_message='Pb allocation tmp_flux_o'
    50        CALL abort_gcm(modname,abort_message,1)
    51     ENDIF
    52 
    53     ALLOCATE(tmp_flux_g(klon), stat = error)
    54     IF (error /= 0) THEN
    55        abort_message='Pb allocation tmp_flux_g'
    56        CALL abort_gcm(modname,abort_message,1)
    57     ENDIF
    5840
    5941! Initialize module cpl_init
    … …  
    7153       p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
    7254       petAcoef, peqAcoef, petBcoef, peqBcoef, &
    73        ps, u1_lay, v1_lay, pctsrf_in, &
     55       ps, u1_lay, v1_lay, &
    7456       radsol, snow, agesno, &
    7557       qsurf, evap, fluxsens, fluxlat, &
    76        tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     58       tsurf_new, dflux_s, dflux_l)
    7759!
    7860! This subroutine treats the "open ocean", all grid points that are not entierly covered
    … …  
    10183    REAL, DIMENSION(klon), INTENT(IN)        :: ps
    10284    REAL, DIMENSION(klon), INTENT(IN)        :: u1_lay, v1_lay
    103     REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf_in
    10485
    10586! In/Output arguments
    … …  
    11596    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    11697    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    117     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_oce
    11898
    11999! Local variables
    … …  
    122102    INTEGER, DIMENSION(1) :: iloc
    123103    REAL, DIMENSION(klon) :: cal, beta, dif_grnd
    124     REAL, DIMENSION(klon) :: zx_sl
    125104    REAL, DIMENSION(klon) :: fder_new
    126105    REAL, DIMENSION(klon) :: tsurf_cpl
    … …  
    134113
    135114!****************************************************************************************
    136 ! Receive sea-surface temperature(tsurf_cpl) and new fraction of ocean surface(pctsrf_oce)
    137 ! from coupler
    138 !
    139 !****************************************************************************************
    140     CALL cpl_receive_ocean_fields(itime, dtime, knon, knindex, pctsrf_in, &
    141          tsurf_cpl, pctsrf_oce)
     115! Receive sea-surface temperature(tsurf_cpl) from coupler
     116!
     117!****************************************************************************************
     118    CALL cpl_receive_ocean_fields(knon, knindex, tsurf_cpl)
    142119
    143120!****************************************************************************************
    … …  
    176153
    177154!****************************************************************************************
    178 ! Flux ocean-atmosphere useful for "slab" ocean but here calculated only for printing
    179 ! usage later in physiq 
    180 !
    181 !****************************************************************************************
    182     tmp_flux_o(:) = 0.0
    183     DO i=1, knon
    184        zx_sl(i) = RLVTT
    185        IF (tsurf_new(i) .LT. RTT) zx_sl(i) = RLSTT
    186        !IM     flux_o(i) = fluxsens(i)-evap(i)*zx_sl(i)
    187        !       flux_o(i) = fluxsens(i) + fluxlat(i)
    188        IF (pctsrf_oce(knindex(i)) .GT. epsfra) THEN
    189           tmp_flux_o(knindex(i)) = fluxsens(i) + fluxlat(i)
    190        ENDIF
    191     ENDDO
    192 
    193 
    194 !****************************************************************************************
    195155! Send and cumulate fields to the coupler
    196156!
    … …  
    213173       p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
    214174       petAcoef, peqAcoef, petBcoef, peqBcoef, &
    215        ps, u1_lay, v1_lay, pctsrf_in, &
     175       ps, u1_lay, v1_lay, pctsrf, &
    216176       radsol, snow, qsurf, &
    217177       alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    218        tsurf_new, dflux_s, dflux_l, pctsrf_sic)
     178       tsurf_new, dflux_s, dflux_l)
    219179!
    220180! This subroutine treats the ocean where there is ice. The subroutine first receives
    … …  
    244204    REAL, DIMENSION(klon), INTENT(IN)        :: ps
    245205    REAL, DIMENSION(klon), INTENT(IN)        :: u1_lay, v1_lay
    246     REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf_in
     206    REAL, DIMENSION(klon,nbsrf), INTENT(IN)  :: pctsrf
    247207
    248208! In/output arguments
    … …  
    258218    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    259219    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    260     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_sic
    261220
    262221! Local variables
    … …  
    277236
    278237!****************************************************************************************
    279 ! Receive ocean temperature(tsurf_cpl), albedo(alb_cpl) and new fraction of
    280 ! seaice(pctsrf_sic) from coupler
     238! Receive ocean temperature(tsurf_cpl) and albedo(alb_new) from coupler
    281239!
    282240!****************************************************************************************
    283241
    284242    CALL cpl_receive_seaice_fields(knon, knindex, &
    285          tsurf_cpl, alb_cpl, pctsrf_sic)
     243         tsurf_cpl, alb_cpl)
    286244
    287245    alb1_new(1:knon) = alb_cpl(1:knon)
    … …  
    308266    CALL calcul_wind_flux(knon, dtime, taux, tauy)
    309267   
    310 !****************************************************************************************
    311 ! Flux ocean-atmosphere useful for "slab" ocean but here calculated only for printing
    312 ! usage later in physiq 
    313 !
    314 !  IM: faire dependre le coefficient de conduction de la glace de mer
    315 !      de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
    316 !      actuel correspond a 3m de glace de mer, cf. L.Li
    317 !
    318 !****************************************************************************************
    319     tmp_flux_g(:) = 0.0
    320     DO i = 1, knon
    321        IF (cal(i) .GT. 1.0e-15 .AND. pctsrf_sic(knindex(i)) .GT. epsfra) &
    322             tmp_flux_g(knindex(i)) = (tsurf_new(i)-t_grnd) * &
    323             dif_grnd(i) * RCPD/cal(i)
    324     ENDDO
    325    
     268   
    326269!****************************************************************************************
    327270! Calculate fder : flux derivative (sensible and latente)
    … …  
    344287
    345288    CALL cpl_send_seaice_fields(itime, dtime, knon, knindex, &
    346        pctsrf_in, lafin, rlon, rlat, &
     289       pctsrf, lafin, rlon, rlat, &
    347290       swnet, lwnet, fluxlat, fluxsens, &
    348291       precip_rain, precip_snow, evap, tsurf_new, fder_new, alb1, taux, tauy)
    … …  
    353296!****************************************************************************************
    354297!
    355   SUBROUTINE ocean_cpl_get_vars(flux_o, flux_g)
    356 
    357 ! This subroutine returns variables private in this module to an external
    358 ! routine (physiq).
    359 
    360     REAL, DIMENSION(klon), INTENT(OUT) :: flux_o
    361     REAL, DIMENSION(klon), INTENT(OUT) :: flux_g
    362 
    363 ! Set the output variables
    364     flux_o(:) = tmp_flux_o(:)
    365     flux_g(:) = tmp_flux_g(:)
    366 
    367   END SUBROUTINE ocean_cpl_get_vars
    368 !
    369 !****************************************************************************************
    370 !
    371298END MODULE ocean_cpl_mod

  • LMDZ4/trunk/libf/phylmd/ocean_forced_mod.F90

    r888 r996  
    1414  IMPLICIT NONE
    1515
    16   REAL, ALLOCATABLE, DIMENSION(:), SAVE, PRIVATE   :: tmp_flux_o, tmp_flux_g
    17   !$OMP THREADPRIVATE(tmp_flux_o,tmp_flux_g)
    18 
    1916CONTAINS
    20 !
    21 !****************************************************************************************
    22 !
    23   SUBROUTINE ocean_forced_init
    24 ! Allocate fields needed for this module
    25 !   
    26     INTEGER              :: error
    27     CHARACTER (len = 80) :: abort_message
    28     CHARACTER (len = 20) :: modname = 'ocean_forced_init'
    29 !****************************************************************************************
    30 
    31     ALLOCATE(tmp_flux_o(1:klon), stat = error)
    32     IF (error /= 0) THEN
    33        abort_message='Pb allocation tmp_flux_o'
    34        CALL abort_gcm(modname,abort_message,1)
    35     ENDIF
    36 
    37     ALLOCATE(tmp_flux_g(1:klon), stat = error)
    38     IF (error /= 0) THEN
    39        abort_message='Pb allocation tmp_flux_g'
    40        CALL abort_gcm(modname,abort_message,1)
    41     ENDIF   
    42 
    43 
    44   END SUBROUTINE ocean_forced_init
    45 !
    46 !****************************************************************************************
    47 !****************************************************************************************
    48 !
    49   SUBROUTINE ocean_forced_final
    50 ! Allocate fields needed for this module
    51 !   
    52     INTEGER              :: error
    53     CHARACTER (len = 80) :: abort_message
    54     CHARACTER (len = 20) :: modname = 'ocean_forced_init'
    55 !****************************************************************************************
    56 
    57     DEALLOCATE(tmp_flux_o)
    58     DEALLOCATE(tmp_flux_g)
    59 
    60 
    61   END SUBROUTINE ocean_forced_final
    6217!
    6318!****************************************************************************************
    … …  
    7126       radsol, snow, agesno, &
    7227       qsurf, evap, fluxsens, fluxlat, &
    73        tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     28       tsurf_new, dflux_s, dflux_l)
    7429!
    7530! This subroutine treats the "open ocean", all grid points that are not entierly covered
    7631! by ice.
    77 ! The routine reads data from climatologie file and does some calculations at the
     32! The routine receives data from climatologie file limit.nc and does some calculations at the
    7833! surface.
    7934!
     35    USE limit_read_mod
    8036    INCLUDE "indicesol.h"
    8137    INCLUDE "YOMCST.h"
    … …  
    10864    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    10965    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l
    110     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_oce
    11166
    11267! Local variables
    … …  
    11671    REAL, DIMENSION(klon)       :: alb_neig, tsurf_lim, zx_sl
    11772    LOGICAL                     :: check=.FALSE.
    118     REAL, DIMENSION(klon,nbsrf) :: pctsrf_lim
    11973
    12074!****************************************************************************************
    … …  
    12579!****************************************************************************************
    12680! 1)   
    127 ! Read from climatologie file SST and fraction of sub-surfaces
    128 !
    129 !****************************************************************************************
    130 ! Get from file tsurf_lim and pctsrf_lim
    131     CALL interfoce_lim(itime, dtime, jour, &
    132          knon, knindex, &
    133          debut, &
    134          tsurf_lim, pctsrf_lim)
    135    
     81! Read sea-surface temperature from file limit.nc
     82!
     83!****************************************************************************************
     84    CALL limit_read_sst(knon,knindex,tsurf_lim)
     85
    13686!****************************************************************************************
    13787! 2)
    … …  
    154104         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
    155105
    156 !****************************************************************************************
    157 ! 3)
    158 ! Calculate tmp_flux_o
    159 !
    160 !****************************************************************************************
    161 !IM: flux ocean-atmosphere utile pour le "slab" ocean
    162 ! The flux are written to hist file
    163     tmp_flux_o(:) = 0.0
    164     DO i=1, knon
    165        zx_sl(i) = RLVTT
    166        IF (tsurf_new(i) .LT. RTT) zx_sl(i) = RLSTT
    167 
    168 !IM     flux_o(i) = fluxsens(i)-evap(i)*zx_sl(i)
    169 !       flux_o(i) = fluxsens(i) + fluxlat(i)
    170        IF (pctsrf_lim(knindex(i),is_oce) .GT. epsfra) THEN
    171           tmp_flux_o(knindex(i)) = fluxsens(i) + fluxlat(i)
    172        ENDIF
    173     ENDDO
    174 
    175 !****************************************************************************************
    176 ! 4)
    177 ! Return the new values for the ocean fraction
    178 !
    179 !****************************************************************************************
    180 
    181     pctsrf_oce(:) = pctsrf_lim(:,is_oce)
    182106
    183107  END SUBROUTINE ocean_forced_noice
    … …  
    192116       radsol, snow, qsol, agesno, tsoil, &
    193117       qsurf, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    194        tsurf_new, dflux_s, dflux_l, pctsrf_sic)
     118       tsurf_new, dflux_s, dflux_l)
    195119!
    196120! This subroutine treats the ocean where there is ice.
    197121! The routine reads data from climatologie file and does flux calculations at the
    198122! surface.
    199 !   
     123!
     124    USE limit_read_mod
     125
    200126    INCLUDE "indicesol.h"
    201127    INCLUDE "dimsoil.h"
    … …  
    219145    REAL, DIMENSION(klon), INTENT(IN)    :: u1_lay, v1_lay
    220146
    221 
    222147! In/Output arguments
    223148!****************************************************************************************
    … …  
    226151    REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
    227152    REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
    228 
    229153
    230154! Output arguments
    … …  
    236160    REAL, DIMENSION(klon), INTENT(OUT)            :: tsurf_new
    237161    REAL, DIMENSION(klon), INTENT(OUT)            :: dflux_s, dflux_l     
    238     REAL, DIMENSION(klon), INTENT(OUT)            :: pctsrf_sic
    239 
    240162
    241163! Local variables
    … …  
    248170    REAL, DIMENSION(klon)       :: alb_neig, tsurf_tmp
    249171    REAL, DIMENSION(klon)       :: soilcap, soilflux
    250     REAL, DIMENSION(klon,nbsrf) :: pctsrf_lim
    251172
    252173!****************************************************************************************
    … …  
    256177
    257178!****************************************************************************************
    258 ! 1)   
    259 ! Read from climatologie file SST and fraction of sub-surfaces
    260 !
    261 !****************************************************************************************
    262     CALL interfoce_lim(itime, dtime, jour, &
    263          knon, knindex, &
    264          debut, &
    265          tsurf_tmp, pctsrf_lim)
    266  
    267     DO i = 1, knon
    268        tsurf_tmp(i) = tsurf_in(i)
    269        IF (pctsrf_lim(knindex(i),is_sic) < EPSFRA) THEN
    270           snow(i) = 0.0
    271           tsurf_tmp(i) = RTT - 1.8
    272           IF (soil_model) tsoil(i,:) = RTT -1.8
    273        ENDIF
    274     ENDDO
    275 
    276 !****************************************************************************************
    277 ! 2)
     179! 1)
    278180! Flux calculation : tsurf_new, evap, fluxlat, fluxsens,
    279181!                    dflux_s, dflux_l and qsurf
    280182!****************************************************************************************
     183    tsurf_tmp(:) = tsurf_in(:)
    281184
    282185! calculate the parameters cal, beta, capsol and dif_grnd
    283 
    284186    CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, capsol, dif_grnd)
     187
    285188   
    286189    IF (soil_model) THEN
    287190! update tsoil and calculate soilcap and soilflux
    288        CALL soil(dtime, is_sic, knon,snow, tsurf_tmp, tsoil,soilcap, soilflux)
     191       CALL soil(dtime, is_sic, knon, snow, tsurf_tmp, tsoil,soilcap, soilflux)
    289192       cal(1:knon) = RCPD / soilcap(1:knon)
    290193       radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)
    … …  
    305208
    306209!****************************************************************************************
    307 ! 3)
     210! 2)
    308211! Calculations due to snow and runoff
    309212!
    … …  
    327230    alb2_new(:) = alb1_new(:)
    328231
    329 !****************************************************************************************
    330 ! 4)
    331 ! Calculate tmp_flux_g
    332 !
    333 !****************************************************************************************
    334 !
    335     tmp_flux_g(:) = 0.0
    336     DO i = 1, knon
    337 !IM: faire dependre le coefficient de conduction de la glace de mer
    338 !    de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
    339 !    actuel correspond a 3m de glace de mer, cf. L.Li
    340 !
    341        IF (cal(i) .GT. 1.0e-15 .AND. pctsrf_lim(knindex(i),is_sic) .GT. epsfra) &
    342             tmp_flux_g(knindex(i)) = (tsurf_new(i)-t_grnd) * dif_grnd(i) *RCPD/cal(i)
    343        
    344     ENDDO
    345 
    346 
    347 !****************************************************************************************
    348 ! 5)
    349 ! Return the new values for the seaice fraction
    350 !
    351 !****************************************************************************************
    352 
    353     pctsrf_sic(:) = pctsrf_lim(:,is_sic)
    354 
    355232  END SUBROUTINE ocean_forced_ice
    356233!
    357234!****************************************************************************************
    358235!
    359   SUBROUTINE ocean_forced_get_vars(flux_o, flux_g)
    360 ! Get some variables from module oceanforced.
    361 ! This subroutine returns variables to a external routine
    362 
    363     REAL, DIMENSION(klon), INTENT(OUT) :: flux_o
    364     REAL, DIMENSION(klon), INTENT(OUT) :: flux_g
    365 
    366 ! Initialize the output variables
    367     flux_o(:) = tmp_flux_o(:)
    368     flux_g(:) = tmp_flux_g(:)
    369 
    370   END SUBROUTINE ocean_forced_get_vars
    371 !
    372 !****************************************************************************************
    373 !
    374236END MODULE ocean_forced_mod
    375237

  • LMDZ4/trunk/libf/phylmd/ocean_slab_mod.F90

    r888 r996  
    77! "ocean=slab".
    88!
    9   USE surface_data,     ONLY : tau_gl, calice, calsno
     9  USE surface_data
    1010  USE fonte_neige_mod,  ONLY : fonte_neige
    1111  USE calcul_fluxs_mod, ONLY : calcul_fluxs
    … …  
    1313 
    1414  IMPLICIT NONE
    15 
    16   INTEGER, PRIVATE, SAVE                           :: lmt_pas, julien, idayvrai
    17   !$OMP THREADPRIVATE(lmt_pas,julien,idayvrai)
    18   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_seaice
    19   !$OMP THREADPRIVATE(tmp_seaice)
    20   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_tslab_loc
    21   !$OMP THREADPRIVATE(tmp_tslab_loc)
    22   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: slab_bils
    23   !$OMP THREADPRIVATE(slab_bils)
    24   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE , SAVE  :: lmt_bils
    25   !$OMP THREADPRIVATE(lmt_bils)
    26   REAL, ALLOCATABLE, DIMENSION(:,:), PRIVATE, SAVE :: tmp_pctsrf_slab
    27   !$OMP THREADPRIVATE(tmp_pctsrf_slab)
    28   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_tslab
    29   !$OMP THREADPRIVATE(tmp_tslab)
    30   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_radsol
    31   !$OMP THREADPRIVATE(tmp_radsol)
    32   REAL, ALLOCATABLE, DIMENSION(:), PRIVATE, SAVE   :: tmp_flux_o, tmp_flux_g
    33   !$OMP THREADPRIVATE(tmp_flux_o,tmp_flux_g)
    34   LOGICAL, PRIVATE, SAVE                           :: check = .FALSE.
    35   !$OMP THREADPRIVATE(check)
     15  PRIVATE
     16  PUBLIC :: ocean_slab_frac, ocean_slab_noice
    3617
    3718CONTAINS
    … …  
    3920!****************************************************************************************
    4021!
    41   SUBROUTINE ocean_slab_init(dtime, tslab_rst, seaice_rst, pctsrf_rst)
     22  SUBROUTINE ocean_slab_frac(itime, dtime, jour, pctsrf, is_modified)
    4223
     24    USE limit_read_mod
    4325    INCLUDE "indicesol.h"
    44     INCLUDE "iniprint.h"
     26!    INCLUDE "clesphys.h"
    4527
    46 ! Input variables
     28! Arguments
    4729!****************************************************************************************
    48     REAL, INTENT(IN)                         :: dtime
    49 ! Variables read from restart file
    50     REAL, DIMENSION(klon), INTENT(IN)        :: tslab_rst         
    51     REAL, DIMENSION(klon), INTENT(IN)        :: seaice_rst
    52     REAL, DIMENSION(klon, nbsrf), INTENT(IN) :: pctsrf_rst
    53 
     30    INTEGER, INTENT(IN)                        :: itime   ! numero du pas de temps courant
     31    INTEGER, INTENT(IN)                        :: jour    ! jour a lire dans l'annee
     32    REAL   , INTENT(IN)                        :: dtime   ! pas de temps de la physique (en s)
     33    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: pctsrf  ! sub-surface fraction
     34    LOGICAL, INTENT(OUT)                       :: is_modified ! true if pctsrf is modified at this time step
    5435
    5536! Local variables
    5637!****************************************************************************************
    57     INTEGER                :: error
    5838    CHARACTER (len = 80)   :: abort_message
    59     CHARACTER (len = 20)   :: modname = 'ocean_slab_intit'
     39    CHARACTER (len = 20)   :: modname = 'ocean_slab_frac'
    6040
    6141
    62     WRITE(lunout,*) '************************'
    63     WRITE(lunout,*) 'SLAB OCEAN est actif, prenez precautions !'
    64     WRITE(lunout,*) '************************'   
     42    IF (version_ocean=='sicOBS') THEN   
     43       CALL limit_read_frac(itime, dtime, jour, pctsrf, is_modified)
     44    ELSE
     45       abort_message='Ocean slab model without forced sea-ice fractions has to be rewritten!!!'
     46       CALL abort_gcm(modname,abort_message,1)
     47! Here should sea-ice/ocean fraction either be calculated or returned if saved as a module varaiable
     48! (in the case the new fractions are calculated in ocean_slab_ice or ocean_slab_noice subroutines). 
     49    END IF
    6550
    66 ! Allocate variables initialize from restart fields
    67     ALLOCATE(tmp_tslab(klon), stat = error)
    68     IF (error /= 0) THEN
    69        abort_message='Pb allocation tmp_tslab'
    70        CALL abort_gcm(modname,abort_message,1)
    71     ENDIF
    72     tmp_tslab(:) = tslab_rst(:)
    73 
    74     ALLOCATE(tmp_tslab_loc(klon), stat = error)
    75     IF (error /= 0) THEN
    76        abort_message='Pb allocation tmp_tslab_loc'
    77        CALL abort_gcm(modname,abort_message,1)
    78     ENDIF
    79     tmp_tslab_loc(:) = tslab_rst(:)
    80 
    81     ALLOCATE(tmp_seaice(klon), stat = error)
    82     IF (error /= 0) THEN
    83        abort_message='Pb allocation tmp_seaice'
    84        CALL abort_gcm(modname,abort_message,1)
    85     ENDIF
    86     tmp_seaice(:) = seaice_rst(:)
    87 
    88     ALLOCATE(tmp_pctsrf_slab(klon,nbsrf), stat = error)
    89     IF (error /= 0) THEN
    90        abort_message='Pb allocation tmp_pctsrf_slab'
    91        CALL abort_gcm(modname,abort_message,1)
    92     ENDIF
    93     tmp_pctsrf_slab(:,:) = pctsrf_rst(:,:)
    94    
    95 ! Allocate some other variables internal in module mod_oceanslab
    96     ALLOCATE(tmp_radsol(klon), stat = error)
    97     IF (error /= 0) THEN
    98        abort_message='Pb allocation tmp_radsol'
    99        CALL abort_gcm(modname,abort_message,1)
    100     ENDIF
    101 
    102     ALLOCATE(tmp_flux_o(klon), stat = error)
    103     IF (error /= 0) THEN
    104        abort_message='Pb allocation tmp_flux_o'
    105        CALL abort_gcm(modname,abort_message,1)
    106     ENDIF
    107    
    108     ALLOCATE(tmp_flux_g(klon), stat = error)
    109     IF (error /= 0) THEN
    110        abort_message='Pb allocation tmp_flux_g'
    111        CALL abort_gcm(modname,abort_message,1)
    112     ENDIF
    113 
    114 ! a mettre un slab_bils aussi en force !!!
    115     ALLOCATE(slab_bils(klon), stat = error)
    116     IF (error /= 0) THEN
    117        abort_message='Pb allocation slab_bils'
    118        CALL abort_gcm(modname,abort_message,1)
    119     ENDIF
    120     slab_bils(:) = 0.0   
    121 
    122     ALLOCATE(lmt_bils(klon), stat = error)
    123     IF (error /= 0) THEN
    124        abort_message='Pb allocation lmt_bils'
    125        CALL abort_gcm(modname,abort_message,1)
    126     ENDIF
    127 
    128 
    129 ! pour une lecture une fois par jour   
    130     lmt_pas = NINT(86400./dtime * 1.0)
    131 
    132   END SUBROUTINE ocean_slab_init
     51  END SUBROUTINE ocean_slab_frac
    13352!
    13453!****************************************************************************************
    13554!
    13655  SUBROUTINE ocean_slab_noice( &
    137        dtime, knon, knindex, &
     56       itime, dtime, jour, knon, knindex, &
    13857       p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
    13958       petAcoef, peqAcoef, petBcoef, peqBcoef, &
    140        ps, u1_lay, v1_lay, &
     59       ps, u1_lay, v1_lay, tsurf_in, &
    14160       radsol, snow, agesno, &
    14261       qsurf, evap, fluxsens, fluxlat, &
    143        tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     62       tsurf_new, dflux_s, dflux_l, lmt_bils)
    14463
    14564    INCLUDE "indicesol.h"
    … …  
    14867! Input arguments
    14968!****************************************************************************************
     69    INTEGER, INTENT(IN)                  :: itime
     70    INTEGER, INTENT(IN)                  :: jour
    15071    INTEGER, INTENT(IN)                  :: knon
    15172    INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
    … …  
    15980    REAL, DIMENSION(klon), INTENT(IN)    :: ps
    16081    REAL, DIMENSION(klon), INTENT(IN)    :: u1_lay, v1_lay
     82    REAL, DIMENSION(klon), INTENT(IN)    :: tsurf_in
    16183
    16284! In/Output arguments
    … …  
    17294    REAL, DIMENSION(klon), INTENT(OUT)   :: tsurf_new
    17395    REAL, DIMENSION(klon), INTENT(OUT)   :: dflux_s, dflux_l     
    174     REAL, DIMENSION(klon), INTENT(OUT)   :: pctsrf_oce
     96    REAL, DIMENSION(klon), INTENT(OUT)   :: lmt_bils
    17597
    17698! Local variables
    17799!****************************************************************************************
    178     INTEGER                :: i
    179     REAL, DIMENSION(klon)  :: cal, beta, dif_grnd
    180     REAL, DIMENSION(klon)  :: tsurf_temp
     100    INTEGER               :: i
     101    REAL, DIMENSION(klon) :: cal, beta, dif_grnd
     102    REAL, DIMENSION(klon) :: lmt_bils_oce, lmt_foce, diff_sst
     103    REAL                  :: calc_bils_oce, deltat
     104    REAL, PARAMETER       :: cyang=50.0 * 4.228e+06 ! capacite calorifique volumetrique de l'eau J/(m2 K)
    181105
    182106!****************************************************************************************
    183     IF (check) WRITE(*,*)' Entering ocean_slab_noice'   
    184 
    185     tsurf_new(1:knon) = tmp_tslab(knindex(1:knon))
    186     pctsrf_oce(:)   = tmp_pctsrf_slab(:,is_oce)
    187    
    188     tsurf_temp(:) = tsurf_new(:)
    189     cal = 0.
    190     beta = 1.
    191     dif_grnd = 0.
    192     agesno(:) = 0.
     107! 1) Flux calculation
     108!
     109!****************************************************************************************
     110    cal(:)      = 0.
     111    beta(:)     = 1.
     112    dif_grnd(:) = 0.
     113    agesno(:)   = 0.
    193114   
    194115    CALL calcul_fluxs(knon, is_oce, dtime, &
    195          tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
     116         tsurf_in, p1lay, cal, beta, tq_cdrag, ps, &
    196117         precip_rain, precip_snow, snow, qsurf,  &
    197118         radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    198119         petAcoef, peqAcoef, petBcoef, peqBcoef, &
    199120         tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
    200    
    201     tmp_flux_o(:) = 0.0
    202     tmp_radsol(:) = 0.0
    203121
    204     DO i=1, knon
    205        tmp_radsol(knindex(i))=radsol(i)
    206        
    207        IF (pctsrf_oce(knindex(i)) .GT. epsfra) &
    208             tmp_flux_o(knindex(i)) = fluxsens(i) + fluxlat(i)
    209     ENDDO
    210    
     122!****************************************************************************************
     123! 2) Get global variables lmt_bils and lmt_foce from file limit_slab.nc
     124!
     125!****************************************************************************************
     126    CALL limit_slab(itime, dtime, jour, lmt_bils, lmt_foce, diff_sst)  ! global pour un processus
     127
     128    lmt_bils_oce(:) = 0.
     129    WHERE (lmt_foce > 0.)
     130       lmt_bils_oce = lmt_bils / lmt_foce ! global
     131    END WHERE
     132
     133!****************************************************************************************
     134! 3) Recalculate new temperature
     135!
     136!****************************************************************************************
     137    DO i = 1, knon
     138       calc_bils_oce = radsol(i) + fluxsens(i) + fluxlat(i)
     139       deltat        = (calc_bils_oce - lmt_bils_oce(knindex(i)))*dtime/cyang +diff_sst(knindex(i))
     140       tsurf_new(i)  = tsurf_in(i) + deltat
     141    END DO
     142
    211143  END SUBROUTINE ocean_slab_noice
    212144!
    213145!****************************************************************************************
    214146!
    215   SUBROUTINE ocean_slab_ice(   &
    216        itime, dtime, jour, knon, knindex, &
    217        debut, &
    218        tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum, &
    219        petAcoef, peqAcoef, petBcoef, peqBcoef, &
    220        ps, u1_lay, v1_lay, &
    221        radsol, snow, qsurf, qsol, agesno, tsoil, &
    222        alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    223        tsurf_new, dflux_s, dflux_l, pctsrf_sic)
    224 
    225     INCLUDE "indicesol.h"
    226     INCLUDE "dimsoil.h"
    227     INCLUDE "YOMCST.h"
    228     INCLUDE "iniprint.h"
    229     INCLUDE "clesphys.h"
    230 
    231 ! Input arguments 
    232 !****************************************************************************************
    233     INTEGER, INTENT(IN)                  :: itime, jour, knon
    234     INTEGER, DIMENSION(klon), INTENT(IN) :: knindex
    235     REAL, INTENT(IN)                     :: dtime
    236     REAL, DIMENSION(klon), INTENT(IN)    :: tsurf
    237     REAL, DIMENSION(klon), INTENT(IN)    :: p1lay
    238     REAL, DIMENSION(klon), INTENT(IN)    :: tq_cdrag
    239     REAL, DIMENSION(klon), INTENT(IN)    :: precip_rain, precip_snow
    240     REAL, DIMENSION(klon), INTENT(IN)    :: temp_air, spechum
    241     REAL, DIMENSION(klon), INTENT(IN)    :: petAcoef, peqAcoef
    242     REAL, DIMENSION(klon), INTENT(IN)    :: petBcoef, peqBcoef
    243     REAL, DIMENSION(klon), INTENT(IN)    :: ps
    244     REAL, DIMENSION(klon), INTENT(IN)    :: u1_lay, v1_lay
    245     LOGICAL, INTENT(IN)                  :: debut
    246 
    247 !In/Output arguments
    248 !****************************************************************************************
    249     REAL, DIMENSION(klon), INTENT(INOUT)          :: radsol
    250     REAL, DIMENSION(klon), INTENT(INOUT)          :: snow, qsol
    251     REAL, DIMENSION(klon), INTENT(INOUT)          :: agesno
    252     REAL, DIMENSION(klon, nsoilmx), INTENT(INOUT) :: tsoil
    253 
    254 ! Output arguments
    255 !****************************************************************************************
    256     REAL, DIMENSION(klon), INTENT(OUT)            :: qsurf
    257     REAL, DIMENSION(klon), INTENT(OUT)            :: alb1_new  ! new albedo in visible SW interval
    258     REAL, DIMENSION(klon), INTENT(OUT)            :: alb2_new  ! new albedo in near IR interval
    259     REAL, DIMENSION(klon), INTENT(OUT)            :: evap, fluxsens, fluxlat
    260     REAL, DIMENSION(klon), INTENT(OUT)            :: tsurf_new
    261     REAL, DIMENSION(klon), INTENT(OUT)            :: dflux_s, dflux_l     
    262     REAL, DIMENSION(klon), INTENT(OUT)            :: pctsrf_sic
    263 
    264 ! Local variables
    265 !****************************************************************************************
    266     INTEGER                              :: i
    267     REAL, DIMENSION(klon)                :: cal, beta, dif_grnd, capsol
    268     REAL, DIMENSION(klon)                :: alb_neig, tsurf_temp
    269     REAL, DIMENSION(klon)                :: soilcap, soilflux
    270     REAL, DIMENSION(klon)                :: zfra
    271     REAL, PARAMETER                      :: t_grnd=271.35
    272     REAL                                 :: amn, amx
    273     REAL, DIMENSION(klon)                :: tslab
    274     REAL, DIMENSION(klon)                :: seaice ! glace de mer (kg/m2)
    275     REAL, DIMENSION(klon,nbsrf)          :: pctsrf_new
    276 
    277 !****************************************************************************************
    278 
    279     IF (check) WRITE(*,*)'Entering surface_seaice, knon=',knon
    280 
    281 ! Initialization of output variables
    282     alb1_new(:) = 0.0
    283 
    284 !****************************************************************************************
    285 !
    286 !
    287 !****************************************************************************************
    288     IF ( ok_slab_sicOBS) THEN   
    289        ! glace de mer observee, lecture conditions limites
    290        CALL interfoce_lim(itime, dtime, jour, &
    291             knon, knindex, &
    292             debut, &
    293             tsurf_new, pctsrf_new)
    294 
    295        tmp_pctsrf_slab(:,:) = pctsrf_new(:,:)
    296        WRITE(lunout,*) 'jour lecture pctsrf_new sic =',jour
    297 
    298     ELSE
    299        pctsrf_new=tmp_pctsrf_slab
    300     ENDIF
    301 
    302     DO i = 1, knon
    303        tsurf_new(i) = tsurf(i)
    304        IF (pctsrf_new(knindex(i),is_sic) < EPSFRA) THEN
    305           snow(i) = 0.0
    306           tsurf_new(i) = RTT - 1.8
    307           IF (soil_model) tsoil(i,:) = RTT -1.8
    308        ENDIF
    309     ENDDO
    310    
    311     CALL calbeta(dtime, is_sic, knon, snow, qsol, beta, capsol, dif_grnd)
    312    
    313     IF (soil_model) THEN
    314        CALL soil(dtime, is_sic, knon, snow, tsurf_new, tsoil, soilcap, soilflux)
    315        cal(1:knon) = RCPD / soilcap(1:knon)
    316        radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)
    317     ELSE
    318        dif_grnd = 1.0 / tau_gl
    319        cal = RCPD * calice
    320        WHERE (snow > 0.0) cal = RCPD * calsno
    321     ENDIF
    322     tsurf_temp = tsurf_new
    323     beta = 1.0
    324 
    325 !****************************************************************************************
    326 !
    327 !
    328 !****************************************************************************************
    329     CALL calcul_fluxs(knon, is_sic, dtime, &
    330          tsurf_temp, p1lay, cal, beta, tq_cdrag, ps, &
    331          precip_rain, precip_snow, snow, qsurf,  &
    332          radsol, dif_grnd, temp_air, spechum, u1_lay, v1_lay, &
    333          petAcoef, peqAcoef, petBcoef, peqBcoef, &
    334          tsurf_new, evap, fluxlat, fluxsens, dflux_s, dflux_l)
    335    
    336     CALL fonte_neige( knon, is_sic, knindex, dtime, &
    337          tsurf_temp, precip_rain, precip_snow, &
    338          snow, qsol, tsurf_new, evap)
    339 
    340 !****************************************************************************************
    341 !     calcul albedo
    342 !
    343 !****************************************************************************************
    344     CALL albsno(klon,knon,dtime,agesno(:),alb_neig(:), precip_snow(:)) 
    345     WHERE (snow(1 : knon) .LT. 0.0001) agesno(1 : knon) = 0.
    346     zfra(1:knon) = MAX(0.0,MIN(1.0,snow(1:knon)/(snow(1:knon)+10.0)))
    347     alb1_new(1 : knon) = alb_neig(1 : knon) *zfra(1:knon) + &
    348          0.6 * (1.0-zfra(1:knon))
    349    
    350     alb2_new(:) = alb1_new(:)
    351 
    352 !
    353 !IM: flux entre l'ocean et la glace de mer pour le "slab" ocean
    354     tmp_flux_g(:) = 0.0
    355     DO i = 1, knon
    356 !
    357 !IM: faire dependre le coefficient de conduction de la glace de mer
    358 !    de l'epaisseur de la glace de mer, dans l'hypothese ou le coeff.
    359 !    actuel correspond a 3m de glace de mer, cf. L.Li
    360 !
    361        IF ((cal(i).GT.1.0e-15) .AND. (pctsrf_new(knindex(i),is_sic) .GT. epsfra)) THEN
    362           tmp_flux_g(knindex(i))=(tsurf_new(i)-t_grnd) &
    363                * dif_grnd(i) *RCPD/cal(i)
    364        ENDIF
    365 !
    366 !IM: Attention: ne pas initialiser le tmp_radsol puisque c'est deja fait sur is_oce;
    367 !IM:            tmp_radsol doit etre le flux solaire qui arrive sur l'ocean
    368 !IM:            et non pas celui qui arrive sur la glace de mer
    369     ENDDO
    370    
    371 !
    372 ! 2eme appel a interfoce pour le cumul et le passage des flux a l'ocean
    373 !
    374 
    375     IF (check) THEN
    376        amn=MIN(tmp_tslab(1),1000.)
    377        amx=MAX(tmp_tslab(1),-1000.)
    378        DO i=2, klon
    379           amn=MIN(tmp_tslab(i),amn)
    380           amx=MAX(tmp_tslab(i),amx)
    381        ENDDO
    382        
    383        WRITE(lunout,*) ' debut avant interfoce_slab min max tmp_tslab',amn,amx
    384     ENDIF !(check) THEN
    385 
    386     tslab = tmp_tslab   
    387 
    388     CALL interfoce_slab(klon, debut, itime, dtime, jour, &
    389          tmp_radsol, tmp_flux_o, tmp_flux_g, tmp_pctsrf_slab, &
    390          tslab, seaice, pctsrf_new)
    391    
    392     tmp_pctsrf_slab=pctsrf_new
    393 
    394     DO i=1, knon
    395        tmp_tslab(knindex(i))=tslab(knindex(i))
    396     ENDDO
    397      
    398        
    399 !****************************************************************************************
    400 ! Return the fraction of sea-ice
    401 ! NB! jg : Peut-etre un probleme, faut-il prend aussi tmp_pctsrf_slab(:,is_oce)???
    402 !****************************************************************************************
    403     pctsrf_sic(:) =  tmp_pctsrf_slab(:,is_sic)
    404 
    405 
    406   END SUBROUTINE ocean_slab_ice
    407 !
    408 !****************************************************************************************
    409 !
    410   SUBROUTINE interfoce_slab(klon, debut, itap, dtime, ijour, &
    411        radsol, fluxo, fluxg, pctsrf, &
    412        tslab, seaice, pctsrf_slab)
    413 !
    414 ! Cette routine calcule la temperature d'un slab ocean, la glace de mer
    415 ! et les pourcentages de la maille couverte par l'ocean libre et/ou
    416 ! la glace de mer pour un "slab" ocean de 50m
    417 !
    418 ! Conception: Laurent Li
    419 ! Re-ecriture + adaptation LMDZ4: I. Musat
    420 !
    421 ! input:
    422 !   klon         nombre total de points de grille
    423 !   debut        logical: 1er appel a la physique
    424 !   itap         numero du pas de temps
    425 !   dtime        pas de temps de la physique (en s)
    426 !   ijour        jour dans l'annee en cours
    427 !   radsol       rayonnement net au sol (LW + SW)
    428 !   fluxo        flux turbulent (sensible + latent) sur les mailles oceaniques
    429 !   fluxg        flux de conduction entre la surface de la glace de mer et l'ocean
    430 !   pctsrf       tableau des pourcentages de surface de chaque maille
    431 ! output:
    432 !   tslab        temperature de l'ocean libre
    433 !   seaice       glace de mer (kg/m2)
    434 !   pctsrf_slab  "pourcentages" (valeurs entre 0. et 1.) surfaces issus du slab
    435 
    436     INCLUDE "indicesol.h"
    437     INCLUDE "YOMCST.h"
    438     INCLUDE "iniprint.h"
    439     INCLUDE "clesphys.h"
    440 
    441 ! Input arguments
    442 !****************************************************************************************
    443     INTEGER, INTENT(IN)                       :: klon
    444     LOGICAL, INTENT(IN)                       :: debut    ! not used
    445     INTEGER, INTENT(IN)                       :: itap
    446     REAL, INTENT(IN)                          :: dtime       ! not used
    447     INTEGER, INTENT(IN)                       :: ijour
    448     REAL, DIMENSION(klon), INTENT(IN)         :: radsol
    449     REAL, DIMENSION(klon), INTENT(IN)         :: fluxo
    450     REAL, DIMENSION(klon), INTENT(IN)         :: fluxg
    451     REAL, DIMENSION(klon, nbsrf), INTENT(IN)  :: pctsrf
    452 
    453 ! Output arguments
    454 !****************************************************************************************
    455     REAL, DIMENSION(klon), INTENT(OUT)        :: tslab
    456     REAL, DIMENSION(klon), INTENT(OUT)        :: seaice ! glace de mer (kg/m2)
    457     REAL, DIMENSION(klon, nbsrf), INTENT(OUT) :: pctsrf_slab
    458 
    459 ! Local variables
    460 !****************************************************************************************
    461     REAL                    :: amn, amx
    462     REAL, PARAMETER         :: unjour=86400.
    463     REAL, PARAMETER         :: cyang=50.0 * 4.228e+06 ! capacite calorifique volumetrique de l'eau J/(m2 K)
    464     REAL, PARAMETER         :: cbing=0.334e+05        ! J/kg
    465     REAL, DIMENSION(klon)   :: siceh !hauteur de la glace de mer (m)
    466     INTEGER                 :: i
    467     REAL                    :: zz, za, zb
    468 !
    469 !****************************************************************************************
    470 !
    471     julien = MOD(ijour,360)
    472 
    473     IF (check ) THEN
    474        amn=MIN(tmp_tslab_loc(1),1000.)
    475        amx=MAX(tmp_tslab_loc(1),-1000.)
    476        DO i=2, klon
    477           amn=MIN(tmp_tslab_loc(i),amn)
    478           amx=MAX(tmp_tslab_loc(i),amx)
    479        ENDDO
    480 
    481        WRITE(lunout,*) ' debut min max tslab',amn,amx
    482        WRITE(lunout,*) ' itap,lmt_pas unjour',itap,lmt_pas,unjour
    483     ENDIF
    484 
    485     pctsrf_slab(1:klon,1:nbsrf) = pctsrf(1:klon,1:nbsrf)
    486 !
    487 ! lecture du bilan au sol lmt_bils issu d'une simulation forcee en debut de journee
    488 !
    489     IF (MOD(itap,lmt_pas) .EQ. 1) THEN
    490        ! 1er pas de temps de la journee
    491        idayvrai = ijour
    492        CALL condsurf(julien,idayvrai, lmt_bils)
    493     ENDIF
    494 
    495     DO i = 1, klon
    496        IF((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
    497             (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
    498 !
    499 ! fabriquer de la glace si congelation atteinte:
    500 !         
    501           IF (tmp_tslab_loc(i).LT.(RTT-1.8)) THEN
    502              zz =  (RTT-1.8)-tmp_tslab_loc(i)
    503              tmp_seaice(i) = tmp_seaice(i) + cyang/cbing * zz
    504              seaice(i) = tmp_seaice(i)
    505              tmp_tslab_loc(i) = RTT-1.8
    506           ENDIF
    507 !
    508 ! faire fondre de la glace si temperature est superieure a 0:
    509 !
    510           IF ((tmp_tslab_loc(i).GT.RTT) .AND. (tmp_seaice(i).GT.0.0)) THEN
    511              zz = cyang/cbing * (tmp_tslab_loc(i)-RTT)
    512              zz = MIN(zz,tmp_seaice(i))
    513              tmp_seaice(i) = tmp_seaice(i) - zz
    514              seaice(i) = tmp_seaice(i)
    515              tmp_tslab_loc(i) = tmp_tslab_loc(i) - zz*cbing/cyang
    516           ENDIF
    517 !
    518 ! limiter la glace de mer a 10 metres (10000 kg/m2)
    519 !
    520           IF(tmp_seaice(i).GT.45.) THEN
    521              tmp_seaice(i) = MIN(tmp_seaice(i),10000.0)
    522           ELSE
    523              tmp_seaice(i) = 0.
    524           ENDIF
    525           seaice(i) = tmp_seaice(i)
    526           siceh(i)=tmp_seaice(i)/1000. !en metres
    527 !
    528 ! determiner la nature du sol (glace de mer ou ocean libre):
    529 !
    530 ! on fait dependre la fraction de seaice "pctsrf(i,is_sic)"
    531 ! de l'epaisseur de seaice :
    532 ! pctsrf(i,is_sic)=1. si l'epaisseur de la glace de mer est >= a 20cm
    533 ! et pctsrf(i,is_sic) croit lineairement avec seaice de 0. a 20cm d'epaisseur
    534 !
    535 
    536           IF(.NOT.ok_slab_sicOBS) THEN
    537              pctsrf_slab(i,is_sic)=MIN(siceh(i)/0.20, &
    538                   1.-(pctsrf_slab(i,is_ter)+pctsrf_slab(i,is_lic)))
    539              pctsrf_slab(i,is_oce)=1.0 - &
    540                   (pctsrf_slab(i,is_ter)+pctsrf_slab(i,is_lic)+pctsrf_slab(i,is_sic))
    541           ELSE
    542              IF (i.EQ.1) WRITE(lunout,*) 'cas ok_slab_sicOBS TRUE : passe sur la modif.'
    543           ENDIF !(.NOT.ok_slab_sicOBS) then
    544        ENDIF !pctsrf
    545     ENDDO
    546 !
    547 ! Calculer le bilan du flux de chaleur au sol :
    548 !
    549     DO i = 1, klon
    550        za = radsol(i) + fluxo(i)
    551        zb = fluxg(i)
    552        IF((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
    553             (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
    554           slab_bils(i)=slab_bils(i)+(za*pctsrf_slab(i,is_oce) &
    555                +zb*pctsrf_slab(i,is_sic))/ FLOAT(lmt_pas)
    556        ENDIF
    557     ENDDO !klon
    558 !
    559 ! calcul tslab
    560 !
    561     IF (MOD(itap,lmt_pas).EQ.0) THEN !fin de journee
    562 !
    563 ! calcul tslab
    564        amn=MIN(tmp_tslab_loc(1),1000.)
    565        amx=MAX(tmp_tslab_loc(1),-1000.)
    566        DO i = 1, klon
    567           IF ((pctsrf_slab(i,is_oce).GT.epsfra).OR. &
    568                (pctsrf_slab(i,is_sic).GT.epsfra)) THEN
    569              tmp_tslab_loc(i) = tmp_tslab_loc(i) + &
    570                   (slab_bils(i)-lmt_bils(i)) &
    571                   /cyang*unjour
    572 
    573 ! on remet l'accumulation a 0
    574              slab_bils(i) = 0.
    575           ENDIF !pctsrf
    576 !
    577           IF (check) THEN
    578              IF(i.EQ.1) THEN 
    579                 WRITE(lunout,*) 'i tmp_tslab_loc MOD(itap,lmt_pas).EQ.0 sicINTER',i,itap, &
    580                      tmp_tslab_loc(i), tmp_seaice(i)
    581              ENDIF
    582              
    583              amn=MIN(tmp_tslab_loc(i),amn)
    584              amx=MAX(tmp_tslab_loc(i),amx)
    585           ENDIF
    586        ENDDO !klon
    587     ENDIF !(MOD(itap,lmt_pas).EQ.0) THEN
    588 
    589     IF ( check ) THEN
    590        WRITE(lunout,*) 'fin min max tslab',amn,amx
    591     ENDIF
    592 
    593     tslab  = tmp_tslab_loc
    594     seaice = tmp_seaice
    595 
    596   END SUBROUTINE interfoce_slab
    597 !
    598 !**************************************************************************************** 
    599 !
    600   SUBROUTINE ocean_slab_final(tslab_rst, seaice_rst)
    601 
    602 ! This subroutine will send to phyredem the variables concerning the slab
    603 ! ocean that should be written to restart file.
    604 
    605 !****************************************************************************************
    606 
    607     REAL, DIMENSION(klon), INTENT(OUT) :: tslab_rst
    608     REAL, DIMENSION(klon), INTENT(OUT) :: seaice_rst
    609 
    610 !****************************************************************************************
    611 ! Set the output variables
    612     tslab_rst(:)  = tmp_tslab(:)
    613 !    tslab_rst(:)  = tmp_tslab_loc(:)
    614     seaice_rst(:) = tmp_seaice(:)
    615 
    616 ! Deallocation of all variables in module
    617     DEALLOCATE(tmp_tslab, tmp_tslab_loc, tmp_pctsrf_slab)
    618     DEALLOCATE(tmp_seaice, tmp_radsol, tmp_flux_o, tmp_flux_g)
    619     DEALLOCATE(slab_bils, lmt_bils)
    620 
    621   END SUBROUTINE ocean_slab_final
    622 !
    623 !****************************************************************************************
    624 !
    625   SUBROUTINE ocean_slab_get_vars(tslab_loc, seaice_loc, flux_o_loc, flux_g_loc)
    626 ! "Get some variables from module ocean_slab_mod"
    627 ! This subroutine prints variables to a external routine
    628 
    629     REAL, DIMENSION(klon), INTENT(OUT) :: tslab_loc
    630     REAL, DIMENSION(klon), INTENT(OUT) :: seaice_loc
    631     REAL, DIMENSION(klon), INTENT(OUT) :: flux_o_loc
    632     REAL, DIMENSION(klon), INTENT(OUT) :: flux_g_loc
    633 
    634 ! Set the output variables
    635     tslab_loc(:)  = tmp_tslab(:)
    636     seaice_loc(:) = tmp_seaice(:)
    637     flux_o_loc(:) = tmp_flux_o(:)
    638     flux_g_loc(:) = tmp_flux_g(:)
    639 
    640   END SUBROUTINE ocean_slab_get_vars
    641 !
    642 !****************************************************************************************
    643 !
    644147END MODULE ocean_slab_mod

  • LMDZ4/trunk/libf/phylmd/pbl_surface_mod.F90

    r987 r996  
    1313  USE mod_phys_lmdz_para,  ONLY : mpi_size
    1414  USE ioipsl
    15   USE surface_data,        ONLY : ocean, ok_veget
     15  USE surface_data,        ONLY : type_ocean, ok_veget
    1616  USE surf_land_mod,       ONLY : surf_land
    1717  USE surf_landice_mod,    ONLY : surf_landice
    … …  
    151151!****************************************************************************************
    152152
    153     IF (ocean /= 'slab  ' .AND. ocean /= 'force ' .AND. ocean /= 'couple') THEN
     153    IF (type_ocean /= 'slab  ' .AND. type_ocean /= 'force ' .AND. type_ocean /= 'couple') THEN
    154154      WRITE(lunout,*)' *** Warning ***'
    155       WRITE(lunout,*)'Option couplage pour l''ocean = ', ocean
     155      WRITE(lunout,*)'Option couplage pour l''ocean = ', type_ocean
    156156      abort_message='option pour l''ocean non valable'
    157157      CALL abort_gcm(modname,abort_message,1)
    … …  
    187187       zxtsol,    zxfluxlat, zt2m,     qsat2m,        &
    188188       d_t,       d_q,       d_u,      d_v,           &
    189        zcoefh,    pctsrf_new,                         &
     189       zcoefh,    slab_wfbils,                        &
    190190       qsol_d,    zq2m,      s_pblh,   s_plcl,        &
    191191       s_capCL,   s_oliqCL,  s_cteiCL, s_pblT,        &
    … …  
    299299    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: u10m    ! u speed at 10m
    300300    REAL, DIMENSION(klon, nbsrf), INTENT(INOUT)     :: v10m    ! v speed at 10m
     301    REAL, DIMENSION(klon, klev+1, nbsrf), INTENT(INOUT) :: tke
    301302
    302303! Output variables
    … …  
    321322    REAL, DIMENSION(klon, klev),  INTENT(OUT)       :: d_v        ! change in v speed
    322323    REAL, DIMENSION(klon, klev),  INTENT(OUT)       :: zcoefh     ! coef for turbulent diffusion of T and Q, mean for each grid point
    323     REAL, DIMENSION(klon, nbsrf), INTENT(OUT)       :: pctsrf_new ! new sub-surface fraction
    324324
    325325! Output only for diagnostics
     326    REAL, DIMENSION(klon),        INTENT(OUT)       :: slab_wfbils! heat balance at surface only for slab at ocean points
    326327    REAL, DIMENSION(klon),        INTENT(OUT)       :: qsol_d     ! water height in the soil (mm)
    327328    REAL, DIMENSION(klon),        INTENT(OUT)       :: zq2m       ! water vapour at 2m, mean for each grid point
    … …  
    367368    REAL, DIMENSION(klon, nbsrf),INTENT(OUT)        :: q2m        ! water vapour at 2 meter height
    368369    REAL, DIMENSION(klon, klev, nbsrf), INTENT(OUT) :: flux_q     ! water vapour flux(latent flux) (kg/m**2/s)
    369 
    370 ! Input/output
    371     REAL, DIMENSION(klon, klev+1, nbsrf), INTENT(INOUT) :: tke
    372370
    373371
    … …  
    431429    REAL, DIMENSION(klon)              :: ypsref
    432430    REAL, DIMENSION(klon)              :: yevap, ytsurf_new, yalb1_new, yalb2_new
    433     REAL, DIMENSION(klon)              :: pctsrf_nsrf
    434431    REAL, DIMENSION(klon)              :: ztsol
    435432    REAL, DIMENSION(klon)              :: alb_m  ! mean albedo for whole SW interval
    … …  
    446443    REAL, DIMENSION(klon,klev+1)       :: ytke
    447444    REAL, DIMENSION(klon,nsoilmx)      :: ytsoil
    448     REAL, DIMENSION(klon,nbsrf)        :: pctsrf_pot
    449445    CHARACTER(len=80)                  :: abort_message
    450446    CHARACTER(len=20)                  :: modname = 'pbl_surface'
    … …  
    470466    REAL, DIMENSION(klon,nbsrf)        :: trmb2        ! inhibition
    471467    REAL, DIMENSION(klon,nbsrf)        :: trmb3        ! point Omega
    472     REAL, DIMENSION(klon,nbsrf)        :: zx_rh2m, zx_qsat2m
    473     REAL, DIMENSION(klon,nbsrf)        :: zx_qs1, zx_t1
    474     REAL, DIMENSION(klon,nbsrf)        :: zdelta1, zcor1
     468    REAL, DIMENSION(klon,nbsrf)        :: zx_t1
    475469    REAL, DIMENSION(klon, nbsrf)       :: alb          ! mean albedo for whole SW interval
    476470    REAL, DIMENSION(klon)              :: ylwdown      ! jg : temporary (ysollwdown)
    477471
    478 
    479 !jg+ temporary test
    480     REAL, DIMENSION(klon,klev)         :: y_flux_u_old, y_flux_v_old
    481     REAL, DIMENSION(klon,klev)         :: y_d_u_old, y_d_v_old
    482 !jg-
    483    
     472    REAL                               :: zx_qs1, zcor1, zdelta1
     473
    484474!****************************************************************************************
    485475! Declarations specifiques pour le 1D. A reprendre
    … …  
    558548    yv1 = 0.0     ; ypaprs = 0.0     ; ypplay = 0.0
    559549    ydelp = 0.0   ; yu = 0.0         ; yv = 0.0        ; yt = 0.0         
    560     yq = 0.0      ; pctsrf_new = 0.0 ; y_dflux_t = 0.0 ; y_dflux_q = 0.0
     550    yq = 0.0      ; y_dflux_t = 0.0 ; y_dflux_q = 0.0
    561551    yrugoro = 0.0 ; yu10mx = 0.0     ; yu10my = 0.0    ; ywindsp = 0.0   
    562552    d_ts = 0.0    ; yfluxlat=0.0     ; flux_t = 0.0    ; flux_q = 0.0     
    … …  
    661651! 4) Loop over different surfaces
    662652!
    663 ! All points with a possibility of the current surface are used. This is done
    664 ! to allow the sea-ice to appear or disappear. It is considered here that the
    665 ! entier domaine of the ocean possibly can contain sea-ice.
     653! Only points containing a fraction of the sub surface will be threated.
    666654!
    667655!****************************************************************************************
    668 
    669     pctsrf_pot = pctsrf
    670     pctsrf_pot(:,is_oce) = 1. - zmasq(:)
    671     pctsrf_pot(:,is_sic) = 1. - zmasq(:)
    672      
     656   
    673657    loop_nbsrf: DO nsrf = 1, nbsrf
    674658
    … …  
    677661       knon  = 0
    678662       DO i = 1, klon
    679           IF (pctsrf_pot(i,nsrf).GT.epsfra) THEN
     663          IF (pctsrf(i,nsrf) > 0.) THEN
    680664             knon = knon + 1
    681665             ni(knon) = i
    … …  
    730714             ypaprs(j,k) = paprs(i,k)
    731715             ypplay(j,k) = pplay(i,k)
    732              ydelp(j,k) = delp(i,k)
    733              ytke(j,k)=tke(i,k,nsrf)
     716             ydelp(j,k)  = delp(i,k)
     717             ytke(j,k)   = tke(i,k,nsrf)
    734718             yu(j,k) = u(i,k)
    735719             yv(j,k) = v(i,k)
    … …  
    762746       CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
    763747            ypaprs, ypplay, yu, yv, yq, yt, yts, yrugos, yqsurf,  &
    764             ycoefm, ycoefh,ytke)
     748            ycoefm, ycoefh, ytke)
    765749       
    766750!****************************************************************************************
    … …  
    817801               ysnow, yqsol, yagesno, ytsoil, &
    818802               yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, &
    819                yqsurf, ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf, &
     803               yqsurf, ytsurf_new, y_dflux_t, y_dflux_q, &
    820804               ylwdown)
    821805     
    … …  
    828812               ysnow, yqsurf, yqsol, yagesno, &
    829813               ytsoil, yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, &
    830                ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf)
     814               ytsurf_new, y_dflux_t, y_dflux_q)
    831815         
    832816       CASE(is_oce)
    833817          CALL surf_ocean(rlon, rlat, ysolsw, ysollw, yalb1, &
    834                yrugos, ywindsp, rmu0, yfder, &
     818               yrugos, ywindsp, rmu0, yfder, yts, &
    835819               itap, dtime, jour, knon, ni, &
    836820               debut, &
    … …  
    840824               ysnow, yqsurf, yagesno, &
    841825               yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, &
    842                ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf)
     826               ytsurf_new, y_dflux_t, y_dflux_q, slab_wfbils)
    843827         
    844828       CASE(is_sic)
    … …  
    852836               ysnow, yqsurf, yqsol, yagesno, ytsoil, &
    853837               yz0_new, yalb1_new, yalb2_new, yevap, yfluxsens, yfluxlat, &
    854                ytsurf_new, y_dflux_t, y_dflux_q, pctsrf_nsrf)
     838               ytsurf_new, y_dflux_t, y_dflux_q)
    855839         
    856840
    … …  
    861845       END SELECT
    862846
    863 !****************************************************************************************
    864 ! Save the fraction of this sub-surface
    865 !
    866 !****************************************************************************************
    867        pctsrf_new(:,nsrf) = pctsrf_nsrf(:)
    868847
    869848!****************************************************************************************
    … …  
    882861!****************************************************************************************
    883862! H and Q
    884 !      print *,'pbl_surface: ok_flux_surf=',ok_flux_surf
    885 !      print *,'pbl_surface: fsens flat RLVTT=',fsens,flat,RLVTT
     863!       print *,'pbl_surface: ok_flux_surf=',ok_flux_surf
     864!       print *,'pbl_surface: fsens flat RLVTT=',fsens,flat,RLVTT
    886865       if (ok_flux_surf) then
    887866          y_flux_t1(:) =  fsens
    … …  
    916895!****************************************************************************************
    917896
    918        tke(:,:,nsrf)=0.
     897       tke(:,:,nsrf) = 0.
    919898       DO k = 1, klev
    920899          DO j = 1, knon
    … …  
    922901             ycoefh(j,k) = ycoefh(j,k) * ypct(j)
    923902             ycoefm(j,k) = ycoefm(j,k) * ypct(j)
    924              y_d_t(j,k) = y_d_t(j,k) * ypct(j)
    925              y_d_q(j,k) = y_d_q(j,k) * ypct(j)
    926              y_d_u(j,k) = y_d_u(j,k) * ypct(j)
    927              y_d_v(j,k) = y_d_v(j,k) * ypct(j)
     903             y_d_t(j,k)  = y_d_t(j,k) * ypct(j)
     904             y_d_q(j,k)  = y_d_q(j,k) * ypct(j)
     905             y_d_u(j,k)  = y_d_u(j,k) * ypct(j)
     906             y_d_v(j,k)  = y_d_v(j,k) * ypct(j)
    928907
    929908             flux_t(i,k,nsrf) = y_flux_t(j,k)
    … …  
    932911             flux_v(i,k,nsrf) = y_flux_v(j,k)
    933912
    934              tke(i,k,nsrf)=ytke(j,k)
     913             tke(i,k,nsrf)    = ytke(j,k)
    935914
    936915          ENDDO
    … …  
    10211000       trmb2(:,nsrf)  = 0.        ! inhibition
    10221001       trmb3(:,nsrf)  = 0.        ! Point Omega
     1002       rh2m(:)        = 0.
     1003       qsat2m(:)      = 0.
    10231004
    10241005#undef T2m     
    10251006#define T2m     
    10261007#ifdef T2m
    1027 ! diagnostic t,q a 2m et u, v a 10m
     1008! Calculations of diagnostic t,q at 2m and u, v at 10m
    10281009
    10291010       DO j=1, knon
    … …  
    10551036          i = ni(j)
    10561037          t2m(i,nsrf)=yt2m(j)
     1038          q2m(i,nsrf)=yq2m(j)
    10571039         
    1058           q2m(i,nsrf)=yq2m(j)
    1059 
    1060 ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman
     1040          ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman
    10611041          u10m(i,nsrf)=(yu10m(j) * uzon(j))/SQRT(uzon(j)**2+vmer(j)**2)
    10621042          v10m(i,nsrf)=(yu10m(j) * vmer(j))/SQRT(uzon(j)**2+vmer(j)**2)
    1063 
    10641043       END DO
    10651044
     1045!IM Calcule de l'humidite relative a 2m (rh2m) pour diagnostique
     1046!IM Ajoute dependance type surface
     1047       IF (thermcep) THEN
     1048          DO j = 1, knon
     1049             i=ni(j)
     1050             zdelta1 = MAX(0.,SIGN(1., rtt-yt2m(j) ))
     1051             zx_qs1  = r2es * FOEEW(yt2m(j),zdelta1)/paprs(i,1)
     1052             zx_qs1  = MIN(0.5,zx_qs1)
     1053             zcor1   = 1./(1.-RETV*zx_qs1)
     1054             zx_qs1  = zx_qs1*zcor1
     1055             
     1056             rh2m(i)   = rh2m(i)   + yq2m(j)/zx_qs1 * pctsrf(i,nsrf)
     1057             qsat2m(i) = qsat2m(i) + zx_qs1  * pctsrf(i,nsrf)
     1058          END DO
     1059       END IF
    10661060
    10671061       CALL HBTM(knon, ypaprs, ypplay, &
    … …  
    10861080       
    10871081#else
    1088 ! not defined T2m
     1082! T2m not defined
    10891083! No calculation
    1090 ! Set output variables to zero
    1091        DO j = 1, knon
    1092           i = ni(j)
    1093           pblh(i,nsrf)   = 0.
    1094           plcl(i,nsrf)   = 0.
    1095           capCL(i,nsrf)  = 0.
    1096           oliqCL(i,nsrf) = 0.
    1097           cteiCL(i,nsrf) = 0.
    1098           pblT(i,nsrf)   = 0.
    1099           therm(i,nsrf)  = 0.
    1100           trmb1(i,nsrf)  = 0.
    1101           trmb2(i,nsrf)  = 0.
    1102           trmb3(i,nsrf)  = 0.
    1103        END DO
    1104        DO j = 1, knon
    1105           i = ni(j)
    1106           t2m(i,nsrf)=0.
    1107           q2m(i,nsrf)=0.
    1108           u10m(i,nsrf)=0.
    1109           v10m(i,nsrf)=0.
    1110        END DO
     1084       PRINT*,' Warning !!! No T2m calculation. Output is set to zero.'
    11111085#endif
    11121086
    … …  
    11201094! 16) Calculate the mean value over all sub-surfaces for som variables
    11211095!
    1122 !     NB!!! jg : Pour garder la convergence numerique j'utilise pctsrf_new comme c'etait
    1123 !     fait dans physiq.F mais ca devrait plutot etre pctsrf???!!!!! A verifier!
    11241096!****************************************************************************************
    11251097   
    … …  
    11291101       DO k = 1, klev
    11301102          DO i = 1, klon
    1131              zxfluxt(i,k) = zxfluxt(i,k) + flux_t(i,k,nsrf) * pctsrf_new(i,nsrf)
    1132              zxfluxq(i,k) = zxfluxq(i,k) + flux_q(i,k,nsrf) * pctsrf_new(i,nsrf)
    1133              zxfluxu(i,k) = zxfluxu(i,k) + flux_u(i,k,nsrf) * pctsrf_new(i,nsrf)
    1134              zxfluxv(i,k) = zxfluxv(i,k) + flux_v(i,k,nsrf) * pctsrf_new(i,nsrf)
     1103             zxfluxt(i,k) = zxfluxt(i,k) + flux_t(i,k,nsrf) * pctsrf(i,nsrf)
     1104             zxfluxq(i,k) = zxfluxq(i,k) + flux_q(i,k,nsrf) * pctsrf(i,nsrf)
     1105             zxfluxu(i,k) = zxfluxu(i,k) + flux_u(i,k,nsrf) * pctsrf(i,nsrf)
     1106             zxfluxv(i,k) = zxfluxv(i,k) + flux_v(i,k,nsrf) * pctsrf(i,nsrf)
    11351107          END DO
    11361108       END DO
    … …  
    11431115    ENDDO
    11441116   
    1145 
    1146     DO i = 1, klon
    1147        IF ( ABS( pctsrf_new(i, is_ter) + pctsrf_new(i, is_lic) + &
    1148             pctsrf_new(i, is_oce) + pctsrf_new(i, is_sic)  - 1.) .GT. EPSFRA) &
    1149             THEN
    1150           WRITE(*,*) 'physiq : pb sous surface au point ', i, &
    1151                pctsrf_new(i, 1 : nbsrf)
    1152        ENDIF
    1153     ENDDO
    1154 
    11551117!
    11561118! Incrementer la temperature du sol
    … …  
    11711133         
    11721134          wfbils(i,nsrf) = ( solsw(i,nsrf) + sollw(i,nsrf) &
    1173                + flux_t(i,1,nsrf) + fluxlat(i,nsrf) ) * pctsrf_new(i,nsrf)
     1135               + flux_t(i,1,nsrf) + fluxlat(i,nsrf) ) * pctsrf(i,nsrf)
    11741136          wfbilo(i,nsrf) = (evap(i,nsrf) - (rain_f(i) + snow_f(i))) * &
    1175                pctsrf_new(i,nsrf)
    1176 
    1177           zxtsol(i)    = zxtsol(i)    + ts(i,nsrf)      * pctsrf_new(i,nsrf)
    1178           zxfluxlat(i) = zxfluxlat(i) + fluxlat(i,nsrf) * pctsrf_new(i,nsrf)
     1137               pctsrf(i,nsrf)
     1138
     1139          zxtsol(i)    = zxtsol(i)    + ts(i,nsrf)      * pctsrf(i,nsrf)
     1140          zxfluxlat(i) = zxfluxlat(i) + fluxlat(i,nsrf) * pctsrf(i,nsrf)
    11791141         
    1180           zt2m(i)  = zt2m(i)  + t2m(i,nsrf)  * pctsrf_new(i,nsrf)
    1181           zq2m(i)  = zq2m(i)  + q2m(i,nsrf)  * pctsrf_new(i,nsrf)
    1182           zu10m(i) = zu10m(i) + u10m(i,nsrf) * pctsrf_new(i,nsrf)
    1183           zv10m(i) = zv10m(i) + v10m(i,nsrf) * pctsrf_new(i,nsrf)
    1184 
    1185           s_pblh(i)   = s_pblh(i)   + pblh(i,nsrf)  * pctsrf_new(i,nsrf)
    1186           s_plcl(i)   = s_plcl(i)   + plcl(i,nsrf)  * pctsrf_new(i,nsrf)
    1187           s_capCL(i)  = s_capCL(i)  + capCL(i,nsrf) * pctsrf_new(i,nsrf)
    1188           s_oliqCL(i) = s_oliqCL(i) + oliqCL(i,nsrf)* pctsrf_new(i,nsrf)
    1189           s_cteiCL(i) = s_cteiCL(i) + cteiCL(i,nsrf)* pctsrf_new(i,nsrf)
    1190           s_pblT(i)   = s_pblT(i)   + pblT(i,nsrf)  * pctsrf_new(i,nsrf)
    1191           s_therm(i)  = s_therm(i)  + therm(i,nsrf) * pctsrf_new(i,nsrf)
    1192           s_trmb1(i)  = s_trmb1(i)  + trmb1(i,nsrf) * pctsrf_new(i,nsrf)
    1193           s_trmb2(i)  = s_trmb2(i)  + trmb2(i,nsrf) * pctsrf_new(i,nsrf)
    1194           s_trmb3(i)  = s_trmb3(i)  + trmb3(i,nsrf) * pctsrf_new(i,nsrf)
     1142          zt2m(i)  = zt2m(i)  + t2m(i,nsrf)  * pctsrf(i,nsrf)
     1143          zq2m(i)  = zq2m(i)  + q2m(i,nsrf)  * pctsrf(i,nsrf)
     1144          zu10m(i) = zu10m(i) + u10m(i,nsrf) * pctsrf(i,nsrf)
     1145          zv10m(i) = zv10m(i) + v10m(i,nsrf) * pctsrf(i,nsrf)
     1146
     1147          s_pblh(i)   = s_pblh(i)   + pblh(i,nsrf)  * pctsrf(i,nsrf)
     1148          s_plcl(i)   = s_plcl(i)   + plcl(i,nsrf)  * pctsrf(i,nsrf)
     1149          s_capCL(i)  = s_capCL(i)  + capCL(i,nsrf) * pctsrf(i,nsrf)
     1150          s_oliqCL(i) = s_oliqCL(i) + oliqCL(i,nsrf)* pctsrf(i,nsrf)
     1151          s_cteiCL(i) = s_cteiCL(i) + cteiCL(i,nsrf)* pctsrf(i,nsrf)
     1152          s_pblT(i)   = s_pblT(i)   + pblT(i,nsrf)  * pctsrf(i,nsrf)
     1153          s_therm(i)  = s_therm(i)  + therm(i,nsrf) * pctsrf(i,nsrf)
     1154          s_trmb1(i)  = s_trmb1(i)  + trmb1(i,nsrf) * pctsrf(i,nsrf)
     1155          s_trmb2(i)  = s_trmb2(i)  + trmb2(i,nsrf) * pctsrf(i,nsrf)
     1156          s_trmb3(i)  = s_trmb3(i)  + trmb3(i,nsrf) * pctsrf(i,nsrf)
    11951157       END DO
    11961158    END DO
    … …  
    12051167       PRINT*,' debut apres d_ts min max ftsol(ts)',itap,amn,amx
    12061168    ENDIF
    1207 !
    1208 ! If a sub-surface does not exsist for a grid point, the mean value for all
    1209 ! sub-surfaces is distributed.
    1210 !
    1211     DO nsrf = 1, nbsrf
    1212        DO i = 1, klon
    1213           IF ((pctsrf_new(i,nsrf) .LT. epsfra) .OR. (t2m(i,nsrf).EQ.0.)) THEN
    1214              ts(i,nsrf)     = zxtsol(i)
    1215              t2m(i,nsrf)    = zt2m(i)
    1216              q2m(i,nsrf)    = zq2m(i)
    1217              u10m(i,nsrf)   = zu10m(i)
    1218              v10m(i,nsrf)   = zv10m(i)
    1219 
    1220 ! Les variables qui suivent sont plus utilise, donc peut-etre pas la peine a les mettre ajour
    1221              pblh(i,nsrf)   = s_pblh(i)
    1222              plcl(i,nsrf)   = s_plcl(i)
    1223              capCL(i,nsrf)  = s_capCL(i)
    1224              oliqCL(i,nsrf) = s_oliqCL(i)
    1225              cteiCL(i,nsrf) = s_cteiCL(i)
    1226              pblT(i,nsrf)   = s_pblT(i)
    1227              therm(i,nsrf)  = s_therm(i)
    1228              trmb1(i,nsrf)  = s_trmb1(i)
    1229              trmb2(i,nsrf)  = s_trmb2(i)
    1230              trmb3(i,nsrf)  = s_trmb3(i)
    1231           ENDIF
    1232        ENDDO
    1233     ENDDO
     1169!!$!
     1170!!$! If a sub-surface does not exsist for a grid point, the mean value for all
     1171!!$! sub-surfaces is distributed.
     1172!!$!
     1173!!$    DO nsrf = 1, nbsrf
     1174!!$       DO i = 1, klon
     1175!!$          IF ((pctsrf_new(i,nsrf) .LT. epsfra) .OR. (t2m(i,nsrf).EQ.0.)) THEN
     1176!!$             ts(i,nsrf)     = zxtsol(i)
     1177!!$             t2m(i,nsrf)    = zt2m(i)
     1178!!$             q2m(i,nsrf)    = zq2m(i)
     1179!!$             u10m(i,nsrf)   = zu10m(i)
     1180!!$             v10m(i,nsrf)   = zv10m(i)
     1181!!$
     1182!!$! Les variables qui suivent sont plus utilise, donc peut-etre pas la peine a les mettre ajour
     1183!!$             pblh(i,nsrf)   = s_pblh(i)
     1184!!$             plcl(i,nsrf)   = s_plcl(i)
     1185!!$             capCL(i,nsrf)  = s_capCL(i)
     1186!!$             oliqCL(i,nsrf) = s_oliqCL(i)
     1187!!$             cteiCL(i,nsrf) = s_cteiCL(i)
     1188!!$             pblT(i,nsrf)   = s_pblT(i)
     1189!!$             therm(i,nsrf)  = s_therm(i)
     1190!!$             trmb1(i,nsrf)  = s_trmb1(i)
     1191!!$             trmb2(i,nsrf)  = s_trmb2(i)
     1192!!$             trmb3(i,nsrf)  = s_trmb3(i)
     1193!!$          ENDIF
     1194!!$       ENDDO
     1195!!$    ENDDO
    12341196
    12351197
    … …  
    12421204    DO nsrf = 1, nbsrf
    12431205       DO i = 1, klon
    1244           zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf_new(i,nsrf)
    1245           zxsnow(i)  = zxsnow(i)  + snow(i,nsrf)  * pctsrf_new(i,nsrf)
     1206          zxqsurf(i) = zxqsurf(i) + qsurf(i,nsrf) * pctsrf(i,nsrf)
     1207          zxsnow(i)  = zxsnow(i)  + snow(i,nsrf)  * pctsrf(i,nsrf)
    12461208       END DO
    12471209    END DO
    12481210
    1249 !
    1250 !IM Calculer l'humidite relative a 2m (rh2m) pour diagnostique
    1251 !IM ajout dependance type surface
    1252     rh2m(:)   = 0.0
    1253     qsat2m(:) = 0.0
    1254 
    1255     DO i = 1, klon
    1256        DO nsrf=1, nbsrf
    1257           zx_t1(i,nsrf) = t2m(i,nsrf)
    1258           IF (thermcep) THEN
    1259              zdelta1(i,nsrf) = MAX(0.,SIGN(1.,rtt-zx_t1(i,nsrf)))
    1260              zx_qs1(i,nsrf)  = r2es * &
    1261                   FOEEW(zx_t1(i,nsrf),zdelta1(i,nsrf))/paprs(i,1)
    1262              zx_qs1(i,nsrf)  = MIN(0.5,zx_qs1(i,nsrf))
    1263              zcor1(i,nsrf)   = 1./(1.-retv*zx_qs1(i,nsrf))
    1264              zx_qs1(i,nsrf)  = zx_qs1(i,nsrf)*zcor1(i,nsrf)
    1265           END IF
    1266           zx_rh2m(i,nsrf) = q2m(i,nsrf)/zx_qs1(i,nsrf)
    1267           zx_qsat2m(i,nsrf)=zx_qs1(i,nsrf)
    1268           rh2m(i) = rh2m(i)+zx_rh2m(i,nsrf)*pctsrf_new(i,nsrf)
    1269           qsat2m(i)=qsat2m(i)+zx_qsat2m(i,nsrf)*pctsrf_new(i,nsrf)
    1270        END DO
    1271     END DO
    12721211
    12731212! Some of the module declared variables are returned for printing in physiq.F
    … …  
    13221261! 
    13231262!****************************************************************************************
     1263!
     1264  SUBROUTINE pbl_surface_newfrac(itime, pctsrf_new, pctsrf_old, tsurf, alb1, alb2, u10m, v10m, tke)
     1265
     1266    ! Give default values where new fraction has appread
     1267
     1268    INCLUDE "indicesol.h"
     1269    INCLUDE "dimsoil.h"
     1270    INCLUDE "clesphys.h"
     1271
     1272! Input variables
     1273!****************************************************************************************
     1274    INTEGER, INTENT(IN)                     :: itime
     1275    REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf_new, pctsrf_old
     1276
     1277! InOutput variables
     1278!****************************************************************************************
     1279    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: tsurf
     1280    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: alb1, alb2
     1281    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT)        :: u10m, v10m
     1282    REAL, DIMENSION(klon,klev+1,nbsrf), INTENT(INOUT) :: tke
     1283
     1284! Local variables
     1285!****************************************************************************************
     1286    INTEGER           :: nsrf, nsrf_comp1, nsrf_comp2, nsrf_comp3, i
     1287    CHARACTER(len=80) :: abort_message
     1288    CHARACTER(len=20) :: modname = 'pbl_surface_newfrac'
     1289    INTEGER, DIMENSION(nbsrf) :: nfois=0, mfois=0, pfois=0
     1290    LOGICAL           :: debug=.FALSE.
     1291!
     1292! All at once !!
     1293!****************************************************************************************
     1294   
     1295    DO nsrf = 1, nbsrf
     1296       ! First decide complement sub-surfaces
     1297       SELECT CASE (nsrf)
     1298       CASE(is_oce)
     1299          nsrf_comp1=is_sic
     1300          nsrf_comp2=is_ter
     1301          nsrf_comp3=is_lic
     1302       CASE(is_sic)
     1303          nsrf_comp1=is_oce
     1304          nsrf_comp2=is_ter
     1305          nsrf_comp3=is_lic
     1306       CASE(is_ter)
     1307          nsrf_comp1=is_lic
     1308          nsrf_comp2=is_oce
     1309          nsrf_comp3=is_sic
     1310       CASE(is_lic)
     1311          nsrf_comp1=is_ter
     1312          nsrf_comp2=is_oce
     1313          nsrf_comp3=is_sic
     1314       END SELECT
     1315
     1316       ! Initialize all new fractions
     1317       DO i=1, klon
     1318          IF (pctsrf_new(i,nsrf) > 0. .AND. pctsrf_old(i,nsrf) == 0.) THEN
     1319
     1320             IF (pctsrf_old(i,nsrf_comp1) > 0.) THEN
     1321                ! Use the complement sub-surface, keeping the continents unchanged
     1322                qsurf(i,nsrf) = qsurf(i,nsrf_comp1)
     1323                evap(i,nsrf)  = evap(i,nsrf_comp1)
     1324                rugos(i,nsrf) = rugos(i,nsrf_comp1)
     1325                tsurf(i,nsrf) = tsurf(i,nsrf_comp1)
     1326                alb1(i,nsrf)  = alb1(i,nsrf_comp1)
     1327                alb2(i,nsrf)  = alb2(i,nsrf_comp1)
     1328                u10m(i,nsrf)  = u10m(i,nsrf_comp1)
     1329                v10m(i,nsrf)  = v10m(i,nsrf_comp1)
     1330                tke(i,:,nsrf) = tke(i,:,nsrf_comp1)
     1331                mfois(nsrf) = mfois(nsrf) + 1
     1332             ELSE
     1333                ! The continents have changed. The new fraction receives the mean sum of the existent fractions
     1334                qsurf(i,nsrf) = qsurf(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + qsurf(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3)
     1335                evap(i,nsrf)  = evap(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + evap(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1336                rugos(i,nsrf) = rugos(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + rugos(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3)
     1337                tsurf(i,nsrf) = tsurf(i,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + tsurf(i,nsrf_comp3)*pctsrf_old(i,nsrf_comp3)
     1338                alb1(i,nsrf)  = alb1(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb1(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1339                alb2(i,nsrf)  = alb2(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + alb2(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1340                u10m(i,nsrf)  = u10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + u10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1341                v10m(i,nsrf)  = v10m(i,nsrf_comp2) *pctsrf_old(i,nsrf_comp2) + v10m(i,nsrf_comp3) *pctsrf_old(i,nsrf_comp3)
     1342                tke(i,:,nsrf) = tke(i,:,nsrf_comp2)*pctsrf_old(i,nsrf_comp2) + tke(i,:,nsrf_comp3)*pctsrf_old(i,nsrf_comp3)
     1343           
     1344                ! Security abort. This option has never been tested. To test, comment the following line.
     1345!                abort_message='The fraction of the continents have changed!'
     1346!                CALL abort_gcm(modname,abort_message,1)
     1347                nfois(nsrf) = nfois(nsrf) + 1
     1348             END IF
     1349             snow(i,nsrf)     = 0.
     1350             agesno(i,nsrf)   = 0.
     1351             ftsoil(i,:,nsrf) = tsurf(i,nsrf)
     1352          ELSE
     1353             pfois(nsrf) = pfois(nsrf)+ 1
     1354          END IF
     1355       END DO
     1356       
     1357    END DO
     1358
     1359    IF (debug) THEN
     1360       print*,'itime=,',itime, 'Pas de nouveau fraction',pfois,'fois'
     1361       print*,'itime=,',itime, 'The fraction of the continents have changed',nfois,'fois'
     1362       print*,'itime=,',itime, 'The fraction ocean-seaice has changed',mfois,'fois'
     1363    END IF
     1364
     1365  END SUBROUTINE pbl_surface_newfrac
     1366
    13241367! 
     1368!****************************************************************************************
     1369! 
    13251370
    13261371END MODULE pbl_surface_mod

  • LMDZ4/trunk/libf/phylmd/phyetat0.F

    r987 r996  
    55c
    66      SUBROUTINE phyetat0 (fichnom,
    7      .           ocean_in, ok_veget_in,
    87     .           clesphy0,
    98     .           tabcntr0)
    … …  
    1312      USE mod_phys_lmdz_para
    1413      USE iophy
    15       USE ocean_slab_mod,   ONLY : ocean_slab_init
    1614      USE ocean_cpl_mod,    ONLY : ocean_cpl_init
    17       USE ocean_forced_mod, ONLY : ocean_forced_init
    1815      USE fonte_neige_mod,  ONLY : fonte_neige_init
    1916      USE pbl_surface_mod,  ONLY : pbl_surface_init
    20       USE surface_data,     ONLY : ocean, ok_veget
     17      USE surface_data,     ONLY : type_ocean
    2118      USE phys_state_var_mod
    2219
    … …  
    6865      REAL wake_fip_glo(klon_glo)
    6966      REAL tsoil_p(klon,nsoilmx,nbsrf)
    70       REAL tslab_p(klon), seaice_p(klon)
    7167      REAL qsurf_p(klon,nbsrf)
    7268      REAL qsol_p(klon)
    … …  
    8379      REAL zmasq_glo(klon_glo)
    8480      REAL tsoil(klon_glo,nsoilmx,nbsrf)
    85 cIM "slab" ocean
    86       REAL tslab(klon_glo), seaice(klon_glo)
    8781      REAL qsurf(klon_glo,nbsrf)
    8882      REAL qsol(klon_glo)
    … …  
    492486      ENDDO
    493487      ENDDO
    494 c
    495 cIM "slab" ocean
    496 c
    497 c Lecture de tslab (pour slab ocean seulement):     
    498 c
    499       IF (ocean_in .eq. 'slab  ') then
    500         ierr = NF_INQ_VARID (nid, "TSLAB", nvarid)
    501         IF (ierr.NE.NF_NOERR) THEN
    502           PRINT*, "phyetat0: Le champ <TSLAB> est absent"
    503           CALL abort
    504         ENDIF
    505 #ifdef NC_DOUBLE
    506         ierr = NF_GET_VAR_DOUBLE(nid, nvarid, tslab)
    507 #else
    508         ierr = NF_GET_VAR_REAL(nid, nvarid, tslab)
    509 #endif
    510         IF (ierr.NE.NF_NOERR) THEN
    511           PRINT*, "phyetat0: Lecture echouee pour <TSLAB>"
    512           CALL abort
    513         ENDIF
    514         xmin = 1.0E+20
    515         xmax = -1.0E+20
    516         DO i = 1, klon_glo
    517           xmin = MIN(tslab(i),xmin)
    518           xmax = MAX(tslab(i),xmax)
    519         ENDDO
    520         PRINT*,'Min, Max tslab (utilise si OCEAN=slab )', xmin, xmax
    521 c
    522 c Lecture de seaice (pour slab ocean seulement):
    523 c
    524         ierr = NF_INQ_VARID (nid, "SEAICE", nvarid)
    525         IF (ierr.NE.NF_NOERR) THEN
    526           PRINT*, "phyetat0: Le champ <SEAICE> est absent"
    527           CALL abort
    528         ENDIF
    529 #ifdef NC_DOUBLE
    530         ierr = NF_GET_VAR_DOUBLE(nid, nvarid, seaice)
    531 #else
    532         ierr = NF_GET_VAR_REAL(nid, nvarid, seaice)
    533 #endif
    534         IF (ierr.NE.NF_NOERR) THEN
    535           PRINT*, "phyetat0: Lecture echouee pour <SEAICE>"
    536           CALL abort
    537         ENDIF
    538         xmin = 1.0E+20
    539         xmax = -1.0E+20
    540         DO i = 1, klon_glo
    541           xmin = MIN(seaice(i),xmin)
    542           xmax = MAX(seaice(i),xmax)
    543         ENDDO
    544         PRINT*,'Masse de la glace de mer (utilise si OCEAN=slab)',
    545      $  xmin, xmax
    546       ELSE
    547         tslab = 0.
    548         seaice = 0.
    549       ENDIF
    550488c
    551489c Lecture de l'humidite de l'air juste au dessus du sol:
    … …  
    17621700      call Scatter( wake_fip_glo, wake_fip)
    17631701      call Scatter( tsoil,tsoil_p)
    1764       call Scatter( tslab,tslab_p)
    1765       call Scatter( seaice,seaice_p)
    17661702      call Scatter( qsurf,qsurf_p)
    17671703      call Scatter( qsol,qsol_p)
    … …  
    17961732
    17971733c
    1798 c Initilalize variables in module surface_data
    1799 c
    1800       ok_veget = ok_veget_in
    1801       ocean    = ocean_in
    1802 c
    18031734c Initialize module pbl_surface_mod
    18041735c
    … …  
    18061737     $     evap_p, frugs_p, agesno_p, tsoil_p)
    18071738
    1808 c Initialize ocean module according to ocean type
    1809       IF ( ocean == 'slab' ) THEN
    1810 c        initilalize module ocean_slab_init
    1811          CALL ocean_slab_init(dtime, tslab_p, seaice_p, pctsrf)
    1812       ELSEIF ( ocean == 'couple' ) THEN
    1813 c        initilalize module ocean_cpl_init
     1739c Initialize module ocean_cpl_mod for the case of coupled ocean
     1740      IF ( type_ocean == 'couple' ) THEN
    18141741         CALL ocean_cpl_init(dtime, rlon, rlat)
    1815       ELSE
    1816 c        initilalize module ocean_forced_init
    1817          CALL ocean_forced_init
    18181742      ENDIF
    18191743c

  • LMDZ4/trunk/libf/phylmd/phyredem.F

    r975 r996  
    88      USE mod_grid_phy_lmdz
    99      USE mod_phys_lmdz_para
    10       USE ocean_slab_mod,   ONLY : ocean_slab_final
    1110      USE fonte_neige_mod,  ONLY : fonte_neige_final
    1211      USE pbl_surface_mod,  ONLY : pbl_surface_final
    13       USE surface_data,     ONLY : ocean, ok_veget
    1412      USE phys_state_var_mod
    1513
    … …  
    6058      REAL wake_fip_glo(klon_glo)
    6159
    62 cIM "slab" ocean
    6360      REAL tsoil_p(klon,nsoilmx,nbsrf)
    64       REAL tslab_p(klon), seaice_p(klon)
    6561      REAL qsurf_p(klon,nbsrf)
    6662      REAL qsol_p(klon)
    … …  
    7369     
    7470      REAL tsoil(klon_glo,nsoilmx,nbsrf)
    75       REAL tslab(klon_glo), seaice(klon_glo)
    7671      REAL qsurf(klon_glo,nbsrf)
    7772      REAL qsol(klon_glo)
    … …  
    10398c Get a variable calculated in module fonte_neige_mod
    10499      CALL fonte_neige_final(run_off_lic_0_p)
    105 
    106 c If slab ocean then get 2 varaibles from module ocean_slab_mod
    107       IF ( ocean == 'slab' ) THEN
    108          CALL ocean_slab_final(tslab_p, seaice_p)
    109       ELSE
    110          tslab_p(:)  = 0.0
    111          seaice_p(:) = 0.0
    112       ENDIF     
    113100
    114101c======================================================================
    … …  
    150137
    151138      call Gather( tsoil_p,tsoil)
    152       call Gather( tslab_p,tslab)
    153       call Gather( seaice_p,seaice)
    154139      call Gather( qsurf_p,qsurf)
    155140      call Gather( qsol_p,qsol)
    … …  
    390375      ENDDO
    391376c
    392 cIM "slab" ocean
    393       ierr = NF_REDEF (nid)
    394 #ifdef NC_DOUBLE
    395       ierr = NF_DEF_VAR (nid, "TSLAB", NF_DOUBLE, 1, idim2,nvarid)
    396 #else
    397       ierr = NF_DEF_VAR (nid, "TSLAB", NF_FLOAT, 1, idim2,nvarid)
    398 #endif
    399       ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 33,
    400      .                        "Ecart de la SST (pour slab-ocean)")
    401       ierr = NF_ENDDEF(nid)
    402 #ifdef NC_DOUBLE
    403       ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tslab)
    404 #else
    405       ierr = NF_PUT_VAR_REAL (nid,nvarid,tslab)
    406 #endif
    407 c
    408       ierr = NF_REDEF (nid)
    409 #ifdef NC_DOUBLE
    410       ierr = NF_DEF_VAR (nid, "SEAICE", NF_DOUBLE, 1, idim2,nvarid)
    411 #else
    412       ierr = NF_DEF_VAR (nid, "SEAICE", NF_FLOAT, 1, idim2,nvarid)
    413 #endif
    414       ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 33,
    415      .                        "Glace de mer kg/m2 (pour slab-ocean)")
    416       ierr = NF_ENDDEF(nid)
    417 #ifdef NC_DOUBLE
    418       ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,seaice)
    419 #else
    420       ierr = NF_PUT_VAR_REAL (nid,nvarid,seaice)
    421 #endif
    422 c
    423377      DO nsrf = 1, nbsrf
    424378        IF (nsrf.LE.99) THEN

  • LMDZ4/trunk/libf/phylmd/phys_output_mod.F90

    r976 r996  
    188188  integer, dimension(nfiles) , save :: flag_philevsSTD   = (/ 1, 1, 1, 10 /)
    189189
    190   integer, dimension(nfiles) , save :: flag_fluxo        = (/ 1, 1, 10, 10 /)
    191   integer, dimension(nfiles) , save :: flag_fluxg        = (/ 1, 1, 10, 10 /)
    192190  integer, dimension(nfiles) , save :: flag_t_oce_sic    = (/ 1, 10, 10, 10 /)
    193 
    194   integer, dimension(nfiles) , save :: flag_lmt_bils     = (/ 1, 1, 10, 10 /)
    195   integer, dimension(nfiles) , save :: flag_tslab        = (/ 1, 1, 10, 10 /)
    196   integer, dimension(nfiles) , save :: flag_seaice       = (/ 1, 1, 10, 10 /)
    197   integer, dimension(nfiles) , save :: flag_siceh        = (/ 1, 1, 10, 10 /)
    198191
    199192  integer, dimension(nfiles) , save :: flag_weakinv      = (/ 10, 1, 10, 10 /)
    … …  
    590583       ENDDO
    591584
    592 !IM diagnostiques flux ocean-atm ou ocean-glace de mer
    593 !IM pour utilisation dans un modele de "slab" ocean
    594  CALL histdef2d(iff,flag_fluxo,"fluxo","Flux turbulents ocean-atmosphere","W/m2")
    595  CALL histdef2d(iff,flag_fluxg,"fluxg","Flux turbulents ocean-glace de mer","W/m2")
    596585 CALL histdef2d(iff,flag_t_oce_sic,"t_oce_sic","Temp mixte oce-sic","K")
    597586
    598      IF (OCEAN.EQ.'force ') THEN
    599  CALL histdef2d(iff,flag_lmt_bils,"lmt_bils","Bilan au sol atmosphere forcee","W/m2")
    600      ELSE IF (OCEAN.EQ.'slab  ') THEN
    601  CALL histdef2d(iff,flag_slab_bils,"slab_bils","Bilan au sol Slab","W/m2")
    602  CALL histdef2d(iff,flag_tslab,"tslab", "Slab SST ", "K")
    603  CALL histdef2d(iff,flag_seaice,"seaice","Slab seaice","kg/m2")
    604  CALL histdef2d(iff,flag_siceh,"siceh","Slab seaice height","m")
    605      ENDIF
     587 IF (ocean=='slab') &
     588      CALL histdef2d(iff,flag_slab_bils, "slab_wbils_oce","Bilan au sol sur ocean slab", "W/m2")
    606589
    607590 CALL histdef2d(iff,flag_ale_bl,"ale_bl","ALE BL","m2/s2")

  • LMDZ4/trunk/libf/phylmd/phys_output_write.h

    r973 r996  
    4141       IF (flag_contfracOR(iff)<=lev_files(iff)) THEN
    4242      CALL histwrite_phy(nid_files(iff),"contfracOR",itau_w,
    43      $                   pctsrf_new(:,is_ter))
     43     $                   pctsrf(:,is_ter))
    4444       ENDIF
    4545
    … …  
    599599       ENDIF !(bb2.EQ."850".OR.bb2.EQ."700".OR.
    600600       ENDDO
    601 
    602 !IM diagnostiques flux ocean-atm ou ocean-glace de mer
    603 !IM pour utilisation dans un modele de "slab" ocean
    604 
    605       IF (flag_fluxo(iff)<=lev_files(iff)) THEN
    606       DO i=1, klon
    607        IF (pctsrf(i,is_oce).GT.epsfra) THEN
    608         zx_tmp_fi2d(i) = fluxo(i)
    609        ELSE
    610         zx_tmp_fi2d(i) = 0.
    611        ENDIF
    612       ENDDO
    613       CALL histwrite_phy(nid_files(iff),"fluxo",itau_w,zx_tmp_fi2d)
    614       ENDIF
    615 
    616       IF (flag_fluxg(iff)<=lev_files(iff)) THEN
    617       DO i=1, klon
    618        IF (pctsrf(i,is_sic).GT.epsfra) THEN
    619         zx_tmp_fi2d(i) = fluxg(i)
    620        ELSE
    621         zx_tmp_fi2d(i) = 0.
    622        ENDIF
    623       ENDDO
    624       CALL histwrite_phy(nid_files(iff),"fluxg",itau_w,zx_tmp_fi2d)
    625       ENDIF
    626601
    627602      IF (flag_t_oce_sic(iff)<=lev_files(iff)) THEN
    … …  
    639614      ENDIF
    640615
    641       IF (OCEAN.EQ.'force ') THEN
    642       IF (flag_lmt_bils(iff)<=lev_files(iff)) THEN
    643       DO i=1, klon
    644       IF((pctsrf(i,is_oce).GT.epsfra).OR.
    645      $   (pctsrf(i,is_sic).GT.epsfra)) THEN
    646        zx_tmp_fi2d(i) = (radsol(i) + fluxo(i))*pctsrf(i,is_oce)+
    647      $                  fluxg(i)*pctsrf(i,is_sic)
    648       ELSE
    649        zx_tmp_fi2d(i) = 1.E+20
    650       ENDIF
    651       ENDDO
    652       CALL histwrite_phy(nid_files(iff),"lmt_bils",itau_w,zx_tmp_fi2d)
    653       ENDIF
     616      IF (type_ocean=='force ') THEN
    654617
    655618      IF (iflag_coupl.EQ.1) THEN
    … …  
    683646       ENDIF
    684647      ENDIF
    685 
    686       ELSE IF (OCEAN.EQ.'slab  ') THEN
    687 
    688       IF (flag_slab_bils(iff)<=lev_files(iff)) THEN
    689       DO i=1, klon
    690       IF((pctsrf(i,is_oce).GT.epsfra).OR.
    691      $   (pctsrf(i,is_sic).GT.epsfra)) THEN
    692        zx_tmp_fi2d(i) = (radsol(i) + fluxo(i))*pctsrf(i,is_oce)+
    693      $                  fluxg(i)*pctsrf(i,is_sic)
    694       ELSE
    695        zx_tmp_fi2d(i) = 1.E+20
    696       ENDIF
    697       ENDDO
    698       CALL histwrite_phy(nid_files(iff),"slab_bils",itau_w,zx_tmp_fi2d)
    699       ENDIF
    700 
    701       IF (flag_tslab(iff)<=lev_files(iff)) THEN
    702       DO i=1, klon
    703        IF(pctsrf(i,is_oce).GT.epsfra.OR.
    704      $    pctsrf(i,is_sic).GT.epsfra) THEN
    705         zx_tmp_fi2d(i)=tslab(i)
    706        ELSE
    707         zx_tmp_fi2d(i) = 1.E+20
    708        ENDIF
    709       ENDDO !i=1, klon
    710       CALL histwrite_phy(nid_files(iff),"tslab",itau_w,zx_tmp_fi2d)
    711       ENDIF
    712 
    713       IF (flag_seaice(iff)<=lev_files(iff)) THEN
    714       CALL histwrite_phy(nid_files(iff),"seaice",itau_w,seaice)
    715       ENDIF
    716 
    717       IF (flag_siceh(iff)<=lev_files(iff)) THEN
    718       CALL histwrite_phy(nid_files(iff),"siceh",itau_w, seaice/1000.)
    719       ENDIF
    720       ENDIF !OCEAN.EQ.force/slab
     648 
     649      ELSE IF (type_ocean=='slab  ') THEN
     650
     651      IF ( flag_slab_bils(iff)<=lev_files(iff))
     652     $     CALL histwrite_phy(
     653     $     nid_files(iff),"slab_wbils_oce",itau_w,slab_wfbils)
     654
     655      ENDIF !type_ocean == force/slab
    721656
    722657      IF (flag_weakinv(iff)<=lev_files(iff)) THEN

  • LMDZ4/trunk/libf/phylmd/physiq.F

    r987 r996  
    2424      USE vampir
    2525      USE pbl_surface_mod, ONLY : pbl_surface
    26       USE surface_data,     ONLY : ocean, ok_veget
     26      USE change_srf_frac_mod
     27      USE surface_data,     ONLY : type_ocean, ok_veget
    2728      USE phys_local_var_mod ! Variables internes non sauvegardees de la physique
    2829      USE phys_state_var_mod ! Variables sauvegardees de la physique
    29 
    30 
    31       USE ocean_slab_mod, ONLY   : ocean_slab_get_vars
    32       USE ocean_cpl_mod, ONLY    : ocean_cpl_get_vars
    33       USE ocean_forced_mod, ONLY : ocean_forced_get_vars
    3430      USE fonte_neige_mod, ONLY  : fonte_neige_get_vars
    3531      USE phys_output_mod
    … …  
    120116      LOGICAL, SAVE :: rnpb=.TRUE.
    121117c$OMP THREADPRIVATE(rnpb)
    122 cIM "slab" ocean
    123       REAL tslab(klon)    !Temperature du slab-ocean
    124       REAL seaice(klon)   !glace de mer (kg/m2)
    125       REAL fluxo(klon)    !flux turbulents ocean-glace de mer
    126       REAL fluxg(klon)    !flux turbulents ocean-atmosphere
    127118      REAL amn, amx
    128119      INTEGER igout
    … …  
    699690      REAL wfbils(klon,nbsrf) ! bilan de chaleur au sol, pour chaque
    700691C                             ! type de sous-surface et pondere par la fraction
     692      REAL slab_wfbils(klon)  ! bilan de chaleur au sol pour le cas de slab, sur les points d'ocean
     693
    701694      REAL fder(klon)         
    702695      REAL ve(klon) ! integr. verticale du transport meri. de l'energie
    … …  
    715708      SAVE lmt_pas                ! frequence de mise a jour
    716709c$OMP THREADPRIVATE(lmt_pas)
    717 cIM
    718       REAL pctsrf_new(klon,nbsrf) !pourcentage surfaces issus d'ORCHIDEE
    719 
    720 cym      SAVE pctsrf                 ! sous-fraction du sol
    721710
    722711cIM sorties
    … …  
    12391228c appel a la lecture du run.def physique
    12401229c
    1241          call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel,
     1230         call conf_phys(ok_journe, ok_mensuel,
    12421231     .                  ok_instan, ok_hf,
    12431232     .                  solarlong0,seuil_inversion,
    … …  
    12791268!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    12801269
    1281          CALL phyetat0 ("startphy.nc",ocean, ok_veget,clesphy0,tabcntr0)
     1270         CALL phyetat0 ("startphy.nc",clesphy0,tabcntr0)
    12821271cIM begin
    12831272          print*,'physiq: clwcon rnebcon ratqs',clwcon(1,1),rnebcon(1,1)
    … …  
    12881277
    12891278!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    1290 
    1291 
    1292        DO i=1,klon
    1293          IF ( abs( pctsrf(i, is_ter) + pctsrf(i, is_lic) +
    1294      $       pctsrf(i, is_oce) + pctsrf(i, is_sic)  - 1.) .GT. EPSFRA)
    1295      $       THEN
    1296             WRITE(*,*)
    1297      $   'physiq apres lecture de restart: pb sous surface au point ',
    1298      $   i, pctsrf(i, 1 : nbsrf)
    1299          ENDIF
    1300       ENDDO
    1301 
    13021279c
    13031280C on remet le calendrier a zero
    … …  
    15181495       call phys_output_open(jjmp1,nqmax,nlevSTD,clevSTD,nbteta,
    15191496     &                        ctetaSTD,dtime,presnivs,ok_veget,
    1520      &                        ocean,iflag_pbl,ok_mensuel,ok_journe,
     1497     &                        type_ocean,iflag_pbl,ok_mensuel,ok_journe,
    15211498     &                        ok_hf,ok_instan,nid_files)
    15221499c$OMP END MASTER
    … …  
    15811558c
    15821559      ENDIF
    1583 c
    1584 c   ****************     Fin  de   IF ( debut  )   ***************
    1585 c
     1560!
     1561!   ****************     Fin  de   IF ( debut  )   ***************
     1562!
     1563!
     1564! Incrementer le compteur de la physique
     1565!
     1566      itap   = itap + 1
     1567      julien = MOD(NINT(xjour),360)
     1568      if (julien .eq. 0) julien = 360
     1569
     1570!
     1571! Update fraction of the sub-surfaces (pctsrf) and
     1572! initialize, where a new fraction has appeared, all variables depending
     1573! on the surface fraction.
     1574!
     1575      CALL change_srf_frac(itap, dtime, julien,
     1576     *     pctsrf, falb1, falb2, ftsol, u10m, v10m, pbl_tke)
     1577
    15861578
    15871579! Tendances bidons pour les processus qui n'affectent pas certaines
    … …  
    17251717cIM END
    17261718c
    1727 c Incrementer le compteur de la physique
    1728 c
    1729       itap   = itap + 1
    1730       julien = MOD(NINT(xjour),360)
    1731       if (julien .eq. 0) julien = 360
    1732 c
    17331719c Mettre en action les conditions aux limites (albedo, sst, etc.).
    17341720c Prescrire l'ozone et calculer l'albedo sur l'ocean.
    … …  
    18361822     s     zxtsol,    zxfluxlat, zt2m,    qsat2m,
    18371823     s     d_t_vdf,   d_q_vdf,   d_u_vdf, d_v_vdf,
    1838      s     ycoefh,    pctsrf_new,               
     1824     s     ycoefh,    slab_wfbils,               
    18391825     d     qsol,      zq2m,      s_pblh,  s_lcl,
    18401826     d     s_capCL,   s_oliqCL,  s_cteiCL,s_pblT,
    … …  
    18471833     -     dsens,     devap,     zxsnow,
    18481834     -     zxfluxt,   zxfluxq,   q2m,     fluxq, pbl_tke )
    1849 c
    1850 c
    1851       pctsrf(:,:) = pctsrf_new(:,:)
     1835
    18521836     
    18531837!-----------------------------------------------------------------------------------------
    … …  
    33323316     .     zxfqcalving, zxfqfonte, zxffonte)
    33333317
    3334       IF (ocean == 'slab') THEN
    3335          ! Get some variables from module ocean_slab_mod
    3336          CALL ocean_slab_get_vars(tslab, seaice, fluxo, fluxg)
    3337       ELSEIF (ocean == 'couple') THEN
    3338          ! Get some variables from module ocean_cpl_mod
    3339          CALL ocean_cpl_get_vars(fluxo, fluxg)
    3340       ELSE
    3341          ! Get some variables from module ocean_forced_mod
    3342          CALL ocean_forced_get_vars(fluxo, fluxg)
    3343       ENDIF
    3344 
    33453318
    33463319c Commente par abderrahmane le 11 2 08

  • LMDZ4/trunk/libf/phylmd/surf_land_bucket_mod.F90

    r888 r996  
    88! This module is used when no external land model is choosen.
    99!
    10   USE fonte_neige_mod
    11   USE calcul_fluxs_mod
    12   USE dimphy
    13   USE mod_grid_phy_lmdz
    14   USE mod_phys_lmdz_para
    15  
    1610  IMPLICIT NONE
    1711
    … …  
    2620       fluxsens, fluxlat, tsurf_new, dflux_s, dflux_l)
    2721
     22    USE limit_read_mod
     23    USE surface_data
     24    USE fonte_neige_mod
     25    USE calcul_fluxs_mod
     26    USE cpl_mod
     27    USE dimphy
     28    USE mod_grid_phy_lmdz
     29    USE mod_phys_lmdz_para
    2830!****************************************************************************************
    2931! Bucket calculations for surface.
    … …  
    7476    REAL, DIMENSION(klon) :: zfra
    7577    REAL, DIMENSION(klon) :: radsol       ! total net radiance at surface
     78    REAL, DIMENSION(klon) :: dummy_riverflow,dummy_coastalflow
    7679    INTEGER               :: i
    7780!
    … …  
    8285!* Read from limit.nc : albedo(alb_lim), length of rugosity(z0_new)
    8386!
    84     CALL interfsur_lim(itime, dtime, jour, &
    85          knon, knindex, debut,  &
    86          alb_lim, z0_new)
    87    
     87    CALL limit_read_rug_alb(itime, dtime, jour,&
     88         knon, knindex, &
     89         z0_new, alb_lim)
    8890!
    8991!* Calcultaion of fluxes
    … …  
    145147!* Calculate the rugosity
    146148!
    147     z0_new = SQRT(z0_new**2+rugoro**2)
    148        
     149    DO i = 1, knon
     150       z0_new(i) = SQRT(z0_new(i)**2+rugoro(i)**2)
     151    END DO
     152
     153!* Send to coupler
     154!  The run-off from river and coast are not calculated in the bucket modele.
     155!  For testing purpose of the coupled modele we put the run-off to zero.
     156    IF (type_ocean=='couple') THEN
     157       dummy_riverflow(:)   = 0.0
     158       dummy_coastalflow(:) = 0.0
     159       CALL cpl_send_land_fields(itime, knon, knindex, &
     160            dummy_riverflow, dummy_coastalflow)
     161    ENDIF
     162
    149163!
    150164!* End
    … …  
    154168!****************************************************************************************
    155169!
    156   SUBROUTINE interfsur_lim(itime, dtime, jour, &
    157        knon, knindex, debut, &
    158        lmt_alb_p, lmt_rug_p)
    159    
    160 ! Cette routine sert d'interface entre le modele atmospherique et un fichier
    161 ! de conditions aux limites
    162 !
    163 ! L. Fairhead 02/2000
    164 !
    165 ! input:
    166 !   itime        numero du pas de temps courant
    167 !   dtime        pas de temps de la physique (en s)
    168 !   jour         jour a lire dans l'annee
    169 !   knon         nombre de points dans le domaine a traiter
    170 !   knindex      index des points de la surface a traiter
    171 !   debut        logical: 1er appel a la physique (initialisation)
    172 !
    173 ! output:
    174 !   lmt_alb_p      Albedo lu
    175 !   lmt_rug_p      longueur de rugosite lue
    176 
    177     INCLUDE "netcdf.inc"
    178 
    179 ! Input variables
    180 !****************************************************************************************
    181     INTEGER, INTENT(IN)                      :: itime
    182     REAL   , INTENT(IN)                      :: dtime
    183     INTEGER, INTENT(IN)                      :: jour
    184     INTEGER, INTENT(IN)                      :: knon
    185     INTEGER, DIMENSION(klon_loc), INTENT(IN) :: knindex
    186     LOGICAL, INTENT(IN)                      :: debut
    187 
    188 ! Output variables
    189 !****************************************************************************************
    190     REAL, INTENT(out), DIMENSION(klon_loc)   :: lmt_alb_p
    191     REAL, INTENT(out), DIMENSION(klon_loc)   :: lmt_rug_p
    192 
    193 ! Local variables with attribute SAVE
    194 !****************************************************************************************
    195     INTEGER,SAVE   :: lmt_pas     ! frequence de lecture des conditions limites
    196                                   ! (en pas de physique)
    197     !$OMP THREADPRIVATE(lmt_pas)
    198     LOGICAL,SAVE   :: deja_lu_sur ! pour indiquer que le jour a lire a deja
    199                                   ! lu pour une surface precedente
    200     !$OMP THREADPRIVATE(deja_lu_sur)
    201     INTEGER,SAVE                           :: jour_lu_sur
    202     !$OMP THREADPRIVATE(jour_lu_sur)
    203     CHARACTER (len = 20),SAVE              :: fich ='limit.nc'
    204     !$OMP THREADPRIVATE(fich)
    205     LOGICAL,SAVE                           :: check = .FALSE.
    206     !$OMP THREADPRIVATE(check)
    207 ! Champs lus dans le fichier de CL
    208     REAL, ALLOCATABLE , SAVE, DIMENSION(:) :: alb_lu_p, rug_lu_p
    209     !$OMP THREADPRIVATE(alb_lu_p, rug_lu_p)
    210 
    211 ! quelques variables pour netcdf
    212     INTEGER ,SAVE                          :: nid, nvarid
    213     !$OMP THREADPRIVATE(nid, nvarid)
    214     INTEGER, DIMENSION(2),SAVE             :: start, epais
    215     !$OMP THREADPRIVATE(start, epais)
    216 
    217 ! Other local variables
    218 !****************************************************************************************
    219     INTEGER                                :: ii, ierr
    220     CHARACTER (len = 20)                   :: modname = 'interfsur_lim'
    221     CHARACTER (len = 80)                   :: abort_message
    222     REAL, DIMENSION(klon_glo)              :: alb_lu
    223     REAL, DIMENSION(klon_glo)              :: rug_lu
    224 
    225 !
    226 ! End delaration
    227 !****************************************************************************************
    228 
    229     IF (debut) THEN
    230        lmt_pas = NINT(86400./dtime * 1.0) ! pour une lecture une fois par jour
    231        jour_lu_sur = jour - 1
    232        ALLOCATE(alb_lu_p(klon_loc))
    233        ALLOCATE(rug_lu_p(klon_loc))
    234     ENDIF
    235    
    236     IF ((jour - jour_lu_sur) /= 0) deja_lu_sur = .FALSE.
    237  
    238     IF (check) WRITE(*,*) modname,':: jour_lu_sur, deja_lu_sur', jour_lu_sur, deja_lu_sur
    239     IF (check) WRITE(*,*) modname,':: itime, lmt_pas', itime, lmt_pas
    240     IF (check) CALL flush(6)
    241 
    242 !   
    243 ! Tester d'abord si c'est le moment de lire le fichier
    244 !
    245     IF (MOD(itime-1, lmt_pas) == 0 .AND. .NOT. deja_lu_sur) THEN
    246 
    247 !
    248 ! Ouverture et lecture du fichier
    249 !
    250 !$OMP MASTER
    251        IF (is_mpi_root) THEN
    252           fich = TRIM(fich)
    253           IF (check) WRITE(*,*)modname,' ouverture fichier ',fich
    254           IF (check) CALL flush(6)
    255           ierr = NF_OPEN (fich, NF_NOWRITE,nid)
    256           IF (ierr.NE.NF_NOERR) THEN
    257              abort_message = 'Pb d''ouverture du fichier de conditions aux limites'
    258              CALL abort_gcm(modname,abort_message,1)
    259           ENDIF
    260 !
    261 ! La tranche de donnees a lire:
    262           start(1) = 1
    263           start(2) = jour
    264           epais(1) = klon_glo
    265           epais(2) = 1
    266 !
    267 ! Lecture albedo
    268           ierr = NF_INQ_VARID(nid, 'ALB', nvarid)
    269           IF (ierr /= NF_NOERR) THEN
    270              abort_message = 'Le champ <ALB> est absent'
    271              CALL abort_gcm(modname,abort_message,1)
    272           ENDIF
    273 #ifdef NC_DOUBLE
    274           ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais, alb_lu)
    275 #else
    276           ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais, alb_lu)
    277 #endif
    278           IF (ierr /= NF_NOERR) THEN
    279              abort_message = 'Lecture echouee pour <ALB>'
    280              CALL abort_gcm(modname,abort_message,1)
    281           ENDIF
    282 !
    283 ! Lecture rugosite!
    284           ierr = NF_INQ_VARID(nid, 'RUG', nvarid)
    285           IF (ierr /= NF_NOERR) THEN
    286              abort_message = 'Le champ <RUG> est absent'
    287              CALL abort_gcm(modname,abort_message,1)
    288           ENDIF
    289 #ifdef NC_DOUBLE
    290           ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais, rug_lu)
    291 #else
    292           ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais, rug_lu)
    293 #endif
    294           IF (ierr /= NF_NOERR) THEN
    295              abort_message = 'Lecture echouee pour <RUG>'
    296              CALL abort_gcm(modname,abort_message,1)
    297           ENDIF
    298 
    299 !
    300 ! Fin de lecture
    301           ierr = NF_CLOSE(nid)
    302 
    303        ENDIF ! is_mpi_root
    304 !$OMP END MASTER
    305 
    306        CALL Scatter(alb_lu,alb_lu_p)
    307        CALL Scatter(rug_lu,rug_lu_p)
    308 
    309        deja_lu_sur = .TRUE.
    310        jour_lu_sur = jour       
    311 
    312     ENDIF
    313  
    314 !
    315 ! Recopie des variables dans les champs de sortie
    316 !
    317     lmt_alb_p(:) = 999999.
    318     lmt_rug_p(:) = 999999.
    319     DO ii = 1, knon
    320        lmt_alb_p(ii) = alb_lu_p(knindex(ii))
    321        lmt_rug_p(ii) = rug_lu_p(knindex(ii))
    322     ENDDO
    323    
    324 
    325   END SUBROUTINE interfsur_lim
    326 !
    327 !****************************************************************************************
    328 !
    329170END MODULE surf_land_bucket_mod

  • LMDZ4/trunk/libf/phylmd/surf_land_mod.F90

    r888 r996  
    2626       snow, qsol, agesno, tsoil, &
    2727       z0_new, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    28        qsurf, tsurf_new, dflux_s, dflux_l, pctsrf_ter, &
     28       qsurf, tsurf_new, dflux_s, dflux_l, &
    2929       lwdown_m)
    3030
    … …  
    7575    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    7676    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    77     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_ter
    7877
    7978! Local variables
    … …  
    154153    ENDIF ! ok_veget
    155154
    156 !****************************************************************************************
    157 ! Return the pourcentage of land in each grid cell, even if not changed in here!
    158 !
    159 !****************************************************************************************
    160     pctsrf_ter(:) = pctsrf(:,is_ter)
    161 
    162 
    163155  END SUBROUTINE surf_land
    164156!

  • LMDZ4/trunk/libf/phylmd/surf_land_orchidee_mod.F90

    r987 r996  
    1515  USE intersurf     ! module d'ORCHIDEE
    1616  USE cpl_mod,      ONLY : cpl_send_land_fields
    17   USE surface_data, ONLY : ocean, ok_veget
     17  USE surface_data, ONLY : type_ocean
    1818  USE comgeomphy,   ONLY : cuphy, cvphy
    1919  USE mod_grid_phy_lmdz
    … …  
    199199
    200200    IF (check) WRITE(lunout,*)'Entree ', modname
    201     IF (check) WRITE(lunout,*)'ok_veget = ',ok_veget
    202201 
    203202! Initialisation
    … …  
    416415!* Send to coupler
    417416!
    418     IF (ocean=='couple') THEN
     417    IF (type_ocean=='couple') THEN
    419418       CALL cpl_send_land_fields(itime, knon, knindex, &
    420419            riverflow, coastalflow)

  • LMDZ4/trunk/libf/phylmd/surf_landice_mod.F90

    r888 r996  
    55 
    66  USE dimphy
    7   USE surface_data,     ONLY : ocean, calice, calsno
     7  USE surface_data,     ONLY : type_ocean, calice, calsno
    88  USE fonte_neige_mod,  ONLY : fonte_neige, run_off_lic
    99  USE cpl_mod,          ONLY : cpl_send_landice_fields
    … …  
    2323       snow, qsurf, qsol, agesno, &
    2424       tsoil, z0_new, alb1, alb2, evap, fluxsens, fluxlat, &
    25        tsurf_new, dflux_s, dflux_l, pctsrf_lic)
     25       tsurf_new, dflux_s, dflux_l)
    2626
    2727    INCLUDE "indicesol.h"
    … …  
    6464    REAL, DIMENSION(klon), INTENT(OUT)            :: tsurf_new
    6565    REAL, DIMENSION(klon), INTENT(OUT)            :: dflux_s, dflux_l     
    66     REAL, DIMENSION(klon), INTENT(OUT)            :: pctsrf_lic
    6766
    6867! Local variables
    … …  
    153152    z0_new(:) = rugoro(:)
    154153
    155 
    156 !****************************************************************************************
    157 ! Return the pourcentage for this sub-surface
    158 !
    159 !****************************************************************************************
    160     pctsrf_lic(:) = pctsrf(:,is_lic)
    161    
    162 
    163154!****************************************************************************************
    164155! Send run-off on land-ice to coupler if coupled ocean.
    … …  
    166157!
    167158!****************************************************************************************
    168     IF (ocean=='couple') THEN
     159    IF (type_ocean=='couple') THEN
    169160       CALL cpl_send_landice_fields(itime, knon, knindex, run_off_lic)
    170161    ENDIF

  • LMDZ4/trunk/libf/phylmd/surf_ocean_mod.F90

    r888 r996  
    55
    66  USE dimphy
    7   USE surface_data, ONLY     : ocean
     7  USE surface_data, ONLY     : type_ocean
    88  USE ocean_forced_mod, ONLY : ocean_forced_noice
    99  USE ocean_slab_mod, ONLY   : ocean_slab_noice
    … …  
    1717!
    1818  SUBROUTINE surf_ocean(rlon, rlat, swnet, lwnet, alb1, &
    19        rugos, windsp, rmu0, fder, &
     19       rugos, windsp, rmu0, fder, tsurf_in, &
    2020       itime, dtime, jour, knon, knindex, &
    2121       debut, &
    … …  
    2525       snow, qsurf, agesno, &
    2626       z0_new, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    27        tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     27       tsurf_new, dflux_s, dflux_l, lmt_bils)
    2828!
    2929! This subroutine will make a call to ocean_XXX_noice according to the ocean mode (force,
    … …  
    4747    REAL, DIMENSION(klon), INTENT(IN)        :: rmu0 
    4848    REAL, DIMENSION(klon), INTENT(IN)        :: fder
     49    REAL, DIMENSION(klon), INTENT(IN)        :: tsurf_in
    4950    REAL, DIMENSION(klon), INTENT(IN)        :: p1lay
    5051    REAL, DIMENSION(klon), INTENT(IN)        :: tq_cdrag
    … …  
    7475    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    7576    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    76     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_oce
    77 
     77    REAL, DIMENSION(klon), INTENT(OUT)       :: lmt_bils
    7878
    7979! Local variables
    … …  
    9797! Switch according to type of ocean (couple, slab or forced)
    9898!****************************************************************************************
    99     SELECT CASE(ocean)
     99    SELECT CASE(type_ocean)
    100100    CASE('couple')
    101101       CALL ocean_cpl_noice( &
    … …  
    106106            p1lay, tq_cdrag, precip_rain, precip_snow,temp_air,spechum,&
    107107            petAcoef, peqAcoef, petBcoef, peqBcoef, &
    108             ps, u1_lay, v1_lay, pctsrf, &
     108            ps, u1_lay, v1_lay, &
    109109            radsol, snow, agesno, &
    110110            qsurf, evap, fluxsens, fluxlat, &
    111             tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     111            tsurf_new, dflux_s, dflux_l)
    112112
    113113    CASE('slab')
    114114       CALL ocean_slab_noice( &
    115             dtime, knon, knindex, &
     115            itime, dtime, jour, knon, knindex, &
    116116            p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
    117117            petAcoef, peqAcoef, petBcoef, peqBcoef, &
    118             ps, u1_lay, v1_lay, &
     118            ps, u1_lay, v1_lay, tsurf_in, &
    119119            radsol, snow, agesno, &
    120120            qsurf, evap, fluxsens, fluxlat, &
    121             tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     121            tsurf_new, dflux_s, dflux_l, lmt_bils)
    122122       
    123123    CASE('force')
    … …  
    131131            radsol, snow, agesno, &
    132132            qsurf, evap, fluxsens, fluxlat, &
    133             tsurf_new, dflux_s, dflux_l, pctsrf_oce)
     133            tsurf_new, dflux_s, dflux_l)
    134134    END SELECT
    135135

  • LMDZ4/trunk/libf/phylmd/surf_seaice_mod.F90

    r888 r996  
    55
    66  USE dimphy
    7   USE surface_data, ONLY     : ocean
    8   USE ocean_slab_mod, ONLY   : ocean_slab_ice
     7  USE surface_data
    98  USE ocean_forced_mod, ONLY : ocean_forced_ice
    109  USE ocean_cpl_mod, ONLY    : ocean_cpl_ice
    … …  
    2423       snow, qsurf, qsol, agesno, tsoil, &
    2524       z0_new, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    26        tsurf_new, dflux_s, dflux_l, pctsrf_sic)
     25       tsurf_new, dflux_s, dflux_l)
    2726!
    2827! This subroutine will make a call to ocean_XXX_ice according to the ocean mode (force,
    … …  
    7069    REAL, DIMENSION(klon), INTENT(OUT)       :: tsurf_new
    7170    REAL, DIMENSION(klon), INTENT(OUT)       :: dflux_s, dflux_l     
    72     REAL, DIMENSION(klon), INTENT(OUT)       :: pctsrf_sic
    7371
    7472! Local arguments
    … …  
    9189!
    9290!****************************************************************************************
    93     SELECT CASE(ocean)
    94     CASE('couple')
     91    IF (type_ocean == 'couple') THEN
     92       
    9593       CALL ocean_cpl_ice( &
    96           rlon, rlat, swnet, lwnet, alb1, &
    97           fder, &
    98           itime, dtime, knon, knindex, &
    99           lafin,&
    100           p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
    101           petAcoef, peqAcoef, petBcoef, peqBcoef, &
    102           ps, u1_lay, v1_lay, pctsrf, &
    103           radsol, snow, qsurf, &
    104           alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    105           tsurf_new, dflux_s, dflux_l, pctsrf_sic)
     94            rlon, rlat, swnet, lwnet, alb1, &
     95            fder, &
     96            itime, dtime, knon, knindex, &
     97            lafin,&
     98            p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
     99            petAcoef, peqAcoef, petBcoef, peqBcoef, &
     100            ps, u1_lay, v1_lay, pctsrf, &
     101            radsol, snow, qsurf, &
     102            alb1_new, alb2_new, evap, fluxsens, fluxlat, &
     103            tsurf_new, dflux_s, dflux_l)
     104       
     105    ELSE IF (type_ocean == 'force' .OR. (type_ocean == 'slab' .AND. version_ocean=='sicOBS')) THEN
     106       CALL ocean_forced_ice(itime, dtime, jour, knon, knindex, &
     107            debut, &
     108            tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
     109            petAcoef, peqAcoef, petBcoef, peqBcoef, &
     110            ps, u1_lay, v1_lay, &
     111            radsol, snow, qsol, agesno, tsoil, &
     112            qsurf, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
     113            tsurf_new, dflux_s, dflux_l)
    106114
    107     CASE('slab')
    108        CALL ocean_slab_ice( &
    109           itime, dtime, jour, knon, knindex, &
    110           debut, &
    111           tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
    112           petAcoef, peqAcoef, petBcoef, peqBcoef, &
    113           ps, u1_lay, v1_lay, &
    114           radsol, snow, qsurf, qsol, agesno, tsoil, &
    115           alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    116           tsurf_new, dflux_s, dflux_l, pctsrf_sic)
    117    
    118     CASE('force')
    119        CALL ocean_forced_ice(itime, dtime, jour, knon, knindex, &
    120           debut, &
    121           tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
    122           petAcoef, peqAcoef, petBcoef, peqBcoef, &
    123           ps, u1_lay, v1_lay, &
    124           radsol, snow, qsol, agesno, tsoil, &
    125           qsurf, alb1_new, alb2_new, evap, fluxsens, fluxlat, &
    126           tsurf_new, dflux_s, dflux_l, pctsrf_sic)
    127     END SELECT
     115    ELSE IF (type_ocean == 'slab') THEN
     116!!$       CALL ocean_slab_ice( &
     117!!$          itime, dtime, jour, knon, knindex, &
     118!!$          debut, &
     119!!$          tsurf, p1lay, tq_cdrag, precip_rain, precip_snow, temp_air, spechum,&
     120!!$          petAcoef, peqAcoef, petBcoef, peqBcoef, &
     121!!$          ps, u1_lay, v1_lay, pctsrf, &
     122!!$          radsol, snow, qsurf, qsol, agesno, tsoil, &
     123!!$          alb1_new, alb2_new, evap, fluxsens, fluxlat, &
     124!!$          tsurf_new, dflux_s, dflux_l)
     125
     126    END IF
    128127
    129128!****************************************************************************************

  • LMDZ4/trunk/libf/phylmd/surface_data.F90

    r803 r996  
    44MODULE surface_data
    55
    6   USE dimphy, ONLY : klon
    7 
    86  REAL, PARAMETER        :: calice=1.0/(5.1444e+06*0.15)
    97  REAL, PARAMETER        :: tau_gl=86400.*5.
    108  REAL, PARAMETER        :: calsno=1./(2.3867e+06*.15)
    119 
    12   LOGICAL, SAVE          :: ok_veget       ! true for use of vegetation model ORCHIDEE
     10  LOGICAL, SAVE          :: ok_veget      ! true for use of vegetation model ORCHIDEE
    1311  !$OMP THREADPRIVATE(ok_veget)
    14   CHARACTER(len=6), SAVE :: ocean          ! force/slab/couple
    15   !$OMP THREADPRIVATE(ocean)
     12
     13  CHARACTER(len=6), SAVE :: type_ocean    ! force/slab/couple
     14  !$OMP THREADPRIVATE(type_ocean)
     15
     16  ! if type_ocean=couple : version_ocean=opa8 ou nemo
     17  ! if type_ocean=slab   : version_ocean=sicOBS
     18  CHARACTER(len=6), SAVE :: version_ocean
     19  !$OMP THREADPRIVATE(version_ocean)
    1620
    1721END MODULE surface_data

  • LMDZ4/trunk/libf/phylmd/write_histday.h

    r895 r996  
    3232      CALL histwrite_phy(nid_day,"contfracATM",itau_w,zx_tmp_fi2d)
    3333c
    34 cym      CALL gr_fi_ecrit(1,klon,iim,jjmp1,pctsrf_new(:,is_ter),zx_tmp_2d)
     34cym      CALL gr_fi_ecrit(1,klon,iim,jjmp1,pctsrf(:,is_ter),zx_tmp_2d)
    3535      CALL histwrite_phy(nid_day,"contfracOR",itau_w,
    36      &                   pctsrf_new(:,is_ter))
     36     &                   pctsrf(:,is_ter))
    3737c
    3838

  • LMDZ4/trunk/libf/phylmd/write_histday_seri.h

    r778 r996  
    191191c
    192192c     ok_msk=.TRUE.
    193 c     msk(1:klon)=pctsrf_new(1:klon,is_ter)
     193c     msk(1:klon)=pctsrf(1:klon,is_ter)
    194194c     CALL moyglo_pondaire(klon, zx_tmp_fi2d, airephy,
    195195c    .                     ok_msk, msk, moyglo)

  • LMDZ4/trunk/libf/phylmd/write_histhf.h

    r888 r996  
    2424      CALL histwrite_phy(nid_hf,"contfracATM",itau_w,zx_tmp_fi2d)
    2525c
    26 cym      CALL gr_fi_ecrit(1,klon,iim,jjmp1,pctsrf_new(:,is_ter),zx_tmp_2d)
     26cym      CALL gr_fi_ecrit(1,klon,iim,jjmp1,pctsrf(:,is_ter),zx_tmp_2d)
    2727      CALL histwrite_phy(nid_hf,"contfracOR",itau_w,
    28      .                   pctsrf_new(:,is_ter))
     28     .                   pctsrf(:,is_ter))
    2929c
    3030cym      CALL gr_fi_ecrit(1, klon,iim,jjmp1, zt2m,zx_tmp_2d)

  • LMDZ4/trunk/libf/phylmd/write_paramLMDZ_phy.h

    r956 r996  
    77c Variables type caractere : plusieurs valeurs possibles
    88c
    9       IF(ocean.EQ.'force ') THEN
    10        zx_tmp_2d(1:iim,1:jjmp1)=1.
    11       ELSE IF(ocean.EQ.'slab  ') THEN
     9      IF(type_ocean.EQ.'force ') THEN
     10       zx_tmp_2d(1:iim,1:jjmp1)=1.
     11      ELSE IF(type_ocean.EQ.'slab  ') THEN
    1212       zx_tmp_2d(1:iim,1:jjmp1)=2.
    13       ELSE IF(ocean.EQ.'couple') THEN
     13      ELSE IF(type_ocean.EQ.'couple') THEN
    1414       zx_tmp_2d(1:iim,1:jjmp1)=3.
    1515      ENDIF
    1616      CALL histwrite(nid_ctesGCM,"ocean",itau_w,
    17      .               zx_tmp_2d,iim*jjmp1,ndex2d)
    18 c
    19       IF(ok_slab_sicOBS) THEN
    20        zx_tmp_2d(1:iim,1:jjmp1)=1.
    21       ELSE
    22        zx_tmp_2d(1:iim,1:jjmp1)=0.
    23       ENDIF
    24       CALL histwrite(nid_ctesGCM,"ok_slab_sicOBS",itau_w,
    2517     .               zx_tmp_2d,iim*jjmp1,ndex2d)
    2618c
Note: See TracChangeset for help on using the changeset viewer.