# Mathematical Control of Space-Based Kinetic Energy Weapons Based on Partial Differential Equations and Evaluation of Their Destructive Effects

Authors: Zhixin Lian; Tiehua Ma
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

This paper presents an in-depth study and analysis of the mathematical control of space-based kinetic energy weapons and the evaluation of the damaging effect by partial differential equations. The spectral element discrete format of the optimal control problem is constructed, the a priori error estimate of the control problem solution is proved theoretically, the posteriori error estimator is constructed, and the adaptive solution algorithm is designed. The posteriori error estimator is used as the encryption criterion to guide the local encryption of the grid so that the distribution of the dissection nodes is denser where the function regularity is poor. In the case of a compartment subjected to a shell attack, the effect of different factors on the structural damage of a single compartment under two explosions is investigated by varying the explosive mass, the location of the blast point, and the interval between the two explosions. In this paper, the problem of implosion in the cabin is studied, and the main factors affecting the response of the implosion structure are analyzed using dimensionality. The data of various simulation conditions are counted, and the dimensionless damage number is fitted with the deformation results of the bulkhead under each damage mode. It is still difficult to obtain the acceleration information of the target. Experimental studies were conducted on the composite honeycomb sandwich structure penetrated by fragmentation at different velocities and angles, and the accuracy of the theoretical model and the fragmentation residual formula was verified based on the test results. The analysis found that the antidamage performance of the composite honeycomb sandwich structure is better than the existing honeycomb sandwich body structure, and the energy absorption per unit volume of the structure is 27.2%–84.2% higher than that of the existing structure in the range of different penetration velocities. The average error between the theoretical calculation results and experimental results of the remaining velocity of the broken piece is within 6%, the error with the numerical simulation is about 8%, and the composite honeycomb sandwich structure with the best energy absorption characteristics and structural parameters is obtained.


---

## Body

## 1. Introduction

After having a high-power laser and a stable and reliable follow-up system, the high-energy laser reaches the target through atmospheric transmission, and its physical and chemical properties are changed by irradiating it for a certain period to finally complete the damage to the target structure or function. In this process, the influence of the atmospheric environment on the high-energy laser transmission cannot be ignored, and the atmospheric transmission characteristics are an important factor influencing the application effect of high-energy laser weapon systems and an important basis for their system testing and evaluation. In addition to the attenuation effect and turbulence effect of the atmosphere on the laser, the atmospheric transmission characteristics must also consider the effect of the thermal corona effect for high-energy lasers [1]. The attenuation effect makes the laser energy decay after transmission, the turbulence effect makes the laser transmission time intensity or phase undulation, beam drift or expansion, the thermal corona effect makes the beam further distorted and quality degradation, and the combined effect of turbulence and thermal corona will have further impact. The fly-net bomb ultra-close range active protection interception system is a new concept of ultra-close range active protection system, which detects the incoming munitions into a certain range in real time, solves the ballistic trajectory of the target munitions, and launches the fly-net bomb with the action. The fly-net projectile captures and induces detonation of the target munition outside a certain range from the protection target so that the protection fortification is changed from possible penetration damage to fragmentation and shock wave damage, thus increasing its survival probability [2]. However, the maneuvering acceleration of the target plays a key role in the terminal interception of the interceptor, so this chapter also designs a disturbance observer based on high-order sliding mode for the terminal guidance stage. Underwater explosion load is very powerful, whether shock wave load or bubble pulsation and bubble jet can cause devastating damage to the ship, especially the ship’s structure in the whole process of near-field explosion in different stages of the combined effect of different forms of load, and it is very easy to lead to the complete loss of life, so in the field of ship impact resistance and underwater attack weapons, research and development are required to carry out in-depth research on this. US Navy combat ships are required to undergo underwater blast impact tests before they are commissioned.The study of inverse problems and mathematical control of stochastic partial differential equations is still in a very preliminary stage. The main difficulties in the theoretical analysis of stochastic equations are that the solutions of differential equations are affected by random noise, the classical derivatives concerning timet no longer exist, some tight embedding theorems for the deterministic case no longer hold in the stochastic case, and the immovable point method for nonlinear partial differential equations is not always feasible for stochastic nonlinear partial differential equations [3]. In general, it is difficult to find the exact solution to these control problems; therefore, it is important to study the efficient numerical solution algorithm to check the correctness of the mathematical model and the successful application of the optimal control problem. At present, the most important numerical method for solving optimal control problems with restricted partial differential equations is the finite element method [4]. In some extreme conditions, the air defense interceptor position may be destroyed by the attacking bomb before it is fully operational. Second, the higher flight speed possessed by hypersonic weapons as the attacker poses a great threat to the interception capability of interceptor bombs. The model based on the target motion equation is used to filter the system state with the augmentation of position parameters, which is equivalent to the identification of the parameters of the motion model. Traditional interceptor systems in the interception of cruise missiles, aircraft and other typical targets, the use of their speed advantage, etc., can more easily achieve interception. However, in the face of hypersonic targets, the speed advantage no longer exists, the inescapable area for interception of hypersonic weapons by interceptors is significantly compressed, and the interceptor’s operational effectiveness is severely challenged.The study of explosions in confined spaces began with the analysis of the factors causing structural damage, and as the research progressed, the study of the destruction of cabin structures by shock wave overpressure, quasi-static overpressure, fragmentation loads, and coupling effects between them attracted extensive attention from scholars. In this paper, we focus on the effects when a hull hits a ship and explodes in the compartment and the change law of compartment wall response when multiple weapons explode in the compartment. This has some application value for the design of ship protection structures and various damage assessment techniques that are compatible with intelligent combat. The main forms of implosion loads are high-speed fragmentation, shock wave overpressure, and quasi-static overpressure. High-speed fragmentation duration is the shortest, mainly resulting in localized damage, and cabin equipment and other obstacles can have a dampening effect on the high-speed fragmentation load. Shock wave overpressure peak is a large, wide range of action, reflected and superimposed at the bulkhead many times, and the load form is complex and related to the specific structure of the cabin. Quasi-static overpressure peak is smaller, the longest duration, the scope of action in the chamber where the explosion occurred, and the adjacent chamber that caused the breach, and the cabin equipment will enhance the destruction of quasi-static overpressure loads. The shock wave overpressure and quasi-static overpressure can be collectively referred to as the “pressure load.”


## 2. Related Works

Goldberg pilots did not have any feedback on the damage of the target, and decision-makers needed to know the damaging effect of the target being struck to reasonably allocate combat resources and develop the next battle plan, so the Air Force pilots started to take pictures of the target during the strike and write postbattle reports [5]. Iohara et al. elaborated on the main differences between system targets and typical targets, gave the basic principles and general steps for the assessment of the damaging effect of system targets according to their characteristics, and discussed the effect decay model and the selection of indicators for the damaging effect accordingly [6]. Boukera et al. proposed the module composition, overall structure, and evaluation process of the general target destruction effect evaluation simulation system, and wrote the program to give the relevant calculation examples [7]. The research on the inverse problem of the stochastic partial differential equation and the problem of mathematical control is still in a very preliminary stage [8]. The main difficulties in the theoretical analysis of stochastic equations are that the solutions of differential equations are affected by random noise, the classical derivatives concerning time t no longer exist, some tight embedding theorems of the deterministic case no longer hold in the stochastic case, and the immovable point method for nonlinear partial differential equations is not always feasible for stochastic nonlinear partial differential equations [9]. These difficulties lead to the fact that many method dealings with the controllability and inverse problems of deterministic equations cannot be applied directly to the stochastic case.Harfash considers a system of stochastic second-order and fourth-order parabolic equations [10]. The authors establish Carleman estimates for the stochastic heat equation and the fourth-order operator, respectively, and sum up to obtain the Carleman estimates for the coupled system of equations. Based on this Carleman estimate, the zero controllability of the coupled system is obtained [11]. Gupta proves the zero controllability of the coupled inverse stochastic heat equation system [12]. Because of the transient occurrence, this process can usually be regarded as a linear problem, so the shock wave load can be accurately calculated by empirical formula, then the energy carried by the shock wave can be obtained by integrating the relationship between the wave energy and the shock wave pressure, and then the energy carried by the shock wave can be obtained. There is no study related to the Carleman estimation of the system of stochastic degenerate equations and its application. The structure was tested for intrusion using a falling hammer experimental machine, the intrusion energy was relatively low compared to Wang’s study, and the main energy dissipation of the structure was still dominated by the fiber tensile fracture of the composite skin [13]. Strength gradient and precompression experiments were conducted, respectively, and it was found that the tandem structure exhibited better mechanical properties, while no peak force appeared in the subsequent deformation phase of the structure under precompression, the transverse resistance had a significant effect on the deformation mode, and filling the honeycomb tube with tandem honeycomb could significantly improve the mechanical properties of such structures [14]. For honeycomb structures, the three commonly used analysis methods are static compression, ballistic damage, and blast impact. Under static compression, the compressive strength is much smaller than that under dynamic loading, and the plastic deformation region of the honeycomb structure is near the region where the sandwich structure is in contact with the intruder, while its damage is mainly concentrated inside the specimen, and the structure is more easily damaged during the intrusion process, which leads to fiber fracture, debonding delamination, or even perforation.However, the propagation of pressure load inside the compartment and the effect law of multi-shot explosives on the destructive effect of the compartment wall when multi-shot explosives explode inside the compartment have not been fully studied. In the actual battle, the ship is likely to be attacked by multiple weapons, so this paper studies the damaging effect of multi-shot explosives on the cabin structure. The main purpose of this paper is to provide a quick prediction of the damage range of the compartment when it is subjected to multiple explosive attacks and to provide a reference for the damage assessment of multiple weapons. In this paper, we use the method of dimensional analysis to clarify the main factors affecting the structural response of the internal explosion and derive the dimensionless damage number of the internal explosion when the explosive charge and the blast distance change, the data of the simulation conditions are counted, and the dimensionless damage number is fitted with the deformation results of the bulkhead in each damage mode.


## 3. Analysis of Mathematical Control Partial Differential Equations for Space-Based Kinetic Energy Weapons

The optimality condition is given, the discrete format of spectral element approximation is constructed, the a priori error estimate of the optimal control problem is proved, the posteriori error estimator is constructed, and the adaptive solution algorithm is designed. Finally, the discrete system is solved using the gradient projection algorithm, and the efficiency of the adaptive solution algorithm is verified by extensive numerical experiments on the indication of the error by the posteriori error estimator. In this section, we discuss the integral control-constrained optimal control problem and consider the following model problem. Define the local velocities of the body as specified, including coupled translational and rotational motions.(1)maxu∈FadJu=12∫Ωy−yc2−12∫Ωu.The basic idea of the Eulerian finite element method is to separate the transport term from the Eulerian equation, thus splitting it into two equations, solving them separately, and then combining the results of the two parts to obtain an approximate solution of the original equation. The method has a clear physical meaning: the equation without the transport term can be considered as fluid motion from the Lagrangian point of view; i.e., the grid follows the fluid mass under the action of the pressure gradient [15]. Since there is no transport term in this phase, it can be solved directly by the explicit finite element method without the need to construct the windward format; the second equation, which contains only the transport term, indicates that at the end of the time step, the fluid mass is kept still and the deformed mesh is moved back to its original position.In this process, there is a relative motion between the cell and the fluid, and the transport of matter between the cells needs to be considered, a phase called the Eulerian phase. Through the combination of the two phases, the fluid matter keeps moving forward with time, while the mesh is fixed to maintain the original mesh mass. The method uses the proven techniques of the finite element method and the finite volume method, which makes the method highly robust and stable, especially for violent bubble motions with relatively high density at the interface. This kind of mesh control is almost exclusively used to solve problems involving rigid body motion, and it is necessary to define the mesh deformation range, the physical parameters of the rigid body, and the force and motion constraints on the rigid body in the dynamic mesh file Dynamometric.(2)aiL=1V−ΔVaiV+K−KiΔV2.In the Lagrangian step, the mesh and the fluid physical quantities are moved forward first, but in large deformation problems, the movement of the mesh can produce distortions and other factors that are detrimental to computational accuracy and stability. Therefore, it is also necessary to move the mesh back to its original position by a Eulerian step after each Lagrangian step. In the Eulerian step, the fluid mass is fixed, the mesh is moved backward to its initial position to eliminate the cell deformation according to the velocity calculated in the Lagrangian step, and the relative motion between the fluid and the mesh occurs.(3)G=Ac∑j=138Rjaj2.With the rapid development of software and hardware technologies, numerical computation methods have become a common means for people. Therefore, an effective mathematical model can simplify the complex problems caused by numerical simulation and other methods. Compared with analytical methods, numerical computation methods lack certain rigor, but they are extremely adaptable and have outstanding advantages when dealing with practical problems. The core of numerical computational methods is discretization, and in general, the discretization in time and space is significantly different, while the discretization in space is relatively complex and has a great impact on the computational results [16]. Another focus of numerical computation method is the discrete format, which is the embodiment of its mathematical ideas and according to which different methods are divided accordingly, and the most common ones are finite difference method, finite volume method, and finite element method.(4)CT˙−KT˙=P,∂T∂tI+Δt=1∇tTt+Δt+Tt−O∇t.The motion law of the vehicle is determined by the external load acting on the vehicle, and the trajectory and motion state of the object can be determined by solving the six-degree-of-freedom motion control equations. The independent part is updated by numerical integration, the nonindependent part is updated by iteratively solving a nonlinear system of equations, and after one step of integration, error checking and integration step size selection are performed. If the tolerance error is satisfied, the state variables are updated and the simulation results are output to the simulation result file at the specified time interval before the next integration step.(5)Pw=2γpSij+wSnmδnmSij.Although it is possible to simulate most of the penetration process and the breakage of the target plate with the help of finite element simulation, the results are always poor compared to experiments. Therefore, an effective mathematical model can simplify the complexity of the problem brought about by numerical simulations, etc. By simplifying the problem, the elements of the material can follow simple stress-strain laws, thus making the main features of the intrusion process obvious.The elastic energy absorption of the underlying structure fills the problem of the lowering of the energy absorption curve caused by the stress drop of the upper core layer during the softening stage, resulting in the phenomenon that the structural energy absorption does not change significantly when the softening stage transitions to the crushing stage. Predictive corrected guidance is a common algorithm for ballistic planning of re-entry vehicles, but the prediction process of predictive corrected guidance is complex and is generally obtained by integrating differential equations, which is computationally intensive and has poor real-time performance. As the target model in this paper, the accuracy requirement for ballistic planning is not the primary consideration, so an energy-based analytical predictive guidance algorithm is used to perform ballistic simulation for the boosted glide hypersonic vehicle, as shown in Figure1.Figure 1
Mathematical control of partial differential equation steps.When making predictions, future ballistic predictions can be made based on the lift-to-resistance ratio at the current moment [17]. Between 3.41 ms and 6.21 ms, it can be observed that the deflection value of measuring point 6 continues to increase, and the increasing trend of deflection value of measuring point 7 is obviously slowed down, which indicates that the deformation of the bulkhead at this stage is mainly concentrated in the area near the center of the bulkhead; after 6.52 ms, the deflection value of each place no longer increases substantially and begins to converge to a certain fixed value. Although the lift-to-resistance ratio changes throughout the ballistic process, this error is acceptable when the distance is far, the range of error in the lift-to-resistance ratio becomes smaller and smaller during the approach to the target point, and the prediction accuracy becomes higher and higher. Although this analytical method has limited accuracy compared with the traditional prediction method, the computing efficiency is greatly improved, which is very suitable for the research object of this paper.In addition, in the end-guidance phase, the general infrared detectors can only obtain the relative angular information of the target, and it is still difficult to obtain the acceleration information of the target even when supplemented with the distance and velocity information from external detection systems. However, the maneuvering acceleration of the target plays a critical role in the implementation of terminal interception by the interceptor. Therefore, in this chapter, a perturbation observer based on a high-order sliding mode is also designed for the terminal guidance phase, and the unknown acceleration information of the target is observed and estimated as a random perturbation. Regarding the convergence of the sliding mode observer, the detailed proof is given in this paper to ensure that it can achieve accurate estimation in finite time while satisfying the assumptions. From the overload characteristics, the relative change in the target’s overload has a large magnitude and a long change duration, and although the absolute value of the overload change is not large in a short period, more accurate modeling of the target’s motion must also be performed to perform a 50–100 s ballistic forecast, as shown in Figure2. With the continuous action of the pressure load in the cabin, the units in the direction of the final explosive connection line are also all separated from the bulkhead.Figure 2
Variation of center-of-mass position.However, in this paper, in addition to the target tracking, it is also necessary to predict the future motion of the target for 50–100 seconds to give the predicted collision point [18]. To solve this problem, a model based on the equation of motion of the target, supplemented with position parameter augmentation, is used to filter the system state, which is equivalent to the identification of the parameters of the motion model. However, when using this type of motion model, many unknown parameters cannot be obtained directly, they need to be assumed using a stochastic process and augmented to the state of the original equation of the target, and the specific process of this method is described below.When performing such state augmentation and filtering, we must admit that the final obtained parameter identification results do not necessarily achieve high accuracy due to the noise of the observations, the degree of approximation of the augmented state model, and other factors, but this part of the study is mainly to make a general forecast of the future state of the target vehicle and does not require a very accurate prediction. Moreover, as the prediction time is shortened, the prediction error will be reduced continuously.


## 4. Design for Mathematically Controlled Destruction Effects of Space-Based Kinetic Energy Weapons

As mentioned earlier, the underwater explosion process includes shock wave and bubble stages, both of which account for an important proportion of the total energy. In the shock wave phase, the blasted product generates a strong intermittent shock wave by compressing the surrounding aqueous medium, which carries the shock wave energy consisting of kinetic and internal energy to propagate outward and eventually dissipate in the distant flow field. Because of the transient occurrence, this process can usually be considered a linear problem, so the shock wave load can be calculated accurately by empirical formulas, and then the energy carried by the shock wave is obtained by integrating the relationship between the wave energy and the shock wave pressure, and then obtains an empirical formula that can meet the engineering needs [19]. Therefore, although the underwater explosion shock wave and bubble two-stage full process numerical simulation method established in this paper can simulate the shock wave stage well, this section does not focus on the energy dissipation in the shock wave stage as the research focus.In the bubble phase, when the underwater explosion bubble moves in the free field, the internal gaseous explosion products and the external flow field constitute a dynamical system, and under the balance of pressure and inertial forces, the bubble generates multiple pulsations and radiates the bubble pulsation load outward. Experimental studies have found that the maximum radius of multiple bubble pulsations decreases one by one, which indicates that the energy of the whole system is continuously dissipated. This energy dissipation is essential for the accurate prediction of bubble pulsation loads. However, bubble motion is a typical nonlinear large deformation process, which is usually accompanied by complex processes such as nonspherical motion, jet development, and tearing and fusion, while the complete compressibility of the flow field needs to be accounted for, which poses great difficulties for accurate calculation of energy dissipation in this phase, as shown in Figure3. The attenuation effect makes the energy attenuate after laser transmission, the turbulence effect makes the light intensity or phase fluctuate, the beam drifts or expands, and the thermal halo effect further distorts the beam and reduces its quality. The combined effect of turbulence and thermal halo will have further influence.Figure 3
Mathematical control model for space-based kinetic weapons.The angle selected in Figure3 can observe the propagation pattern of the shock wave in the center of the chamber and the area near the bulkhead. The phase represented by 0.7∼0.8 ms indicates that the shock wave propagates along the bulkhead boundary to the corner of the chamber and converges at the corner; the phase represented by 1.1∼1.6 ms indicates that the shock wave converges at the corner and then reflects along the diagonal of the chamber and bulkhead toward the center of the chamber and bulkhead, and the high-pressure area is formed again at the center of the chamber with the highest pressure inside the chamber [20]. From the center of the chamber, the shock wave starts to propagate to the bulkhead once again. The difference with the first incident phase is that before the subsequent shock wave reaches the other areas of the chamber, there is already a certain value of overpressure, and the shock wave propagation path is more “dispersed” and no longer concentrated at the boundary of the bulkhead. By simplifying the problem, the elements of the material can follow a simple stress-strain law, so that the main characteristics of the penetration process become obvious.The mesh motion is calculated based on the pressure on these boundaries, and in turn, the solver provides feedback to the fluid simulation. It changes the velocity boundary conditions on the included boundaries to specify the local velocity of the defined body, including coupled translational and rotational motions. This mesh control is used almost exclusively to solve problems involving rigid body motion and requires information on the extent of mesh deformation, physical parameters of the rigid body, and force and motion constraints on the rigid body to be defined in the dynamic mesh file DynamicMeshDict. In particular, under the combined action of different forms of loads at different stages in the entire near-field explosion process, the ship structure can easily lead to the complete loss of vitality.The density, pressure, and velocity values of the gas change drastically under the interference of the excitation waves, and the pressure at the excitation wave action increases significantly [21]. Then, during 0.0002 s∼0.0015 s, due to the blockage of the sub-bullet body, the excitation waves are repeatedly superimposed, causing multiple reflections of excitation waves and mutual intervention, resulting in the gas pressure in this region which has been increasing; due to the blockage of the support body surface, a stronger beam of excitation waves is reflected toward the bullet, and the range of the low-pressure region formed in the lower side of the bullet increases significantly; the tail of the bullet is complicated by the excitation wave disturbance. The complex structure of the flow field at the tail of the projectile due to the complex surge disturbance and the obstruction at the bottom of the support body is shown in Figure 4.Figure 4
Dimensionless parameters at different sound velocity working conditions.The two sets of numerical simulations start to have errors with the experimental data at a later stage, and this error accumulates with time. The distance between the same node of the experiment and the simulation is getting larger at the same moment. Therefore, it is necessary to carry out in-depth research in the field of ship shock resistance and the research and development of underwater attack weapons. The above set of experiments and two sets of simulations are briefly analyzed and evaluated, and the following conclusions are drawn: both sets of simulation models can simulate the motion attitude of the bead chain in the experiment to some extent, the simulation effect is good in the early stage, and the error gradually increases in the later stage.After this, the upper layer of the composite honeycomb sandwich structure is in the softening stage, while the deformation of the middle skin is gradually increasing the load applied to the lower honeycomb core layer. Due to the inertia effect, the deformation trend of the upper layer structure is larger than that of the lower layer structure, the upper layer structure is the first to enter the compacting stage, and the energy absorption curve of the structure can be seen to increase significantly. After 150μs, the bottom skin is penetrated completely and the energy absorption curve reaches the limit value.


## 5. Analysis of Results

### 5.1. Analysis of Algorithm Performance Results

The consistent and adaptive grids are used to compute Example 2, and the control, state, and conjugate state variables are approximated by the piecewise higher-order elements. In Table1, two different sets of results are presented, where M in the consistent and adaptive grids denotes the overall number of cells dissected, except that the consistent grid is uniformly dissected each time, while the adaptive grid is more densely dissected at each cell to be encrypted and the rest of the cells are sparsely dissected. CS denotes the number of encryptions of the cells to be encrypted in the adaptive grid, and N denotes the polynomial number.Table 1
Numerical calculation results forN  = 7.
TypeCSN149.323.577.271.279215.124.556.934.622.9357.880.713.774.234.3480.614.950.552.262.8554.41487.455.252.4Numerical results are given for each cell with seven and eight polynomial approximations. The results show that for the same number of cells, the adaptive grid converges significantly faster than the consistent grid, greatly reduces the number of heavy operations, the errors are within the desired range of controllability, and the computational efficiency is significantly improved.Theoretical proofs of the a priori and a posteriori error estimates of the numerical solutions are presented. The hN adaptive spectral element algorithm is designed, and two numerical examples are used to verify the correctness of the theoretical results. In the first example, the solution regularity of the elliptic equation is good, and a consistent grid is used for the calculation. The numerical results show that the adaptive grid converges faster than the consistent grid, which greatly reduces the computational effort and improves the computational efficiency. As an attacker, the higher flight speed of hypersonic weapons poses a huge threat to the interception capability of the interceptor.The adaptive spectral element method for the control-constrained elliptic optimal control problem is studied. The spectral element discrete format of the optimal control problem is constructed, the a priori error estimate of the solution of the control problem is proved theoretically, the a posteriori error estimator is constructed, the adaptive solution algorithm is designed, and the discrete system is solved by combining the gradient projection algorithm and finally verified using numerical experiments. The solution regularity of the first example is good, the computation is performed using a consistent grid, and the numerical results show that the indicator is efficient; the solution regularity of the second example is not good, the computation is performed using a consistent grid and an adaptive grid, the numerical results again show that the adaptive grid converges faster than the consistent grid, and the singularity of the solution can be captured quickly, which greatly reduces the number of heavy operations and the error. The error is also within the ideal range of controllability, as shown in Figure5. Results are always far from satisfactory when compared to experiments.Figure 5
Location of measurement points and the corresponding deflection time course curves.It can be seen from Figure5 that the deformation trend of the three measurement points at the boundary is the same before 1.35 ms, and the deflection value increases simultaneously in the region inside the bulkhead boundary and at the boundary when the bulkhead boundary starts to deform and the middle region of the bulkhead deforms together with the boundary. After 1.35 ms to 3.41 ms, the slope of measurement points 7 and 8 changed less, and the slope of measurement point 8 decreased significantly, indicating that the deflection value at the boundary increased at a significantly slower rate, and the deflection value in the middle region of the bulkhead continued to increase at a higher rate synchronously with the center of the bulkhead. The deformation in the bulkhead boundary region is significantly reduced because of the applied solid support boundary condition. Between 3.41 ms and 6.21 ms, it can be observed that the deflection value of measurement point 6 continues to increase, while the deflection value of measurement point 7 tends to slow down significantly, which indicates that the deformation of the bulkhead is mainly concentrated in the area near the bulkhead center during this period, and the deflection value of each area stops increasing significantly after 6.52 ms and starts to converge to a fixed value.All other boundaries are assumed to be far enough away from the bubble except for the horizontal wall located near the bubble. Due to the pressure difference between the bubble interfaces, an initial Riemann problem is formed. Then, the excitation waves propagate outward in the surrounding water, while the rarefaction waves propagate inward in the inner gas. The main forms of implosion loads are high-speed fragmentation, shock wave overpressure, and quasi-static overpressure. At the same time, the bubble interface expands with the intermittent motion of the contact. If the effects of rigid walls and buoyancy are not considered, the bubble will oscillate in a spherical shape until all the energy is dissipated by the emitted pressure waves. If the direction normal to the wall is parallel to the direction of gravity, the problem can be considered axisymmetric.However, one significant drawback of strategy iteration is that the strategy must be evaluated in each iteration, and strategy evaluation requires traversing each state, which is itself a computationally intensive step. To reduce the computational time consumption of strategy iteration, it is not necessary to use a fully convergent value function for strategy evaluation. In many cases, too many iterations of the value function do not affect the strategy update, so a simplification of the strategy evaluation can be considered, and a more effective simplification method is to perform only one round of iterations of the strategy evaluation and then go directly to the strategy update.


### 5.2. Space-Based Kinetic Weapon Mathematical Control Destruction Effect Assessment Design

In the case of interception using the extended proportional guidance law, the initial overload of the interceptor is higher than that of the proportional guidance law; i.e., the target is intercepted with a more advanced posture, so the overload command of the interceptor does not increase rapidly and saturates quickly compared with that of the proportional guidance law. However, as the target with differential countermeasure guidance law performed a reverse maneuver at 3.5 seconds, the interceptor overload demand changed rapidly. During the final seconds of the interceptor maneuver, the interceptor continues to approach the target with the maximum saturation maneuver, but it is still not enough to achieve the interception of the target.Target vulnerability refers to the ease with which a target can be damaged by enemy fire. In an operational scenario, a weapon system killing a target is a complex random event that can be viewed as two mutually independent events in sequence. The equation without the transport term can be regarded as the fluid motion in the Lagrangian view; that is, the grid moves along with the fluid particles under the action of the pressure gradient. The first random event is the probability of the combat section hitting the target, a process that affects many factors and is usually described in terms of tactical vulnerability; the second random event is the probability of causing corresponding damage to the target after the combat section successfully hits the target, which is a conditional probability and is usually described in terms of structural vulnerability. Therefore, to cause a certain degree of damage to the target, two independent events must be guaranteed to occur simultaneously; i.e., the final kill probability is the product of the probability of hitting the target and the probability of damage to the target after the hit. In this study, considering the specific research focus and concern, the analysis of the structural vulnerability, i.e., the probability of the target being damaged under the hit condition, is mainly carried out for the convenience of the study.The above situation arises because of the convergence of the shock wave loads reflected in the cabin several times, and the form of the shock wave loads in the enclosed space is very complex and difficult to derive directly from the theory. Next, the propagation of shock waves will be illustrated by analyzing the pressure load distribution in some air domains inside the chamber during offset detonation to further explain the above situation where the highest value of deflection of the surrounding bulkhead appears to shift to the left, as shown in Figure6.Figure 6
Destruction situation and conditional destruction probability.The evaluation process can be divided into two major parts: the first part is the calculation of the thermal ablation effect of the laser on the target in the atmosphere, whose specific mathematical model and numerical calculation methods have been described more completely in the previous section. Discrete in space is relatively complex and has a great influence on the calculation result. The second part is the conversion of the thermal ablation effect to the degree of the target function and combat capability destruction, which is closely related to the target characteristics. Taking the typical target missile of laser weapon as the object of study, the analysis of its target vulnerability is carried out, including the determination of its structural composition and key components, the selection of damage criteria, and the classification of damage level, taking into account the differences in structure and material of different components, their interaction mechanism with laser and damage effect is also different, and the laser-to-target energy density is no longer used as the judgment criterion. Instead of using the laser-to-target energy density as a criterion, the components were subdivided according to their characteristics, and the component temperature, component strength, and component ablation volume were selected as the indicators of the damaging effect for the missile’s combat wing, engine, and tail fin, respectively, as shown in Figure7.Figure 7
Deformation process of bulkhead.As can be seen from the figure, after the multiple reflection and convergence of the shock wave, the plastic deformation at the center of the bulkhead is the largest, and with the increase of time the center of the bulkhead begins to appear a breach, in which the damage along the direction of the explosive arrangement line is greater than the vertical direction. Another key point of the numerical calculation method is the discrete format. In the breach with the explosives connected to the vertical direction of the unit first all failed to delete, while the unit in the direction of the explosive line did not all fail; with the continued action of the cabin pressure load, the final explosive line of direction of the unit also all separated from the bulkhead. For delayed initiation, because the right side of the explosive is first detonated, the shock wave converges on the left half of the compartment area, resulting in the maximum deflection of the bulkhead to the left and the crack area to the left.


## 5.1. Analysis of Algorithm Performance Results

The consistent and adaptive grids are used to compute Example 2, and the control, state, and conjugate state variables are approximated by the piecewise higher-order elements. In Table1, two different sets of results are presented, where M in the consistent and adaptive grids denotes the overall number of cells dissected, except that the consistent grid is uniformly dissected each time, while the adaptive grid is more densely dissected at each cell to be encrypted and the rest of the cells are sparsely dissected. CS denotes the number of encryptions of the cells to be encrypted in the adaptive grid, and N denotes the polynomial number.Table 1
Numerical calculation results forN  = 7.
TypeCSN149.323.577.271.279215.124.556.934.622.9357.880.713.774.234.3480.614.950.552.262.8554.41487.455.252.4Numerical results are given for each cell with seven and eight polynomial approximations. The results show that for the same number of cells, the adaptive grid converges significantly faster than the consistent grid, greatly reduces the number of heavy operations, the errors are within the desired range of controllability, and the computational efficiency is significantly improved.Theoretical proofs of the a priori and a posteriori error estimates of the numerical solutions are presented. The hN adaptive spectral element algorithm is designed, and two numerical examples are used to verify the correctness of the theoretical results. In the first example, the solution regularity of the elliptic equation is good, and a consistent grid is used for the calculation. The numerical results show that the adaptive grid converges faster than the consistent grid, which greatly reduces the computational effort and improves the computational efficiency. As an attacker, the higher flight speed of hypersonic weapons poses a huge threat to the interception capability of the interceptor.The adaptive spectral element method for the control-constrained elliptic optimal control problem is studied. The spectral element discrete format of the optimal control problem is constructed, the a priori error estimate of the solution of the control problem is proved theoretically, the a posteriori error estimator is constructed, the adaptive solution algorithm is designed, and the discrete system is solved by combining the gradient projection algorithm and finally verified using numerical experiments. The solution regularity of the first example is good, the computation is performed using a consistent grid, and the numerical results show that the indicator is efficient; the solution regularity of the second example is not good, the computation is performed using a consistent grid and an adaptive grid, the numerical results again show that the adaptive grid converges faster than the consistent grid, and the singularity of the solution can be captured quickly, which greatly reduces the number of heavy operations and the error. The error is also within the ideal range of controllability, as shown in Figure5. Results are always far from satisfactory when compared to experiments.Figure 5
Location of measurement points and the corresponding deflection time course curves.It can be seen from Figure5 that the deformation trend of the three measurement points at the boundary is the same before 1.35 ms, and the deflection value increases simultaneously in the region inside the bulkhead boundary and at the boundary when the bulkhead boundary starts to deform and the middle region of the bulkhead deforms together with the boundary. After 1.35 ms to 3.41 ms, the slope of measurement points 7 and 8 changed less, and the slope of measurement point 8 decreased significantly, indicating that the deflection value at the boundary increased at a significantly slower rate, and the deflection value in the middle region of the bulkhead continued to increase at a higher rate synchronously with the center of the bulkhead. The deformation in the bulkhead boundary region is significantly reduced because of the applied solid support boundary condition. Between 3.41 ms and 6.21 ms, it can be observed that the deflection value of measurement point 6 continues to increase, while the deflection value of measurement point 7 tends to slow down significantly, which indicates that the deformation of the bulkhead is mainly concentrated in the area near the bulkhead center during this period, and the deflection value of each area stops increasing significantly after 6.52 ms and starts to converge to a fixed value.All other boundaries are assumed to be far enough away from the bubble except for the horizontal wall located near the bubble. Due to the pressure difference between the bubble interfaces, an initial Riemann problem is formed. Then, the excitation waves propagate outward in the surrounding water, while the rarefaction waves propagate inward in the inner gas. The main forms of implosion loads are high-speed fragmentation, shock wave overpressure, and quasi-static overpressure. At the same time, the bubble interface expands with the intermittent motion of the contact. If the effects of rigid walls and buoyancy are not considered, the bubble will oscillate in a spherical shape until all the energy is dissipated by the emitted pressure waves. If the direction normal to the wall is parallel to the direction of gravity, the problem can be considered axisymmetric.However, one significant drawback of strategy iteration is that the strategy must be evaluated in each iteration, and strategy evaluation requires traversing each state, which is itself a computationally intensive step. To reduce the computational time consumption of strategy iteration, it is not necessary to use a fully convergent value function for strategy evaluation. In many cases, too many iterations of the value function do not affect the strategy update, so a simplification of the strategy evaluation can be considered, and a more effective simplification method is to perform only one round of iterations of the strategy evaluation and then go directly to the strategy update.


## 5.2. Space-Based Kinetic Weapon Mathematical Control Destruction Effect Assessment Design

In the case of interception using the extended proportional guidance law, the initial overload of the interceptor is higher than that of the proportional guidance law; i.e., the target is intercepted with a more advanced posture, so the overload command of the interceptor does not increase rapidly and saturates quickly compared with that of the proportional guidance law. However, as the target with differential countermeasure guidance law performed a reverse maneuver at 3.5 seconds, the interceptor overload demand changed rapidly. During the final seconds of the interceptor maneuver, the interceptor continues to approach the target with the maximum saturation maneuver, but it is still not enough to achieve the interception of the target.Target vulnerability refers to the ease with which a target can be damaged by enemy fire. In an operational scenario, a weapon system killing a target is a complex random event that can be viewed as two mutually independent events in sequence. The equation without the transport term can be regarded as the fluid motion in the Lagrangian view; that is, the grid moves along with the fluid particles under the action of the pressure gradient. The first random event is the probability of the combat section hitting the target, a process that affects many factors and is usually described in terms of tactical vulnerability; the second random event is the probability of causing corresponding damage to the target after the combat section successfully hits the target, which is a conditional probability and is usually described in terms of structural vulnerability. Therefore, to cause a certain degree of damage to the target, two independent events must be guaranteed to occur simultaneously; i.e., the final kill probability is the product of the probability of hitting the target and the probability of damage to the target after the hit. In this study, considering the specific research focus and concern, the analysis of the structural vulnerability, i.e., the probability of the target being damaged under the hit condition, is mainly carried out for the convenience of the study.The above situation arises because of the convergence of the shock wave loads reflected in the cabin several times, and the form of the shock wave loads in the enclosed space is very complex and difficult to derive directly from the theory. Next, the propagation of shock waves will be illustrated by analyzing the pressure load distribution in some air domains inside the chamber during offset detonation to further explain the above situation where the highest value of deflection of the surrounding bulkhead appears to shift to the left, as shown in Figure6.Figure 6
Destruction situation and conditional destruction probability.The evaluation process can be divided into two major parts: the first part is the calculation of the thermal ablation effect of the laser on the target in the atmosphere, whose specific mathematical model and numerical calculation methods have been described more completely in the previous section. Discrete in space is relatively complex and has a great influence on the calculation result. The second part is the conversion of the thermal ablation effect to the degree of the target function and combat capability destruction, which is closely related to the target characteristics. Taking the typical target missile of laser weapon as the object of study, the analysis of its target vulnerability is carried out, including the determination of its structural composition and key components, the selection of damage criteria, and the classification of damage level, taking into account the differences in structure and material of different components, their interaction mechanism with laser and damage effect is also different, and the laser-to-target energy density is no longer used as the judgment criterion. Instead of using the laser-to-target energy density as a criterion, the components were subdivided according to their characteristics, and the component temperature, component strength, and component ablation volume were selected as the indicators of the damaging effect for the missile’s combat wing, engine, and tail fin, respectively, as shown in Figure7.Figure 7
Deformation process of bulkhead.As can be seen from the figure, after the multiple reflection and convergence of the shock wave, the plastic deformation at the center of the bulkhead is the largest, and with the increase of time the center of the bulkhead begins to appear a breach, in which the damage along the direction of the explosive arrangement line is greater than the vertical direction. Another key point of the numerical calculation method is the discrete format. In the breach with the explosives connected to the vertical direction of the unit first all failed to delete, while the unit in the direction of the explosive line did not all fail; with the continued action of the cabin pressure load, the final explosive line of direction of the unit also all separated from the bulkhead. For delayed initiation, because the right side of the explosive is first detonated, the shock wave converges on the left half of the compartment area, resulting in the maximum deflection of the bulkhead to the left and the crack area to the left.


## 6. Conclusion

In this paper, the adaptive spectral element calculation of elliptic partial differential equations and their optimal control problems are systematically studied. In addition, to accurately describe the dynamics of the target motion and better evaluate the actual operational effectiveness of the interceptor system, a target motion trajectory planning algorithm based on the prediction correction method is also given, which enables the hypersonic boost-glide vehicle to plan the attack trajectory to the landing area according to its lift-to-drag ratio characteristics. The discrete format is the embodiment of its mathematical ideas, and different methods are divided according to its difference. To analyze the antidamage performance of the composite honeycomb sandwich structure, an experimental study was conducted on the fragmentation of the composite honeycomb sandwich structure at different speeds and angles, and the accuracy of the theoretical model and the fragmentation residual formula was verified. It is found that the accuracy of the residual velocity calculation formula of the fragmentation penetration composite honeycomb sandwich structure is very high, and the theoretical value is always between the numerical simulation and the experimental value with an average error of 6%. The destruction of the circumferential bulkhead was greater than that of both bulkheads under each condition, the destruction of the left bulkhead was greater than that of the right bulkhead during the delayed detonation, and the maximum deflection value of the bulkhead and the area of the breach tended to increase with the increase of the detonation time interval. Comparison of simulation results in the cabin pressure load, velocity curve, the load process, and the classical “three-wave method” is consistent. The difference between different targets also plays a crucial role in the damaging effect. The smaller the target thickness, the lower the thermal conductivity of the material and the lower the curvature of the irradiated surface.



---

*Source: 3420088-2022-09-06.xml*