:py:mod:`opendrift.models.physics_methods` ========================================== .. py:module:: opendrift.models.physics_methods Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: opendrift.models.physics_methods.PhysicsMethods Functions ~~~~~~~~~ .. autoapisummary:: opendrift.models.physics_methods.wind_drift_factor_from_trajectory opendrift.models.physics_methods.distance_between_trajectories opendrift.models.physics_methods.distance_along_trajectory opendrift.models.physics_methods.skillscore_darpa opendrift.models.physics_methods.skillscore_liu_weissberg opendrift.models.physics_methods.plot_wind_drift_factor opendrift.models.physics_methods.oil_wave_entrainment_rate_li2017 opendrift.models.physics_methods.significant_wave_height_from_wind_neumann_pierson opendrift.models.physics_methods.wave_breaking_fraction_from_wind opendrift.models.physics_methods.wave_period_from_wind opendrift.models.physics_methods.verticaldiffusivity_Sundby1983 opendrift.models.physics_methods.verticaldiffusivity_Large1994 opendrift.models.physics_methods.verticaldiffusivity_stepfunction opendrift.models.physics_methods.gls_tke opendrift.models.physics_methods.plot_stokes_profile opendrift.models.physics_methods.stokes_transport_monochromatic opendrift.models.physics_methods.stokes_drift_profile_monochromatic opendrift.models.physics_methods.stokes_drift_profile_exponential opendrift.models.physics_methods.stokes_drift_profile_phillips opendrift.models.physics_methods.stokes_drift_profile_windsea_swell opendrift.models.physics_methods.ftle opendrift.models.physics_methods.wave_stokes_drift_parameterised opendrift.models.physics_methods.wave_significant_height_parameterised opendrift.models.physics_methods.wind_drag_coefficient opendrift.models.physics_methods.windspeed_from_stress_polyfit opendrift.models.physics_methods.declination opendrift.models.physics_methods.equation_of_time opendrift.models.physics_methods.hour_angle opendrift.models.physics_methods.solar_elevation Attributes ~~~~~~~~~~ .. autoapisummary:: opendrift.models.physics_methods.logger opendrift.models.physics_methods.__stokes_coefficients__ opendrift.models.physics_methods.__stokes_hs_coefficients__ .. py:data:: logger .. py:function:: wind_drift_factor_from_trajectory(trajectory_dict, min_period=None) Estimate wind_drift_fator based on wind and current along given trajectory trajectory_dict: dictionary with arrays of same length of the following variables: time, lon, lat, x_wind, y_wind, x_sea_water_velocity, y_sea_water_velocity Returns array of same length minus one of the fitted wind_drift_factor .. py:function:: distance_between_trajectories(lon1, lat1, lon2, lat2) Calculate the distances [m] between two trajectories .. py:function:: distance_along_trajectory(lon, lat) Calculate the distances [m] between points along a trajectory .. py:function:: skillscore_darpa(lon1, lat1, lon2, lat2) Scoring algorithm used for DARPA float challenge 2021. Copied from implementation made by Jean Rabault. Assuming 6 positions, separated by 2 days, where first pos (start) are identical. .. py:function:: skillscore_liu_weissberg(lon_obs, lat_obs, lon_model, lat_model, tolerance_threshold=1) calculate skill score from normalized cumulative seperation distance. Liu and Weisberg 2011. Returns the skill score bewteen 0. and 1. .. py:function:: plot_wind_drift_factor(wdf, azimuth, wmax_plot=None) Polar plot of array of wind drift factor, with associated azimuthal offset .. py:function:: oil_wave_entrainment_rate_li2017(dynamic_viscosity, oil_density, interfacial_tension, significant_wave_height=None, wave_breaking_fraction=None, wind_speed=None, sea_water_density=1028.0) .. py:function:: significant_wave_height_from_wind_neumann_pierson(wind_speed) .. py:function:: wave_breaking_fraction_from_wind(wind_speed, wave_period=None) .. py:function:: wave_period_from_wind(wind_speed) .. py:function:: verticaldiffusivity_Sundby1983(windspeed, depth, mixedlayerdepth=50, background_diffusivity=0) Vertical diffusivity from Sundby (1983) S. Sundby (1983): A one-dimensional model for the vertical distribution of pelagic fish eggs in the mixed layer Deep Sea Research (30) pp. 645-661 .. py:function:: verticaldiffusivity_Large1994(windspeed, depth, mixedlayerdepth=50, background_diffusivity=0) Vertical diffusivity from Large et al. (1994) Depending on windspeed, depth and mixed layer depth (default 50m). .. py:function:: verticaldiffusivity_stepfunction(depth, MLD=20, k_above=0.1, k_below=0.02) eddy diffusivity with discontinuity for testing of mixing scheme .. py:function:: gls_tke(windstress, depth, sea_water_density, tke, generic_length_scale, gls_parameters=None) From LADIM model. .. py:function:: plot_stokes_profile(profiles, view=['vertical', 'birdseye', 'u', 'v'], filename=None) Plot vertical profile of Stokes drift Args: list of dictionary with: u, v: components of Stokes drift vector z: depth in meters, 0 at surface and negative below (e.g. -5 = 5m depth) kwargs: forwarded to plot method view: - 'vertical': magnitude versus depth - 'birdseye': viewed from above - or list of both above with one axis for each case (default) .. py:function:: stokes_transport_monochromatic(mean_wave_period, significant_wave_height) .. py:function:: stokes_drift_profile_monochromatic(stokes_u_surface, stokes_v_surface, significant_wave_height, mean_wave_period, z) Vertical Stokes drift profile assuming a single monochromatic wave Breivik, O., Janssen, P., Bidlot, J., 2014. Approximate stokes drift profiles in deep water. J. Phys. Oceanogr. 44, 2433-2445. doi:10.1175/JPO-D-14-0020.1. .. py:function:: stokes_drift_profile_exponential(stokes_u_surface, stokes_v_surface, significant_wave_height, mean_wave_period, z) Breivik, O., Janssen, P., Bidlot, J., 2014. Approximate stokes drift profiles in deep water. J. Phys. Oceanogr. 44, 2433-2445. doi:10.1175/JPO-D-14-0020.1. .. py:function:: stokes_drift_profile_phillips(stokes_u_surface, stokes_v_surface, significant_wave_height, mean_wave_period, z) Calculate vertical Stokes drift profile from Breivik et al. 2016, A Stokes drift approximation based on the Phillips spectrum, Ocean Mod. 100 .. py:function:: stokes_drift_profile_windsea_swell(stokes_u_surface, stokes_v_surface, swell_mean_direction_to, swell_mean_period, swell_height, wind_sea_mean_direction_to, wind_sea_mean_period, wind_sea_height, z) Calculate vertical Stokes drift profile from Breivik, O., and K. H. Christensen, 2020: A Combined Stokes Drift Profile under Swell and Wind Sea. J. Phys. Oceanogr., 50, 2819-2833, https://doi.org/10.1175/JPO-D-20-0087.1. .. py:function:: ftle(X, Y, delta, duration) Calculate Finite Time Lyapunov Exponents .. py:data:: __stokes_coefficients__ .. py:function:: wave_stokes_drift_parameterised(wind, fetch) Parameterise stokes drift based on pre calculated tables and fetch. .. py:data:: __stokes_hs_coefficients__ .. py:function:: wave_significant_height_parameterised(wind, fetch) Parameterise significant wave height based on pre calculated tables and fetch. .. py:class:: PhysicsMethods Physics methods to be inherited by OpenDriftSimulation class .. py:method:: sea_water_density(T=10.0, S=35.0) :staticmethod: The function gives the density of seawater at one atmosphere pressure as given in : N.P. Fofonoff and R.C. Millard Jr.,1983, Unesco technical papers in marine science no. 44. S = Salinity in promille of the seawater T = Temperature of the seawater in degrees Celsius .. py:method:: skillscore_trajectory(lon_obs, lat_obs, time_obs, duration=None, max_time_offset=timedelta(seconds=60), method='liu-weissberg', **kwargs) Calculate skill-score comparing simulated and observed trajectories A skill score is calculated for each of the modelled trajectories, starting from the closest time/obs of the given modelled trajectory Parameters ---------- lon_obs : array_like lat_obs : array_like The longitude and latitudes of the observed trajectory time_obs : array of datetime The time corresponding to the observed positions duration : timedelta If None, skill score is calculated for the full trajectory, from the start Otherwise, only this sub-part of the trajectory is used. max_time_offset : timedelta The maximum allowed time difference between observed and modelled time Default: 1 minute method : string Currently available: 'liu-wessberg' To be implemented: 'darpa' **kwargs : further arguments for skillscore method (e.g. tolerance_threshold) Returns ------- skillscore : array_like One value for each trajectory .. py:method:: advect_ocean_current(factor=1) .. py:method:: advect_with_sea_ice(factor=1) .. py:method:: advect_wind(factor=1) .. py:method:: stokes_drift(factor=1) .. py:method:: resurface_elements(minimum_depth) .. py:method:: calculate_missing_environment_variables() .. py:method:: wind_speed() .. py:method:: current_speed() .. py:method:: significant_wave_height() .. py:method:: _wave_frequency() .. py:method:: wave_period() .. py:method:: wave_energy() .. py:method:: wave_energy_dissipation() .. py:method:: wave_damping_coefficient() .. py:method:: sea_surface_wave_breaking_fraction() .. py:method:: air_density() .. py:method:: windspeed_from_stress() .. py:method:: solar_elevation() Solar elevation at present time and position of active elements. .. py:method:: sea_floor_depth() Sea floor depth (positive) for presently active elements .. py:method:: sea_surface_height() Sea surface height (positive/negative) for presently active elements .. py:function:: wind_drag_coefficient(windspeed) Large and Pond (1981), J. Phys. Oceanog., 11, 324-336. .. py:function:: windspeed_from_stress_polyfit(wind_stress) Inverting Large and Pond (1981) using polyfit .. py:function:: declination(time) Solar declination in degrees. .. py:function:: equation_of_time(time) Equation of time in minutes. .. py:function:: hour_angle(time, longitude) Solar hour angle in degrees. .. py:function:: solar_elevation(time, longitude, latitude) Solar elevation in degrees.