A Modular Arbitrary-Order Ocean-Atmosphere Model -- Stochastic implementation
rk4_tl_ad_integrator.f90
Go to the documentation of this file.
1 
2 ! tl_ad_integrator.f90
3 !
4 !> Tangent Linear (TL) and Adjoint (AD) model versions of MAOOAM.
5 !> Integrators module.
6 !
7 !> @copyright
8 !> 2016 Lesley De Cruz, Jonathan Demaeyer & Sebastian Schubert.
9 !> See LICENSE.txt for license information.
10 !
11 !---------------------------------------------------------------------------!
12 !
13 !> @remark
14 !> This module actually contains the RK4 algorithm routines.
15 !> The user can modify it according to its preferred integration scheme.
16 !> For higher-order schemes, additional bufers will probably have to be defined.
17 !
18 !---------------------------------------------------------------------------
19 
20 
21 
22 MODULE tl_ad_integrator
23 
24  USE params, only: ndim
25  USE aotensor_def, only: aotensor
26 
27  USE tl_ad_tensor , only: ad,tl
28  IMPLICIT NONE
29 
30  PRIVATE
31 
32  REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: buf_y1 !< Buffer to hold the intermediate position of the tangent linear model
33  REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: buf_ka !< Buffer to hold tendencies in the RK4 scheme for the tangent linear model
34  REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: buf_kb !< Buffer to hold tendencies in the RK4 scheme for the tangent linear model
35 
36 
38 
39 CONTAINS
40 
41  !> Routine to initialise the TL-AD integration bufers.
42  SUBROUTINE init_tl_ad_integrator
43  INTEGER :: allocstat
44  ALLOCATE(buf_y1(0:ndim),buf_ka(0:ndim),buf_kb(0:ndim),stat=allocstat)
45  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
46  END SUBROUTINE init_tl_ad_integrator
47 
48 
49  !-----------------------------------------------------!
50  ! !
51  ! Adjoint model integrator !
52  ! !
53  !-----------------------------------------------------!
54 
55 
56 
57 
58  !> Routine to perform an integration step (RK4 algorithm) of the adjoint model. The incremented time is returned.
59  !> @param y Initial point.
60  !> @param ystar Adjoint model at the point ystar.
61  !> @param t Actual integration time
62  !> @param dt Integration timestep.
63  !> @param res Final point after the step.
64  SUBROUTINE ad_step(y,ystar,t,dt,res)
65  REAL(KIND=8), DIMENSION(0:ndim), INTENT(IN) :: y,ystar
66  REAL(KIND=8), INTENT(INOUT) :: t
67  REAL(KIND=8), INTENT(IN) :: dt
68  REAL(KIND=8), DIMENSION(0:ndim), INTENT(OUT) :: res
69 
70  CALL ad(t,ystar,y,buf_ka)
71  buf_y1 = y+0.5*dt*buf_ka
72  CALL ad(t+0.5*dt,ystar,buf_y1,buf_kb)
73  buf_y1 = y+0.5*dt*buf_kb
74  buf_ka = buf_ka+2*buf_kb
75  CALL ad(t+0.5*dt,ystar,buf_y1,buf_kb)
76  buf_y1 = y+0.5*dt*buf_kb
77  buf_ka = buf_ka+2*buf_kb
78  CALL ad(t+dt,ystar,buf_y1,buf_kb)
80  res=y+buf_ka*dt/6
81  t=t+dt
82  END SUBROUTINE ad_step
83 
84 
85  !-----------------------------------------------------!
86  ! !
87  ! Tangent linear model integrator !
88  ! !
89  !-----------------------------------------------------!
90 
91  !> Routine to perform an integration step (RK4 algorithm) of the tangent linear model. The incremented time is returned.
92  !> @param y Initial point.
93  !> @param ystar Adjoint model at the point ystar.
94  !> @param t Actual integration time
95  !> @param dt Integration timestep.
96  !> @param res Final point after the step.
97  SUBROUTINE tl_step(y,ystar,t,dt,res)
98  REAL(KIND=8), DIMENSION(0:ndim), INTENT(IN) :: y,ystar
99  REAL(KIND=8), INTENT(INOUT) :: t
100  REAL(KIND=8), INTENT(IN) :: dt
101  REAL(KIND=8), DIMENSION(0:ndim), INTENT(OUT) :: res
102 
103  CALL tl(t,ystar,y,buf_ka)
104  buf_y1 = y+0.5*dt*buf_ka
105  CALL tl(t+0.5*dt,ystar,buf_y1,buf_kb)
106  buf_y1 = y+0.5*dt*buf_kb
107  buf_ka = buf_ka+2*buf_kb
108  CALL tl(t+0.5*dt,ystar,buf_y1,buf_kb)
109  buf_y1 = y+0.5*dt*buf_kb
110  buf_ka = buf_ka+2*buf_kb
111  CALL tl(t+dt,ystar,buf_y1,buf_kb)
113  res=y+buf_ka*dt/6
114  t=t+dt
115  END SUBROUTINE tl_step
116 
117 END MODULE tl_ad_integrator
integer ndim
Number of variables (dimension of the model)
Definition: params.f90:85
Statistics accumulators.
Definition: stat.f90:14
Tangent Linear (TL) and Adjoint (AD) model versions of MAOOAM. Integrators module.
subroutine, public ad_step(y, ystar, t, dt, res)
Routine to perform an integration step (Heun algorithm) of the adjoint model. The incremented time is...
The equation tensor for the coupled ocean-atmosphere model with temperature which allows for an exten...
subroutine, public ad(t, ystar, deltay, buf)
Tendencies for the AD of MAOOAM in point ystar for perturbation deltay.
subroutine, public tl(t, ystar, deltay, buf)
Tendencies for the TL of MAOOAM in point ystar for perturbation deltay.
The model parameters module.
Definition: params.f90:18
Tangent Linear (TL) and Adjoint (AD) model versions of MAOOAM. Tensors definition module...
real(kind=8), dimension(:), allocatable buf_y1
Buffer to hold the intermediate position (Heun algorithm) of the tangent linear model.
subroutine, public init_tl_ad_integrator
Routine to initialise the integration buffers.
type(coolist), dimension(:), allocatable, public aotensor
- Tensor representation of the tendencies.
real(kind=8), dimension(:), allocatable buf_ka
Buffer to hold tendencies in the RK4 scheme for the tangent linear model.
subroutine, public tl_step(y, ystar, t, dt, res)
Routine to perform an integration step (Heun algorithm) of the tangent linear model. The incremented time is returned.
real(kind=8), dimension(:), allocatable buf_kb
Buffer to hold tendencies in the RK4 scheme for the tangent linear model.