# Theoretical and Numerical Studies of Dynamic Scaling of  a Six-Degree-of-Freedom Laser Propulsion Vehicle

Authors: Masayuki Takahashi; Naofumi Ohnishi
Journal: International Journal of Aerospace Engineering
 (2015)

---

## Abstract

To estimate the flight reactions of a full-scale vehicle from reduced-scale tests, we constructed a scaling theory for the vehicle size, input energy, moment of inertia, and pulse frequency needed to maintain dynamic equivalence between a laboratory-scale and full-scale launch of a laser propulsion vehicle. The dynamic scaling law for a single pulse was constructed using translational and angular equations of motion. The analytical scaling was confirmed for a single-pulse incident using a fluid-orbit coupling simulator for the interaction between the blast wave and the vehicle. Motion equivalence was maintained for multiple pulses by adjusting the repetition frequency of the pulse incident to correct for the effect of aerodynamic drag during the free flight of the pulse-to-pulse interval. The flight of a full-scale vehicle can be estimated for single- and multiple-pulse operations from the flight data for a small-scale vehicle using the proposed scaling theory, which provides correlations between the characteristics of small-scale and large-scale flight systems. Small-scale tests were shown to be useful in estimating the flight of a full-scale vehicle using dynamic scaling theory.


---

## Body

## 1. Introduction

In recent years, the business of space travel has gradually come into being as a result of the opening of the Space Ship One Program by Scaled Composites and the Space Ship mission organized by Virgin Galactic, which proposes transporting passengers to an altitude of 100 km above the Earth for a fee of twenty million dollars. The Space Exploration Technologies Corporation has proposed a commercial orbit transportation service using the Falcon rocket and the Dragon spacecraft. In addition to commercial plans for space travel, the demand for space utilization for academic research and national security purposes has been growing. The costs of launching payloads should be reduced by using different launching techniques for different types of missions because high launch costs interfere with mission execution. Laser propulsion is a novel launch system involving reception of high-power repetitive beams transmitted from a ground base with little or no fuel loading on a vehicle [1–14]. The vehicle design is simplified because a complex gas-burning engine is not required on the vehicle to provide propulsive force. This simplified design can improve the reusability of the vehicle. It is expensive to construct a huge laser oscillator on Earth, but the total launching cost can be reduced by recouping the initial cost of the laser oscillator using a simply designed vehicle with high reusability and repeatability [15]. In aerodriving mode, the vehicle receives an impulsive thrust through interactions between the vehicle surface and the positive pressure of a strong blast wave, induced by laser irradiation in the atmosphere acting as a driving gas. A high-power repetitive-pulse laser is irradiated to the vehicle from a ground laser oscillator and is focused by a parabolic mirror installed on the vehicle. The intense field of the focused laser beam causes gas breakdown, and a blast wave is generated as propulsive power to the vehicle. However, the generated thrust decreases if the vehicle cannot capture the laser beam because of incident misalignment with the parabolic mirror of the vehicle. If this occurs, the vehicle will lose thrust and flight control as a result of failure to capture the beam, and it will not be possible to maintain the beam-driven flight.Deviation from the laser beam line must be avoided to obtain continuous thrust in a beam-driven flight. The degree of success achieved in capturing the laser beam line, which is referred to as the beam-riding performance, is crucial to maintaining stable flight without deviation from the laser beam [16–24]. The beam-riding performance is totally assessed from the feedback performance with respect to the translational and angular misalignments of the incident laser. When the laser is irradiated with lateral error with respect to the body axis, the vehicle should cancel out the lateral offset to conform the vehicle axis to the laser axis (centering performance). On the other hand, an angular restoring moment must be generated to reduce the angular offset when the laser beam is not irradiated parallel to the vehicle axis (tipping performance). The beam-riding flight becomes unstable by losing posture control due to rapid increases in the lateral and angular offsets if an insufficient restoring force and moment are generated for the lateral and angular misalignments.The vehicle shape must be designed to achieve good beam-riding performance. The “lightcraft” vehicle was proposed by Myrabo [10] as the best beam-riding vehicle among the existing laser propulsion vehicles. The lightcraft is constructed with a front body, a rear body containing a parabolic mirror, and an annular shroud to improve the centering performance. During the flight of the lightcraft, a laser beam is reflected by the rear mirror and is focused inside the shroud, and the interaction between the blast wave and the shroud pushes the vehicle back toward the center of the laser beam if the laser incident is laterally misaligned with respect to the vehicle axis [21]. High-speed axial spin is provided before the launch to obtain gyroscope stabilization for the angular motion, in addition to centering feedback provided by the shock wave interaction with the shroud. A flight test of the lightcraft was conducted in the USA in 2001 using high-power repetitive pulses. A type-200 spinning lightcraft weighing 50 g was launched and achieved an altitude of 71 m [10]. In previous studies of the lightcraft, experiments for the beam-riding physics were conducted to examine the recentering and angular impulses for lateral offset only with a single pulse [18, 19]. In addition, a flight orbit was simulated for multiple pulses using experimental impulse data for lateral offset only [16, 19–21]. However, it is difficult to predict the flight trajectory from the lateral offset alone because the laser incident condition is complicated by the angular offset and the combined lateral-angular offset during the actual flight. State-of-the-art flight simulation is required to achieve more precise prediction of the beaming flight by characterizing the vehicle motion with respect to the angular and combined lateral-angular offsets.In our previous studies, we developed a three-dimensional computational fluid dynamics (CFD) code coupled with a six-degree-of-freedom (6-DOF) orbital calculation code to assess the lateral and angular flight reactions driven by the blast wave interaction [25–30]. A single-pulse reaction was assessed using the code developed, and a vehicle reaction data map was constructed for the lateral, angular, and combined lateral-angular offsets to simulate multiple-pulse flight of the lightcraft. Based on the vehicle reaction data map for a single pulse, the computed flight altitude of the lightcraft was found to agree well with the experimental flight data obtained in a previous flight demonstration for a 32.5 g vehicle using a multiple-pulse laser. The deviation mechanism from the laser beam line was identified from the flight simulation results for the multiple pulses, and an increase in the angular offset was found to cause flight instability. An active laser control concept was proposed based on a genetic algorithm (GA) to suppress any angular offset during the flight. A vehicle weighing 32.5 g can fly at altitudes on the order of kilometers with keeping the beam-riding strategy by introducing a GA-controlled beam incident. However, it is necessary to examine the feasibility of such a flight control technique for a scaled-up vehicle with a high-power laser beam in expectation of the launch of an actual satellite, which would have a larger and heavier body, using a more intense laser beam. A dynamic scaling theory should be constructed to correlate the characteristics of small-scale flight of a g-order vehicle and the launch of a larger vehicle. A laboratory-scale flight experiment can be meaningful in describing an actual satellite launch using active laser control if dynamic equivalence of the translational and rotational flight motions is maintained between the laboratory test and the full-scale flight.The objective of this study is to construct a scaling theory to establish the dynamic equivalence of laboratory-scale and full-scale 6-DOF flight motion. The vehicle mass, vehicle size, moment of inertia, and input energy of the laser beam are theoretically determined to maintain dynamic similarity for a single-pulse reaction based on the translational and rotational equations of motion. The theoretical scaling of the single-pulse reaction is numerically confirmed through a fluid-orbit coupling simulation in which the interactions between a blast wave and a lightcraft are reproduced. The dynamic equivalence of small-scale and larger-scale flight reactions can also be maintained for multiple pulses using the proposed scaling theory.


## 2. Numerical Methods

### 2.1. Outline of Numerical Simulation

The three-dimensional Navier-Stokes equation is numerically solved by integrating the 6-DOF equations of motion to estimate the vehicle reaction for a single-pulse incident. In the single-pulse simulation, the ray-tracing method is used to reproduce beam irradiation with lateral and angular incident offsets, beam reflection by the parabolic mirror, and laser focusing. Using a fluid-orbit coupling simulator, a vehicle reaction data map is constructed with respect to the lateral, angular, and combined lateral-angular offsets to reduce the computational cost of the orbital prediction for multiple pulses. The impulsive reaction driven by the blast wave is interpolated from the impulse data map at the time of the laser incident in the multiple-pulse flight. During free flight of pulse-to-pulse intervals, the 6-DOF equations of motion are only integrated using aerodynamic coefficients only to predict the laser offset for the next pulse. The simulation code developed is used to validate the analytical scaling theory for single and multiple pulses.


### 2.2. Flowfield Calculation

An annular blast wave is generated by the gas breakdown that occurs when a laser beam is irradiated on the parabolic mirror of the lightcraft. A shock wave and an expansion wave are repeatedly reflected by the parabolic mirror and the shroud of the vehicle, generating axial thrust, lateral force, and angular momentum. The flowfield around the vehicle should be solved using the simulation code because the wave reflection complicates the forces generated. In this study, the three-dimensional Navier-Stokes equation, shown below, is numerically integrated to predict the forces and moments acting on the high-speed spinning lightcraft:(1)∂Q∂t+∂(E-E~)∂x+∂(F-F~)∂y+∂(G-G~)∂z=0.Discretization is performed in a cell-centered finite-volume manner [31]. The advection upstream splitting method (AUSM-DV) [32] is employed for numerical flux using the second-order monotone upstream-centered schemes for conservation laws (MUSCL) method [33]. Characteristic variables are chosen as interpolation quantities in the MUSCL method and a minmod limiter is used to fulfill the total variation diminishing (TVD) condition. The numerical flux is calculated with a local rotation matrix to introduce a general curvilinear coordinate [34, 35]. The viscous flux is estimated using the second-order central difference method, and the viscous coefficient is determined using Sutherland’s formula [31]. The Spalart-Allmaras model [36] is employed, assuming fully developed turbulence due to the axial spin of the lightcraft, to introduce the effect of turbulent flow. The initial spin rate is set to 10,000 rpm, based on the results of a previous study [21]. Before the blast wave propagation is calculated, the lower-upper successive overrelaxation (LU-SSOR) method [37, 38] is utilized to obtain a steady solution as an initial flow for the spinning lightcraft to reduce computational time. The laser beam is irradiated to the vehicle using the ray-tracing method when the rotating fluid reaches a steady state, and the unsteady dynamics of the blast wave propagation is predicted. The inviscid eigenvalues are modified using eigenvalues corresponding to the viscous terms to maintain the stability of the scheme because, in the LU-SSOR method, ignoring viscous flux Jacobians results in instability in viscous flows [39]. The stable computation can be conducted with Courant-Friedrichs-Lewy (CFL) numbers up to 20 using the LU-SSOR because of the modification of the viscous-term eigenvalues. The time integration for the unsteady blast wave propagation is performed using the first-order explicit Euler method. A comparison of the fourth-order Runge-Kutta method and the first-order Euler method revealed no obvious discrepancies between them.During the blast wave propagation, the flowfield calculation is combined with the 6-DOF motion calculation of the rigid body of the vehicle because drag is induced by the instantaneous vehicle movement and rotation in the fluid. The translational and angular velocities of the vehicle estimated using the equations of motion are subtracted from the flowfield to account for the effect of motion of the vehicle. The dynamic pressure induced by the rapid vehicle motion acts in the form of lateral and angular drag by coupling the flowfield with the 6-DOF vehicle motion.


### 2.3. Orbital Calculation

The vehicle’s flight motion is determined by solving 6-DOF equations of motion using aerodynamic forces and moments, estimated from the three-dimensional fluid calculation, that occur when the strong interaction with the blast wave pushes and rotates the lightcraft. The conservation of linear and angular momenta for a rigid body results in the 6-DOF equations of motion expressed by the body-fixed coordinate system. The 6-DOF equations of motion to trace the vehicle trajectory are as follows:(2)m(U˙+QW-RV)=Xg+Xa,m(V˙+RU-PW)=Yg+Ya,m(W˙+PV-QU)=Zg+Za,IxxP˙+(Izz-Iyy)QR=Lg+La,IyyQ˙+(Ixx-Izz)RP=Mg+Ma,IzzR˙+(Iyy-Ixx)PQ=Ng+Na.The 6-DOF flight system of the laser propulsion vehicle is illustrated in Figure1. Gravity forces Fg=[Xg,Yg,Zg]T and gravity gradient torques Mg=[Lg,Mg,Ng]T are introduced as external gravity forces and moments acting on the vehicle [40]. When the blast wave strongly interacts strongly with the vehicle surface, the aerodynamic forces Fa=[Xa,Ya,Za]T and the aerodynamic moments Ma=[La,Ma,Na]T should be estimated by integrating the surface pressure distribution and the viscous stress from the flowfield solution in the body-fixed coordinate system. On the other hand, the 6-DOF equations of motion are integrated using only aerodynamic coefficients estimated from the experiments if the interaction between the blast wave and the vehicle surface can be ignored for the free flight in the multiple-pulse simulation. The time integration of the equations of motion is performed using the fourth-order Runge-Kutta method.Figure 1
Velocities and angular velocities of the vehicle in the body-fixed coordinate system.


### 2.4. Fluid-Orbit Coupling Method

The effect of the vehicle motion on the flowfield computation is introduced because accelerative motion of the lightcraft, driven by the strong blast wave, generates additional drag in the translational and rotational directions. The 6-DOF flight calculation is coupled with a flowfield simulation in the body-fixed coordinate system to consider the effect of motion of the vehicle. The time variations of the translational and rotational velocities of the vehicle are obtained by solving the 6-DOF equations of motion, and these variations are fed back to the flowfield calculation at each time step except for an axial spin:(3)unew=u-ΔU-ΔΩ×(rc-rg).A rotational boundary condition for the vehicle surface is employed to model the effects of a rotating fluid due to axial high-speed spin [41]. A flowfield-orbit coupling calculation is employed until the forces from the blast wave can be ignored. For a 32.5 g vehicle, the coupling calculation is terminated at 1,000 μs because the positive and negative pressures of the blast wave pass completely through the vehicle surface, and the generated thrust becomes 0 at 1,000 μs. Negative thrust by an expansion region of the blast wave is included in this calculation, but the net impulse becomes positive for little or no laser offset. The blast wave weakens if the laser incident has lateral or angular offsets, and this weakening reduces or eliminates the thrust.


### 2.5. Energy Deposition Model for Flowfield Calculation

The focal position and gathered energy are not clear if the laser beam is irradiated on the parabolic mirror with lateral or angular offsets. We predict the position of the point at which the reflected laser is focused and the amount of energy to be deposited at the focal point to reproduce off-axis and/or inclined laser incidence. The ray-tracing method is employed to determine the initial energy deposition [42–45].A laser beam of 420 J is irradiated to a vehicle weighing 32.5 g based on the single-pulse calculation to construct an impulse data map of lateral, angular, and combined lateral-angular offsets for comparison with data from previous experiments [21, 25]. The laser beam is numerically divided into ~3,000 rays for the purpose of the energy deposition, and reflections of each ray are geometrically traced to estimate the focal point position in the three-dimensional simulation domain. The laser intensity is set to be spatially uniform within a circular cross-section with a 7.2 cm radius [21]. A pulse duration of 18 μs was used in the flight experiment conducted using the 420-J laser beam with the 32.5 g vehicle. However, the pulse duration can be ignored in the ray tracing. The laser energy is immediately deposited on the flowfield because the pulse duration is much shorter than the time scale of the shock wave interaction. The maximum energy density point is searched as the focal point in a circular plane around the vehicle axis using a three-dimensional computation. The laser energy is input to the source term of the Navier-Stokes equation to drive the blast wave so that the ratio of the internal energy to the kinetic energy is 5 to 2 [46]. A denser grid is used for the ray tracing to improve the accuracy of the energy deposition, and the deposited distribution is interpolated to a coarser grid for the flowfield computation. The efficiency of the laser energy absorption by the plasma is assumed to be 60% to adjust the generated impulse to the experimental data [21, 25].


### 2.6. Vehicle Reaction Data Map and Multiple-Pulse Calculation

The computational cost becomes enormous for multiple-pulse operation if the fluid-orbit coupling calculation is conducted for each pulse incidence with offsets. The coupling simulation is therefore omitted to reduce the computational load by introducing the impulse reaction data map driven by the blast wave [25]. A data map with respect to the parameters Yl-, Zl-, Y^l-, and Z^l- and angular offsets at a reference altitude of 1.194 km [21] was obtained via the fluid-orbit coupling calculation. Figure 2 shows an example of the impulsive reaction data map for laser offsets with an angular offset of 3 degrees. A total of 166 impulse data points were assessed to converge the numerical results of the trajectory simulation for multiple pulses [25]. To obtain the impulse data map, a 420-J laser beam and a 32.5 g vehicle were used for comparison with the experimental flight data [21]. In our previous study [25], the experimental flight of a 32.5 g lightcraft was successfully reproduced using our simulation code based on the impulse data map.Figure 2
Impulse data map against lateral offsets for angular offset of 3 degrees [30].During the multiple-pulse flight, the impulsive reaction of the vehicle is linearly interpolated from the data map at the pulse incident timing. Impulsive increments of the translational and angular velocities induced by the blast wave are added to the current vehicle velocitiesUec and angular velocities Ωec in the Earth-fixed coordinate system, as shown below:(4)Ue=Uec+ΔUe0,Ωe=Ωec+ΔΩe0.However, the axial spin rate is maintained because the axial spin of the vehicle is not affected by the blast wave. After the pulse incidence, for the free-flight interval in which the interaction between the blast wave and the vehicle can be ignored, the 6-DOF equations of motion are temporally integrated using the aerodynamic coefficients to estimate the lateral and angular offsets at the next laser pulse. Based on the results of the experiments, the aerodynamic coefficients of the nonspinning lightcraft are approximated by the following forms [21]:(5)CD=0.3928ζ3-0.8463ζ2+0.0398ζ+0.731,CL=0.1993ζ5-0.3317ζ4-0.2724ζ3+0.5531ζ2-0.0316ζ,CM=0.1598ζ6-0.6619ζ5+0.8794ζ4-0.3585ζ3-0.0583ζ2+0.080ζ.The flight trajectory for the multiple pulses is calculated by repeating the above process.The data map was constructed based on the experimental conditions of the 32.5 g vehicle with the 420-J laser beam. However, the conventional data map cannot be used if the vehicle mass, size, and beam power are increased to correspond to the flight of an actual future satellite. The impulsive data map must be scaled up to predict the flight characteristics of a larger-scale vehicle while maintaining dynamic equivalence for the single-pulse reaction between the conventional scale of the 32.5 g vehicle and the larger scale.


### 2.7. Thrust Decay Model by Beam Expansion and Density Decrease

The thrust decrease caused by beam expansion and density decrease cannot be ignored in multiple-pulse operation because the flight reaches altitudes of tens of kilometers. Analytical modeling of the thrust decrease is conducted for beam expansion and air rarefaction at higher altitudes. The beam radius is modeled using the Gaussian beam profile [47] with a beam wavelength of 1 μm and an initial beam radius of 7.2 cm for the 32.5 g vehicle [15, 21]:(6)w(h)=w01+λhπw0221/2.Because the blast wave generating the thrust obeys the spherically symmetric Sedov solution [46], the translational and angular impulse dependencies on the incident energy decrease are analytically modeled using the shock wave pressure of the Sedov solution [29, 30, 48]. From the Gaussian beam profile and the impulse modeling based on the Sedov solution, (4) are rewritten to describe the decreases in the translational and angular velocities at higher altitudes:(7)Ue=Uec+ρhρ01/2w0whΔUe0,Ωe=Ωec+ρhρ01/2w0whΔΩe0.The ambient density ρ decreases with flight height h, based on the US standard atmosphere 1976 model, in the multiple-pulse simulation.


### 2.8. Laser Position Control Using Optimization Method

Active control of the incident laser is introduced to maintain the stability of the flight because the beam-riding flight becomes unstable if no flight control is employed during multiple-pulse operation. The position and strength of the generated blast wave are controlled by controlling the incident position of the laser beam to adjust the vehicle posture and the translational position in keeping the beam-riding strategy. The optimal position of the pulse incident is obtained by monitoring the flight trajectory and posture. The incident laser position is determined using a genetic algorithm (GA) that has been found to perform well in avoiding locally optimized solutions, in comparison to other optimization schemes. A real-coded-type algorithm [49] is used to code the GA to obtain better solutions than those obtainable using a binary-coded GA. The BLX-α model [49] is employed for crossover in the use of the GA. The minimal generation gap (MGG) model [50] is employed for generation alternation. Mutation is not employed in this study because variety is maintained when the MGG model is combined with the BLX-α model.The GA procedure employed in this study is as follows. Initially, twenty laser position vectors are randomly created. Parents are selected from this group, and the crossover of the BLX-α model is operated twenty times to generate forty children. Then, using the MGG model, two individuals are selected from the group that includes the original parents and children, and parents are replaced by these individuals. This process is repeated until the value of an evaluation function becomes sufficiently small or until the GA loops 100 times. For the GA optimization, the evaluation function Π is defined as follows:(8)Π=R0ifRpast<REth,i,R0≥REth,i,αotherwise.REth,i (i=1,2,…,69) is defined in 69 stages, as shown below:(9)REth,i=(1+0.125(i-1))×0.1.The angular offset α is usually regarded as the evaluation function except when the vehicle passes across the threshold radii REth,i. The centering feedback force is generated by switching the evaluation function to control the translational position of the vehicle if the translational position crosses the threshold radii. To estimate the value of the evaluation function, the 6-DOF equations of motion are integrated for forty laser position vectors at the pulse incident timing. Parallelization by the message-passing interface (MPI) library is conducted in the optimization module to accelerate calculation of the value of the evaluation function.


## 2.1. Outline of Numerical Simulation

The three-dimensional Navier-Stokes equation is numerically solved by integrating the 6-DOF equations of motion to estimate the vehicle reaction for a single-pulse incident. In the single-pulse simulation, the ray-tracing method is used to reproduce beam irradiation with lateral and angular incident offsets, beam reflection by the parabolic mirror, and laser focusing. Using a fluid-orbit coupling simulator, a vehicle reaction data map is constructed with respect to the lateral, angular, and combined lateral-angular offsets to reduce the computational cost of the orbital prediction for multiple pulses. The impulsive reaction driven by the blast wave is interpolated from the impulse data map at the time of the laser incident in the multiple-pulse flight. During free flight of pulse-to-pulse intervals, the 6-DOF equations of motion are only integrated using aerodynamic coefficients only to predict the laser offset for the next pulse. The simulation code developed is used to validate the analytical scaling theory for single and multiple pulses.


## 2.2. Flowfield Calculation

An annular blast wave is generated by the gas breakdown that occurs when a laser beam is irradiated on the parabolic mirror of the lightcraft. A shock wave and an expansion wave are repeatedly reflected by the parabolic mirror and the shroud of the vehicle, generating axial thrust, lateral force, and angular momentum. The flowfield around the vehicle should be solved using the simulation code because the wave reflection complicates the forces generated. In this study, the three-dimensional Navier-Stokes equation, shown below, is numerically integrated to predict the forces and moments acting on the high-speed spinning lightcraft:(1)∂Q∂t+∂(E-E~)∂x+∂(F-F~)∂y+∂(G-G~)∂z=0.Discretization is performed in a cell-centered finite-volume manner [31]. The advection upstream splitting method (AUSM-DV) [32] is employed for numerical flux using the second-order monotone upstream-centered schemes for conservation laws (MUSCL) method [33]. Characteristic variables are chosen as interpolation quantities in the MUSCL method and a minmod limiter is used to fulfill the total variation diminishing (TVD) condition. The numerical flux is calculated with a local rotation matrix to introduce a general curvilinear coordinate [34, 35]. The viscous flux is estimated using the second-order central difference method, and the viscous coefficient is determined using Sutherland’s formula [31]. The Spalart-Allmaras model [36] is employed, assuming fully developed turbulence due to the axial spin of the lightcraft, to introduce the effect of turbulent flow. The initial spin rate is set to 10,000 rpm, based on the results of a previous study [21]. Before the blast wave propagation is calculated, the lower-upper successive overrelaxation (LU-SSOR) method [37, 38] is utilized to obtain a steady solution as an initial flow for the spinning lightcraft to reduce computational time. The laser beam is irradiated to the vehicle using the ray-tracing method when the rotating fluid reaches a steady state, and the unsteady dynamics of the blast wave propagation is predicted. The inviscid eigenvalues are modified using eigenvalues corresponding to the viscous terms to maintain the stability of the scheme because, in the LU-SSOR method, ignoring viscous flux Jacobians results in instability in viscous flows [39]. The stable computation can be conducted with Courant-Friedrichs-Lewy (CFL) numbers up to 20 using the LU-SSOR because of the modification of the viscous-term eigenvalues. The time integration for the unsteady blast wave propagation is performed using the first-order explicit Euler method. A comparison of the fourth-order Runge-Kutta method and the first-order Euler method revealed no obvious discrepancies between them.During the blast wave propagation, the flowfield calculation is combined with the 6-DOF motion calculation of the rigid body of the vehicle because drag is induced by the instantaneous vehicle movement and rotation in the fluid. The translational and angular velocities of the vehicle estimated using the equations of motion are subtracted from the flowfield to account for the effect of motion of the vehicle. The dynamic pressure induced by the rapid vehicle motion acts in the form of lateral and angular drag by coupling the flowfield with the 6-DOF vehicle motion.


## 2.3. Orbital Calculation

The vehicle’s flight motion is determined by solving 6-DOF equations of motion using aerodynamic forces and moments, estimated from the three-dimensional fluid calculation, that occur when the strong interaction with the blast wave pushes and rotates the lightcraft. The conservation of linear and angular momenta for a rigid body results in the 6-DOF equations of motion expressed by the body-fixed coordinate system. The 6-DOF equations of motion to trace the vehicle trajectory are as follows:(2)m(U˙+QW-RV)=Xg+Xa,m(V˙+RU-PW)=Yg+Ya,m(W˙+PV-QU)=Zg+Za,IxxP˙+(Izz-Iyy)QR=Lg+La,IyyQ˙+(Ixx-Izz)RP=Mg+Ma,IzzR˙+(Iyy-Ixx)PQ=Ng+Na.The 6-DOF flight system of the laser propulsion vehicle is illustrated in Figure1. Gravity forces Fg=[Xg,Yg,Zg]T and gravity gradient torques Mg=[Lg,Mg,Ng]T are introduced as external gravity forces and moments acting on the vehicle [40]. When the blast wave strongly interacts strongly with the vehicle surface, the aerodynamic forces Fa=[Xa,Ya,Za]T and the aerodynamic moments Ma=[La,Ma,Na]T should be estimated by integrating the surface pressure distribution and the viscous stress from the flowfield solution in the body-fixed coordinate system. On the other hand, the 6-DOF equations of motion are integrated using only aerodynamic coefficients estimated from the experiments if the interaction between the blast wave and the vehicle surface can be ignored for the free flight in the multiple-pulse simulation. The time integration of the equations of motion is performed using the fourth-order Runge-Kutta method.Figure 1
Velocities and angular velocities of the vehicle in the body-fixed coordinate system.


## 2.4. Fluid-Orbit Coupling Method

The effect of the vehicle motion on the flowfield computation is introduced because accelerative motion of the lightcraft, driven by the strong blast wave, generates additional drag in the translational and rotational directions. The 6-DOF flight calculation is coupled with a flowfield simulation in the body-fixed coordinate system to consider the effect of motion of the vehicle. The time variations of the translational and rotational velocities of the vehicle are obtained by solving the 6-DOF equations of motion, and these variations are fed back to the flowfield calculation at each time step except for an axial spin:(3)unew=u-ΔU-ΔΩ×(rc-rg).A rotational boundary condition for the vehicle surface is employed to model the effects of a rotating fluid due to axial high-speed spin [41]. A flowfield-orbit coupling calculation is employed until the forces from the blast wave can be ignored. For a 32.5 g vehicle, the coupling calculation is terminated at 1,000 μs because the positive and negative pressures of the blast wave pass completely through the vehicle surface, and the generated thrust becomes 0 at 1,000 μs. Negative thrust by an expansion region of the blast wave is included in this calculation, but the net impulse becomes positive for little or no laser offset. The blast wave weakens if the laser incident has lateral or angular offsets, and this weakening reduces or eliminates the thrust.


## 2.5. Energy Deposition Model for Flowfield Calculation

The focal position and gathered energy are not clear if the laser beam is irradiated on the parabolic mirror with lateral or angular offsets. We predict the position of the point at which the reflected laser is focused and the amount of energy to be deposited at the focal point to reproduce off-axis and/or inclined laser incidence. The ray-tracing method is employed to determine the initial energy deposition [42–45].A laser beam of 420 J is irradiated to a vehicle weighing 32.5 g based on the single-pulse calculation to construct an impulse data map of lateral, angular, and combined lateral-angular offsets for comparison with data from previous experiments [21, 25]. The laser beam is numerically divided into ~3,000 rays for the purpose of the energy deposition, and reflections of each ray are geometrically traced to estimate the focal point position in the three-dimensional simulation domain. The laser intensity is set to be spatially uniform within a circular cross-section with a 7.2 cm radius [21]. A pulse duration of 18 μs was used in the flight experiment conducted using the 420-J laser beam with the 32.5 g vehicle. However, the pulse duration can be ignored in the ray tracing. The laser energy is immediately deposited on the flowfield because the pulse duration is much shorter than the time scale of the shock wave interaction. The maximum energy density point is searched as the focal point in a circular plane around the vehicle axis using a three-dimensional computation. The laser energy is input to the source term of the Navier-Stokes equation to drive the blast wave so that the ratio of the internal energy to the kinetic energy is 5 to 2 [46]. A denser grid is used for the ray tracing to improve the accuracy of the energy deposition, and the deposited distribution is interpolated to a coarser grid for the flowfield computation. The efficiency of the laser energy absorption by the plasma is assumed to be 60% to adjust the generated impulse to the experimental data [21, 25].


## 2.6. Vehicle Reaction Data Map and Multiple-Pulse Calculation

The computational cost becomes enormous for multiple-pulse operation if the fluid-orbit coupling calculation is conducted for each pulse incidence with offsets. The coupling simulation is therefore omitted to reduce the computational load by introducing the impulse reaction data map driven by the blast wave [25]. A data map with respect to the parameters Yl-, Zl-, Y^l-, and Z^l- and angular offsets at a reference altitude of 1.194 km [21] was obtained via the fluid-orbit coupling calculation. Figure 2 shows an example of the impulsive reaction data map for laser offsets with an angular offset of 3 degrees. A total of 166 impulse data points were assessed to converge the numerical results of the trajectory simulation for multiple pulses [25]. To obtain the impulse data map, a 420-J laser beam and a 32.5 g vehicle were used for comparison with the experimental flight data [21]. In our previous study [25], the experimental flight of a 32.5 g lightcraft was successfully reproduced using our simulation code based on the impulse data map.Figure 2
Impulse data map against lateral offsets for angular offset of 3 degrees [30].During the multiple-pulse flight, the impulsive reaction of the vehicle is linearly interpolated from the data map at the pulse incident timing. Impulsive increments of the translational and angular velocities induced by the blast wave are added to the current vehicle velocitiesUec and angular velocities Ωec in the Earth-fixed coordinate system, as shown below:(4)Ue=Uec+ΔUe0,Ωe=Ωec+ΔΩe0.However, the axial spin rate is maintained because the axial spin of the vehicle is not affected by the blast wave. After the pulse incidence, for the free-flight interval in which the interaction between the blast wave and the vehicle can be ignored, the 6-DOF equations of motion are temporally integrated using the aerodynamic coefficients to estimate the lateral and angular offsets at the next laser pulse. Based on the results of the experiments, the aerodynamic coefficients of the nonspinning lightcraft are approximated by the following forms [21]:(5)CD=0.3928ζ3-0.8463ζ2+0.0398ζ+0.731,CL=0.1993ζ5-0.3317ζ4-0.2724ζ3+0.5531ζ2-0.0316ζ,CM=0.1598ζ6-0.6619ζ5+0.8794ζ4-0.3585ζ3-0.0583ζ2+0.080ζ.The flight trajectory for the multiple pulses is calculated by repeating the above process.The data map was constructed based on the experimental conditions of the 32.5 g vehicle with the 420-J laser beam. However, the conventional data map cannot be used if the vehicle mass, size, and beam power are increased to correspond to the flight of an actual future satellite. The impulsive data map must be scaled up to predict the flight characteristics of a larger-scale vehicle while maintaining dynamic equivalence for the single-pulse reaction between the conventional scale of the 32.5 g vehicle and the larger scale.


## 2.7. Thrust Decay Model by Beam Expansion and Density Decrease

The thrust decrease caused by beam expansion and density decrease cannot be ignored in multiple-pulse operation because the flight reaches altitudes of tens of kilometers. Analytical modeling of the thrust decrease is conducted for beam expansion and air rarefaction at higher altitudes. The beam radius is modeled using the Gaussian beam profile [47] with a beam wavelength of 1 μm and an initial beam radius of 7.2 cm for the 32.5 g vehicle [15, 21]:(6)w(h)=w01+λhπw0221/2.Because the blast wave generating the thrust obeys the spherically symmetric Sedov solution [46], the translational and angular impulse dependencies on the incident energy decrease are analytically modeled using the shock wave pressure of the Sedov solution [29, 30, 48]. From the Gaussian beam profile and the impulse modeling based on the Sedov solution, (4) are rewritten to describe the decreases in the translational and angular velocities at higher altitudes:(7)Ue=Uec+ρhρ01/2w0whΔUe0,Ωe=Ωec+ρhρ01/2w0whΔΩe0.The ambient density ρ decreases with flight height h, based on the US standard atmosphere 1976 model, in the multiple-pulse simulation.


## 2.8. Laser Position Control Using Optimization Method

Active control of the incident laser is introduced to maintain the stability of the flight because the beam-riding flight becomes unstable if no flight control is employed during multiple-pulse operation. The position and strength of the generated blast wave are controlled by controlling the incident position of the laser beam to adjust the vehicle posture and the translational position in keeping the beam-riding strategy. The optimal position of the pulse incident is obtained by monitoring the flight trajectory and posture. The incident laser position is determined using a genetic algorithm (GA) that has been found to perform well in avoiding locally optimized solutions, in comparison to other optimization schemes. A real-coded-type algorithm [49] is used to code the GA to obtain better solutions than those obtainable using a binary-coded GA. The BLX-α model [49] is employed for crossover in the use of the GA. The minimal generation gap (MGG) model [50] is employed for generation alternation. Mutation is not employed in this study because variety is maintained when the MGG model is combined with the BLX-α model.The GA procedure employed in this study is as follows. Initially, twenty laser position vectors are randomly created. Parents are selected from this group, and the crossover of the BLX-α model is operated twenty times to generate forty children. Then, using the MGG model, two individuals are selected from the group that includes the original parents and children, and parents are replaced by these individuals. This process is repeated until the value of an evaluation function becomes sufficiently small or until the GA loops 100 times. For the GA optimization, the evaluation function Π is defined as follows:(8)Π=R0ifRpast<REth,i,R0≥REth,i,αotherwise.REth,i (i=1,2,…,69) is defined in 69 stages, as shown below:(9)REth,i=(1+0.125(i-1))×0.1.The angular offset α is usually regarded as the evaluation function except when the vehicle passes across the threshold radii REth,i. The centering feedback force is generated by switching the evaluation function to control the translational position of the vehicle if the translational position crosses the threshold radii. To estimate the value of the evaluation function, the 6-DOF equations of motion are integrated for forty laser position vectors at the pulse incident timing. Parallelization by the message-passing interface (MPI) library is conducted in the optimization module to accelerate calculation of the value of the evaluation function.


## 3. Results

### 3.1. Dynamic Scaling Theory for Single Pulse

In our previous studies [25–27] based on experiments with the lightcraft [10, 21], beam-riding flight of a 32.5 g lightcraft was achieved using active laser control by GA optimization. To examine the feasibility of using this technique to control the flight of a larger vehicle in an actual satellite launch, it is necessary to scale up the impulse data map while maintaining the dynamic equivalence between the small scale of the 32.5 g vehicle and the scale of a larger vehicle. The dynamic equivalence theory can correlate the laboratory scale to a full-scale satellite launch, and we can provide more weight to the flight test at the laboratory scale to achieve a full-scale launch. When the vehicle mass and size are increased, theoretical scaling is conducted to determine the input energy and the moment of inertia required to maintain dynamic similarity for a single pulse, based on the equations of motion for the translational and angular directions.A lightcraft vehicle of the 200-5/6 type [10] is scaled up to 1 kg and 10 kg masses to construct the dynamic scaling law. Depending on the mass increase, the vehicle size is enlarged while the mass density and geometric similarity are maintained, and the length scale and area scale are expanded in proportion to the mass using the following equations:(10)L=mm01/3L0,S=mm02/3S0,where the reference flight system of the 32.5 g vehicle is indexed by 0. The scaling law for the moment of inertia is as follows:(11)Ii=mm0ξIi0,where ξ is the scaling parameter for the moment of inertia and i=xx,yy, and zz. The reference moments of inertia Ixx0, Iyy0, and Izz0 are set to 465 g·cm2 for the x-axis, 603 g·cm2 for the y-axis, and 603 g·cm2 for the z-axis, based on previous works [21, 25].During the interaction between the vehicle and the blast wave after the pulse incident, the translational and angular equations of motion are given by the following equations:(12)dudt=-g+Φt,dωdt=Φa,where u and ω are the translational and angular velocities of the vehicle. The thrust term Φt and the angular momentum term Φa driven by the blast wave are defined using the thrust Fs and the external moment Ms:(13)Φt≡Fsm,Φa≡MsIi.We can ignore effects of aerodynamic drag because the aerodynamic drag is smaller than the impulsive thrust during the blast wave interaction. To maintain dynamic equivalence between the scaled-up system and the reference system, the laser energy at one pulse and the moment of inertia should be determined to satisfy the following equations:(14)ΦtΦt0=1,ΦaΦa0=1.The thrust Fs given by the blast wave is modeled using the three-dimensional Sedov solution:(15)Fs=1Δts∫0Δts‍ps-p0rs2πdt,where rs, ps, and Δts are the front position of the shock wave, the shock wave pressure, and the interaction interval between the vehicle body and the shock wave, respectively, as given by the following equations:(16)rs=ξ0μEρ01/5t2/5,ps=825ρ0γ+1ξ02μEρ02/5t-6/5,Δts=trs=L=ξ0-5/2μEρ0-1/2L5/2.rs refers to the propagation distance of the shock wave from an explosion source. The value of rs is analytically determined as a function of the input energy and the time, based on the three-dimensional Sedov solution. We assume that the area surrounded by rs causes the propulsive thrust of the vehicle in (15). Any dependence of Δts on the pulse duration can be ignored because the pulse duration of the laser beam is much shorter than the interaction time scale of Δts. The pulse duration is several tens of μs, which is much shorter than Δts, which is on the order of several hundreds or thousands of μs. By integrating (15) using (16), the thrust generated by the blast wave can be expressed in the following form:(17)Fs=η1EL-1-η2L2,where(18)η1≡815πξ05ρ0γ+1μρ0,η2≡59p0π.The external moment Ms driven by the blast wave is modeled using the thrust modeling of (17):(19)Ms=kmsFsL.By substituting (11), (17), and (19) into (14), the dynamic equivalence conditions of (14) can be rewritten for the translational and angular motions as follows:(20)ΦtΦt0=η1EL-1-η2L2η1E0L0-1-η2L02m0m=1,ΦaΦa0=η1EL-1-η2L2η1E0L0-1-η2L02m0mmm04/3-ξ=1.Because the input energy of the single pulse E and the multiplier for the moment of inertia ξ are calculated from (20), E and Ii are obtained from the following equations to maintain the dynamic equivalence:(21)E=mm04/3E0+L03η2η1mm01-mm01/3,(22)Ii=Ii0mm04/3.For the 1 kg vehicle, the input energy E and the moments of inertia Ixx, Iyy, and Izz are determined to maintain the dynamic equivalence, based on (21) and (22), as follows:(23)E=1.3×104J,Ixx=4.5×10-3kg·m2,Iyy=Izz=5.8×10-3kg·m2.For the 10 kg vehicle, E, Ixx, Iyy, and Izz are calculated as follows:(24)E=1.4×105J,Ixx=9.7×10-2kg·m2,Iyy=Izz=1.3×10-1kg·m2.To calculate the values of (23) and (24), we use the actual parameter values for the single-pulse flight: m0=0.0325 kg, ρ0=1.09 kg/m3, E0=420 J, μ=0.6, L0=0.11 m, ξ0=1.03, p0=87,782 Pa, γ=1.4, Ixx0=4.653×10-5 kg·m2, and Iyy0=6.027×10-5 kg·m2. Dynamic scaling for the pulse duration is not necessary because the pulse duration is much shorter than the time scale of the shock wave interaction. The pulse duration remains very short in comparison with the time scale of the shock wave interaction, even if the pulse duration increases with increasing laser power because the interaction time scale also increases with the system scale-up.


### 3.2. Validation of Scaling Law by Numerical Simulation

Based on the theoretical scaling for the single pulse, equivalent motions can be obtained for other scales for the translational and rotational directions, if the input energy and the moment of inertia are determined using (21) and (22). The dynamic equivalence given using the scaling theory described by (21) and (22) should be validated by single-pulse simulation using our fluid-orbit coupling code. In the three-dimensional simulation, the vehicle size and the moment of inertia are determined using (10) and (22) to examine the dynamic equivalence for the 1 kg and 10 kg vehicles. The vehicle shape is maintained for all of the cases regardless of the vehicle size. A single pulse with the energy given by (21) is irradiated on the vehicle with a lateral y-offset of (m/m0)1/3×(-5) mm using the ray-tracing method, and the blast wave propagation is simulated for the 1 kg and 10 kg vehicles (Figures 3(a) and 3(b)).Figure 3
Pressure distribution of blast wave driven by laser incident at 200μs. The laser beam is irradiated using the ray-tracing method with an initial offset of (m/m0)1/3×(-5) mm for 1 kg and 10 kg vehicles, and the fluid-orbit coupling calculation is conducted.
(a)
1 kg vehicle with laser beam of1.3×104 J
(b)
10 kg vehicle with laser beam of1.4×105 JThe translational and rotational vehicle motions are traced by combining the orbit calculation with the fluid calculation. For the coupling calculation, the time scale is normalized by the characteristic time of the self-similar solution of the blast wave [46]:(25)t*=mm01/3L05/2ρ0μE1/2.Dynamic equivalence is successfully maintained for the 32.5 g, 1 kg, and 10 kg cases in the translational and angular directions (Figures 4(a) and 4(b)). The numerical simulation results show the correctness of the analytical scaling theory described by (21) and (22). The larger-scale flight can be predicted for a single pulse based on the laboratory-scale flight test using the scaling theory described by (21) and (22) because motion equivalence is maintained for the different scales. By conducting the small-scale flight, we can examine the flight of the larger-scale system using the scaling theory for a single pulse.Figure 4
Time evolutions of translational and angular velocities during blast wave interaction for 32.5 g, 1 kg, and 10 kg vehicles. The vehicle size, laser energy, and moment of inertia are determined using (10), (21), and (22) for 1 kg and 10 kg vehicles in the fluid-orbit coupling simulation.
(a)
Translational velocities
(b)
Angular velocities


### 3.3. Dynamic Equivalence for Multiple Pulses

The scaling theory was constructed to maintain dynamic equivalence for a single pulse, but there is no guarantee that dynamic equivalence will be maintained for a multiple-pulse flight. Actually, the dynamic equivalence is broken down for multiple pulses by the aerodynamic drag during the free flight of the pulse-to-pulse interval because the effective aerodynamic drag decreases as a result of the greater inertia of the heavier vehicle [29]. If the pulse frequency is held constant at 100 Hz, the vehicle velocity increases with increasing mass because the greater inertia of a heavier vehicle accelerates the vehicle during free flight in the pulse-to-pulse intervals. The pulse frequency should be decreased depending on the mass increase to maintain the dynamic equivalence. By correcting the pulse frequency of the repetitive operation, we can maintain the dynamic equivalence of the smaller and larger systems for multiple pulses. The following 1-DOF equation of motion applies to the reference system during multiple-pulse operation:(26)m0dudt=-m0g-12ρ0u2S0CD+m0ΔU0f0.Similarly, the 1-DOF equation of motion for the larger system is as follows:(27)mdvdt=-mg-12ρ0v2SCD+mΔUf.In a multiple-pulse flight, the vehicle velocity reaches a terminal velocity of uterm because of the balance between the gravity force, the aerodynamic drag, and the thrust. The dynamic equivalence conditions for multiple pulses can therefore be expressed as follows:(28)dudt=dvdt=0,u=v=uterm.Because dynamic equivalence is maintained for a single pulse if the incident energy obeys (21), the impulsive velocity for a single pulse can be expressed as follows:(29)ΔU=ΔU0.Using (26)–(29) and the area scaling S=(m/m)2/3S0 given by (10), the pulse frequency is determined to obtain dynamic equivalence for multiple pulses:(30)f=1m0ΔU012ρ0uterm2S0CDmm0-1/3-1+f0.Using (30), values of f=35.7 Hz and f=19.6 Hz that maintain dynamic equivalence are obtained for the 1 kg and 10 kg vehicles, respectively. We used the following parameter values for the multiple-pulse flight: m0=0.0325 kg, S0=0.0095 m2, CD=0.731, uterm=40 m/s, f0=100 Hz, ρ0=1.09 kg/m3, and ΔU0=1.93 m/s.The 6-DOF multiple-pulse simulations are conducted for the 32.5 g, 1 kg, and 10 kg vehicles with an initial offset of(m/m0)1/3×5 mm to confirm the dynamic equivalence of the multiple pulses. The conventional data map for the impulsive motion can be employed to simulate the multiple-pulse trajectory of the 1 kg and 10 kg vehicles because dynamic equivalence is maintained for each case by following the scaling theory described by (21) and (22). The pulse frequency is determined using (30) for the larger-scale flight. Perturbations (e.g., wind, system delay for the active control, and positioning error of the laser incident) are ignored for the sake of focusing on the motion equivalence for multiple-pulse flight. The stability of the flight of the laser-boosted vehicle is maintained using the active laser control technique based on the GA optimization [26, 27]. The centering threshold REth,i of the evaluation function Π is scaled up for the larger vehicle to maintain the equivalence of the posture control: (31)REth,i=1+0.125i-1×0.1f0f.Dynamic equivalence is maintained in the multiple-pulse operation because identical velocities are obtained for the 32.5 g, 1 kg, and 10 kg cases by adjusting the pulse frequency based on (30) (Figure 5(a)). The vehicle velocity decreases when centering feedback is provided by the active laser control because the propulsive force is used as the lateral feedback force to control the translational position of the vehicle. This results in velocity fluctuations during the flight. The velocity fluctuations increase with increasing vehicle mass because the centering control requires more time for a larger vehicle because of the decrease in the pulse frequency. Flights at km-order altitudes are successfully achieved for the three cases with keeping the beam-riding strategy by suppressing the angular offset using active laser control (Figure 5(b)). Trajectory fluctuations also increase with increasing vehicle mass because the pulse frequency decreases and the centering threshold of the active control REth,i increases with system scale-up. The spatial and temporal scales of the flight trajectories increase with system scale-up, but the flight dynamics of the scaled-up vehicle is approximately equivalent. Without conducting a flight experiment or a flight simulation with a larger vehicle, the multiple-pulse flight reaction of the larger vehicle can be estimated by conducting a flight test with a small-scale vehicle because dynamic equivalence is maintained between the larger and smaller systems by applying the scaling theory to the input energy, the geometric scale of the vehicle, and the pulse frequency. We can assign more weight to the small-scale flight using the scaling theory, which correlates smaller and larger flight systems, because the results of small-scale flight tests are equivalent to those for larger-scale flights as scaled up using the scaling theory.Figure 5
Flight trajectories and vehicle velocities using scaled-up data.
(a)
Time evolution of vehicle velocities
(b)
6-DOF flight trajectories


## 3.1. Dynamic Scaling Theory for Single Pulse

In our previous studies [25–27] based on experiments with the lightcraft [10, 21], beam-riding flight of a 32.5 g lightcraft was achieved using active laser control by GA optimization. To examine the feasibility of using this technique to control the flight of a larger vehicle in an actual satellite launch, it is necessary to scale up the impulse data map while maintaining the dynamic equivalence between the small scale of the 32.5 g vehicle and the scale of a larger vehicle. The dynamic equivalence theory can correlate the laboratory scale to a full-scale satellite launch, and we can provide more weight to the flight test at the laboratory scale to achieve a full-scale launch. When the vehicle mass and size are increased, theoretical scaling is conducted to determine the input energy and the moment of inertia required to maintain dynamic similarity for a single pulse, based on the equations of motion for the translational and angular directions.A lightcraft vehicle of the 200-5/6 type [10] is scaled up to 1 kg and 10 kg masses to construct the dynamic scaling law. Depending on the mass increase, the vehicle size is enlarged while the mass density and geometric similarity are maintained, and the length scale and area scale are expanded in proportion to the mass using the following equations:(10)L=mm01/3L0,S=mm02/3S0,where the reference flight system of the 32.5 g vehicle is indexed by 0. The scaling law for the moment of inertia is as follows:(11)Ii=mm0ξIi0,where ξ is the scaling parameter for the moment of inertia and i=xx,yy, and zz. The reference moments of inertia Ixx0, Iyy0, and Izz0 are set to 465 g·cm2 for the x-axis, 603 g·cm2 for the y-axis, and 603 g·cm2 for the z-axis, based on previous works [21, 25].During the interaction between the vehicle and the blast wave after the pulse incident, the translational and angular equations of motion are given by the following equations:(12)dudt=-g+Φt,dωdt=Φa,where u and ω are the translational and angular velocities of the vehicle. The thrust term Φt and the angular momentum term Φa driven by the blast wave are defined using the thrust Fs and the external moment Ms:(13)Φt≡Fsm,Φa≡MsIi.We can ignore effects of aerodynamic drag because the aerodynamic drag is smaller than the impulsive thrust during the blast wave interaction. To maintain dynamic equivalence between the scaled-up system and the reference system, the laser energy at one pulse and the moment of inertia should be determined to satisfy the following equations:(14)ΦtΦt0=1,ΦaΦa0=1.The thrust Fs given by the blast wave is modeled using the three-dimensional Sedov solution:(15)Fs=1Δts∫0Δts‍ps-p0rs2πdt,where rs, ps, and Δts are the front position of the shock wave, the shock wave pressure, and the interaction interval between the vehicle body and the shock wave, respectively, as given by the following equations:(16)rs=ξ0μEρ01/5t2/5,ps=825ρ0γ+1ξ02μEρ02/5t-6/5,Δts=trs=L=ξ0-5/2μEρ0-1/2L5/2.rs refers to the propagation distance of the shock wave from an explosion source. The value of rs is analytically determined as a function of the input energy and the time, based on the three-dimensional Sedov solution. We assume that the area surrounded by rs causes the propulsive thrust of the vehicle in (15). Any dependence of Δts on the pulse duration can be ignored because the pulse duration of the laser beam is much shorter than the interaction time scale of Δts. The pulse duration is several tens of μs, which is much shorter than Δts, which is on the order of several hundreds or thousands of μs. By integrating (15) using (16), the thrust generated by the blast wave can be expressed in the following form:(17)Fs=η1EL-1-η2L2,where(18)η1≡815πξ05ρ0γ+1μρ0,η2≡59p0π.The external moment Ms driven by the blast wave is modeled using the thrust modeling of (17):(19)Ms=kmsFsL.By substituting (11), (17), and (19) into (14), the dynamic equivalence conditions of (14) can be rewritten for the translational and angular motions as follows:(20)ΦtΦt0=η1EL-1-η2L2η1E0L0-1-η2L02m0m=1,ΦaΦa0=η1EL-1-η2L2η1E0L0-1-η2L02m0mmm04/3-ξ=1.Because the input energy of the single pulse E and the multiplier for the moment of inertia ξ are calculated from (20), E and Ii are obtained from the following equations to maintain the dynamic equivalence:(21)E=mm04/3E0+L03η2η1mm01-mm01/3,(22)Ii=Ii0mm04/3.For the 1 kg vehicle, the input energy E and the moments of inertia Ixx, Iyy, and Izz are determined to maintain the dynamic equivalence, based on (21) and (22), as follows:(23)E=1.3×104J,Ixx=4.5×10-3kg·m2,Iyy=Izz=5.8×10-3kg·m2.For the 10 kg vehicle, E, Ixx, Iyy, and Izz are calculated as follows:(24)E=1.4×105J,Ixx=9.7×10-2kg·m2,Iyy=Izz=1.3×10-1kg·m2.To calculate the values of (23) and (24), we use the actual parameter values for the single-pulse flight: m0=0.0325 kg, ρ0=1.09 kg/m3, E0=420 J, μ=0.6, L0=0.11 m, ξ0=1.03, p0=87,782 Pa, γ=1.4, Ixx0=4.653×10-5 kg·m2, and Iyy0=6.027×10-5 kg·m2. Dynamic scaling for the pulse duration is not necessary because the pulse duration is much shorter than the time scale of the shock wave interaction. The pulse duration remains very short in comparison with the time scale of the shock wave interaction, even if the pulse duration increases with increasing laser power because the interaction time scale also increases with the system scale-up.


## 3.2. Validation of Scaling Law by Numerical Simulation

Based on the theoretical scaling for the single pulse, equivalent motions can be obtained for other scales for the translational and rotational directions, if the input energy and the moment of inertia are determined using (21) and (22). The dynamic equivalence given using the scaling theory described by (21) and (22) should be validated by single-pulse simulation using our fluid-orbit coupling code. In the three-dimensional simulation, the vehicle size and the moment of inertia are determined using (10) and (22) to examine the dynamic equivalence for the 1 kg and 10 kg vehicles. The vehicle shape is maintained for all of the cases regardless of the vehicle size. A single pulse with the energy given by (21) is irradiated on the vehicle with a lateral y-offset of (m/m0)1/3×(-5) mm using the ray-tracing method, and the blast wave propagation is simulated for the 1 kg and 10 kg vehicles (Figures 3(a) and 3(b)).Figure 3
Pressure distribution of blast wave driven by laser incident at 200μs. The laser beam is irradiated using the ray-tracing method with an initial offset of (m/m0)1/3×(-5) mm for 1 kg and 10 kg vehicles, and the fluid-orbit coupling calculation is conducted.
(a)
1 kg vehicle with laser beam of1.3×104 J
(b)
10 kg vehicle with laser beam of1.4×105 JThe translational and rotational vehicle motions are traced by combining the orbit calculation with the fluid calculation. For the coupling calculation, the time scale is normalized by the characteristic time of the self-similar solution of the blast wave [46]:(25)t*=mm01/3L05/2ρ0μE1/2.Dynamic equivalence is successfully maintained for the 32.5 g, 1 kg, and 10 kg cases in the translational and angular directions (Figures 4(a) and 4(b)). The numerical simulation results show the correctness of the analytical scaling theory described by (21) and (22). The larger-scale flight can be predicted for a single pulse based on the laboratory-scale flight test using the scaling theory described by (21) and (22) because motion equivalence is maintained for the different scales. By conducting the small-scale flight, we can examine the flight of the larger-scale system using the scaling theory for a single pulse.Figure 4
Time evolutions of translational and angular velocities during blast wave interaction for 32.5 g, 1 kg, and 10 kg vehicles. The vehicle size, laser energy, and moment of inertia are determined using (10), (21), and (22) for 1 kg and 10 kg vehicles in the fluid-orbit coupling simulation.
(a)
Translational velocities
(b)
Angular velocities


## 3.3. Dynamic Equivalence for Multiple Pulses

The scaling theory was constructed to maintain dynamic equivalence for a single pulse, but there is no guarantee that dynamic equivalence will be maintained for a multiple-pulse flight. Actually, the dynamic equivalence is broken down for multiple pulses by the aerodynamic drag during the free flight of the pulse-to-pulse interval because the effective aerodynamic drag decreases as a result of the greater inertia of the heavier vehicle [29]. If the pulse frequency is held constant at 100 Hz, the vehicle velocity increases with increasing mass because the greater inertia of a heavier vehicle accelerates the vehicle during free flight in the pulse-to-pulse intervals. The pulse frequency should be decreased depending on the mass increase to maintain the dynamic equivalence. By correcting the pulse frequency of the repetitive operation, we can maintain the dynamic equivalence of the smaller and larger systems for multiple pulses. The following 1-DOF equation of motion applies to the reference system during multiple-pulse operation:(26)m0dudt=-m0g-12ρ0u2S0CD+m0ΔU0f0.Similarly, the 1-DOF equation of motion for the larger system is as follows:(27)mdvdt=-mg-12ρ0v2SCD+mΔUf.In a multiple-pulse flight, the vehicle velocity reaches a terminal velocity of uterm because of the balance between the gravity force, the aerodynamic drag, and the thrust. The dynamic equivalence conditions for multiple pulses can therefore be expressed as follows:(28)dudt=dvdt=0,u=v=uterm.Because dynamic equivalence is maintained for a single pulse if the incident energy obeys (21), the impulsive velocity for a single pulse can be expressed as follows:(29)ΔU=ΔU0.Using (26)–(29) and the area scaling S=(m/m)2/3S0 given by (10), the pulse frequency is determined to obtain dynamic equivalence for multiple pulses:(30)f=1m0ΔU012ρ0uterm2S0CDmm0-1/3-1+f0.Using (30), values of f=35.7 Hz and f=19.6 Hz that maintain dynamic equivalence are obtained for the 1 kg and 10 kg vehicles, respectively. We used the following parameter values for the multiple-pulse flight: m0=0.0325 kg, S0=0.0095 m2, CD=0.731, uterm=40 m/s, f0=100 Hz, ρ0=1.09 kg/m3, and ΔU0=1.93 m/s.The 6-DOF multiple-pulse simulations are conducted for the 32.5 g, 1 kg, and 10 kg vehicles with an initial offset of(m/m0)1/3×5 mm to confirm the dynamic equivalence of the multiple pulses. The conventional data map for the impulsive motion can be employed to simulate the multiple-pulse trajectory of the 1 kg and 10 kg vehicles because dynamic equivalence is maintained for each case by following the scaling theory described by (21) and (22). The pulse frequency is determined using (30) for the larger-scale flight. Perturbations (e.g., wind, system delay for the active control, and positioning error of the laser incident) are ignored for the sake of focusing on the motion equivalence for multiple-pulse flight. The stability of the flight of the laser-boosted vehicle is maintained using the active laser control technique based on the GA optimization [26, 27]. The centering threshold REth,i of the evaluation function Π is scaled up for the larger vehicle to maintain the equivalence of the posture control: (31)REth,i=1+0.125i-1×0.1f0f.Dynamic equivalence is maintained in the multiple-pulse operation because identical velocities are obtained for the 32.5 g, 1 kg, and 10 kg cases by adjusting the pulse frequency based on (30) (Figure 5(a)). The vehicle velocity decreases when centering feedback is provided by the active laser control because the propulsive force is used as the lateral feedback force to control the translational position of the vehicle. This results in velocity fluctuations during the flight. The velocity fluctuations increase with increasing vehicle mass because the centering control requires more time for a larger vehicle because of the decrease in the pulse frequency. Flights at km-order altitudes are successfully achieved for the three cases with keeping the beam-riding strategy by suppressing the angular offset using active laser control (Figure 5(b)). Trajectory fluctuations also increase with increasing vehicle mass because the pulse frequency decreases and the centering threshold of the active control REth,i increases with system scale-up. The spatial and temporal scales of the flight trajectories increase with system scale-up, but the flight dynamics of the scaled-up vehicle is approximately equivalent. Without conducting a flight experiment or a flight simulation with a larger vehicle, the multiple-pulse flight reaction of the larger vehicle can be estimated by conducting a flight test with a small-scale vehicle because dynamic equivalence is maintained between the larger and smaller systems by applying the scaling theory to the input energy, the geometric scale of the vehicle, and the pulse frequency. We can assign more weight to the small-scale flight using the scaling theory, which correlates smaller and larger flight systems, because the results of small-scale flight tests are equivalent to those for larger-scale flights as scaled up using the scaling theory.Figure 5
Flight trajectories and vehicle velocities using scaled-up data.
(a)
Time evolution of vehicle velocities
(b)
6-DOF flight trajectories


## 4. Conclusion

To predict the flight reactions of full-scale vehicles, a scaling theory was constructed to maintain the dynamic equivalence of laboratory-scale and full-scale laser propulsion flights. Translational and angular equations of motion are used to construct a dynamic scaling law for a single-pulse incident during blast wave interaction. The theoretical scaling law was validated by numerical simulation of the blast wave interaction for a single pulse using a fluid-orbital coupling simulator, and equivalent vehicle motions for the smaller- and larger-scale flights were obtained for the translational and angular directions. Dynamic equivalence was also maintained for multiple pulses by adjusting the repetition frequency of the incident laser to correct for the effect of aerodynamic drag during the free flight of the pulse-to-pulse interval. The flight of a larger vehicle can be estimated for single- and multiple-pulse operations using the flight data for a small-scale vehicle, using the scaling theory to maintain dynamic equivalence between the two systems. A laboratory-scale flight test can be used to estimate the flight characteristics of a full-scale vehicle using dynamic scaling because the small-scale flight is equivalent to the larger-scale flight when expanded using the proposed scaling theory.



---

*Source: 801371-2015-03-08.xml*