source: LMDZ6/branches/Ocean_skin/libf/phylmd/cpl_mod.F90 @ 5453

Last change on this file since 5453 was 4368, checked in by lguez, 2 years ago

Sync latest trunk changes to Ocean_skin

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 67.1 KB
Line 
1!
2MODULE cpl_mod
3!
4! This module excahanges and transforms all fields that should be recieved or sent to
5! coupler. The transformation of the fields are done from the grid 1D-array in phylmd
6! to the regular 2D grid accepted by the coupler. Cumulation of the fields for each
7! timestep is done in here.
8!
9! Each type of surface that recevie fields from the coupler have a subroutine named
10! cpl_receive_XXX_fields and each surface that have fields to be sent to the coupler
11! have a subroutine named cpl_send_XXX_fields.
12!
13!*************************************************************************************
14
15! Use statements
16!*************************************************************************************
17  USE dimphy, ONLY : klon
18  USE mod_phys_lmdz_para
19  USE ioipsl
20  USE iophy
21
22! The module oasis is always used. Without the cpp key CPP_COUPLE only the parameters
23! in the module are compiled and not the subroutines.
24  USE oasis
25  USE write_field_phy
26  USE time_phylmdz_mod, ONLY: day_step_phy
27 
28! Global attributes
29!*************************************************************************************
30  IMPLICIT NONE
31  PRIVATE
32
33  ! All subroutine are public except cpl_send_all
34  PUBLIC :: cpl_init, cpl_receive_frac, cpl_receive_ocean_fields, cpl_receive_seaice_fields, &
35       cpl_send_ocean_fields, cpl_send_seaice_fields, cpl_send_land_fields, &
36       cpl_send_landice_fields, gath2cpl
37 
38
39! Declaration of module variables
40!*************************************************************************************
41! variable for coupling period
42  INTEGER, SAVE :: nexca
43  !$OMP THREADPRIVATE(nexca)
44
45! variables for cumulating fields during a coupling periode :
46  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_sols, cpl_nsol, cpl_rain
47  !$OMP THREADPRIVATE(cpl_sols,cpl_nsol,cpl_rain)
48  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_snow, cpl_evap, cpl_tsol
49  !$OMP THREADPRIVATE(cpl_snow,cpl_evap,cpl_tsol)
50
51  REAL, ALLOCATABLE, SAVE:: cpl_delta_sst(:), cpl_delta_sal(:), cpl_dter(:), &
52       cpl_dser(:), cpl_dt_ds(:)
53  !$OMP THREADPRIVATE(cpl_delta_sst, cpl_delta_sal, cpl_dter, cpl_dser)
54  !$OMP THREADPRIVATE(cpl_dt_ds)
55 
56  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_fder, cpl_albe, cpl_taux, cpl_tauy
57  !$OMP THREADPRIVATE(cpl_fder,cpl_albe,cpl_taux,cpl_tauy)
58  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_windsp
59  !$OMP THREADPRIVATE(cpl_windsp)
60  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_sens_rain, cpl_sens_snow
61  !$OMP THREADPRIVATE(cpl_sens_rain, cpl_sens_snow)
62  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_taumod
63  !$OMP THREADPRIVATE(cpl_taumod)
64  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_atm_co2
65  !$OMP THREADPRIVATE(cpl_atm_co2)
66  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_rriv2D, cpl_rcoa2D, cpl_rlic2D
67  !$OMP THREADPRIVATE(cpl_rriv2D,cpl_rcoa2D,cpl_rlic2D)
68  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: rlic_in_frac2D  ! fraction for continental ice
69  !$OMP THREADPRIVATE(rlic_in_frac2D)
70
71! variables read from coupler :
72  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_sst     ! sea surface temperature
73  !$OMP THREADPRIVATE(read_sst)
74  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_sit     ! sea ice temperature
75  !$OMP THREADPRIVATE(read_sit)
76
77  REAL, ALLOCATABLE, SAVE:: read_sss(:, :)
78  ! bulk salinity of the surface layer of the ocean, in ppt
79  !$OMP THREADPRIVATE(read_sss)
80
81  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_sic     ! sea ice fraction
82  !$OMP THREADPRIVATE(read_sic)
83  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_alb_sic ! albedo at sea ice
84  !$OMP THREADPRIVATE(read_alb_sic)
85  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_u0, read_v0 ! ocean surface current
86  !$OMP THREADPRIVATE(read_u0,read_v0)
87  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: read_co2     ! ocean co2 flux
88  !$OMP THREADPRIVATE(read_co2)
89  INTEGER, ALLOCATABLE, DIMENSION(:), SAVE  :: unity
90  !$OMP THREADPRIVATE(unity)
91  INTEGER, SAVE                             :: nidct, nidcs
92  !$OMP THREADPRIVATE(nidct,nidcs)
93
94! variables to be sent to the coupler
95  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: cpl_sols2D, cpl_nsol2D, cpl_rain2D
96  !$OMP THREADPRIVATE(cpl_sols2D, cpl_nsol2D, cpl_rain2D)
97  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: cpl_snow2D, cpl_evap2D, cpl_tsol2D
98  !$OMP THREADPRIVATE(cpl_snow2D, cpl_evap2D, cpl_tsol2D)
99
100  REAL, ALLOCATABLE, SAVE:: cpl_delta_sst_2D(:, :), cpl_delta_sal_2D(:, :), &
101       cpl_dter_2D(:, :), cpl_dser_2D(:, :), cpl_dt_ds_2D(:, :)
102  !$OMP THREADPRIVATE(cpl_delta_sst_2D, cpl_delta_sal_2D)
103  !$OMP THREADPRIVATE(cpl_dter_2D, cpl_dser_2D, cpl_dt_ds_2D)
104
105  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: cpl_fder2D, cpl_albe2D
106  !$OMP THREADPRIVATE(cpl_fder2D, cpl_albe2D)
107  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: cpl_taux2D, cpl_tauy2D
108  !$OMP THREADPRIVATE(cpl_taux2D, cpl_tauy2D)
109  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE :: cpl_taumod2D
110  !$OMP THREADPRIVATE(cpl_taumod2D)
111  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_windsp2D
112  !$OMP THREADPRIVATE(cpl_windsp2D)
113  REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE   :: cpl_sens_rain2D, cpl_sens_snow2D
114  !$OMP THREADPRIVATE(cpl_sens_rain2D, cpl_sens_snow2D)
115  REAL, ALLOCATABLE, DIMENSION(:,:), SAVE   :: cpl_atm_co22D
116  !$OMP THREADPRIVATE(cpl_atm_co22D)
117
118!!!!!!!!!! variable for calving
119  INTEGER, PARAMETER :: nb_zone_calving = 3
120  REAL,ALLOCATABLE, DIMENSION(:,:,:),SAVE :: area_calving
121  !$OMP THREADPRIVATE(area_calving)
122  REAL,ALLOCATABLE, DIMENSION(:,:),SAVE :: cell_area2D
123  !$OMP THREADPRIVATE(cell_area2D)
124  INTEGER, SAVE :: ind_calving(nb_zone_calving)
125  !$OMP THREADPRIVATE(ind_calving)
126
127  LOGICAL,SAVE :: cpl_old_calving
128  !$OMP THREADPRIVATE(cpl_old_calving)
129 
130CONTAINS
131!
132!************************************************************************************
133!
134  SUBROUTINE cpl_init(dtime, rlon, rlat)
135    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, fco2_ocn_day
136    USE surface_data
137    USE indice_sol_mod
138    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat, grid1dTo2d_glo, klon_glo, grid_type, unstructured, regular_lonlat
139    USE time_phylmdz_mod, ONLY: annee_ref, day_ini, itau_phy, itaufin_phy
140    USE print_control_mod, ONLY: lunout
141    USE geometry_mod, ONLY : longitude_deg, latitude_deg, ind_cell_glo, cell_area
142    USE ioipsl_getin_p_mod, ONLY: getin_p
143    use config_ocean_skin_m, only: activate_ocean_skin
144
145! Input arguments
146!*************************************************************************************
147    REAL, INTENT(IN)                  :: dtime
148    REAL, DIMENSION(klon), INTENT(IN) :: rlon, rlat
149
150! Local variables
151!*************************************************************************************
152    INTEGER                           :: error, sum_error, ig, i
153    INTEGER                           :: jf, nhoridct
154    INTEGER                           :: nhoridcs
155    INTEGER                           :: idtime
156    INTEGER                           :: idayref
157    INTEGER                           :: npas ! only for OASIS2
158    REAL                              :: zjulian
159    REAL, DIMENSION(nbp_lon,nbp_lat)  :: zx_lon, zx_lat
160    CHARACTER(len = 20)               :: modname = 'cpl_init'
161    CHARACTER(len = 80)               :: abort_message
162    CHARACTER(len=80)                 :: clintocplnam, clfromcplnam
163    REAL, DIMENSION(klon_mpi)         :: rlon_mpi, rlat_mpi, cell_area_mpi
164    INTEGER, DIMENSION(klon_mpi)           :: ind_cell_glo_mpi
165    REAL, DIMENSION(nbp_lon,jj_nb)         :: lon2D, lat2D
166    INTEGER :: mask_calving(nbp_lon,jj_nb,nb_zone_calving)
167    REAL :: pos
168
169!***************************************
170! Use old calving or not (default new calving method)
171! New calving method should be used with DYNAMICO and when using new coupling
172! weights.
173    cpl_old_calving=.FALSE.
174    CALL getin_p("cpl_old_calving",cpl_old_calving)
175    WRITE(lunout,*)' cpl_old_calving = ', cpl_old_calving
176
177
178!*************************************************************************************
179! Calculate coupling period
180!
181!*************************************************************************************
182     
183    npas = itaufin_phy
184!    nexca = 86400 / dtime
185    nexca = t_coupl / dtime
186    WRITE(lunout,*)' ##### Ocean couple #####'
187    WRITE(lunout,*)' Valeurs des pas de temps'
188    WRITE(lunout,*)' npas = ', npas
189    WRITE(lunout,*)' nexca = ', nexca
190   
191!*************************************************************************************
192! Allocate variables
193!
194!*************************************************************************************
195    error = 0
196    sum_error = 0
197
198    ALLOCATE(unity(klon), stat = error)
199    sum_error = sum_error + error
200    ALLOCATE(cpl_sols(klon,2), stat = error)
201    sum_error = sum_error + error
202    ALLOCATE(cpl_nsol(klon,2), stat = error)
203    sum_error = sum_error + error
204    ALLOCATE(cpl_rain(klon,2), stat = error)
205    sum_error = sum_error + error
206    ALLOCATE(cpl_snow(klon,2), stat = error)
207    sum_error = sum_error + error
208    ALLOCATE(cpl_evap(klon,2), stat = error)
209    sum_error = sum_error + error
210    ALLOCATE(cpl_tsol(klon,2), stat = error)
211    sum_error = sum_error + error
212    ALLOCATE(cpl_fder(klon,2), stat = error)
213    sum_error = sum_error + error
214    ALLOCATE(cpl_albe(klon,2), stat = error)
215    sum_error = sum_error + error
216    ALLOCATE(cpl_taux(klon,2), stat = error)
217    sum_error = sum_error + error
218    ALLOCATE(cpl_tauy(klon,2), stat = error)
219    sum_error = sum_error + error
220    ALLOCATE(cpl_windsp(klon,2), stat = error)
221    sum_error = sum_error + error
222    ALLOCATE(cpl_taumod(klon,2), stat = error)
223    sum_error = sum_error + error
224    ALLOCATE(cpl_sens_rain(klon,2), stat = error)
225    sum_error = sum_error + error
226    ALLOCATE(cpl_sens_snow(klon,2), stat = error)
227    sum_error = sum_error + error
228    ALLOCATE(cpl_rriv2D(nbp_lon,jj_nb), stat=error)
229    sum_error = sum_error + error
230    ALLOCATE(cpl_rcoa2D(nbp_lon,jj_nb), stat=error)
231    sum_error = sum_error + error
232    ALLOCATE(cpl_rlic2D(nbp_lon,jj_nb), stat=error)
233    sum_error = sum_error + error
234    ALLOCATE(rlic_in_frac2D(nbp_lon,jj_nb), stat=error)
235    sum_error = sum_error + error
236    ALLOCATE(read_sst(nbp_lon, jj_nb), stat = error)
237    sum_error = sum_error + error
238    ALLOCATE(read_sic(nbp_lon, jj_nb), stat = error)
239    sum_error = sum_error + error
240    ALLOCATE(read_sit(nbp_lon, jj_nb), stat = error)
241    sum_error = sum_error + error
242
243    if (activate_ocean_skin >= 1) then
244       ALLOCATE(read_sss(nbp_lon, jj_nb), stat = error)
245       sum_error = sum_error + error
246   
247       if (activate_ocean_skin == 2) then
248          ALLOCATE(cpl_delta_sst(klon), cpl_delta_sal(klon), cpl_dter(klon), &
249               cpl_dser(klon), cpl_dt_ds(klon), stat = error)
250          sum_error = sum_error + error
251       end if
252    end if
253
254    ALLOCATE(read_alb_sic(nbp_lon, jj_nb), stat = error)
255    sum_error = sum_error + error
256    ALLOCATE(read_u0(nbp_lon, jj_nb), stat = error)
257    sum_error = sum_error + error
258    ALLOCATE(read_v0(nbp_lon, jj_nb), stat = error)
259    sum_error = sum_error + error
260
261    IF (carbon_cycle_cpl) THEN
262       ALLOCATE(read_co2(nbp_lon, jj_nb), stat = error)
263       sum_error = sum_error + error
264       ALLOCATE(cpl_atm_co2(klon,2), stat = error)
265       sum_error = sum_error + error
266
267! Allocate variable in carbon_cycle_mod
268       IF (.NOT.ALLOCATED(fco2_ocn_day)) ALLOCATE(fco2_ocn_day(klon), stat = error)
269       sum_error = sum_error + error
270    ENDIF
271
272! calving initialization
273    ALLOCATE(area_calving(nbp_lon, jj_nb, nb_zone_calving), stat = error)
274    sum_error = sum_error + error
275    ALLOCATE(cell_area2D(nbp_lon, jj_nb), stat = error)   
276    sum_error = sum_error + error
277
278    CALL gather_omp(longitude_deg,rlon_mpi)
279    CALL gather_omp(latitude_deg,rlat_mpi)
280    CALL gather_omp(ind_cell_glo,ind_cell_glo_mpi)
281    CALL gather_omp(cell_area,cell_area_mpi)
282     
283    IF (is_omp_master) THEN
284      CALL Grid1DTo2D_mpi(rlon_mpi,lon2D)
285      CALL Grid1DTo2D_mpi(rlat_mpi,lat2D)
286      CALL Grid1DTo2D_mpi(cell_area_mpi,cell_area2D)
287      !--the next line is required for lat-lon grid and should have no impact
288      !--for an unstructured grid for which nbp_lon=1
289      !--if north pole in process mpi then divide cell area of pole cell by number of replicates
290      IF (is_north_pole_dyn) cell_area2D(:,1)=cell_area2D(:,1)/FLOAT(nbp_lon)
291      !--if south pole in process mpi then divide cell area of pole cell by number of replicates
292      IF (is_south_pole_dyn) cell_area2D(:,jj_nb)=cell_area2D(:,jj_nb)/FLOAT(nbp_lon)
293      mask_calving(:,:,:) = 0
294      WHERE ( lat2D >= 40) mask_calving(:,:,1) = 1
295      WHERE ( lat2D < 40 .AND. lat2D > -50) mask_calving(:,:,2) = 1
296      WHERE ( lat2D <= -50) mask_calving(:,:,3) = 1
297   
298   
299      DO i=1,nb_zone_calving
300        area_calving(:,:,i)=mask_calving(:,:,i)*cell_area2D(:,:)
301        pos=1
302        IF (i>1) pos = 1 + ((nbp_lon*nbp_lat-1)*(i-1))/(nb_zone_calving-1)
303     
304        ind_calving(i)=0
305        IF (grid_type==unstructured) THEN
306
307          DO ig=1,klon_mpi
308            IF (ind_cell_glo_mpi(ig)==pos) ind_calving(i)=ig
309          ENDDO
310
311        ELSE IF (grid_type==regular_lonlat) THEN
312          IF ((ij_begin<=pos .AND. ij_end>=pos) .OR. (ij_begin<=pos .AND. is_south_pole_dyn )) THEN
313            ind_calving(i)=pos-(jj_begin-1)*nbp_lon
314          ENDIF
315        ENDIF
316     
317      ENDDO
318    ENDIF
319   
320    IF (sum_error /= 0) THEN
321       abort_message='Pb allocation variables couplees'
322       CALL abort_physic(modname,abort_message,1)
323    ENDIF
324!*************************************************************************************
325! Initialize the allocated varaibles
326!
327!*************************************************************************************
328    DO ig = 1, klon
329       unity(ig) = ig
330    ENDDO
331
332!*************************************************************************************
333! Initialize coupling
334!
335!*************************************************************************************
336    idtime = INT(dtime)
337#ifdef CPP_COUPLE
338    CALL inicma
339#endif
340
341!*************************************************************************************
342! initialize NetCDF output
343!
344!*************************************************************************************
345    IF (is_sequential) THEN
346       idayref = day_ini
347       CALL ymds2ju(annee_ref, 1, idayref, 0.0, zjulian)
348       CALL grid1dTo2d_glo(rlon,zx_lon)
349       DO i = 1, nbp_lon
350          zx_lon(i,1) = rlon(i+1)
351          zx_lon(i,nbp_lat) = rlon(i+1)
352       ENDDO
353       CALL grid1dTo2d_glo(rlat,zx_lat)
354       clintocplnam="cpl_atm_tauflx"
355       CALL histbeg(clintocplnam,nbp_lon,zx_lon(:,1),nbp_lat,zx_lat(1,:),&
356            1,nbp_lon,1,nbp_lat, itau_phy,zjulian,dtime,nhoridct,nidct)
357! no vertical axis
358       CALL histdef(nidct, 'tauxe','tauxe', &
359            "-",nbp_lon,nbp_lat, nhoridct, 1, 1, 1, -99, 32, "inst", dtime,dtime)
360       CALL histdef(nidct, 'tauyn','tauyn', &
361            "-",nbp_lon,nbp_lat, nhoridct, 1, 1, 1, -99, 32, "inst", dtime,dtime)
362       CALL histdef(nidct, 'tmp_lon','tmp_lon', &
363            "-",nbp_lon,nbp_lat, nhoridct, 1, 1, 1, -99, 32, "inst", dtime,dtime)
364       CALL histdef(nidct, 'tmp_lat','tmp_lat', &
365            "-",nbp_lon,nbp_lat, nhoridct, 1, 1, 1, -99, 32, "inst", dtime,dtime)
366       DO jf=1,maxsend
367         IF (infosend(i)%action) THEN
368             CALL histdef(nidct, infosend(i)%name ,infosend(i)%name , &
369                "-",nbp_lon,nbp_lat,nhoridct,1,1,1,-99,32,"inst",dtime,dtime)
370         ENDIF
371       ENDDO
372       CALL histend(nidct)
373       CALL histsync(nidct)
374       
375       clfromcplnam="cpl_atm_sst"
376       CALL histbeg(clfromcplnam,nbp_lon,zx_lon(:,1),nbp_lat,zx_lat(1,:),1,nbp_lon,1,nbp_lat, &
377            0,zjulian,dtime,nhoridcs,nidcs)
378! no vertical axis
379       DO jf=1,maxrecv
380         IF (inforecv(i)%action) THEN
381             CALL histdef(nidcs,inforecv(i)%name ,inforecv(i)%name , &
382                "-",nbp_lon,nbp_lat,nhoridcs,1,1,1,-99,32,"inst",dtime,dtime)
383         ENDIF
384       ENDDO
385       CALL histend(nidcs)
386       CALL histsync(nidcs)
387
388    ENDIF    ! is_sequential
389   
390!*************************************************************************************
391! compatibility test
392!
393!*************************************************************************************
394    IF (carbon_cycle_cpl .AND. version_ocean=='opa8') THEN
395       abort_message='carbon_cycle_cpl does not work with opa8'
396       CALL abort_physic(modname,abort_message,1)
397    ENDIF
398
399  END SUBROUTINE cpl_init
400 
401!
402!*************************************************************************************
403!
404 
405  SUBROUTINE cpl_receive_frac(itime, dtime, pctsrf, is_modified)
406! This subroutine receives from coupler for both ocean and seaice
407! 4 fields : read_sst, read_sic, read_sit and read_alb_sic.
408! The new sea-ice-land-landice fraction is returned. The others fields
409! are stored in this module.
410    USE surface_data
411    USE geometry_mod, ONLY : longitude_deg, latitude_deg
412    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl
413    USE indice_sol_mod
414    USE time_phylmdz_mod, ONLY: start_time, itau_phy
415    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
416    use config_ocean_skin_m, only: activate_ocean_skin
417
418    INCLUDE "YOMCST.h"
419
420! Arguments
421!************************************************************************************
422    INTEGER, INTENT(IN)                        :: itime
423    REAL, INTENT(IN)                           :: dtime
424    REAL, DIMENSION(klon,nbsrf), INTENT(INOUT) :: pctsrf
425    LOGICAL, INTENT(OUT)                       :: is_modified
426
427! Local variables
428!************************************************************************************
429    INTEGER                                 :: j, i, time_sec
430    INTEGER                                 :: itau_w
431    INTEGER, DIMENSION(nbp_lon*nbp_lat)     :: ndexcs
432    CHARACTER(len = 20)                     :: modname = 'cpl_receive_frac'
433    CHARACTER(len = 80)                     :: abort_message
434    REAL, DIMENSION(klon)                   :: read_sic1D
435    REAL, DIMENSION(nbp_lon,jj_nb,maxrecv)      :: tab_read_flds
436    REAL, DIMENSION(klon,nbsrf)             :: pctsrf_old
437    REAL, DIMENSION(klon_mpi)               :: rlon_mpi, rlat_mpi
438    REAL, DIMENSION(nbp_lon,jj_nb)             :: tmp_lon, tmp_lat
439    REAL, DIMENSION(nbp_lon,jj_nb)             :: tmp_r0
440
441!*************************************************************************************
442! Start calculation
443! Get fields from coupler
444!
445!*************************************************************************************
446
447    is_modified=.FALSE.
448
449! Check if right moment to receive from coupler
450    IF (MOD(itime, nexca) == 1) THEN
451       is_modified=.TRUE.
452 
453       time_sec=(itime-1)*dtime
454#ifdef CPP_COUPLE
455!$OMP MASTER
456    CALL fromcpl(time_sec, tab_read_flds)
457!$OMP END MASTER
458#endif
459   
460! NetCDF output of received fields
461       IF (is_sequential) THEN
462          ndexcs(:) = 0
463          itau_w = itau_phy + itime + start_time * day_step_phy
464          DO i = 1, maxrecv
465            IF (inforecv(i)%action) THEN
466                CALL histwrite(nidcs,inforecv(i)%name,itau_w,tab_read_flds(:,:,i),nbp_lon*(nbp_lat),ndexcs)
467            ENDIF
468          ENDDO
469       ENDIF
470
471
472! Save each field in a 2D array.
473!$OMP MASTER
474       read_sst(:,:)     = tab_read_flds(:,:,idr_sisutw)  ! Sea surface temperature
475       read_sic(:,:)     = tab_read_flds(:,:,idr_icecov)  ! Sea ice concentration
476       read_alb_sic(:,:) = tab_read_flds(:,:,idr_icealw)  ! Albedo at sea ice
477       read_sit(:,:)     = tab_read_flds(:,:,idr_icetem)  ! Sea ice temperature
478       if (activate_ocean_skin >= 1) read_sss(:,:) = tab_read_flds(:,:,idr_sss)
479!$OMP END MASTER
480
481       IF (cpl_current) THEN
482
483! Transform the longitudes and latitudes on 2D arrays
484          CALL gather_omp(longitude_deg,rlon_mpi)
485          CALL gather_omp(latitude_deg,rlat_mpi)
486!$OMP MASTER
487          CALL Grid1DTo2D_mpi(rlon_mpi,tmp_lon)
488          CALL Grid1DTo2D_mpi(rlat_mpi,tmp_lat)
489
490! Transform the currents from cartesian to spheric coordinates
491! tmp_r0 should be zero
492          CALL geo2atm(nbp_lon, jj_nb, tab_read_flds(:,:,idr_curenx), &
493             tab_read_flds(:,:,idr_cureny), tab_read_flds(:,:,idr_curenz), &
494               tmp_lon, tmp_lat, &
495               read_u0(:,:), read_v0(:,:), tmp_r0(:,:))
496!$OMP END MASTER
497
498      ELSE
499          read_u0(:,:) = 0.
500          read_v0(:,:) = 0.
501      ENDIF
502
503       IF (carbon_cycle_cpl) THEN
504!$OMP MASTER
505           read_co2(:,:) = tab_read_flds(:,:,idr_oceco2) ! CO2 flux
506!$OMP END MASTER
507       ENDIF
508
509!*************************************************************************************
510!  Transform seaice fraction (read_sic : ocean-seaice mask) into global
511!  fraction (pctsrf : ocean-seaice-land-landice mask)
512!
513!*************************************************************************************
514       CALL cpl2gath(read_sic, read_sic1D, klon, unity)
515
516       pctsrf_old(:,:) = pctsrf(:,:)
517       DO i = 1, klon
518          ! treatment only of points with ocean and/or seaice
519          ! old land-ocean mask can not be changed
520          IF (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic) > 0.) THEN
521             pctsrf(i,is_sic) = (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic)) &
522                  * read_sic1D(i)
523             pctsrf(i,is_oce) = (pctsrf_old(i,is_oce) + pctsrf_old(i,is_sic)) &
524                  - pctsrf(i,is_sic)
525          ENDIF
526       ENDDO
527
528    ENDIF ! if time to receive
529
530  END SUBROUTINE cpl_receive_frac
531
532!
533!*************************************************************************************
534!
535
536  SUBROUTINE cpl_receive_ocean_fields(knon, knindex, tsurf_new, u0_new, &
537       v0_new, sss)
538!
539! This routine returns the field for the ocean that has been read from the coupler
540! (done earlier with cpl_receive_frac). The field is the temperature.
541! The temperature is transformed into 1D array with valid points from index 1 to knon.
542!
543    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, fco2_ocn_day
544    USE indice_sol_mod
545    use config_ocean_skin_m, only: activate_ocean_skin
546
547! Input arguments
548!*************************************************************************************
549    INTEGER, INTENT(IN)                     :: knon
550    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
551
552! Output arguments
553!*************************************************************************************
554    REAL, DIMENSION(klon), INTENT(OUT)      :: tsurf_new
555
556    REAL, INTENT(OUT):: sss(:) ! (klon)
557    ! bulk salinity of the surface layer of the ocean, in ppt
558
559    REAL, DIMENSION(klon), INTENT(OUT)      :: u0_new
560    REAL, DIMENSION(klon), INTENT(OUT)      :: v0_new
561
562! Local variables
563!*************************************************************************************
564    INTEGER                  :: i
565    INTEGER, DIMENSION(klon) :: index
566    REAL, DIMENSION(klon)    :: sic_new
567
568!*************************************************************************************
569! Transform read_sst into compressed 1D variable tsurf_new
570!
571!*************************************************************************************
572    CALL cpl2gath(read_sst, tsurf_new, knon, knindex)
573    if (activate_ocean_skin >= 1) CALL cpl2gath(read_sss, sss, knon, knindex)
574    CALL cpl2gath(read_sic, sic_new, knon, knindex)
575    CALL cpl2gath(read_u0, u0_new, knon, knindex)
576    CALL cpl2gath(read_v0, v0_new, knon, knindex)
577
578!*************************************************************************************
579! Transform read_co2 into uncompressed 1D variable fco2_ocn_day added directly in
580! the module carbon_cycle_mod
581!
582!*************************************************************************************
583    IF (carbon_cycle_cpl) THEN
584       DO i=1,klon
585          index(i)=i
586       ENDDO
587       CALL cpl2gath(read_co2, fco2_ocn_day, klon, index)
588    ENDIF
589
590!*************************************************************************************
591! The fields received from the coupler have to be weighted with the fraction of ocean
592! in relation to the total sea-ice+ocean
593!
594!*************************************************************************************
595    DO i=1, knon
596       tsurf_new(i) = tsurf_new(i)/(1. - sic_new(i))
597    ENDDO
598
599  END SUBROUTINE cpl_receive_ocean_fields
600
601!
602!*************************************************************************************
603!
604
605  SUBROUTINE cpl_receive_seaice_fields(knon, knindex, &
606       tsurf_new, alb_new, u0_new, v0_new)
607!
608! This routine returns the fields for the seaice that have been read from the coupler
609! (done earlier with cpl_receive_frac). These fields are the temperature and
610! albedo at sea ice surface and fraction of sea ice.
611! The fields are transformed into 1D arrays with valid points from index 1 to knon.
612!
613
614! Input arguments
615!*************************************************************************************
616    INTEGER, INTENT(IN)                     :: knon
617    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
618
619! Output arguments
620!*************************************************************************************
621    REAL, DIMENSION(klon), INTENT(OUT)      :: tsurf_new
622    REAL, DIMENSION(klon), INTENT(OUT)      :: alb_new
623    REAL, DIMENSION(klon), INTENT(OUT)      :: u0_new
624    REAL, DIMENSION(klon), INTENT(OUT)      :: v0_new
625
626! Local variables
627!*************************************************************************************
628    INTEGER               :: i
629    REAL, DIMENSION(klon) :: sic_new
630
631!*************************************************************************************
632! Transform fields read from coupler from 2D into compressed 1D variables
633!
634!*************************************************************************************
635    CALL cpl2gath(read_sit, tsurf_new, knon, knindex)
636    CALL cpl2gath(read_alb_sic, alb_new, knon, knindex)
637    CALL cpl2gath(read_sic, sic_new, knon, knindex)
638    CALL cpl2gath(read_u0, u0_new, knon, knindex)
639    CALL cpl2gath(read_v0, v0_new, knon, knindex)
640
641!*************************************************************************************
642! The fields received from the coupler have to be weighted with the sea-ice
643! concentration (in relation to the total sea-ice + ocean).
644!
645!*************************************************************************************
646    DO i= 1, knon
647       tsurf_new(i) = tsurf_new(i) / sic_new(i)
648       alb_new(i)   = alb_new(i)   / sic_new(i)
649    ENDDO
650
651  END SUBROUTINE cpl_receive_seaice_fields
652
653!
654!*************************************************************************************
655!
656
657  SUBROUTINE cpl_send_ocean_fields(itime, knon, knindex, &
658       swdown, lwdown, fluxlat, fluxsens, &
659       precip_rain, precip_snow, evap, tsurf, fder, albsol, taux, tauy, windsp,&
660       sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol, delta_sst, &
661       delta_sal, dTer, dSer, dt_ds)
662
663    ! This subroutine cumulates some fields for each time-step during
664    ! a coupling period. At last time-step in a coupling period the
665    ! fields are transformed to the grid accepted by the coupler. No
666    ! sending to the coupler will be done from here (it is done in
667    ! cpl_send_seaice_fields). Crucial hypothesis is that the surface
668    ! fractions do not change between coupling time-steps.
669
670    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl, co2_send
671    USE indice_sol_mod
672    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
673    use config_ocean_skin_m, only: activate_ocean_skin
674
675! Input arguments
676!*************************************************************************************
677    INTEGER, INTENT(IN)                     :: itime
678    INTEGER, INTENT(IN)                     :: knon
679    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
680    REAL, DIMENSION(klon), INTENT(IN)       :: swdown, lwdown
681    REAL, DIMENSION(klon), INTENT(IN)       :: fluxlat, fluxsens
682    REAL, DIMENSION(klon), INTENT(IN)       :: precip_rain, precip_snow
683    REAL, DIMENSION(klon), INTENT(IN)       :: evap, tsurf, fder, albsol
684    REAL, DIMENSION(klon), INTENT(IN)       :: taux, tauy, windsp
685    REAL, INTENT(IN):: sens_prec_liq(:), sens_prec_sol(:) ! (knon)
686    REAL, DIMENSION(klon), INTENT(IN)       :: lat_prec_liq, lat_prec_sol
687   
688    REAL, intent(in):: delta_sst(:) ! (knon)
689    ! Ocean-air interface temperature minus bulk SST, in
690    ! K. Defined only if activate_ocean_skin >= 1.
691
692    real, intent(in):: delta_sal(:) ! (knon)
693    ! Ocean-air interface salinity minus bulk salinity, in ppt.
694
695    REAL, intent(in):: dter(:) ! (knon)
696    ! Temperature variation in the diffusive microlayer, that is
697    ! ocean-air interface temperature minus subskin temperature. In
698    ! K.
699
700    REAL, intent(in):: dser(:) ! (knon)
701    ! Salinity variation in the diffusive microlayer, that is
702    ! ocean-air interface salinity minus subskin salinity. In ppt.
703
704    real, intent(in):: dt_ds(:) ! (knon)
705    ! (tks / tkt) * dTer, in K
706
707! Local variables
708!*************************************************************************************
709    INTEGER                                 :: cpl_index, ig
710    INTEGER                                 :: error, sum_error
711    CHARACTER(len = 25)                     :: modname = 'cpl_send_ocean_fields'
712    CHARACTER(len = 80)                     :: abort_message
713
714!*************************************************************************************
715! Start calculation
716! The ocean points are saved with second array index=1
717!
718!*************************************************************************************
719    cpl_index = 1
720
721!*************************************************************************************
722! Reset fields to zero in the beginning of a new coupling period
723!
724!*************************************************************************************
725    IF (MOD(itime, nexca) == 1) THEN
726       cpl_sols(1:knon,cpl_index) = 0.0
727       cpl_nsol(1:knon,cpl_index) = 0.0
728       cpl_rain(1:knon,cpl_index) = 0.0
729       cpl_snow(1:knon,cpl_index) = 0.0
730       cpl_evap(1:knon,cpl_index) = 0.0
731       cpl_tsol(1:knon,cpl_index) = 0.0
732       cpl_fder(1:knon,cpl_index) = 0.0
733       cpl_albe(1:knon,cpl_index) = 0.0
734       cpl_taux(1:knon,cpl_index) = 0.0
735       cpl_tauy(1:knon,cpl_index) = 0.0
736       cpl_windsp(1:knon,cpl_index) = 0.0
737       cpl_sens_rain(1:knon,cpl_index) = 0.0
738       cpl_sens_snow(1:knon,cpl_index) = 0.0
739       cpl_taumod(1:knon,cpl_index) = 0.0
740       IF (carbon_cycle_cpl) cpl_atm_co2(1:knon,cpl_index) = 0.0
741
742       if (activate_ocean_skin == 2) then
743          cpl_delta_sst = 0.
744          cpl_delta_sal = 0.
745          cpl_dter = 0.
746          cpl_dser = 0.
747          cpl_dt_ds = 0.
748       end if
749    ENDIF
750       
751!*************************************************************************************
752! Cumulate at each time-step
753!
754!*************************************************************************************   
755    DO ig = 1, knon
756       cpl_sols(ig,cpl_index) = cpl_sols(ig,cpl_index) + &
757            swdown(ig)      / REAL(nexca)
758       cpl_nsol(ig,cpl_index) = cpl_nsol(ig,cpl_index) + &
759            (lwdown(ig) + fluxlat(ig) +fluxsens(ig)) / REAL(nexca)
760       cpl_rain(ig,cpl_index) = cpl_rain(ig,cpl_index) + &
761            precip_rain(ig) / REAL(nexca)
762       cpl_snow(ig,cpl_index) = cpl_snow(ig,cpl_index) + &
763            precip_snow(ig) / REAL(nexca)
764       cpl_evap(ig,cpl_index) = cpl_evap(ig,cpl_index) + &
765            evap(ig)        / REAL(nexca)
766       cpl_tsol(ig,cpl_index) = cpl_tsol(ig,cpl_index) + &
767            tsurf(ig)       / REAL(nexca)
768       cpl_fder(ig,cpl_index) = cpl_fder(ig,cpl_index) + &
769            fder(ig)        / REAL(nexca)
770       cpl_albe(ig,cpl_index) = cpl_albe(ig,cpl_index) + &
771            albsol(ig)      / REAL(nexca)
772       cpl_taux(ig,cpl_index) = cpl_taux(ig,cpl_index) + &
773            taux(ig)        / REAL(nexca)
774       cpl_tauy(ig,cpl_index) = cpl_tauy(ig,cpl_index) + &
775            tauy(ig)        / REAL(nexca)     
776       cpl_windsp(ig,cpl_index) = cpl_windsp(ig,cpl_index) + &
777            windsp(ig)      / REAL(nexca)
778       cpl_sens_rain(ig,cpl_index) = cpl_sens_rain(ig,cpl_index) + &
779            sens_prec_liq(ig)      / REAL(nexca)
780       cpl_sens_snow(ig,cpl_index) = cpl_sens_snow(ig,cpl_index) + &
781            sens_prec_sol(ig)      / REAL(nexca)
782       cpl_taumod(ig,cpl_index) =   cpl_taumod(ig,cpl_index) + &
783          SQRT ( taux(ig)*taux(ig)+tauy(ig)*tauy(ig) ) / REAL (nexca)
784
785       IF (carbon_cycle_cpl) THEN
786          cpl_atm_co2(ig,cpl_index) = cpl_atm_co2(ig,cpl_index) + &
787               co2_send(knindex(ig))/ REAL(nexca)
788!!---OB: this is correct but why knindex ??
789       ENDIF
790
791       if (activate_ocean_skin == 2) then
792          cpl_delta_sst(ig) = cpl_delta_sst(ig) + delta_sst(ig) / REAL(nexca)
793          cpl_delta_sal(ig) = cpl_delta_sal(ig) + delta_sal(ig) / REAL(nexca)
794          cpl_dter(ig) = cpl_dter(ig) + dter(ig) / REAL(nexca)
795          cpl_dser(ig) = cpl_dser(ig) + dser(ig) / REAL(nexca)
796          cpl_dt_ds(ig) = cpl_dt_ds(ig) + dt_ds(ig) / REAL(nexca)
797       end if
798     ENDDO
799
800!*************************************************************************************
801! If the time-step corresponds to the end of coupling period the
802! fields are transformed to the 2D grid.
803! No sending to the coupler (it is done from cpl_send_seaice_fields).
804!
805!*************************************************************************************
806    IF (MOD(itime, nexca) == 0) THEN
807
808       IF (.NOT. ALLOCATED(cpl_sols2D)) THEN
809          sum_error = 0
810          ALLOCATE(cpl_sols2D(nbp_lon,jj_nb,2), stat=error)
811          sum_error = sum_error + error
812          ALLOCATE(cpl_nsol2D(nbp_lon,jj_nb,2), stat=error)
813          sum_error = sum_error + error
814          ALLOCATE(cpl_rain2D(nbp_lon,jj_nb,2), stat=error)
815          sum_error = sum_error + error
816          ALLOCATE(cpl_snow2D(nbp_lon,jj_nb,2), stat=error)
817          sum_error = sum_error + error
818          ALLOCATE(cpl_evap2D(nbp_lon,jj_nb,2), stat=error)
819          sum_error = sum_error + error
820          ALLOCATE(cpl_tsol2D(nbp_lon,jj_nb,2), stat=error)
821          sum_error = sum_error + error
822          ALLOCATE(cpl_fder2D(nbp_lon,jj_nb,2), stat=error)
823          sum_error = sum_error + error
824          ALLOCATE(cpl_albe2D(nbp_lon,jj_nb,2), stat=error)
825          sum_error = sum_error + error
826          ALLOCATE(cpl_taux2D(nbp_lon,jj_nb,2), stat=error)
827          sum_error = sum_error + error
828          ALLOCATE(cpl_tauy2D(nbp_lon,jj_nb,2), stat=error)
829          sum_error = sum_error + error
830          ALLOCATE(cpl_windsp2D(nbp_lon,jj_nb), stat=error)
831          sum_error = sum_error + error
832          ALLOCATE(cpl_sens_rain2D(nbp_lon,jj_nb,2), stat=error)
833          sum_error = sum_error + error
834          ALLOCATE(cpl_sens_snow2D(nbp_lon,jj_nb,2), stat=error)
835          sum_error = sum_error + error
836          ALLOCATE(cpl_taumod2D(nbp_lon,jj_nb,2), stat=error)
837          sum_error = sum_error + error
838         
839          IF (carbon_cycle_cpl) THEN
840             ALLOCATE(cpl_atm_co22D(nbp_lon,jj_nb), stat=error)
841             sum_error = sum_error + error
842          ENDIF
843
844          if (activate_ocean_skin == 2) then
845             ALLOCATE(cpl_delta_sst_2D(nbp_lon, jj_nb), &
846                  cpl_delta_sal_2D(nbp_lon, jj_nb), &
847                  cpl_dter_2D(nbp_lon, jj_nb), cpl_dser_2D(nbp_lon, jj_nb), &
848                  cpl_dt_ds_2D(nbp_lon, jj_nb), stat = error)
849             sum_error = sum_error + error
850          end if
851
852          IF (sum_error /= 0) THEN
853             abort_message='Pb allocation variables couplees pour l''ecriture'
854             CALL abort_physic(modname,abort_message,1)
855          ENDIF
856       ENDIF
857       
858
859       CALL gath2cpl(cpl_sols(:,cpl_index), cpl_sols2D(:,:,cpl_index), &
860            knon, knindex)
861
862       CALL gath2cpl(cpl_nsol(:,cpl_index), cpl_nsol2D(:,:,cpl_index), &
863            knon, knindex)
864
865       CALL gath2cpl(cpl_rain(:,cpl_index), cpl_rain2D(:,:,cpl_index), &
866            knon, knindex)
867
868       CALL gath2cpl(cpl_snow(:,cpl_index), cpl_snow2D(:,:,cpl_index), &
869            knon, knindex)
870
871       CALL gath2cpl(cpl_evap(:,cpl_index), cpl_evap2D(:,:,cpl_index), &
872            knon, knindex)
873
874! cpl_tsol2D(:,:,:) not used!
875       CALL gath2cpl(cpl_tsol(:,cpl_index), cpl_tsol2D(:,:, cpl_index), &
876            knon, knindex)
877
878! cpl_fder2D(:,:,1) not used, only cpl_fder(:,:,2)!
879       CALL gath2cpl(cpl_fder(:,cpl_index), cpl_fder2D(:,:,cpl_index), &
880            knon, knindex)
881
882! cpl_albe2D(:,:,:) not used!
883       CALL gath2cpl(cpl_albe(:,cpl_index), cpl_albe2D(:,:,cpl_index), &
884            knon, knindex)
885
886       CALL gath2cpl(cpl_taux(:,cpl_index), cpl_taux2D(:,:,cpl_index), &
887            knon, knindex)
888
889       CALL gath2cpl(cpl_tauy(:,cpl_index), cpl_tauy2D(:,:,cpl_index), &
890            knon, knindex)
891
892       CALL gath2cpl(cpl_windsp(:,cpl_index), cpl_windsp2D(:,:), &
893            knon, knindex)
894
895       CALL gath2cpl(cpl_sens_rain(:,cpl_index), cpl_sens_rain2D(:,:,cpl_index), &
896            knon, knindex)
897
898       CALL gath2cpl(cpl_sens_snow(:,cpl_index), cpl_sens_snow2D(:,:,cpl_index), &
899            knon, knindex)
900
901       CALL gath2cpl(cpl_taumod(:,cpl_index), cpl_taumod2D(:,:,cpl_index), &
902            knon, knindex)
903
904       IF (carbon_cycle_cpl) &
905            CALL gath2cpl(cpl_atm_co2(:,cpl_index), cpl_atm_co22D(:,:), knon, knindex)
906       if (activate_ocean_skin == 2) then
907          CALL gath2cpl(cpl_delta_sst, cpl_delta_sst_2D, knon, knindex)
908          CALL gath2cpl(cpl_delta_sal, cpl_delta_sal_2D, knon, knindex)
909          CALL gath2cpl(cpl_dter, cpl_dter_2D, knon, knindex)
910          CALL gath2cpl(cpl_dser, cpl_dser_2D, knon, knindex)
911          CALL gath2cpl(cpl_dt_ds, cpl_dt_ds_2D, knon, knindex)
912       end if
913    ENDIF
914
915  END SUBROUTINE cpl_send_ocean_fields
916
917!
918!*************************************************************************************
919!
920
921  SUBROUTINE cpl_send_seaice_fields(itime, dtime, knon, knindex, &
922       pctsrf, lafin, rlon, rlat, &
923       swdown, lwdown, fluxlat, fluxsens, &
924       precip_rain, precip_snow, evap, tsurf, fder, albsol, taux, tauy,&
925       sens_prec_liq, sens_prec_sol, lat_prec_liq, lat_prec_sol)
926!
927! This subroutine cumulates some fields for each time-step during a coupling
928! period. At last time-step in a coupling period the fields are transformed to the
929! grid accepted by the coupler. All fields for all types of surfaces are sent to
930! the coupler.
931!
932    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl
933    USE indice_sol_mod
934    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
935
936! Input arguments
937!*************************************************************************************
938    INTEGER, INTENT(IN)                     :: itime
939    INTEGER, INTENT(IN)                     :: knon
940    INTEGER, DIMENSION(klon), INTENT(IN)    :: knindex
941    REAL, INTENT(IN)                        :: dtime
942    REAL, DIMENSION(klon), INTENT(IN)       :: rlon, rlat
943    REAL, DIMENSION(klon), INTENT(IN)       :: swdown, lwdown
944    REAL, DIMENSION(klon), INTENT(IN)       :: fluxlat, fluxsens
945    REAL, DIMENSION(klon), INTENT(IN)       :: precip_rain, precip_snow
946    REAL, DIMENSION(klon), INTENT(IN)       :: evap, tsurf, fder
947    REAL, DIMENSION(klon), INTENT(IN)       :: albsol, taux, tauy
948    REAL, DIMENSION(klon,nbsrf), INTENT(IN) :: pctsrf
949    REAL, INTENT(IN):: sens_prec_liq(:), sens_prec_sol(:) ! (knon)
950    REAL, DIMENSION(klon), INTENT(IN)       :: lat_prec_liq, lat_prec_sol
951    LOGICAL, INTENT(IN)                     :: lafin
952
953! Local variables
954!*************************************************************************************
955    INTEGER                                 :: cpl_index, ig
956    INTEGER                                 :: error, sum_error
957    CHARACTER(len = 25)                     :: modname = 'cpl_send_seaice_fields'
958    CHARACTER(len = 80)                     :: abort_message
959    REAL, DIMENSION(klon)                   :: cpl_fder_tmp
960
961!*************************************************************************************
962! Start calulation
963! The sea-ice points are saved with second array index=2
964!
965!*************************************************************************************
966    cpl_index = 2
967
968!*************************************************************************************
969! Reset fields to zero in the beginning of a new coupling period
970!
971!*************************************************************************************
972    IF (MOD(itime, nexca) == 1) THEN
973       cpl_sols(1:knon,cpl_index) = 0.0
974       cpl_nsol(1:knon,cpl_index) = 0.0
975       cpl_rain(1:knon,cpl_index) = 0.0
976       cpl_snow(1:knon,cpl_index) = 0.0
977       cpl_evap(1:knon,cpl_index) = 0.0
978       cpl_tsol(1:knon,cpl_index) = 0.0
979       cpl_fder(1:knon,cpl_index) = 0.0
980       cpl_albe(1:knon,cpl_index) = 0.0
981       cpl_taux(1:knon,cpl_index) = 0.0
982       cpl_tauy(1:knon,cpl_index) = 0.0
983       cpl_sens_rain(1:knon,cpl_index) = 0.0
984       cpl_sens_snow(1:knon,cpl_index) = 0.0
985       cpl_taumod(1:knon,cpl_index) = 0.0
986    ENDIF
987       
988!*************************************************************************************
989! Cumulate at each time-step
990!
991!*************************************************************************************   
992    DO ig = 1, knon
993       cpl_sols(ig,cpl_index) = cpl_sols(ig,cpl_index) + &
994            swdown(ig)      / REAL(nexca)
995       cpl_nsol(ig,cpl_index) = cpl_nsol(ig,cpl_index) + &
996            (lwdown(ig) + fluxlat(ig) +fluxsens(ig)) / REAL(nexca)
997       cpl_rain(ig,cpl_index) = cpl_rain(ig,cpl_index) + &
998            precip_rain(ig) / REAL(nexca)
999       cpl_snow(ig,cpl_index) = cpl_snow(ig,cpl_index) + &
1000            precip_snow(ig) / REAL(nexca)
1001       cpl_evap(ig,cpl_index) = cpl_evap(ig,cpl_index) + &
1002            evap(ig)        / REAL(nexca)
1003       cpl_tsol(ig,cpl_index) = cpl_tsol(ig,cpl_index) + &
1004            tsurf(ig)       / REAL(nexca)
1005       cpl_fder(ig,cpl_index) = cpl_fder(ig,cpl_index) + &
1006            fder(ig)        / REAL(nexca)
1007       cpl_albe(ig,cpl_index) = cpl_albe(ig,cpl_index) + &
1008            albsol(ig)      / REAL(nexca)
1009       cpl_taux(ig,cpl_index) = cpl_taux(ig,cpl_index) + &
1010            taux(ig)        / REAL(nexca)
1011       cpl_tauy(ig,cpl_index) = cpl_tauy(ig,cpl_index) + &
1012            tauy(ig)        / REAL(nexca)     
1013       cpl_sens_rain(ig,cpl_index) = cpl_sens_rain(ig,cpl_index) + &
1014            sens_prec_liq(ig)      / REAL(nexca)
1015       cpl_sens_snow(ig,cpl_index) = cpl_sens_snow(ig,cpl_index) + &
1016            sens_prec_sol(ig)      / REAL(nexca)
1017       cpl_taumod(ig,cpl_index) = cpl_taumod(ig,cpl_index) + &
1018            SQRT ( taux(ig)*taux(ig)+tauy(ig)*tauy(ig) ) / REAL(nexca)
1019    ENDDO
1020
1021!*************************************************************************************
1022! If the time-step corresponds to the end of coupling period the
1023! fields are transformed to the 2D grid and all fields are sent to coupler.
1024!
1025!*************************************************************************************
1026    IF (MOD(itime, nexca) == 0) THEN
1027       IF (.NOT. ALLOCATED(cpl_sols2D)) THEN
1028          sum_error = 0
1029          ALLOCATE(cpl_sols2D(nbp_lon,jj_nb,2), stat=error)
1030          sum_error = sum_error + error
1031          ALLOCATE(cpl_nsol2D(nbp_lon,jj_nb,2), stat=error)
1032          sum_error = sum_error + error
1033          ALLOCATE(cpl_rain2D(nbp_lon,jj_nb,2), stat=error)
1034          sum_error = sum_error + error
1035          ALLOCATE(cpl_snow2D(nbp_lon,jj_nb,2), stat=error)
1036          sum_error = sum_error + error
1037          ALLOCATE(cpl_evap2D(nbp_lon,jj_nb,2), stat=error)
1038          sum_error = sum_error + error
1039          ALLOCATE(cpl_tsol2D(nbp_lon,jj_nb,2), stat=error)
1040          sum_error = sum_error + error
1041          ALLOCATE(cpl_fder2D(nbp_lon,jj_nb,2), stat=error)
1042          sum_error = sum_error + error
1043          ALLOCATE(cpl_albe2D(nbp_lon,jj_nb,2), stat=error)
1044          sum_error = sum_error + error
1045          ALLOCATE(cpl_taux2D(nbp_lon,jj_nb,2), stat=error)
1046          sum_error = sum_error + error
1047          ALLOCATE(cpl_tauy2D(nbp_lon,jj_nb,2), stat=error)
1048          sum_error = sum_error + error
1049          ALLOCATE(cpl_windsp2D(nbp_lon,jj_nb), stat=error)
1050          sum_error = sum_error + error
1051          ALLOCATE(cpl_sens_rain2D(nbp_lon,jj_nb,2), stat=error)
1052          sum_error = sum_error + error
1053          ALLOCATE(cpl_sens_snow2D(nbp_lon,jj_nb,2), stat=error)
1054          sum_error = sum_error + error
1055          ALLOCATE(cpl_taumod2D(nbp_lon,jj_nb,2), stat=error)
1056          sum_error = sum_error + error
1057
1058          IF (carbon_cycle_cpl) THEN
1059             ALLOCATE(cpl_atm_co22D(nbp_lon,jj_nb), stat=error)
1060             sum_error = sum_error + error
1061          ENDIF
1062
1063          IF (sum_error /= 0) THEN
1064             abort_message='Pb allocation variables couplees pour l''ecriture'
1065             CALL abort_physic(modname,abort_message,1)
1066          ENDIF
1067       ENDIF
1068
1069       CALL gath2cpl(cpl_sols(:,cpl_index), cpl_sols2D(:,:,cpl_index), &
1070            knon, knindex)
1071
1072       CALL gath2cpl(cpl_nsol(:,cpl_index), cpl_nsol2D(:,:,cpl_index), &
1073            knon, knindex)
1074
1075       CALL gath2cpl(cpl_rain(:,cpl_index), cpl_rain2D(:,:,cpl_index), &
1076            knon, knindex)
1077
1078       CALL gath2cpl(cpl_snow(:,cpl_index), cpl_snow2D(:,:,cpl_index), &
1079            knon, knindex)
1080
1081       CALL gath2cpl(cpl_evap(:,cpl_index), cpl_evap2D(:,:,cpl_index), &
1082            knon, knindex)
1083
1084! cpl_tsol2D(:,:,:) not used!
1085       CALL gath2cpl(cpl_tsol(:,cpl_index), cpl_tsol2D(:,:, cpl_index), &
1086            knon, knindex)
1087
1088       ! Set default value and decompress before gath2cpl
1089       cpl_fder_tmp(:) = -20.
1090       DO ig = 1, knon
1091          cpl_fder_tmp(knindex(ig))=cpl_fder(ig,cpl_index)
1092       ENDDO
1093       CALL gath2cpl(cpl_fder_tmp(:), cpl_fder2D(:,:,cpl_index), &
1094            klon, unity)
1095
1096! cpl_albe2D(:,:,:) not used!
1097       CALL gath2cpl(cpl_albe(:,cpl_index), cpl_albe2D(:,:,cpl_index), &
1098            knon, knindex)
1099
1100       CALL gath2cpl(cpl_taux(:,cpl_index), cpl_taux2D(:,:,cpl_index), &
1101            knon, knindex)
1102
1103       CALL gath2cpl(cpl_tauy(:,cpl_index), cpl_tauy2D(:,:,cpl_index), &
1104            knon, knindex)
1105
1106       CALL gath2cpl(cpl_sens_rain(:,cpl_index), cpl_sens_rain2D(:,:,cpl_index), &
1107            knon, knindex)
1108
1109       CALL gath2cpl(cpl_sens_snow(:,cpl_index), cpl_sens_snow2D(:,:,cpl_index), &
1110            knon, knindex)
1111
1112       CALL gath2cpl(cpl_taumod(:,cpl_index), cpl_taumod2D(:,:,cpl_index), &
1113            knon, knindex)
1114
1115       ! Send all fields
1116       CALL cpl_send_all(itime, dtime, pctsrf, lafin, rlon, rlat)
1117    ENDIF
1118
1119  END SUBROUTINE cpl_send_seaice_fields
1120
1121!
1122!*************************************************************************************
1123!
1124
1125  SUBROUTINE cpl_send_land_fields(itime, knon, knindex, rriv_in, rcoa_in)
1126!
1127! This subroutine cumulates some fields for each time-step during a coupling
1128! period. At last time-step in a coupling period the fields are transformed to the
1129! grid accepted by the coupler. No sending to the coupler will be done from here
1130! (it is done in cpl_send_seaice_fields).
1131!
1132    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
1133
1134! Input arguments
1135!*************************************************************************************
1136    INTEGER, INTENT(IN)                       :: itime
1137    INTEGER, INTENT(IN)                       :: knon
1138    INTEGER, DIMENSION(klon), INTENT(IN)      :: knindex
1139    REAL, DIMENSION(klon), INTENT(IN)         :: rriv_in
1140    REAL, DIMENSION(klon), INTENT(IN)         :: rcoa_in
1141
1142! Local variables
1143!*************************************************************************************
1144    REAL, DIMENSION(nbp_lon,jj_nb)             :: rriv2D
1145    REAL, DIMENSION(nbp_lon,jj_nb)             :: rcoa2D
1146
1147!*************************************************************************************
1148! Rearrange fields in 2D variables
1149! First initialize to zero to avoid unvalid points causing problems
1150!
1151!*************************************************************************************
1152!$OMP MASTER
1153    rriv2D(:,:) = 0.0
1154    rcoa2D(:,:) = 0.0
1155!$OMP END MASTER
1156    CALL gath2cpl(rriv_in, rriv2D, knon, knindex)
1157    CALL gath2cpl(rcoa_in, rcoa2D, knon, knindex)
1158
1159!*************************************************************************************
1160! Reset cumulated fields to zero in the beginning of a new coupling period
1161!
1162!*************************************************************************************
1163    IF (MOD(itime, nexca) == 1) THEN
1164!$OMP MASTER
1165       cpl_rriv2D(:,:) = 0.0
1166       cpl_rcoa2D(:,:) = 0.0
1167!$OMP END MASTER
1168    ENDIF
1169
1170!*************************************************************************************
1171! Cumulate : Following fields should be cumulated at each time-step
1172!
1173!*************************************************************************************   
1174!$OMP MASTER
1175    cpl_rriv2D(:,:) = cpl_rriv2D(:,:) + rriv2D(:,:) / REAL(nexca)
1176    cpl_rcoa2D(:,:) = cpl_rcoa2D(:,:) + rcoa2D(:,:) / REAL(nexca)
1177!$OMP END MASTER
1178
1179  END SUBROUTINE cpl_send_land_fields
1180
1181!
1182!*************************************************************************************
1183!
1184
1185  SUBROUTINE cpl_send_landice_fields(itime, knon, knindex, rlic_in, rlic_in_frac)
1186! This subroutine cumulates the field for melting ice for each time-step
1187! during a coupling period. This routine will not send to coupler. Sending
1188! will be done in cpl_send_seaice_fields.
1189!
1190
1191    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
1192
1193! Input varibales
1194!*************************************************************************************
1195    INTEGER, INTENT(IN)                       :: itime
1196    INTEGER, INTENT(IN)                       :: knon
1197    INTEGER, DIMENSION(klon), INTENT(IN)      :: knindex
1198    REAL, DIMENSION(klon), INTENT(IN)         :: rlic_in
1199    REAL, DIMENSION(klon), INTENT(IN)         :: rlic_in_frac  ! Fraction for continental ice, can be equal to
1200                                                               ! pctsrf(:,is_lic) or not, depending on landice_opt
1201   
1202
1203! Local varibales
1204!*************************************************************************************
1205    REAL, DIMENSION(nbp_lon,jj_nb)             :: rlic2D
1206
1207!*************************************************************************************
1208! Rearrange field in a 2D variable
1209! First initialize to zero to avoid unvalid points causing problems
1210!
1211!*************************************************************************************
1212!$OMP MASTER
1213    rlic2D(:,:) = 0.0
1214!$OMP END MASTER
1215    CALL gath2cpl(rlic_in, rlic2D, knon, knindex)
1216    CALL gath2cpl(rlic_in_frac(:), rlic_in_frac2D(:,:), knon, knindex)
1217!*************************************************************************************
1218! Reset field to zero in the beginning of a new coupling period
1219!
1220!*************************************************************************************
1221    IF (MOD(itime, nexca) == 1) THEN
1222!$OMP MASTER
1223       cpl_rlic2D(:,:) = 0.0
1224!$OMP END MASTER
1225    ENDIF
1226
1227!*************************************************************************************
1228! Cumulate : Melting ice should be cumulated at each time-step
1229!
1230!*************************************************************************************   
1231!$OMP MASTER
1232    cpl_rlic2D(:,:) = cpl_rlic2D(:,:) + rlic2D(:,:) / REAL(nexca)
1233!$OMP END MASTER
1234
1235  END SUBROUTINE cpl_send_landice_fields
1236
1237!
1238!*************************************************************************************
1239!
1240
1241  SUBROUTINE cpl_send_all(itime, dtime, pctsrf, lafin, rlon, rlat)
1242! This routine will send fields for all different surfaces to the coupler.
1243! This subroutine should be executed after calculations by the last surface(sea-ice),
1244! all calculations at the different surfaces have to be done before.
1245!   
1246    USE surface_data
1247    USE carbon_cycle_mod, ONLY : carbon_cycle_cpl
1248    USE indice_sol_mod
1249    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
1250    USE time_phylmdz_mod, ONLY: start_time, itau_phy
1251    use config_ocean_skin_m, only: activate_ocean_skin
1252! Some includes
1253!   
1254! Input arguments
1255!*************************************************************************************
1256    INTEGER, INTENT(IN)                                  :: itime
1257    REAL, INTENT(IN)                                     :: dtime
1258    REAL, DIMENSION(klon), INTENT(IN)                    :: rlon, rlat
1259    REAL, DIMENSION(klon,nbsrf), INTENT(IN)              :: pctsrf
1260    LOGICAL, INTENT(IN)                                  :: lafin
1261   
1262! Local variables
1263!*************************************************************************************
1264    INTEGER                                              :: error, sum_error, i,j,k
1265    INTEGER                                              :: itau_w
1266    INTEGER                                              :: time_sec
1267    INTEGER, DIMENSION(nbp_lon*(nbp_lat))                      :: ndexct
1268    REAL                                                 :: Up, Down
1269    REAL, DIMENSION(nbp_lon, jj_nb)                          :: tmp_lon, tmp_lat
1270    REAL, DIMENSION(nbp_lon, jj_nb, 4)                       :: pctsrf2D
1271    REAL, DIMENSION(nbp_lon, jj_nb)                          :: deno
1272    CHARACTER(len = 20)                                  :: modname = 'cpl_send_all'
1273    CHARACTER(len = 80)                                  :: abort_message
1274   
1275! Variables with fields to coupler
1276    REAL, DIMENSION(nbp_lon, jj_nb)                          :: tmp_taux
1277    REAL, DIMENSION(nbp_lon, jj_nb)                          :: tmp_tauy
1278    REAL, DIMENSION(nbp_lon, jj_nb)                          :: tmp_calv
1279! Table with all fields to send to coupler
1280    REAL, DIMENSION(nbp_lon, jj_nb, maxsend)                 :: tab_flds
1281    REAL, DIMENSION(klon_mpi)                                :: rlon_mpi, rlat_mpi
1282    REAL  :: calving(nb_zone_calving)
1283    REAL  :: calving_glo(nb_zone_calving)
1284   
1285#ifdef CPP_MPI
1286    INCLUDE 'mpif.h'
1287    INTEGER, DIMENSION(MPI_STATUS_SIZE)                  :: status
1288#endif
1289
1290! End definitions
1291!*************************************************************************************
1292   
1293
1294
1295!*************************************************************************************
1296! All fields are stored in a table tab_flds(:,:,:)
1297! First store the fields which are already on the right format
1298!
1299!*************************************************************************************
1300!$OMP MASTER
1301    tab_flds(:,:,ids_windsp) = cpl_windsp2D(:,:)
1302    tab_flds(:,:,ids_shfice) = cpl_sols2D(:,:,2)
1303    tab_flds(:,:,ids_nsfice) = cpl_nsol2D(:,:,2)
1304    tab_flds(:,:,ids_dflxdt) = cpl_fder2D(:,:,2)
1305    tab_flds(:,:,ids_qraioc) = cpl_sens_rain2D(:,:,1)
1306    tab_flds(:,:,ids_qsnooc) = cpl_sens_snow2D(:,:,1)
1307    tab_flds(:,:,ids_qraiic) = cpl_sens_rain2D(:,:,2)
1308    tab_flds(:,:,ids_qsnoic) = cpl_sens_snow2D(:,:,2)
1309
1310    if (activate_ocean_skin == 2) then
1311       tab_flds(:, :, ids_delta_sst) = cpl_delta_sst_2D
1312       tab_flds(:, :, ids_delta_sal) = cpl_delta_sal_2D
1313       tab_flds(:, :, ids_dter) = cpl_dter_2D
1314       tab_flds(:, :, ids_dser) = cpl_dser_2D
1315       tab_flds(:, :, ids_dt_ds) = cpl_dt_ds_2D
1316    end if
1317   
1318    IF (version_ocean=='nemo') THEN
1319       tab_flds(:,:,ids_liqrun) = (cpl_rriv2D(:,:) + cpl_rcoa2D(:,:))
1320       IF (carbon_cycle_cpl) tab_flds(:,:,ids_atmco2)=cpl_atm_co22D(:,:)
1321    ELSE IF (version_ocean=='opa8') THEN
1322       tab_flds(:,:,ids_shfoce) = cpl_sols2D(:,:,1)
1323       tab_flds(:,:,ids_nsfoce) = cpl_nsol2D(:,:,1)
1324       tab_flds(:,:,ids_icevap) = cpl_evap2D(:,:,2)
1325       tab_flds(:,:,ids_ocevap) = cpl_evap2D(:,:,1)
1326       tab_flds(:,:,ids_runcoa) = cpl_rcoa2D(:,:)
1327       tab_flds(:,:,ids_rivflu) = cpl_rriv2D(:,:)
1328    ENDIF
1329
1330!*************************************************************************************
1331! Transform the fraction of sub-surfaces from 1D to 2D array
1332!
1333!*************************************************************************************
1334    pctsrf2D(:,:,:) = 0.
1335!$OMP END MASTER
1336    CALL gath2cpl(pctsrf(:,is_oce), pctsrf2D(:,:,is_oce), klon, unity)
1337    CALL gath2cpl(pctsrf(:,is_sic), pctsrf2D(:,:,is_sic), klon, unity)
1338
1339
1340
1341
1342!*************************************************************************************
1343! Calculate the average calving per latitude
1344! Store calving in tab_flds(:,:,19)
1345!
1346!*************************************************************************************     
1347    IF (is_omp_root) THEN
1348
1349      IF (cpl_old_calving) THEN   ! use old calving
1350
1351        DO j = 1, jj_nb
1352           tmp_calv(:,j) = DOT_PRODUCT (cpl_rlic2D(1:nbp_lon,j), &
1353                rlic_in_frac2D(1:nbp_lon,j)) / REAL(nbp_lon)
1354        ENDDO
1355   
1356   
1357        IF (is_parallel) THEN
1358           IF (.NOT. is_north_pole_dyn) THEN
1359#ifdef CPP_MPI
1360              CALL MPI_RECV(Up,1,MPI_REAL_LMDZ,mpi_rank-1,1234,COMM_LMDZ_PHY,status,error)
1361              CALL MPI_SEND(tmp_calv(1,1),1,MPI_REAL_LMDZ,mpi_rank-1,1234,COMM_LMDZ_PHY,error)
1362#endif
1363           ENDIF
1364       
1365           IF (.NOT. is_south_pole_dyn) THEN
1366#ifdef CPP_MPI
1367              CALL MPI_SEND(tmp_calv(1,jj_nb),1,MPI_REAL_LMDZ,mpi_rank+1,1234,COMM_LMDZ_PHY,error)
1368              CALL MPI_RECV(down,1,MPI_REAL_LMDZ,mpi_rank+1,1234,COMM_LMDZ_PHY,status,error)
1369#endif
1370           ENDIF
1371         
1372           IF (.NOT. is_north_pole_dyn .AND. ii_begin /=1) THEN
1373              Up=Up+tmp_calv(nbp_lon,1)
1374              tmp_calv(:,1)=Up
1375           ENDIF
1376           
1377           IF (.NOT. is_south_pole_dyn .AND. ii_end /= nbp_lon) THEN
1378              Down=Down+tmp_calv(1,jj_nb)
1379              tmp_calv(:,jj_nb)=Down
1380           ENDIF
1381        ENDIF
1382        tab_flds(:,:,ids_calvin) = tmp_calv(:,:)
1383
1384      ELSE
1385         ! cpl_old_calving=FALSE
1386         ! To be used with new method for calculation of coupling weights
1387         DO k=1,nb_zone_calving
1388            calving(k)=0
1389            DO j = 1, jj_nb
1390               calving(k)= calving(k)+DOT_PRODUCT(cpl_rlic2D(:,j)*area_calving(:,j,k),rlic_in_frac2D(:,j))
1391            ENDDO
1392         ENDDO
1393         
1394#ifdef CPP_MPI
1395         CALL MPI_ALLREDUCE(calving, calving_glo, nb_zone_calving, MPI_REAL_LMDZ, MPI_SUM, COMM_LMDZ_PHY, error)
1396#endif
1397         
1398         tab_flds(:,:,ids_calvin) = 0
1399         DO k=1,nb_zone_calving
1400            IF (ind_calving(k)>0 ) THEN
1401               j=(ind_calving(k)-1)/nbp_lon + 1
1402               i=MOD(ind_calving(k)-1,nbp_lon)+1
1403               tab_flds(i,j,ids_calvin) = calving_glo(k)
1404            ENDIF
1405         ENDDO
1406         
1407      ENDIF
1408     
1409!*************************************************************************************
1410! Calculate total flux for snow, rain and wind with weighted addition using the
1411! fractions of ocean and seaice.
1412!
1413!*************************************************************************************   
1414       ! fraction oce+seaice
1415       deno =  pctsrf2D(:,:,is_oce) + pctsrf2D(:,:,is_sic)
1416
1417       IF (version_ocean=='nemo') THEN
1418          tab_flds(:,:,ids_shftot)  = 0.0
1419          tab_flds(:,:,ids_nsftot) = 0.0
1420          tab_flds(:,:,ids_totrai) = 0.0
1421          tab_flds(:,:,ids_totsno) = 0.0
1422          tab_flds(:,:,ids_toteva) = 0.0
1423          tab_flds(:,:,ids_taumod) = 0.0
1424 
1425          tmp_taux(:,:)    = 0.0
1426          tmp_tauy(:,:)    = 0.0
1427          ! For all valid grid cells containing some fraction of ocean or sea-ice
1428          WHERE ( deno(:,:) /= 0 )
1429             tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1430                  cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1431             tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1432                  cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1433
1434             tab_flds(:,:,ids_shftot) = cpl_sols2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1435                  cpl_sols2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1436             tab_flds(:,:,ids_nsftot) = cpl_nsol2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1437                  cpl_nsol2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1438             tab_flds(:,:,ids_totrai) = cpl_rain2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1439                  cpl_rain2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1440             tab_flds(:,:,ids_totsno) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1441                  cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1442             tab_flds(:,:,ids_toteva) = cpl_evap2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1443                  cpl_evap2D(:,:,2)  * pctsrf2D(:,:,is_sic) / deno(:,:)
1444             tab_flds(:,:,ids_taumod) = cpl_taumod2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1445                  cpl_taumod2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1446             
1447         ENDWHERE
1448
1449          tab_flds(:,:,ids_icevap) = cpl_evap2D(:,:,2)
1450         
1451       ELSE IF (version_ocean=='opa8') THEN
1452          ! Store fields for rain and snow in tab_flds(:,:,15) and tab_flds(:,:,16)
1453          tab_flds(:,:,ids_totrai) = 0.0
1454          tab_flds(:,:,ids_totsno) = 0.0
1455          tmp_taux(:,:)    = 0.0
1456          tmp_tauy(:,:)    = 0.0
1457          ! For all valid grid cells containing some fraction of ocean or sea-ice
1458          WHERE ( deno(:,:) /= 0 )
1459             tab_flds(:,:,ids_totrai) = cpl_rain2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1460                  cpl_rain2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1461             tab_flds(:,:,ids_totsno) = cpl_snow2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1462                  cpl_snow2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1463             
1464             tmp_taux = cpl_taux2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1465                  cpl_taux2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1466             tmp_tauy = cpl_tauy2D(:,:,1) * pctsrf2D(:,:,is_oce) / deno(:,:) +    &
1467                  cpl_tauy2D(:,:,2) * pctsrf2D(:,:,is_sic) / deno(:,:)
1468          ENDWHERE
1469       ENDIF
1470
1471    ENDIF ! is_omp_root
1472 
1473!*************************************************************************************
1474! Transform the wind components from local atmospheric 2D coordinates to geocentric
1475! 3D coordinates.
1476! Store the resulting wind components in tab_flds(:,:,1:6)
1477!*************************************************************************************
1478
1479! Transform the longitudes and latitudes on 2D arrays
1480   
1481    CALL gather_omp(rlon,rlon_mpi)
1482    CALL gather_omp(rlat,rlat_mpi)
1483!$OMP MASTER
1484    CALL Grid1DTo2D_mpi(rlon_mpi,tmp_lon)
1485    CALL Grid1DTo2D_mpi(rlat_mpi,tmp_lat)
1486!$OMP END MASTER   
1487
1488    IF (is_sequential) THEN
1489       IF (is_north_pole_dyn) tmp_lon(:,1)     = tmp_lon(:,2)
1490       IF (is_south_pole_dyn) tmp_lon(:,nbp_lat) = tmp_lon(:,nbp_lat-1)
1491    ENDIF
1492     
1493! NetCDF output of the wind before transformation of coordinate system
1494    IF (is_sequential) THEN
1495       ndexct(:) = 0
1496       itau_w = itau_phy + itime + start_time * day_step_phy
1497       CALL histwrite(nidct,'tauxe',itau_w,tmp_taux,nbp_lon*(nbp_lat),ndexct)
1498       CALL histwrite(nidct,'tauyn',itau_w,tmp_tauy,nbp_lon*(nbp_lat),ndexct)
1499       CALL histwrite(nidct,'tmp_lon',itau_w,tmp_lon,nbp_lon*(nbp_lat),ndexct)
1500       CALL histwrite(nidct,'tmp_lat',itau_w,tmp_lat,nbp_lon*(nbp_lat),ndexct)
1501    ENDIF
1502
1503! Transform the wind from spherical atmospheric 2D coordinates to geocentric
1504! cartesian 3D coordinates
1505!$OMP MASTER
1506    CALL atm2geo (nbp_lon, jj_nb, tmp_taux, tmp_tauy, tmp_lon, tmp_lat, &
1507         tab_flds(:,:,ids_tauxxu), tab_flds(:,:,ids_tauyyu), tab_flds(:,:,ids_tauzzu) )
1508   
1509    tab_flds(:,:,ids_tauxxv)  = tab_flds(:,:,ids_tauxxu)
1510    tab_flds(:,:,ids_tauyyv)  = tab_flds(:,:,ids_tauyyu)
1511    tab_flds(:,:,ids_tauzzv)  = tab_flds(:,:,ids_tauzzu)
1512!$OMP END MASTER
1513
1514!*************************************************************************************
1515! NetCDF output of all fields just before sending to coupler.
1516!
1517!*************************************************************************************
1518    IF (is_sequential) THEN
1519        DO j=1,maxsend
1520          IF (infosend(j)%action) CALL histwrite(nidct,infosend(j)%name, itau_w, &
1521             tab_flds(:,:,j),nbp_lon*(nbp_lat),ndexct)
1522        ENDDO
1523    ENDIF
1524!*************************************************************************************
1525! Send the table of all fields
1526!
1527!*************************************************************************************
1528    time_sec=(itime-1)*dtime
1529#ifdef CPP_COUPLE
1530!$OMP MASTER
1531    CALL intocpl(time_sec, lafin, tab_flds(:,:,:))
1532!$OMP END MASTER
1533#endif
1534
1535!*************************************************************************************
1536! Finish with some dellocate
1537!
1538!************************************************************************************* 
1539    sum_error=0
1540    DEALLOCATE(cpl_sols2D, cpl_nsol2D, cpl_rain2D, cpl_snow2D, stat=error )
1541    sum_error = sum_error + error
1542    DEALLOCATE(cpl_evap2D, cpl_tsol2D, cpl_fder2D, cpl_albe2D, stat=error )
1543    sum_error = sum_error + error
1544    DEALLOCATE(cpl_taux2D, cpl_tauy2D, cpl_windsp2D, cpl_taumod2D, stat=error )
1545    sum_error = sum_error + error
1546    DEALLOCATE(cpl_sens_rain2D, cpl_sens_snow2D, stat=error)
1547    sum_error = sum_error + error
1548
1549   
1550    IF (carbon_cycle_cpl) THEN
1551       DEALLOCATE(cpl_atm_co22D, stat=error )
1552       sum_error = sum_error + error
1553    ENDIF
1554
1555    if (activate_ocean_skin == 2) deallocate(cpl_delta_sst_2d, &
1556         cpl_delta_sal_2d, cpl_dter_2d, cpl_dser_2d, cpl_dt_ds_2d)
1557
1558    IF (sum_error /= 0) THEN
1559       abort_message='Pb in deallocation of cpl_xxxx2D coupling variables'
1560       CALL abort_physic(modname,abort_message,1)
1561    ENDIF
1562   
1563  END SUBROUTINE cpl_send_all
1564!
1565!*************************************************************************************
1566!
1567  SUBROUTINE cpl2gath(champ_in, champ_out, knon, knindex)
1568  USE mod_phys_lmdz_para
1569! Cette routine transforme un champs de la grille 2D recu du coupleur sur la grille
1570! 'gathered' (la grille physiq comprime).
1571!
1572!
1573! input:         
1574!   champ_in     champ sur la grille 2D
1575!   knon         nombre de points dans le domaine a traiter
1576!   knindex      index des points de la surface a traiter
1577!
1578! output:
1579!   champ_out    champ sur la grille 'gatherd'
1580!
1581    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
1582
1583! Input
1584    INTEGER, INTENT(IN)                       :: knon
1585    REAL, DIMENSION(nbp_lon,jj_nb), INTENT(IN)    :: champ_in
1586    INTEGER, DIMENSION(klon), INTENT(IN)      :: knindex
1587
1588! Output
1589    REAL, DIMENSION(klon_mpi), INTENT(OUT)        :: champ_out
1590
1591! Local
1592    INTEGER                                   :: i, ig
1593    REAL, DIMENSION(klon_mpi)                 :: temp_mpi
1594    REAL, DIMENSION(klon)                     :: temp_omp
1595
1596!*************************************************************************************
1597!
1598   
1599
1600! Transform from 2 dimensions (nbp_lon,jj_nb) to 1 dimension (klon)
1601!$OMP MASTER
1602    CALL Grid2Dto1D_mpi(champ_in,temp_mpi)
1603!$OMP END MASTER
1604
1605    CALL scatter_omp(temp_mpi,temp_omp)
1606   
1607! Compress from klon to knon
1608    DO i = 1, knon
1609       ig = knindex(i)
1610       champ_out(i) = temp_omp(ig)
1611    ENDDO
1612
1613  END SUBROUTINE cpl2gath
1614!
1615!*************************************************************************************
1616!
1617  SUBROUTINE gath2cpl(champ_in, champ_out, knon, knindex)
1618  USE mod_phys_lmdz_para
1619! Cette routine ecrit un champ 'gathered' sur la grille 2D pour le passer
1620! au coupleur.
1621!
1622! input:         
1623!   champ_in     champ sur la grille gathere       
1624!   knon         nombre de points dans le domaine a traiter
1625!   knindex      index des points de la surface a traiter
1626!
1627! output:
1628!   champ_out    champ sur la grille 2D
1629!
1630    USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
1631   
1632! Input arguments
1633!*************************************************************************************
1634    INTEGER, INTENT(IN)                    :: knon
1635    REAL, DIMENSION(klon), INTENT(IN)      :: champ_in
1636    INTEGER, DIMENSION(klon), INTENT(IN)   :: knindex
1637
1638! Output arguments
1639!*************************************************************************************
1640    REAL, DIMENSION(nbp_lon,jj_nb), INTENT(OUT) :: champ_out
1641
1642! Local variables
1643!*************************************************************************************
1644    INTEGER                                :: i, ig
1645    REAL, DIMENSION(klon)                  :: temp_omp
1646    REAL, DIMENSION(klon_mpi)              :: temp_mpi
1647!*************************************************************************************
1648
1649! Decompress from knon to klon
1650    temp_omp = 0.
1651    DO i = 1, knon
1652       ig = knindex(i)
1653       temp_omp(ig) = champ_in(i)
1654    ENDDO
1655
1656! Transform from 1 dimension (klon) to 2 dimensions (nbp_lon,jj_nb)
1657    CALL gather_omp(temp_omp,temp_mpi)
1658
1659!$OMP MASTER   
1660    CALL Grid1Dto2D_mpi(temp_mpi,champ_out)
1661   
1662    IF (is_north_pole_dyn) champ_out(:,1)=temp_mpi(1)
1663    IF (is_south_pole_dyn) champ_out(:,jj_nb)=temp_mpi(klon)
1664!$OMP END MASTER
1665   
1666  END SUBROUTINE gath2cpl
1667!
1668!*************************************************************************************
1669!
1670END MODULE cpl_mod
1671
Note: See TracBrowser for help on using the repository browser.