A Modular Arbitrary-Order Ocean-Atmosphere Model -- Stochastic implementation
inprod_analytic.f90
Go to the documentation of this file.
1 
2 ! inprod_analytic.f90
3 !
4 !> Inner products between the truncated set of basis functions for the
5 !> ocean and atmosphere streamfunction fields.
6 !> These are partly calculated using the analytical expressions from
7 !> Cehelsky, P., & Tung, K. K. : Theories of multiple equilibria and
8 !> weather regimes-A critical reexamination. Part II: Baroclinic two-layer
9 !> models. Journal of the atmospheric sciences, 44(21), 3282-3303, 1987.
10 !
11 !> @copyright
12 !> 2015 Lesley De Cruz & Jonathan Demaeyer.
13 !> See LICENSE.txt for license information.
14 !
15 !---------------------------------------------------------------------------!
16 !
17 !> @remark
18 !> Generated Fortran90/95 code
19 !> from inprod_analytic.lua
20 !
21 !---------------------------------------------------------------------------!
22 
24 
25  !-----------------------------------------------------!
26  ! !
27  ! Preamble and variables declaration !
28  ! !
29  !-----------------------------------------------------!
30 
31  USE params, only: nbatm, nboc, natm, noc, n, oms, ams, pi
32  USE util, only: isin,piksrt
33  IMPLICIT NONE
34 
35  PRIVATE
36 
37  PUBLIC :: init_inprod
38 
39  !> Atmospheric bloc specification type
40  TYPE :: atm_wavenum
41  CHARACTER :: typ
42  INTEGER :: m=0,p=0,h=0
43  REAL(KIND=8) :: nx=0.,ny=0.
44  END TYPE atm_wavenum
45 
46  !> Oceanic bloc specification type
47  TYPE :: ocean_wavenum
48  INTEGER :: p,h
49  REAL(KIND=8) :: nx,ny
50  END TYPE ocean_wavenum
51 
52  !> Type holding the atmospheric inner products tensors
53  TYPE :: atm_tensors
54  PROCEDURE(calculate_a), POINTER, NOPASS :: a
55  PROCEDURE(calculate_b), POINTER, NOPASS :: b
56  PROCEDURE(calculate_c_atm), POINTER, NOPASS :: c
57  PROCEDURE(calculate_d), POINTER, NOPASS :: d
58  PROCEDURE(calculate_g), POINTER, NOPASS :: g
59  PROCEDURE(calculate_s), POINTER, NOPASS :: s
60  END TYPE atm_tensors
61 
62  !> Type holding the oceanic inner products tensors
63  TYPE :: ocean_tensors
64  PROCEDURE(calculate_k), POINTER, NOPASS :: k
65  PROCEDURE(calculate_m), POINTER, NOPASS :: m
66  PROCEDURE(calculate_c_oc), POINTER, NOPASS :: c
67  PROCEDURE(calculate_n), POINTER, NOPASS :: n
68  PROCEDURE(calculate_o), POINTER, NOPASS :: o
69  PROCEDURE(calculate_w), POINTER, NOPASS :: w
70  END TYPE ocean_tensors
71 
72  !> Atmospheric blocs specification
73  TYPE(atm_wavenum), DIMENSION(:), ALLOCATABLE, PUBLIC :: awavenum
74  !> Oceanic blocs specification
75  TYPE(ocean_wavenum), DIMENSION(:), ALLOCATABLE, PUBLIC :: owavenum
76 
77  !> Atmospheric tensors
78  TYPE(atm_tensors), PUBLIC :: atmos
79  !> Oceanic tensors
80  TYPE(ocean_tensors), PUBLIC :: ocean
81 
82 
83  !-----------------------------------------------------!
84  ! !
85  ! End of preamble !
86  ! !
87  !-----------------------------------------------------!
88 
89 CONTAINS
90 
91  !-----------------------------------------------------!
92  ! !
93  ! Definition of the Helper functions from Cehelsky !
94  ! & Tung !
95  ! !
96  !-----------------------------------------------------!
97 
98  !> Cehelsky & Tung Helper functions
99  REAL(KIND=8) FUNCTION b1(Pi, Pj, Pk)
100  INTEGER :: Pi,Pj,Pk
101  b1 = (pk + pj) / REAL(pi)
102  END FUNCTION b1
103 
104  !> Cehelsky & Tung Helper functions
105  REAL(KIND=8) FUNCTION b2(Pi, Pj, Pk)
106  INTEGER :: Pi,Pj,Pk
107  b2 = (pk - pj) / REAL(pi)
108  END FUNCTION b2
109 
110  !> Integer Dirac delta function
111  REAL(KIND=8) FUNCTION delta(r)
112  INTEGER :: r
113  IF (r==0) THEN
114  delta = 1.d0
115  ELSE
116  delta = 0.d0
117  ENDIF
118  END FUNCTION delta
119 
120  !> "Odd or even" function
121  REAL(KIND=8) FUNCTION flambda(r)
122  INTEGER :: r
123  IF (mod(r,2)==0) THEN
124  flambda = 0.d0
125  ELSE
126  flambda = 1.d0
127  ENDIF
128  END FUNCTION flambda
129 
130  !> Cehelsky & Tung Helper functions
131  REAL(KIND=8) FUNCTION s1(Pj, Pk, Mj, Hk)
132  INTEGER :: Pk,Pj,Mj,Hk
133  s1 = -((pk * mj + pj * hk)) / 2.d0
134  END FUNCTION s1
135 
136  !> Cehelsky & Tung Helper functions
137  REAL(KIND=8) FUNCTION s2(Pj, Pk, Mj, Hk)
138  INTEGER :: Pk,Pj,Mj,Hk
139  s2 = (pk * mj - pj * hk) / 2.d0
140  END FUNCTION s2
141 
142  !> Cehelsky & Tung Helper functions
143  REAL(KIND=8) FUNCTION s3(Pj, Pk, Hj, Hk)
144  INTEGER :: Pj,Pk,Hj,Hk
145  s3 = (pk * hj + pj * hk) / 2.d0
146  END FUNCTION s3
147 
148  !> Cehelsky & Tung Helper functions
149  REAL(KIND=8) FUNCTION s4(Pj, Pk, Hj, Hk)
150  INTEGER :: Pj,Pk,Hj,Hk
151  s4 = (pk * hj - pj * hk) / 2.d0
152  END FUNCTION s4
153 
154  !-----------------------------------------------------!
155  ! Inner products definition routines !
156  !--------------------------------------------------------!
157  ! 1. Inner products in the equations for the atmosphere !
158  !--------------------------------------------------------!
159 
160  !> Eigenvalues of the Laplacian (atmospheric)
161  !>
162  !> \f$ a_{i,j} = (F_i, \nabla^2 F_j)\f$ .
163  REAL(KIND=8) FUNCTION calculate_a(i,j)
164  INTEGER, INTENT(IN) :: i,j
165  TYPE(atm_wavenum) :: Ti
166 
167  calculate_a = 0.d0
168  IF (i==j) THEN
169  ti = awavenum(i)
170  calculate_a = -(n**2) * ti%Nx**2 - ti%Ny**2
171  END IF
172  END FUNCTION calculate_a
173 
174  !> Streamfunction advection terms (atmospheric)
175  !>
176  !> \f$ b_{i,j,k} = (F_i, J(F_j, \nabla^2 F_k))\f$ .
177  REAL(KIND=8) FUNCTION calculate_b(i,j,k)
178  INTEGER, INTENT(IN) :: i,j,k
179 
180  calculate_b = calculate_a(k,k) * calculate_g(i,j,k)
181 
182  END FUNCTION calculate_b
183 
184  !> Beta term for the atmosphere
185  !>
186  !> \f$ c_{i,j} = (F_i, \partial_x F_j)\f$ .
187  REAL(KIND=8) FUNCTION calculate_c_atm(i,j)
188  INTEGER, INTENT(IN) :: i,j
189  TYPE(atm_wavenum) :: Ti, Tj
190 
191  ti = awavenum(i)
192  tj = awavenum(j)
193  calculate_c_atm = 0.d0
194  IF ((ti%typ == "K") .AND. (tj%typ == "L")) THEN
195  calculate_c_atm = n * ti%M * delta(ti%M - tj%H) * delta(ti%P - tj%P)
196  ELSE IF ((ti%typ == "L") .AND. (tj%typ == "K")) THEN
197  ti = awavenum(j)
198  tj = awavenum(i)
199  calculate_c_atm = - n * ti%M * delta(ti%M - tj%H) * delta(ti%P - tj%P)
200  END IF
201 
202  END FUNCTION calculate_c_atm
203 
204  !> Forcing of the ocean on the atmosphere.
205  !>
206  !> \f$ d_{i,j} = (F_i, \nabla^2 \eta_j)\f$ .
207  REAL(KIND=8) FUNCTION calculate_d(i,j)
208  INTEGER, INTENT(IN) :: i,j
209 
211 
212  END FUNCTION calculate_d
213 
214  !> Temperature advection terms (atmospheric)
215  !>
216  !> \f$ g_{i,j,k} = (F_i, J(F_j, F_k))\f$ .
217  REAL(KIND=8) FUNCTION calculate_g(i,j,k)
218  INTEGER, INTENT(IN) :: i,j,k
219  TYPE(atm_wavenum) :: Ti,Tj,Tk
220  REAL(KIND=8) :: val,vb1, vb2, vs1, vs2, vs3, vs4
221  INTEGER, DIMENSION(3) :: a,b
222  INTEGER, DIMENSION(3,3) :: w
223  CHARACTER, DIMENSION(3) :: s
224  INTEGER :: par
225 
226  ti = awavenum(i)
227  tj = awavenum(j)
228  tk = awavenum(k)
229 
230  a(1)=i
231  a(2)=j
232  a(3)=k
233 
234  val=0.d0
235 
236  IF ((ti%typ == "L") .AND. (tj%typ == "L") .AND. (tk%typ == "L")) THEN
237 
238  CALL piksrt(3,a,par)
239 
240  ti = awavenum(a(1))
241  tj = awavenum(a(2))
242  tk = awavenum(a(3))
243 
244  vs3 = s3(tj%P,tk%P,tj%H,tk%H)
245  vs4 = s4(tj%P,tk%P,tj%H,tk%H)
246  val = vs3 * ((delta(tk%H - tj%H - ti%H) - delta(tk%H &
247  &- tj%H + ti%H)) * delta(tk%P + tj%P - ti%P) +&
248  & delta(tk%H + tj%H - ti%H) * (delta(tk%P - tj%P&
249  & + ti%P) - delta(tk%P - tj%P - ti%P))) + vs4 *&
250  & ((delta(tk%H + tj%H - ti%H) * delta(tk%P - tj&
251  &%P - ti%P)) + (delta(tk%H - tj%H + ti%H) -&
252  & delta(tk%H - tj%H - ti%H)) * (delta(tk%P - tj&
253  &%P - ti%P) - delta(tk%P - tj%P + ti%P)))
254  ELSE
255 
256  s(1)=ti%typ
257  s(2)=tj%typ
258  s(3)=tk%typ
259 
260  w(1,:)=isin("A",s)
261  w(2,:)=isin("K",s)
262  w(3,:)=isin("L",s)
263 
264  IF (any(w(1,:)/=0) .AND. any(w(2,:)/=0) .AND. any(w(3,:)/=0)) THEN
265  b=w(:,1)
266  ti = awavenum(a(b(1)))
267  tj = awavenum(a(b(2)))
268  tk = awavenum(a(b(3)))
269  call piksrt(3,b,par)
270  vb1 = b1(ti%P,tj%P,tk%P)
271  vb2 = b2(ti%P,tj%P,tk%P)
272  val = -2 * sqrt(2.) / pi * tj%M * delta(tj%M - tk%H) * flambda(ti%P + tj%P + tk%P)
273  IF (val /= 0.d0) val = val * (vb1**2 / (vb1**2 - 1) - vb2**2 / (vb2**2 - 1))
274  ELSEIF ((w(2,2)/=0) .AND. (w(2,3)==0) .AND. any(w(3,:)/=0)) THEN
275  ti = awavenum(a(w(2,1)))
276  tj = awavenum(a(w(2,2)))
277  tk = awavenum(a(w(3,1)))
278  b(1)=w(2,1)
279  b(2)=w(2,2)
280  b(3)=w(3,1)
281  call piksrt(3,b,par)
282  vs1 = s1(tj%P,tk%P,tj%M,tk%H)
283  vs2 = s2(tj%P,tk%P,tj%M,tk%H)
284  val = vs1 * (delta(ti%M - tk%H - tj%M) * delta(ti%P -&
285  & tk%P + tj%P) - delta(ti%M- tk%H - tj%M) *&
286  & delta(ti%P + tk%P - tj%P) + (delta(tk%H - tj%M&
287  & + ti%M) + delta(tk%H - tj%M - ti%M)) *&
288  & delta(tk%P + tj%P - ti%P)) + vs2 * (delta(ti%M&
289  & - tk%H - tj%M) * delta(ti%P - tk%P - tj%P) +&
290  & (delta(tk%H - tj%M - ti%M) + delta(ti%M + tk%H&
291  & - tj%M)) * (delta(ti%P - tk%P + tj%P) -&
292  & delta(tk%P - tj%P + ti%P)))
293  ENDIF
294  ENDIF
295  calculate_g=par*val*n
296 
297  END FUNCTION calculate_g
298 
299  !> Forcing (thermal) of the ocean on the atmosphere.
300  !>
301  !> \f$ s_{i,j} = (F_i, \eta_j)\f$ .
302  REAL(KIND=8) FUNCTION calculate_s(i,j)
303  INTEGER, INTENT(IN) :: i,j
304  TYPE(atm_wavenum) :: Ti
305  TYPE(ocean_wavenum) :: Dj
306  REAL(KIND=8) :: val
307 
308  ti = awavenum(i)
309  dj = owavenum(j)
310  val=0.d0
311  IF (ti%typ == "A") THEN
312  val = flambda(dj%H) * flambda(dj%P + ti%P)
313  IF (val /= 0.d0) THEN
314  val = val*8*sqrt(2.)*dj%P/(pi**2 * (dj%P**2 - ti%P**2) * dj%H)
315  END IF
316  ELSEIF (ti%typ == "K") THEN
317  val = flambda(2 * ti%M + dj%H) * delta(dj%P - ti%P)
318  IF (val /= 0.d0) THEN
319  val = val*4*dj%H/(pi * (-4 * ti%M**2 + dj%H**2))
320  END IF
321  ELSEIF (ti%typ == "L") THEN
322  val = delta(dj%P - ti%P) * delta(2 * ti%H - dj%H)
323  END IF
324  calculate_s=val
325 
326  END FUNCTION calculate_s
327 
328  !--------------------------------------------------------!
329  ! 2. Inner products in the equations for the ocean !
330  !--------------------------------------------------------!
331 
332  !> Forcing of the atmosphere on the ocean.
333  !>
334  !> \f$ K_{i,j} = (\eta_i, \nabla^2 F_j)\f$ .
335  REAL(KIND=8) FUNCTION calculate_k(i,j)
336  INTEGER, INTENT(IN) :: i,j
337 
338  calculate_k = calculate_s(j,i) * calculate_a(j,j)
339  END FUNCTION calculate_k
340 
341  !> Forcing of the ocean fields on the ocean.
342  !>
343  !> \f$ M_{i,j} = (eta_i, \nabla^2 \eta_j)\f$ .
344  REAL(KIND=8) FUNCTION calculate_m(i,j)
345  INTEGER, INTENT(IN) :: i,j
346  TYPE(ocean_wavenum) :: Di
347 
348  calculate_m=0.d0
349  IF (i==j) THEN
350  di = owavenum(i)
351  calculate_m = -(n**2) * di%Nx**2 - di%Ny**2
352  END IF
353  END FUNCTION calculate_m
354 
355  !> Beta term for the ocean
356  !>
357  !> \f$ N_{i,j} = (\eta_i, \partial_x \eta_j) \f$.
358  REAL(KIND=8) FUNCTION calculate_n(i,j)
359  INTEGER, INTENT(IN) :: i,j
360  TYPE(ocean_wavenum) :: Di,Dj
361  REAL(KIND=8) :: val
362 
363  di = owavenum(i)
364  dj = owavenum(j)
365  calculate_n = 0.d0
366  IF (dj%H/=di%H) THEN
367  val = delta(di%P - dj%P) * flambda(di%H + dj%H)
368  calculate_n = val * (-2) * dj%H * di%H * n / ((dj%H**2 - di%H**2) * pi)
369  ENDIF
370 
371  END FUNCTION calculate_n
372 
373  !> Temperature advection term (passive scalar)
374  !>
375  !> \f$ O_{i,j,k} = (\eta_i, J(\eta_j, \eta_k))\f$ .
376  REAL(KIND=8) FUNCTION calculate_o(i,j,k)
377  INTEGER, INTENT(IN) :: i,j,k
378  TYPE(ocean_wavenum) :: Di,Dj,Dk
379  REAL(KIND=8) :: vs3,vs4,val
380  INTEGER, DIMENSION(3) :: a
381  INTEGER :: par
382 
383  val=0.d0
384 
385  a(1)=i
386  a(2)=j
387  a(3)=k
388 
389  CALL piksrt(3,a,par)
390 
391  di = owavenum(a(1))
392  dj = owavenum(a(2))
393  dk = owavenum(a(3))
394 
395  vs3 = s3(dj%P,dk%P,dj%H,dk%H)
396  vs4 = s4(dj%P,dk%P,dj%H,dk%H)
397  val = vs3*((delta(dk%H - dj%H - di%H) - delta(dk%H - dj&
398  &%H + di%H)) * delta(dk%P + dj%P - di%P) + delta(dk&
399  &%H + dj%H - di%H) * (delta(dk%P - dj%P + di%P) -&
400  & delta(dk%P - dj%P - di%P))) + vs4 * ((delta(dk%H &
401  &+ dj%H - di%H) * delta(dk%P - dj%P - di%P)) +&
402  & (delta(dk%H - dj%H + di%H) - delta(dk%H - dj%H -&
403  & di%H)) * (delta(dk%P - dj%P - di%P) - delta(dk%P &
404  &- dj%P + di%P)))
405  calculate_o = par * val * n / 2
406  END FUNCTION calculate_o
407 
408  !> Streamfunction advection terms (oceanic)
409  !>
410  !> \f$ C_{i,j,k} = (\eta_i, J(\eta_j,\nabla^2 \eta_k))\f$ .
411  REAL(KIND=8) FUNCTION calculate_c_oc(i,j,k)
412  INTEGER, INTENT(IN) :: i,j,k
413 
414  calculate_c_oc = calculate_m(k,k) * calculate_o(i,j,k)
415 
416  END FUNCTION calculate_c_oc
417 
418  !> Short-wave radiative forcing of the ocean
419  !>
420  !> \f$ W_{i,j} = (\eta_i, F_j)\f$ .
421  REAL(KIND=8) FUNCTION calculate_w(i,j)
422  INTEGER, INTENT(IN) :: i,j
423 
424  calculate_w = calculate_s(j,i)
425 
426  END FUNCTION calculate_w
427 
428  !-----------------------------------------------------!
429  ! !
430  ! Initialisation routine !
431  ! !
432  !-----------------------------------------------------!
433 
434  !> Initialisation of the inner product
435  SUBROUTINE init_inprod
436  INTEGER :: i,j
437  INTEGER :: AllocStat
438 
439  IF (natm == 0 ) THEN
440  stop "*** Problem : natm==0 ! ***"
441  ELSEIF (noc == 0) then
442  stop "*** Problem : noc==0 ! ***"
443  END IF
444 
445 
446  ! Definition of the types and wave numbers tables
447 
448  ALLOCATE(owavenum(noc),awavenum(natm), stat=allocstat)
449  IF (allocstat /= 0) stop "*** Not enough memory ! ***"
450 
451  j=0
452  DO i=1,nbatm
453  IF (ams(i,1)==1) THEN
454  awavenum(j+1)%typ='A'
455  awavenum(j+2)%typ='K'
456  awavenum(j+3)%typ='L'
457 
458  awavenum(j+1)%P=ams(i,2)
459  awavenum(j+2)%M=ams(i,1)
460  awavenum(j+2)%P=ams(i,2)
461  awavenum(j+3)%H=ams(i,1)
462  awavenum(j+3)%P=ams(i,2)
463 
464  awavenum(j+1)%Ny=REAL(ams(i,2))
465  awavenum(j+2)%Nx=REAL(ams(i,1))
466  awavenum(j+2)%Ny=REAL(ams(i,2))
467  awavenum(j+3)%Nx=REAL(ams(i,1))
468  awavenum(j+3)%Ny=REAL(ams(i,2))
469 
470  j=j+3
471  ELSE
472  awavenum(j+1)%typ='K'
473  awavenum(j+2)%typ='L'
474 
475  awavenum(j+1)%M=ams(i,1)
476  awavenum(j+1)%P=ams(i,2)
477  awavenum(j+2)%H=ams(i,1)
478  awavenum(j+2)%P=ams(i,2)
479 
480  awavenum(j+1)%Nx=REAL(ams(i,1))
481  awavenum(j+1)%Ny=REAL(ams(i,2))
482  awavenum(j+2)%Nx=REAL(ams(i,1))
483  awavenum(j+2)%Ny=REAL(ams(i,2))
484 
485  j=j+2
486 
487  ENDIF
488  ENDDO
489 
490  DO i=1,noc
491  owavenum(i)%H=oms(i,1)
492  owavenum(i)%P=oms(i,2)
493 
494  owavenum(i)%Nx=oms(i,1)/2.d0
495  owavenum(i)%Ny=oms(i,2)
496 
497  ENDDO
498 
499  ! Pointing to the atmospheric inner products functions
500 
501  atmos%a => calculate_a
502  atmos%g => calculate_g
503  atmos%s => calculate_s
504  atmos%b => calculate_b
505  atmos%d => calculate_d
507 
508  ! Pointing to the oceanic inner products functions
509 
510  ocean%M => calculate_m
511  ocean%N => calculate_n
512  ocean%O => calculate_o
513  ocean%C => calculate_c_oc
514  ocean%W => calculate_w
515  ocean%K => calculate_k
516 
517  END SUBROUTINE init_inprod
518 
519 
520 END MODULE inprod_analytic
521 
integer noc
Number of oceanic basis functions.
Definition: params.f90:84
Utility module.
Definition: util.f90:12
real(kind=8) function calculate_g(i, j, k)
Temperature advection terms (atmospheric)
Atmospheric bloc specification type.
subroutine, public init_inprod
Initialisation of the inner product.
integer nbatm
Number of oceanic blocks.
Definition: params.f90:82
real(kind=8) function calculate_d(i, j)
Forcing of the ocean on the atmosphere.
Inner products between the truncated set of basis functions for the ocean and atmosphere streamfuncti...
Statistics accumulators.
Definition: stat.f90:14
subroutine, public piksrt(k, arr, par)
Simple card player sorting function.
Definition: util.f90:118
real(kind=8) function b1(Pi, Pj, Pk)
Cehelsky & Tung Helper functions.
real(kind=8) function s4(Pj, Pk, Hj, Hk)
Cehelsky & Tung Helper functions.
real(kind=8) pi
Definition: params.f90:48
real(kind=8) function s2(Pj, Pk, Mj, Hk)
Cehelsky & Tung Helper functions.
real(kind=8) function b2(Pi, Pj, Pk)
Cehelsky & Tung Helper functions.
real(kind=8) function delta(r)
Integer Dirac delta function.
real(kind=8) function calculate_c_oc(i, j, k)
Streamfunction advection terms (oceanic)
real(kind=8) function calculate_o(i, j, k)
Temperature advection term (passive scalar)
integer function, dimension(size(s)), public isin(c, s)
Determine if a character is in a string and where.
Definition: util.f90:47
type(ocean_wavenum), dimension(:), allocatable, public owavenum
Oceanic blocs specification.
Type holding the oceanic inner products tensors.
real(kind=8) function calculate_s(i, j)
Forcing (thermal) of the ocean on the atmosphere.
real(kind=8) function calculate_k(i, j)
Forcing of the atmosphere on the ocean.
type(atm_wavenum), dimension(:), allocatable, public awavenum
Atmospheric blocs specification.
Oceanic bloc specification type.
real(kind=8) function s1(Pj, Pk, Mj, Hk)
Cehelsky & Tung Helper functions.
real(kind=8) function s3(Pj, Pk, Hj, Hk)
Cehelsky & Tung Helper functions.
integer natm
Number of atmospheric basis functions.
Definition: params.f90:83
type(atm_tensors), public atmos
Atmospheric tensors.
real(kind=8) function calculate_m(i, j)
Forcing of the ocean fields on the ocean.
real(kind=8) function calculate_w(i, j)
Short-wave radiative forcing of the ocean.
integer nboc
Number of atmospheric blocks.
Definition: params.f90:81
real(kind=8) function flambda(r)
"Odd or even" function
The model parameters module.
Definition: params.f90:18
real(kind=8) function calculate_b(i, j, k)
Streamfunction advection terms (atmospheric)
real(kind=8) function calculate_c_atm(i, j)
Beta term for the atmosphere.
Type holding the atmospheric inner products tensors.
integer, dimension(:,:), allocatable ams
Atmospheric mode selection array.
Definition: params.f90:87
real(kind=8) function calculate_n(i, j)
Beta term for the ocean.
type(ocean_tensors), public ocean
Oceanic tensors.
integer, dimension(:,:), allocatable oms
Ocean mode selection array.
Definition: params.f90:86
real(kind=8) function calculate_a(i, j)
Eigenvalues of the Laplacian (atmospheric)
real(kind=8) n
- Aspect ratio
Definition: params.f90:24