FNFT
Files | Enumerations | Functions
PRIVATE: Internals related to the nonlinear Schroedinger equation

Files

file  fnft__nse_discretization.h
 Properties of the discretizations for the nonlinear Schroedinger equation.
 
file  fnft__nse_finvscatter.h
 Recovers the signal from a scattering matrix.
 
file  fnft__nse_fscatter.h
 Computes the polynomial approximation of the combined scattering matrix.
 
file  fnft__nse_scatter.h
 Slow forward scattering.
 
file  fnft__nsep_testcases.h
 Provides test cases for the tests of fnft_nsep.
 
file  fnft__nsev_testcases.h
 Provides test cases for the tests of fnft_nsev.
 

Enumerations

enum  fnft__nsep_testcases_t { fnft__nsep_testcases_PLANE_WAVE_FOCUSING, fnft__nsep_testcases_CONSTANT_DEFOCUSING }
 
enum  fnft__nsev_testcases_t { fnft__nsev_testcases_SECH_FOCUSING, fnft__nsev_testcases_SECH_DEFOCUSING, fnft__nsev_testcases_TRUNCATED_SOLITON }
 

Functions

FNFT_UINT fnft__nse_discretization_degree (fnft_nse_discretization_t nse_discretization)
 This routine returns the max degree d of the polynomials in a single scattering matrix or zero if the discretization is unknown. More...
 
FNFT_REAL fnft__nse_discretization_boundary_coeff (fnft_nse_discretization_t nse_discretization)
 This routine returns the boundary coefficient based on the discretization. More...
 
FNFT_INT fnft__nse_discretization_to_akns_discretization (fnft_nse_discretization_t nse_discretization, fnft__akns_discretization_t *const akns_discretization)
 This routine returns akns discretization related to the given nse discretization. More...
 
FNFT_INT fnft__nse_lambda_to_z (const FNFT_UINT n, const FNFT_REAL eps_t, FNFT_COMPLEX *const vals, fnft_nse_discretization_t discretization)
 This routine maps lambda from continuous-time domain to z in the discrete-time domain based on the discretization. More...
 
FNFT_INT fnft__nse_z_to_lambda (const FNFT_UINT n, const FNFT_REAL eps_t, FNFT_COMPLEX *const vals, fnft_nse_discretization_t discretization)
 This routine maps z from the discrete-time domain to lambda in the continuous-time domain based on the discretization. More...
 
FNFT_INT fnft__nse_phase_factor_rho (const FNFT_REAL eps_t, const FNFT_REAL T1, FNFT_REAL *const phase_factor_rho, fnft_nse_discretization_t nse_discretization)
 This routine returns the phase factor for reflection coefficient (rho). It is required for applying boundary conditions to the transfer_matrix based on the discretization. More...
 
FNFT_INT fnft__nse_phase_factor_a (const FNFT_REAL eps_t, const FNFT_UINT D, FNFT_REAL const *const T, FNFT_REAL *const phase_factor_a, fnft_nse_discretization_t nse_discretization)
 This routine returns the phase factor for a coefficient. It is required for applying boundary conditions to the transfer_matrix based on the discretization. More...
 
FNFT_INT fnft__nse_phase_factor_b (const FNFT_REAL eps_t, const FNFT_UINT D, FNFT_REAL const *const T, FNFT_REAL *const phase_factor_b, fnft_nse_discretization_t nse_discretization)
 This routine returns the phase factor for b coefficient. It is required for applying boundary conditions to the transfer_matrix based on the discretization. More...
 
FNFT_INT fnft__nse_finvscatter (const FNFT_UINT deg, FNFT_COMPLEX *const transfer_matrix, FNFT_COMPLEX *const q, const FNFT_REAL eps_t, const FNFT_INT kappa, const fnft_nse_discretization_t discretization)
 Recovers the samples that corresponding to a transfer matrix fast. More...
 
FNFT_UINT fnft__nse_fscatter_numel (FNFT_UINT D, fnft_nse_discretization_t discretization)
 Returns the length of transfer_matrix to be allocated based on the number of samples and discretization. More...
 
FNFT_INT fnft__nse_fscatter (const FNFT_UINT D, FNFT_COMPLEX const *const q, const FNFT_REAL eps_t, const FNFT_INT kappa, FNFT_COMPLEX *const result, FNFT_UINT *const deg_ptr, FNFT_INT *const W_ptr, fnft_nse_discretization_t discretization)
 Fast computation of polynomial approximation of the combined scattering matrix. More...
 
FNFT_INT fnft__nse_scatter_bound_states (const FNFT_UINT D, FNFT_COMPLEX const *const q, FNFT_REAL const *const T, FNFT_UINT *trunc_index_ptr, FNFT_UINT K, FNFT_COMPLEX *bound_states, FNFT_COMPLEX *a_vals, FNFT_COMPLEX *aprime_vals, FNFT_COMPLEX *b, fnft_nse_discretization_t discretization)
 Computes \(a(\lambda)\), \( a'(\lambda) = \frac{\partial a(\lambda)}{\partial \lambda}\) and \(b(\lambda)\) for complex values \(\lambda\) assuming that they are very close to the true bound-states. More...
 
FNFT_INT fnft__nse_scatter_matrix (const FNFT_UINT D, FNFT_COMPLEX const *const q, const FNFT_REAL eps_t, const FNFT_INT kappa, const FNFT_UINT K, FNFT_COMPLEX const *const lambda, FNFT_COMPLEX *const result, fnft_nse_discretization_t discretization)
 Computes the scattering matrix and its derivative. More...
 
FNFT_INT fnft__nsep_testcases_test_fnft (fnft__nsep_testcases_t tc, FNFT_UINT D, FNFT_REAL error_bounds[3], fnft_nsep_opts_t *const opts)
 This routine is used by the tests for fnft_nsep.
More...
 
FNFT_INT fnft__nsev_testcases_test_fnft (fnft__nsev_testcases_t tc, FNFT_UINT D, const FNFT_REAL eb[6], fnft_nsev_opts_t *const opts)
 Routine to run tests for fnft_nsev.
. More...
 

Detailed Description

Enumeration Type Documentation

◆ fnft__nsep_testcases_t

List of currently implemented test cases for the NSE with periodic boundary conditions.

fnft__nsep_testcases_PLANE_WAVE_FOCUSING : Test for focusing NSE (kappa = +1) with a plane wave signal.

fnft__nsep_testcases_CONSTANT_DEFOCUSING : Test for defocusing NSE (kappa = -1) with a constant signal.

◆ fnft__nsev_testcases_t

List of currently implemented test cases for the NSE with vanishing boundary conditions.

fnft__nsev_testcases_SECH_FOCUSING : Test for focusing NSE (kappa = +1) with a sech potential.

fnft__nsev_testcases_SECH_DEFOCUSING : Test for defocusing NSE (kappa = -1) with a sech potential.

fnft__nsev_testcases_TRUNCATED_SOLITON : Test for focusing NSE (kappa = +1) with a truncated solitonic potential.

Function Documentation

◆ fnft__nse_discretization_boundary_coeff()

FNFT_REAL fnft__nse_discretization_boundary_coeff ( fnft_nse_discretization_t  nse_discretization)

This routine returns the boundary coefficient based on the discretization.

The boundary coefficient is the fraction of the step size that a discretized potential extends beyond the last sample. This routine returns this value based on the discretization of type fnft_nse_discretization_t.

Parameters
[in]nse_discretizationThe type of discretization to be used. Should be of type fnft_nse_discretization_t.
Returns
the boundary coefficient, or NAN for discretizations not supported by fnft__nse_fscatter.

◆ fnft__nse_discretization_degree()

FNFT_UINT fnft__nse_discretization_degree ( fnft_nse_discretization_t  nse_discretization)

This routine returns the max degree d of the polynomials in a single scattering matrix or zero if the discretization is unknown.

Parameters
[in]nse_discretizationThe type of discretization to be used. Should be of type fnft_nse_discretization_t.
Returns
polynomial degree, or 0 for discretizations not supported by fnft__nse_fscatter.

◆ fnft__nse_discretization_to_akns_discretization()

FNFT_INT fnft__nse_discretization_to_akns_discretization ( fnft_nse_discretization_t  nse_discretization,
fnft__akns_discretization_t *const  akns_discretization 
)

This routine returns akns discretization related to the given nse discretization.

The function is used by nse specific functions to convert discretization type from fnft_nse_discretization_t to fnft__akns_discretization_t.

Parameters
[in]nse_discretizationThe type of nse discretization. Should be of type fnft_nse_discretization_t.
[out]akns_discretizationThe pointer to the converted discretization of type fnft__akns_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_finvscatter()

FNFT_INT fnft__nse_finvscatter ( const FNFT_UINT  deg,
FNFT_COMPLEX *const  transfer_matrix,
FNFT_COMPLEX *const  q,
const FNFT_REAL  eps_t,
const FNFT_INT  kappa,
const fnft_nse_discretization_t  discretization 
)

Recovers the samples that corresponding to a transfer matrix fast.

Recovers samples q[0], q[1], ..., q[D-1] that, when applied to fnft__nse_fscatter, the transfer matrix provided by the user is reconstructed. Note that the transfer matrix has to be (at least close to) realizable, i.e., such samples have to exist. Also, the transfer matrix has to be sufficiently nice to keep numerical errors tolerable, which in my experience means that the q[n] it corresponds to are "smooth enough" and such that |eps_t*q[n]|<<1 for all n.

More information about the algorithm used here can be found in Wahls and Poor (Proc. IEEE ISIT 2015) and McClary (Geophysics 48(10), 1983).

Parameters
degDegree of the polynomials in the transfer matrix.
[in]transfer_matrixA transfer matrix in the same format as used by fnft__nse_fscatter.
[out]qArray with D=deg/base_deg entries in which the samples are stored, where base_deg is the output of fnft__nse_discretization_degree.
[in]eps_tSee fnft__nse_fscatter.
[in]kappaSee fnft__nse_fscatter.
[in]discretizationSee fnft__nse_fscatter. Currently, only the 2SPLIT2_MODAL and 2SPLIT2A discretizations are supported.

◆ fnft__nse_fscatter()

FNFT_INT fnft__nse_fscatter ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
const FNFT_REAL  eps_t,
const FNFT_INT  kappa,
FNFT_COMPLEX *const  result,
FNFT_UINT *const  deg_ptr,
FNFT_INT *const  W_ptr,
fnft_nse_discretization_t  discretization 
)

Fast computation of polynomial approximation of the combined scattering matrix.

This routine computes the polynomial approximation of the combined scattering matrix by multipying together individual scattering matrices.
Individual scattering matrices depend on the chosen discretization.
The main reference is Wahls and Poor (Proc. ICASSP 2013 ).

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]eps_tStep-size, eps_t \(= (T[1]-T[0])/(D-1) \).
[in]kappa=+1 for the focusing nonlinear Schroedinger equation, =-1 for the defocusing one
[out]resultarray of length nse_fscatter_numel(D,discretization), will contain the combined scattering matrix. Result needs to be pre-allocated with malloc(nse_fscatter_numel(D,discretization)*sizeof(COMPLEX)).
[out]deg_ptrPointer to variable containing degree of the discretization. Determined based on discretization by fnft__nse_discretization_degree.
[in]W_ptrPointer to normalization flag fnft_nsev_opts_t::normalization_flag.
[in]discretizationThe type of discretization to be used. Should be of type fnft_nse_discretization_t. Not all fnft_nse_discretization_t discretizations are supported. Check fnft_nse_discretization_t for list of supported types.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_fscatter_numel()

FNFT_UINT fnft__nse_fscatter_numel ( FNFT_UINT  D,
fnft_nse_discretization_t  discretization 
)

Returns the length of transfer_matrix to be allocated based on the number of samples and discretization.

This routine returns the length 4*D*(nse_discretization_degree(discretization) + 1) to be allocated based on the number of samples and discretization of type discretization.

Parameters
[in]DNumber of samples.
[in]discretizationType of discretization from fnft_nse_discretization_t.
Returns
Returns the length to be allocated. Returns 0 for unknown discretizations.

◆ fnft__nse_lambda_to_z()

FNFT_INT fnft__nse_lambda_to_z ( const FNFT_UINT  n,
const FNFT_REAL  eps_t,
FNFT_COMPLEX *const  vals,
fnft_nse_discretization_t  discretization 
)

This routine maps lambda from continuous-time domain to z in the discrete-time domain based on the discretization.

This routine maps continuous-time domain value lambda to discrete-time domain value z = exp(2i*lambda*eps_t/degree1step), where degree1step is based on the discretization of type fnft_nse_discretization_t. Changes discretization to fnft__akns_discretization_t type and calls fnft__akns_lambda_to_z.

Parameters
[in]nNumber of values to be mapped.
[in]eps_tReal-valued discretization step-size.
[in,out]valsPointer to location of first element of array containing complex-valued continuous-time domain spectral parameter lambda. The values are replaced with discrete-time domain values z.
[in]discretizationDiscretization of type fnft_nse_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_phase_factor_a()

FNFT_INT fnft__nse_phase_factor_a ( const FNFT_REAL  eps_t,
const FNFT_UINT  D,
FNFT_REAL const *const  T,
FNFT_REAL *const  phase_factor_a,
fnft_nse_discretization_t  nse_discretization 
)

This routine returns the phase factor for a coefficient. It is required for applying boundary conditions to the transfer_matrix based on the discretization.

This routine computes the phase correction factor for the computation of the a coefficient from the transfer_matrix. phase_factor_a = -eps_t*D + (T[1]+eps_t*boundary_coeff) - (T[0]-eps_t*boundary_coeff), where eps_t is the step-size, D is the number of samples used to build the transfer_matrix, T is the 2-element time vector defining the signal support and boundary_coeff is based on the discretization of type fnft_nse_discretization_t.

Parameters
[in]eps_tReal-valued discretization step-size.
[in]DPositive interger number of samples used to build transfer_matrix.
[in]TReal-valued 2-element time vector defining the signal support.
[in,out]phase_factor_aPointer to real-valued variable where the computed phase factor will be stored.
[in]nse_discretizationDiscretization of type fnft_nse_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_phase_factor_b()

FNFT_INT fnft__nse_phase_factor_b ( const FNFT_REAL  eps_t,
const FNFT_UINT  D,
FNFT_REAL const *const  T,
FNFT_REAL *const  phase_factor_b,
fnft_nse_discretization_t  nse_discretization 
)

This routine returns the phase factor for b coefficient. It is required for applying boundary conditions to the transfer_matrix based on the discretization.

This routine computes the phase correction factor for the computation of the a coefficient from the transfer_matrix. phase_factor_b = -eps_t*D - (T[1]+eps_t*boundary_coeff) - (T[0]-eps_t*boundary_coeff), where eps_t is the step-size, D is the number of samples used to build the transfer_matrix, T is the 2-element time vector defining the signal support and boundary_coeff is based on the discretization of type fnft_nse_discretization_t. Only for fnft_nse_discretization_2SPLIT2A and fnft_nse_discretization_2SPLIT2_MODAL, phase_factor_b = -eps_t*D - (T[1]+eps_t*boundary_coeff) - (T[0]-eps_t*boundary_coeff)+ eps_t/degree1step where degree1step is based on the discretization of type fnft_nse_discretization_t.

Parameters
[in]eps_tReal-valued discretization step-size.
[in]DPositive interger number of samples used to build transfer_matrix.
[in]TReal-valued 2-element time vector defining the signal support.
[in,out]phase_factor_bPointer to real-valued variable where the computed phase factor will be stored.
[in]nse_discretizationDiscretization of type fnft_nse_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_phase_factor_rho()

FNFT_INT fnft__nse_phase_factor_rho ( const FNFT_REAL  eps_t,
const FNFT_REAL  T1,
FNFT_REAL *const  phase_factor_rho,
fnft_nse_discretization_t  nse_discretization 
)

This routine returns the phase factor for reflection coefficient (rho). It is required for applying boundary conditions to the transfer_matrix based on the discretization.

This routine computes the phase correction factor for the computation of the reflection coefficient from the transfer_matrix. phase_factor_rho = -2.0*(T1 + eps_t*boundary_coeff), where eps_t is the step-size, T1 is the time at the right-boundary and boundary_coeff is based on the discretization of type fnft_nse_discretization_t. Only for fnft_nse_discretization_2SPLIT2A and fnft_nse_discretization_2SPLIT2_MODAL, phase_factor_rho = -2.0*(T1 + eps_t*boundary_coeff)+eps_t/degree1step where degree1step is based on the discretization of type fnft_nse_discretization_t.

Parameters
[in]eps_tReal-valued discretization step-size.
[in]T1Real-valued time at the right-boundary.
[in,out]phase_factor_rhoPointer to real-valued variable where the computed phase factor will be stored.
[in]nse_discretizationDiscretization of type fnft_nse_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_scatter_bound_states()

FNFT_INT fnft__nse_scatter_bound_states ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
FNFT_REAL const *const  T,
FNFT_UINT trunc_index_ptr,
FNFT_UINT  K,
FNFT_COMPLEX bound_states,
FNFT_COMPLEX a_vals,
FNFT_COMPLEX aprime_vals,
FNFT_COMPLEX b,
fnft_nse_discretization_t  discretization 
)

Computes \(a(\lambda)\), \( a'(\lambda) = \frac{\partial a(\lambda)}{\partial \lambda}\) and \(b(\lambda)\) for complex values \(\lambda\) assuming that they are very close to the true bound-states.

The function performs slow direct scattering and is primarily based on the reference Boffetta and Osborne (J. Comput. Physics 1992 ). A forward-backward scheme as mentioned by Aref in (Unpublished) is used to compute the norming constants \(b(\lambda)\).

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,out]trunc_index_ptrPointer containing sample location where the signal will be split to compute the norming constants using the forward-backward scheme. The value should be between 0 and D-1. If the value is D then L1-norm is used to compute the best sample location to split the signal.
[in]KNumber of bound-states.
[in]bound_statesArray of length K, contains the bound-states \(\lambda\).
[out]a_valsArray of length K, contains the values of \(a(\lambda)\).
[out]aprime_valsArray of length K, contains the values of \( a'(\lambda) = \frac{\partial a(\lambda)}{\partial \lambda}\).
[out]bArray of length K, contains the values of \(b(\lambda)\).
[in]discretizationThe type of discretization to be used. Should be of type fnft_nse_discretization_t. Not all nse_discretization_t discretizations are supported. Check fnft_nse_discretization_t for list of supported types.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_scatter_matrix()

FNFT_INT fnft__nse_scatter_matrix ( const FNFT_UINT  D,
FNFT_COMPLEX const *const  q,
const FNFT_REAL  eps_t,
const FNFT_INT  kappa,
const FNFT_UINT  K,
FNFT_COMPLEX const *const  lambda,
FNFT_COMPLEX *const  result,
fnft_nse_discretization_t  discretization 
)

Computes the scattering matrix and its derivative.

The function computes the scattering matrix and the derivative of the scattering matrix with respect to \(\lambda\). The function performs slow direct scattering and is primarily based on the reference Boffetta and Osborne (J. Comput. Physics 1992 ).

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]eps_tStep-size, eps_t \(= (T[1]-T[0])/(D-1) \).
[in]kappa=+1 for the focusing nonlinear Schroedinger equation, =-1 for the defocusing one
[in]KNumber of values of \(\lambda\).
[in]lambdaArray of length K, contains the values of \(\lambda\).
[out]resultArray of length 8*K, contains the values [S11 S12 S21 S22 S11' S12' S21' S22'] where S = [S11, S12; S21, S22] is the scattering matrix.
[in]discretizationThe type of discretization to be used. Should be of type fnft_nse_discretization_t. Not all nse_discretization_t discretizations are supported. Check fnft_nse_discretization_t for list of supported types.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nse_z_to_lambda()

FNFT_INT fnft__nse_z_to_lambda ( const FNFT_UINT  n,
const FNFT_REAL  eps_t,
FNFT_COMPLEX *const  vals,
fnft_nse_discretization_t  discretization 
)

This routine maps z from the discrete-time domain to lambda in the continuous-time domain based on the discretization.

This routine maps discrete-time domain value z to continuous-time domain value lambda = degree1step*log(z)/(2i*eps_t), where degree1step is based on the discretization of type fnft_nse_discretization_t. Changes discretization to fnft__akns_discretization_t type and calls fnft__akns_z_to_lambda.

Parameters
[in]nNumber of values to be mapped.
[in]eps_tReal-valued discretization step-size.
[in,out]valsPointer to location of first element of array containing complex-valued discrete-time domain spectral parameter z. The values are replaced with continuous-time domain values lambda.
[in]discretizationDiscretization of type fnft_nse_discretization_t.
Returns
FNFT_SUCCESS or one of the FNFT_EC_... error codes defined in fnft_errwarn.h.

◆ fnft__nsep_testcases_test_fnft()

FNFT_INT fnft__nsep_testcases_test_fnft ( fnft__nsep_testcases_t  tc,
FNFT_UINT  D,
FNFT_REAL  error_bounds[3],
fnft_nsep_opts_t *const  opts 
)

This routine is used by the tests for fnft_nsep.

It runs the specified test case tc with the specificed number of samples D and the options opts, and tests if several errors stay below the provided error bounds in error_bounds.

Parameters
[in]tcType of test case.
[in]DNumber of samples.
[in]error_boundsReal valued array with 3 elements corresponding to various error bounds.
[in]optsfnft_nsep_opts_t options for the tests.
Returns
If all errors stay below bounds the routine FNFT_SUCCESS. Otherwise, it returns an error code (normally, FNFT_EC_TEST_FAILED).

◆ fnft__nsev_testcases_test_fnft()

FNFT_INT fnft__nsev_testcases_test_fnft ( fnft__nsev_testcases_t  tc,
FNFT_UINT  D,
const FNFT_REAL  eb[6],
fnft_nsev_opts_t *const  opts 
)

Routine to run tests for fnft_nsev.
.

This routine is used by the tests for fnft_nsev. It runs the specified test case tc with the specificed number of samples D and the options opts, and tests if several errors stay below the provided error bounds in eb.

Parameters
[in]tcType of test case.
[in]DNumber of samples.
[in]ebReal valued array with 6 elements corresponding to various error bounds.
[in]optsfnft_nsev_opts_t options for the tests.
Returns
If all errors stay below bounds the routine FNFT_SUCCESS. Otherwise, it returns an error code (normally, FNFT_EC_TEST_FAILED).