Choice of numerical parameters in HOS-ocean
===========================================

This section describes the main elements driving the choice of the numerical parameters in HOS-ocean. The model is based on pseudo-spectral formalism, which decompose the different surface spatial quantities on the normal modes of the domain.

Spatial/Modal discretization
----------------------------

Horizontal: ``n1`` and ``n2``
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
It is possible to provide general recommendations in the choice of the discretization in x-direction ``n1``. However, as any numerical parameter, a convergence study should be carried out to test the influence of changing ``n1``. One should take into account the following points:

- With a pseudo-spectral approach, the discretization defines the largest wavenumber that can be solved k\ :sub:`max` \=(n1-1)*2π/xlen
- One wants to simulate a wavefield with a given peak wavenumber k\ :sub:`p` (or eq. regular wavenumber k\ :sub:`0`)
- **In HOS-ocean one should use:** k\ :sub:`max` \≈ 8 k\ :sub:`p`

For the discretization along y-direction ``n2``, the choice should be done carefully after a convergence study. General recommendation cannot be provided since it is highly dependent on the directional spreading associated with the simulated wave field.

Orders of nonlinearity ``mHOS`` and ``M_l``
-------------------------------------------
The HOS method relies on an iterative procedure that is controlled by the so-called HOS order of nonlinearity ``mHOS`` (and ``M_l`` for varying bathymetry solved with ``i_method=2``). These parameters can be seen as pure numerical parameters and consequently their effect can be characterized by a classical convergence study. However, regarding the HOS formulation, it is also possible to relate ``mHOS`` to physical processes at play during wave propagation. Indeed, this order of nonlinearity corresponds to the nonlinear wave interactions between wave components that you take into account.

- Flat bottom (original HOS method)

    - ``mHOS=3`` corresponds to third-order of nonlinearity, or equivalently the so-called four-wave interactions. With this set-up, the HOS model is equivalent to Zakharov equation for water waves. **This is the minimum value suggested**, if one focuses on nonlinear wave properties since it includes the most important nonlinear features at play during wave propagation
    - ``mHOS=1`` corresponds to a pure linear solution
    - For practical applications, this parameter is set in the range [3,8] and ``mHOS=5`` is a standard value

- Varying bottom

    - ``i_method=1``: the same order of nonlinearity ``mHOS`` is used for free surface nonlinearities and influence of bathymetry on wave propagation. Convergence study has to be carried out to ensure that the bottom variation is correctly taken into account.
    - ``i_method=2``: the choice of ``mHOS`` is the same than for flat bottom configuration and a convergence study with respect to ``M_l`` needs to be carried out to ensure that the bottom variation is correctly taken into account.


Dealiasing ``p1`` and ``p2``
----------------------------
Pseudo-spectral methods face the possible issue of aliasing. A zero-padding procedure is set-up as dealiasing methodology in HOS-ocean, which is controlled thanks to the parameters ``p1`` and ``p2``. The following recommendations can be made

- Increased accuracy is achieved with total dealiasing ``p1=mHOS``. This is at the expense of computational effort as well as numerical stability
- A typical optimal value for this parameter (robustness/accuracy) is ``p1=3``
- For unidirectional configuration ``n2=1``, one should use ``p2=1``
- For multidirectional waves, one should use ``p2=p1``

Temporal integration
--------------------
The accuracy of the time integration scheme is controlled with the tolerance parameter ``toler``. This is a numerical parameter, which convergence can be studied for a detailed study conducted with HOS-ocean. **A value of** ``toler=10^(-6)`` **is suggested**.

.. note:: In pseudo-spectral models, without dissipation terms or breaking models, the stability of the numerical solution is highly dependent on the choice of the numerical parameters. This has been studied in details in :cite:p:`ducrozet2017applicability`