Assessment of Ice Volume Changes in the Cryosphere via Simplified Heat Transport Model ()
1. Introduction
The use of Global Climate Models (or Global Circulation Models, GCMs) appears to be becoming a main stream course in the analysis of climate change effects in order to develop public policy, despite the mixed reactions documented in recent literature regarding the success of their use (e.g., [1] -[3] ). The trend towards the inclusion of additional components of the heat budget and its several transport mechanisms in GCM model development has highlighted the need to better address uncertainty problems, as shown in the focused modeling efforts of some individual transport mechanisms (for example, see [4] ). There are, however, some aspects of the climate which, depending on the time scale involved and the primary heat budget elements may be more predictable. For example, some efforts such as [3] [5] [6] suggest that due to the dominating effects of the latent heat of phase change in freezing and thawing, a model of the heat budget associated with the cyrosphere may be more tractable. In this paper, model parameters from a differential equation based solution representing temperature and ice volume during the Pleistocene are used to approximate effects of climate change.
2. Mathematical Model
The proposed conceptual model will directly relate estimates of temperature found in the literature (as derived from proxy data, such as the ratio dD; for example, see [7] -[13] ) to estimates of ice volume (as derived from proxy data, such as the ratio d18O; see [7] [8] [14] ). Similar to the underpinnings of the Stefan model (e.g., [15] ) for freezing and thawing soils or other medium, the heat budget relationship can be simplified due to the dominating effects of freezing and thawing of water in thecryosphere. From the resulting mathematical analog of temperature and ice volume, an exact solution is obtained. Model calibration is achieved with normalized data of frozen material, I(t), for the Pleistocene time period. To simplify the calibration process, the governing differential equation is normalized to reduce the number of parameters.
The heat transport effects that impact the problem domain include solar and internal heat within the planet, circulation of the atmosphere and oceans with corresponding interface and convective heat transport effects, among others. All of these components interface and integrate together and are represented by the term [Heat] . (The notation [Heat] includes all the components of heat transport that effect phase change in the cryosphere. The total heat budget can be found in [16] , among others.) [Heat] includes all the thermal balance elements of the heat budget, with magnitude, H(t). The function H(t) is determined by calibration to the historic temperature T(t) function developed from historic data. In the current approach, H(t) is developed by examining what H(t) did in the past with respect the T(t) function. This modeling approach of building H(t) based upon what H(t) did in the past may sidestep many complexities that are involved in assembling approximations of the various sub- processes and feedbacks that form the climate and associated boundary conditions. Furthermore, such an H(t) can still be modified to represent changes in the environment.
3. Data and Use of Proxy Relationships
Details regarding conceptual model parameters, data types used and corresponding proxy relationships are presented in [4] but are summarized here for the reader’s convenience. The solution to the initial basic conceptual model relates an input forcing function, H(t), to the output I(t) function. The input function is assumed to follow the trends of reported data, such as ice core temperature estimates ([8] [13] [17] ). Under assumed “natural conditions”, the output function, I(t) , is assumed to also follow the data trends that relate frozen material measure to proxy data such as d18O.
Because normalized data are used, the forcing function, H(t), is linearly related to the selected temperature data proxy for the time domain of interest. That is, for some constants a1 and a2, and assuming H(t) = a1 + a2T(t), then when transformed into normalized N(0,1) form, both functions H(t) and T(t) plot identically. Consequently, one can work directly with the N(0,1) normalized form of a selected T(t) realization. Furthermore, the estimated temperature T(t) is typically given as a linear function of the proxy data, dD (for example, see [7] -[11] , among others). Consequently, the model formulation necessarily is using a linear transform of dD as the forcing function, and therefore the N(0,1) transform of the available dD source data serves to describe the forcing function.
Some papers suggest that ice volume estimates may be written as a linear function of d18O (e.g., [7] [8] [18] ) and other papers suggest that ice volume estimates may be developed using a nonlinear function of d18O (e.g., [19] ). Assuming that the relationship is linear for the range of conditions contemplated in the modeling formulation, the N(0,1) transform of the available d18O source data describes the ice volume trends during the Pleistocene time period. (If the relationship is nonlinear, then the N(0,1) transform of the d18O source data may require further considerations such as including other statistical moments or using another choice of transformation.) In the current paper, d18O data [20] [21] is used as a proxy for representing the trends of ice volume.
The general solution to the conceptual model operates upon estimates of the magnitude of heat, H(t) , assumed to be a linear function of temperature. In turn, these are assumed to be a linear function of dD and produce an output of ice measure that is calibrated to estimates of ice volume. Ice volume is assumed to be a linear function of d18O (as stated previously, the above development can consider a nonlinear transformation of d18O, or other proxy data). Hence, the model solution is a mapping (i.e., a relationship between the two measured variables considered) from dD to d180.
The proxy data types considered in this paper include those listed in Table 1. Each type of proxy shown in bold in Table 2 was collected and processed by partitioning the study time domain into 1000-year intervals (the selection of a 1000-year partition is consistent with the partition size used in [22] -[24] in a study of phase-space modeling for Pleistocene ice volume). The resulting partition was modified to better fit the source data relative to the maximum and minimum data in order to capture the respective peak and minimal values. The modified partitioned data were then transformed using the standard normal transformation N(0,1) transformation. In such a normalized form, calibration of the forcing function H(t) to a selected proxy and calibration of the governing PDE solution I(t) is made simpler as fewer parameters are involved.
4. Mathematical Analog
In [4] (Equation (1) and following), a mathematical analog of changes in volumetric ice in the cryosphere is developed given by,
(1)
where
= total heat into the system affecting phase change in the cryosphere;
= total volumetric ice in the cryosphere including glaciers, tundra ice, snow and other forms of water subject to freeze/thaw effects;
= an initial condition of
, contemplated as a local maximum value such as during a glacial period;
= lumped parameter representing that portion of heat returned due to phase change effects, such as reflected heat due to snow, among other factors, where
;
= lumped parameter representing the area-averaged latent heat of fusion as averaged throughout the
![]()
Table 2. Climate change definitions.
cryosphere; and
= model time.
For details regarding the underpinnings of Equation (1), the reader is referred to [4] .
The ratio 
in Equation (1) represents the impact of changing geometry during phase change effects,
such as growing glaciers, or changing boundary lengths with respect to interface with ocean currents, among other effects.
It is straightforward to extend 
with an exponent, p, such that 
is considered in Equation (1). For the current work, however, 
is retained.
Dividing 
by model initial condition, ![]()
, gives
![]()
(2)
where the initial condition for![]()
, or![]()
. The literature indicates [21] that the minimum value of ![]()
during the Pleistocene period may have a value of about 0.30. For pure water, the latent heat of fusion has value 80 cal/gm which suggests that ![]()
may be an appropriate range for the parameter,![]()
.
Multiplying by![]()
, setting ![]()
and setting H0 as a reference value for heat (see [4] ), then taking the exponentiation gives,
![]()
(3)
or
![]()
(4)
where G(t) is the total heat applied to the cryosphere between model times s = 0 and s = t (see [4] ).
5. Cumulative Departure Model (CDM)
Various fields of geoscience utilize the Cumulative Departure Model (CDM) when looking at long term trends in data where the data under examination (usually annual data) must be considered in the context of previous data influences (for example, groundwater accumulation with respect to annual rainfall). As such, the CDM method involves three steps: 1) calculate the mean for the complete data set; 2) subtract the mean from each data point to determine its relationship or departure from the mean (i.e., surplus or deficit); and 3) add the departures together to create a cumulative departure at each time step. In Equation (4), the second exponential term, G(t), indicates a long term trend of aggradation or degradation of![]()
. If there is no evidence of such long term aggradation or degradation of![]()
, then the argument of the second exponential term may be assumed zero (i.e.,![]()
) and Equation (4) simplifies to
![]()
(5)
where ![]()
is the outcome of the cumulative departure method applied to the reference temperature (it is noted that the reference temperature may be offset by an arbitrary constant and yet ![]()
is not changed in value).
Carrying forward the formulation of Equation (5) (indicating no long term aggradation or degradation of![]()
), the exponential function may be expanded into a Taylor Series, where neglecting terms of order 2 and higher, gives
![]()
(6)
![]()
(7)
where ![]()
is a cluster of parameters given by
![]()
(8)
It is this formulation that is used to approximate climate change impacts upon our cryosphere based upon the above simplified mathematical model.
6. Modeling Effects on the Cryosphere from Climate Change
An opportunity is afforded, in assessing climate change impacts on the cryosphere, by using the simplified modeling equation formulation of Equations (6) through (8). Specifically, the key model D parameter of Equation (8) is a composite of five parameters which results in the variation in the diffusion parameters analogous to the variation in the components of the parameter cluster. For example, a climate change impact that has the effect of increasing the area averaged latent heat of fusion parameter would be represented by a change in the “r” parameter, which inversely impacts the D parameter. Other such influences follow accordingly.
To assess the impact from such parameter changes, the D parameter is varied by 5%, 10% and 25%. The model of (6) through (8) is then run in reverse time for 435,000 years. Plots of the revised volumetric ice contents versus time (expressed in years before present time) are shown in Figures 1-3.
7. Discussion
From Figures 1-3, the simplified mathematical model formulation simulation results indicate that for changes in the D model parameter of 5% primarily shows more thawing of the cryosphere occurring during interglacials but
![]()
Figure 1. Modeling results with 5% variance of D.
![]()
Figure 2. Modeling results with 10% variance of D.
![]()
Figure 3. Modeling results with 25% variance of D.
that recovery substantially occurs during the next freeze cycle. However, at a 25% variation in the D parameter, significant impacts are indicated in the cryosphere ice volume. Even a 10% variation in D indicates significant change in volumetric ice content predictions.
From Equation (8), a variation in D is achieved by variation of several heat budget parameters or a larger variation by only one. For example, changes in the climate due to increased retention of heat, such as contemplated in some climate change model scenarios by others, can be approximately modeled by a decrease in the model return heat parameter, k. From Equation (8), increasing k implies decreasing values of D. Figures 1-3 demonstrate this simplified model estimate of the effect of both increased and decreased values of the D parameter.
From the figures, even for a 5% parameter variation in D, the indicated minimum volumetric ice content is seen to be substantially lower when compared with respect to the recent experience of volumetric ice content. That is, compared to recent measures, which are occurring during a thaw cycle, the predicted ice quantities would be significantly less. However, the next freezing cycle recovers the ice volume substantially. In contrast, a 25% variation in D indicates a far more substantial impact for both the thawing and freezing cycles. It is noted that in all of these simulations, the entire freeze/thaw records of several hundred thousand years is simulated. Simulating only a portion of the record would indicate less modeling impacts.
8. Conclusion
In this paper, the simple heat budget model formulated in [4] is used to develop estimates of cryosphere ice volume impacts due to small changes in modeling parameters caused by hypothesized climate change effects. The resulting mathematical model is a differential equation that can be solved exactly. The mathematical solution to the heat budget model is shown to be a two-parameter model that can be simplified further into a single parameter cluster, D, which is the product and quotients of the several heat budget parameters. Variations in a heat budget parameter result in a similar variation in the single parameter cluster, D. The mathematical model is shown to be an application of the Cumulative Departure Method used in other fields of geoscience. For the problems considered in this paper, computational power needed is provided by a spreadsheet program. Three scenarios of heat parameter effects are considered, of 5%, 10% and 25% variation. The 5% and even the 10% parameter variations indicate significant impacts to the cryosphere ice volume during thaw periods, but it is indicated that these impacts are partially explained by the smaller amounts of ice volume involved in the cryosphere.
Acknowledgements
Acknowledgements are paid to the United States Military Academy, West Point, NY, Department of Mathematical Sciences, for their support to the authors during this research. Also acknowledged are the several individuals who have participated in particular tasks in developing this paper including but by no means limited to Rene Perez, Laura Hromadka, Michael Barton, T.V. Hromadka III, among others.