Hamiltonian system is a mechanical system used to describe physical processes and phenomena. The mathematical framework of Hamiltonian mechanics is symplectic geometry. The study of symplectic geometry can be traced back to Hamilton, the great British mathematician and physicist in the 19th century, its form is:

H(z) is Hamilton function and J is 2n order antisymmetric matrix. In the study of Newtonian mechanics, he introduced generalized coordinates and generalized momentum, namely Hamiltonian function, to represent the energy of the system. Symplectic structure is a basic characteristic of Hamiltonian system. Its solution is symplectic transformation and the canonical equation of Hamiltonian system is invariant under symplectic transformation. From this point of view, the calculation method of symplectic preserving structure of discrete equation of Hamilton system is studied. Hamilton system is the starting point of physics research because of its unique formal symmetry. Hamilton equation includes finite dimension and infinite dimension. They are special forms of ordinary differential equations or partial differential equations, which play a very important role in the research of modern physics and mechanics. Therefore, many scientists have great interest in this kind of equation and are committed to studying its numerical solution. However, traditional numerical methods often destroy the symplectic structure of Hamiltonian system, so that the numerical simulation cannot remain stable for a long time. In some practical problems, such as the calculation and Simulation of long-term propagation of some geophysical wave fields (such as earth self-excited oscillation, seismic wave, etc.), weather prediction, complex hydrodynamics, electromagnetism have become hot spots, which requires the selection of more efficient and stable long-time numerical algorithms with simulation ability and the calculation scale is becoming larger and larger, The corresponding calculation time is longer and longer and the requirements for numerical calculation methods are also gradually improved. Therefore, it is of great practical significance to study the algorithm to maintain the symplectic structure of Hamiltonian system. The algorithm that maintains the symplectic structure of Hamiltonian system is called symplectic algorithm. Symplectic algorithm is a difference method based on the basic principle of Hamiltonian mechanics [1-2].

In the early 1980s, Chinese scholar Feng Kang and American scientist Ruth proposed a numerical integration algorithm that can maintain the symplectic structure of Hamiltonian phase flow almost at the same time, which solves the problem that the traditional algorithm cannot maintain the energy conservation of the system due to long-term time integration. The former mainly uses the implicit midpoint method to construct the high-order implicit symplectic algorithm for the inseparable Hamiltonian system and the latter mainly establishes the explicit symplectic algorithm for the Hamiltonian system which can be decomposed into kinetic energy T and potential energy V. If an algorithm can maintain one or more physical or geometric structures, such as first integral, symmetry, symplectic structure, phase space volume, etc., the algorithm is called structure preserving algorithm. In addition, Mr. Feng Kang systematically developed a set of symplectic geometric algorithm (symplectic algorithm) for solving Hamiltonian system for the first time, systematically described the generating function method of constructing any order symplectic difference scheme from generating function, constructively gave the symplectic difference scheme to obtain any order accuracy and extended the theory and method of one-dimensional regular Hamiltonian system to infinite dimensional case [3].

Because the traditional algorithms are not symplectic except for a few, most of them will produce artificial dissipation However, many algorithms proposed by Feng Kang can maintain the architecture and are obviously superior to traditional algorithms in many aspects. It can be said that the birth and development of symplectic algorithm has great theoretical and practical significance for the whole field of computational mathematics. It promotes the establishment of the concept of more general structure preserving algorithm and the research in related application fields. In fact, in addition to the general geometric structure (symplectic structure, unitary structure, contact structure), the structure preserving algorithm proposed by Feng Kang also covers other conservative properties of conservative systems, such as physical conservation laws (first integration, momentum, energy, etc.) and algebraic properties (system eigenvalues, Lie algebraic structure, differential complex, etc.).

Differential equations can be used to describe some physical phenomena in science and engineering. However, for most partial differential equations, it is very difficult to obtain their exact solutions. Therefore, numerical solution of partial differential equations is a very important method and it is particularly important to pursue more efficient and stable numerical methods. In recent decades, many scholars at home and abroad have carried out algorithm theory and application research in this field, proposed a series of algorithm construction methods and theories and achieved a lot of fruitful research results. The research shows that the architecture and conservation law of the discrete Hashimoto algorithm are completely parallel, highly close to the Hashimoto prototype and has theoretically infinite long-term tracking ability, which fundamentally ensures and explains the uniqueness and superiority of the Hashimoto algorithm. How can discretization approach the Hartmann prototype highly? A correct way of calculation is to hamiltonize the equation first and then solve it by Hamiltonian algorithm, so that the basic characteristics of the problem prototype should be maintained as much as possible after discretization. Therefore, discretization should be carried out in the same formal framework of the problem prototype as far as possible, so the difference scheme used for discretization should be symplectic, that is, the mapping from the previous time to the next time is symplectic [4].

The initial stage of symplectic algorithm is mainly for ordinary differential equations. The practical problems are generally complex and are mostly simulated by nonlinear partial differential equations. Therefore, compared with ordinary differential systems, infinite dimensional partial differential systems have a wider range of applications. The purpose of this paper is to review the research of structure preserving algorithms and summarize the existing methods and conclusions. The next section is devoted to the summary of related research problems.

Summary of Related Research Problems

Example 1 (Seismic wave equation): Seismic wave numerical simulation is an effective means to study the propagation characteristics of seismic wave, which is helpful to understand the propagation law of seismic wave in complex media. It is an important part of seismic exploration and seismology. Seismic wave is an elastic wave propagating in rock stratum. In essence, it is the propagation process of energy. What affects the energy propagation is only the nature of the medium. The process of seismic wave propagation is essentially a process in which energy is gradually lost until it is exhausted. In practical application, it is often described by elastic wave equation or scalar wave equation under the assumption of no energy loss. In the Hamiltonian system, the seismic wave propagation process is the evolution process of an infinite dimensional Hamiltonian system with time. If the energy loss is not considered, the wave field evolution process is essentially a single parameter continuous symplectic transformation, so the corresponding numerical algorithm should be symplectic geometry algorithm.

 The process of seismic wave propagation is very complex and its analytical solution is almost impossible. Generally, it can only be realized by numerical calculation. Among various seismic wave numerical simulation methods, finite difference method and finite element method are the earliest and most widely used methods. The finite difference method is used to calculate the elastic wave equation in homogeneous media. Subsequently, the method has been greatly developed in calculation speed, accuracy and applicability and has been widely used to calculate wave equations and elastic wave equations in complex media, including seismic waves. At present, it has become a relatively mature algorithm. The seismic wave is described in the Hamiltonian system. Under the condition of no energy loss, the seismic wave propagation process described by the wave equation is essentially a single parameter symplectic transformation, so the symplectic algorithm is suitable for numerical calculation.

The propagation of seismic wave in medium can be described by the following three-dimensional scalar wave equation:

In which V represents the velocity of seismic wave propagation in the medium and u represents the seismic wave field. After discretizing the wave field in time and phase space, some symplectic schemes for wave field calculation are derived in Luo et al. [5] and the similarities and differences between finite difference scheme and symplectic scheme are compared and discussed. Through wave field numerical calculation, it can be seen that the results of symplectic algorithm and finite difference method are almost the same when the time step is small. In Luo et al. [6], the spectral factorization method under the helical boundary condition is adopted. Under the helical boundary condition, the inverse matrix needs to be transformed into a banded matrix and the position and size of non-zero elements in each column are very similar. Therefore, the spectral factorization method can be used to realize the fast decomposition. In this paper, the implicit frog leaping symplectic scheme and spectral factorization method with second-order accuracy are adopted, the wave field in the layered medium model with constant velocity is calculated. The calculation shows that the implicit symplectic algorithm is a good method for wave field calculation. Based on the spectral element method for spatial domain discretization and the newly derived third-order symplectic algorithm for time domain discretization, a new algorithm with space-time structure preserving characteristics is constructed in Wang et al. [7]. The comparison results of several groups of numerical experiments given in this paper show that the algorithm has good performance in terms of memory consumption, stability, computing time and long-term tracking ability; In addition, the numerical examples of multi-layer medium model on undulating surface show the effectiveness of the algorithm in dealing with complex geometry and complex medium. The development of the multi symplectic structure spectral element method will provide a more extensive and effective choice for the calculation and Simulation of long-term seismic wave propagation.

Example 2 (Duffing Equation) 

Under the action of periodic external force, Duffing equation can be written as:

Where, F and are the amplitude and frequency of external force respectively.

Aiming at the problem of increasing the amount of calculation due to matrix inversion in the implementation of symplectic precise integration method, the authors in
Du and Hou [8] uses the idea of piecewise differential to approximate homogenize the non-homogeneous differential equation, so that the inversion of coefficient matrix of unsteady term in symplectic precise integration scheme does not contain time variables and reduces the amount of calculation of matrix inversion, This method is applied to the numerical analysis of undamped Duffing equation. The numerical results show that the symplectic precise integration method after approximate homogenization is obviously superior to the classical Runge Kutta method in numerical accuracy, calculation time and maintaining system energy; At the same time, the existing precise integration method is significantly improved in maintaining the system energy. In Gang et al. [9], the dissipation term and forcing term of Duffing equation are regarded as the perturbation of nonlinear Hamiltonian equation and the perturbation effect is applied to the parameters ε express. Firstly, based on the gradient Hamiltonian decomposition theory of vector field, the structure preserving algorithm of Duffing equation is constructed by using the splitting method. Then, according to Liouville formula, it is proved that the determinant of Jacobian matrix of structure preserving algorithm is equal to the determinant of exact flow of Duffing equation. However, considering the explicit Runge Kutta method, it is found that there is a p+1 order error term in the determinant of Jacobian matrix. The volume evolution laws of different algorithms in a given region in phase space are discussed respectively. The results show that the sum of Lyapunov exponents is completely unchanged for the structure preserving algorithm proposed in this paper. Finally, through numerical experiments, the relative norm error and absolute energy error of phase trajectory of structure preserving algorithm and Heun method (second-order Runge Kutta method) are compared. The calculation results show that when ε When less than or equal to zero, the structure preserving algorithm is obviously better than Heun method.

Example 3 (Structure-Preserving Neural Networks)

With the emergence of the so-called fourth paradigm of science, people are more and more interested in machine learning of scientific laws. A large number of methods have been developed. By using the technology from classical regression to the most complex deep learning methods, these methods can more or less accurately predict the response of physical systems in previously invisible situations. The purpose in Hernández et al. [10] is to develop a new structure preserving neural network structure, which can predict the time evolution of the system based on the experimental observation of the system without knowing its governing equations in advance and is suitable for conservative and dissipative systems. The key idea is to combine the proven computing power of neural network in highly nonlinear physics with the data-driven algorithm consistent with thermodynamics. The resulting method, as will be seen, is a powerful neural network architecture, conceptually very simple, based on the standard feedforward method, uses the correct thermodynamic structure of the system revealed from the experimental data and produces interpretable results.

Example 4 (Sine-Gordon Equation)

Sine Gordon equation is a partial differential equation discovered in the 19th century. It has attracted much attention because of its many soliton solutions. In Jiang et al. [11], a new linear implicit local energy preserving scheme for Sine-Gordon equation is proposed. The basic idea is to construct the energy stable scheme of gradient system, namely energy equilibrium, from the invariant energy square method. Taking Sine-Gordon equation as an example, it shows that the invariant energy quadrature method is also an effective method to construct linear implicit and local energy conservation schemes for energy conservation systems. Using the invariant energy quadrature method, the Sine-Gordon equation is transformed into an equivalent system. Then the new system is discretized by the traditional finite difference method and a semi discrete system is obtained, which can maintain the semi discrete local energy conservation law. Then, the linear implicit structure preserving method is applied to the obtained semi discrete system to obtain a fully discretized scheme. It is proved that the obtained scheme can accurately maintain the discrete local energy conservation law. The theoretical analysis is verified by various numerical examples and the superiority of the new scheme over the existing local structure preserving scheme is proved.

In Cai et al. [12], two strategies for constructing homogeneous Neumann boundary condition structure preserving algorithms for sine-Gordon equation are proposed and most of the existing structure preserving algorithms are only suitable for zero or periodic boundary conditions. The first strategy is based on the traditional second-order central difference quotient but uses the cell centered grid, while the other strategy is based on the regular grid but uses some operators to sum. Both methods can provide different forms of Hamiltonian structure and conservative semi discretization of discrete energy. However, using the existing SBP formula, the format obtained by the second strategy can directly obtain high-order accuracy, while the format based on the element center grid is not easy to improve the accuracy. Further, combining the implicit midpoint method and Scalar Auxiliary Variable (SAV) method, we construct the symplectic integrator and linear implicit energy preserving scheme for two-dimensional sine Gordon equation, respectively. A large number of numerical experiments show their effectiveness under homogeneous Neumann boundary conditions.

Example 5 (Fractional Equation)

 The last two decades have witnessed a great progress in fractional calculus and fractional-order dynamical systems. It has been found that fractional calculus is a mathematical tool that works adequately for “anomalous” social and physical systems with non-local, frequency- and history-dependent properties and for intermediate states such as “soft” materials, which are neither idea solid nor idea fluid. Differential equations with fractional-order derivatives/integrals are called fractional differential equations and they have found many successful applications in viscoelasticity, heat conduction, electromagnetic wave, diffusion wave, control theory and so on. Except for a few special cases, it is impossible to find a closed-form solution for a fractional differential equation. Therefore, conditions that govern the existence of some kind of solutions are very important in understanding real systems described by fractional differential equations [13-16]. 

Based on the newly developed block average vector field method, a structure preserving numerical scheme for two-dimensional fractional Klein-Gordon-Schrödinger equation is constructed in Fu et al. [17]. Firstly, the authors derive an equivalent equation and transform the equation into an infinite dimensional regular Hamiltonian system by using the fractional Laplace variational derivative of the functional. Then, the Fourier pseudo spectral method is used to discretize the equation in the spatial direction and a semi discrete conservative system is obtained, which can be transformed into a finite dimensional regular Hamiltonian system. The block average vector field method is further applied to semi discrete systems and a class of fully discrete schemes which can accurately maintain mass and energy are obtained. Finally, the theoretical analysis results are verified by numerical examples.

In Xu et al. [18], taking the spatial fractional nonlinear Schrodinger equation as an example, the convergence error analysis of the conservative Fourier pseudospectral method is established. The authors introduces a new fractional Sobolev norm to construct discrete fractional Sobolev spaces and proves some important lemmas of the new fractional Sobolev norm. Based on these lemmas and energy methods, a priori error estimation of the method can be established. Then, it is proved that the Fourier pseudospectral method is unconditionally convergent. Numerical examples verify the correctness of the theoretical analysis.

In Fu et al. [19], a new method for constructing explicit structure preserving schemes for Nonlinear Fractional wave equations is proposed. Firstly, the equation is reformulated as a regular Hamiltonian system [20]. Then, the equation is discretized in space by using the fourth-order fractional central difference formula and a conservative semi discrete system is obtained. Then, an explicit fully discrete energy preserving scheme for the equation is developed by using the proposed method. Finally, a numerical example is given to verify the theoretical analysis.

This paper reviews the research status of structure preserving algorithm, summarizes some existing methods and conclusions and many meaningful results have not been summarized in this paper. Although many fruitful research results have been achieved in the structure preserving algorithm, there are still some research difficulties and bottlenecks to be solved. It is still necessary to further excavate and develop the structure preserving algorithm suitable for practical scientific and engineering problems.

Acknowledgment

Project Supported by the Major Projects of Natural Sciences of University in Jiangsu Province of China (No. 18KJA110003) and the Qianfan Plan of Jiangsu Maritime Institute (No. 014060).

1. Yang, Tao. Structure-Preserving Algorithms for Beam Vibration Equations. Master’s thesis, Nanjing Normal University, 2015.

2. Peng, Yingying. Local Structure-Preserving Algorithms for the Damped KdV Equation. Master’s thesis, Nanjing Normal University, 2018.

3. Xu, Qian. Structure-Preserving Algorithms for Some Partial Differential Equations. Doctoral dissertation, National University of Defense Science and Technology, 2014.

4. Sun, Lang et al. "Classification and Development of Symplectic Algorithm." Progress in Astronomy, vol. 39, 2021, pp. 211-233.

5. Luo, Ming-Qiu et al. "Hamiltonian Description and Symplectic Method of Seismic Wave Propagation." Chinese Journal of Geophysics, vol. 44, 2001, pp. 120-128.

6. Luo, Ming-Qiu et al. "Seismic Wave Modeling with Implicit Symplectic Method Based on Spectral Factorization on Helix." Chinese Journal of Geophysics, vol. 44, 2001, pp. 379-388.

7. Wang, W.S. et al. "Structure-Preserving Modeling for Seismic Wavefields Based upon a Multisymplectic Spectral Element Method." Chinese Journal of Geophysics, vol. 55, 2012, pp. 3427-3439.

8. Du, Lin and Pinglan Hou. "Symplectic Precise Integration Method for Duffing Equation." Journal of Dynamics and Control, vol. 15, 2017, pp. 1-5.

9. Gang, Tie-Qiang et al. "Structure-Preserving Algorithms for the Duffing Equation." Chinese Physics B, vol. 17, 2008, pp. 3623-3628.

10. Hernández, Quercus et al. "Structure-Preserving Neural Networks." Journal of Computational Physics, vol. 426, 2021, 109950.

11. Jiang, Chaolong et al. "A Linearly Implicit and Local Energy-Preserving Scheme for the Sine-Gordon Equation Based on the Invariant Energy Quadratization Approach." Journal of Scientific Computing, vol. 80, 2019, pp. 1629-1655.

12. Cai, Wenjun et al. "Structure-Preserving Algorithms for the Two-Dimensional Sine-Gordon Equation with Neumann Boundary Conditions." Journal of Computational Physics, vol. 395, 2019, pp. 166-185.

13. Podlubny, I. Fractional Differential Equations. Vol. 198, Mathematics in Science and Engineering, Academic Press, 1999.

14. Hilfer, R. Applications of Fractional Calculus in Physics. World Scientific, 2000.

15. Xu, M.Y. and W. Tan. "Intermediate Processes and Critical Phenomena: Theory, Method and Progress of Fractional Operators and Their Applications to Modern Mechanics." Science in China (Series G), vol. 49, no. 3, 2006, pp. 257-272.

16. Kilbas, A.A. et al. Theory and Applications of Fractional Differential Equations. Vol. 204, North-Holland Mathematics Studies, Elsevier, 2006.

17. Fu, Yayun et al. "Structure-Preserving Algorithms for the Two-Dimensional Fractional Klein-Gordon-Schrödinger Equation." Applied Numerical Mathematics, vol. 156, 2020, pp. 77-93.

18. Xu, Zhuangzhi et al. "On the Convergence of a Conservative Fourier Pseudo-Spectral Method for the Space Fractional Nonlinear Schrödinger Equation." Numerical Methods for Partial Differential Equations, 2020, pp. 1-21.

19. Fu, Yayun et al. "An Explicit Structure-Preserving Algorithm for the Nonlinear Fractional Hamiltonian Wave Equation." Applied Mathematics Letters, vol. 102, 2020, 106123.

20. Fu, Yayun et al. "A Linearly Implicit Structure-Preserving Scheme for the Fractional Sine-Gordon Equation Based on the IEQ Approach." Applied Numerical Mathematics, vol. 160, 2021, pp. 368-385.