shelf_density Function

  pure function shelf_density(this) result(density)
    !* Author: Christopher MacMackin
    !  Date: April 2016
    !
    ! Returns the density of the ice in the shelf, which is assumed to be
    ! uniform across its domain.
    !
    ! @NOTE Based on my approach to non-dimensionalisation, I'm pretty
    ! sure the density should always be 1, making this method
    ! unneccessary.
    !
    class(ice_shelf), intent(in) :: this
    real(r8)                     :: density !! The ice density.
    density = 1.0_r8/1.12_r8 !TODO: Will probably want to change this at some point
#ifdef DEBUG
    call logger%debug('ice_shelf%density','Ice shelf has density '// &
                      trim(str(density))//'.')
#endif
  end function shelf_density