Abstract
Based on the measured data and characteristics of sea ice temperature distribution in space and time, this study is intended to consider a parabolic partial differential equation of the thermodynamic field of sea ice (coupled by snow, ice, and sea water layers) with a time-dependent domain and its parameter identification problem. An optimal model with state constraints is presented with the thicknesses of snow and sea ice as parametric variables and the deviation between the calculated and measured sea ice temperatures as the performance criterion. The unique existence of the weak solution of the thermodynamic system is proved. The properties of the identification problem and the existence of the optimal parameter are discussed, and the one-order necessary condition is derived. Finally, based on the nonoverlapping domain decomposition method and semi-implicit difference scheme, an optimization algorithm is proposed for the numerical simulation. Results show that the simulated temperature of sea ice fit well with the measured data, and the better fit is corresponding to the deeper sea ice.
1. Introduction
Arctic sea ice cover affects the exchange of heat, energy, mass, and momentum between the atmosphere and the Arctic ocean, and it is a major component of the Arctic environment. Parkinson and Cavalieri [1] suggested that sea ice in the Arctic is a key indicator of climate change. Sea ice has been the principal threat to the development of the Arctic cruise which is cheap, reliable, and capable of year-around operation. Extremely formidable sea ice cover can restrict the operation of surface ships to a few months each summer, and even then there are areas where icebreaker assistance is still required. Meanwhile, sea ice is also a significant threat for the safe design of offshore structures and harbor facilities. Therefore, the research on sea ice has aroused interests of the scientists in many fields. Unfortunately, the physical parameters of sea ice have been obtained mainly by field measurement until now, which are spare and unsatisfactory due to the difficulties in situ, especially during the polar winter. The parameter identification method can reinforce the field data for improved understanding of the physical mechanism and variabilities of parameters of sea ice because sea ice temperature can be measured continuously and automatically. Moreover, it can help to forecast the variabilities of sea ice and climate [2–4].
Many researchers have been devoted to the establishment and improvement of the thermodynamic model of sea ice system. Maykut and Untersteinter [5] first formally proposed a comprehensive one-dimensional sea ice thermodynamic model named MU model, but the corresponding fluxes of atmosphere and ocean were not involved. Semtner [6] simplified MU model, his model has been widely used in climate simulations, while lacks universal and effectiveness. Parkinson and Washington [7] made an improvement of Semtner's model [6], extended the thermodynamic model to a three-dimensional model on a large scale, and considered especially the impact of the leads, which are stretches of open water within fields of sea ice. Hibler [8] established a thermal-dynamic model and carried out first the numerical simulation using finite difference method. Numerical simulation and parameter identification of sea ice thermodynamic system have been paid increasing attention over the past several years, and the consideration of the thermodynamic process of sea ice has been more detailed and comprehensive [9–12]. The thickness of sea ice is one of the most important parameters in thermodynamic and dynamic models of sea ice, and it describes the vertical scale of sea ice. Unfortunately, owing to the randomness and variability of the distribution of sea ice, the continuous and accurate measurement of the thickness is very difficult. Meanwhile, the continuous changes of sea ice thickness in the numerical simulation of large-scale and long process may lead to comparatively large deviation between the simulated and measured sea ice temperatures. The introduction of the change of sea ice thickness in the thermodynamic system thus cannot only describe the thermodynamic behavior of sea ice better but also make the numerical simulation more accurately, while there is few consideration of it at present [13].
The thermodynamic process of sea ice is a heat conductivity process in nature which can be described by a parabolic partial differential equation generally. Therefore, we describe the variability of sea ice temperature in space and time by a parabolic distributed parameter system, with sea ice temperature as the state variable. A one-dimensional three-layer thermodynamic system coupled by snow, ice, and sea water is considered, and the thicknesses of snow and sea ice are identified. A distributed parameter system in a time-dependent domain is established and transferred to an abstract parabolic evolution system by some mathematical methods, and the existence and uniqueness of the weak solution of the evolution equation are proved. By taking the thicknesses of snow and sea ice as identified variables and the deviation between the calculated and measured temperature of sea ice as the performance criterion, an optimal model with state constraints is presented. Then, the existence of the optimal parameter is discussed, and the one-order necessary condition of the optimality is given. Finally, based on the nonoverlapping domain decomposition method and semi-implicit difference schemes, an optimization algorithm is constructed and applied to study the thickness and temperature of the Arctic sea ice. Through numerical simulation, we obtained the change characteristics of snow and ice thicknesses and temperature distribution in sea ice.
This paper is organized as follows. In Section 2, we present a thermodynamic system coupled by snow, ice, and sea water. Section 3 derives some important properties of a bilinear functional and proves the existence and uniqueness of the weak solution of the evolution equation. Section 4 establishes an optimal model with state constraints and proves the strong continuity of the weak solution with respect to the parameters and the existence of the optimal parameter. Moreover, the one-order necessary condition of the optimality is derived. An numerical optimization algorithm is constructed and applied to numerical simulation in Section 5, and numerical results are presented in Section 6. Some conclusions are presented in the final section.
2. The Coupled Thermodynamic System of Sea Ice
Consider the thermodynamic system of sea ice coupled by snow, ice, and sea water layers (Figure 1), denoted by snow-ice-water system. Since the gradient variation of sea ice temperature in the vertical direction is far greater than that in the horizontal direction, we only consider the heat flux in the vertical direction. Let denote the time (unit: s), the final time, and the observation period. Set .
Let axis be the depth direction of the snow-ice-water system (unit: m), and the point on snow surface at the initial measured time the coordinate origin (Figure 1). Let , , and be the thicknesses of snow, ice, and sea water layers at the time , respectively, then the total thickness at the time can be expressed as
Obviously, the thickness of sea water at the time can be expressed as According to the physical properties of sea ice, we have , .
Let and be the initial thicknesses of snow and sea ice, respectively, which are known constants, then the initial thickness of sea water layer is
Let the initial space domains of snow, ice, and sea water layers be , and , respectively, then the initial space domain of the system can be described as obviously, is a nonempty, bounded, and connected set.
Additionally, let , , and be the space domains of snow, ice, and sea water layers at the time , respectively, then the space domain of the snow-ice-water system at the time can be expressed as
Apparently, the measured space-time domain denoted by is opened and bounded.
By the energy conversation and the Fourier law, the coupled thermodynamic sea ice model is thus described by the following parabolic partial differential dynamic system where is the density, the specific heat and the thermal conductivity, the temperature of the snow-ice-water system at the depth and the time (unit: Kelvin), the temperature at the depth and the time , and the heat source term. , , denotes the snow, ice, and sea water layer, respectively. and are both known functions on the time .
The system (2.6) is piecewise smooth owing to the difference of the densities, specific heats, and thermal conductivities in the three layers. By the physical properties of sea ice, we give the following assumptions.(A1) The thermodynamic parameters (, , and ) remain unchanged during the measured period. (A2) The upper and lower boundaries remain unmoved, namely, the total thickness of the system is a constant (). (A3) is continuous on . (A4), are bounded functions with positive values, and there exist and satisfying , such that .
Additionally, functions in the system (2.6) should be differential on or in some sense; thus, we assume that (A5), , , , .
According to the assumptions (A1) and (A2), we set , and , then the system (2.6) can be described as the following homogeneous boundary value problem denoted by (HBP): where
3. Properties of the Coupled Thermodynamic System
For convenience, we let , , , and is measurable, is the admissible parameter set, is a known nonempty, bounded, closed convex subset on . And let be the spatial domain of the distributed parameter system (HBP) at the time , where .
For any and , let be a separable Hilbert space, , with norms denoted by and , respectively, then . Let and denote the inner product on and , respectively, and the dual spaces of and , respectively, and the dual paring between and , then is a Gelf-triple space with , namely, the embedding is continuous, and is dense in . The inner product and the induced norm on are defined by The inner product and the induced norm on are defined by where and are the derivatives of and on the depth , respectively.
The dual paring between and is expressed as
Let , we multiply both sides of (2.7) by and integrate them over
Using the method of integration by parts and the homogeneous boundary conditions, (3.4) can be described as
For any , we define a functional on as
The properties of the functional are presented in the following.
Property 1. Suppose that assumptions (A1)–(A5) hold, then(1) defined by (3.6) is a bilinear functional on .(2)For all , , is a measurable function on , and there exists , such that .
Proof. (1) For all , , for all , , , , we have
So is a bilinear functional on .
(2) Obviously, for all , , is a measurable function on according to (3.6), and
Taking , then .
Property 2. Suppose that assumptions (A1)–(A5) hold, then, for the bilinear functional defined by (3.6), there exist and , such that , for all .
Proof. The conclusion follows directly from (3.6) and Garding inequality [14].
From the Properties 1 and 2 and Lax-Milgram theorem [15], there exists a bounded continuous linear operator , such that
According to (3.9), the Properties 1 and 2, we can get the following property.
Property 3. Suppose that assumptions (A1)–(A5) hold, then for all , for all , and for all , there exist and , such that
For any , , and , set , then is a continuous linear functional on . According to Riesz representation theorem [14], there exists a unique element , such that
From the above, for each , the system (HBP) can be written as a parabolic evolution equation denoted by (PEE):
where , , and , with the inner product defined as
Obviously, the above process is reversible, the system (PEE) is thus equivalent to (HBP). The definition of weak solution of the system (PEE) is given in the following.
Definition 3.1. A function is a weak solution of (PEE), if satisfies where and .
Theorem 3.2. Suppose that assumptions (A1)–(A5), Properties 1 and 3 hold. For all , let is bounded, and there exists a countable base on , such that , . Then, the system (PEE) has a unique weak solution which depends continuously on and .
Proof. Let , here is the outward normal vector of at . By the definition of , we know that satisfies the Dirichlet condition (, ) for any given . The system (PEE) thus has a unique weak solution which depends continuously on and [15].
From the above theorem, we can also know that the system (HBP) has a unique solution .
4. Parameter Identification of the Coupled System
4.1. Identification Model
We define the performance criterion as with where is the temperature function obtained from the system (2.6), is the observed temperature. Let , then (4.1) is expressed as
The goal of this study is to make the temperature obtained from the system (2.6) fits the measured data as far as possible. Then, the identification model of the system (HBP) is expressed as where is the solution of system (HBP) corresponding to , and is defined by (4.1).
Obviously, from (4.1)–(4.2), we can conclude that the problem (SIP) is equivalent to the following problem where is the solution of system (2.6) corresponding to .
4.2. Strong Continuity of the Weak Solution on Parameters
Theorem 4.1. Suppose that assumptions (A1)–(A5) hold, then the mapping is strongly continuous.
Proof. For any given parameter , let be a feasible parameter sequence, such that as , and and the solutions of (PEE) corresponding to and , respectively. Set , , , then we obtain the following system
Multiply both sides of the first equation of (4.4) by and integrate them over , then we have
For the first term of (4.5) on the left side, we have
For the second term of (4.5) on the left side, according to Property 2, we have
where and are the same with that of Property 2.
Using the elementary inequality and taking , then, for the right terms of (4.5), we have
Let . Substituting (4.6)–(4.8) into (4.5), then we have
Using Gronwall inequality, we obtain
From all the above and , we get , the mapping is thus strongly continuous.
4.3. Existence of Optimal Parameter
Theorem 4.2. Suppose that assumptions (A1)–(A5) hold, then there exists at least one optimal parameter satisfies the identification problem (SIP).
Proof. Obviously, is continuous on , and, from (4.3), . By Theorem 4.1, the mapping is continuous. Hence, the mapping is continuous on . Since is a nonempty, bounded, and closed set, there exists , such that for all , , e.g., is an optimal parameter, thus the desired result is obtained.
4.4. Necessary Optimality Condition
Let be the optimal parameter of the problem (SIP), following the convex compactness of and the continuity of on , we can prove that is Gateaux differentiable at and its Gateaux derivative exists. Hence, we can get the following necessary conditions for optimality.
Theorem 4.3. Suppose that assumptions (A1)–(A5) are valid, and let be the optimal parameter of the system (SIP), then satisfies the following inequality:
5. Optimization and Numerical Algorithm
5.1. Optimization of the Parameter Identification Model
In this section, we aim at constructing an optimal algorithm to solve the parameter identification problem (SIPO). On the domain , we decompose (2.6) as with the initial condition and the boundary and penetrating conditions
Denote the subsystem composed by (5.1), (5.2), and (5.3) as , by (5.1), (5.2), and (5.4) as , and (5.1), (5.2), and (5.5) as .
According to the partial differential theory, we can obtain the following theorem.
Theorem 5.1. Suppose that assumptions (A1)–(A5) are valid, then, for any , the subsystem has a unique weak solution which depends continuously on .
Let , be the numbers of the measured spatial and temporal spots, respectively, the measured depth, the measured time, the calculated temperature of sea ice at depth and time from the system (2.6), the measured temperature data at depth and time , then the practical identification problem can be expressed as
Let is the solution of the subsystem corresponding to , and for the numerical implementation, the problem (SIPOD) can be decomposed to the following three connected subproblems
Obviously, the performance function of (SIPOD) can be described as
From Theorem 4.2, we can get that there also exists at least one optimal parameter satisfying the subproblem .
Set then, the identification problem is expressed as and the subidentification problem is
Theorem 5.2. If is the optimal parameter of SIPOD, then is also that of .
5.2. Numerical Optimization Algorithm
In this section, the semi-implicit finite difference scheme [15] is employed to discretize, and the nonoverlapping Schwarz alternating direction method to solve (5.1). And related to the measured temperature data , a numerical optimization algorithm is constructed with the following steps.
Step 1. Select starting points , directions , , starting step length , , , acceleration factor , accuracy , and the maximum iteration number , set , ;
Step 2. Calculate the numerical solution , by nonoverlapping Schwarz alternating direction method and semi-implicit difference scheme.
Step 3. If , set ; for the case , if , set , else, set .
Step 4. If , set , go to Step 3; for the case , if , go to Step 5, else, go to Step 6.
Step 5. Set , , , and , go to Step 2.
Step 6. If , stop, set , else, set , , , and , go to Step 2.
6. Numerical Results
Using the above optimal algorithm, we can obtain the temperature distribution of sea ice in the Arctic. The data set was measured by a monitor buoy installed on floe ice in the Arctic.
The comparison between the calculated and measured sea ice temperatures from November 1, 2003 to February 29, 2004 was shown in Figure 2, with the test temperatures data being smoothed for the unomitted transparent errors in some of the original data. The measured and simulated data in Figure 2 both indicate a general increasing trend of the temperature of sea ice, while a decreasing increment with increasing depth of sea ice. Meanwhile, because the data had been measured during the freezing period of sea ice, the temperature generally decreases with increasing time. The calculated temperature fits well with the measured data; furthermore, the better fit between the calculated and measured temperatures is corresponding to the deeper sea ice. The temperature deviation between the simulated and measured temperatures is , which is a much small error for the temperature of sea ice. It can be concluded that the numerical results by our method reflect the actual variation of the sea ice temperature distribution and approach the measured data well.
7. Conclusions
In this study, we aim to simulate the temperature of sea ice with thicknesses of snow and sea ice as identified parameters. We established a parabolic distributed system of the temperature field of sea ice with a time-dependent domain in the Arctic sea ice and proposed its parameter identification problem with thicknesses of snow and sea ice as identified variables and the deviation between the calculated and measured sea ice temperatures as the performance criterion. And we proved the unique existence of the weak solution of the thermodynamic system, discussed the existence of the optimal parameter, and derived the one-order necessary condition. Finally, based on the nonoverlapping domain decomposition method and semi-implicit difference scheme, we constructed an optimization algorithm to simulate the sea ice temperature in the Arctic. The results (Figure 2) show that the temperature of sea ice increases with increasing depth of sea ice, while, because the measurement had been made during the freezing period of sea ice, it decreases with increasing time.
The simulated temperatures of sea ice fitted well with the measured data, and the deviation between the simulated and measured temperature of sea ice was slight. It indicates that the numerical results by our method reflect the actual variation of the sea ice temperature distribution in space and time. Meanwhile, it is illuminated that the algorithm proposed in the present study is valid and more suitable for the deeper sea ice.
Nevertheless, it is important to note that our mathematical framework was verified only using in situ data from the Arctic, and the variabilities of the physical parameters of the coupled sea ice system, such as density, specific heat, were not considered, more measurements would be required for a true verification of this method in the future work. However, we believe that our method is a promising approach worthy of further studies under different environmental conditions in the Arctic as well as in the Antarctic.
Acknowledgments
The authors are grateful to Chinese Arctic and Antarctic Administration for the logistic support during the field investigation in the Arctic Ocean. This study was supported by the National Nature Science Foundation of China (Grants no. 40930848, 50921001, and 10871033), and Norwegian Research Council Project AMORA (no. 193592/S30).