init_temp_water_age_m.F90

Go to the documentation of this file.

!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!
!  Module :  i n i t _ t e m p _ w a t e r _ a g e _ m
!
!> @file
!!
!! Initial temperature, water content and age.
!!
!! @section Copyright
!!
!! Copyright 2009-2017 Ralf Greve, Thorben Dunse
!!
!! @section License
!!
!! This file is part of SICOPOLIS.
!!
!! SICOPOLIS is free software: you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation, either version 3 of the License, or
!! (at your option) any later version.
!!
!! SICOPOLIS is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with SICOPOLIS.  If not, see <http://www.gnu.org/licenses/>.
!<
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age.
!<------------------------------------------------------------------------------
module init_temp_water_age_m

  use sico_types_m
  use sico_variables_m
  use sico_vars_m

  implicit none

  private
  public :: init_temp_water_age_1_1, init_temp_water_age_1_2
  public :: init_temp_water_age_1_3, init_temp_water_age_1_4
  public :: init_temp_water_age_2

contains

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age
!! (case ANF_DAT==1, TEMP_INIT==1:
!! present-day initial topography, isothermal conditions).
!<------------------------------------------------------------------------------
  subroutine init_temp_water_age_1_1()

  implicit none

  integer(i4b) :: i, j, kc

!-------- Initial ice temperature --------

  do i=0, imax
  do j=0, jmax

#if (defined(NMARS) || defined(SMARS))   /* Polar caps of Mars */

     do kc=0, kcmax
        temp_c(kc,j,i) = -100.0_dp
     end do

#else   /* all other domains */

     do kc=0, kcmax
        temp_c(kc,j,i) =  -10.0_dp
     end do

#endif

  end do
  end do

!-------- Initial lithosphere temperature, water content and age --------

  call init_temp_r()
  call init_water()
  call init_age()

  end subroutine init_temp_water_age_1_1

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age
!! (case ANF_DAT==1, TEMP_INIT==2:
!! present-day initial topography,
!! ice temperature equal to local surface temperature).
!<------------------------------------------------------------------------------
  subroutine init_temp_water_age_1_2()

  implicit none

  integer(i4b) :: i, j, kc

!-------- Initial ice temperature --------

  do i=0, imax
  do j=0, jmax

     do kc=0, kcmax
        temp_c(kc,j,i) = temp_s(j,i)
     end do

  end do
  end do

!-------- Initial lithosphere temperature, water content and age --------

  call init_temp_r()
  call init_water()
  call init_age()

  end subroutine init_temp_water_age_1_2

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age
!! (case ANF_DAT==1, TEMP_INIT==3:
!! present-day initial topography,
!! ice temperature linearly increasing with depth).
!<------------------------------------------------------------------------------
  subroutine init_temp_water_age_1_3()

  use ice_material_properties_m, only : kappa_val

  implicit none

  integer(i4b) :: i, j, kc
  real(dp)     :: kappa_const_val
  real(dp)     :: temp_ice_base

!-------- Initial ice temperature --------

#if (defined(NMARS) || defined(SMARS))   /* Polar caps of Mars */
  kappa_const_val = kappa_val(-100.0_dp)
#else   /* all other domains */
  kappa_const_val = kappa_val(-10.0_dp)
#endif

  do i=0, imax
  do j=0, jmax

     if (maske(j,i)<=2_i2b) then

        do kc=0, kcmax

           temp_c(kc,j,i) = temp_s(j,i) &
                            + (q_geo(j,i)/kappa_const_val) &
                              *h_c(j,i)*(1.0_dp-eaz_c_quotient(kc))
                            ! linear temperature distribution according to the
                            ! geothermal heat flux
        end do

        if (temp_c(0,j,i) >= -beta*h_c(j,i)) then

           temp_ice_base = -beta*h_c(j,i)

           do kc=0, kcmax
              temp_c(kc,j,i) = temp_s(j,i) &
                               + (temp_ice_base-temp_s(j,i)) &
                                 *(1.0_dp-eaz_c_quotient(kc))
           end do

        end if

     else   ! maske(j,i)==3_i2b, floating ice

        temp_ice_base = -beta*h_c(j,i) - delta_tm_sw

        do kc=0, kcmax
           temp_c(kc,j,i) = temp_s(j,i) &
                            + (temp_ice_base-temp_s(j,i)) &
                              *(1.0_dp-eaz_c_quotient(kc))
        end do

     end if

  end do
  end do

!-------- Initial lithosphere temperature, water content and age --------

  call init_temp_r()
  call init_water()
  call init_age()

  end subroutine init_temp_water_age_1_3

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age
!! (case ANF_DAT==1, TEMP_INIT==4:
!! present-day initial topography, ice temperature from Robin (1955) solution).
!<------------------------------------------------------------------------------
  subroutine init_temp_water_age_1_4()

  use ice_material_properties_m, only : kappa_val, c_val

  implicit none

  integer(i4b) :: i, j, kc
  real(dp)     :: kappa_const_val, c_const_val
  real(dp)     :: as_val, H_val, qgeo_val, K, z_above_base
  real(dp)     :: erf_val_1, erf_val_2
  real(dp)     :: temp_ice_base, temp_scale_factor

!-------- Initial ice temperature --------

#if (defined(NMARS) || defined(SMARS))   /* Polar caps of Mars */
  kappa_const_val = kappa_val(-100.0_dp)
      c_const_val =     c_val(-100.0_dp)
#else   /* all other domains */
  kappa_const_val = kappa_val(-10.0_dp)
      c_const_val =     c_val(-10.0_dp)
#endif

  do i=0, imax
  do j=0, jmax

     if (maske(j,i)<=2_i2b) then
        as_val = max(as_perp(j,i), epsi)
     else   ! maske(j,i)==3_i2b, floating ice
        as_val = epsi   ! this will produce an almost linear temperature profile
     end if

     h_val    = max(h_c(j,i)  , eps)
     qgeo_val = max(q_geo(j,i), eps)

     k = sqrt( (kappa_const_val/(rho*c_const_val)) * (h_val/as_val) )

     erf_val_1 = erf(h_c(j,i)/(sqrt(2.0_dp)*k))

     do kc=0, kcmax
        z_above_base   = h_c(j,i)*eaz_c_quotient(kc)
        erf_val_2      = erf(z_above_base/(sqrt(2.0_dp)*k))
        temp_c(kc,j,i) = temp_s(j,i) &
                          + (qgeo_val/kappa_const_val) &
                            * sqrt(0.5_dp*pi)*k*(erf_val_1-erf_val_2)
     end do

     if ( (maske(j,i) <= 2_i2b).and.(temp_c(0,j,i) >= -beta*h_c(j,i)) ) then
        temp_ice_base     = -beta*h_c(j,i)
        temp_scale_factor = (temp_ice_base-temp_s(j,i)) &
                                  /(temp_c(0,j,i)-temp_s(j,i))
     else if (maske(j,i) == 3_i2b) then
        temp_ice_base     = -beta*h_c(j,i)-delta_tm_sw
        temp_scale_factor = (temp_ice_base-temp_s(j,i)) &
                                /(temp_c(0,j,i)-temp_s(j,i))
     else
        temp_scale_factor = 1.0_dp
     end if

     do kc=0, kcmax
        temp_c(kc,j,i) = temp_s(j,i) &
                         + temp_scale_factor*(temp_c(kc,j,i)-temp_s(j,i))
     end do

  end do
  end do

!-------- Initial lithosphere temperature, water content and age --------

  call init_temp_r()
  call init_water()
  call init_age()

  end subroutine init_temp_water_age_1_4

!-------------------------------------------------------------------------------
!> Initial temperature, water content and age
!! (case ANF_DAT==2: ice-free conditions with relaxed bedrock).
!<------------------------------------------------------------------------------
  subroutine init_temp_water_age_2()

  implicit none

  integer(i4b) :: i, j, kc

!-------- Initial ice temperature --------

  do i=0, imax
  do j=0, jmax

     do kc=0, kcmax
        temp_c(kc,j,i) = temp_s(j,i)
     end do

  end do
  end do

!-------- Initial lithosphere temperature, water content and age --------

  call init_temp_r()
  call init_water()
  call init_age()

  end subroutine init_temp_water_age_2

!-------------------------------------------------------------------------------
!> Initial lithosphere temperature.
!<------------------------------------------------------------------------------
  subroutine init_temp_r()

  implicit none

  integer(i4b) :: i, j, kr

  do i=0, imax
  do j=0, jmax

     do kr=0, krmax
        temp_r(kr,j,i) = temp_c(0,j,i) &
                         + (q_geo(j,i)/kappa_r) &
                           *h_r*(1.0_dp-zeta_r(kr))
              ! linear temperature distribution according to the
              ! geothermal heat flux
     end do

  end do
  end do

  end subroutine init_temp_r

!-------------------------------------------------------------------------------
!> Initial water content.
!<------------------------------------------------------------------------------
  subroutine init_water()

  implicit none

  omega_c = 0.0_dp   ! only required for the enthalpy method
  omega_t = 0.0_dp

  end subroutine init_water

!-------------------------------------------------------------------------------
!> Initial age.
!<------------------------------------------------------------------------------
  subroutine init_age()

  implicit none

#if (defined(ASF))   /* Austfonna */

  age_c =  3500.0_dp*year_sec
  age_t =  3500.0_dp*year_sec

#elif (defined(NMARS) || defined(SMARS))   /* Polar caps of Mars */

  age_c =  1.0e+06_dp*year_sec
  age_t =  1.0e+06_dp*year_sec

#else   /* all other domains */

  age_c = 15000.0_dp*year_sec
  age_t = 15000.0_dp*year_sec

#endif

  end subroutine init_age

!-------------------------------------------------------------------------------

end module init_temp_water_age_m
!