source: LMDZ6/trunk/libf/phylmd/carbon_cycle_mod.F90 @ 3420

Last change on this file since 3420 was 3391, checked in by oboucher, 6 years ago

Fields can now be passed to ORCHIDEE through surf_land_orchidee.
Fields are compressed in pbl_surface_mod and uncompressed in surf_land_orchidee.
Cosmetic changes, including in surf_land_mod.

  • 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
File size: 33.0 KB
RevLine 
[1227]1MODULE carbon_cycle_mod
[3387]2!=======================================================================
3!   Authors: Patricia Cadule and Laurent Fairhead
4!            base sur un travail anterieur mene par Patricia Cadule et Josefine Ghattas
5!
6!  Purpose and description:
7!  -----------------------
8! Control module for the carbon CO2 tracers :
[1454]9!   - Identification
10!   - Get concentrations comming from coupled model or read from file to tracers
11!   - Calculate new RCO2 for radiation scheme
12!   - Calculate new carbon flux for sending to coupled models (PISCES and ORCHIDEE)
13!
[3387]14! Module permettant de mettre a jour les champs (puits et sources) pour le
15! transport de CO2 en online (IPSL-CM et LMDZOR) et offline (lecture de carte)
16!
17! Le cas online/offline est defini par le flag carbon_cycle_cpl (y/n)
18! Le transport du traceur CO2 est defini par le flag carbon_cycle_tr (y/n)
19! la provenance des champs (termes de puits) est denini par le flag level_coupling_esm
20!
21! level_coupling_esm : level of coupling of the biogeochemical fields between
22! LMDZ, ORCHIDEE and NEMO
23! Definitions of level_coupling_esm in physiq.def
24! level_coupling_esm = 0  ! No field exchange between LMDZ and ORCHIDEE models
25!                         ! No field exchange between LMDZ and NEMO
26! level_coupling_esm = 1  ! Field exchange between LMDZ and ORCHIDEE models
27!                         ! No field exchange between LMDZ and NEMO models
28! level_coupling_esm = 2  ! No field exchange between LMDZ and ORCHIDEE models
29!                         ! Field exchange between LMDZ and NEMO models
30! level_coupling_esm = 3  ! Field exchange between LMDZ and ORCHIDEE models
31!                         ! Field exchange between LMDZ and NEMO models
32!=======================================================================
[1227]33
34  IMPLICIT NONE
35  SAVE
36  PRIVATE
[3387]37  PUBLIC :: carbon_cycle_init, carbon_cycle, infocfields_init
[1227]38
39! Variables read from parmeter file physiq.def
[1250]40  LOGICAL, PUBLIC :: carbon_cycle_tr        ! 3D transport of CO2 in the atmosphere, parameter read in conf_phys
[1249]41!$OMP THREADPRIVATE(carbon_cycle_tr)
[1250]42  LOGICAL, PUBLIC :: carbon_cycle_cpl       ! Coupling of CO2 fluxes between LMDZ/ORCHIDEE and LMDZ/OCEAN(PISCES)
[1249]43!$OMP THREADPRIVATE(carbon_cycle_cpl)
[3387]44  INTEGER, PUBLIC :: level_coupling_esm ! Level of coupling for the ESM - 0, 1, 2, 3
[3384]45!$OMP THREADPRIVATE(level_coupling_esm)
[1454]46
[1759]47  LOGICAL :: carbon_cycle_emis_comp_omp=.FALSE.
[1227]48  LOGICAL :: carbon_cycle_emis_comp=.FALSE. ! Calculation of emission compatible
[1454]49!$OMP THREADPRIVATE(carbon_cycle_emis_comp)
[1227]50
[1759]51  LOGICAL :: RCO2_inter_omp
[1454]52  LOGICAL :: RCO2_inter  ! RCO2 interactive : if true calculate new value RCO2 for the radiation scheme
53!$OMP THREADPRIVATE(RCO2_inter)
54
[3385]55! Scalare values when no transport, from physiq.def
56  REAL :: fos_fuel_s_omp
57  REAL :: fos_fuel_s  ! carbon_cycle_fos_fuel dans physiq.def
58!$OMP THREADPRIVATE(fos_fuel_s)
[1227]59  REAL :: emis_land_s ! not yet implemented
[1250]60!$OMP THREADPRIVATE(emis_land_s)
[1227]61
[1454]62  REAL :: airetot     ! Total area of the earth surface
63!$OMP THREADPRIVATE(airetot)
[1250]64
[3385]65  INTEGER :: ntr_co2  ! Number of tracers concerning the carbon cycle
[1454]66!$OMP THREADPRIVATE(ntr_co2)
67
[3385]68! fco2_ocn_day : flux CO2 from ocean for 1 day (cumulated) [gC/m2/d]. Allocation and initalization done in cpl_mod
69  REAL, DIMENSION(:), ALLOCATABLE, PUBLIC :: fco2_ocn_day
70!$OMP THREADPRIVATE(fco2_ocn_day)
[1454]71
[3385]72  REAL, DIMENSION(:), ALLOCATABLE :: fco2_land_day   ! flux CO2 from land for 1 day (cumulated)  [gC/m2/d]
73!$OMP THREADPRIVATE(fco2_land_day)
74  REAL, DIMENSION(:), ALLOCATABLE :: fco2_lu_day     ! Emission from land use change for 1 day (cumulated) [gC/m2/d]
75!$OMP THREADPRIVATE(fco2_lu_day)
[1227]76
[1454]77  REAL, DIMENSION(:,:), ALLOCATABLE :: dtr_add       ! Tracer concentration to be injected
78!$OMP THREADPRIVATE(dtr_add)
79
[3385]80! Following 2 fields will be allocated and initialized in surf_land_orchidee
81  REAL, DIMENSION(:), ALLOCATABLE, PUBLIC :: fco2_land_inst  ! flux CO2 from land at one time step
82!$OMP THREADPRIVATE(fco2_land_inst)
83  REAL, DIMENSION(:), ALLOCATABLE, PUBLIC :: fco2_lu_inst    ! Emission from land use change at one time step
84!$OMP THREADPRIVATE(fco2_lu_inst)
[1227]85
[3385]86! Calculated co2 field to be send to the ocean via the coupler and to ORCHIDEE
[1454]87  REAL, DIMENSION(:), ALLOCATABLE, PUBLIC :: co2_send ! Field allocated in phyetat0
[1250]88!$OMP THREADPRIVATE(co2_send)
[1227]89
[3387]90! nbfields : total number of fields
91  INTEGER, PUBLIC :: nbcf
92!$OMP THREADPRIVATE(nbcf)
[1454]93
[3387]94! nbcf_in : number of fields IN
95  INTEGER, PUBLIC  :: nbcf_in
96!$OMP THREADPRIVATE(nbcf_in)
97
98! nbcf_in_orc : number of fields IN
99  INTEGER, PUBLIC  :: nbcf_in_orc
100!$OMP THREADPRIVATE(nbcf_in_orc)
101
102! nbcf_in_inca : number of fields IN (from INCA)
103  INTEGER, PUBLIC  :: nbcf_in_inca
104!$OMP THREADPRIVATE(nbcf_in_inca)
105
106! nbcf_in_nemo : number of fields IN (from nemo)
107  INTEGER, PUBLIC  :: nbcf_in_nemo
108!$OMP THREADPRIVATE(nbcf_in_nemo)
109
110! nbcf_in_ant : number of fields IN (from anthropogenic sources)
111  INTEGER, PUBLIC  :: nbcf_in_ant
112!$OMP THREADPRIVATE(nbcf_in_ant)
113
114! nbcf_out : number of fields OUT
115  INTEGER, PUBLIC :: nbcf_out
116!$OMP THREADPRIVATE(nbcf_out)
117
118! Name of variables
119  CHARACTER(len=25), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfname     ! coupling field short name for restart (?) and diagnostics
120!$OMP THREADPRIVATE(cfname)
121
122  CHARACTER(len=25), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfname_in  ! coupling field short name for restart (?) and diagnostics
123!$OMP THREADPRIVATE(cfname_in)
124
125  CHARACTER(len=25), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfname_out ! coupling field short name for restart (?) and diagnostics
126!$OMP THREADPRIVATE(cfname_out)
127
128  CHARACTER(len=15), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfunits_in  !  coupling field units for diagnostics
129!$OMP THREADPRIVATE(cfunits_in)
130
131  CHARACTER(len=15), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfunits_out !  coupling field units for diagnostics
132!$OMP THREADPRIVATE(cfunits_out)
133
134  CHARACTER(len=120), ALLOCATABLE, DIMENSION(:), PUBLIC :: cftext_in  ! coupling field long name for diagnostics
135!$OMP THREADPRIVATE(cftext_in)
136
137  CHARACTER(len=120), ALLOCATABLE, DIMENSION(:), PUBLIC :: cftext_out ! coupling field long name for diagnostics
138!$OMP THREADPRIVATE(cftext_out)
139
140  CHARACTER(len=5), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfmod1 ! model 1 (rreference) : LMDz
141!$OMP THREADPRIVATE(cfmod1)
142
143  CHARACTER(len=5), ALLOCATABLE, DIMENSION(:), PUBLIC :: cfmod2 ! model 2
144!$OMP THREADPRIVATE(cfmod2)
145
[3391]146  CHARACTER(LEN=20), ALLOCATABLE, DIMENSION(:), PUBLIC :: field_out_names
147!$OMP THREADPRIVATE(field_out_names)
[3387]148
[3391]149  CHARACTER(LEN=20), ALLOCATABLE, DIMENSION(:), PUBLIC :: field_in_names
150!$OMP THREADPRIVATE(field_in_names)
[3387]151
[3391]152  REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC :: fields_in   !  klon,nbcf_in
153!$OMP THREADPRIVATE(fields_in)
154
155  REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC :: yfields_in  !  knon,nbcf_in
156!$OMP THREADPRIVATE(yfields_in)
157
158  REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC :: fields_out  !  klon,nbcf_out
159!$OMP THREADPRIVATE(fields_out)
160
161  REAL, ALLOCATABLE, DIMENSION(:,:), PUBLIC :: yfields_out !  knon,nbcf_out
162!$OMP THREADPRIVATE(yfields_out)
163
[3390]164  TYPE, PUBLIC :: co2_trac_type
[1454]165     CHARACTER(len = 8) :: name       ! Tracer name in tracer.def
166     INTEGER            :: id         ! Index in total tracer list, tr_seri
167     CHARACTER(len=30)  :: file       ! File name
[3385]168     LOGICAL            :: cpl        ! True if this tracers is coupled from ORCHIDEE or PISCES.
[1454]169                                      ! False if read from file.
170     INTEGER            :: updatefreq ! Frequence to inject in second
171     INTEGER            :: readstep   ! Actual time step to read in file
172     LOGICAL            :: updatenow  ! True if this tracer should be updated this time step
173  END TYPE co2_trac_type
174  INTEGER,PARAMETER :: maxco2trac=5  ! Maximum number of different CO2 fluxes
175  TYPE(co2_trac_type), DIMENSION(maxco2trac) :: co2trac
176
[1227]177CONTAINS
178 
[1454]179  SUBROUTINE carbon_cycle_init(tr_seri, pdtphys, aerosol, radio)
180! This subroutine is called from traclmdz_init, only at first timestep.
181! - Read controle parameters from .def input file
182! - Search for carbon tracers and set default values
183! - Allocate variables
184! - Test for compatibility
185
[1227]186    USE dimphy
[2351]187    USE geometry_mod, ONLY : cell_area
[1454]188    USE mod_phys_lmdz_transfert_para
[2320]189    USE infotrac_phy, ONLY: nbtr, nqo, niadv, tname
[1227]190    USE IOIPSL
191    USE surface_data, ONLY : ok_veget, type_ocean
[1454]192    USE phys_cal_mod, ONLY : mth_len
[2311]193    USE print_control_mod, ONLY: lunout
[1227]194
195    IMPLICIT NONE
196    INCLUDE "clesphys.h"
197 
198! Input argument
199    REAL,DIMENSION(klon,klev,nbtr),INTENT(IN) :: tr_seri ! Concentration Traceur [U/KgA] 
[1454]200    REAL,INTENT(IN)                           :: pdtphys ! length of time step in physiq (sec)
[1227]201
202! InOutput arguments
203    LOGICAL,DIMENSION(nbtr), INTENT(INOUT) :: aerosol
204    LOGICAL,DIMENSION(nbtr), INTENT(INOUT) :: radio
205
206! Local variables
[1454]207    INTEGER               :: ierr, it, iiq, itc
208    INTEGER               :: teststop
[1227]209
[1454]210! 1) Read controle parameters from .def input file
211! ------------------------------------------------
212    ! Read fosil fuel value if no transport
213    IF (.NOT. carbon_cycle_tr) THEN
[1759]214!$OMP MASTER
[3385]215       fos_fuel_s_omp = 0.
216       CALL getin ('carbon_cycle_fos_fuel',fos_fuel_s_omp)
[1759]217!$OMP END MASTER
218!$OMP BARRIER
[3385]219       fos_fuel_s=fos_fuel_s_omp
220       WRITE(lunout,*) 'carbon_cycle_fos_fuel = ', fos_fuel_s
[1227]221    END IF
222
[1454]223    ! Read parmeter for calculation compatible emission
224    IF (.NOT. carbon_cycle_tr) THEN
[1759]225!$OMP MASTER
226       carbon_cycle_emis_comp_omp=.FALSE.
227       CALL getin('carbon_cycle_emis_comp',carbon_cycle_emis_comp_omp)
228!$OMP END MASTER
229!$OMP BARRIER
230       carbon_cycle_emis_comp=carbon_cycle_emis_comp_omp
[1454]231       WRITE(lunout,*) 'carbon_cycle_emis_comp = ',carbon_cycle_emis_comp
232       IF (carbon_cycle_emis_comp) THEN
[2311]233          CALL abort_physic('carbon_cycle_init', 'carbon_cycle_emis_comp option not yet implemented!!',1)
[1454]234       END IF
235    END IF
[1227]236
[1454]237    ! Read parameter for interactive calculation of the CO2 value for the radiation scheme
[1759]238!$OMP MASTER
239    RCO2_inter_omp=.FALSE.
240    CALL getin('RCO2_inter',RCO2_inter_omp)
241!$OMP END MASTER
242!$OMP BARRIER
243    RCO2_inter=RCO2_inter_omp
[1454]244    WRITE(lunout,*) 'RCO2_inter = ', RCO2_inter
245    IF (RCO2_inter) THEN
246       WRITE(lunout,*) 'RCO2 will be recalculated once a day'
247       WRITE(lunout,*) 'RCO2 initial = ', RCO2
248    END IF
[1227]249
[1454]250
[3385]251! 2) Search for carbon tracers and set default values
[1454]252! ---------------------------------------------------
253    itc=0
[3385]254    DO it=1,nbtr
255!!       iiq=niadv(it+2)                                                            ! jyg
256       iiq=niadv(it+nqo)                                                            ! jyg
257       
258       SELECT CASE(tname(iiq))
259       CASE("fCO2_ocn")
[1454]260          itc = itc + 1
[3385]261          co2trac(itc)%name='fCO2_ocn'
262          co2trac(itc)%id=it
263          co2trac(itc)%file='fl_co2_ocean.nc'
264          IF (carbon_cycle_cpl .AND. type_ocean=='couple') THEN
265             co2trac(itc)%cpl=.TRUE.
266             co2trac(itc)%updatefreq = 86400 ! Once a day as the coupling with OASIS/PISCES
267          ELSE
268             co2trac(itc)%cpl=.FALSE.
269             co2trac(itc)%updatefreq = 86400*mth_len ! Once a month
[1227]270          END IF
[3385]271       CASE("fCO2_land")
[1454]272          itc = itc + 1
[3385]273          co2trac(itc)%name='fCO2_land'
274          co2trac(itc)%id=it
275          co2trac(itc)%file='fl_co2_land.nc'
[1454]276          IF (carbon_cycle_cpl .AND. ok_veget) THEN
[3385]277             co2trac(itc)%cpl=.TRUE.
278             co2trac(itc)%updatefreq = INT(pdtphys) ! Each timestep as the coupling with ORCHIDEE
[1454]279          ELSE
[3385]280             co2trac(itc)%cpl=.FALSE.
281!             co2trac(itc)%updatefreq = 10800   ! 10800sec = 3H
282             co2trac(itc)%updatefreq = 86400*mth_len ! Once a month
[1454]283          END IF
[3385]284       CASE("fCO2_land_use")
[1454]285          itc = itc + 1
[3385]286          co2trac(itc)%name='fCO2_land_use'
287          co2trac(itc)%id=it
288          co2trac(itc)%file='fl_co2_land_use.nc'
[1454]289          IF (carbon_cycle_cpl .AND. ok_veget) THEN
[3385]290             co2trac(it)%cpl=.TRUE.
291             co2trac(itc)%updatefreq = INT(pdtphys) ! Each timestep as the coupling with ORCHIDEE
[1454]292          ELSE
[3385]293             co2trac(itc)%cpl=.FALSE.
294             co2trac(itc)%updatefreq = 10800   ! 10800sec = 3H
[1454]295          END IF
[3385]296       CASE("fCO2_fos_fuel")
[1454]297          itc = itc + 1
[3385]298          co2trac(itc)%name='fCO2_fos_fuel'
299          co2trac(itc)%id=it
300          co2trac(itc)%file='fossil_fuel.nc'
301          co2trac(itc)%cpl=.FALSE.       ! This tracer always read from file
302!         co2trac(itc)%updatefreq = 86400  ! 86400sec = 24H Cadule case
303          co2trac(itc)%updatefreq = 86400*mth_len ! Once a month
304       CASE("fCO2_bbg")
[1454]305          itc = itc + 1
[3385]306          co2trac(itc)%name='fCO2_bbg'
307          co2trac(itc)%id=it
308          co2trac(itc)%file='fl_co2_bbg.nc'
309          co2trac(itc)%cpl=.FALSE.       ! This tracer always read from file
310          co2trac(itc)%updatefreq = 86400*mth_len ! Once a month
[1454]311       CASE("fCO2")
312          ! fCO2 : One tracer transporting the total CO2 flux
313          itc = itc + 1
[3385]314          co2trac(itc)%name='fCO2'
315          co2trac(itc)%id=it
316          co2trac(itc)%file='fl_co2.nc'
[1454]317          IF (carbon_cycle_cpl) THEN
[3385]318             co2trac(itc)%cpl=.TRUE.
[1454]319          ELSE
[3385]320             co2trac(itc)%cpl=.FALSE.
[1454]321          END IF
[3385]322          co2trac(itc)%updatefreq = 86400
[1454]323          ! DOES THIS WORK ???? Problematic due to implementation of the coupled fluxes...
[2311]324          CALL abort_physic('carbon_cycle_init','transport of total CO2 has to be implemented and tested',1)
[1454]325       END SELECT
326    END DO
[1227]327
[3385]328    ! Total number of carbon CO2 tracers
329    ntr_co2 = itc
330   
[1454]331    ! Definition of control varaiables for the tracers
332    DO it=1,ntr_co2
333       aerosol(co2trac(it)%id) = .FALSE.
334       radio(co2trac(it)%id)   = .FALSE.
335    END DO
336   
337    ! Vector indicating which timestep to read for each tracer
338    ! Always start read in the beginning of the file
339    co2trac(:)%readstep = 0
340   
[1227]341
[1454]342! 3) Allocate variables
343! ---------------------
344    ! Allocate vector for storing fluxes to inject
[3385]345    ALLOCATE(dtr_add(klon,maxco2trac), stat=ierr)
[2311]346    IF (ierr /= 0) CALL abort_physic('carbon_cycle_init', 'pb in allocation 11',1)       
[3385]347   
[1454]348    ! Allocate variables for cumulating fluxes from ORCHIDEE
[3385]349    IF (RCO2_inter) THEN
350       IF (.NOT. carbon_cycle_tr .AND. carbon_cycle_cpl) THEN
351          ALLOCATE(fco2_land_day(klon), stat=ierr)
[2311]352          IF (ierr /= 0) CALL abort_physic('carbon_cycle_init', 'pb in allocation 2',1)
[3385]353          fco2_land_day(1:klon) = 0.
[1454]354         
[3385]355          ALLOCATE(fco2_lu_day(klon), stat=ierr)
[2311]356          IF (ierr /= 0) CALL abort_physic('carbon_cycle_init', 'pb in allocation 3',1)
[3385]357          fco2_lu_day(1:klon)   = 0.
[1227]358       END IF
[3385]359    END IF
[1227]360
[1454]361
362! 4) Test for compatibility
363! -------------------------
364!    IF (carbon_cycle_cpl .AND. type_ocean/='couple') THEN
365!       WRITE(lunout,*) 'Coupling with ocean model is needed for carbon_cycle_cpl'
[2311]366!       CALL abort_physic('carbon_cycle_init', 'coupled ocean is needed for carbon_cycle_cpl',1)
[1454]367!    END IF
368!
369!    IF (carbon_cycle_cpl .AND..NOT. ok_veget) THEN
370!       WRITE(lunout,*) 'Coupling with surface land model ORCHDIEE is needed for carbon_cycle_cpl'
[2311]371!       CALL abort_physic('carbon_cycle_init', 'ok_veget is needed for carbon_cycle_cpl',1)
[1454]372!    END IF
373
374    ! Compiler test : following should never happen
375    teststop=0
376    DO it=1,teststop
[2311]377       CALL abort_physic('carbon_cycle_init', 'Entering loop from 1 to 0',1)
[1454]378    END DO
379
[3385]380    IF (ntr_co2==0) THEN
[1454]381       ! No carbon tracers found in tracer.def. It is not possible to do carbon cycle
[3385]382       WRITE(lunout,*) 'No carbon tracers found in tracer.def. Not ok with carbon_cycle_tr and/or carbon_cycle_cp'
[2311]383       CALL abort_physic('carbon_cycle_init', 'No carbon tracers found in tracer.def',1)
[1227]384    END IF
[1454]385   
386! 5) Calculate total area of the earth surface
387! --------------------------------------------
[2351]388    CALL reduce_sum(SUM(cell_area),airetot)
[1454]389    CALL bcast(airetot)
[1227]390
391  END SUBROUTINE carbon_cycle_init
392
[1454]393  SUBROUTINE carbon_cycle(nstep, pdtphys, pctsrf, tr_seri, source)
394! Subroutine for injection of co2 in the tracers
[1227]395!
[1454]396! - Find out if it is time to update
397! - Get tracer from coupled model or from file
398! - Calculate new RCO2 value for the radiation scheme
399! - Calculate CO2 flux to send to ocean and land models (PISCES and ORCHIDEE)
[1227]400
[2320]401    USE infotrac_phy, ONLY: nbtr
[1227]402    USE dimphy
[1454]403    USE mod_phys_lmdz_transfert_para
[1227]404    USE phys_cal_mod, ONLY : mth_cur, mth_len
405    USE phys_cal_mod, ONLY : day_cur
[1785]406    USE indice_sol_mod
[2311]407    USE print_control_mod, ONLY: lunout
[2351]408    USE geometry_mod, ONLY : cell_area
[1227]409
410    IMPLICIT NONE
411
412    INCLUDE "clesphys.h"
[1454]413    INCLUDE "YOMCST.h"
[1227]414
415! In/Output arguments
416    INTEGER,INTENT(IN) :: nstep      ! time step in physiq
417    REAL,INTENT(IN)    :: pdtphys    ! length of time step in physiq (sec)
[1454]418    REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf            ! Surface fraction
419    REAL, DIMENSION(klon,klev,nbtr), INTENT(INOUT)  :: tr_seri ! All tracers
420    REAL, DIMENSION(klon,nbtr), INTENT(INOUT)       :: source  ! Source for all tracers
[1227]421
422! Local variables
[3385]423    INTEGER :: it
[1227]424    LOGICAL :: newmonth ! indicates if a new month just started
425    LOGICAL :: newday   ! indicates if a new day just started
426    LOGICAL :: endday   ! indicated if last time step in a day
427
[3385]428    REAL, PARAMETER :: fact=1.E-15/2.12  ! transformation factor from gC/m2/day => ppm/m2/day
[1454]429    REAL, DIMENSION(klon) :: fco2_tmp
[1227]430    REAL :: sumtmp
431    REAL :: delta_co2_ppm
432   
433
[1454]434! 1) Calculate logicals indicating if it is a new month, new day or the last time step in a day (end day)
435! -------------------------------------------------------------------------------------------------------
436
[3385]437    newday = .FALSE.; endday = .FALSE.; newmonth = .FALSE.
[1227]438
439    IF (MOD(nstep,INT(86400./pdtphys))==1) newday=.TRUE.
440    IF (MOD(nstep,INT(86400./pdtphys))==0) endday=.TRUE.
441    IF (newday .AND. day_cur==1) newmonth=.TRUE.
[1454]442
443! 2)  For each carbon tracer find out if it is time to inject (update)
444! --------------------------------------------------------------------
[3385]445    DO it = 1, ntr_co2
446       IF ( MOD(nstep,INT(co2trac(it)%updatefreq/pdtphys)) == 1 ) THEN
447          co2trac(it)%updatenow = .TRUE.
[1454]448       ELSE
[3385]449          co2trac(it)%updatenow = .FALSE.
[1227]450       END IF
[1454]451    END DO
[1227]452
[1454]453! 3) Get tracer update
454! --------------------------------------
[3385]455    DO it = 1, ntr_co2
456       IF ( co2trac(it)%updatenow ) THEN
457          IF ( co2trac(it)%cpl ) THEN
[1454]458             ! Get tracer from coupled model
[3385]459             SELECT CASE(co2trac(it)%name)
460             CASE('fCO2_land')     ! from ORCHIDEE
461                dtr_add(:,it) = fco2_land_inst(:)*pctsrf(:,is_ter)*fact ! [ppm/m2/day]
462             CASE('fCO2_land_use') ! from ORCHIDEE
463                dtr_add(:,it) = fco2_lu_inst(:)  *pctsrf(:,is_ter)*fact ! [ppm/m2/day]
464             CASE('fCO2_ocn')      ! from PISCES
465                dtr_add(:,it) = fco2_ocn_day(:)  *pctsrf(:,is_oce)*fact ! [ppm/m2/day]
[1454]466             CASE DEFAULT
[3385]467                WRITE(lunout,*) 'Error with tracer ',co2trac(it)%name
[2311]468                CALL abort_physic('carbon_cycle', 'No coupling implemented for this tracer',1)
[1454]469             END SELECT
470          ELSE
471             ! Read tracer from file
[3385]472             co2trac(it)%readstep = co2trac(it)%readstep + 1 ! increment time step in file
473! Patricia   CALL read_map2D(co2trac(it)%file,'fco2',co2trac(it)%readstep,.FALSE.,dtr_add(:,it))
474             CALL read_map2D(co2trac(it)%file,'fco2',co2trac(it)%readstep,.TRUE.,dtr_add(:,it))
[1227]475
[3385]476             ! Converte from kgC/m2/h to kgC/m2/s
477             dtr_add(:,it) = dtr_add(:,it)/3600
[1454]478             ! Add individual treatment of values read from file
[3385]479             SELECT CASE(co2trac(it)%name)
480             CASE('fCO2_land')
[1454]481                dtr_add(:,it) = dtr_add(:,it) *pctsrf(:,is_ter)
[3385]482             CASE('fCO2_land_use')
[1454]483                dtr_add(:,it) = dtr_add(:,it) *pctsrf(:,is_ter)
[3385]484             CASE('fCO2_ocn')
[1454]485                dtr_add(:,it) = dtr_add(:,it) *pctsrf(:,is_oce)
486! Patricia :
487!             CASE('fCO2_fos_fuel')
488!                dtr_add(:,it) = dtr_add(:,it)/mth_len
[3385]489!                co2trac(it)%readstep = 0 ! Always read same value for fossil fuel(Cadule case)
[1454]490             END SELECT
[1227]491          END IF
[1454]492       END IF
493    END DO
[1227]494
[1454]495! 4) Update co2 tracers :
496!    Loop over all carbon tracers and add source
497! ------------------------------------------------------------------
498    IF (carbon_cycle_tr) THEN
499       DO it = 1, ntr_co2
[3385]500          IF (.FALSE.) THEN
501             tr_seri(1:klon,1,co2trac(it)%id) = tr_seri(1:klon,1,co2trac(it)%id) + dtr_add(1:klon,it)
502             source(1:klon,co2trac(it)%id) = 0.
503          ELSE
504             source(1:klon,co2trac(it)%id) = dtr_add(1:klon,it)
505          END IF
[1454]506       END DO
[1227]507    END IF
508
509
[1454]510! 5) Calculations for new CO2 value for the radiation scheme(instead of reading value from .def)
511! ----------------------------------------------------------------------------------------------
[3385]512    IF (RCO2_inter) THEN
513       ! Cumulate fluxes from ORCHIDEE at each timestep
514       IF (.NOT. carbon_cycle_tr .AND. carbon_cycle_cpl) THEN
515          IF (newday) THEN ! Reset cumulative variables once a day
516             fco2_land_day(1:klon) = 0.
517             fco2_lu_day(1:klon)   = 0.
518          END IF
519          fco2_land_day(1:klon) = fco2_land_day(1:klon) + fco2_land_inst(1:klon) ![gC/m2/day]
520          fco2_lu_day(1:klon)   = fco2_lu_day(1:klon)   + fco2_lu_inst(1:klon)   ![gC/m2/day]
521       END IF
[1227]522
[3385]523       ! At the end of a new day, calculate a mean scalare value of CO2
524       ! JG : Ici on utilise uniquement le traceur du premier couche du modele. Est-ce que c'est correcte ?
525       IF (endday) THEN
[1227]526
[1454]527          IF (carbon_cycle_tr) THEN
528             ! Sum all co2 tracers to get the total delta CO2 flux
529             fco2_tmp(:) = 0.
530             DO it = 1, ntr_co2
531                fco2_tmp(1:klon) = fco2_tmp(1:klon) + tr_seri(1:klon,1,co2trac(it)%id)
532             END DO
[3385]533             
534          ELSE IF (carbon_cycle_cpl) THEN ! no carbon_cycle_tr
[1454]535             ! Sum co2 fluxes comming from coupled models and parameter for fossil fuel
[3385]536             fco2_tmp(1:klon) = fos_fuel_s + ((fco2_lu_day(1:klon) + fco2_land_day(1:klon))*pctsrf(1:klon,is_ter) &
537                  + fco2_ocn_day(:)*pctsrf(:,is_oce)) * fact
[1227]538          END IF
539
[1454]540          ! Calculate a global mean value of delta CO2 flux
[3385]541          fco2_tmp(1:klon) = fco2_tmp(1:klon) * cell_area(1:klon)
[1227]542          CALL reduce_sum(SUM(fco2_tmp),sumtmp)
[1454]543          CALL bcast(sumtmp)
[1227]544          delta_co2_ppm = sumtmp/airetot
545         
[1454]546          ! Add initial value for co2_ppm and delta value
547          co2_ppm = co2_ppm0 + delta_co2_ppm
[3385]548         
[1454]549          ! Transformation of atmospheric CO2 concentration for the radiation code
550          RCO2 = co2_ppm * 1.0e-06  * 44.011/28.97
551         
[3385]552          WRITE(lunout,*) 'RCO2 is now updated! RCO2 = ', RCO2
553       END IF ! endday
[1227]554
[3385]555    END IF ! RCO2_inter
[1227]556
557
[1454]558! 6) Calculate CO2 flux to send to ocean and land models : PISCES and ORCHIDEE         
559! ----------------------------------------------------------------------------
[3385]560    IF (carbon_cycle_cpl) THEN
[1454]561
[3385]562       IF (carbon_cycle_tr) THEN
563          ! Sum all co2 tracers to get the total delta CO2 flux at first model layer
564          fco2_tmp(:) = 0.
[1454]565          DO it = 1, ntr_co2
[3385]566             fco2_tmp(1:klon) = fco2_tmp(1:klon) + tr_seri(1:klon,1,co2trac(it)%id)
[1454]567          END DO
568          co2_send(1:klon) = fco2_tmp(1:klon) + co2_ppm0
[1227]569       ELSE
[3385]570          ! Send a scalare value in 2D variable to ocean and land model (PISCES and ORCHIDEE)
[1454]571          co2_send(1:klon) = co2_ppm
[1227]572       END IF
573
[3385]574    END IF
[1227]575
576  END SUBROUTINE carbon_cycle
577 
[3387]578  SUBROUTINE infocfields_init
579
580    USE control_mod, ONLY: planet_type
581    USE phys_cal_mod, ONLY : mth_cur
582    USE mod_synchro_omp
583    USE mod_phys_lmdz_para, ONLY: is_mpi_root, is_omp_root
584    USE mod_phys_lmdz_transfert_para
585    USE mod_phys_lmdz_omp_transfert
586    USE dimphy, ONLY: klon
587
588    IMPLICIT NONE
589
590!=======================================================================
591!
592!   Authors: Patricia Cadule and Laurent Fairhead 
593!   -------
594!
595!  Purpose and description:
596!  -----------------------
597!
598! Infofields
599! this routine enables to define the field exchanges in both directions between
600! the atmospheric circulation model (LMDZ) and ORCHIDEE. In the future this
601! routing might apply to other models (e.g., NEMO, INCA, ...).
602! Therefore, currently with this routine, it is possible to define the coupling
603! fields only between LMDZ and ORCHIDEE.
604! The coupling_fields.def file enables to define the name of the exchanged
605! fields at the coupling interface.
606! field_in_names : the set of names of the exchanged fields in input to ORCHIDEE
607! (LMDZ to ORCHIDEE)
608! field_out_names : the set of names of the exchanged fields in output of
609! ORCHIDEE (ORCHIDEE to LMDZ)
610! n : the number of exchanged fields at th coupling interface
611! nb_fields_in : number of inputs fields to ORCHIDEE (LMDZ to ORCHIDEE)
612! nb_fields_out : number of ouput fields of ORCHIDEE (ORCHIDEE to LMDZ)
613!
614! The syntax for coupling_fields.def is as follows:
615! IMPORTANT: each column entry must be separated from the previous one by 3
616! spaces and only that
617! field name         coupling          model 1         model 2         long_name
618!                    direction
619!   10char  -3spaces-  3char  -3spaces- 4char -3spaces- 4char -3spaces-  30char
620!
621! n
622! FIELD1 IN LMDZ ORC
623! ....
624! FIELD(j) IN LMDZ ORC
625! FIELD(j+1) OUT LMDZ ORC
626! ...
627! FIELDn OUT LMDZ ORC   
628!
629!=======================================================================
630!   ... 22/12/2017 ....
631!-----------------------------------------------------------------------
632! Declarations
633
634  INCLUDE "clesphys.h"
635  INCLUDE "dimensions.h"
636  INCLUDE "iniprint.h"
637
638! Local variables
639
640  INTEGER :: iq,  ierr, stat, error
641
642  CHARACTER(LEN=20), ALLOCATABLE, DIMENSION(:), SAVE  :: cfname_root
643  CHARACTER(LEN=120), ALLOCATABLE, DIMENSION(:), SAVE :: cftext_root
644  CHARACTER(LEN=15), ALLOCATABLE, DIMENSION(:), SAVE  :: cfunits_root
645
646  CHARACTER(len=3), ALLOCATABLE, DIMENSION(:) :: cfintent_root
647  CHARACTER(len=5), ALLOCATABLE, DIMENSION(:) :: cfmod1_root
648  CHARACTER(len=5), ALLOCATABLE, DIMENSION(:) :: cfmod2_root
649
650  LOGICAL, ALLOCATABLE, DIMENSION(:), SAVE :: mask_in_root
651  LOGICAL, ALLOCATABLE, DIMENSION(:), SAVE :: mask_out_root
652
653  CHARACTER(len=*),parameter :: modname="infocfields"
654
655!-----------------------------------------------------------------------
656
657nbcf=0
658nbcf_in=0
659nbcf_out=0
660
661IF (planet_type=='earth') THEN
662
663    IF (is_mpi_root .AND. is_omp_root) THEN
664
665       IF (level_coupling_esm.GT.0) THEN
666
667          OPEN(200,file='coupling_fields.def',form='formatted',status='old', iostat=ierr)
668
669          IF (ierr.EQ.0) THEN
670
671             WRITE(lunout,*) trim(modname),': Open coupling_fields.def : ok'
672             READ(200,*) nbcf
673             WRITE(lunout,*) 'infocfields_mod.F90 --- nbcf=',nbcf
[3390]674             ALLOCATE(cfname_root(nbcf))
675             ALLOCATE(cfintent_root(nbcf))
676             ALLOCATE(cfmod1_root(nbcf))
677             ALLOCATE(cfmod2_root(nbcf))
678             ALLOCATE(cftext_root(nbcf))
679             ALLOCATE(cfunits_root(nbcf))
680             ALLOCATE(mask_in_root(nbcf))
681             ALLOCATE(mask_out_root(nbcf))
[3387]682
683             nbcf_in=0
684             nbcf_out=0
[3390]685
[3387]686             DO iq=1,nbcf
687                WRITE(lunout,*) 'infofields : field=',iq
688                READ(200,'(A15,3X,A3,3X,A5,3X,A5,3X,A120,3X,A15)',IOSTAT=ierr) &
689                   cfname_root(iq),cfintent_root(iq),cfmod1_root(iq),cfmod2_root(iq),cftext_root(iq),cfunits_root(iq)
690                cfname_root(iq)=TRIM(cfname_root(iq))
691                cfintent_root(iq)=TRIM(cfintent_root(iq))
692                cfmod1_root(iq)=TRIM(cfmod1_root(iq))
693                cfmod2_root(iq)=TRIM(cfmod2_root(iq))
694                cftext_root(iq)=TRIM(cftext_root(iq))
695                cfunits_root(iq)=TRIM(cfunits_root(iq))
696                WRITE(lunout,*) 'coupling field: ',cfname_root(iq), & 
697                               ', number: ',iq,', INTENT: ',cfintent_root(iq)
698                WRITE(lunout,*) 'coupling field: ',cfname_root(iq), &
699                               ', number: ',iq,', model 1 (ref): ',cfmod1_root(iq),', model 2: ',cfmod2_root(iq)
700                WRITE(lunout,*) 'coupling field: ',cfname_root(iq), &
701                               ', number: ',iq,', long name: ',cftext_root(iq),', units ',cfunits_root(iq)
[3390]702                IF (nbcf_in+nbcf_out.LT.nbcf) THEN
703                  IF (cfintent_root(iq).NE.'OUT') THEN
704                    nbcf_in=nbcf_in+1
705                    mask_in_root(iq)=.TRUE.
706                    mask_out_root(iq)=.FALSE.
707                  ELSE IF (cfintent_root(iq).EQ.'OUT') THEN
708                    nbcf_out=nbcf_out+1
709                    mask_in_root(iq)=.FALSE.
710                    mask_out_root(iq)=.TRUE.
711                  ENDIF
[3387]712                ELSE
[3390]713                  WRITE(lunout,*) 'abort_gcm --- nbcf    : ',nbcf
714                  WRITE(lunout,*) 'abort_gcm --- nbcf_in : ',nbcf_in
[3387]715                  WRITE(lunout,*) 'abort_gcm --- nbcf_out: ',nbcf_out
716                  CALL abort_gcm('infocfields_init','Problem in the definition of the coupling fields',1)
717               ENDIF
718             ENDDO !DO iq=1,nbcf
719          ELSE
720             WRITE(lunout,*) trim(modname),': infocfields_mod.F90 --- Problem in opening coupling_fields.def'
721             WRITE(lunout,*) trim(modname),': infocfields_mod.F90 --- WARNING using defaut values'
722          ENDIF ! ierr
723          CLOSE(200)
724
725       ENDIF ! level_coupling_esm
726
727    ENDIF !   (is_mpi_root .AND. is_omp_root)
728!$OMP BARRIER
729
730    CALL bcast(nbcf)
731    CALL bcast(nbcf_in)
732    CALL bcast(nbcf_out)
733
[3390]734    WRITE(lunout,*) 'infocfields_mod.F90 --- nbcf    =',nbcf
735    WRITE(lunout,*) 'infocfields_mod.F90 --- nbcf_in =',nbcf_in
[3387]736    WRITE(lunout,*) 'infocfields_mod.F90 --- nbcf_out=',nbcf_out
737
[3390]738    ALLOCATE(cfname(nbcf))
739    ALLOCATE(cfname_in(nbcf_in))
740    ALLOCATE(cftext_in(nbcf_in))
741    ALLOCATE(cfname_out(nbcf_out))
742    ALLOCATE(cftext_out(nbcf_out))
743    ALLOCATE(cfmod1(nbcf))
744    ALLOCATE(cfmod2(nbcf))
745    ALLOCATE(cfunits_in(nbcf_in))
746    ALLOCATE(cfunits_out(nbcf_out))
[3387]747       
748    IF (is_mpi_root .AND. is_omp_root) THEN
749
[3390]750        IF (nbcf.GT.0)     cfname=cfname_root
751        IF (nbcf_in.GT.0)  cfname_in=PACK(cfname_root,mask_in_root)
[3387]752        IF (nbcf_out.GT.0) cfname_out=PACK(cfname_root,mask_out_root)
[3390]753        IF (nbcf_in.GT.0)  cftext_in=PACK(cftext_root,mask_in_root)
[3387]754        IF (nbcf_out.GT.0) cftext_out=PACK(cftext_root,mask_out_root)
[3390]755        IF (nbcf.GT.0)     cfmod1=cfmod1_root
756        IF (nbcf.GT.0)     cfmod2=cfmod2_root
757        IF (nbcf_in.GT.0)  cfunits_in=PACK(cfunits_root,mask_in_root)
[3387]758        IF (nbcf_out.GT.0) cfunits_out=PACK(cfunits_root,mask_out_root)
759
760        nbcf_in_orc=0
761        nbcf_in_nemo=0
762        nbcf_in_inca=0
763        nbcf_in_ant=0
764
765        DO iq=1,nbcf
[3390]766            IF (cfmod1(iq) == "ORC")  nbcf_in_orc  = nbcf_in_orc  + 1 
767            IF (cfmod1(iq) == "NEMO") nbcf_in_nemo = nbcf_in_nemo + 1 
768            IF (cfmod1(iq) == "INCA") nbcf_in_inca = nbcf_in_inca + 1
769            IF (cfmod1(iq) == "ALL")  nbcf_in_orc  = nbcf_in_orc  + 1  ! ALL = ORC/NEMO/INCA
770            IF (cfmod1(iq) == "ALL")  nbcf_in_nemo = nbcf_in_nemo + 1  ! ALL = ORC/NEMO/INCA
771            IF (cfmod1(iq) == "ALL")  nbcf_in_inca = nbcf_in_inca + 1  ! ALL = ORC/NEMO/INCA
772            IF (cfmod1(iq) == "ANT")  nbcf_in_ant  = nbcf_in_ant  + 1 
[3387]773        ENDDO
774
775    ENDIF !   (is_mpi_root .AND. is_omp_root)
776!$OMP BARRIER
777
778    CALL bcast(nbcf_in_orc)
779    CALL bcast(nbcf_in_nemo)
780    CALL bcast(nbcf_in_inca)
781    CALL bcast(nbcf_in_ant)
782
[3390]783    WRITE(lunout,*) 'nbcf_in_orc  =',nbcf_in_orc
784    WRITE(lunout,*) 'nbcf_in_nemo =',nbcf_in_nemo
785    WRITE(lunout,*) 'nbcf_in_inca =',nbcf_in_inca
786    WRITE(lunout,*) 'nbcf_in_ant  =',nbcf_in_ant
[3387]787
788    IF (nbcf_in.GT.0) THEN
789        DO iq=1,nbcf_in
790          CALL bcast(cfname_in(iq))
791          CALL bcast(cftext_in(iq))
792          CALL bcast(cfunits_in(iq))
793        ENDDO
794    ENDIF
795
796    IF (nbcf_out.GT.0) THEN
797        DO iq=1,nbcf_out
798          CALL bcast(cfname_out(iq))
799          CALL bcast(cftext_out(iq))
800          CALL bcast(cfunits_out(iq))
801        ENDDO
802    ENDIF
803
804    IF (nbcf.GT.0) THEN
805        DO iq=1,nbcf
806          CALL bcast(cfmod1(iq))
807          CALL bcast(cfmod2(iq))
808        ENDDO
809    ENDIF
810
[3390]811    IF (nbcf_in.GT.0)  WRITE(lunout,*)'infocfields_mod --- cfname_in: ',cfname_in
[3387]812    IF (nbcf_out.GT.0) WRITE(lunout,*)'infocfields_mod --- cfname_out: ',cfname_out
813
[3390]814    IF (nbcf_in.GT.0)  WRITE(lunout,*)'infocfields_mod --- cftext_in: ',cftext_in
[3387]815    IF (nbcf_out.GT.0) WRITE(lunout,*)'infocfields_mod --- cftext_out: ',cftext_out
816
817    IF (nbcf.GT.0) WRITE(lunout,*)'infocfields_mod --- cfmod1: ',cfmod1
818    IF (nbcf.GT.0) WRITE(lunout,*)'infocfields_mod --- cfmod2: ',cfmod2
819
[3390]820    IF (nbcf_in.GT.0)  WRITE(lunout,*)'infocfunits_mod --- cfunits_in: ',cfunits_in
[3387]821    IF (nbcf_out.GT.0) WRITE(lunout,*)'infocfunits_mod --- cfunits_out: ',cfunits_out
822
[3390]823    IF (nbcf_in.GT.0)  WRITE(*,*)'infocfields_init --- number of fields in to LMDZ: ',nbcf_in
[3387]824    IF (nbcf_out.GT.0) WRITE(*,*)'infocfields_init --- number of fields out of LMDZ: ',nbcf_out
825
826 ELSE
827 ! Default values for other planets
828    nbcf=0
829    nbcf_in=0
830    nbcf_out=0
831 ENDIF ! planet_type
832
[3391]833 ALLOCATE(fields_in(klon,nbcf_in),stat=error)
834 IF (error /= 0)  CALL abort_gcm(modname,'Pb in allocation fields_in',1)
835 ALLOCATE(yfields_in(klon,nbcf_in),stat=error)
836 IF (error /= 0)  CALL abort_gcm(modname,'Pb in allocation yfields_in',1)
837 ALLOCATE(fields_out(klon,nbcf_out),stat=error)
838 IF (error /= 0)  CALL abort_gcm(modname,'Pb in allocation fields_out',1)
839 ALLOCATE(yfields_out(klon,nbcf_out),stat=error)
840 IF (error /= 0)  CALL abort_gcm(modname,'Pb in allocation yfields_out',1)
[3387]841
842END SUBROUTINE infocfields_init
843
[1227]844END MODULE carbon_cycle_mod
Note: See TracBrowser for help on using the repository browser.