source: LMDZ4/trunk/libf/phylmd/cpl_mod.F90 @ 1351

Last change on this file since 1351 was 1279, checked in by Laurent Fairhead, 15 years ago

Merged LMDZ4-dev branch changes r1241:1278 into the trunk
Running trunk and LMDZ4-dev in LMDZOR configuration on local
machine (sequential) and SX8 (4-proc) yields identical results
(restart and restartphy are identical binarily)
Log history from r1241 to r1278 is available by switching to
source:LMDZ4/branches/LMDZ4-dev-20091210

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