6 Ground surface

The time integration of 1D thermal conduction equation of grand surface (equation (55) in Part I is performed by the Crank-Nicolson scheme. The space differencing is evaluated by the second order centered scheme. The grand temperature and depth are evaluated on the grid point and the heat flux is evaluated on the half grid point. The number of vertical grid point is $J'$ and the suffix of the lowest grid point is $j = 1$. The $T_{i,J'}$ is assumed to the surface temperature $T_{sfc,i}$. The finite difference 1D thermal conduction equation is represented as follows.

$$\frac{T_{i,j}^{n+1}-T_{i,j}^{n}}{\Delta t} = \frac{\kappa}{4\Delta z_{j}} \left( \frac{T_{i,j+1}^{n+1} - T_{i,... ... - \frac{T_{i,j}^{n+1} - T_{i,j-1}^{n+1}} {\Delta z_{j}+\Delta z_{j-1}} \right.$$

$$+ \left. \frac{T_{i,j+1}^{n} - T_{i,j+1}^{n}} {\Delta z_{j+1}+\De... ...j}} - \frac{T_{i,j}^{n} - T_{i,j-1}^{n}} {\Delta z_{j}+\Delta z_{j-1}} \right).$$ (66)

or,

$$- \frac{\kappa \Delta t}{\Delta z_{j}} \frac{T_{i,j-1}^{n+1}} {\o... ...a t}{\Delta z_{j}} \frac{T_{i,j+1}^{n+1}} {\overline{\Delta z}_{j+\frac{1}{2}}}$$

$$= + \frac{\kappa \Delta t}{\Delta z_{j}} \frac{T_{i,j-1}^{n}} {\ove... ... t}{\Delta z_{j}} \frac{T_{i,j+1}^{n+1}} {\overline{\Delta z}_{j+\frac{1}{2}}}.$$ (67)

where $\kappa = k_{g}/\rho _{g}c_{p,g}$. This equation can be represented in matrix form as follows.

$$\Dvect{A}\cdot \Dvect{T}^{n+1} = \Dvect{B}\cdot \Dvect{T}^{n},$$ (68)

where $\Dvect{T}^{n}=(..., T_{i,j}^{n}, T_{i,j+1}^{n}, T_{i,j+2}^{n}, ...)^{T}$. The elements of $\Dvect{A}$, $\Dvect{B}$ are represented as follows.

$$\begin{eqnarray*} && A_{jj}=4 + \frac{\kappa \Delta t}{\Delta z_{j}} \left( \f... ...t}{\Delta z_{i}} \frac{1}{\overline{\Delta z}_{i-\frac{1}{2}}}, \end{eqnarray*}$$

Considering the boundary condition of upper and lower boundaries, (68) is modified as follows.

$$\Dvect{A}\cdot \Dvect{T}^{n+1} = \Dvect{B}\cdot \Dvect{T}^{n} + \Dvect{S}$$ (69)

Therefore, the grand temperature is given by the solution of the following matrix equation.

$$\Dvect{T}^{n+1} = \Dvect{A}^{-1} \cdot (\Dvect{B}\cdot \Dvect{T}^{n} + \Dvect{S}),$$ (70)

where the elements of $\Dvect{A}$ and $\Dvect{B}$ are modified as follows.

$$\begin{eqnarray*} && A_{11}=4 + \frac{\kappa \Delta t}{\Delta z_{1}} \left( \f... ... \left( \frac{1}{\overline{\Delta z}_{J'-\frac{1}{2}}} \right), \end{eqnarray*}$$

$\Dvect{S}$ is a column vector whose dimension is $J'$ are represented as follows.

$$S_{j} = \left\{ \begin{array}{ll} \frac{\Delta t}{\rho _{... ...F_{IR,net} + H], & j=J' \\ 0, & j\neq J' \end{array}\right.$$