A simple energy balance model of global warming

(relaxation to steady state response temperature)

General problem

One important feature of the anthropogenic climate change (global warming) is that the rise of temperature lags the rise of carbon dioxide concentration in the atmosphere. It is not a simple delay (time shift), but relationship between stock and flow (CO2 concentration makes the flow of energy, and temperature is related to the stock of energy). Thus, it involves an integral in time.

The delay can schematically (not very realistically) be expressed in a zero-dimensional energy balance model, where the stock of energy (per unit area of the earth's surface) is equal to surface temperature T multiplied by constant heat capacity (per unit area) C.

If we tentatively assume that energy input to the climate system (per unit area per unit time) is F, the system can be written in the following equation.

C dT(t)/dt = F(t)

If integrated in time,

T(t) = T(0) + (1/C) çt'=0t F(t') dt' .

Steady-state response

We assume that the pre-industrial climate is in a steady state, with CO2 concentration qpreind and temperature Tpreind.

[Note: I just tentatively use "q" for CO2 concentration. It does not mean specific humidity (concentration of water vapor) in this context.]

We assume that the steady-state response of temperature to CO2 doubling is ƒÂTssr,double. [If CO2 concentration would be 2 qpreind for an indefinitely long time, the average temperature would be Tpreind + ƒÂTssr,double.

From experience of modeling studies, steady-state response is approximately proportional to logarithm of CO2 concentration. Thus, it can be written like this.

Tssr(t) = Tpreind + ƒÂTssr,double E{ (log(q(t)) - log(qpreind) }/ log(2)

(Either log-base-10 (common logarithm) or log-base-e (natural logarithm) can be used as far as the same base is used consistently.)

Simple formulation of energy input: relaxation towards steady-state

As for a simple example, we assume that F is proportional to the difference between the actual temperature and the steady-state response temperature.

F = k E { Tssr(q(t)) - T(t) }

If Tssr is constant, T approaches Tssr according to an exponential function of t (with a negative coefficient). In applied physics, this type of behavior is called "relaxation" or "damping".

Specification of a problem with some CO2 concentration scenarios

The RCP (reference concentration pathways) scenarios have carbon dioxide concentration from year 2010 to 2100 by 10 year intervals. It is shown in a table whose value is taken from the Table A II.4.1 of IPCC AR5 (WG1). For convenience, CO2 concentration at year 2000, known based on observations (Table A II.1.1a), is also included in our table.

Pre-industrial level of CO2 concentration qpreind is assumed to be 280 ppm. (According to IPCC AR5 Table A II.1.1a, CO2 concentration took this value in year 1790.)

Pre-industrial level of temperature is not precisely known. IPCC AR5 Table A II.3 contains historical decadal mean global surface air temperature. According to "HadCRUT4_Median" (compiled by UKMO Hadley Centre and Univ. East Anglia's Climatic Research Unit), the temperature in 1850 was -0.320 K and that in 2000 was +0.453 K, relative to 1961-1990 average. We (tentatively) assume that T in 1850 as Tpreind, and show all temperature values as deviation from this pre-industrial level. Temperature in 2000 is 0.773 K, in this convention.

Steady-state response to CO2 doubling ƒÂTssr,double is (tentatively) assumed to be 3.0 K.

Heat capacity is assumed to be that of a homogeneous layer of water. Heat capacity per unit area C is the product of the depth of layer (tentatively assumed to be 500 m), density of water (1 g/cm3 = 1000 kg/m3), and specific heat capacity of water (approx. 4.2 J/(g K) = 4200 J/(kg K) ).

The coefficient of relaxation k cannot be determined by physics. I tentatively give the value 1.5 W/(m2 K), based on empirical matching of results of test runs reproducing the history of the 20th century with observations.

As temperature is expressed in deviation from the pre-industrial level, the model equations are as follows:

C dT(t)/dt = k { Tssr(q(t)) - T(t) }

Tssr(q(t)) = ƒÂTssr, doubleE{ log (q(t)) - log (qpreind) } / log(2)

And the differential equation is approximated by a finite difference formula like this.

T(t + ƒ¢t) = T(t) + ƒ¢tEkE{ Tssr(q(t)) - T(t) } / C

Programing in R

I have made an example program in R language.

R software is designed for statistical data analysis. It is also good for making graphs. It is not specifically designed for time step simulations, but it can be used.

The implementation software of R is available free from CRAN (https://cran.r-project.org ) or its mirror sites. There are packages for MS Windows, Mac, and Linux.

Download the following file of R script in the same directory as the above data. Open it by R, and run.

By one run, it makes one time series of temperature according to one of the RCP scenarios chosen by the variable "kscenario".

The part for output is unfinished. It just dumps values of array variables. (Also note that the value of "tssr" for the year 2110 is invalid.)

[This paragraph is added on 2019-09-01.] A version of R script which displays the results in a graph is added. (This version does not show numbers on the console or files, however.)

[This part is added on 2019-09-16.] A version of R script which writes the results into a file is added.

As the scripts is run without modification, the results will be stored in a file like this.

Another R script which reads this file and makes the same graph as "co2warming_graph.R" is here.

Formatted output will be added later (but I do not promise I finish by what time).

Programing in spreadsheets

The calculation can be done in spreadsheets of MS Excel or Libreoffice Calc. It is done by putting formulas of calculation in the cells. Users must care "absolute reference" and "relative reference" of cells. When a formula is copied, absolute reference (denoted with "$") does not change, but relative reference (denoted without "$") shifts.

(I do not recommend using spreadsheets if you can use R, Python or Ruby.)

1A simple energy balance model of global warming (relaxation to steady state response tempearture) 2019-07-16 Masuda
2constantsyearCO2 conc. [ppm]Steady state response [K]Temp. (base year 1850) [K]
3CO2 pre-industrial [ppm]280RCP2.6RCP4.5RCP6.0RCP8.5RCP2.6RCP4.5RCP6.0RCP8.5RCP2.6RCP4.5RCP6.0RCP8.5
4SSR to CO2 doubling [K]32000368.7=$E4=$E4=$E4=$B$4*(LOG(E4)-LOG($B$3))/LOG(2)0.773=$O4=$O4=$O4
6relaxation coeff. [W/(m2 K)]1.52020412.1411.1409.4415.8«¨¨¨«¨¨¨
8depth of water [m]5002040............«¨¨¨«¨¨¨
9density of water [kg/m3]10002050............«¨¨¨«¨¨¨
10specific heat of water [J/(kg K)]42002060............«¨¨¨«¨¨¨
11heat capacity / area [J/(m2 K)]=$B$8*$B$9*$B$102070............«¨¨¨«¨¨¨
13time step [year]102090............«¨¨¨«¨¨¨
14day / year365.252100............«¨¨¨«¨¨¨
15second / day864002110«¨¨¨
16time step [s]=$B$13*$B$14*$B$15

Cells with "..." should be filled with data values from RCP (shown above). (Perhaps you can open another sheet and take data from file there, and then copy the values and paste onto the worksheet for simulation.)

"«" and "¨" means that the formula shown at the upper left corner should be copied to all cells in 4 columns and 11 rows.


2019-Jul-16, 2019-Jul-17, 2019-Sep-01, 2019-Sep-16
[Index page of "Climate Change 2019"]