-
Notifications
You must be signed in to change notification settings - Fork 186
/
closure_operators.jl
150 lines (112 loc) · 5.67 KB
/
closure_operators.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
using Oceananigans.Operators
"""
ν_σᶜᶜᶜ(i, j, k, grid, ν, σᶜᶜᶜ, u, v, w)
Multiply the array `ν` located at `ᶜᶜᶜ` by a function
`σᶜᶜᶜ(i, j, k, grid, u, v, w)`
at index `i, j, k` and location `ᶜᶜᶜ`.
"""
@inline ν_σᶜᶜᶜ(i, j, k, grid, ν::TN, σᶜᶜᶜ::TS, u, v, w) where {TN<:AbstractArray, TS} =
@inbounds ν[i, j, k] * σᶜᶜᶜ(i, j, k, grid, u, v, w)
"""
ν_σᶠᶠᶜ(i, j, k, grid, ν, σᶠᶠᶜ, u, v, w)
Multiply the array `ν` located at `ᶜᶜᶜ` by a function
`σᶠᶠᶜ(i, j, k, grid, u, v, w)`
at index `i, j, k` and location `ᶠᶠᶜ`.
"""
@inline ν_σᶠᶠᶜ(i, j, k, grid, ν::TN, σᶠᶠᶜ::TS, u, v, w) where {TN<:AbstractArray, TS} =
@inbounds ℑxyᶠᶠᵃ(i, j, k, grid, ν) * σᶠᶠᶜ(i, j, k, grid, u, v, w)
# These functions are analogous to the two above, but for different locations:
@inline ν_σᶠᶜᶠ(i, j, k, grid, ν::TN, σᶠᶜᶠ::TS, u, v, w) where {TN<:AbstractArray, TS} =
@inbounds ℑxzᶠᵃᶠ(i, j, k, grid, ν) * σᶠᶜᶠ(i, j, k, grid, u, v, w)
@inline ν_σᶜᶠᶠ(i, j, k, grid, ν::TN, σᶜᶠᶠ::TS, u, v, w) where {TN<:AbstractArray, TS} =
@inbounds ℑyzᵃᶠᶠ(i, j, k, grid, ν) * σᶜᶠᶠ(i, j, k, grid, u, v, w)
#####
##### Stress divergences
#####
# At fcc
@inline ∂x_2ν_Σ₁₁(i, j, k, grid, closure, U, diffusivities) =
2 * ∂xᶠᵃᵃ(i, j, k, grid, ν_σᶜᶜᶜ, diffusivities.νₑ, Σ₁₁, U.u, U.v, U.w)
@inline ∂y_2ν_Σ₁₂(i, j, k, grid, closure, U, diffusivities) =
2 * ∂yᵃᶜᵃ(i, j, k, grid, ν_σᶠᶠᶜ, diffusivities.νₑ, Σ₁₂, U.u, U.v, U.w)
@inline ∂z_2ν_Σ₁₃(i, j, k, grid, closure, U, diffusivities) =
2 * ∂zᵃᵃᶜ(i, j, k, grid, ν_σᶠᶜᶠ, diffusivities.νₑ, Σ₁₃, U.u, U.v, U.w)
# At cfc
@inline ∂x_2ν_Σ₂₁(i, j, k, grid, closure, U, diffusivities) =
2 * ∂xᶜᵃᵃ(i, j, k, grid, ν_σᶠᶠᶜ, diffusivities.νₑ, Σ₂₁, U.u, U.v, U.w)
@inline ∂y_2ν_Σ₂₂(i, j, k, grid, closure, U, diffusivities) =
2 * ∂yᵃᶠᵃ(i, j, k, grid, ν_σᶜᶜᶜ, diffusivities.νₑ, Σ₂₂, U.u, U.v, U.w)
@inline ∂z_2ν_Σ₂₃(i, j, k, grid, closure, U, diffusivities) =
2 * ∂zᵃᵃᶜ(i, j, k, grid, ν_σᶜᶠᶠ, diffusivities.νₑ, Σ₂₃, U.u, U.v, U.w)
# At ccf
@inline ∂x_2ν_Σ₃₁(i, j, k, grid, closure, U, diffusivities) =
2 * ∂xᶜᵃᵃ(i, j, k, grid, ν_σᶠᶜᶠ, diffusivities.νₑ, Σ₃₁, U.u, U.v, U.w)
@inline ∂y_2ν_Σ₃₂(i, j, k, grid, closure, U, diffusivities) =
2 * ∂yᵃᶜᵃ(i, j, k, grid, ν_σᶜᶠᶠ, diffusivities.νₑ, Σ₃₂, U.u, U.v, U.w)
@inline ∂z_2ν_Σ₃₃(i, j, k, grid, closure, U, diffusivities) =
2 * ∂zᵃᵃᶠ(i, j, k, grid, ν_σᶜᶜᶜ, diffusivities.νₑ, Σ₃₃, U.u, U.v, U.w)
"""
κ_∂x_c(i, j, k, grid, c, κ, closure, args...)
Return `κ ∂x c`, where `κ` is an array or function that computes
diffusivity at cell centers (location `ccc`), and `c` is an array of scalar
data located at cell centers.
"""
@inline function κ_∂x_c(i, j, k, grid, κ, c, closure, args...)
κ = ℑxᶠᵃᵃ(i, j, k, grid, κ, closure, args...)
∂x_c = ∂xᶠᵃᵃ(i, j, k, grid, c)
return κ * ∂x_c
end
"""
κ_∂y_c(i, j, k, grid, c, κ, closure, args...)
Return `κ ∂y c`, where `κ` is an array or function that computes
diffusivity at cell centers (location `ccc`), and `c` is an array of scalar
data located at cell centers.
"""
@inline function κ_∂y_c(i, j, k, grid, κ, c, closure, args...)
κ = ℑyᵃᶠᵃ(i, j, k, grid, κ, closure, args...)
∂y_c = ∂yᵃᶠᵃ(i, j, k, grid, c)
return κ * ∂y_c
end
"""
κ_∂z_c(i, j, k, grid, c, κ, closure, buoyancy, u, v, w, T, S)
Return `κ ∂z c`, where `κ` is an array or function that computes
diffusivity at cell centers (location `ccc`), and `c` is an array of scalar
data located at cell centers.
"""
@inline function κ_∂z_c(i, j, k, grid, κ, c, closure, args...)
κ = ℑzᵃᵃᶠ(i, j, k, grid, κ, closure, args...)
∂z_c = ∂zᵃᵃᶠ(i, j, k, grid, c)
return κ * ∂z_c
end
"""
∂ⱼ_2ν_Σ₁ⱼ(i, j, k, grid, closure, U, diffusivities)
Return the ``x``-component of the turbulent diffusive flux divergence:
`∂x(2 ν Σ₁₁) + ∂y(2 ν Σ₁₁) + ∂z(2 ν Σ₁₁)`
at the location `fcc`.
"""
@inline ∂ⱼ_2ν_Σ₁ⱼ(i, j, k, grid, closure::IsotropicViscosity, U, diffusivities) = (
∂x_2ν_Σ₁₁(i, j, k, grid, closure, U, diffusivities)
+ ∂y_2ν_Σ₁₂(i, j, k, grid, closure, U, diffusivities)
+ ∂z_2ν_Σ₁₃(i, j, k, grid, closure, U, diffusivities)
)
"""
∂ⱼ_2ν_Σ₂ⱼ(i, j, k, grid, closure, U, diffusivities)
Return the ``y``-component of the turbulent diffusive flux divergence:
`∂x(2 ν Σ₂₁) + ∂y(2 ν Σ₂₂) + ∂z(2 ν Σ₂₂)`
at the location `ccf`.
"""
@inline ∂ⱼ_2ν_Σ₂ⱼ(i, j, k, grid, closure::IsotropicViscosity, U, diffusivities) = (
∂x_2ν_Σ₂₁(i, j, k, grid, closure, U, diffusivities)
+ ∂y_2ν_Σ₂₂(i, j, k, grid, closure, U, diffusivities)
+ ∂z_2ν_Σ₂₃(i, j, k, grid, closure, U, diffusivities)
)
"""
∂ⱼ_2ν_Σ₃ⱼ(i, j, k, grid, closure, diffusivities)
Return the ``z``-component of the turbulent diffusive flux divergence:
`∂x(2 ν Σ₃₁) + ∂y(2 ν Σ₃₂) + ∂z(2 ν Σ₃₃)`
at the location `ccf`.
"""
@inline ∂ⱼ_2ν_Σ₃ⱼ(i, j, k, grid, closure::IsotropicViscosity, U, diffusivities) = (
∂x_2ν_Σ₃₁(i, j, k, grid, closure, U, diffusivities)
+ ∂y_2ν_Σ₃₂(i, j, k, grid, closure, U, diffusivities)
+ ∂z_2ν_Σ₃₃(i, j, k, grid, closure, U, diffusivities)
)