FNFT
Loading...
Searching...
No Matches
Files | Classes | Functions
Fast nonlinear Fourier transforms
Collaboration diagram for Fast nonlinear Fourier transforms:

Files

file  fnft_kdv_discretization_t.h
 Lists discretizations for the Korteweg-de Vries equation.
 
file  fnft_kdvv.h
 Fast nonlinear Fourier transform for the vanishing nonlinear Schroedinger equation.
 
file  fnft_manakov_discretization_t.h
 Lists discretizations for the Manakov Equation.
 
file  fnft_manakovv.h
 Fast nonlinear Fourier transform for the vanishing Manakov equation.
 
file  fnft_nse_discretization_t.h
 Lists discretizations for the nonlinear Schroedinger equation.
 
file  fnft_nsep.h
 Fast nonlinear Fourier transform for the periodic nonlinear Schroedinger equation.
 
file  fnft_nsev.h
 Fast nonlinear Fourier transform for the vanishing nonlinear Schroedinger equation.
 
file  fnft_version.h
 Provides a function to query the version of the FNFT library.
 

Classes

struct  fnft_kdvv_opts_t
 Stores additional options for the routine fnft_kdvv. More...
 
struct  fnft_manakovv_opts_t
 Stores additional options for the routine fnft_manakovv. More...
 
struct  fnft_nsep_opts_t
 Stores additional options for the routine fnft_nsep. More...
 
struct  fnft_nsev_opts_t
 Stores additional options for the routine fnft_nsev. More...
 

Functions

fnft_kdvv_opts_t fnft_kdvv_default_opts ()
 Creates a new options variable for fnft_kdvv with default settings.
 
FNFT_INT fnft_kdvv (const FNFT_UINT D, FNFT_COMPLEX const *const q, FNFT_REAL const *const T, const FNFT_UINT M, FNFT_COMPLEX *const contspec, FNFT_REAL const *const XI, FNFT_UINT *const K_ptr, FNFT_COMPLEX *const bound_states, FNFT_COMPLEX *const normconsts_or_residues, fnft_kdvv_opts_t *opts)
 Nonlinear Fourier transform for the Korteweg-de Vries equation with vanishing boundary conditions.
 
fnft_manakovv_opts_t fnft_manakovv_default_opts ()
 Creates a new options variable for fnft_manakovv with default settings.
 
FNFT_UINT fnft_manakovv_max_K (const FNFT_UINT D, fnft_manakovv_opts_t const *const opts)
 Returns the maximum number of bound states that can be detected by fnft_manakovv.
 
FNFT_INT fnft_manakovv (const FNFT_UINT D, FNFT_COMPLEX const *const q1, FNFT_COMPLEX const *const q2, FNFT_REAL const *const T, const FNFT_UINT M, FNFT_COMPLEX *const contspec, FNFT_REAL const *const XI, FNFT_UINT *const K_ptr, FNFT_COMPLEX *const bound_states, FNFT_COMPLEX *const normconsts_or_residues, const FNFT_INT kappa, fnft_manakovv_opts_t *opts)
 Nonlinear Fourier transform for the Manakov equation with vanishing boundary conditions. Fast algorithms are used if the discretization supports it.
 
fnft_nsep_opts_t fnft_nsep_default_opts ()
 Creates a new options variable for fnft_nsep with default settings.
 
FNFT_INT fnft_nsep (const FNFT_UINT D, FNFT_COMPLEX const *const q, FNFT_REAL const *const T, FNFT_REAL const phase_shift, FNFT_UINT *const K_ptr, FNFT_COMPLEX *const main_spec, FNFT_UINT *const M_ptr, FNFT_COMPLEX *const aux_spec, FNFT_REAL *const sheet_indices, const FNFT_INT kappa, fnft_nsep_opts_t *opts)
 Fast nonlinear Fourier transform for the nonlinear Schroedinger equation with (quasi-)periodic boundary conditions.
 
fnft_nsev_opts_t fnft_nsev_default_opts ()
 Creates a new options variable for fnft_nsev with default settings.
 
FNFT_UINT fnft_nsev_max_K (const FNFT_UINT D, fnft_nsev_opts_t const *const opts)
 Returns the maximum number of bound states that can be detected by fnft_nsev.
 
FNFT_INT fnft_nsev (const FNFT_UINT D, FNFT_COMPLEX const *const q, FNFT_REAL const *const T, const FNFT_UINT M, FNFT_COMPLEX *const contspec, FNFT_REAL const *const XI, FNFT_UINT *const K_ptr, FNFT_COMPLEX *const bound_states, FNFT_COMPLEX *const normconsts_or_residues, const FNFT_INT kappa, fnft_nsev_opts_t *opts)
 Nonlinear Fourier transform for the nonlinear Schroedinger equation with vanishing boundary conditions. Fast algorithms are used if the discretization supports it.
 
FNFT_INT fnft_version (FNFT_UINT *const major, FNFT_UINT *const minor, FNFT_UINT *const patch, char suffix[FNFT_VERSION_SUFFIX_MAXLEN+1])
 

Detailed Description

Function Documentation

◆ fnft_kdvv()

FNFT_INT fnft_kdvv ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
FNFT_REAL const *const  T,
const FNFT_UINT  M,
FNFT_COMPLEX *const  contspec,
FNFT_REAL const *const  XI,
FNFT_UINT *const  K_ptr,
FNFT_COMPLEX *const  bound_states,
FNFT_COMPLEX *const  normconsts_or_residues,
fnft_kdvv_opts_t opts 
)

Nonlinear Fourier transform for the Korteweg-de Vries equation with vanishing boundary conditions.

This routine computes the nonlinear Fourier transform for the Korteweg-de Vries equation

\[ q_x + 6qq_{t} + q_{ttt}=0, \quad q=q(x,t), \]

of Gardner et al. (Phys. Rev. Lett., 1967) for initial conditions with vanishing boundaries

\[ \lim_{t\to \pm \infty }q(x_0,t) = 0 \text{ sufficiently rapidly.} \]

Fast algorithms are used if the discretization supports it.
The main references are:

The routine also utilizes ideas from the following papers:

The routine supports all discretizations of type fnft_kdv_discretization_t. The following discretizations use fast algorithms which have a computational complexity of \( \mathcal{O}(D\log^2 D)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_kdv_discretization_2SPLIT1A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT1B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT2A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT2B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT2S(_VANILLA)
  • fnft_kdv_discretization_2SPLIT2_MODAL(_VANILLA)
  • fnft_kdv_discretization_2SPLIT3A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT3B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT3S(_VANILLA)
  • fnft_kdv_discretization_2SPLIT4A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT4B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT5A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT5B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT6A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT6B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT7A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT7B(_VANILLA)
  • fnft_kdv_discretization_2SPLIT8A(_VANILLA)
  • fnft_kdv_discretization_2SPLIT8B(_VANILLA)
  • fnft_kdv_discretization_4SPLIT4A(_VANILLA)
  • fnft_kdv_discretization_4SPLIT4B(_VANILLA)

The discretizations nft_kdv_discretization_2SPLIT2_MODAL(_VANILLA) require that q[i] * eps_t^2 is unequal to -1 for all samples i.

The following discretizations use classical algorithms which have a computational complexity of \( \mathcal{O}(D^2)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_kdv_discretization_BO(_VANILLA)
  • fnft_kdv_discretization_CF4_2(_VANILLA)
  • fnft_kdv_discretization_CF4_3(_VANILLA)
  • fnft_kdv_discretization_CF5_3(_VANILLA)
  • fnft_kdv_discretization_CF6_4(_VANILLA)
  • fnft_kdv_discretization_ES4(_VANILLA)
  • fnft_kdv_discretization_TES4(_VANILLA)

The accuray of the computed quantities for a given signal depends primarily on the number of samples \( D\) and the numerical method. When the exact spectrum is is know, the accuracy can be quantified by defining a suitable error. The error usually decreases with increasing \( D\) assuming everthing else remains the same. The rate at which the error decreases with increase in \( D\) is known as the order of the method. The orders of the various discretizations can be found at fnft_kdv_discretization_t. The orders of the discretizations which use exponential splitting schemes should be the same as their base methods but can deviate when accuracy of the splitting scheme is low.

Application of one step Richardson extrapolation should in theory increase the order by one. However, it can deviate both ways. In case of some smooth signals the order may increase by two instead of one. On the other hand for discontinuous signals it maybe deterimental to apply Richardson extrapolation. It is observed that Richardson extrapolation has no effect on the order of the residues and norming constants, except for the following case:

  • CF5_3, residues and norming constants have order six instead of five.
Parameters
[in]DNumber of samples
[in]qArray of length D, contains samples \( q(t_n)=q(x_0, t_n) \), where \( t_n = T[0] + n(T[1]-T[0])/(D-1) \) and \(n=0,1,\dots,D-1\), of the to-be-transformed signal in ascending order (i.e., \( q(t_0), q(t_1), \dots, q(t_{D-1}) \)). The imaginary part of every sample must be 0.
[in]TArray of length 2, contains the position in time of the first and of the last sample. It should be \(T[0]<T[1]\).
[in]MNumber of points at which the continuous spectrum (aka reflection coefficient) should be computed.
[out]contspecArray of length M in which the routine will store the desired samples \( r(\xi_m) \) of the continuous spectrum (aka reflection coefficient) in ascending order, where \( \xi_m = XI[0]+m(XI[1]-XI[0])/(M-1) \) and \(m=0,1,\dots,M-1\). Has to be preallocated by the user. If NULL is passed instead, the continuous spectrum will not be computed. By changing the options, it is also possible to compute the values of \( a(\xi) \) and \( b(\xi) \) instead. In that case, twice the amount of memory has to be allocated.
[in]XIArray of length 2, contains the position of the first and the last sample of the continuous spectrum. It should be \(XI[0]<XI[1]\). The grid that is defined by \( \xi_m = XI[0]+m(XI[1]-XI[0])/(M-1) \) and \(m=0,1,\dots,M-1\) is not allowed to contain \( \xi_m=0\), because \( \xi=0\) is a singularity of the scattering problem for the KdV equation. \( XI\) can also be NULL if contspec==NULL.
[in,out]K_ptrUpon entry, *K_ptr should contain the length of the array bound_states. Upon return, *K_ptr contains the number of actually detected bound states. If the length of the array bound_states was not sufficient to store all of the detected bound states, a warning is printed and as many bound states as possible are returned instead. Note that in order to skip the computation of the bound states completely, it is not sufficient to pass *K_ptr==0. Instead, one needs to pass bound_states==NULL.
[in,out]bound_statesArray. Upon entry, only if the option bound_state_localization=fnft_kdvv_bsloc_NEWTON is passed, this array should contain initial guesses for the bound states. (These initial guesses should be positive imaginary numbers.) Otherwise the input values are ignored. Upon return, the routine has stored the detected bound states (aka eigenvalues) in the first *K_ptr entries of this array. If NULL is passed instead, the discrete spectrum will not be computed. Has to be preallocated by the user. The user can choose an arbitrary length. Typically, D is a good choice.
[out]normconsts_or_residuesArray of the same length as bound_states. Upon return, the routine has stored the residues (aka spectral amplitudes) \(\rho_k = b_k\big/ \frac{da(\lambda_k)}{d\lambda}\) in the first *K_ptr entries of this array. By passing a proper opts, it is also possible to store the norming constants (the ' \( b_k \)') or both. Has to be pre-allocated by the user. If NULL is passed instead, the residues will not be computed.
[in]optsPointer to a fnft_kdvv_opts_t object. The object can be used to modify the behavior of the routine. Use the routine fnft_kdvv_default_opts to generate such an object and modify as desired. It is also possible to pass NULL, in which case the routine will use the default options. The user is reponsible to freeing the object after the routine has returned.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft_kdvv_default_opts()

fnft_kdvv_opts_t fnft_kdvv_default_opts ( )

Creates a new options variable for fnft_kdvv with default settings.

Returns
A fnft_kdvv_opts_t object with the following options.
bound_state_filtering = fnft_kdvv_bsfilt_FULL
bound_state_localization = fnft_kdvv_bsloc_SUBSAMPLE_AND_REFINE
niter = 10
discspec_type = fnft_kdvv_dstype_NORMING_CONSTANTS
contspec_type = fnft_kdvv_cstype_REFLECTION_COEFFICIENT
normalization_flag = 1
discretization = fnft_kdv_discretization_2SPLIT4B
richardson_extrapolation_flag = 0

◆ fnft_manakovv()

FNFT_INT fnft_manakovv ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q1,
FNFT_COMPLEX const *const  q2,
FNFT_REAL const *const  T,
const FNFT_UINT  M,
FNFT_COMPLEX *const  contspec,
FNFT_REAL const *const  XI,
FNFT_UINT *const  K_ptr,
FNFT_COMPLEX *const  bound_states,
FNFT_COMPLEX *const  normconsts_or_residues,
const FNFT_INT  kappa,
fnft_manakovv_opts_t opts 
)

Nonlinear Fourier transform for the Manakov equation with vanishing boundary conditions. Fast algorithms are used if the discretization supports it.

This routine computes the nonlinear Fourier transform for the Manakov equation (vector nonlinear Schrödinger equation)

\[ iq_x + q_{tt} \pm 2q|q|^2=0, \quad q=[q1(x,t); q2(x,t)] \]

of S.V. Manakov (Soviet Physics-JETP 38.2 (1974): 248-253) for initial conditions with vanishing boundaries

\[ \lim_{t\to \pm \infty }q(x_0,t) = 0 \text{ sufficiently rapidly.} \]


The main references are:

The routine also utilizes ideas from the following papers:

The routine supports all discretizations of type fnft_manakov_discretization_t. The following discretizations use fast algorithms which have a computational complexity of \( \mathcal{O}(D\log^2 D)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_manakov_discretization_2SPLIT3A
  • fnft_manakov_discretization_2SPLIT3B
  • fnft_manakov_discretization_2SPLIT4A
  • fnft_manakov_discretization_2SPLIT4B
  • fnft_manakov_discretization_2SPLIT6B
  • fnft_manakov_discretization_4SPLIT4A
  • fnft_manakov_discretization_4SPLIT4B
    • fnft_manakov_discretization_4SPLIT6B
    • fnft_manakov_discretization_FTES4_4A
    • fnft_manakov_discretization_FTES4_suzuki

The following discretizations use classical algorithms which have a computational complexity of \( \mathcal{O}(D^2)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_manakov_discretization_BO
  • fnft_manakov_discretization_CF4_2

The accuray of the computed quantities for a given signal depends primarily on the number of samples \( D\) and the numerical method. When the exact spectrum is is know, the accuracy can be quantified by defining a suitable error. The error usually decreases with increasing \( D\) assuming everthing else remains the same. The rate at which the error decreases with increase in \( D\) is known as the order of the method. The orders of the various discretizations can be found at fnft_manakov_discretization_t. The orders of the discretizations which use exponential splitting schemes should be the same as their base methods but can deviate when accuracy of the splitting scheme is low.

Application of one step Richardson extrapolation should in theory increase the order by one. However, it can deviate both ways. In case of some smooth signals the order may increase by two instead of one. On the other hand for discontinuous signals it maybe detrimental to apply Richardson extrapolation. In the following case the error with Richardson extrapolation has been observed to be higher than without:

  • Focusing case 4split6B
Parameters
[in]DNumber of samples
[in]q1First element of the potential function. Array of length D, contains samples \( q1(t_n)=q1(x_0, t_n) \), where \( t_n = T[0] + n(T[1]-T[0])/(D-1) \) and \(n=0,1,\dots,D-1\), of the to-be-transformed signal in ascending order (i.e., \( q1(t_0), q1(t_1), \dots, q1(t_{D-1}) \))
[in]q2Second element of the potential function, analoguous to q1
[in]TArray of length 2, contains the position in time of the first and of the last sample. It should be \(T[0]<T[1]\).
[in]MNumber of points at which the continuous spectrum (aka reflection coefficient and/or NFT coefficient a, b1, b2) should be computed.
[out]contspecArray of length 2*M, 3*M or 5*M dependent on desired contspec of type fnft_manakovv_cstype_t in which the routine will store the reflection coefficients \( \rho \) and/or the NFT coefficients a, b1, b2 in ascending order: \( \rho_1(\xi_m) = b1(\xi_m)/a(\xi_m) \) followed by \( \rho_1(\xi_m) = b1(\xi_m)/a(\xi_m) \) if only the reflection coefficients are requested, where \( \xi_m = XI[0]+m(XI[1]-XI[0])/(M-1) \) and \(m=0,1,\dots,M-1\). \( a(\xi_m) \), \( b_1(\xi_m) \) followed by \( b_2(\xi_m) \) if the NFT coefficients are requested and both these arrays if both are requested. If NULL is passed instead, the continuous spectrum will not be computed.
[in]XIArray of length 2, contains the position of the first and the last sample of the continuous spectrum. It should be \(XI[0]<XI[1]\). Can also be NULL if contspec==NULL.
[in,out]K_ptrUpon entry, *K_ptr should contain the length of the array bound_states. Upon return, *K_ptr contains the number of actually detected bound states. If the length of the array bound_states was not sufficient to store all of the detected bound states, a warning is printed and as many bound states as possible are returned instead. Note that in order to skip_ the computation of the bound states completely, it is not sufficient to pass *K_ptr==0. Instead, one needs to pass bound_states==NULL.
[in,out]bound_statesArray. Upon entry, only if the option bound_state_localization=fnft_manakovvv_bsloc_NEWTON is passed, this array should contain initial guesses for the bound states. Otherwise the input values are ignored. Upon return, the routine has stored the detected bound states (aka eigenvalues) in the first *K_ptr entries of this array. If NULL is passed instead, the discrete spectrum will not be computed. Has to be preallocated by the user. The user can choose an arbitrary length. Typically, D is a good choice.
[out]normconsts_or_residuesArray of the same length as bound_states. Upon return, the routine has stored the residues (aka spectral amplitudes) \(\rho_k = b_k\big/ \frac{da(\lambda_k)}{d\lambda}\) in the first *K_ptr entries of this array. By passing a proper opts, it is also possible to store the norming constants (the ' \( b_k \)') or both. Has to be pre-allocated by the user. If NULL is passed instead, the residues will not be computed.
[in]kappa=+1 for the focusing Manakov equation, =-1 for the defocusing one.
[in]optsPointer to a fnft_manakovv_opts_t object. The object can be used to modify the behavior of the routine. Use the routine fnft_manakovv_default_opts to generate such an object and modify as desired. It is also possible to pass NULL, in which case the routine will use the default options. The user is reponsible to freeing the object after the routine has returned.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft_manakovv_default_opts()

fnft_manakovv_opts_t fnft_manakovv_default_opts ( )

Creates a new options variable for fnft_manakovv with default settings.

Returns
A fnft_manakovv_opts_t object with the following options.
bound_state_filtering = fnft_manakovv_bsfilt_FULL
bound_state_localization = fnft_manakovv_bsloc_SUBSAMPLE_AND_REFINE
niter = 10
Dsub = 0
discspec_type = fnft_manakovv_dstype_NORMING_CONSTANTS
contspec_type = fnft_manakovv_cstype_REFLECTION_COEFFICIENT
normalization_flag = 1
discretization = fnft_manakov_discretization_2SPLIT4B
richardson_extrapolation_flag = 0

◆ fnft_manakovv_max_K()

FNFT_UINT fnft_manakovv_max_K ( const FNFT_UINT  D,
fnft_manakovv_opts_t const *const  opts 
)

Returns the maximum number of bound states that can be detected by fnft_manakovv.

Parameters
[in]Dnumber of samples
[in]optsoptions for the fnft, used for getting the discretization method. Should be of type fnft_manakovv_opts_t
Returns
Returns the maximum number of bound states that can be detected by fnft_manakovv, based on the degree of the total transition matrix

◆ fnft_nsep()

FNFT_INT fnft_nsep ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
FNFT_REAL const *const  T,
FNFT_REAL const  phase_shift,
FNFT_UINT *const  K_ptr,
FNFT_COMPLEX *const  main_spec,
FNFT_UINT *const  M_ptr,
FNFT_COMPLEX *const  aux_spec,
FNFT_REAL *const  sheet_indices,
const FNFT_INT  kappa,
fnft_nsep_opts_t opts 
)

Fast nonlinear Fourier transform for the nonlinear Schroedinger equation with (quasi-)periodic boundary conditions.



This routine computes the nonlinear Fourier transform for the nonlinear Schroedinger equation

\[ iq_x + q_{tt} \pm 2q|q|^2=0, \quad q=q(x,t), \]

of Kotlyarov and Its (Problems of Mathemetical Physics and Functional Analysis, 1976) for initial conditions with (quasi-)periodic boundary conditions,

\[ q(x_0,t + L) = mq(x_0,t),|m|=1. \]


The main references for the numerical algorithm are:

The routine supports the following discretizations of type fnft_nse_discretization_t:

  • fnft_nse_discretization_2SPLIT1A
  • fnft_nse_discretization_2SPLIT1B
  • fnft_nse_discretization_2SPLIT2A
  • fnft_nse_discretization_2SPLIT2B
  • fnft_nse_discretization_2SPLIT2S
  • fnft_nse_discretization_2SPLIT2_MODAL
  • fnft_nse_discretization_2SPLIT3A
  • fnft_nse_discretization_2SPLIT3B
  • fnft_nse_discretization_2SPLIT3S
  • fnft_nse_discretization_2SPLIT4A
  • fnft_nse_discretization_2SPLIT4B
  • fnft_nse_discretization_2SPLIT5A
  • fnft_nse_discretization_2SPLIT5B
  • fnft_nse_discretization_2SPLIT6A
  • fnft_nse_discretization_2SPLIT6B
  • fnft_nse_discretization_2SPLIT7A
  • fnft_nse_discretization_2SPLIT7B
  • fnft_nse_discretization_2SPLIT8A
  • fnft_nse_discretization_2SPLIT8B
  • fnft_nse_discretization_4SPLIT4A

For the Newton refinement mode only, one of the following should be used instead:

  • fnft_nse_discretization_BO
  • fnft_nse_discretization_CF4_2
Parameters
[in]DNumber of samples. Should be power of two and \( D\geq 2\).
[in]qArray of length D, contains samples \( q(t_n)=q(x_0, t_n) \), where \( t_n = T[0] + nL/D \), where \(L=T[1]-T[0]\) is the period and \(n=0,1,\dots,D-1\), of the to-be-transformed signal in ascending order (i.e., \( q(t_0), q(t_1), \dots, q(t_{D-1}) \))
[in]TArray of length 2. T[0] is the position in time of the first sample. T[1] is the beginning of the next period. (The location of the last sample is thus \(t_{D-1}=T[1]-L/D\).) It should be \(T[0]<T[1]\).
[in]phase_shiftReal scalar constant. It is the change in the phase over one quasi-period, \( \angle(q(t+L)/q(t))\). For periodic signals it will be 0.
[in,out]K_ptrUpon entry, *K_ptr should contain the length of the array main_spec (except for Newton localization, see below). Upon return, *K_ptr contains the number of actually detected points in the main spectrum. If the length of the array main_spec was not sufficient to store all of the detected points in the main spectrum, a warning is printed and as many points as possible are returned instead. Note that in order to skip the computation of the main spectrum completely, it is not sufficient to pass *K_ptr==NULL. Instead, pass main_spec==NULL. EXCEPTION: If opts->localization==fnft_nsep_loc_NEWTON, then *K_ptr should contain the number of initial guess per spine point (see opts->points_per_spine).
[out]main_specArray. Upon return, the routine has stored the detected main specrum points (i.e., the points for which the trace of the monodromy matrix is either +2 or -2) in the first *K_ptr entries of this array. If NULL is passed instead, the main spectrum will not be computed. Has to be preallocated by the user. The user can choose an arbitrary length. Typically, D is a good choice. If opts->localization==fnft_nsep_loc_NEWTON, then main_spec should be of the form [g1_1, ..., gK_1, g1_2, ..., gK_2, ..., g1_P, ..., gK_P] when the routine is called, where gk_n is the k-th initial guess for the n-th spine point, K=*K_ptr, and P=opts->points_per_spine.
[in,out]M_ptrUpon entry, *M_ptr should contain the length of the array aux_spec. Upon return, *M_ptr contains the number of actually detected points in the auxiliary spectrum. If the length of the array aux_spec was not sufficient to store all of the detected points in the auxiliary spectrum, a warning is printed and as many points as possible are returned instead. Note that in order to skip the computation of the auxiliary spectrum completely, it is not sufficient to pass *M_ptr==NULL. Instead, pass aux_spec==NULL.
[out]aux_specArray. Upon return, the routine has stored the detected auxiliary specrum points (i.e., the points for which the upper right element of the monodromy matrix is zero) in the first *M_ptr entries of this array. If NULL is passed instead, the auxiliary spectrum will not be computed. Has to be preallocated by the user. The user can choose an arbitrary length. Typically, D is a good choice.
[in]sheet_indicesNot yet implemented. Pass NULL.
[in]kappa=+1 for the focusing nonlinear Schroedinger equation, =-1 for the defocusing one.
[in]optsPointer to a fnft_nsep_opts_t object. The object can be used to modify the behavior of the routine. Use the routine fnft_nsep_default_opts to generate such an object and modify as desired. It is also possible to pass NULL, in which case the routine will use the default options. The user is reponsible to freeing the object after the routine has returned.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft_nsep_default_opts()

fnft_nsep_opts_t fnft_nsep_default_opts ( )

Creates a new options variable for fnft_nsep with default settings.

Returns
A fnft_nsep_opts_t object with the following options.
localization = fnft_nsep_loc_MIXED
filtering = fnft_nsep_filt_AUTO
max_evals = 20
bounding_box[0] = -FNFT_INF
bounding_box[1] = FNFT_INF
bounding_box[2] = -FNFT_INF
bounding_box[3] = FNFT_INF
normalization_flag = 1
discretization = fnft_nse_discretization_2SPLIT2A
floquet_range = {-1, 1}
floquet_nvals = 2
Dsub = 0

◆ fnft_nsev()

FNFT_INT fnft_nsev ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
FNFT_REAL const *const  T,
const FNFT_UINT  M,
FNFT_COMPLEX *const  contspec,
FNFT_REAL const *const  XI,
FNFT_UINT *const  K_ptr,
FNFT_COMPLEX *const  bound_states,
FNFT_COMPLEX *const  normconsts_or_residues,
const FNFT_INT  kappa,
fnft_nsev_opts_t opts 
)

Nonlinear Fourier transform for the nonlinear Schroedinger equation with vanishing boundary conditions. Fast algorithms are used if the discretization supports it.

This routine computes the nonlinear Fourier transform for the nonlinear Schroedinger equation

\[ iq_x + q_{tt} \pm 2q|q|^2=0, \quad q=q(x,t), \]

of Zakharov and Shabat (Soviet. Phys. JTEP 31(1), 1972) for initial conditions with vanishing boundaries

\[ \lim_{t\to \pm \infty }q(x_0,t) = 0 \text{ sufficiently rapidly.} \]


The main references are:

The routine also utilizes ideas from the following papers:

The routine supports all discretizations of type fnft_nse_discretization_t. The following discretizations use fast algorithms which have a computational complexity of \( \mathcal{O}(D\log^2 D)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_nse_discretization_2SPLIT1A
  • fnft_nse_discretization_2SPLIT1B
  • fnft_nse_discretization_2SPLIT2A
  • fnft_nse_discretization_2SPLIT2B
  • fnft_nse_discretization_2SPLIT2S
  • fnft_nse_discretization_2SPLIT2_MODAL
  • fnft_nse_discretization_2SPLIT3A
  • fnft_nse_discretization_2SPLIT3B
  • fnft_nse_discretization_2SPLIT3S
  • fnft_nse_discretization_2SPLIT4A
  • fnft_nse_discretization_2SPLIT4B
  • fnft_nse_discretization_2SPLIT5A
  • fnft_nse_discretization_2SPLIT5B
  • fnft_nse_discretization_2SPLIT6A
  • fnft_nse_discretization_2SPLIT6B
  • fnft_nse_discretization_2SPLIT7A
  • fnft_nse_discretization_2SPLIT7B
  • fnft_nse_discretization_2SPLIT8A
  • fnft_nse_discretization_2SPLIT8B
  • fnft_nse_discretization_4SPLIT4A
  • fnft_nse_discretization_4SPLIT4B

The following discretizations use classical algorithms which have a computational complexity of \( \mathcal{O}(D^2)\) for \( D\) point continuous spectrum given \( D\) samples:

  • fnft_nse_discretization_BO
  • fnft_nse_discretization_CF4_2
  • fnft_nse_discretization_CF4_3
  • fnft_nse_discretization_CF5_3
  • fnft_nse_discretization_CF6_4
  • fnft_nse_discretization_ES4
  • fnft_nse_discretization_TES4

The accuray of the computed quantities for a given signal depends primarily on the number of samples \( D\) and the numerical method. When the exact spectrum is is know, the accuracy can be quantified by defining a suitable error. The error usually decreases with increasing \( D\) assuming everthing else remains the same. The rate at which the error decreases with increase in \( D\) is known as the order of the method. The orders of the various discretizations can be found at fnft_nse_discretization_t. The orders of the discretizations which use exponential splitting schemes should be the same as their base methods but can deviate when accuracy of the splitting scheme is low.

Application of one step Richardson extrapolation should in theory increase the order by one. However, it can deviate both ways. In case of some smooth signals the order may increase by two instead of one. On the other hand for discontinuous signals it maybe deterimental to apply Richardson extrapolation. In the following case the order of the method has been observed to be less than expected:

  • Focusing case CF5_3, residues and norming constants have order five instead of six.
Parameters
[in]DNumber of samples
[in]qArray of length D, contains samples \( q(t_n)=q(x_0, t_n) \), where \( t_n = T[0] + n(T[1]-T[0])/(D-1) \) and \(n=0,1,\dots,D-1\), of the to-be-transformed signal in ascending order (i.e., \( q(t_0), q(t_1), \dots, q(t_{D-1}) \))
[in]TArray of length 2, contains the position in time of the first and of the last sample. It should be \(T[0]<T[1]\).
[in]MNumber of points at which the continuous spectrum (aka reflection coefficient) should be computed.
[out]contspecArray of length M in which the routine will store the desired samples \( r(\xi_m) \) of the continuous spectrum (aka reflection coefficient) in ascending order, where \( \xi_m = XI[0]+m(XI[1]-XI[0])/(M-1) \) and \(m=0,1,\dots,M-1\). Has to be preallocated by the user. If NULL is passed instead, the continuous spectrum will not be computed. By changing the options, it is also possible to compute the values of \( a(\xi) \) and \( b(\xi) \) instead. In that case, twice the amount of memory has to be allocated.
[in]XIArray of length 2, contains the position of the first and the last sample of the continuous spectrum. It should be \(XI[0]<XI[1]\). Can also be NULL if contspec==NULL.
[in,out]K_ptrUpon entry, *K_ptr should contain the length of the array bound_states. Upon return, *K_ptr contains the number of actually detected bound states. If the length of the array bound_states was not sufficient to store all of the detected bound states, a warning is printed and as many bound states as possible are returned instead. Note that in order to skip the computation of the bound states completely, it is not sufficient to pass *K_ptr==0. Instead, one needs to pass bound_states==NULL.
[in,out]bound_statesArray. Upon entry, only if the option bound_state_localization=fnft_nsev_bsloc_NEWTON is passed, this array should contain initial guesses for the bound states. Otherwise the input values are ignored. Upon return, the routine has stored the detected bound states (aka eigenvalues) in the first *K_ptr entries of this array. If NULL is passed instead, the discrete spectrum will not be computed. Has to be preallocated by the user. The user can choose an arbitrary length. Typically, D is a good choice.
[out]normconsts_or_residuesArray of the same length as bound_states. Upon return, the routine has stored the residues (aka spectral amplitudes) \(\rho_k = b_k\big/ \frac{da(\lambda_k)}{d\lambda}\) in the first *K_ptr entries of this array. By passing a proper opts, it is also possible to store the norming constants (the ' \( b_k \)') or both. Has to be pre-allocated by the user. If NULL is passed instead, the residues will not be computed.
[in]kappa=+1 for the focusing nonlinear Schroedinger equation, =-1 for the defocusing one.
[in]optsPointer to a fnft_nsev_opts_t object. The object can be used to modify the behavior of the routine. Use the routine fnft_nsev_default_opts to generate such an object and modify as desired. It is also possible to pass NULL, in which case the routine will use the default options. The user is reponsible to freeing the object after the routine has returned.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft_nsev_default_opts()

fnft_nsev_opts_t fnft_nsev_default_opts ( )

Creates a new options variable for fnft_nsev with default settings.

Returns
A fnft_nsev_opts_t object with the following options.
bound_state_filtering = fnft_nsev_bsfilt_FULL
bound_state_localization = fnft_nsev_bsloc_SUBSAMPLE_AND_REFINE
niter = 100
tol = -1.0
discspec_type = fnft_nsev_dstype_NORMING_CONSTANTS
contspec_type = fnft_nsev_cstype_REFLECTION_COEFFICIENT
normalization_flag = 1
discretization = fnft_nse_discretization_2SPLIT4B
richardson_extrapolation_flag = 0
bounding_box = {NAN, NAN, NAN, NAN}

◆ fnft_nsev_max_K()

FNFT_UINT fnft_nsev_max_K ( const FNFT_UINT  D,
fnft_nsev_opts_t const *const  opts 
)

Returns the maximum number of bound states that can be detected by fnft_nsev.

Parameters
[in]DNumber of samples that will be passed to fnft_nsev. Should be larger than zero.
[in]optsOptions that will be passed to fnft_nsev. If NULL is passed, the default options will be used.
Returns
Returns the maximum number of bound states or zero on error.

◆ fnft_version()

FNFT_INT fnft_version ( FNFT_UINT *const  major,
FNFT_UINT *const  minor,
FNFT_UINT *const  patch,
char  suffix[FNFT_VERSION_SUFFIX_MAXLEN+1] 
)

Provides the version of the FNFT library.

The version is of the form $major.$minor.$patch$suffix, e.g., "0.2.1-al". The function simply returns the corresponding values defined fnft_config.h.

Parameters
[out]major*major is overwritten with the value of FNFT_VERSION_MAJOR
[out]minor*minor is overwritten with the value of FNFT_VERSION_MINOR
[out]patch*patch is overwritten with the value of FNFT_VERSION_PATCH
[out]suffixsuffix is overwritten with the value of FNFT_VERSION_SUFFIX
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.