1
+ device = Oceananigans. Architectures. device
2
+
3
+ @inline convert_to_0_360 (x) = ((x % 360 ) + 360 ) % 360
4
+
1
5
""" a structure to represent a tripolar grid on a spherical shell """
2
6
struct Tripolar{N, F, S}
3
7
north_poles_latitude :: N
@@ -56,21 +60,21 @@ The north singularities are located at
56
60
`i = 1, j = Nφ` and `i = Nλ ÷ 2 + 1, j = Nλ`
57
61
"""
58
62
function TripolarGrid (arch = CPU (), FT:: DataType = Float64;
59
- size,
63
+ size,
60
64
southernmost_latitude = - 80 , # The southermost `Center` latitude of the grid
61
- halo = (4 , 4 , 4 ),
62
- radius = R_Earth,
65
+ halo = (4 , 4 , 4 ),
66
+ radius = R_Earth,
63
67
z = (0 , 1 ),
64
68
north_poles_latitude = 55 ,
65
69
first_pole_longitude = 70 ) # The second pole is at `λ = first_pole_longitude + 180ᵒ`
66
70
67
- # TODO : change a couple of allocations here and there to be able
71
+ # TODO : change a couple of allocations here and there to be able
68
72
# to construct the grid on the GPU. This is not a huge problem as
69
73
# grid generation is quite fast, but it might become for sub-kilometer grids
70
74
71
75
latitude = (southernmost_latitude, 90 )
72
76
longitude = (- 180 , 180 )
73
-
77
+
74
78
focal_distance = tand ((90 - north_poles_latitude) / 2 )
75
79
76
80
Nλ, Nφ, Nz = size
@@ -106,7 +110,7 @@ function TripolarGrid(arch = CPU(), FT::DataType = Float64;
106
110
φCC = zeros (Nλ, Nφ)
107
111
108
112
loop! = _compute_tripolar_coordinates! (device (CPU ()), (16 , 16 ), (Nλ, Nφ))
109
-
113
+
110
114
loop! (λFF, φFF, λFC, φFC, λCF, φCF, λCC, φCC,
111
115
λᶠᵃᵃ, λᶜᵃᵃ, φᵃᶠᵃ, φᵃᶜᵃ,
112
116
first_pole_longitude,
@@ -127,7 +131,7 @@ function TripolarGrid(arch = CPU(), FT::DataType = Float64;
127
131
128
132
Nx = Nλ
129
133
Ny = Nφ
130
-
134
+
131
135
# return λFF, φFF, λFC, φFC, λCF, φCF, λCC, φCC
132
136
# Helper grid to fill halo
133
137
grid = RectilinearGrid (; size = (Nx, Ny), halo = (Hλ, Hφ), topology = (Periodic, RightConnected, Flat), x = (0 , 1 ), y = (0 , 1 ))
@@ -148,7 +152,7 @@ function TripolarGrid(arch = CPU(), FT::DataType = Float64;
148
152
149
153
lFC = Field ((Face, Center, Center), grid; boundary_conditions = default_boundary_conditions)
150
154
pFC = Field ((Face, Center, Center), grid; boundary_conditions = default_boundary_conditions)
151
-
155
+
152
156
lCF = Field ((Center, Face, Center), grid; boundary_conditions = default_boundary_conditions)
153
157
pCF = Field ((Center, Face, Center), grid; boundary_conditions = default_boundary_conditions)
154
158
0 commit comments