Add frequency (Hz) coordinate to datasets

[ryancoe]

I often do the following after running the Capytaine solver:

ds = ds.assign_coords(freq=("omega",ds['omega'].values/2/np.pi))
ds['freq'].attrs['long_name'] = 'Frequency'
ds['freq'].attrs['units'] = 'Hz'

which allows me to more easily make plots with the x_label being Frequency [Hz] by using da.plot(x='freq', ...), e.g.,

Fexc = (ds['diffraction_force'] + ds['Froude_Krylov_force']).sel(wave_direction=0)

fig, ax = plt.subplots(nrows=2,
                       sharex=True)

np.abs(Fexc).plot(x='freq', ax=ax[0], hue='influenced_dof')
ax[0].set_ylabel('$|F_e|$')

np.arctan(np.real(Fexc)/np.imag(Fexc)).plot(x='freq', ax=ax[1], hue='influenced_dof', add_legend=False)
ax[1].set_ylabel('$\\angle{ F_e } $')

for axi in ax:
    axi.set_title('')
    axi.label_outer()
    axi.grid()
    axi.autoscale(enable=True, axis='x', tight=True)
    axi.spines[['right', 'top']].set_visible(False)
    axi.set_xscale('log')

Looking at the source (for the first time in a while...), I think this could be handled cleanly in xarray.py in two ways:

1. (more complete) Add frequency to the frequency_keys -- this would allow the user to formulate problems using frequency in addition to the existing options.

capytaine/capytaine/io/xarray.py

Line 77 in f5f0b9f

frequency_keys = keys_in_dataset & {'omega', 'period', 'wavelength', 'wavenumber'}

2. (easier) Add my lines above to the assemble_dataset function to add the frequency coordinate to the results of any problem.

I can definitely do 2, and probably do 1 with some guidance. Any thoughts @mancellin?

[mancellin]

We could completely support frequency in the same way that we support omega, period, wavenumber and wavelength.
We'd need to:

• Add a frequency attribute to LinearPotentialFlowProblem and update accordingly
  

  capytaine/capytaine/bem/problems_and_results.py

  Line 77 in f5f0b9f

  def _get_frequencies(self, omega, period, wavenumber, wavelength):

  
  and
  

  capytaine/capytaine/bem/problems_and_results.py

  Line 182 in f5f0b9f

  def _asdict(self):

• Update cpt.io.xarray.problems_from_dataset to pass the frequency data to the initialization of LinearPotentialFlow
• Add an output in cpt.io.xarray.assemble_dataset following the pattern of

  capytaine/capytaine/io/xarray.py

  Line 368 in f5f0b9f

  if period and main_freq_type != "period":

This fifth case is making the code a bit longer but not more complicated.

But I don't know if it is worth it, given that the conversion from omega is straightforward.
I don't remember seeing any paper using the frequency in Hz to plot hydrodynamic data. Do you have any example of usage in the literature that would justify adding it to the code?

[ryancoe]

Thanks @mancellin! As far as how useful this is, it might be good to let others weigh in. I definitely do present plots like this with frequency in Hz because I find it easier to interpret (no need to divide by $2\pi$ to understand what the dynamics are at, e.g., 5s). It also aligns with the more common presentation of ocean spectra, which I think are most typically in $S(f)$. Maybe @dav-og @cmichelenstrofer @RubendeBruin @salhus and others can say if they feel similarly?

[dav-og]

I agree with you @ryancoe

Prof. Salter once chastised me for plotting hydrodynamic coefficients as a function of rad/s..."engineers have a hard time thinking in rad/s...its much easier for them to understand if you plot these in Hz or S"

I don't remember seeing any paper using the frequency in Hz to plot hydrodynamic data. Do you have any example of usage in the literature that would justify adding it to the code?

@mancellin it used to be more common back in the day - the Edinburgh Wave Power Project used seconds extensively (e.g. Pizer - Numerical Modelling of Wave Energy Absorbers:

[RubendeBruin]

I guess it really depends on who you are working with. For practical reports and presentations, I always use seconds (like above).
When programming or presenting to scientists I usually use rad/s because you are working with radians for phase angles as well. But Hz may be the best (or worst) of both worlds: "easy" to convert to seconds while having a good view of infinite periods (f->0).

[mancellin]

@mancellin it used to be more common back in the day - the Edinburgh Wave Power Project used seconds extensively (e.g. Pizer - Numerical Modelling of Wave Energy Absorbers

Well, that is an example of using periods in seconds (already implemented in Capytaine). My question was about frequencies in Hertz.

But if several of you find useful to add that, I won't refuse a PR.

[salhus]

I think the consensus seems to be towards Hz being the best of both worlds..., but I still personally find rad/s quite useful. The dynamic response shows up as a well distributed lobe, it is easier to interpret phases and how it relates to excitations, easy to map it to diff equations, and to Laplace domain, especially for applications where I want to see the resonant peaks, such as controls. Whereas the periods way, plots it the other way round... slow/large waves on the left and fast/short waves on the right....

Again a matter of preference and application...

@mancellin BEM Rosetta seems to switch between these worlds quite quickly and seamlessly. A tighter integration with Rosetta, might be an expedient route.