Changeset 3857 for LMDZ6/trunk/libf/phylmd/tracco2i_mod.F90

Timestamp:

Feb 26, 2021, 8:37:59 PM (3 years ago)

Author:

oboucher

Message:

Additions for the interactive carbon cycle

File:

• LMDZ6/trunk/libf/phylmd/tracco2i_mod.F90 (modified) (8 diffs)

LMDZ6/trunk/libf/phylmd/tracco2i_mod.F90

@@ -34 +34 @@
     USE carbon_cycle_mod, ONLY: id_CO2, nbcf_in, fields_in, cfname_in
     USE carbon_cycle_mod, ONLY: fco2_ocn_day, fco2_ff, fco2_bb, fco2_land, fco2_ocean
+    USE carbon_cycle_mod, ONLY: read_fco2_ocean_cor,var_fco2_ocean_cor,fco2_ocean_cor 
+    USE carbon_cycle_mod, ONLY: read_fco2_land_cor,var_fco2_land_cor,fco2_land_cor 
+    USE carbon_cycle_mod, ONLY: co2_send
     USE carbon_cycle_mod, ONLY: fco2_land_nbp, fco2_land_nep, fco2_land_fLuc
     USE carbon_cycle_mod, ONLY: fco2_land_fwoodharvest, fco2_land_fHarvest
     USE carbon_cycle_mod, ONLY: carbon_cycle_cpl, carbon_cycle_tr, carbon_cycle_rad, RCO2_glo, RCO2_tot
+    USE carbon_cycle_mod, ONLY: ocean_area_tot
+    USE carbon_cycle_mod, ONLY: land_area_tot
     USE mod_grid_phy_lmdz
     USE mod_phys_lmdz_mpi_data, ONLY: is_mpi_root
     USE mod_phys_lmdz_para, ONLY: gather, bcast, scatter
+    USE mod_phys_lmdz_omp_data, ONLY: is_omp_root
     USE phys_cal_mod
     USE phys_state_var_mod, ONLY: pctsrf
@@ -75 +81 @@
     REAL, DIMENSION(klon_glo,klev) :: co2_glo   ! variable temporaire sur la grille global
     REAL, DIMENSION(klon_glo,klev) :: m_air_glo ! variable temporaire sur la grille global
+    REAL, DIMENSION(klon_glo,nbsrf):: pctsrf_glo      !--fractions de maille sur la grille globale
+    REAL, DIMENSION(klon_glo)      :: pctsrf_ter_glo
+    REAL, DIMENSION(klon_glo)      :: pctsrf_oce_glo
+    REAL, DIMENSION(klon_glo)      :: pctsrf_sic_glo
+    REAL, DIMENSION(klon_glo)      :: cell_area_glo   !--aire des mailles sur la grille globale
 
     LOGICAL, SAVE :: check_fCO2_nbp_in_cfname
@@ -80 +91 @@
     INTEGER, SAVE :: day_pre=-1
 !$OMP THREADPRIVATE(day_pre)
+
+    REAL, PARAMETER :: secinday=86400.
 
     IF (is_mpi_root) THEN
@@ -100 +113 @@
         IF (cfname_in(nb)=="fCO2_nbp") check_fCO2_nbp_in_cfname=.TRUE.
       ENDDO
+
+      CALL gather(pctsrf,pctsrf_glo)
+      CALL gather(pctsrf(:,is_ter),pctsrf_ter_glo)
+      CALL gather(pctsrf(:,is_oce),pctsrf_oce_glo)
+      CALL gather(pctsrf(:,is_sic),pctsrf_sic_glo)
+      CALL gather(cell_area(:),cell_area_glo)
 
     ENDIF
@@ -146 +165 @@
     ENDDO
 
+    PRINT *, 'tracco2i_mod.F90 --- read_fco2_ocean_cor ',read_fco2_ocean_cor
+    PRINT *, 'tracco2i_mod.F90 --- read_fco2_land_cor ',read_fco2_land_cor
+
+IF (debutphy) THEN
+
+    IF (read_fco2_ocean_cor) THEN
+!$OMP MASTER
+       IF (is_mpi_root .AND. is_omp_root) THEN
+          ocean_area_tot=0.
+          PRINT *, 'tracco2i_mod.F90 --- var_fco2_ocean_cor (PgC/yr) ',var_fco2_ocean_cor
+          DO i=1, klon_glo
+             ocean_area_tot = ocean_area_tot + (pctsrf_oce_glo(i)+pctsrf_sic_glo(i))*cell_area_glo(i)
+          ENDDO
+      ENDIF !--is_mpi_root and is_omp_root 
+!$OMP END MASTER
+      CALL bcast(ocean_area_tot)
+     PRINT *, 'tracco2i_mod.F90 --- ocean_area_tot (debutphy) ',ocean_area_tot
+    ENDIF
+
+    IF (read_fco2_land_cor) THEN
+!$OMP MASTER
+       IF (is_mpi_root .AND. is_omp_root) THEN
+          land_area_tot=0.
+          PRINT *, 'tracco2i_mod.F90 --- var_fco2_land_cor (PgC/yr) ',var_fco2_land_cor
+          DO i=1, klon_glo
+             land_area_tot = land_area_tot + pctsrf_ter_glo(i)*cell_area_glo(i) 
+          ENDDO
+      ENDIF !--is_mpi_root and is_omp_root 
+!$OMP END MASTER
+      CALL bcast(land_area_tot)
+     PRINT *, 'tracco2i_mod.F90 --- land_area_tot (debutphy) ',land_area_tot 
+ENDIF
+
+    ENDIF !-- debutphy  
+
+    PRINT *, 'tracco2i_mod.F90 --- ocean_area_tot (m2) ',ocean_area_tot
+    PRINT *, 'tracco2i_mod.F90 --- land_area_tot (m2) ',land_area_tot
+
+    IF (read_fco2_ocean_cor) THEN
+! var_fco2_ocean_cor: correction of the surface downward CO2 flux into the ocean fgco2 (PgC/yr)
+! This is the correction of the the net air to ocean carbon flux. Positive flux is into the ocean.
+!    PRINT *, 'tracco2i_mod.F90 --- var_fco2_ocean_cor (PgC/yr) ',var_fco2_ocean_cor
+
+!var_fco2_ocean_cor: correction of the net air to ocean carbon flux (input data is a scalar in PgC/yr and must be converted in kg CO2 m-2 s-1)
+
+! Factors for carbon and carbon dioxide
+! 1 mole CO2 = 44.009 g CO2 = 12.011 g C
+! 1 ppm by volume of atmosphere CO2 = 2.13 Gt C
+! 1 gC = 44.009/12.011 gCO2
+
+! ocean_area_tot: ocean area (m2)
+
+! year_len: year length (in days)
+
+! conversion: PgC/yr --> kg CO2 m-2 s-1
+! fco2_ocean_cor  / (86400.*year_len): PgC/yr to PgC/s
+! fco2_ocean_cor  / (86400.*year_len)*(pctsrf(i,is_oce)+pctsrf(i,is_sic))/ocean_area_tot: PgC/s to PgC/s/m2
+! (fco2_ocean_cor / (86400.*year_len)*(pctsrf(i,is_oce)+pctsrf(i,is_sic))/ocean_area_tot) *1e12: PgC/s/m2 to kgC/s/m2
+! (fco2_ocean_cor / (86400.*year_len)*(pctsrf(i,is_oce)+pctsrf(i,is_sic))/ocean_area_tot) * 1e12 * (RMCO2/RMC): kgC/s/m2 to kgCO2/s/m2
+
+      DO i=1, klon  
+          fco2_ocean_cor(i)=(var_fco2_ocean_cor*(RMCO2/RMC)*(pctsrf(i,is_oce)+pctsrf(i,is_sic))/ocean_area_tot/(secinday*year_len))*1.e12
+      ENDDO
+
+      PRINT *, 'tracco2i_mod.F90 --- MINVAL(fco2_ocean_cor) ',MINVAL(fco2_ocean_cor)
+      PRINT *, 'tracco2i_mod.F90 --- MAXVAL(fco2_ocean_cor) ',MAXVAL(fco2_ocean_cor)
+
+    ELSE
+    fco2_ocean_cor(:)=0.
+    ENDIF
+
+    IF (read_fco2_land_cor) THEN
+! var_fco2_land_cor: correction of the carbon Mass Flux out of Atmosphere Due to Net Biospheric Production on Land  (PgC/yr)
+! This is the correction of the net mass flux of carbon between land and atmosphere calculated as
+! photosynthesis MINUS the sum of plant and soil respiration, carbon fluxes from
+! fire, harvest, grazing and land use change. Positive flux is into the land.
+!    PRINT *, 'tracco2i_mod.F90 --- var_fco2_land_cor (m2) ',var_fco2_land_cor
+
+!var_fco2_land_cor: correction of the et air to land carbon flux (input data is a scalar in PgC/yr and must be converted in kg CO2 m-2 s-1)
+
+! Factors for carbon and carbon dioxide
+! 1 mole CO2 = 44.009 g CO2 = 12.011 g C
+! 1 ppm by volume of atmosphere CO2 = 2.13 Gt C
+! 1 gC = 44.009/12.011 gCO2
+
+! land_area_tot: land area (m2)
+
+! year_len: year length (in days)
+
+! conversion: PgC/yr --> kg CO2 m-2 s-1
+! fco2_land_cor  / (86400.*year_len): PgC/yr to PgC/s
+! fco2_land_cor  / (86400.*year_len)*pctsrf(i,is_ter)/land_area_tot: PgC/s to PgC/s/m2
+! (fco2_land_cor / (86400.*year_len)*pctsrf(i,is_ter)/land_area_tot) *1e12: PgC/s/m2 to kgC/s/m2
+! (fco2_land_cor / (86400.*year_len)*pctsrf(i,is_ter)/land_area_tot) * 1e12 * (RMCO2/RMC): kgC/s/m2 to kgCO2/s/m2
+
+      DO i=1, klon
+         fco2_land_cor(i)=var_fco2_land_cor*RMCO2/RMC*pctsrf(i,is_ter)/land_area_tot/(secinday*year_len)*1.e12
+      ENDDO
+
+      PRINT *, 'tracco2i_mod.F90 --- MINVAL(fco2_land_cor) ',MINVAL(fco2_land_cor)
+      PRINT *, 'tracco2i_mod.F90 --- MAXVAL(fco2_land_cor) ',MAXVAL(fco2_land_cor)
+
+    ELSE
+      fco2_land_cor(:)=0.
+    ENDIF
+
 !--if fCO2_nbp is transferred we use it, otherwise we use the sum of what has been passed from ORCHIDEE
     IF (check_fCO2_nbp_in_cfname)  THEN
@@ -167 +292 @@
 !
 !--build final source term for CO2
-    source(:,id_CO2)=fco2_ff(:)+fco2_bb(:)+fco2_land(:)+fco2_ocean(:)
+    source(:,id_CO2)=fco2_ff(:)+fco2_bb(:)+fco2_land(:)+fco2_ocean(:)-fco2_ocean_cor(:)-fco2_land_cor(:)
 
 !--computing global mean CO2 for radiation
@@ -195 +320 @@
     ENDIF
 
+    PRINT *, 'tracco2i_mod.F90 --- MINVAL(tr_seri(:,1,id_CO2)*1.e6*RMD/RMCO2): L1: ',MINVAL(tr_seri(:,1,id_CO2)*1.e6*RMD/RMCO2)
+    PRINT *, 'tracco2i_mod.F90 --- MAXVAL(tr_seri(:,1,id_CO2)*1.e6*RMD/RMCO2): L1: ',MAXVAL(tr_seri(:,1,id_CO2)*1.e6*RMD/RMCO2)
+
+    PRINT *, 'tracco2i_mod.F90 --- MINVAL(tr_seri(:,79,id_CO2)*1.e6*RMD/RMCO2): L79: ',MINVAL(tr_seri(:,79,id_CO2)*1.e6*RMD/RMCO2)
+    PRINT *, 'tracco2i_mod.F90 --- MAXVAL(tr_seri(:,79,id_CO2)*1.e6*RMD/RMCO2): L79: ',MAXVAL(tr_seri(:,79,id_CO2)*1.e6*RMD/RMCO2)
+
+    co2_send(:) = tr_seri(:,1,id_CO2)*1.e6*RMD/RMCO2
+
+    PRINT *, 'tracco2i_mod.F90 --- MINVAL(co2_send) ',MINVAL(co2_send)
+    PRINT *, 'tracco2i_mod.F90 --- MAXVAL(co2_send) ',MAXVAL(co2_send)
+
   END SUBROUTINE tracco2i
 
@@ -252 +388 @@
       IF (readco2ff) THEN 
 
-        ! ... Open the COZff file
+        ! ... Open the CO2ff file
         CALL nf95_open("sflx_lmdz_co2_ff.nc", nf90_nowrite, ncid_in)