Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/field_def_bgc.xml (revision 9351)
@@ -86,40 +86,53 @@
PHN * e3t
+
PHD * e3t
+
ZMI * e3t
+
ZME * e3t
+
DIN * e3t
+
SIL * e3t
+
FER * e3t
+
DET * e3t
+
PDS * e3t
+
DTC * e3t
+
DiC * e3t
+
ALK * e3t
+
OXY * e3t
+
@@ -1646,8 +1659,6 @@
-
-
@@ -1676,6 +1687,23 @@
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
@@ -1730,5 +1758,5 @@
-
+
@@ -1781,33 +1809,18 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
@@ -1909,33 +1922,18 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
@@ -1956,33 +1954,18 @@
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
@@ -2003,13 +1986,11 @@
-
+
-
-
@@ -2019,18 +2000,15 @@
-
-
+
-
-
@@ -2040,18 +2018,15 @@
-
-
+
-
-
@@ -2061,18 +2036,15 @@
-
-
+
-
-
@@ -2082,18 +2054,15 @@
-
-
+
-
-
@@ -2103,18 +2072,15 @@
-
-
+
-
-
@@ -2124,5 +2090,4 @@
-
@@ -2161,49 +2126,73 @@
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-
-
+
-
-
-
+
-
-
-
+
-
-
-
+
-
-
-
+
-
-
@@ -2421,4 +2410,17 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
@@ -2432,15 +2434,5 @@
-
-
-
-
-
-
-
-
-
-
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/xco2.atm
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/xco2.atm (revision 9351)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/CONFIG/SHARED/xco2.atm (revision 9351)
@@ -0,0 +1,166 @@
+% Data source is Meinshausen et al., GMD, doi:10.5194/gmd-10-2057-2017, 2017
+%
+% Year xCO2
+% (ppm)
+1849.5 284.3170000
+1850.5 284.3169999
+1851.5 284.4509998
+1852.5 284.5979998
+1853.5 284.7309998
+1854.5 284.8459998
+1855.5 284.9409998
+1856.5 285.0489998
+1857.5 285.2039998
+1858.5 285.3689998
+1859.5 285.5449998
+1860.5 285.7389997
+1861.5 285.9329997
+1862.5 286.0999997
+1863.5 286.2709997
+1864.5 286.4419997
+1865.5 286.6139996
+1866.5 286.7809996
+1867.5 286.9549996
+1868.5 287.1049996
+1869.5 287.2249996
+1870.5 287.3549996
+1871.5 287.4939995
+1872.5 287.6639995
+1873.5 287.8599995
+1874.5 288.0609995
+1875.5 288.2909994
+1876.5 288.5199994
+1877.5 288.7519994
+1878.5 288.9929994
+1879.5 289.2209994
+1880.5 289.4699993
+1881.5 289.7369993
+1882.5 290.0189993
+1883.5 290.2629992
+1884.5 290.5119992
+1885.5 290.7969992
+1886.5 291.0999992
+1887.5 291.4139991
+1888.5 291.7629991
+1889.5 292.1129991
+1890.5 292.4579990
+1891.5 292.8159990
+1892.5 293.1669989
+1893.5 293.4769989
+1894.5 293.7909989
+1895.5 294.0789989
+1896.5 294.3649988
+1897.5 294.6459988
+1898.5 294.9539987
+1899.5 295.2999986
+1900.5 295.6749984
+1901.5 296.0069980
+1902.5 296.3249974
+1903.5 296.6539968
+1904.5 296.9539963
+1905.5 297.2889957
+1906.5 297.6619951
+1907.5 298.0979945
+1908.5 298.5179940
+1909.5 298.9359935
+1910.5 299.3769929
+1911.5 299.8289923
+1912.5 300.3529917
+1913.5 300.9099911
+1914.5 301.4189908
+1915.5 301.9369903
+1916.5 302.4849897
+1917.5 303.0109891
+1918.5 303.4489886
+1919.5 303.8139883
+1920.5 304.2459877
+1921.5 304.5999873
+1922.5 304.9449868
+1923.5 305.2709861
+1924.5 305.6299855
+1925.5 305.8129850
+1926.5 305.9539845
+1927.5 306.1769839
+1928.5 306.3289833
+1929.5 306.4949827
+1930.5 306.6199824
+1931.5 306.8219821
+1932.5 307.0929818
+1933.5 307.4019813
+1934.5 307.7849806
+1935.5 308.2269800
+1936.5 309.0119793
+1937.5 309.7639786
+1938.5 310.2939782
+1939.5 310.8509776
+1940.5 311.3569769
+1941.5 311.8109763
+1942.5 312.1719757
+1943.5 312.3899752
+1944.5 312.4129748
+1945.5 312.3849746
+1946.5 312.3899740
+1947.5 312.4859732
+1948.5 312.5209725
+1949.5 312.6319721
+1950.5 312.8209713
+1951.5 313.0139706
+1952.5 313.3419700
+1953.5 313.7299694
+1954.5 314.0949688
+1955.5 314.4149680
+1956.5 314.6979672
+1957.5 314.9919665
+1958.5 315.3449659
+1959.5 315.8069651
+1960.5 316.6249644
+1961.5 317.2989638
+1962.5 318.0439631
+1963.5 318.6499622
+1964.5 319.3329614
+1965.5 319.8159605
+1966.5 320.8799597
+1967.5 321.4799589
+1968.5 322.3889580
+1969.5 323.2509570
+1970.5 324.7829559
+1971.5 325.3999550
+1972.5 327.3489541
+1973.5 329.9089531
+1974.5 330.7559525
+1975.5 330.8269520
+1976.5 331.5449510
+1977.5 333.3529501
+1978.5 335.0099494
+1979.5 336.6049484
+1980.5 338.7049479
+1981.5 340.0589478
+1982.5 340.6439474
+1983.5 342.2659473
+1984.5 344.0079495
+1985.5 345.4589538
+1986.5 346.9029481
+1987.5 348.7749474
+1988.5 351.2759329
+1989.5 352.8939214
+1990.5 354.0729301
+1991.5 355.3529349
+1992.5 356.2289489
+1993.5 356.9249589
+1994.5 358.2539622
+1995.5 360.2389564
+1996.5 362.0049615
+1997.5 363.2519618
+1998.5 365.9329580
+1999.5 367.8449748
+2000.5 369.1249800
+2001.5 370.6729901
+2002.5 372.8349943
+2003.5 375.4109950
+2004.5 376.9870025
+2005.5 378.9070050
+2006.5 381.0100070
+2007.5 382.6030247
+2008.5 384.7390186
+2009.5 386.2800190
+2010.5 388.7170293
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90 (revision 9351)
@@ -153,14 +153,11 @@
IF (tmask(ji,jj,1) == 1) then
IF (lk_oasis) THEN
- !! use 2D atm xCO2 from atm coupling
- f_xco2a(ji,jj) = PCO2a_in_cpl(ji,jj)
- !!!
!!! Jpalm test on atm xCO2
- IF ( (f_xco2a(ji,jj) > 1500.0 ).OR.(f_xco2a(ji,jj) < 100.0 ) ) THEN
+ IF ( (f_xco2a(ji,jj) .GT. 10000.0 ).OR.(f_xco2a(ji,jj) .LE. 0.0 ) ) THEN
IF(lwp) THEN
WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj), &
' -- ji =', mig(ji),' jj = ', mjg(jj)
ENDIF
- CALL ctl_stop( 'MEDUSA - Air-Sea :', 'unrealistic atm xCO2 ' )
+ CALL ctl_stop( 'MEDUSA - Air-Sea :', 'unrealistic coupled atm xCO2 ' )
ENDIF
ENDIF
@@ -261,28 +258,36 @@
!! must be associated to air-sea flux or air xCO2...
!! Check MOCSY outputs
- IF ( (f_pco2w(ji,jj) > 1.E4 ).OR.(f_pco2w(ji,jj) < 0.0 ) .OR. &
- (f_fco2w(ji,jj) > 1.E4 ).OR.(f_fco2w(ji,jj) < 0.0 ) .OR. &
- (f_fco2atm(ji,jj) > 1.E4 ).OR.(f_fco2atm(ji,jj) < 0.0 ) .OR. &
- (f_co2flux(ji,jj) > 1.E-2 ).OR.(f_co2flux(ji,jj) < -1.E-2 ) .OR. &
- (f_dpco2(ji,jj) > 1.E4 ).OR.(f_dpco2(ji,jj) < -1.E4 ) ) THEN
- IF(lwp) THEN
- WRITE(numout,*) ' surface T = ',ztmp(ji,jj)
- WRITE(numout,*) ' surface S = ',zsal(ji,jj)
- WRITE(numout,*) ' surface ALK = ',zalk(ji,jj)
- WRITE(numout,*) ' surface DIC = ',zdic(ji,jj)
- WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj)
- WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)
- WRITE(numout,*) ' surface pco2w = ',f_pco2w(ji,jj)
- WRITE(numout,*) ' surface fco2w = ',f_fco2w(ji,jj)
- WRITE(numout,*) ' surface fco2a = ',f_fco2atm(ji,jj)
- WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj)
- WRITE(numout,*) ' surface dpco2 = ',f_dpco2(ji,jj)
- WRITE(numout,*) ' MOCSY output: ji =', mig(ji),' jj = ', mjg(jj), &
- ' kt = ', kt
- WRITE(numout,*) 'MEDUSA - Air-Sea OUTPUT: unrealistic surface Carb. Chemistry'
- ENDIF
- CALL ctl_stop( 'MEDUSA - Air-Sea OUTPUT: ', &
- 'unrealistic surface Carb. Chemistry -- OUTPUTS' )
- ENDIF
+ !!===================
+ !! Jpalm -- 19-02-2018 -- remove the cap - only check MOCSY inputs
+ !! because of specific area in arabic sea where strangely
+ !! with core 2 forcing, ALK is lower than DIC and result in
+ !! Enormous dpco2 - even if all carb chem caract are OK.
+ !! and this check stops the model.
+ !! ==> let's run the model without it and see how it goes.
+ !!===================
+ !!IF ( (f_pco2w(ji,jj) > 1.E4 ).OR.(f_pco2w(ji,jj) < 0.0 ) .OR. &
+ !! (f_fco2w(ji,jj) > 1.E4 ).OR.(f_fco2w(ji,jj) < 0.0 ) .OR. &
+ !! (f_fco2atm(ji,jj) > 1.E4 ).OR.(f_fco2atm(ji,jj) < 0.0 ) .OR. &
+ !! (f_co2flux(ji,jj) > 1.E-1 ).OR.(f_co2flux(ji,jj) < -1.E-1 ) .OR. &
+ !! (f_dpco2(ji,jj) > 1.E4 ).OR.(f_dpco2(ji,jj) < -1.E4 ) ) THEN
+ !! IF(lwp) THEN
+ !! WRITE(numout,*) ' surface T = ',ztmp(ji,jj)
+ !! WRITE(numout,*) ' surface S = ',zsal(ji,jj)
+ !! WRITE(numout,*) ' surface ALK = ',zalk(ji,jj)
+ !! WRITE(numout,*) ' surface DIC = ',zdic(ji,jj)
+ !! WRITE(numout,*) ' KW660 = ',f_kw660(ji,jj)
+ !! WRITE(numout,*) ' atm xCO2 = ',f_xco2a(ji,jj)
+ !! WRITE(numout,*) ' surface pco2w = ',f_pco2w(ji,jj)
+ !! WRITE(numout,*) ' surface fco2w = ',f_fco2w(ji,jj)
+ !! WRITE(numout,*) ' surface fco2a = ',f_fco2atm(ji,jj)
+ !! WRITE(numout,*) ' surface co2flx = ',f_co2flux(ji,jj)
+ !! WRITE(numout,*) ' surface dpco2 = ',f_dpco2(ji,jj)
+ !! WRITE(numout,*) ' MOCSY output: ji =', mig(ji),' jj = ', mjg(jj), &
+ !! ' kt = ', kt
+ !! WRITE(numout,*) 'MEDUSA - Air-Sea OUTPUT: unrealistic surface Carb. Chemistry'
+ !! ENDIF
+ !! CALL ctl_stop( 'MEDUSA - Air-Sea OUTPUT: ', &
+ !! 'unrealistic surface Carb. Chemistry -- OUTPUTS' )
+ !!ENDIF
ENDIF
ENDDO
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/carb_chem.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/carb_chem.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/carb_chem.F90 (revision 9351)
@@ -103,8 +103,4 @@
!! OPEN wet point IF..THEN loop
IF (tmask(ji,jj,jk).eq.1) THEN
- IF (lk_oasis) THEN
- !! use 2D atm xCO2 from atm coupling
- f_xco2a(ji,jj) = PCO2a_in_cpl(ji,jj)
- ENDIF
!! Do carbonate chemistry
!!
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/sms_medusa.F90 (revision 9351)
@@ -283,112 +283,20 @@
#if defined key_roam
!!----------------------------------------------------------------------
-!! Atmospheric pCO2 data (1859 to 2100 inclusive)
-!!----------------------------------------------------------------------
-!!
- REAL(wp), DIMENSION(242) :: hist_pco2 !: pCO2
-
-# if defined key_rcp26
- !! UKMO, run AJKKH + KAAEC, RCP 2.6, pCO2 time evolution
- DATA hist_pco2 / 286.0230, 286.1730, 286.3230, 286.4480, 286.5730, &
- & 286.7230, 286.8480, 286.9480, 287.0480, 287.1730, &
- & 287.3230, 287.4730, 287.6480, 287.8480, 288.0730, &
- & 288.3480, 288.6480, 288.9730, 289.3470, 289.7470, &
- & 290.1730, 290.6470, 291.1470, 291.6220, 292.0720, &
- & 292.5220, 292.9220, 293.2470, 293.5220, 293.7470, &
- & 293.9470, 294.1220, 294.2720, 294.4220, 294.5470, &
- & 294.6470, 294.7470, 294.8470, 294.9710, 295.1710, &
- & 295.4460, 295.7470, 296.0720, 296.4210, 296.7710, &
- & 297.1460, 297.5710, 298.0210, 298.4460, 298.8460, &
- & 299.2460, 299.6450, 300.0210, 300.3710, 300.7200, &
- & 301.0450, 301.3460, 301.6710, 302.0200, 302.3450, &
- & 302.6450, 302.9700, 303.3450, 303.7200, 304.0700, &
- & 304.4700, 304.9200, 305.3440, 305.7700, 306.2450, &
- & 306.7190, 307.1700, 307.6440, 308.1190, 308.5440, &
- & 308.9440, 309.3440, 309.6940, 309.9440, 310.1190, &
- & 310.2440, 310.3190, 310.3190, 310.2440, 310.1440, &
- & 310.0690, 310.0440, 310.0690, 310.1440, 310.2690, &
- & 310.4440, 310.6940, 311.0430, 311.4440, 311.8690, &
- & 312.3680, 312.9430, 313.5430, 314.1680, 314.7900, &
- & 315.4430, 316.2150, 317.0170, 317.7370, 318.3400, &
- & 318.8680, 319.5900, 320.5890, 321.5470, 322.5770, &
- & 323.8440, 324.9260, 325.7960, 327.0810, 328.6180, &
- & 329.6830, 330.5250, 331.6880, 333.2120, 334.7870, &
- & 336.4640, 338.2990, 339.6660, 340.7310, 342.1360, &
- & 343.7200, 345.2200, 346.7350, 348.5820, 350.6740, &
- & 352.4230, 353.7910, 354.9530, 355.8210, 356.7130, &
- & 358.0630, 359.7720, 361.3970, 363.0900, 365.2560, &
- & 367.2810, 368.7980, 370.4000, 372.4550, 374.6920, &
- & 376.7440, 378.7440, 380.7580, 382.7080, 384.7300, &
- & 386.9310, 389.2150, 391.4910, 393.7710, 396.0460, &
- & 398.3240, 400.6080, 402.8950, 405.1780, 407.4550, &
- & 409.7260, 411.9930, 414.2500, 416.4410, 418.5280, &
- & 420.5250, 422.4390, 424.2720, 426.0200, 427.6750, &
- & 429.2360, 430.7050, 432.0850, 433.3580, 434.5140, &
- & 435.5740, 436.5490, 437.4420, 438.2550, 438.9810, &
- & 439.6110, 440.1430, 440.5770, 440.9450, 441.2660, &
- & 441.5410, 441.7840, 442.0050, 442.2040, 442.3780, &
- & 442.5210, 442.6200, 442.6720, 442.6810, 442.6540, &
- & 442.5830, 442.4670, 442.3270, 442.1680, 441.9960, &
- & 441.8060, 441.5930, 441.3440, 441.0540, 440.7230, &
- & 440.3510, 439.9300, 439.4650, 438.9730, 438.4630, &
- & 437.9400, 437.4020, 436.8400, 436.2640, 435.6850, &
- & 435.1030, 434.5160, 433.9170, 433.3060, 432.7010, &
- & 432.1110, 431.5380, 430.9810, 430.4320, 429.8860, &
- & 429.3370, 428.7810, 428.2220, 427.6490, 427.0660, &
- & 426.4890, 425.9270, 425.3840, 424.8610, 424.3540, &
- & 423.8540, 423.3540, 422.8530, 422.3510, 421.8410, &
- & 421.3250, 420.8190 /
-# else
- !! UKMO, run AJKKH + KAAEF, RCP 8.5, pCO2 time evolution
- DATA hist_pco2 / 286.0230, 286.1730, 286.3230, 286.4480, 286.5730, &
- & 286.7230, 286.8480, 286.9480, 287.0480, 287.1730, &
- & 287.3230, 287.4730, 287.6480, 287.8480, 288.0730, &
- & 288.3480, 288.6480, 288.9730, 289.3470, 289.7470, &
- & 290.1730, 290.6470, 291.1470, 291.6220, 292.0720, &
- & 292.5220, 292.9220, 293.2470, 293.5220, 293.7470, &
- & 293.9470, 294.1220, 294.2720, 294.4220, 294.5470, &
- & 294.6470, 294.7470, 294.8470, 294.9710, 295.1710, &
- & 295.4460, 295.7470, 296.0720, 296.4210, 296.7710, &
- & 297.1460, 297.5710, 298.0210, 298.4460, 298.8460, &
- & 299.2460, 299.6450, 300.0210, 300.3710, 300.7200, &
- & 301.0450, 301.3460, 301.6710, 302.0200, 302.3450, &
- & 302.6450, 302.9700, 303.3450, 303.7200, 304.0700, &
- & 304.4700, 304.9200, 305.3440, 305.7700, 306.2450, &
- & 306.7190, 307.1700, 307.6440, 308.1190, 308.5440, &
- & 308.9440, 309.3440, 309.6940, 309.9440, 310.1190, &
- & 310.2440, 310.3190, 310.3190, 310.2440, 310.1440, &
- & 310.0690, 310.0440, 310.0690, 310.1440, 310.2690, &
- & 310.4440, 310.6940, 311.0430, 311.4440, 311.8690, &
- & 312.3680, 312.9430, 313.5430, 314.1680, 314.7900, &
- & 315.4430, 316.2150, 317.0170, 317.7370, 318.3400, &
- & 318.8680, 319.5900, 320.5890, 321.5470, 322.5770, &
- & 323.8440, 324.9260, 325.7960, 327.0810, 328.6180, &
- & 329.6830, 330.5250, 331.6880, 333.2120, 334.7870, &
- & 336.4640, 338.2990, 339.6660, 340.7310, 342.1360, &
- & 343.7200, 345.2200, 346.7350, 348.5820, 350.6740, &
- & 352.4230, 353.7910, 354.9530, 355.8210, 356.7130, &
- & 358.0630, 359.7720, 361.3970, 363.0900, 365.2560, &
- & 367.2810, 368.7980, 370.4000, 372.4550, 374.6920, &
- & 376.7440, 378.7440, 380.7580, 382.7080, 384.7300, &
- & 386.9420, 389.2540, 391.5670, 393.9370, 396.3920, &
- & 398.9320, 401.5550, 404.2550, 407.0220, 409.8530, &
- & 412.7470, 415.7050, 418.7210, 421.7880, 424.9180, &
- & 428.1200, 431.3970, 434.7470, 438.1650, 441.6410, &
- & 445.1700, 448.7530, 452.3920, 456.0950, 459.8810, &
- & 463.7680, 467.7660, 471.8750, 476.0960, 480.4210, &
- & 484.8390, 489.3470, 493.9430, 498.6400, 503.4380, &
- & 508.3410, 513.3630, 518.5160, 523.8050, 529.2290, &
- & 534.7780, 540.4450, 546.2230, 552.1120, 558.1110, &
- & 564.2110, 570.4130, 576.7390, 583.1990, 589.7980, &
- & 596.5390, 603.4110, 610.4060, 617.4940, 624.6500, &
- & 631.8800, 639.1750, 646.5360, 653.9800, 661.5230, &
- & 669.1840, 676.9570, 684.8290, 692.7790, 700.7690, &
- & 708.8050, 716.8870, 725.0020, 733.1770, 741.3900, &
- & 749.6700, 758.0480, 766.5050, 775.0350, 783.6110, &
- & 792.2200, 800.8740, 809.5680, 818.2760, 827.0090, &
- & 835.8020, 844.6550, 853.5730, 862.5690, 871.6190, &
- & 880.7020, 889.8240, 898.9590, 908.1270, 917.3080, &
- & 926.4960, 935.7040 /
-# endif
+!! JPALM -- change hist_pco2 init
+!!----------------------------------------------------------------------
+!!
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: hist_pco2 !: pCO2
+ INTEGER :: co2_rec
+ REAL(wp) :: co2_yinit, co2_yend
+#endif
+
+!!----------------------------------------------------------------------
+!! JPALM -- PI CO2 key
+!!----------------------------------------------------------------------
+!!
+#if defined key_axy_pi_co2
+ LOGICAL , PUBLIC :: lk_pi_co2 = .TRUE. !: PI xCO2 used
+#else
+ LOGICAL , PUBLIC :: lk_pi_co2 = .FALSE. !: PI xCO2 unused
#endif
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90 (revision 9351)
@@ -81,5 +81,7 @@
gdept_0, gdept_n, &
gdepw_0, gdepw_n, &
- nday_year, nsec_day, nyear, &
+ nday_year, nsec_day, &
+ nyear, nyear_len, ndastp, &
+ nsec_month, &
rdt, tmask, mig, mjg, nimpp, &
njmpp
@@ -92,5 +94,5 @@
mpp_min, mpp_minloc, &
ctl_stop, ctl_warn, lk_mpp
- USE oce, ONLY: tsb, tsn
+ USE oce, ONLY: tsb, tsn, PCO2a_in_cpl
USE par_kind, ONLY: wp
USE par_medusa, ONLY: jpalk, jpchd, jpchn, jpdet, &
@@ -102,5 +104,6 @@
!! JPALM (27-06-2016): add lk_oasis for CO2 and DMS coupling with atm
USE sbc_oce, ONLY: lk_oasis
- USE sms_medusa, ONLY: hist_pco2
+ USE sms_medusa, ONLY: hist_pco2, co2_yinit, co2_yend, &
+ lk_pi_co2
USE trc, ONLY: ln_rsttr, nittrc000, trn
USE bio_medusa_init_mod, ONLY: bio_medusa_init
@@ -114,4 +117,5 @@
USE bio_medusa_diag_slice_mod, ONLY: bio_medusa_diag_slice
USE bio_medusa_fin_mod, ONLY: bio_medusa_fin
+ USE trcstat, ONLY: trc_rst_dia_stat
IMPLICIT NONE
@@ -181,5 +185,6 @@
!!
!! temporary variables
- REAL(wp) :: fq0,fq1,fq2,fq3,fq4
+ REAL(wp) :: fq3,fq4
+ REAL(wp) :: this_y, this_d, this_s, fyear
!!
!! T and S check temporary variable
@@ -293,44 +298,64 @@
!!------------------------------------------------------------------
!!
- !! what's atmospheric pCO2 doing? (data start in 1859)
- iyr1 = nyear - 1859 + 1
- iyr2 = iyr1 + 1
- if (iyr1 .le. 1) then
- !! before 1860
- f_xco2a(:,:) = hist_pco2(1)
- elseif (iyr2 .ge. 242) then
- !! after 2099
- f_xco2a(:,:) = hist_pco2(242)
- else
- !! just right
- fq0 = hist_pco2(iyr1)
- fq1 = hist_pco2(iyr2)
- fq2 = real(nsec_day) / (60.0 * 60.0 * 24.0)
- !! AXY (14/06/12): tweaked to make more sense (and be correct)
-# if defined key_bs_axy_yrlen
- !! bugfix: for 360d year with HadGEM2-ES forcing
- fq3 = (real(nday_year) - 1.0 + fq2) / 360.0
-# else
- !! original use of 365 days (not accounting for leap year or
- !! 360d year)
- fq3 = (real(nday_year) - 1.0 + fq2) / 365.0
-# endif
- fq4 = (fq0 * (1.0 - fq3)) + (fq1 * fq3)
- f_xco2a(:,:) = fq4
- endif
-# if defined key_axy_pi_co2
- !! OCMIP pre-industrial pCO2
- !! f_xco2a(:,:) = 284.725 !! CMIP5 pre-industrial pCO2
- f_xco2a = 284.317 !! CMIP6 pre-industrial pCO2
-# endif
- !! IF(lwp) WRITE(numout,*) ' MEDUSA nyear =', nyear
- !! IF(lwp) WRITE(numout,*) ' MEDUSA nsec_day =', real(nsec_day)
- !! IF(lwp) WRITE(numout,*) ' MEDUSA nday_year =', real(nday_year)
- !! AXY (29/01/14): remove surplus diagnostics
- !! IF(lwp) WRITE(numout,*) ' MEDUSA fq0 =', fq0
- !! IF(lwp) WRITE(numout,*) ' MEDUSA fq1 =', fq1
- !! IF(lwp) WRITE(numout,*) ' MEDUSA fq2 =', fq2
- !! IF(lwp) WRITE(numout,*) ' MEDUSA fq3 =', fq3
- IF(lwp) WRITE(numout,*) ' MEDUSA atm pCO2 =', f_xco2a(1,1)
+ IF (lk_oasis) THEN
+ !! xCO2 from coupled
+ IF ( ( kt == nittrc000 ) .AND. lwp ) &
+ WRITE(numout,*) '** MEDUSA Atm xCO2 given by the UM **'
+ f_xco2a(:,:) = PCO2a_in_cpl(:,:)
+ ELSEIF (lk_pi_co2) THEN
+ !! OCMIP pre-industrial xCO2
+ IF ( ( kt == nittrc000 ) .AND. lwp ) &
+ WRITE(numout,*) '** MEDUSA Atm xCO2 fixed to pre-industrial value **'
+ !! f_xco2a(:,:) = 284.725 !! CMIP5 pre-industrial pCO2
+ f_xco2a(:,:) = 284.317 !! CMIP6 pre-industrial pCO2
+ ELSE
+ !! xCO2 from file
+ !! AXY - JPALM new interpolation scheme usinf nyear_len
+ this_y = real(nyear)
+ this_d = real(nday_year)
+ this_s = real(nsec_day)
+ !!
+ fyear = this_y + ((this_d - 1) + (this_s / (60. * 60. * 24.))) / real(nyear_len(1))
+ !!
+ IF ( ( kt == nittrc000 ) .AND. lwp ) THEN
+ WRITE(numout,*) '** MEDUSA Atm xCO2 from file **'
+ WRITE(numout,*) ' MEDUSA year =', this_y
+ WRITE(numout,*) ' Year length =', real(nyear_len(1))
+ WRITE(numout,*) ' MEDUSA nday_year =', this_d
+ WRITE(numout,*) ' MEDUSA nsec_day =', this_s
+ ENDIF
+ !!
+ !! different case test
+ IF (fyear .LE. co2_yinit) THEN
+ !! before first record -- pre-industrial value
+ f_xco2a(:,:) = hist_pco2(1)
+ ELSEIF (fyear .GE. co2_yend) THEN
+ !! after last record - continue to use the last value
+ f_xco2a(:,:) = hist_pco2(int(co2_yend - co2_yinit + 1.) )
+ ELSE
+ !! just right
+ iyr1 = int(fyear - co2_yinit) + 1
+ iyr2 = iyr1 + 1
+ fq3 = fyear - real(iyr1) - co2_yinit + 1.
+ fq4 = ((1 - fq3) * hist_pco2(iyr1)) + (fq3 * hist_pco2(iyr2))
+ f_xco2a(:,:) = fq4
+ !!
+ IF ( ( kt == nittrc000 ) .AND. lwp ) THEN
+ WRITE(numout,*) ' MEDUSA year1 =', iyr1
+ WRITE(numout,*) ' MEDUSA year2 =', iyr2
+ WRITE(numout,*) ' xCO2 year1 =', hist_pco2(iyr1)
+ WRITE(numout,*) ' xCO2 year2 =', hist_pco2(iyr2)
+ WRITE(numout,*) ' Year2 weight =', fq3
+ ENDIF
+ ENDIF
+ ENDIF
+
+ !! Writing xCO2 in output on start and on the 1st tsp of each month
+ IF ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
+ ( nsec_month .LE. INT(rdt) ) ) THEN
+ IF ( lwp ) WRITE(numout,*) ' *** Atm xCO2 *** -- kt:', kt, &
+ '; current date:', ndastp
+ call trc_rst_dia_stat(f_xco2a(:,:), 'atm xCO2')
+ ENDIF
# endif
@@ -358,6 +383,23 @@
!! x * 30d + 1*rdt(i.e: mod = rdt)
!! ++ need to pass carb-chem output var through restarts
- If ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
- ((86400*mod(nn_date0,100) + mod(kt*rdt,2592000.)) == rdt) ) THEN
+ !!If ( (kt == nitt8rc000 .AND. .NOT.ln_rsttr) .OR. &
+ !! ( (mod(kt*rdt,2592000.)) == rdt) THEN
+ !!=============================
+ !! (Jpalm -- updated for restartability issues)
+ !! We want this to be start of month or if starting afresh from
+ !! climatology - marc 20/6/17
+ !!If ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
+ !! ((86400*mod(nn_date0,100) + mod(kt*rdt,2592000.)) == rdt) ) THEN
+ !!=============================
+ !! Jpalm -- 15-02-2018 -- need to change 3D carb-chem call freq again.
+ !! previous call did not work, probably the (86400*mod(nn_date0,100) part
+ !! should not be in...
+ !! now use the NEMO calendar tool : nsec_month to be sure to call
+ !! at the beginning of a new month .
+ IF ( (kt == nittrc000 .AND. .NOT.ln_rsttr) .OR. &
+ ( nsec_month .LE. INT(rdt) ) ) THEN
+ IF ( lwp ) WRITE(numout,*) &
+ ' *** 3D carb chem call *** -- kt:', kt, &
+ '; current date:', ndastp
!!---------------------------------------------------------------
!! Calculate the carbonate chemistry for the whole ocean on the first
@@ -502,5 +544,5 @@
!! Exceptionnal value did exist
!!
- Call trc_bio_check(kt)
+ Call trc_bio_check(kt, jk)
!!================================================================
@@ -810,5 +852,5 @@
END SUBROUTINE trc_bio_exceptionnal_fix
- SUBROUTINE trc_bio_check(kt)
+ SUBROUTINE trc_bio_check(kt, jk)
!!-----------------------------------
!! JPALM -- 14-12-2017 -- Still dealing with this micro-boil/carb failure
@@ -824,5 +866,5 @@
INTEGER :: ii,ij ! temporary scalars
INTEGER, DIMENSION(2) :: ilocs !
- INTEGER, INTENT( in ) :: kt
+ INTEGER, INTENT( in ) :: kt, jk
!!
!!==========================
@@ -852,6 +894,6 @@
IF(lwp) THEN
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
- WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC > 4000 '
- WRITE(numout,9600) kt, zmax, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: DIC concentration > 4000 '
+ WRITE(numout,9600) kt, zmax, ii, ij, jk
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
ENDIF
@@ -869,6 +911,6 @@
IF(lwp) THEN
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
- WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface DIC <= 0 '
- WRITE(numout,9700) kt, zmin, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: DIC concentration <= 0 '
+ WRITE(numout,9700) kt, zmin, ii, ij, jk
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
ENDIF
@@ -901,6 +943,6 @@
IF(lwp) THEN
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
- WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity > 4000 '
- WRITE(numout,9800) kt, zmax, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ALK concentration > 4000 '
+ WRITE(numout,9800) kt, zmax, ii, ij, jk
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
ENDIF
@@ -918,6 +960,6 @@
IF(lwp) THEN
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** WARNING *****'
- WRITE(numout,*) 'trc_bio:tracer anomaly: sea surface Alkalinity <= 0 '
- WRITE(numout,9900) kt, zmin, ii, ij
+ WRITE(numout,*) 'trc_bio:tracer anomaly: ALK concentration <= 0 '
+ WRITE(numout,9900) kt, zmin, ii, ij, jk
WRITE(numout,*) 'trc_bio:tracer anomaly: ***** END OF WARNING *****'
ENDIF
@@ -925,8 +967,8 @@
-9600 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j: ',2i5)
-9700 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j: ',2i5)
-9800 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j: ',2i5)
-9900 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j: ',2i5)
+9600 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max DIC: ',f16.10,', i j k: ',3i5)
+9700 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min DIC: ',f16.10,', i j k: ',3i5)
+9800 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' max ALK: ',f16.10,', i j k: ',3i5)
+9900 FORMAT ('trc_bio:tracer anomaly: kt=',i6,' min ALK: ',f16.10,', i j k: ',3i5)
END SUBROUTINE trc_bio_check
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcini_medusa.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcini_medusa.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcini_medusa.F90 (revision 9351)
@@ -324,4 +324,14 @@
IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
IF(lwp) CALL flush(numout)
+ !!
+ !!----------------------------------------------------------------------
+ !! JPALM (23-01-2017): new way to initialize CO2-atm for cmip6
+ !! initially done in trcsms_medusa
+ !!----------------------------------------------------------------------
+ !!
+ IF( ( .NOT.lk_oasis ) .AND. ( .NOT.lk_pi_co2 ) ) THEN
+ IF(lwp) WRITE(numout,*) ' trc_ini_medusa: initialisating atm CO2 record'
+ CALL trc_ini_medusa_co2atm
+ ENDIF
END SUBROUTINE trc_ini_medusa
@@ -480,4 +490,75 @@
END SUBROUTINE trc_ini_medusa_river
+ SUBROUTINE trc_ini_medusa_co2atm
+ !!----------------------------------------------------------------------
+ !! *** trc_ini_medusa_co2atm ***
+ !!
+ !! ** Purpose : initialization atmospheric co2 record
+ !!
+ !! ** Method : - Read the xco2 file
+ !!----------------------------------------------------------------------
+ INTEGER :: jn, jm, io, ierr, inum
+ INTEGER :: iskip = 4 ! number of 1st descriptor lines
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: zyy !: xCO2 record years
+ CHARACTER (len=10) :: clname = 'xco2.atm' !! atm CO2 record file
+ !!----------------------------------------------------------------------
+
+ IF(lwp) WRITE(numout,*)
+ IF(lwp) WRITE(numout,*) ' trc_ini_medusa_co2atm: initialisation of atm CO2 historical record'
+ IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~'
+
+
+ IF(lwp) WRITE(numout,*) 'read of formatted file xco2.atm'
+
+ CALL ctl_opn( inum, clname, 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE. )
+ REWIND(inum)
+
+ ! compute the number of year in the file
+ ! file starts in 1849 do jn represent the record number in the file.
+ ! the year is jn + yinit - 1
+ jn = 1
+ DO
+ READ(inum,'(1x)',END=100)
+ jn = jn + 1
+ END DO
+ 100 co2_rec = jn - 1 - iskip
+ IF ( lwp) WRITE(numout,*) ' ', co2_rec ,' years read in the file'
+ ! ! Allocate CO2 hist arrays
+ ierr = 0
+ ALLOCATE( hist_pco2(co2_rec),zyy(co2_rec), STAT=ierr )
+ IF( ierr > 0 ) THEN
+ CALL ctl_stop( 'trc_ini_medusa_co2atm: unable to allocate array' ) ; RETURN
+ ENDIF
+
+ REWIND(inum)
+
+ DO jm = 1, iskip ! Skip over 1st six descriptor lines
+ READ(inum,'(1x)')
+ END DO
+ ! file starts in 1931 do jn represent the year in the century.jhh
+ ! Read file till the end
+ ! allocate start and end year of the file
+ jn = 1
+ DO
+ READ(inum,*, IOSTAT=io) zyy(jn), hist_pco2(jn)
+ IF( io < 0 ) exit
+ IF(jn==1) co2_yinit = zyy(jn)
+ jn = jn + 1
+ END DO
+ co2_yend = co2_yinit + real(co2_rec) - 1.
+
+ IF(lwp) THEN ! Control print
+ WRITE(numout,*)
+ WRITE(numout,*) 'CO2 hist start year: ', co2_yinit
+ WRITE(numout,*) 'CO2 hist end year: ', co2_yend
+ WRITE(numout,*) ' Year xCO2 atm '
+ DO jn = 1, co2_rec
+ WRITE(numout, '( 5F7.1, 6F9.2)') zyy(jn), hist_pco2(jn)
+ END DO
+ ENDIF
+
+ END SUBROUTINE trc_ini_medusa_co2atm
+
+
#else
!!----------------------------------------------------------------------
Index: /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/trcrst.F90
===================================================================
--- /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/trcrst.F90 (revision 9350)
+++ /branches/NERC/dev_r5518_GO6_rev9312_xCO2/NEMOGCM/NEMO/TOP_SRC/trcrst.F90 (revision 9351)
@@ -378,5 +378,5 @@
!!
IF(lwp) WRITE(numout,*) ' CFC averaged properties absent - setting to zero ...'
- qint_cfc(:,:,jn) = 0.0 !! CHN
+ qint_cfc(:,:,jl) = 0.0 !! CHN
ENDIF
!!