A Modular Arbitrary-Order Ocean-Atmosphere Model -- Stochastic implementation
Functions/Subroutines | Variables
mtv_int_tensor Module Reference

The MTV tensors used to integrate the MTV model. More...

Functions/Subroutines

subroutine, public init_mtv_int_tensor
 Subroutine to initialise the MTV tensor. More...
 

Variables

real(kind=8), dimension(:), allocatable, public h1
 First constant vector. More...
 
real(kind=8), dimension(:), allocatable, public h2
 Second constant vector. More...
 
real(kind=8), dimension(:), allocatable, public h3
 Third constant vector. More...
 
real(kind=8), dimension(:), allocatable, public htot
 Total constant vector. More...
 
type(coolist), dimension(:), allocatable, public l1
 First linear tensor. More...
 
type(coolist), dimension(:), allocatable, public l2
 Second linear tensor. More...
 
type(coolist), dimension(:), allocatable, public l3
 Third linear tensor. More...
 
type(coolist), dimension(:), allocatable, public ltot
 Total linear tensor. More...
 
type(coolist), dimension(:), allocatable, public b1
 First quadratic tensor. More...
 
type(coolist), dimension(:), allocatable, public b2
 Second quadratic tensor. More...
 
type(coolist), dimension(:), allocatable, public btot
 Total quadratic tensor. More...
 
type(coolist4), dimension(:), allocatable, public mtot
 Tensor for the cubic terms. More...
 
real(kind=8), dimension(:,:), allocatable, public q1
 Constant terms for the state-dependent noise covariance matrix. More...
 
real(kind=8), dimension(:,:), allocatable, public q2
 Constant terms for the state-independent noise covariance matrix. More...
 
type(coolist), dimension(:), allocatable, public utot
 Linear terms for the state-dependent noise covariance matrix. More...
 
type(coolist4), dimension(:), allocatable, public vtot
 Quadratic terms for the state-dependent noise covariance matrix. More...
 
real(kind=8), dimension(:), allocatable dumb_vec
 Dummy vector. More...
 
real(kind=8), dimension(:,:), allocatable dumb_mat1
 Dummy matrix. More...
 
real(kind=8), dimension(:,:), allocatable dumb_mat2
 Dummy matrix. More...
 
real(kind=8), dimension(:,:), allocatable dumb_mat3
 Dummy matrix. More...
 
real(kind=8), dimension(:,:), allocatable dumb_mat4
 Dummy matrix. More...
 

Detailed Description

The MTV tensors used to integrate the MTV model.

Remarks
See : Franzke, C., Majda, A. J., & Vanden-Eijnden, E. (2005). Low-order stochastic mode reduction for a realistic barotropic model climate. Journal of the atmospheric sciences, 62(6), 1722-1745.

Function/Subroutine Documentation

subroutine, public mtv_int_tensor::init_mtv_int_tensor ( )

Subroutine to initialise the MTV tensor.

Definition at line 89 of file MTV_int_tensor.f90.

89  INTEGER :: allocstat,i,j,k,l
90 
91  print*, 'Initializing the decomposition tensors...'
92  CALL init_dec_tensor
93  print*, "Initializing the correlation matrices and tensors..."
94  CALL init_corrint
95  print*, "Computing the correlation integrated matrices and tensors..."
96  CALL comp_corrint
97 
98  !H part
99  print*, "Computing the H term..."
100 
101  ALLOCATE(h1(0:ndim), h2(0:ndim), h3(0:ndim), htot(0:ndim), stat=allocstat)
102  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
103 
104  ALLOCATE(dumb_mat1(ndim,ndim), dumb_mat2(ndim,ndim), stat=allocstat)
105  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
106 
107  ALLOCATE(dumb_mat3(ndim,ndim), dumb_mat4(ndim,ndim), stat=allocstat)
108  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
109 
110  !H1
111  CALL coo_to_mat_ik(lxy,dumb_mat1)
112  dumb_mat2=matmul(dumb_mat1,corrint)
113  CALL sparse_mul3_with_mat(bxxy,dumb_mat2,h1)
114 
115  ! H2
116  h2=0.d0
117  IF (mode.ne.'ures') THEN
118  CALL coo_to_mat_ik(lyy,dumb_mat1)
119  dumb_mat1=matmul(inv_corr_i_full,dumb_mat1)
120 
121  DO i=1,ndim
122  CALL coo_to_mat_i(i,bxyy,dumb_mat2)
123  CALL sparse_mul4_with_mat_jl(corr2int,dumb_mat2,dumb_mat3)
124  CALL sparse_mul4_with_mat_jl(corr2int,transpose(dumb_mat2),dumb_mat4)
125  dumb_mat3=dumb_mat3+dumb_mat4
126  h2(i)=mat_contract(dumb_mat1,dumb_mat3)
127  ENDDO
128  ENDIF
129 
130  !H3
131  h3=0.d0
132  CALL sparse_mul3_with_mat(bxyy,corr_i_full,h3)
133 
134  !Htot
135  htot=0.d0
136  htot=h1+h2+h3
137 
138  print*, "Computing the L terms..."
139  ALLOCATE(l1(ndim), l2(ndim), l3(ndim), stat=allocstat)
140  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
141 
142  !L1
143  CALL coo_to_mat_ik(lyx,dumb_mat1)
144  CALL coo_to_mat_ik(lxy,dumb_mat2)
145  dumb_mat3=matmul(inv_corr_i_full,corrint)
146  dumb_mat4=matmul(dumb_mat2,matmul(transpose(dumb_mat3),dumb_mat1))
147  CALL matc_to_coo(dumb_mat4,l1)
148 
149  !L2
150  dumb_mat4=0.d0
151  DO i=1,ndim
152  DO j=1,ndim
153  CALL coo_to_mat_i(i,bxyy,dumb_mat1)
154  CALL sparse_mul4_with_mat_jl(corr2int,dumb_mat1+transpose(dumb_mat1),dumb_mat2)
155 
156  CALL coo_to_mat_j(j,byxy,dumb_mat1)
157  dumb_mat1=matmul(inv_corr_i_full,dumb_mat1)
158  dumb_mat4(i,j)=mat_contract(dumb_mat1,dumb_mat2)
159  END DO
160  END DO
161  CALL matc_to_coo(dumb_mat4,l2)
162 
163  !L3
164  dumb_mat4=0.d0
165  DO i=1,ndim
166  DO j=1,ndim
167  CALL coo_to_mat_j(j,bxxy,dumb_mat1)
168  CALL coo_to_mat_i(i,bxxy,dumb_mat2)
169  dumb_mat4(i,j)=mat_trace(matmul(dumb_mat1,matmul(corrint,transpose(dumb_mat2))))
170  ENDDO
171  END DO
172  CALL matc_to_coo(dumb_mat4,l3)
173 
174  !Ltot
175 
176  ALLOCATE(ltot(ndim), stat=allocstat)
177  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
178 
179  CALL add_to_tensor(l1,ltot)
180  CALL add_to_tensor(l2,ltot)
181  CALL add_to_tensor(l3,ltot)
182 
183  print*, "Computing the B terms..."
184  ALLOCATE(b1(ndim), b2(ndim), btot(ndim), stat=allocstat)
185  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
186  ALLOCATE(dumb_vec(ndim), stat=allocstat)
187  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
188 
189  ! B1
190  CALL coo_to_mat_ik(lxy,dumb_mat1)
191  dumb_mat2=matmul(inv_corr_i_full,corrint)
192 
193  dumb_mat3=matmul(dumb_mat1,transpose(dumb_mat2))
194  DO j=1,ndim
195  DO k=1,ndim
196  CALL coo_to_vec_jk(j,k,byxx,dumb_vec)
197  dumb_vec=matmul(dumb_mat3,dumb_vec)
198  CALL add_vec_jk_to_tensor(j,k,dumb_vec,b1)
199  ENDDO
200  END DO
201 
202  ! B2
203  CALL coo_to_mat_ik(lyx,dumb_mat3)
204  dumb_mat2=matmul(inv_corr_i_full,corrint)
205 
206  dumb_mat4=matmul(transpose(dumb_mat2),dumb_mat3)
207  DO i=1,ndim
208  CALL coo_to_mat_i(i,bxxy,dumb_mat1)
209  dumb_mat2=matmul(dumb_mat1,dumb_mat4)
210  CALL add_matc_to_tensor(i,dumb_mat2,b2)
211  ENDDO
212 
213  CALL add_to_tensor(b1,btot)
214  CALL add_to_tensor(b2,btot)
215 
216  !M
217 
218  print*, "Computing the M term..."
219 
220  ALLOCATE(mtot(ndim), stat=allocstat)
221  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
222 
223  dumb_mat2=matmul(inv_corr_i_full,corrint)
224 
225  DO i=1,ndim
226  CALL coo_to_mat_i(i,bxxy,dumb_mat1)
227  dumb_mat3=matmul(dumb_mat1,transpose(dumb_mat2))
228  DO k=1,ndim
229  DO l=1,ndim
230  CALL coo_to_vec_jk(k,l,byxx,dumb_vec)
231  dumb_vec=matmul(dumb_mat3,dumb_vec)
232  CALL add_vec_ikl_to_tensor4(i,k,l,dumb_vec,mtot)
233  ENDDO
234  END DO
235  END DO
236 
237  !Q
238 
239  print*, "Computing the Q terms..."
240  ALLOCATE(q1(ndim,ndim), q2(ndim,ndim), stat=allocstat)
241  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
242 
243  !Q1
244 
245  CALL coo_to_mat_ik(lxy,dumb_mat1)
246  q1=matmul(dumb_mat1,matmul(corrint,transpose(dumb_mat1)))
247 
248  !Q2
249 
250  DO i=1,ndim
251  DO j=1,ndim
252  CALL coo_to_mat_i(i,bxyy,dumb_mat1)
253  CALL coo_to_mat_i(j,bxyy,dumb_mat2)
254  CALL sparse_mul4_with_mat_jl(corr2int,dumb_mat2,dumb_mat3)
255  CALL sparse_mul4_with_mat_jl(corr2int,transpose(dumb_mat2),dumb_mat4)
256  dumb_mat2=dumb_mat3+dumb_mat4
257  q2(i,j)=mat_contract(dumb_mat1,dumb_mat2)
258  END DO
259  END DO
260 
261  !U
262 
263  ALLOCATE(utot(ndim), stat=allocstat)
264  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
265 
266  CALL coo_to_mat_ik(lxy,dumb_mat1)
267  DO i=1,ndim
268  CALL coo_to_mat_i(i,bxxy,dumb_mat2)
269  dumb_mat3=matmul(dumb_mat1,matmul(corrint,transpose(dumb_mat2)))
270  CALL add_matc_to_tensor(i,dumb_mat3,utot)
271  ENDDO
272 
273  DO j=1,ndim
274  CALL coo_to_mat_i(j,bxxy,dumb_mat2)
275  dumb_mat3=matmul(dumb_mat1,matmul(corrint,transpose(dumb_mat2)))
276  DO k=1,ndim
277  CALL add_vec_jk_to_tensor(j,k,dumb_mat3(:,k),utot)
278  ENDDO
279  ENDDO
280 
281  !V
282 
283  ALLOCATE(vtot(ndim), stat=allocstat)
284  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
285 
286  DO i=1,ndim
287  DO j=1,ndim
288  CALL coo_to_mat_i(i,bxxy,dumb_mat1)
289  CALL coo_to_mat_i(j,bxxy,dumb_mat2)
290  dumb_mat3=matmul(dumb_mat1,matmul(corrint,transpose(dumb_mat2)))
291  CALL add_matc_to_tensor4(j,i,dumb_mat3,vtot)
292  ENDDO
293  ENDDO
294 
295  DEALLOCATE(dumb_mat1, dumb_mat2, stat=allocstat)
296  IF (allocstat /= 0) stop "*** Problem to deallocate ! ***"
297 
298  DEALLOCATE(dumb_mat3, dumb_mat4, stat=allocstat)
299  IF (allocstat /= 0) stop "*** Problem to deallocate ! ***"
300 
301  DEALLOCATE(dumb_vec, stat=allocstat)
302  IF (allocstat /= 0) stop "*** Problem to deallocate ! ***"
303 
304 
integer ndim
Number of variables (dimension of the model)
Definition: params.f90:85
Statistics accumulators.
Definition: stat.f90:14
subroutine, public add_to_tensor(src, dst)
Routine to add a rank-3 tensor to another one.
Definition: tensor.f90:350
character(len=4) mode
Stochastic mode parameter.

Variable Documentation

type(coolist), dimension(:), allocatable, public mtv_int_tensor::b1

First quadratic tensor.

Definition at line 54 of file MTV_int_tensor.f90.

54  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: b1 !< First quadratic tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist), dimension(:), allocatable, public mtv_int_tensor::b2

Second quadratic tensor.

Definition at line 55 of file MTV_int_tensor.f90.

55  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: b2 !< Second quadratic tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist), dimension(:), allocatable, public mtv_int_tensor::btot

Total quadratic tensor.

Definition at line 56 of file MTV_int_tensor.f90.

56  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: btot !< Total quadratic tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
real(kind=8), dimension(:,:), allocatable mtv_int_tensor::dumb_mat1
private

Dummy matrix.

Definition at line 67 of file MTV_int_tensor.f90.

67  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: dumb_mat1 !< Dummy matrix
real(kind=8), dimension(:,:), allocatable mtv_int_tensor::dumb_mat2
private

Dummy matrix.

Definition at line 68 of file MTV_int_tensor.f90.

68  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: dumb_mat2 !< Dummy matrix
real(kind=8), dimension(:,:), allocatable mtv_int_tensor::dumb_mat3
private

Dummy matrix.

Definition at line 69 of file MTV_int_tensor.f90.

69  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: dumb_mat3 !< Dummy matrix
real(kind=8), dimension(:,:), allocatable mtv_int_tensor::dumb_mat4
private

Dummy matrix.

Definition at line 70 of file MTV_int_tensor.f90.

70  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE :: dumb_mat4 !< Dummy matrix
real(kind=8), dimension(:), allocatable mtv_int_tensor::dumb_vec
private

Dummy vector.

Definition at line 66 of file MTV_int_tensor.f90.

66  REAL(KIND=8), DIMENSION(:), ALLOCATABLE :: dumb_vec !< Dummy vector
real(kind=8), dimension(:), allocatable, public mtv_int_tensor::h1

First constant vector.

Definition at line 42 of file MTV_int_tensor.f90.

42  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: h1 !< First constant vector
real(kind=8), dimension(:), allocatable, public mtv_int_tensor::h2

Second constant vector.

Definition at line 43 of file MTV_int_tensor.f90.

43  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: h2 !< Second constant vector
real(kind=8), dimension(:), allocatable, public mtv_int_tensor::h3

Third constant vector.

Definition at line 44 of file MTV_int_tensor.f90.

44  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: h3 !< Third constant vector
real(kind=8), dimension(:), allocatable, public mtv_int_tensor::htot

Total constant vector.

Definition at line 45 of file MTV_int_tensor.f90.

45  REAL(KIND=8), DIMENSION(:), ALLOCATABLE, PUBLIC :: htot !< Total constant vector
type(coolist), dimension(:), allocatable, public mtv_int_tensor::l1

First linear tensor.

Definition at line 48 of file MTV_int_tensor.f90.

48  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: l1 !< First linear tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist), dimension(:), allocatable, public mtv_int_tensor::l2

Second linear tensor.

Definition at line 49 of file MTV_int_tensor.f90.

49  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: l2 !< Second linear tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist), dimension(:), allocatable, public mtv_int_tensor::l3

Third linear tensor.

Definition at line 50 of file MTV_int_tensor.f90.

50  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: l3 !< Third linear tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist), dimension(:), allocatable, public mtv_int_tensor::ltot

Total linear tensor.

Definition at line 51 of file MTV_int_tensor.f90.

51  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: ltot !< Total linear tensor
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist4), dimension(:), allocatable, public mtv_int_tensor::mtot

Tensor for the cubic terms.

Definition at line 58 of file MTV_int_tensor.f90.

58  TYPE(coolist4), DIMENSION(:), ALLOCATABLE, PUBLIC :: mtot !< Tensor for the cubic terms
4d coordinate list. Type used to represent the rank-4 sparse tensor.
Definition: tensor.f90:44
real(kind=8), dimension(:,:), allocatable, public mtv_int_tensor::q1

Constant terms for the state-dependent noise covariance matrix.

Definition at line 61 of file MTV_int_tensor.f90.

61  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: q1 !< Constant terms for the state-dependent noise covariance matrix
real(kind=8), dimension(:,:), allocatable, public mtv_int_tensor::q2

Constant terms for the state-independent noise covariance matrix.

Definition at line 62 of file MTV_int_tensor.f90.

62  REAL(KIND=8), DIMENSION(:,:), ALLOCATABLE, PUBLIC :: q2 !< Constant terms for the state-independent noise covariance matrix
type(coolist), dimension(:), allocatable, public mtv_int_tensor::utot

Linear terms for the state-dependent noise covariance matrix.

Definition at line 63 of file MTV_int_tensor.f90.

63  TYPE(coolist), DIMENSION(:), ALLOCATABLE, PUBLIC :: utot !< Linear terms for the state-dependent noise covariance matrix
Coordinate list. Type used to represent the sparse tensor.
Definition: tensor.f90:38
type(coolist4), dimension(:), allocatable, public mtv_int_tensor::vtot

Quadratic terms for the state-dependent noise covariance matrix.

Definition at line 64 of file MTV_int_tensor.f90.

64  TYPE(coolist4), DIMENSION(:), ALLOCATABLE, PUBLIC :: vtot !< Quadratic terms for the state-dependent noise covariance matrix
4d coordinate list. Type used to represent the rank-4 sparse tensor.
Definition: tensor.f90:44