Hamiltonian Module

The Hamiltonian module provides comprehensive tools for constructing polynomial normal forms of Hamiltonian systems around equilibrium points in the Circular Restricted Three-Body Problem (CR3BP). It implements the complete pipeline from physical Hamiltonians through coordinate transformations to normal forms and center manifold reductions.

The package supports analysis around all five libration points (L1-L5) with appropriate handling of hyperbolic directions at collinear points and elliptic directions at triangular points.

hamiltonian.py

The hamiltonian module provides polynomial Hamiltonian construction using Chebyshev and Legendre expansions.

hiten.algorithms.hamiltonian.hamiltonian._build_T_polynomials(poly_x, poly_y, poly_z, max_deg, psi_table, clmo_table, encode_dict_list)[source]

Compute Chebyshev polynomials of the first kind for collinear point expansions.

Generate the sequence T_n(r) where r = x / sqrt(x^2 + y^2 + z^2) using the classical Chebyshev recurrence relation adapted for Cartesian coordinates in nondimensional CR3BP units.

Parameters:
  • poly_x (List[ndarray]) – Polynomial representations of Cartesian coordinates x, y, z in nondimensional distance units.

  • poly_y (List[ndarray]) – Polynomial representations of Cartesian coordinates x, y, z in nondimensional distance units.

  • poly_z (List[ndarray]) – Polynomial representations of Cartesian coordinates x, y, z in nondimensional distance units.

  • max_deg (int) – Maximum polynomial degree n such that T_n is computed.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo_table (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup tables mapping multi-indices to coefficient positions.

Returns:

Numba typed list where element i contains coefficients of T_i(r).

Return type:

List[List[ndarray]]

Notes

Implements the Chebyshev recurrence relation: T_0 = 1, T_1 = r, T_n = 2*r*T_{n-1} - T_{n-2}

Adapted for Cartesian coordinates as: T_n = ((2n-1)/n) * x * T_{n-1} - ((n-1)/n) * (x^2 + y^2 + z^2) * T_{n-2}

Used for inverse distance expansions around collinear equilibrium points (L1, L2, L3) in the CR3BP where the gravitational potential has the form U = -sum_{n>=2} c_n * T_n(r).

See also

_build_R_polynomials()

Auxiliary polynomials for Lindstedt-Poincare formulation.

_build_potential_U()

Gravitational potential construction using Chebyshev polynomials.

References

Szebehely, V. (1967). Theory of Orbits. Academic Press, Chapter 7.

Compute Chebyshev polynomials of the first kind for collinear point expansions. Generates the sequence T_n(r) where r = x / sqrt(x^2 + y^2 + z^2) using the classical Chebyshev recurrence relation adapted for Cartesian coordinates in nondimensional CR3BP units.

hiten.algorithms.hamiltonian.hamiltonian._build_R_polynomials(poly_x, poly_y, poly_z, poly_T, max_deg, psi_table, clmo_table, encode_dict_list)[source]

Generate auxiliary R_n polynomials for Lindstedt-Poincare formulation.

Compute the sequence R_n required for the Lindstedt-Poincare right-hand side polynomials in collinear point normal form computations.

Parameters:
  • poly_x (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • poly_y (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • poly_z (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • poly_T (List[List[ndarray]]) – Chebyshev polynomials from _build_T_polynomials().

  • max_deg (int) – Maximum polynomial degree for R_n computation.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo_table (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup tables mapping multi-indices to coefficient positions.

Returns:

Auxiliary polynomials R_0, R_1, …, R_max_deg.

Return type:

List[List[ndarray]]

Notes

Implements the recurrence relation: R_0 = -1 R_1 = -3*x R_n = ((2n+3)/(n+2)) * x * R_{n-1} - ((2n+2)/(n+2)) * T_n - ((n+1)/(n+2)) * rho^2 * R_{n-2}

where rho^2 = x^2 + y^2 + z^2.

These polynomials appear in the Lindstedt-Poincare method for constructing periodic solutions around collinear equilibrium points in the CR3BP.

See also

_build_T_polynomials()

Chebyshev polynomials used in this computation.

_build_lindstedt_poincare_rhs_polynomials()

Final right-hand side construction using R_n polynomials.

References

Jorba, A., Masdemont, J. (1999). Dynamics in the center manifold of the collinear points. Physica D, 132(1-2), 189-213.

Generate auxiliary R_n polynomials for Lindstedt-Poincare formulation. Computes the sequence R_n required for the Lindstedt-Poincare right-hand side polynomials in collinear point normal form computations.

hiten.algorithms.hamiltonian.hamiltonian._build_potential_U(poly_T, point, max_deg, psi_table)[source]

Assemble gravitational potential expansion for collinear points.

Construct the effective potential U = -sum_{n>=2} c_n * T_n(r) using Chebyshev polynomial expansions and libration point coefficients.

Parameters:
  • poly_T (List[List[ndarray]]) – Chebyshev polynomials T_n from _build_T_polynomials().

  • point (object) – Libration point object with method cn(k) returning the potential coefficient c_k (dimensionless).

  • max_deg (int) – Maximum polynomial degree for potential truncation.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

Returns:

Polynomial representation of the gravitational potential U in nondimensional energy units.

Return type:

List[ndarray]

Notes

The gravitational potential around collinear equilibrium points has the series expansion U = -sum_{n>=2} c_n * T_n(r) where T_n are Chebyshev polynomials and c_n are coefficients determined by the mass parameter and equilibrium point location.

The expansion starts from n=2 since the linear terms (n=0,1) are absorbed into the equilibrium point definition and coordinate transformation.

See also

_build_T_polynomials()

Chebyshev polynomials used in this expansion.

_build_physical_hamiltonian_collinear()

Complete Hamiltonian construction using this potential.

References

Szebehely, V. (1967). Theory of Orbits. Academic Press, Section 7.2.

Assemble gravitational potential expansion for collinear points. Constructs the effective potential U = -sum_{n>=2} c_n * T_n(r) using Chebyshev polynomial expansions and libration point coefficients.

hiten.algorithms.hamiltonian.hamiltonian._build_kinetic_energy_terms(poly_px, poly_py, poly_pz, max_deg, psi_table, clmo_table, encode_dict_list)[source]

Build kinetic energy polynomial T = (1/2) * (px^2 + py^2 + pz^2).

Construct the kinetic energy contribution to the Hamiltonian using canonical momentum polynomials in nondimensional CR3BP units.

Parameters:
  • poly_px (List[ndarray]) – Polynomial representations of canonical momenta px, py, pz in nondimensional momentum units.

  • poly_py (List[ndarray]) – Polynomial representations of canonical momenta px, py, pz in nondimensional momentum units.

  • poly_pz (List[ndarray]) – Polynomial representations of canonical momenta px, py, pz in nondimensional momentum units.

  • max_deg (int) – Maximum polynomial degree for kinetic energy computation.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo_table (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup tables mapping multi-indices to coefficient positions.

Returns:

Polynomial representation of kinetic energy T in nondimensional energy units.

Return type:

List[ndarray]

Notes

The kinetic energy in the rotating frame has the standard form T = (1/2) * (px^2 + py^2 + pz^2) where px, py, pz are the canonical momenta conjugate to the position coordinates x, y, z.

In nondimensional CR3BP units, energy is normalized by mu * n^2 * a^2 where mu is the total mass, n is the mean motion, and a is the primary-secondary separation.

See also

_build_rotational_terms()

Coriolis terms that couple with kinetic energy in rotating frame.

_build_physical_hamiltonian_collinear()

Complete Hamiltonian assembly including kinetic energy.

Build kinetic energy polynomial T = (1/2) * (px^2 + py^2 + pz^2). Constructs the kinetic energy contribution to the Hamiltonian using canonical momentum polynomials in nondimensional CR3BP units.

hiten.algorithms.hamiltonian.hamiltonian._build_rotational_terms(poly_x, poly_y, poly_px, poly_py, max_deg, psi_table, clmo_table, encode_dict_list)[source]

Construct Coriolis (rotational) terms C = y*px - x*py for rotating frame.

Build the Coriolis contribution to the Hamiltonian that arises from the transformation to the rotating synodic reference frame in the CR3BP.

Parameters:
  • poly_x (List[ndarray]) – Polynomial representations of position coordinates x, y in nondimensional distance units.

  • poly_y (List[ndarray]) – Polynomial representations of position coordinates x, y in nondimensional distance units.

  • poly_px (List[ndarray]) – Polynomial representations of canonical momenta px, py in nondimensional momentum units.

  • poly_py (List[ndarray]) – Polynomial representations of canonical momenta px, py in nondimensional momentum units.

  • max_deg (int) – Maximum polynomial degree for Coriolis term computation.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo_table (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup tables mapping multi-indices to coefficient positions.

Returns:

Polynomial representation of Coriolis terms C in nondimensional energy units.

Return type:

List[ndarray]

Notes

The Coriolis terms C = y*px - x*py arise from the transformation to the rotating synodic frame where the two primaries remain fixed. These terms couple position and momentum and are essential for capturing the dynamics in the rotating reference frame.

The Coriolis terms vanish in the inertial frame but become important in the rotating frame where they represent fictitious forces due to the frame rotation.

See also

_build_kinetic_energy_terms()

Kinetic energy terms that couple with Coriolis terms.

_build_physical_hamiltonian_collinear()

Complete Hamiltonian including Coriolis terms.

References

Szebehely, V. (1967). Theory of Orbits. Academic Press, Section 2.3.

Construct Coriolis (rotational) terms C = y*px - x*py for rotating frame. Builds the Coriolis contribution to the Hamiltonian that arises from the transformation to the rotating synodic reference frame in the CR3BP.

hiten.algorithms.hamiltonian.hamiltonian._build_physical_hamiltonian_collinear(point, max_deg)[source]

Build complete rotating-frame Hamiltonian H = T + U + C for collinear points.

Construct the full CR3BP Hamiltonian around collinear equilibrium points by combining kinetic energy, gravitational potential, and Coriolis terms.

Parameters:
  • point (object) – Libration point object with method cn(k) returning the potential coefficient c_k (dimensionless) for the gravitational expansion.

  • max_deg (int) – Maximum polynomial degree for Hamiltonian truncation.

Returns:

Complete Hamiltonian polynomial coefficients up to max_deg in nondimensional energy units.

Return type:

List[ndarray]

Notes

The rotating-frame Hamiltonian has the form: H = T + U + C = (1/2)*(px^2 + py^2 + pz^2) + U(x,y,z) + (y*px - x*py)

where: - T: kinetic energy in canonical coordinates - U: gravitational potential expanded in Chebyshev polynomials - C: Coriolis terms from rotating frame transformation

The constant term (equilibrium potential energy) is removed since we work with perturbations around the equilibrium point.

All coordinates are in nondimensional CR3BP units: - Positions: primary-secondary separation = 1 - Momenta: normalized by sqrt(mu * n^2 * a^3) - Energy: normalized by mu * n^2 * a^2

See also

_build_kinetic_energy_terms()

Kinetic energy polynomial construction.

_build_potential_U()

Gravitational potential using Chebyshev expansion.

_build_rotational_terms()

Coriolis terms for rotating frame.

_build_physical_hamiltonian_triangular()

Hamiltonian construction for triangular points.

Examples

>>> # Construct Hamiltonian for L1 point with degree 6 truncation
>>> H = _build_physical_hamiltonian_collinear(l1_point, max_deg=6)

Build complete rotating-frame Hamiltonian H = T + U + C for collinear points. Constructs the full CR3BP Hamiltonian around collinear equilibrium points by combining kinetic energy, gravitational potential, and Coriolis terms.

hiten.algorithms.hamiltonian.hamiltonian._build_A_polynomials(poly_x, poly_y, poly_z, d_x, d_y, max_deg, psi_table, clmo_table, encode_dict_list)[source]

Generate Legendre-type polynomials A_n for triangular point expansions.

Compute the sequence A_n used for inverse distance expansions around triangular equilibrium points (L4/L5) in the CR3BP using a Legendre-type recurrence relation.

Parameters:
  • poly_x (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • poly_y (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • poly_z (List[ndarray]) – Polynomial representations of coordinates x, y, z in nondimensional distance units.

  • d_x (float) – Components of primary offset vector in nondimensional distance units. For L4/L5, these represent the displacement from the equilibrium to the primary masses.

  • d_y (float) – Components of primary offset vector in nondimensional distance units. For L4/L5, these represent the displacement from the equilibrium to the primary masses.

  • max_deg (int) – Maximum polynomial degree for A_n computation.

  • psi_table (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo_table (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup tables mapping multi-indices to coefficient positions.

Returns:

Legendre-type polynomials A_0, A_1, …, A_max_deg for inverse distance expansions.

Return type:

List[List[ndarray]]

Notes

Implements the recurrence relation: A_0 = 1 A_1 = d_x*x + d_y*y A_{n+1} = ((2n+1)/(n+1)) * (d_dot_r) * A_n - (n/(n+1)) * rho^2 * A_{n-1}

where d_dot_r = d_x*x + d_y*y and rho^2 = x^2 + y^2 + z^2.

These polynomials are used to expand the inverse distances 1/r_PS and 1/r_PJ in the triangular point Hamiltonian, where PS and PJ denote distances to the primary and secondary masses respectively.

See also

_build_physical_hamiltonian_triangular()

Triangular point Hamiltonian using these polynomials.

_build_T_polynomials()

Analogous Chebyshev polynomials for collinear points.

References

Gomez, G., et al. (2001). Dynamics and Mission Design Near Libration Points. World Scientific, Chapter 3.

Generate Legendre-type polynomials A_n for triangular point expansions. Computes the sequence A_n used for inverse distance expansions around triangular equilibrium points (L4/L5) in the CR3BP using a Legendre-type recurrence relation.

hiten.algorithms.hamiltonian.hamiltonian._build_physical_hamiltonian_triangular(point, max_deg)[source]

Build rotating-frame Hamiltonian for triangular points (L4/L5).

Construct the complete CR3BP Hamiltonian around triangular equilibria using Legendre-type polynomial expansions for the inverse distance terms to the primary masses.

Parameters:
  • point (object) – Triangular libration point object (L4 or L5) with attributes mu (mass parameter) and sign (+1 for L4, -1 for L5).

  • max_deg (int) – Maximum polynomial degree for Hamiltonian truncation.

Returns:

Complete triangular point Hamiltonian coefficients up to max_deg in nondimensional energy units.

Return type:

List[ndarray]

Notes

The triangular point Hamiltonian has the form: H = (1/2)*(px^2 + py^2 + pz^2) + y*px - x*py + ((1/2) - mu)*x + s*sqrt(3)/2*y - (1-mu)/r_PS - mu/r_PJ

where: - (px, py, pz): canonical momenta in nondimensional momentum units - (x, y, z): local coordinates around the triangular equilibrium - mu: mass parameter mu = m2/(m1 + m2) - s: sign parameter (+1 for L4, -1 for L5) - r_PS, r_PJ: distances to primary and secondary masses

The distances to the primaries in local coordinates are: r_PS^2 = (x - x_S)^2 + (y - y_S)^2 + z^2 r_PJ^2 = (x - x_J)^2 + (y - y_J)^2 + z^2

with offset vectors: (x_S, y_S) = (+1/2, s*sqrt(3)/2) (x_J, y_J) = (-1/2, s*sqrt(3)/2)

At equilibrium, both distances equal 1, allowing expansion of inverse distances using Legendre-type polynomials via the binomial series: 1/r = (1 + u)^(-1/2) = sum_{k>=0} C(-1/2,k) * u^k where u = r^2 - 1 and C(alpha,k) are binomial coefficients.

The constant term (equilibrium potential energy = -1) is removed.

See also

_build_A_polynomials()

Legendre-type polynomials for inverse distance expansions.

_build_physical_hamiltonian_collinear()

Analogous Hamiltonian for collinear points.

References

Gomez, G., Llibre, J., Martinez, R., Simo, C. (2001). Dynamics and Mission Design Near Libration Points. World Scientific, Equation 3.57.

Build rotating-frame Hamiltonian for triangular points (L4/L5). Constructs the complete CR3BP Hamiltonian around triangular equilibria using Legendre-type polynomial expansions for the inverse distance terms to the primary masses.

hiten.algorithms.hamiltonian.hamiltonian._build_lindstedt_poincare_rhs_polynomials(point, max_deg)[source]

Compute right-hand side polynomials for Lindstedt-Poincare method.

Generate the polynomial right-hand sides for the x, y, z equations in the first iteration of the Lindstedt-Poincare method for constructing periodic solutions around collinear equilibrium points.

Parameters:
  • point (object) – Libration point object with method cn(k) returning the potential coefficient c_k (dimensionless) for the gravitational expansion.

  • max_deg (int) – Maximum polynomial degree for right-hand side truncation.

Returns:

Polynomial coefficients for the x-, y-, and z-equation right-hand sides respectively, in nondimensional acceleration units.

Return type:

tuple of (List[ndarray], List[ndarray], List[ndarray])

Notes

The Lindstedt-Poincare method constructs periodic solutions by expanding both the solution and the frequency in powers of a small parameter. The right-hand sides have the form:

RHS_x = sum_{n>=2} c_{n+1} * (n+1) * T_n(r) RHS_y = y * sum_{n>=2} c_{n+1} * R_{n-1}(r) RHS_z = z * sum_{n>=2} c_{n+1} * R_{n-1}(r)

where T_n are Chebyshev polynomials, R_n are auxiliary polynomials, and c_k are the potential coefficients from the libration point.

These polynomials appear in the differential equations for the first-order corrections to the linear center manifold dynamics around collinear points.

See also

_build_T_polynomials()

Chebyshev polynomials used in x-equation.

_build_R_polynomials()

Auxiliary polynomials used in y,z-equations.

Examples

>>> # Generate RHS polynomials for L1 point with degree 6
>>> rhs_x, rhs_y, rhs_z = _build_lindstedt_poincare_rhs_polynomials(l1_point, 6)

References

Jorba, A., Masdemont, J. (1999). Dynamics in the center manifold of the collinear points. Physica D, 132(1-2), 189-213.

Compute right-hand side polynomials for Lindstedt-Poincare method. Generates the polynomial right-hand sides for the x, y, z equations in the first iteration of the Lindstedt-Poincare method for constructing periodic solutions around collinear equilibrium points.

lie.py

The lie module provides Lie series transformations for canonical coordinate changes.

hiten.algorithms.hamiltonian.lie._solve_homological_equation(p_elim, n, eta, clmo)[source]

Solve homological equation for polynomial generating function coefficients.

Compute the generating function coefficients that eliminate specified terms from a polynomial Hamiltonian through canonical transformation. The homological equation relates the generating function to the terms being eliminated.

Parameters:
  • p_elim (ndarray) – Coefficient array of polynomial terms to be eliminated, in nondimensional energy units. Each element corresponds to a monomial term.

  • n (int) – Degree of the homogeneous polynomial terms being processed.

  • eta (ndarray, shape (3,)) – Complex eigenvalue array [lambda, i*omega1, i*omega2] where lambda is the real eigenvalue and omega1, omega2 are imaginary frequencies, all in nondimensional frequency units.

  • clmo (ndarray) – Packed multi-index table from _init_index_tables().

Returns:

Coefficient array for the generating function of degree n, in nondimensional action units (energy * time).

Return type:

ndarray

Notes

The homological equation for canonical transformations has the form: {H_0, G} + H_1 = N

where H_0 is the quadratic normal form, G is the generating function, H_1 contains terms to eliminate, and N is the new normal form.

For each monomial with multi-index k = [k0, k1, k2, k3, k4, k5], the generating function coefficient is: g_k = -h_k / ((k3-k0)*lambda + (k4-k1)*i*omega1 + (k5-k2)*i*omega2)

Small divisor terms (denominator < 1e-14) are skipped to avoid numerical instability, corresponding to near-resonant terms that remain in the normal form.

See also

_apply_poly_transform()

Apply the computed generating function via Lie series.

_decode_multiindex()

Decode packed multi-indices used in coefficient computation.

References

Meyer, K.R., Hall, G.R. (1992). Introduction to Hamiltonian Dynamical Systems and the N-Body Problem. Springer-Verlag, Chapter 4.

Solve homological equation for polynomial generating function coefficients. Computes the generating function coefficients that eliminate specified terms from a polynomial Hamiltonian through canonical transformation.

hiten.algorithms.hamiltonian.lie._apply_poly_transform(poly_H, p_G_n, deg_G, N_max, psi, clmo, encode_dict_list, tol)[source]

Apply Lie series transformation to polynomial Hamiltonian.

Transform a polynomial Hamiltonian using a canonical transformation generated by a polynomial generating function via the Lie series method. This implements the exponential Lie operator through iterated Poisson brackets.

Parameters:
  • poly_H (List[ndarray]) – Original Hamiltonian polynomial coefficients organized by degree, in nondimensional energy units.

  • p_G_n (ndarray) – Coefficient array for the generating function of degree deg_G, in nondimensional action units (energy * time).

  • deg_G (int) – Degree of the generating function polynomial.

  • N_max (int) – Maximum polynomial degree to retain in the transformed Hamiltonian.

  • psi (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo (ndarray) – Packed multi-index table from _init_index_tables().

  • encode_dict_list (List[dict]) – Lookup dictionaries mapping multi-indices to coefficient positions.

  • tol (float) – Numerical tolerance for cleaning small coefficients to prevent accumulation of round-off errors.

Returns:

Transformed Hamiltonian polynomial coefficients organized by degree, in nondimensional energy units.

Return type:

List[ndarray]

Notes

The Lie series transformation implements the canonical transformation: H’ = exp(L_G) H = H + {H,G} + (1/2!)*{{H,G},G} + (1/3!)*{{{H,G},G},G} + …

where: - L_G is the Lie operator associated with generating function G - {A,B} denotes the Poisson bracket of polynomials A and B - The series is truncated when terms exceed the maximum degree N_max

The transformation preserves the Hamiltonian structure and implements a canonical coordinate change that eliminates specific polynomial terms while maintaining the symplectic nature of the flow.

The number of terms in the Lie series is determined by the degree constraints: for deg_G > 2, at most (N_max - deg_G)/(deg_G - 2) + 1 terms can contribute to the final result within the degree bound.

See also

_solve_homological_equation()

Compute generating function coefficients for term elimination.

_polynomial_poisson_bracket()

Poisson bracket computation used in the Lie series.

References

Deprit, A. (1969). Canonical transformations depending on a small parameter. Celestial Mechanics, 1(1), 12-30.

Hori, G. (1966). Theory of general perturbations with unspecified canonical variables. Publications of the Astronomical Society of Japan, 18, 287-296.

Apply Lie series transformation to polynomial Hamiltonian. Transforms a polynomial Hamiltonian using a canonical transformation generated by a polynomial generating function via the Lie series method.

transforms.py

The transforms module provides coordinate system transformations and complexification utilities.

hiten.algorithms.hamiltonian.transforms._build_complexification_matrix(mix_indices)[source]

Build complexification transformation matrix for canonical coordinate pairs.

The transformation is given by: q_j(complex) = (q_j(real) + i*p_j(real)) / sqrt(2) p_j(complex) = (i*q_j(real) + p_j(real)) / sqrt(2) where j is the index of the canonical coordinate pair.

Parameters:

mix_indices (tuple of int) – Indices of canonical coordinate pairs to complexify.

Returns:

Complexification transformation matrix.

Return type:

np.ndarray

References

Jorba, A. (1999). A methodology for the numerical computation of normal forms, centre manifolds and first integrals of Hamiltonian systems. Experimental Mathematics, 8(2), 155-195.

Build complexification transformation matrix for coordinate pairs. Creates the matrix that transforms real canonical coordinate pairs to complex coordinates.

hiten.algorithms.hamiltonian.transforms._M(mix_pairs=(1, 2))[source]

Return complexification transformation matrix for canonical coordinate pairs.

Generate the unitary transformation matrix that converts real canonical coordinates to complex coordinates for specified coordinate pairs, typically used to complexify elliptic directions while leaving hyperbolic directions real.

Parameters:

mix_pairs (tuple of int, default (1, 2)) – Indices of canonical coordinate pairs to complexify. Each index j corresponds to the pair (q_{j+1}, p_{j+1}) where j=0,1,2 maps to pairs (q1,p1), (q2,p2), (q3,p3) respectively.

Returns:

Unitary complexification matrix M where complex_coords = M @ real_coords. The matrix is structured as 2x2 blocks for each canonical pair.

Return type:

ndarray, shape (6, 6), complex

Notes

For each coordinate pair j in mix_pairs, the transformation is: q_j(complex) = (q_j(real) + i*p_j(real)) / sqrt(2) p_j(complex) = (i*q_j(real) + p_j(real)) / sqrt(2)

This complexification is particularly useful for collinear libration points where (q1, p1) corresponds to the hyperbolic direction (left real) while (q2, p2) and (q3, p3) correspond to elliptic directions (complexified).

The transformation preserves the canonical structure and enables efficient normal form computations in the complex domain where elliptic motion becomes purely rotational.

See also

_M_inv()

Inverse transformation from complex to real coordinates.

_substitute_complex()

Apply complexification to polynomial expressions.

References

Meyer, K.R., Hall, G.R. (1992). Introduction to Hamiltonian Dynamical Systems and the N-Body Problem. Springer-Verlag, Section 4.3.

Return complexification transformation matrix for canonical coordinate pairs. Provides the matrix M that transforms real coordinates to complex coordinates for elliptic directions in the CR3BP.

hiten.algorithms.hamiltonian.transforms._M_inv(mix_pairs=(1, 2))[source]

Return inverse complexification matrix for real coordinate recovery.

Compute the inverse of the complexification matrix M, which transforms complex canonical coordinates back to real coordinates. Since M is unitary, the inverse is computed as the conjugate transpose.

Parameters:

mix_pairs (tuple of int, default (1, 2)) – Indices of canonical coordinate pairs that were complexified. Must match the mix_pairs used in the forward transformation.

Returns:

Inverse complexification matrix M_inv where real_coords = M_inv @ complex_coords. Computed as the conjugate transpose of M due to unitarity.

Return type:

ndarray, shape (6, 6), complex

Notes

The inverse transformation recovers real coordinates from complex coordinates: q_j(real) = Re(q_j(complex)) * sqrt(2) p_j(real) = Im(q_j(complex)) * sqrt(2)

This is used to convert results from complex normal form computations back to real coordinates for physical interpretation and integration.

See also

_M()

Forward complexification transformation.

_substitute_real()

Apply inverse complexification to polynomial expressions.

Return inverse complexification matrix for real coordinate recovery. Provides the inverse matrix M^{-1} that transforms complex coordinates back to real coordinates.

hiten.algorithms.hamiltonian.transforms._substitute_complex(poly_rn, max_deg, psi, clmo, tol=1e-12, *, mix_pairs=(1, 2))[source]

Transform polynomial from real to complex coordinates.

Apply complexification transformation to a polynomial expressed in real normal form coordinates, converting it to complex coordinates suitable for elliptic normal form analysis.

Parameters:
  • poly_rn (List[ndarray]) – Polynomial coefficients in real normal form coordinates, organized by degree, in nondimensional energy units.

  • max_deg (int) – Maximum polynomial degree to retain in the transformation.

  • psi (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • tol (float, default 1e-12) – Numerical tolerance for cleaning small coefficients after transformation.

  • mix_pairs (tuple of int, default (1, 2)) – Canonical coordinate pairs to complexify.

Returns:

Polynomial coefficients in complex coordinates, organized by degree, in nondimensional energy units.

Return type:

List[ndarray]

Notes

The transformation uses the complexification matrix M to convert polynomial expressions from real to complex coordinates. This enables efficient computation of normal forms around elliptic equilibria where the complex representation naturally captures rotational motion.

The transformation preserves polynomial structure and energy scaling, with complex coordinates maintaining the same physical units as their real counterparts.

See also

_substitute_real()

Inverse transformation from complex to real coordinates.

_M()

Complexification matrix used in the transformation.

Transform polynomial from real to complex coordinates. Converts polynomial coefficients from real modal coordinates to complex coordinates using the complexification transformation.

hiten.algorithms.hamiltonian.transforms._substitute_real(poly_cn, max_deg, psi, clmo, tol=1e-12, *, mix_pairs=(1, 2))[source]

Transform polynomial from complex to real coordinates.

Apply inverse complexification transformation to convert a polynomial from complex coordinates back to real normal form coordinates for physical interpretation and further analysis.

Parameters:
  • poly_cn (List[ndarray]) – Polynomial coefficients in complex coordinates, organized by degree, in nondimensional energy units.

  • max_deg (int) – Maximum polynomial degree to retain in the transformation.

  • psi (ndarray) – Combinatorial index table from _init_index_tables().

  • clmo (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • tol (float, default 1e-12) – Numerical tolerance for cleaning small coefficients after transformation.

  • mix_pairs (tuple of int, default (1, 2)) – Canonical coordinate pairs that were complexified.

Returns:

Polynomial coefficients in real coordinates, organized by degree, in nondimensional energy units.

Return type:

List[ndarray]

Notes

This is the inverse of _substitute_complex(), using the inverse complexification matrix M_inv to recover real polynomial expressions from their complex representations.

The transformation is essential for converting normal form results back to real coordinates for physical interpretation, trajectory integration, and comparison with numerical simulations.

See also

_substitute_complex()

Forward transformation from real to complex coordinates.

_M_inv()

Inverse complexification matrix used in the transformation.

Transform polynomial from complex to real coordinates. Converts polynomial coefficients from complex coordinates back to real modal coordinates using the inverse complexification transformation.

hiten.algorithms.hamiltonian.transforms._solve_complex(real_coords, tol=1e-30, *, mix_pairs=(1, 2))[source]

Transform real coordinates to complex coordinates.

Convert a real coordinate vector to complex coordinates using the complexification transformation, typically for elliptic normal form analysis.

Parameters:
  • real_coords (ndarray, shape (6,)) – Real canonical coordinates [q1, q2, q3, p1, p2, p3] in nondimensional position and momentum units.

  • tol (float, default 1e-30) – Tolerance for cleaning small imaginary parts in the result.

  • mix_pairs (tuple of int, default (1, 2)) – Canonical coordinate pairs to complexify.

Returns:

Complex coordinates [q1c, q2c, q3c, p1c, p2c, p3c] where complexified pairs are in complex units and unmixed pairs remain real.

Return type:

ndarray, shape (6,), complex

Notes

This transformation is used to convert initial conditions or intermediate results from real coordinates to the complex representation needed for normal form analysis around elliptic equilibria.

The transformation preserves the canonical structure and maintains appropriate scaling for further computations.

See also

_solve_real()

Inverse transformation from complex to real coordinates.

_M_inv()

Complexification matrix used in this transformation.

Transform real coordinates to complex coordinates. Applies the complexification transformation to coordinate vectors.

hiten.algorithms.hamiltonian.transforms._solve_real(real_coords, tol=1e-30, *, mix_pairs=(1, 2))[source]

Return real coordinates given complex coordinates using the map M.

Parameters:
  • real_coords (np.ndarray) – Real coordinates [q1, q2, q3, p1, p2, p3]

  • tol (float)

  • mix_pairs (tuple[int, ...])

Returns:

Real coordinates [q1r, q2r, q3r, p1r, p2r, p3r]

Return type:

np.ndarray

Transform complex coordinates to real coordinates. Applies the inverse complexification transformation to coordinate vectors.

hiten.algorithms.hamiltonian.transforms._polylocal2realmodal(point, poly_local, max_deg, psi, clmo, tol=1e-12)[source]

Transform a polynomial from local frame to real modal frame.

Parameters:
  • point (object) – An object with a normal_form_transform method that returns the transformation matrix

  • poly_phys (List[numpy.ndarray]) – Polynomial in physical coordinates

  • max_deg (int) – Maximum degree for polynomial representations

  • psi (numpy.ndarray) – Combinatorial table from _init_index_tables

  • clmo (numba.typed.List) – List of arrays containing packed multi-indices

  • poly_local (List[ndarray])

Returns:

Polynomial in real modal coordinates

Return type:

List[numpy.ndarray]

Notes

This function transforms a polynomial from local coordinates to real modal coordinates using the transformation matrix obtained from the point object.

Transform polynomial from local frame to real modal frame. Converts polynomial coefficients from local coordinates centered at the equilibrium point to real modal coordinates aligned with eigenvectors.

hiten.algorithms.hamiltonian.transforms._polyrealmodal2local(point, poly_realmodal, max_deg, psi, clmo, tol=1e-12)[source]

Transform a polynomial from real modal frame to local frame.

Parameters:
  • point (object) – An object with a normal_form_transform method that returns the transformation matrix

  • poly_realmodal (List[numpy.ndarray]) – Polynomial in real modal coordinates

  • max_deg (int) – Maximum degree for polynomial representations

  • psi (numpy.ndarray) – Combinatorial table from _init_index_tables

  • clmo (numba.typed.List) – List of arrays containing packed multi-indices

Returns:

Polynomial in local coordinates

Return type:

List[numpy.ndarray]

Notes

This function transforms a polynomial from real modal coordinates to local coordinates using the inverse of the transformation matrix obtained from the point object.

Transform polynomial from real modal frame to local frame. Converts polynomial coefficients from real modal coordinates back to local coordinates centered at the equilibrium point.

hiten.algorithms.hamiltonian.transforms._coordrealmodal2local(point, modal_coords, tol=1e-30)[source]

Transform coordinates from real modal to local frame.

Parameters:
  • point (object) – An object with a normal_form_transform method that returns the transformation matrix

  • modal_coords (np.ndarray) – Coordinates in real modal frame

Returns:

Coordinates in local frame

Return type:

np.ndarray

Notes

  • Modal coordinates are ordered as [q1, q2, q3, px1, px2, px3].

  • Local coordinates are ordered as [x1, x2, x3, px1, px2, px3].

Transform coordinates from real modal to local frame. Converts coordinate vectors from real modal coordinates to local coordinates.

hiten.algorithms.hamiltonian.transforms._coordlocal2realmodal(point, local_coords, tol=1e-30)[source]

Transform coordinates from local to real modal frame.

Parameters:
  • point (object) – An object with a normal_form_transform method that returns the transformation matrix

  • local_coords (np.ndarray) – Coordinates in local frame

Returns:

Coordinates in real modal frame

Return type:

np.ndarray

Notes

  • Local coordinates are ordered as [x1, x2, x3, px1, px2, px3].

  • Modal coordinates are ordered as [q1, q2, q3, px1, px2, px3].

Transform coordinates from local to real modal frame. Converts coordinate vectors from local coordinates to real modal coordinates.

hiten.algorithms.hamiltonian.transforms._local2synodic_collinear(point, local_coords, tol=1e-14)[source]

Transform coordinates from local to synodic frame for collinear points.

Convert coordinates from the local frame centered at a collinear equilibrium point to the standard CR3BP synodic rotating frame coordinates.

Parameters:
  • point (CollinearPoint) – Collinear libration point object providing geometric parameters gamma, mu, sign, and a for the coordinate transformation.

  • local_coords (ndarray, shape (6,)) – Local coordinates [x1, x2, x3, px1, px2, px3] in nondimensional units where positions are in distance units and momenta are canonical.

  • tol (float, default 1e-14) – Tolerance for detecting non-negligible imaginary parts in input coordinates.

Returns:

Synodic coordinates [X, Y, Z, Vx, Vy, Vz] in nondimensional CR3BP units where positions are distances and velocities are time derivatives.

Return type:

ndarray, shape (6,)

Raises:

ValueError – If local_coords is not a flat array of length 6 or contains imaginary parts larger than the specified tolerance.

Notes

The transformation implements the mapping from local coordinates centered at the equilibrium point to the standard synodic frame where the primaries are located at (-mu, 0, 0) and (1-mu, 0, 0).

The coordinate transformation includes: - Position scaling by the characteristic length gamma - Translation by the equilibrium point offset (mu + a) - Momentum to velocity conversion with Coriolis terms - Sign convention adjustments for NASA/Szebehely compatibility

Local frame: centered at equilibrium, canonical coordinates Synodic frame: standard CR3BP rotating frame, Cartesian coordinates

See also

_synodic2local_collinear()

Inverse transformation from synodic to local coordinates.

CollinearPoint

Collinear point class providing transformation parameters.

References

Szebehely, V. (1967). Theory of Orbits. Academic Press, Chapter 7.

Transform coordinates from local to synodic frame for collinear points. Converts coordinates from the local frame centered at a collinear libration point to the synodic rotating frame coordinates.

hiten.algorithms.hamiltonian.transforms._synodic2local_collinear(point, synodic_coords, tol=1e-14)[source]

Transform coordinates from synodic to local frame for the collinear points.

This is the exact inverse of _local2synodic_collinear().

Parameters:
  • point (CollinearPoint) – Collinear libration point providing the geometric parameters gamma, mu, sign and a.

  • synodic_coords (np.ndarray) – Coordinates in synodic frame [X, Y, Z, Vx, Vy, Vz].

Returns:

Coordinates in local frame [x1, x2, x3, px1, px2, px3].

Return type:

np.ndarray

Raises:

ValueError – If synodic_coords is not a flat array of length 6 or contains an imaginary part larger than the tolerance (1e-16).

Transform coordinates from synodic to local frame for collinear points. Converts coordinates from the synodic rotating frame to the local frame centered at a collinear libration point.

hiten.algorithms.hamiltonian.transforms._local2synodic_triangular(point, local_coords, tol=1e-14)[source]

Transform coordinates from local to synodic frame for the equilateral points.

Parameters:
  • point (object) – An object with a normal_form_transform method that returns the transformation matrix

  • local_coords (np.ndarray) – Coordinates in local frame

Returns:

Coordinates in synodic frame

Return type:

np.ndarray

Notes

  • Local coordinates are ordered as [x1, x2, x3, px1, px2, px3].

  • Synodic coordinates are ordered as [X, Y, Z, Vx, Vy, Vz].

Raises:

ValueError – If local_coords is not a flat array of length 6 or contains an imaginary part larger than the tolerance (1e-16).

Parameters:
Return type:

np.ndarray

Transform coordinates from local to synodic frame for triangular points. Converts coordinates from the local frame centered at a triangular libration point to the synodic rotating frame coordinates.

hiten.algorithms.hamiltonian.transforms._synodic2local_triangular(point, synodic_coords, tol=1e-14)[source]

Transform coordinates from synodic to local frame for the triangular (equilateral) points.

This is the exact inverse of _local2synodic_triangular().

Parameters:
  • point (TriangularPoint) – Triangular libration point providing the geometric parameters mu and sign.

  • synodic_coords (np.ndarray) – Coordinates in synodic frame [X, Y, Z, Vx, Vy, Vz].

Returns:

Coordinates in local frame [x1, x2, x3, px1, px2, px3].

Return type:

np.ndarray

Raises:

ValueError – If synodic_coords is not a flat array of length 6 or contains an imaginary part larger than the tolerance (1e-16).

Transform coordinates from synodic to local frame for triangular points. Converts coordinates from the synodic rotating frame to the local frame centered at a triangular libration point.

hiten.algorithms.hamiltonian.transforms._restrict_poly_to_center_manifold(point, poly_H, clmo, tol=1e-14)[source]

Restrict polynomial Hamiltonian to center manifold by eliminating hyperbolic terms.

Project a polynomial Hamiltonian onto the center manifold by removing all terms that depend on hyperbolic variables, retaining only the dynamics within the center-stable/center-unstable subspace.

Parameters:
  • point (CollinearPoint or TriangularPoint) – Libration point object determining the manifold structure.

  • poly_H (List[ndarray]) – Polynomial Hamiltonian coefficients organized by degree, in nondimensional energy units.

  • clmo (List[ndarray]) – Packed multi-index table from _init_index_tables().

  • tol (float, default 1e-14) – Tolerance for zeroing small coefficients during restriction.

Returns:

Restricted polynomial with hyperbolic terms eliminated, maintaining the same degree structure as the input.

Return type:

List[ndarray]

Notes

For collinear points, the first canonical pair (q1, p1) corresponds to the hyperbolic direction, so all terms with non-zero exponents in these variables are eliminated. This leaves only the center manifold dynamics in the (q2, p2, q3, p3) subspace.

For triangular points, all directions are elliptic (center-type), so the function returns a copy of the original polynomial without restriction.

This restriction is fundamental to center manifold theory, which reduces the dimensionality of the dynamics by focusing on the neutrally stable directions while eliminating the exponentially growing/decaying modes.

See also

_decode_multiindex()

Multi-index decoding used to identify hyperbolic terms.

CollinearPoint

Collinear points with hyperbolic directions.

TriangularPoint

Triangular points with all elliptic directions.

References

Carr, J. (1981). Applications of Centre Manifold Theory. Springer-Verlag.

Jorba, A., Masdemont, J. (1999). Dynamics in the center manifold of the collinear points. Physica D, 132(1-2), 189-213.

Restrict polynomial Hamiltonian to center manifold by eliminating hyperbolic terms. Removes all terms that depend on hyperbolic variables, effectively restricting the Hamiltonian to the center-stable manifold subspace.

wrappers.py

The wrappers module provides a registry-based conversion system for automatic Hamiltonian transformations.

hiten.algorithms.hamiltonian.wrappers.register_conversion(src_name, dst, required_context=None, default_params=None)[source]

Decorator to register Hamiltonian coordinate system conversions.

Register a conversion function in the global conversion registry, enabling automatic transformation between different Hamiltonian representations in the CR3BP normal form pipeline.

Parameters:
  • src_name (str) – Name identifier of the source Hamiltonian representation (e.g., “physical”, “real_modal”, “complex_modal”).

  • dst (type[Hamiltonian] or str) – Either the destination Hamiltonian class or its name string identifier. Using strings allows registration before class definition.

  • required_context (list, optional) – List of required context keys that must be provided during conversion. Common requirements include [“point”] for libration point information.

  • default_params (dict, optional) – Default parameter values for the conversion function (e.g., tolerances, numerical settings).

Returns:

Decorator function that registers the conversion and returns the original function unchanged.

Return type:

callable

Notes

The conversion registry enables automatic transformation between coordinate systems by maintaining a lookup table of (source, destination) to function mappings. This supports the normal form computation pipeline where Hamiltonians must be transformed through multiple coordinate systems.

Registered functions should follow the signature: conversion_func(ham: Hamiltonian, kwargs) returns Hamiltonian or tuple

The kwargs typically include context (like libration points) and numerical parameters (like tolerances). Some conversions return tuples containing both the transformed Hamiltonian and auxiliary data like generating functions.

See also

Hamiltonian

Base Hamiltonian class used in conversions.

_CONVERSION_REGISTRY()

Global registry storing conversion functions.

Examples

>>> @register_conversion("physical", "real_modal",
...                     required_context=["point"],
...                     default_params={"tol": 1e-12})
... def physical_to_modal(ham: Hamiltonian, kwargs) -> Hamiltonian:
...     point = kwargs["point"]
...     tol = kwargs.get("tol", 1e-12)
...     # Transform polynomial using point's eigenvectors
...     return transformed_hamiltonian

Register a conversion function in the Hamiltonian conversion registry. Decorator function that registers conversion functions between different Hamiltonian representations in the automatic conversion system.

hiten.algorithms.hamiltonian.wrappers._physical_to_real_modal(ham, **kwargs)[source]

Transform Hamiltonian from physical to real modal coordinates.

Convert a polynomial Hamiltonian from local physical coordinates centered at an equilibrium point to real modal coordinates aligned with the linear stability eigenvectors of the equilibrium.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in physical coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point providing the modal transformation matrix.

    • tolfloat, default 1e-12

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in real modal coordinates with name “real_modal”.

Return type:

Hamiltonian

Notes

The transformation uses the eigenvector matrix from the libration point’s linearization to convert from local coordinates to modal coordinates where each coordinate pair corresponds to a specific eigenvalue/eigenvector pair of the linearized dynamics.

This transformation is the first step in the normal form pipeline, enabling subsequent complexification and normal form reductions by working in the natural coordinate system defined by the linear stability properties.

See also

_polylocal2realmodal()

Underlying transformation function.

_real_modal_to_physical()

Inverse transformation back to physical coordinates.

Transform Hamiltonian from physical to real modal coordinates. Converts the physical Hamiltonian from local coordinates to real modal coordinates aligned with the linear stability eigenvectors.

hiten.algorithms.hamiltonian.wrappers._real_modal_to_physical(ham, **kwargs)[source]

Transform Hamiltonian from real modal to physical coordinates.

Convert a polynomial Hamiltonian from real modal coordinates aligned with linear stability eigenvectors back to local physical coordinates centered at the equilibrium point.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in real modal coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point providing the inverse modal transformation matrix.

    • tolfloat, default 1e-12

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in physical coordinates with name “physical”.

Return type:

Hamiltonian

Notes

This transformation is the inverse of the physical to real modal conversion, using the inverse of the eigenvector matrix from the libration point’s linearization to convert back from modal coordinates to local coordinates.

The transformation restores the original physical coordinate system where coordinates represent displacements from the equilibrium point in the standard CR3BP reference frame.

See also

_polyrealmodal2local()

Underlying transformation function.

_physical_to_real_modal()

Forward transformation to real modal coordinates.

Transform Hamiltonian from real modal to physical coordinates. Converts the Hamiltonian from real modal coordinates back to physical local coordinates centered at the equilibrium point.

hiten.algorithms.hamiltonian.wrappers._real_modal_to_complex_modal(ham, **kwargs)[source]

Transform Hamiltonian from real modal to complex modal coordinates.

Convert a polynomial Hamiltonian from real modal coordinates to complex modal coordinates by complexifying elliptic coordinate pairs. This enables the use of complex normal form techniques for analyzing the dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in real modal coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to complexify.

    • tolfloat, default 1e-12

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in complex modal coordinates with name “complex_modal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are complexified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are complexified since all directions are elliptic.

The complexification process introduces complex variables that simplify the analysis of elliptic dynamics and enable the application of complex normal form theory. The transformation preserves the Hamiltonian structure while providing a more convenient representation for subsequent normal form computations.

See also

_substitute_complex()

Underlying complexification function.

_complex_modal_to_real_modal()

Inverse transformation back to real modal coordinates.

Transform Hamiltonian from real modal to complex modal coordinates. Converts the Hamiltonian from real modal coordinates to complex modal coordinates for elliptic directions using complexification.

hiten.algorithms.hamiltonian.wrappers._complex_modal_to_real_modal(ham, **kwargs)[source]

Transform Hamiltonian from complex modal to real modal coordinates.

Convert a polynomial Hamiltonian from complex modal coordinates back to real modal coordinates by realifying complexified elliptic coordinate pairs. This provides the inverse of the complexification process.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex modal coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to realify.

    • tolfloat, default 1e-12

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in real modal coordinates with name “real_modal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are realified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are realified since all directions are elliptic.

The realification process converts complex variables back to real variables while preserving the Hamiltonian structure. This transformation is typically used after complex normal form computations to return to a real representation suitable for physical interpretation or further real-coordinate analysis.

See also

_substitute_real()

Underlying realification function.

_real_modal_to_complex_modal()

Forward transformation to complex modal coordinates.

Transform Hamiltonian from complex modal to real modal coordinates. Converts the Hamiltonian from complex modal coordinates back to real modal coordinates using inverse complexification.

hiten.algorithms.hamiltonian.wrappers._complex_modal_to_complex_partial_normal(ham, **kwargs)[source]

Transform Hamiltonian to partial normal form via Lie series method.

Apply partial normal form transformation to eliminate non-resonant terms from a complex modal Hamiltonian, retaining only terms that contribute to the center manifold dynamics and specific resonances.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex modal coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point providing eigenvalue information for resonance analysis.

    • tol_liefloat, default 1e-30

      Numerical tolerance for Lie series computations and coefficient cleaning.

Returns:

  • Transformed Hamiltonian in partial normal form with name “complex_partial_normal”

  • Generating functions used in the transformation, containing both the total generating function and eliminated terms

Return type:

tuple of (Hamiltonian, LieGeneratingFunction)

Notes

The partial normal form eliminates terms that do not contribute to the center manifold dynamics while preserving resonant terms and the structure needed for center manifold analysis. This is achieved through a sequence of canonical transformations generated by polynomial generating functions.

The Lie series method implements transformations of the form: H’ = exp(L_G) H where L_G is the Lie operator associated with generating function G.

The returned generating functions contain the complete transformation history and can be used for coordinate transformations or analysis of the eliminated dynamics.

See also

_lie_transform()

Underlying partial normal form computation.

_complex_modal_to_complex_full_normal()

Complete normal form transformation.

LieGeneratingFunction

Container for generating function data.

Transform Hamiltonian to partial normal form via Lie series method. Applies partial Lie series normalization to eliminate non-resonant terms while preserving resonant terms for center manifold analysis.

hiten.algorithms.hamiltonian.wrappers._complex_partial_normal_to_real_partial_normal(ham, **kwargs)[source]

Transform Hamiltonian from complex partial normal to real partial normal form.

Convert a polynomial Hamiltonian from complex partial normal form to real partial normal form by realifying complexified elliptic coordinate pairs. This provides a real representation of the partially normalized dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex partial normal form, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to realify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in real partial normal form with name “real_partial_normal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are realified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are realified since all directions are elliptic.

The partial normal form retains only terms that contribute to center manifold dynamics while eliminating non-resonant terms. The realification process converts the complex representation back to real variables while preserving the normalized structure, making the dynamics more interpretable in physical terms.

This transformation is typically used after complex partial normal form computations to obtain a real representation suitable for center manifold analysis or further real-coordinate processing.

See also

_substitute_real()

Underlying realification function.

_real_partial_normal_to_complex_partial_normal()

Inverse transformation back to complex partial normal form.

Transform Hamiltonian from complex partial normal to real partial normal form. Converts the complex partial normal form back to real coordinates while maintaining the partial normalization structure.

hiten.algorithms.hamiltonian.wrappers._real_partial_normal_to_complex_partial_normal(ham, **kwargs)[source]

Transform Hamiltonian from real partial normal to complex partial normal form.

Convert a polynomial Hamiltonian from real partial normal form to complex partial normal form by complexifying elliptic coordinate pairs. This enables the use of complex analysis techniques on the partially normalized dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in real partial normal form, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to complexify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in complex partial normal form with name “complex_partial_normal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are complexified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are complexified since all directions are elliptic.

The partial normal form retains only terms that contribute to center manifold dynamics while eliminating non-resonant terms. The complexification process converts real variables to complex variables while preserving the normalized structure, enabling the application of complex analysis techniques.

This transformation is typically used when complex analysis methods are needed on the partially normalized dynamics, such as for further normal form computations or complex center manifold analysis.

See also

_substitute_complex()

Underlying complexification function.

_complex_partial_normal_to_real_partial_normal()

Inverse transformation back to real partial normal form.

Transform Hamiltonian from real partial normal to complex partial normal form. Converts the real partial normal form to complex coordinates for further processing in the normalization pipeline.

hiten.algorithms.hamiltonian.wrappers._complex_partial_normal_to_center_manifold_complex(ham, **kwargs)[source]

Restrict Hamiltonian to center manifold by eliminating hyperbolic terms.

Project a partial normal form Hamiltonian onto the center manifold by removing all terms that depend on hyperbolic variables, retaining only the dynamics within the center-stable/center-unstable subspace.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex partial normal form, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining the manifold structure and hyperbolic directions.

    • tolfloat, default 1e-14

      Numerical tolerance for zeroing small coefficients during restriction.

Returns:

Restricted Hamiltonian on center manifold with name “center_manifold_complex”.

Return type:

Hamiltonian

Notes

For collinear points, the first canonical pair corresponds to the hyperbolic direction, so all terms with non-zero exponents in these variables are eliminated. This reduces the 6D phase space to a 4D center manifold.

For triangular points, all directions are elliptic (center-type), so the function returns the original Hamiltonian without restriction.

This restriction is fundamental to center manifold theory, which enables dimensional reduction by focusing on neutrally stable directions while eliminating exponentially growing/decaying hyperbolic modes.

The resulting center manifold Hamiltonian captures the long-term dynamics and is used for constructing periodic orbits, invariant tori, and other dynamical structures that persist in the full system.

See also

_restrict_poly_to_center_manifold()

Underlying center manifold restriction function.

_center_manifold_complex_to_center_manifold_real()

Conversion to real center manifold coordinates.

References

Carr, J. (1981). Applications of Centre Manifold Theory. Springer-Verlag.

Restrict Hamiltonian to center manifold by eliminating hyperbolic terms. Removes all terms that depend on hyperbolic variables, effectively restricting the Hamiltonian to the center-stable manifold subspace.

hiten.algorithms.hamiltonian.wrappers._center_manifold_complex_to_center_manifold_real(ham, **kwargs)[source]

Transform Hamiltonian from complex center manifold to real center manifold.

Convert a polynomial Hamiltonian from complex center manifold representation to real center manifold representation by realifying complexified elliptic coordinate pairs. This provides a real representation of the center manifold dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex center manifold representation, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to realify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in real center manifold representation with name “center_manifold_real”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are realified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are realified since all directions are elliptic.

The center manifold Hamiltonian captures the long-term dynamics by eliminating hyperbolic terms and focusing on neutrally stable directions. The realification process converts complex variables back to real variables while preserving the center manifold structure, making the dynamics more interpretable in physical terms.

This transformation is typically used after complex center manifold computations to obtain a real representation suitable for physical interpretation or further real-coordinate analysis of the center manifold dynamics.

See also

_substitute_real()

Underlying realification function.

_center_manifold_real_to_center_manifold_complex()

Inverse transformation back to complex center manifold representation.

Transform center manifold Hamiltonian from complex to real coordinates. Converts the center manifold Hamiltonian from complex coordinates back to real coordinates while maintaining the center manifold restriction.

hiten.algorithms.hamiltonian.wrappers._center_manifold_real_to_center_manifold_complex(ham, **kwargs)[source]

Transform Hamiltonian from real center manifold to complex center manifold.

Convert a polynomial Hamiltonian from real center manifold representation to complex center manifold representation by complexifying elliptic coordinate pairs. This enables the use of complex analysis techniques on center manifold dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in real center manifold representation, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to complexify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in complex center manifold representation with name “center_manifold_complex”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are complexified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are complexified since all directions are elliptic.

The center manifold Hamiltonian captures the long-term dynamics by eliminating hyperbolic terms and focusing on neutrally stable directions. The complexification process converts real variables to complex variables while preserving the center manifold structure, enabling the application of complex analysis techniques.

This transformation is typically used when complex analysis methods are needed on the center manifold dynamics, such as for further normal form computations or complex center manifold analysis.

See also

_substitute_complex()

Underlying complexification function.

_center_manifold_complex_to_center_manifold_real()

Inverse transformation back to real center manifold representation.

Transform center manifold Hamiltonian from real to complex coordinates. Converts the center manifold Hamiltonian from real coordinates to complex coordinates for further analysis or processing.

hiten.algorithms.hamiltonian.wrappers._complex_modal_to_complex_full_normal(ham, **kwargs)[source]

Transform Hamiltonian to full normal form via Lie series method.

Apply complete normal form transformation to eliminate all non-resonant terms from a complex modal Hamiltonian, achieving the maximum possible simplification of the dynamics through canonical transformations.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex modal coordinates, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or

      TriangularPoint Libration point providing eigenvalue information for resonance analysis.

    • tol_liefloat, default 1e-30

      Numerical tolerance for Lie series computations and coefficient cleaning.

    • resonance_tolfloat, default 1e-14

      Tolerance for detecting resonance conditions between eigenvalues.

Returns:

  • Transformed Hamiltonian in full normal form with name “complex_full_normal”

  • Generating functions used in the transformation, containing both the total generating function and eliminated terms

Return type:

tuple of (Hamiltonian, LieGeneratingFunction)

Notes

The full normal form eliminates all non-resonant terms through a sequence of canonical transformations, achieving the maximum possible simplification of the Hamiltonian while preserving the essential dynamics. This is the most complete form of normalization possible.

The Lie series method implements transformations of the form: H’ = exp(L_G) H where L_G is the Lie operator associated with generating function G.

Resonance conditions are detected based on the eigenvalues of the libration point, and only resonant terms are preserved in the final normal form. The resulting Hamiltonian contains only terms that cannot be eliminated due to resonance conditions.

The returned generating functions contain the complete transformation history and can be used for coordinate transformations or analysis of the eliminated dynamics.

See also

_lie_transform()

Underlying full normal form computation.

_complex_modal_to_complex_partial_normal()

Partial normal form transformation.

LieGeneratingFunction

Container for generating function data.

Transform Hamiltonian to full normal form via complete Lie series method. Applies full Lie series normalization to eliminate all non-resonant terms, producing the maximally simplified canonical form.

hiten.algorithms.hamiltonian.wrappers._complex_full_normal_to_real_full_normal(ham, **kwargs)[source]

Transform Hamiltonian from complex full normal to real full normal form.

Convert a polynomial Hamiltonian from complex full normal form to real full normal form by realifying complexified elliptic coordinate pairs. This provides a real representation of the completely normalized dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in complex full normal form, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or

      TriangularPoint Libration point determining which coordinate pairs to realify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in real full normal form with name “real_full_normal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are realified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are realified since all directions are elliptic.

The full normal form represents the maximum possible simplification of the Hamiltonian through canonical transformations, eliminating all non-resonant terms. The realification process converts complex variables back to real variables while preserving the normalized structure, making the dynamics more interpretable in physical terms.

This transformation is typically used after complex full normal form computations to obtain a real representation suitable for physical interpretation or further real-coordinate analysis of the completely normalized dynamics.

See also

_substitute_real()

Underlying realification function.

_real_full_normal_to_complex_full_normal()

Inverse transformation back to complex full normal form.

Transform Hamiltonian from complex full normal to real full normal form. Converts the complex full normal form back to real coordinates while maintaining the complete normalization structure.

hiten.algorithms.hamiltonian.wrappers._real_full_normal_to_complex_full_normal(ham, **kwargs)[source]

Transform Hamiltonian from real full normal to complex full normal form.

Convert a polynomial Hamiltonian from real full normal form to complex full normal form by complexifying elliptic coordinate pairs. This enables the use of complex analysis techniques on the completely normalized dynamics.

Parameters:
  • ham (Hamiltonian) – Source Hamiltonian in real full normal form, with polynomial coefficients in nondimensional energy units.

  • kwargs

    Conversion context and parameters:

    • pointCollinearPoint or TriangularPoint

      Libration point determining which coordinate pairs to complexify.

    • tolfloat, default 1e-14

      Numerical tolerance for cleaning small coefficients.

Returns:

Transformed Hamiltonian in complex full normal form with name “complex_full_normal”.

Return type:

Hamiltonian

Notes

For collinear points, coordinate pairs (1, 2) are complexified, corresponding to the elliptic directions. For triangular points, all coordinate pairs (0, 1, 2) are complexified since all directions are elliptic.

The full normal form represents the maximum possible simplification of the Hamiltonian through canonical transformations, eliminating all non-resonant terms. The complexification process converts real variables to complex variables while preserving the normalized structure, enabling the application of complex analysis techniques.

This transformation is typically used when complex analysis methods are needed on the completely normalized dynamics, such as for further analysis or complex coordinate transformations.

See also

_substitute_complex()

Underlying complexification function.

_complex_full_normal_to_real_full_normal()

Inverse transformation back to real full normal form.

Transform Hamiltonian from real full normal to complex full normal form. Converts the real full normal form to complex coordinates for further analysis or processing in the normalization pipeline.

center/

The center module provides partial normal form computations for center manifold analysis.

normal/

The normal module provides full normal form computations for complete dynamical reduction.