Mountain glacier tangent linear model

From Interactive System for Ice sheet Simulation
Jump to: navigation, search
!                           DISCLAIMER
! Removed by P.H.
!
!   This file was generated by TAF version 1.9.54
!
subroutine g_kees1( xx, g_xx, fc, g_fc )
!******************************************************************
!******************************************************************
!** This routine was generated by Automatic differentiation.     **
!** FastOpt: Transformation of Algorithm in Fortran, TAF 1.9.54  **
!******************************************************************
!******************************************************************
!==============================================
! referencing used modules
!==============================================
use kees1_variables
 
!==============================================
! all entries are defined explicitly
!==============================================
implicit none
 
!==============================================
! declare parameters
!==============================================
real(kind=8) a
parameter ( a = 1.e-16 )
real(kind=8) g
parameter ( g = 9.2 )
real(kind=8) n
parameter ( n = 3 )
real(kind=8) rho
parameter ( rho = 920 )
real(kind=8) c
parameter ( c = 2*a*(rho*g)**n/(n+2)*1000.**n )
real(kind=8) dt
parameter ( dt = 1/12. )
real(kind=8) tend
parameter ( tend = 5000 )
 
!==============================================
! declare arguments
!==============================================
real(kind=8), intent(inout) :: fc
real(kind=8), intent(inout) :: g_fc
real(kind=8), intent(inout) :: g_xx(int(xend/dx)+1)
real(kind=8), intent(inout) :: xx(int(xend/dx)+1)
 
!==============================================
! declare local variables
!==============================================
real(kind=8) b(int(xend/dx)+1)
real(kind=8) diff(int(xend/dx))
real(kind=8) diffh(:)
allocatable diffh
real(kind=8) g_diff(int(xend/dx))
real(kind=8) g_h(int(xend/dx)+1)
real(kind=8) g_mb(int(xend/dx)+1)
real(kind=8) g_s(int(tend/dt)+1,int(xend/dx)+1)
real(kind=8) h(int(xend/dx)+1)
integer i
real(kind=8) mb(int(xend/dx)+1)
integer :: nt = int(tend/dt)
integer :: nx = int(xend/dx)
real(kind=8) s(int(tend/dt)+1,int(xend/dx)+1)
integer t
real(kind=8) :: xarr(int(xend/dx)+1) = (/((i-1)*dx,i=1,int(xend/dx)+1)/)
 
!----------------------------------------------
! TANGENT LINEAR AND FUNCTION STATEMENTS
!----------------------------------------------
allocate( diffh(2:nx+1) )
b = 1-0.0001*xarr
g_mb = 0._8
mb = 0.004-xarr*0.0002
g_s(1,:) = 0._8
s(1,:) = b
g_mb = g_mb+g_xx
mb = mb+xx
do t = 2, nt+1
  g_s(t,1) = 0._8
  s(t,1) = 0
  g_h = g_s(t-1,:)
  h = s(t-1,:)-b
  diffh = (s(t-1,2:nx+1)-s(t-1,1:nx))/dx
  g_diff = (g_s(t-1,2:nx+1)/dx-g_s(t-1,1:nx)/dx)*c*((h(1:nx)+h(2:nx+1))/float(2))**(2+n)*(n-1)*abs(diffh)**((-2)+n)*sign(1._8,&
&diffh)+(g_h(2:nx+1)+g_h(1:nx))*c/float(2)*(2+n)*((h(1:nx)+h(2:nx+1))/float(2))**(1+n)*abs(diffh)**(n-1)
  diff = c*((h(1:nx)+h(2:nx+1))/2)**(n+2)*abs(diffh)**(n-1)
  g_s(t,2:nx) = g_diff(2:nx)*dt/dx**2*(s(t-1,3:nx+1)-s(t-1,2:nx))-g_diff(1:nx-1)*dt/dx**2*(s(t-1,2:nx)-s(t-1,1:nx-1))+g_mb(2:nx)*&
&dt+g_s(t-1,3:nx+1)*dt/dx**2*diff(2:nx)+g_s(t-1,2:nx)*(1-dt/dx**2*(diff(2:nx)+diff(1:nx-1)))+g_s(t-1,1:nx-1)*dt/dx**2*diff(1:nx-1)
  s(t,2:nx) = s(t-1,2:nx)+dt/dx**2*(diff(2:nx)*(s(t-1,3:nx+1)-s(t-1,2:nx))-diff(1:nx-1)*(s(t-1,2:nx)-s(t-1,1:nx-1)))+mb(2:nx)*dt
  where (s(t,:) .lt. b)
    g_s(t,:) = 0._8
    s(t,:) = b
  end where
end do
g_fc = 0._8
fc = 0.
do i = 1, size(h)
  g_fc = g_fc+g_h(i)*dx
  fc = fc+h(i)*dx
end do
deallocate( diffh )
 
end subroutine g_kees1