next up previous contents
Next: Restarting Runs Up: Specifying initial and boundary Previous: No-slip boundary conditions   Contents


Specifying heat flux

To specify source terms in the temperature equation, you will need to edit the HeatSource() function in suntans/main/sources.c. The implementation assumes adiabatic surface and bottom boundary conditions, and so all heat fluxes are implemented as source terms in the temperature equation of the form

\begin{eqnarray*}
\frac{\partial{T}}{\partial{t}} + \mathbf{u}\cdot\nabla T &=& ...
...ppa_{tv}\right) \frac{\partial{T}}{\partial{z}}\right] + Q(T) ,
\end{eqnarray*}

where $\kappa_T$ is the molecular diffusivity (defined in suntans.dat) and $\kappa_{tv}$ is the vertical turbulent eddy-diffusivity, computed in turbulence.c. terms to simplify the discussion, the temperature equation is given by
\begin{displaymath}
\frac{\partial{T}}{\partial{t}} = Q(T) .
\end{displaymath} (1)

The source term can be linearized about some temperature $T_0$ so that, defining $T=T_0 + \Delta T$, we have

\begin{eqnarray*}
Q(T) &=& Q(T_0 + \Delta T) ,\\
&=& Q(T_0) + \left(\frac{\pa...
... + O(\Delta T^2) ,\\
&=& A(T_0) + B(T_0) T + O(\Delta T^2) ,
\end{eqnarray*}

where

\begin{eqnarray*}
A(T_0) = Q(T_0) - \left(\frac{\partial{Q}}{\partial{T}}\right)...
...
B(T_0) = \left(\frac{\partial{Q}}{\partial{T}}\right)_{T_0} .
\end{eqnarray*}

The heat equation is discretized with the theta method so that equation (1) becomes

\begin{eqnarray*}
\frac{T^{n+1} - T^n}{\Delta t} &=& A(T_0) + B(T_0)\left[\theta T^{n+1}+(1-\theta)T^n\right] ,
\end{eqnarray*}

where $\theta$ is specified by the value of theta in suntans.dat. After rearranging, the temperature at the new time step is obtained with

\begin{eqnarray*}
T^{n+1} = \left[1 - \theta \Delta t B(T_0)\right]^{-1}\left[T^n + \Delta t A(T_0)
+ (1-\theta) \Delta t B(T_0)T^n\right] .
\end{eqnarray*}

Typically the best linearization is to linearize about the temperature at time step $n$ so that $T_0=T^n$. A simple implementation of heat flux is provided in the windstress example in suntans/main/examples/windstress and the code is in the HeatSource() function in suntans/main/examples/windstress/sources.c.
next up previous contents
Next: Restarting Runs Up: Specifying initial and boundary Previous: No-slip boundary conditions   Contents
2014-08-06