do_temperature - temperature equation

INTERFACE:

    subroutine do_temperature(n)
DESCRIPTION:

Here, one time step for the temperature equation is performed. First, preparations for the call to the advection schemes are made, i.e. calculating the necessary metric coefficients. After the call to the advection schemes, which actually perform the advection (and horizontal diffusion) step as an operational split step, the solar radiation at the interfaces (rad(k)) is calculated from given surface radiation (swr_loc) by means of a double exponential approach, see equation (102) on page [*]). An option to reflect part of the short wave radiation that reaches the bottom has been implemented. In very shallow waters - or with very clear waters - a significant part of the incoming radiation will reach the bottom. Setting swr_bot_refl_frac to a value between 0 and 1 will reflect this fraction of what ever the value of SWR is at the bottom. The default value of swr_bot_refl_frac is 0. The reflection is only done if the ratio between the surface and bottom values of SWR is greater than swr_min_bot_frac (default 0.01). Furthermore, the surface heat flux sfl_loc is given a value. The sea surface temperature is limited by the freezing point temperature (as a most primitive sea ice model). The next step is to set up the tri-diagonal matrix for calculating the new temperature by means of a semi-implicit central scheme for the vertical diffusion. Source terms which appear on the right hand sides are due to the divergence of the solar radiation at the interfaces. The subroutine is completed by solving the tri-diagonal linear equation by means of a tri-diagonal solver. USES:

    use time, only: month,timestr
    use advection_3d, only: do_advection_3d
    use variables_3d, only: dt,cnpar,hn,ho,nuh,uu,vv,ww,hun,hvn,S
    use domain,       only: imin,imax,jmin,jmax,kmax,az
    use meteo,        only: swr,shf
    use parameters,   only: rho_0,cp
    use parameters, only: avmolt
    use getm_timers, only: tic, toc, TIM_TEMP, TIM_MIXANALYSIS
    use variables_3d, only: do_numerical_analyses
    use variables_3d, only: nummix3d_T,nummix2d_T
    use variables_3d, only: phymix3d_T,phymix2d_T
  $ use omp_lib
    IMPLICIT NONE
INPUT PARAMETERS:
   integer, intent(in) :: n
LOCAL VARIABLES:
    integer                   :: i,j,k,rc
    REALTYPE                  :: T2(I3DFIELD)
   OMP-NOTE: The pointer declarations is to allow each omp thread to
     have its own work storage (over a vertical).
    REALTYPE, POINTER         :: Res(:)
    REALTYPE, POINTER         :: auxn(:),auxo(:)
    REALTYPE, POINTER         :: a1(:),a2(:),a3(:),a4(:)
    REALTYPE, POINTER         :: rad1d(:)
    REALTYPE                  :: zz,swr_loc,shf_loc
    REALTYPE                  :: swr_refl
    REALTYPE                  :: rho_0_cpi
    integer                   :: status