Skip to content

Commit a8634a5

Browse files
knutfrodeknutfrode
authored
Merge pull request #1760 from knutfrode/dev
[run-ex] If wave direction is needed but not provided, it is now set to wind d…
2 parents c7b9a0d + 7494e61 commit a8634a5

3 files changed

Lines changed: 90 additions & 9 deletions

File tree

‎opendrift/models/openberg.py‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -133,6 +133,8 @@ def wind_force(iceb_vel, wind_vel, Ava, Aha, wind_form_drag_coef, wind_skin_drag
133133
wind_form_drag_coef : Ca is the air form drag coefficient
134134
wind_skin_drag_coef : Cda is the air skin drag coefficient
135135
"""
136+
137+
# TODO: wind_force seems to be more than 10 times smaller than wave_force. This should be checked.
136138
vxa, vya = wind_vel[0], wind_vel[1]
137139
x_vel, y_vel = iceb_vel[0], iceb_vel[1]
138140
rel_wind_x_vel = vxa - x_vel
@@ -452,7 +454,7 @@ def advect_iceberg(self):
452454
Ai = sea_ice_thickness * length
453455

454456
mass = width * (Ava + Avo) * rho_iceb * weight_coef
455-
k = (rho_air * wind_form_drag_coef * Ava / (rho_water * water_form_drag_coef * Avo))
457+
k = (rho_air * wind_form_drag_coef * Ava) / (rho_water * water_form_drag_coef * Avo)
456458
f = np.sqrt(k) / (1 + np.sqrt(k)) # (f is the wind drift factor, only used in the no acceleration model)
457459

458460
wave_rad = self.get_config('drift:wave_rad')

‎opendrift/models/physics_methods.py‎

Lines changed: 11 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -814,6 +814,7 @@ def resurface_elements(self, minimum_depth):
814814

815815
def calculate_missing_environment_variables(self):
816816

817+
# TODO: we need a better mechanism to detect missing variables
817818
# Missing significant wave height
818819
if hasattr(self.environment,
819820
'sea_surface_wave_significant_height') and \
@@ -822,7 +823,16 @@ def calculate_missing_environment_variables(self):
822823
logger.debug('Calculating Hs from wind, min: %f, mean: %f, max: %f' %
823824
(Hs.min(), Hs.mean(), Hs.max()))
824825

825-
# Missing wave periode
826+
# Missing wave direction related to previously calculated significant wave height.
827+
# Set equal to wind direction. (Andrea Gierisch)
828+
if hasattr(self.environment, 'sea_surface_wave_from_direction') and \
829+
self.environment.sea_surface_wave_from_direction.max() == 0:
830+
wave_direction = np.rad2deg(np.arctan2(self.environment.x_wind, self.environment.y_wind))
831+
self.environment.sea_surface_wave_from_direction = wave_direction
832+
logger.warning('Setting wave direction equal to wind direction, min: %f, mean: %f, max: %f' %
833+
(wave_direction.min(), wave_direction.mean(), wave_direction.max()))
834+
835+
# Missing wave period
826836
if hasattr(self.environment,
827837
'sea_surface_wave_mean_period_from_variance_spectral_density_second_frequency_moment') and \
828838
self.environment.sea_surface_wave_mean_period_from_variance_spectral_density_second_frequency_moment.max() == 0:

‎tests/models/test_openberg.py‎

Lines changed: 76 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -40,8 +40,8 @@ def test_openberg_constant_forcing():
4040
o.set_config('drift:coriolis', False)
4141
o.seed_elements(4, 60, time=datetime.now())
4242
o.run(steps=2)
43-
np.testing.assert_almost_equal(o.result.lon[0][1], 4.0, 3)
44-
np.testing.assert_almost_equal(o.result.lat[0][1], 60.031, 3)
43+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4.0, 3)
44+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.063, 3)
4545

4646
# Northwards current, no wind, with Coriolis
4747
o = OpenBerg(loglevel=50)
@@ -53,11 +53,11 @@ def test_openberg_constant_forcing():
5353
o.set_config('drift:coriolis', True)
5454
o.seed_elements(4, 60, time=datetime.now())
5555
o.run(steps=2)
56-
np.testing.assert_almost_equal(o.result.lon[0][1], 4.017, 3)
57-
np.testing.assert_almost_equal(o.result.lat[0][1], 60.031, 3)
56+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4.034, 3)
57+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.060, 3)
5858

5959
# Northwards current, eastwards wind, with Coriolis
60-
o = OpenBerg(loglevel=0)
60+
o = OpenBerg(loglevel=50)
6161
o.set_config('environment:constant:x_sea_water_velocity', 0)
6262
o.set_config('environment:constant:y_sea_water_velocity', 1)
6363
o.set_config('environment:constant:x_wind', 10)
@@ -66,8 +66,77 @@ def test_openberg_constant_forcing():
6666
o.set_config('drift:coriolis', True)
6767
o.seed_elements(4, 60, time=datetime.now())
6868
o.run(steps=2)
69-
np.testing.assert_almost_equal(o.result.lon[0][1], 4.018, 3)
70-
np.testing.assert_almost_equal(o.result.lat[0][1], 60.057, 3)
69+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4.114, 3)
70+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.058, 3)
71+
72+
# No current, eastwards wind
73+
o = OpenBerg(loglevel=0)
74+
o.set_config('environment:constant:x_sea_water_velocity', 0)
75+
o.set_config('environment:constant:y_sea_water_velocity', 0)
76+
o.set_config('environment:constant:x_wind', 10)
77+
o.set_config('environment:constant:y_wind', 0)
78+
o.set_config('drift:horizontal_diffusivity', 0)
79+
o.set_config('drift:coriolis', False)
80+
o.seed_elements(4, 60, time=datetime.now())
81+
o.run(steps=2)
82+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4.107, 3)
83+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.000, 3)
84+
85+
# No current, weaker eastwards wind
86+
o = OpenBerg(loglevel=0)
87+
o.set_config('environment:constant:x_sea_water_velocity', 0)
88+
o.set_config('environment:constant:y_sea_water_velocity', 0)
89+
o.set_config('environment:constant:x_wind', 5)
90+
o.set_config('environment:constant:y_wind', 0)
91+
o.set_config('drift:horizontal_diffusivity', 0)
92+
o.set_config('drift:coriolis', False)
93+
o.seed_elements(4, 60, time=datetime.now())
94+
o.run(steps=2)
95+
# TODO: Half wind gives only 25% of eastwards movement. Should be checked.
96+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4.028, 3)
97+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.000, 3)
98+
99+
# No current, westwards wind
100+
o = OpenBerg(loglevel=0)
101+
o.set_config('environment:constant:x_sea_water_velocity', 0)
102+
o.set_config('environment:constant:y_sea_water_velocity', 0)
103+
o.set_config('environment:constant:x_wind', -10)
104+
o.set_config('environment:constant:y_wind', 0)
105+
o.set_config('drift:horizontal_diffusivity', 0)
106+
o.set_config('drift:coriolis', False)
107+
o.seed_elements(4, 60, time=datetime.now())
108+
o.run(steps=2)
109+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 3.892, 3)
110+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.000, 3)
111+
112+
# No current, northwards wind
113+
o = OpenBerg(loglevel=0)
114+
o.set_config('environment:constant:x_sea_water_velocity', 0)
115+
o.set_config('environment:constant:y_sea_water_velocity', 0)
116+
o.set_config('environment:constant:x_wind', 0)
117+
o.set_config('environment:constant:y_wind', 10)
118+
o.set_config('drift:horizontal_diffusivity', 0)
119+
o.set_config('drift:coriolis', False)
120+
o.seed_elements(4, 60, time=datetime.now())
121+
o.run(steps=2)
122+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4, 3)
123+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 60.054, 3)
124+
125+
# No current, southwards wind
126+
o = OpenBerg(loglevel=0)
127+
o.set_config('environment:constant:x_sea_water_velocity', 0)
128+
o.set_config('environment:constant:y_sea_water_velocity', 0)
129+
o.set_config('environment:constant:x_wind', 0)
130+
o.set_config('environment:constant:y_wind', -10)
131+
o.set_config('drift:horizontal_diffusivity', 0)
132+
o.set_config('drift:coriolis', False)
133+
o.seed_elements(4, 60, time=datetime.now())
134+
o.run(steps=2)
135+
np.testing.assert_almost_equal(o.result.lon.isel(time=-1), 4, 3)
136+
np.testing.assert_almost_equal(o.result.lat.isel(time=-1), 59.946, 3)
137+
138+
139+
71140

72141

73142
def test_openberg_norkyst():

0 commit comments

Comments
 (0)