source: LMDZ5/trunk/libf/phylmd/cpl_mod.F90 @ 2539

Last change on this file since 2539 was 2538, checked in by Laurent Fairhead, 8 years ago

Computation of heat fluxes associated with solid and liquid precipitations
over ocean and seaice. Quantities are sent to the coupler
LF

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