The Method for Optimum Design of Water Rocket Flight Stability Performance Conditions Using CAE with T Method and Robust Parameter Design

Eiji Toma, Yoshihiro Ito

National Institute of Technology, Tomakomai College, Tomakomai, Japan.

DOI: 10.4236/jamp.2021.911172                

Abstract

A water rocket is a rocket system that obtains thrust by injecting water with compressed air of up to about 8 atmospheres. It is usually manufactured using a pressure-resistant PET bottle. The mechanical elements and principles contained in the water rocket have much in common with the actual small rocket system, and are suitable as educational and research teaching materials in the field of mechanics. Especially in the field of disaster prevention and mitigation, the use of water rockets is being researched and developed as a rescue tool in the event of a flood or earthquake as a disaster countermeasure. However, since the water rocket is a flying object based on the mechanical principle, it is important to ensure the accuracy and stability of the flight path. In this paper, a mechanical simulator is developed with a numerical calculation program based on the mechanical consideration of water rocket flight performance. In addition, the correlation between the flight distance obtained in the simulation and the estimated flight distance is analyzed by applying a multivariate analysis method and verifying the validity of the flight distance calculated from the result. Based on the verification results, we will apply a statistical optimization method to approach the optimization of flight stability performance conditions for water rockets.

Keywords

Flying Principle, Multivariate Analysis, T Method, Robust Parameter Design, Flight Stability, Energetic SN Ratio

1. Introduction

A water rocket is a rocket system that obtains thrust by injecting water with compressed air of up to about 8 atmospheres. Normally, pressure-resistant PET bottles are used as the material, they are small and easy to manufacture (Figure 1). Unlike model rockets that use gunpowder, there is no danger of fire, and as an example of practical research in the industrial world, there is a method for wire drawing in mountainous areas by the “water rocket wire drawing method” [1]. In addition, the mechanical elements and principles contained in the water rocket have much in common with the actual rocket system, making it an appropriate system for educational and research teaching materials in the field of mechanics.

However, since the water rocket is a flying object based on the mechanical principle, it is important to ensure the accuracy and stability of the flight path. When manufacturing and launching in a short time, such as at an event, a water rocket with a simple configuration is manufactured, so it is necessary to consider an airframe configuration that will fly stably no matter who makes it. For that purpose, it is important to perform numerical analysis of the flight trajectory and quantitative analysis of flight stability by the mechanical simulation to evaluate the flight characteristics in advance. However, at present, basic aerodynamic data is indispensable when analyzing the flight trajectory. There are many empirical factors in the technical guidelines for aerodynamic design and structural design which are important for achieving stable flight. In addition, many experiments are required under various parameter conditions such as the amount of water to be put in the bottle and the pressure of compressed air in order to obtain the flight distance. In reality, most of the optimum values are the data obtained by repeating the launch, and the scientific basis for this has not been clarified [2]. In previous studies, examples of optimization research on the flight distance of water rockets have been reported [3] [4]. However, no verification cases have been reported regarding the estimation of flight distance and the optimization of flight stability performance conditions using a mechanical simulator.

Figure 1. Main structural model of water rocket. A water rocket consists of a fuel tank, a nose cone, and an injection section (including fins, injection nozzles, and blades).

Therefore, in this research, we will develop the original simulator that mechanically considers the flight performance of water rockets. The correlation between the flight distance obtained by the simulation and the estimated flight distance is analyzed using a multivariate analysis method (T Method) to verify the validity of the flight distance.

Based on the verification results, important parameters that affect the flight performance of the water rocket are extracted by a statistical optimization method (Robust Parameter Design), and the optimum flight conditions for obtaining the flight stability of the water rocket are clarified [5].

2. Mechanical Consideration of Water Rocket

2.1. Equation of Motion

The equation of motion that determines the motion of water rocket thrust flight is expressed as follows.

$m\frac{\text{d}v}{\text{d}t}=F-k{\left(v\right)}^2-mg\sin \left(\theta \right)$ (1)

$v\frac{\text{d}\theta }{\text{d}t}=-g\cos \left(\theta \right),\frac{\text{d}m}{\text{d}t}=-\beta$ (2)

$\frac{\text{d}x}{\text{d}t}=v\cos \left(\theta \right),\frac{\text{d}y}{\text{d}t}=v\sin \left(\theta \right)$ (3)

where:

m: Mass of rocket (including mass of water); θ: Orbital angle of the rocket;

F: Thrust force; v: Speed of the rocket; R: Air resistance force;

β: Water discharge per unit time; g: Gravitational acceleration.

2.2. Force Acting on the Water Rocket and Flight Trajectory

Various forces act when a water rocket flies (Figure 2). Even if the rocket is launched upward at a certain angle, the orbit gradually becomes downward. This phenomenon is called “Gravity turn”. The direction of the gravity turn is determined by the resultant force of the rocket’s thrust, air resistance, and the three main forces of gravity acting on the rocket, so the orbit gradually changes downward. In reality, in addition to the three resultant forces in the flight motion of the rocket, complex forces such as lift due to the progress of the rocket, wind force, and rotation motion of the rocket generated by the blades act.

2.3. Thrust by Water Injection

Unlike rockets that use general fuel, the thrust of water rockets is due to the ejection of water and the thrust of air that occurs after the water runs out. This thrust is generated by the acceleration of water due to the pressure difference when the pressure P₀ in the bottle is greater than the atmospheric pressure P_a (Figure 3).

First, consider the thrust obtained by the ejection of water. When water (density: ρ) is ejected from the mouth of cross-sectional area A to the rocket at a

Figure 2. Thrust model. F₁ : Thrust, F₂ : Air resistance, F₃ : Gravity, F₄ : Component force in direction of gravity thrust, F₅ : Gravity thrust and vertical component force, F₆ : Difference between thrust and reverse force, F₇ : All combined forces.

Figure 3. Thrust model. F : Thrust, P₀ : Bottle internal pressure, P_a : Atmospheric pressure, A₀ : Tank cross-sectional area, u₀ : Flow velocity, A : Injection port cross-sectional area, u: Injection flow velocity.

speed u, the mass dm of the water ejected in a minute time dt is given by the following equation.

$dm=\rho Au\cdot dt$ (4)

The impulse $F\cdot dt$ received by the rocket is calculated by ejecting the mass dm backward from the rocket at a speed u,

$F\cdot dt=dm\cdot u$ (5)

By substituting Equation (4) into this and rearranging it, the equation for the thrust F of the rocket can be obtained.

$F=\rho A{u}^2$ (6)

Assuming that the cross-sectional area in the tank is A₀, the flow velocity is u₀, and that at the injection port is A, u, the following equation holds from Bernoulli’s equation and the continuity equation, which are the basic equations of fluid mechanics.

$\frac{1}{2}\rho u^2+P_a=\frac{1}{2}\rho u_0^2+\rho gh+P_0$ (7)

$u\cdot A=u_0\cdot A_0$ (8)

where, P_a is the atmospheric pressure, and P₀ is the pressure inside the tank, g is the gravitational acceleration and h is the height from the injection port to the water surface in the tank. From Equations (7) and (8), the following equation can be obtained.

$u=\frac{1}{\sqrt{1-{\left(\frac{A}{A_0}\right)}^2}}\sqrt{\frac{2\left(P_0-P_a\right)}{\rho }+2gh}$ (9)

By substituting u in Equation (9) into the thrust F equation, the thrust due to water injection becomes the following equation.

$F=\frac{2A}{1-{\left(\frac{A}{A_0}\right)}^2}\left(P_0-P_a+\rho gh\right)$ (10)

The pressure P₀ in the tank decreases rapidly with the injection of water, but this change can be considered as an adiabatic change.

2.4. Thrust by Air Injection

Compressed air still remains in the tank at the end of the water injection. By ejecting the remaining compressed air, the rocket gains thrust and is further accelerated. Here, let us consider the case where the pressure and density of air are P₀ and ρ₀, respectively, and the air flows out of the tank at the pressure P_a. If the pressure and density of gas at the injection port are P₁ and ρ₁, respectively, the outflow velocity u₁ can be obtained from the equation for compressible fluid.

$u_1=\sqrt{\frac{2\gamma }{\gamma -1}\cdot \frac{P_0}{{\rho }_0}\left\{1-{\left(\frac{P_1}{P_0}\right)}^{\frac{\gamma -1}{\gamma }}\right\}}$ (11)

where, $\gamma$ is the specific heat ratio of air, and if $\gamma =1.4$ is substituted and arranged, the outflow velocity u₁ can be given by the following equation.

$u_1=\sqrt{7\cdot \frac{P_0}{{\rho }_0}\left\{1-{\left(\frac{P_1}{P_0}\right)}^{\frac{2}{7}}\right\}}$ (12)

where, the pressure P₁ at the injection port is not always equal to the external pressure P_a. When the flow velocity u₁ of the injected gas does not exceed the speed of sound, the injection port pressure P₁ can be considered to be almost equal to the external pressure P_a. However, as the pressure difference between the inside and outside becomes large and the sound velocity cannot be exceeded even if the flow velocity increases, the injection port pressure P₁ does not fall to the external pressure P_a.

$P_1={\left(\frac{2}{\gamma +1}\right)}^{\frac{\gamma }{\gamma -1}}\cdot P_0=0.528P_0>P_a$ (13)

Such a state is called choking (blockage of flow). The thrust F_a which is the reaction of momentum given to the blast air is ${\rho }_{1}A{u}_1^2$ as in the case of water injection, and the following equation is obtained by using the adiabatic expansion equation.

$F_a={\rho }_{1}A{u}_1^2=7A{P}_1\left\{{\left(\frac{P_0}{P_1}\right)}^{\frac{2}{7}}-1\right\}$ (14)

when the air injection port pressure is higher than atmospheric pressure, thrust F_p (called pressure thrust) due to the pressure difference is added to this as shown in the following equation.

$F_p=A\left(P_1-P_a\right)$ (15)

In summary, the total thrust F of the water rocket is obtained by the following equation.

$F=F_a+F_p=7A{P}_1\left\{{\left(\frac{P_0}{P_1}\right)}^{\frac{2}{7}}-1\right\}+A\left(P_1-P_a\right)$ (16)

where the following conditional expression can be substituted for P₁.

$0.528\times P_0\ge P_a\to P_1=0.528\times P_0$ (17)

$0.528\times P_0<P_a\to P_1=P_a$ (18)

3. Development of Mechanical Simulator

Solve the equations of motion and the relations of thrust in the previous section by Euler method, and develop a simulator using the EXCEL macro function numerical calculation program (Table 1). The simulator is applied as a tool for extracting the optimum flight conditions for obtaining the maximum flight distance and flight stability of the water rocket (Figure 4) [6] [7].

Input the initial parameter values from the program start mode and calculate the horizontal reach of the water rocket until the vertical coordinate (y axis) becomes 0 on the numerical calculation program (Figure 5).

4. Correlation Verification of Estimated Flight Distance by Multivariate Analysis

In this study, we apply a multivariate analysis method to the correlation between

Table 1. Basic formula. Solve the equations of motion and the relations of thrust in the previous section by Euler method, and develop a simulator using the EXCEL macro function numerical calculation program.

Figure 4. Simulator operation flowchart. The simulator is applied as a tool for extracting the optimum flight conditions for obtaining the maximum flight distance and flight stability of the water rocket.

Figure 5. Simulator operation screen. Input the initial parameter values from the program start mode and calculate the horizontal reach of the water rocket until the vertical coordinate (y axis) becomes 0 on the numerical calculation program.

the flight distance obtained in the simulation and the estimated flight distance, and verify the validity of the flight distance calculated from the analysis results [8]. For the estimation of the flight distance, the “T method” which is a quality engineering method that can estimate the output value in time series from multivariate data is adopted [9] [10]. The T method is one of the multivariate analysis methods and is a mathematical method that estimates one objective variable (output value) from multiple item variables (continuous variables) in the same way as multiple regression analysis. Comprehensive estimation is possible from the relationship between each item and the output value [11] [12].

4.1. Computation Formula for the T Method

The T Method defines the Unit Space where the output value is in the medium position and homogeneous (densely populated). The computation procedure of the T Method is explained below [13] [14].

4.1.1. Definition of the Unit Space and Computation of the Average of Relevant Items and Outputs

Let’s suppose that n number of data have been obtained for the Unit Space (Table 2). All the items of the data must be in the same dimension as image density or

Table 2. Data for the unit space and averages values of the items and outputs. All the items of the data must be in same dimension as image density or must be no dimension data.

must be no dimension data. From the n number of samples in the Unit Space, we find average values $\overline{x_1},\overline{x_2},\cdots ,\overline{x_k}$ and average output value $\bar{y}=M_0$ for all items. Accordingly, the average values work out as follows:

$\overline{x_j}=\frac{1}{n}\left(x_{1j}+x_{2j}+\cdots +x_{nj}\right)\left(j=1,2,\cdots ,k\right)$ (19)

$\bar{y}=M_0=\frac{1}{n}\left(y_1+y_2+\cdots +y_n\right)$ (20)

Table 2 also shows these average values. One of the average values obtained from the n members of the Unit Space is the center of the Unit Space.

4.1.2. Definition of Signal Data

All data items marked l, left unselected for the Unit Space are treated as Signal data (Table 3). “Signal data” refers to all data used for finding the proportional coefficient β and SN ratio η, which will be discussed in Section 4.1.4.

4.1.3. Normalization of Signal Data

Signal data is normalized using the average values of items and the output values of samples in the Unit Space (Table 4). Normalization is performed by subtracting the average value $\overline{x_j}$ of item j in the Unit Space from value ${x^{\prime }}_{ij}$ of item j of the i-th Signal data.

$X_{ij}={x^{\prime }}_{ij}-\overline{x_j}\left(i=1,2,\cdots ,l;j=1,2,\cdots ,k\right)$ (21)

Likewise, normalization is performed by subtracting average value M₀ of the output from the Unit Space from output value ${y^{\prime }}_i$ of thei-th Signal data.

4.1.4. Computation of Proportional Coefficient β and SN Ratio η (in Duplicate Ratio) for All Items

We will next compute proportional coefficient β and SN ratio η for all items. How the computation is performed is explained with item 1 as an example:

The SN ratio is expressed by the ratio of the signal amount S and the noise N, and is calculated by the signal (active component of the input energy)/noise (harmful component that did not work effectively as an output).

Proportional coefficient;

Table 3. Signal data. “Signal data” refers to all data used for finding the proportional coefficient β and SN ratio η.

Table 4. Normalized Signal data. Signal data is normalized using the average values of items and the output values of samples in the unit space.

${\beta }_1=\frac{M_1{X}_{11}+M_2{X}_{21}+\cdots +M_l{X}_{l1}}{r}$ (23)

SN ratio;

${\eta }_1=\{\begin{array}{l}\frac{\frac{1}{r}\left(S_{\beta 1}-V_{e1}\right)}{V_{e1}}\left(\text{when}{S}_{\beta 1}>V_{e1}\right)\\ 0\left(\text{when}{S}_{\beta 1}\le V_{e1}\right)\end{array}$ (24)

where:

Effective divider;

$r=M_1^2+M_2^2+\cdots +M_l^2$ (25)

Total variation;

$S_{T1}=X_{11}^2+X_{21}^2+\cdots +X_{l1}^2\left(f=l\right)$ (26)

Variation of Proportional term;

$S_{\beta 1}=\frac{{\left(M_1{X}_{11}+M_2{X}_{21}+\cdots +M_l{X}_{l1}\right)}^2}{r}\left(f=1\right)$ (27)

Error variation;

$S_{e1}=S_{T1}-S_{\beta 1}\left(f=l-1\right)$ (28)

Error variance;

$V_{e1}=\frac{S_{e1}}{l-1}$ (29)

From item 2 up to itemk, we will likewise find proportional coefficient β and SN ratio η. This operation yields the results that are shown in Table 5.

4.1.5. Computation, Signal by Signal, of Integrated Estimate Value $\hat{M}$ of Output

An item-by-item estimated value is found for each piece of Signal data using the proportional coefficient β and SN ratio η, item by item. The estimated value of the output of item 1 for the i-th Signal Data is:

${\hat{M}}_{i1}=\frac{X_{i1}}{{\beta }_1}$ (30)

An estimation is likewise made of item 2 through item l for the i-th Signal data. And finally, integration of the result is performed by weighting it with ${\eta }_1,{\eta }_2,\cdots ,{\eta }_k$, which is the estimated measure of precision of each item. Thus, the integrated estimate value of the output ${\hat{M}}_i$ of the i-th Signal data becomes:

${\hat{M}}_i=\frac{{\eta }_1\times \frac{X_{i1}}{{\beta }_1}+{\eta }_2\times \frac{X_{i2}}{{\beta }_2}+\cdots +{\eta }_k\times \frac{X_{ik}}{{\beta }_k}}{{\eta }_1+{\eta }_2+\cdots +{\eta }_k}\left(i=1,2,\cdots ,l\right)$ (31)

Table 6 shows the real values (measured values) of the Signal data $M_1,M_2,\cdots ,M_l$ and the integrated estimate values ${\hat{M}}_1,{\hat{M}}_2,\cdots ,{\hat{M}}_l$.

4.1.6. Computation of Integrated Estimate SN Ratio

The Integrated Estimate SN Ratio is computed using the following equation based on Table 6 [15].

Integrated Estimate SN Ratio;

$\eta =10{\log }_{10}\left[\frac{\frac{1}{r}\left(S_{\beta }-V_e\right)}{V_e}\right]\left(\text{db}\right)$ (32)

Table 5. Proportional coefficient β and SN ratio η, item by item.

Table 6. Measured values and integrated estimate values of signal data.

where:

Linear equation;

$L=M_1{\hat{M}}_1+M_2{\hat{M}}_2+\cdots +M_l{\hat{M}}_l$ (33)

Effective divider;

$r=M_1^2+M_2^2+\cdots +M_l^2$ (34)

Total variation;

$S_T={\hat{M}}_1^2+{\hat{M}}_2^2+\cdots +{\hat{M}}_l^2\left(f=l\right)$ (35)

Variation of proportional term;

$S_{\beta }=\frac{L^2}{r}\left(f=1\right)$ (36)

Error variation;

$S_e=S_T-S_{\beta }\left(f=l-1\right)$ (37)

Error variance;

$V_{e1}=\frac{S_e}{l-1}$ (38)

4.2. Correlation Verification of Estimated Flight Distance

In this section, we verify the validity of the flight distance calculated from the correlation between the flight distance obtained on the simulation and the estimated flight distance by applying the T method. Table 7 shows the basic data of the item variables for the flight distance obtained in the simulation. With the objective variable as the flight distance, the item variables are bottle capacity (No. 1), guide rail length (No. 2), launch angle (No. 3), water volume ratio (No. 4), and air pressure (No. 5), Air resistance (No. 6), and 6 items in total.

Table 8 and Figure 6 show the calculation results of the flight distance (measured value) obtained on the simulation and the estimated flight distance (estimated value) by the T method, and their correlation distribution. Table 9 shows the calculation results of various coefficients that evaluate the correlation between the measured value and the estimated value. From the calculation results of the correlation coefficient (=0.96) and the coefficient of determination (=0.93), a high correlation can be confirmed in the measured/estimated value of the flight distance. In the T method, it is possible to make a comprehensive estimate from the relationship between each item and the output value, and analyze the contribution to the estimated value for each item [16].

Table 10 shows the calculation process of the integrated estimate SN ratio in the calculation of the estimated flight distance. Table 11 shows the diagnostic results of the contribution to the estimated value for each item extracted when calculating the flight distance obtained in the simulation. The contribution of each item is the result of evaluation by the integrated estimate SN ratio [17] [18].

The contribution of the item is evaluated by how bad (low) is the integrated

Table 7. Basic data of item variables for flight distance obtained on simulation.

Table 8. Flying distance obtained on simulation (measured value) and flying distance by T method (estimated value).

Table 9. Calculation results of correlation coefficient of measured value/estimated value.

Table 10. Calculation process of the integrated estimate SN ratio.

Figure 6. Correlation distribution diagram of measured/estimated values.

Table 11. Contribution of items in calculation of estimated flight distance.

estimate SN ratio, so a two-level orthogonal array is used from the statistical empirical rule. By using an orthogonal array, reliability can be improved by comparing which items are important for estimation accuracy and how important they are by the integrated estimate SN ratio [19] [20].

Figure 7 and Figure 8 show the percentage of contribution and the factor-effect diagram of valid items and ineffective items calculated from the SN ratio difference in the calculation of the estimated flight distance. From the analysis results in Table 11, Figure 7 and Figure 8, the items that have a strong influence on the estimated flight distance are the bottle capacity (No. 1), the water volume ratio (No. 4), and the air pressure (No. 5). Since a water rocket generally manufactured uses a pressure resistant PET bottle, it is necessary to consider the pressure resistance strength. Therefore, as shown in Table 11, the flight conditions for obtaining the flight distance of the water rocket by numerical calculation simulation can be estimated.

5. Analysis of Flight Stability by “Robust Parameter Design”

By adopting the “dynamic simulation system” for the flight characteristics of the

Figure 7. Contribution to the estimated flight distance. Vertical axis; SN ratio difference, horizontal axis; item.

Figure 8. Factor effect diagram of items. Vertical axis; SN ratio, horizontal axis; item level.

water rocket developed in Chapter 3, this chapter approaches the optimization of robust flight stability performance conditions by robust parameter design (Hereafter referred to as RPD) which is the central method of robust design. For the evaluation characteristics of flight stability performance, dynamic characteristic evaluation (Energetic SN ratio; see Section 5.4) is adopted. RPD is a design method based on statistical theory that investigates the relationship between design parameters and the quality of an object and minimizes its influence [21]. It can be said that this is an extremely effective methodology for evaluating the flight stability performance of water rockets that have various parameters such as bottle capacity and water volume ratio.

In this study, we extract factors that affect the flight performance of water rockets and verify the optimization of parameter levels to obtain flight stability.

5.1. Conceptual Diagram of RPD

RPD is a method of reducing variations in characteristics and functions by designing control factors so that they are robust against noise factors such as usage environment conditions. It is a general but proven development philosophy focused on improving the reliability of a process or product. Improving reliability requires that RPD principles be an early and integral part of the development cycle. The objective is to make the end-product immune to factors that could adversely affect reliability. As shown in Figure 9, a general RPD methodology requires that four factors be considered in the design process: Signal factors, Response factors, Noise factors, and Control factors. The basic idea of RPD is to attenuate the effects of noise factors by finding effective control and noise factor interactions in the design. This is an economical and effective method to reduce variability simply by changing the level of control factors [22] [23]. Robustness is a concept used in a wide range of fields such as information engineering such as computer programs, statistics, economics, and biology, in addition to design and manufacturing. Robustness means ensuring “robust stability” in the performance of systems and products and improving reliability. In the RPD experiment, the noise factor intentionally causes variation in the combination of system parameters selected as the control factor, and the robustness is improved by optimizing the level of the strong control factor that can counter the variation.

5.2. Objective Function and Quality Characteristics

In RPD, it is important to understand the objective function of the system. The objective function is the role that the system should do, and the quality characteristics are divided into the following four types. By picking up the objective function, the design information can be utilized in other similar systems, and the versatility of the technology is expanded.

1) Preferably small characteristic (nonnegative and smaller the better).

Figure 9. Conceptual diagram of RPD.

2) Preferably large characteristic (nonnegative and larger the better).

3) Preferably target characteristic (with target value).

4) Dynamic characteristic (Correlation between input and output, output stability).

In this study, the flight stability of the water rocket is the objective function. We apply the “Dynamic characteristics (energetic SN ratio)” with the quality characteristics for extracting the optimum parameter conditions to obtain it as the flight distance [24]. The following data analysis tools are used for quality characteristic evaluation.

Ÿ JUSE. StatWorks/V5 Quality Engineering (Nikka Giken co.,Ltd. in Japan).

5.3. Functionality Evaluation

The ideal function of many system technologies is that the output changes linearly with the change of the input. The total output characteristic value (total work) of the system must be proportional to the energy or work. Especially when a water rocket flies, it requires input energy to generate thrust, and there is energy such as air resistance that disturbs its function.

Therefore, in this study, we focus on the fact that the ideal function that represents the input-output relationship has a linear relationship with energy conversion, and define the ideal function state as shown in Figure 10. The purpose of functionality evaluation is to efficiently extract parameters that bring the function of the water rocket closer to the ideal state and reduce variations and fluctuations. The difference from the ideal state can be considered as some kind of energy loss. That is, the function of the water rocket is consumed in addition to the intended output, and the flight stability is affected as a result. In RPD, it is important to define the function, especially the functionality based on energy conversion, and evaluate its stability efficiently [25].

Figure 10. Ideal function state. Vertical axis; flight distance, horizontal axis; energy parameters.

5.4. Formula of Energetic SN Ratio

The SN ratio (Signal-to-Noise ratio) which has been used in RPD and represents a measure of variation, may change depending on the size (range) of the input signal and the number of data. There is a problem that engineers must keep this in mind when acquiring and analyzing data. The idea of the SN ratio in RPD is to evaluate the variation and non-linearity from the ideal functional state (y = βM) as the poor functional stability. The “Energetic SN ratio” adopted in this study is a new SN ratio announced in June 2008 by a research member of the Kansai Quality Engineering Society which is a certified regional study group of the Quality Engineering Society in japan [26].

All technologies involve energy conversion and transmission, and energy perspective is important in technological development and evaluation of their technical levels. Energetic SN ratio is a technical quality evaluation measurement based on how various energies are used. Technical quality means that more energy input to the system can be stably used for a stable output indefinitely. In the ideal functional state when the water rocket flies, energy is decomposed into effective components (S_β) and harmful components (S_N), and the energetic SN ratio is calculated by the ratio of these. In the calculation of the energetic SN ratio, it can be expressed in decibel value as 10 times the common logarithm like the conventional type SN ratio.

The basic calculation formula of the energetic SN ratio η_E is shown below. As shown in Table 12, it is assumed that the output y_ij at the noise factor level $i\left(i=1,2,\cdots ,n\right)$ and the signal factor level $j\left(j=1,2,\cdots ,k\right)$ is obtained.

1) Total variation component: S_T

This is the sum of the squares of nk pieces of data y_ij and shows the variation from y = 0.

$S_T={\displaystyle \underset{i=1}{\overset{n}{\sum }}{\displaystyle \underset{j=1}{\overset{k}{\sum }}{y}_{ij}^2}}$ (39)

2) Average slope size: β_N0

For some dynamic characteristics of the input signal, the slope of the noise factor N_i level β_Ni is obtained by considering the effective component as the average slope.

${\beta }_{Ni}=\frac{{\displaystyle {\sum }_{j=1}^k{M}_j{y}_{ij}}}{{\displaystyle {\sum }_{j=1}^k{M}_j^2}}=\frac{L_{Ni}}{r}$ (40)

Table 12. Correspondence table for SN ratio calculation.

where:

r: Effective divider, L_Ni: Linear format

$r\equiv {\displaystyle \underset{j=1}{\overset{k}{\sum }}{M}_j^2}$ (41)

$L_{Ni}\equiv {\displaystyle \underset{j=1}{\overset{k}{\sum }}{M}_j{y}_j}$ (42)

The sum of squares of the signal, which is the denominator of β in Equation (40), is called the effective divisor r, and the sum of products of the numerator signal and output is called the linear form L_Ni. From the slope β_Ni for each noise factor level, the average slope β_N0 is given by the following equation.

${\beta }_{N0}=\frac{{\displaystyle {\sum }_{i=1}^n{\beta }_{Ni}}}{n}$ (43)

3) Average slope fluctuation (active ingredient):S_β

The variation of the average slope is expressed as the sum of squares of the size $y={\beta }_{N0}M$.

$S_{\beta }=n\underset{j=1}{\overset{k}{{\displaystyle \sum }}}{\left({\beta }_{N0}{M}_j\right)}^2=n{\beta }_{N0}^2\underset{j=1}{\overset{k}{{\displaystyle \sum }}}{\left(M_j\right)}^2=nr{\beta }_{N0}^2$ (44)

4) Harmful component: S_N

The harmful component represents the variation of the difference between each data and the proportional equation of the average slope, and is obtained by subtracting the effective component from the total variation component.

$S_N={\displaystyle \underset{i=1}{\overset{n}{\sum }}{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\left(y_{ij}-{\beta }_{N0}{M}_j\right)}^2}}=S_T-S_{\beta }$ (45)

5) Energetic SN ratio: η_E

${\eta }_E=10{\log }_{10}\left(\frac{S_{\beta }}{S_T-S_{\beta }}\right)=10{\log }_{10}\left(\frac{S_{\beta }}{S_N}\right)\left(\text{db}\right)$ (46)

5.5. Analysis of Flight Stability by Dynamic Characteristic Evaluation

In this section, we will verify by “dynamic characteristic evaluation (Energetic SN ratio)” with the purpose of extracting the optimum flight conditions for the flight stability of the water rocket. The dynamic characteristic is a characteristic that examines the output by changing the input. It is possible to evaluate the condition that minimizes the two-dimensional variation between input and output among various combinations of factors. It is desirable that the input and output signals are in a proportional relationship even under various usage conditions and environments. It is characterized by a more robustness improvement effect in terms of reproducibility, advancement, and versatility [25].

5.5.1. Determination of Level Tableand Various Factors

Table 13 shows the level table of the control factors in this experiment. The

Table 13. Control factors.

control factor sets four parameters: “BOTTLE capacity”, “guide RAIL length”, “LUNCH angle”, and “WATER volume ratio”. Table 14 shows the noise factors extracted in this experiment, and the air resistance of the rocket body is set to no wind and head wind. Table 15 shows the experimental design and data table of this experiment, and sets the signal factors (initial air pressure; 3atm, 5atm, 7atm) as inputs.

5.5.2. Calculation of SN Ratio and Slopeβ

In RPD, the SN ratio of dynamic characteristics is used for functionality evaluation, and it is ideal that stable output can be obtained even if there are various errors that affect functional characteristics. To evaluate the functionality, change the input signal and check the linearity of the output [27]. Figure 11 shows an input/output relationship graph of initial air pressure (horizontal axis: atm.) and flight distance (vertical axis: m) by simulation for each experiment number. From the graph, it can be seen that the relationship between the input/output initial air pressure and the flight distance is a characteristic factor in the evaluation of the presence or absence of linearity related to flight stability. In this experiment, “Energetic SN ratio based on zero point proportionality” is adopted.

Table 16 shows the values of the SN ratio and the slope β from the output data. Here, the slope β is an index of the degree of variation with respect to the magnitude of the slope of the output average, and is represented by an antilogarithm.

5.5.3. Analysis of Factor Effect and Variance

Figure 12 shows the factor-effect diagram of the SN ratio and slope β of the control factors. The factor-effect table (Table 17) shows the effect or combinations of factors on characteristic values, and the diagram showing it is the factor-effect diagram [28].

The meaning of this diagram is that factors that have a large effect on the SN ratio and do not affect the output slope β are extracted, and it can be judged that the level is highly significant for the flight stability performance of the water rocket. In Figure 12, the optimum levels are the factors “BOTTLE capacity: 1.0 L”, “LAUNCH angle: 55”, “guide RAIL length: 0.38 m”, and “WATER volume ratio: 30%” in the circled parts.

As for the bottle capacity, it can be estimated that a small diameter with less

Table 14. Noise factors.

Table 15. L9 orthogonal table and experimental data.

Table 16. SN ratio and slope β.

Table 17. Factor effect table.

Figure 11. Relationship between initial pressure and flight distance.

Figure 12. Factor effect diagram.

air resistance is more advantageous as a condition for obtaining flight stability. In addition, the obtained characteristic values vary from the average value to each level. This variation scale is evaluated using the analysis of variance table (Table 18). From this table, it can be inferred that the bottle capacity and water volume ratio of factors with a large dispersion ratio are factors that have a high degree of influence on flight stability performance [29]. Table 19 shows the optimal combination of control factor levels for obtaining flight stability in the experimental range with respect to the flight performance of the water rocket.

5.5.4. Reliability Evaluation of L9 Orthogonal Array Experiment

As shown in Table 19, the estimated value of the optimum condition is higher than that of the benchmark condition (BM condition), and the gain is secured. Based on the results, we confirm whether the L9 orthogonal array experimental results are reliable (additive).

The evaluation method is to obtain the maximum and minimum values of the SN ratio obtained from the L9 orthogonal table from statistical empirical rules, and the difference from the SN ratio of the BM condition estimated from the value of 10% of the difference, if it is within the range of 10%, it can be judged to be reliable [30]. Table 20 shows the reliability evaluation results. By evaluating the variation in the proportional relationship between the flight distance and the air pressure, it is a desirable condition for the stability of the flight trajectory, and it can be judged that the L9 orthogonal array experiment is reliable.

5.5.5. Reproducibility Evaluation by Confirmation Experiment

The purpose of the confirmation experiment is to confirm whether the output is stable even if the time, place, environment, etc. change compared to the orthogonal array experiment. Table 21 shows the evaluation results of the confirmation experiment data. In the confirmation experiment, the estimated value of the

Table 18. Variance table.

Table 19. Optimization results.

Table 20. Reliability check of L9 experiment.

Table 21. Confirmation experiment data.

Table 22. Reproducibility of confirmation experiment.

SN ratio of the selected optimum condition and the BM condition is obtained, and the same two types of experiments as the orthogonal array experiment are performed again to obtain the SN ratio. Next, the difference between the two estimates and the gain difference in the confirmation experiment are compared to see if the experiment is reproducible.

Table 22 shows the evaluation results of reproducibility in the confirmation experiment. In RPD, if the gain reproducibility (%) is within the range of 70% to 130%, or if the gain difference (db) is within ±3 db, it can be judged that there is reproducibility [30]. In this experimental result, the gain reproducibility is 87% and the gain difference is 0.56 db, so it can be judged that there is “reproducibility”. Furthermore, the improvement rate of flight stability with respect to the BM condition is about 40% (39%), indicating that the extracted optimum conditions are properly selected and the water rocket is much closer to the ideal functional state.

Figure 13. Comparison of flight stability and linearity.

Figure 13 is a graph comparing the linear relationship between the input and output of the optimum condition and the BM condition in the confirmation experiment. It shows the change in output (flying distance) with respect to input (initial air pressure), and it can be seen that the optimum condition has a more stable linear proportional relationship than the BM condition.

6. Conclusions

In this study, we developed a simulator with a numerical calculation program based on a mechanical consideration of the flight performance of a water rocket and calculated the changes over time in the injection speed of internal pressure, water, and compressed air. We considered how the flight distance changes when the flight conditions are changed, and estimated the flight conditions to obtain the flight distance by solving the equation of motion. Next, the flight conditions that can be set in the simulation of the water rocket flight trajectory were extracted. We considered how the flight distance changes when the flight conditions are changed, and estimated the flight conditions to obtain the flight distance by solving the equation of motion [31] [32]. As a result, it was found that the flight conditions for obtaining the flight distance are the initial air pressure is about 5 atm. to 7 atm., the water volume ratio is about 30% to 40% of the bottle capacity (1.0 L, 1.5 L), and the aircraft weight is 100 g to 130 g, the launch angle is about 45˚ to 55˚.

In Chapter 5, based on the analysis results, we verified the optimum flight conditions for obtaining the flight distance and stability obtained in the simulation of a water rocket by applying RPD. Comparing the analysis results of flight characteristics by simulation and the evaluation results of optimum conditions by RPD, we obtained results with considerable correlation (Table 23). However, since various flight conditions can be set for water rockets, these conditions are not independent of each other but are closely related to each other. Therefore, it was clarified that it is necessary to set the flight conditions to achieve the best balance. Furthermore, from the results of analysis of variance in the evaluation of flight conditions by RPD, it was also verified that it is necessary to add another

Table 23. Summary of optimization analysis results.

parameter level for optimization to obtain maximum flight distance and flight stability. For example, parameters not covered in this study include the size and shape of the tail, its mounting position, and the length of the fuselage. If these are changed, in addition to the weight of the aircraft, the lift, drag, center of gravity, and aerodynamic center position will change, which can be presumed to affect flight stability [33].

As a future research subject, we plan to add parameter conditions related to flight and proceed with a more detailed analysis of the error between the analysis result by simulation and the demonstration experiment result. Furthermore, the effect of the behavior of water in the bottle of water rocket propellant on flight stability will be analyzed. In addition, we plan to acquire hydrodynamic empirical experiment data related to flight performance through wind tunnel experiments, etc., and further, explore research on optimization of flight characteristics [34] [35].

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Sakai, T. and Nakamura, K. (2000) Flight of Aqueduct for Wire Extension. Proceedings of the 2000 Annual Meeting of the Japan Society of Mechanical Engineers (IV), 531-532. https://doi.org/10.1299/jsmemecjo.2000.4.0_531
[2]	Tomita, N., et al. (2003) Development of Science and Technology Education System Using Water Rocket. Proceedings of the 2003 Annual Meeting of the Japan Society of Mechanical Engineers (V), 393-394. https://doi.org/10.1299/jsmemecjo.2003.5.0_393
[3]	Ota, T. and Umemura, A. (2001) Study on Optimal Flight Conditions for Water Rockets. Proceedings of the Japan Society for Aeronautics and Astronautics, 49, 382-387. https://doi.org/10.2322/jjsass.49.382
[4]	Watanabe, R., et al. (2000) Thrust Characteristics of Water Rocket and Its Improvement by Modificating of Nozzle. Proceedings of the 2000 Annual Meeting of the Japan Society of Mechanical Engineers (IV), 525-526. https://doi.org/10.1299/jsmemecjo.2000.4.0_525
[5]	Watanabe, R., et al. (2004) Aerodynamic Characteristics of Water Rocket and Stabilization of Flight Trajectory. Proceedings of the Japan Society for Aeronautics and Astronautics, 51, 449-455. https://doi.org/10.2322/jjsass.52.449
[6]	Kamiashi, F. (2018) Scientific and Technological Calculations with Excel. Maruzen Publishing, Tokyo, 25-158.
[7]	Odaka, T. (2018) Numerical Calculation and Simulation by Python. Ohmsha, Tokyo, 52-81.
[8]	Wakui, Y. and Wakui, S. (2014) Understanding Multivariate Analysis. Gijutsu Hyoronsha, Tokyo, 190-200.
[9]	Oguro, R., et al. (2019) Customer Number Forecasting by Using T-Method in Quality Engineering. Summary of National Research Presentation Conference.
[10]	Soga, M. (2008) A Comparison of Estimate Accuracy with T Method and Multiple Regression Analysis. Quality Engineering Research Presentation Conference Proceedings = QEF/Quality Engineering Research Presentation Conference Executive Committee Edition, 16, 430-433.
[11]	Kawada, H. and Nagata, Y. (2015) Studies on the Item Selection in Taguchi’s T-Method. Journal of the Japanese Society for Quality Control, 45, 179-193.
[12]	Kikuchi, T. (2010) A Research of Estimate Accuracy with One-Side T Method. Proceedings of the Quality Engineering Society, Quality Engineering Research Presentation Conference, 18, 34-37.
[13]	Inabu, J., et al. (2012) Prediction Accuracies of Improved Taguchi’s T Methods Compared to Those of Multiple Regression Analysis. Journal of the Japanese Society for Quality Control, 42, 265-277.
[14]	Tatebayashi, K., et al. (2012) Introduction MT System. Union of Japanese Scientists and Engineers, Tokyo, 133-174.
[15]	Hirotsu, C., et al. (1997) Algorithm for p Value, Power and Sample Size Determination for Max T Method. Society of Applied Statistics, 26, 1-16. https://doi.org/10.5023/jappstat.26.1
[16]	Taguchi, G. (2005) Purpose and Basic Functions (6) Comprehensive Prediction by T Method. Quality Engineering Society, 13, 5-10.
[17]	Taguchi, G. (2006) Purpose and Basic Functions (11) T Method for Recognition. Quality Engineering Society, 14, 5-9.
[18]	Maeda, M. (2017) New Regression Method Based on the Idea of T-Method (1). Journal of the Japanese Society for Quality Control, 47, 185-194.
[19]	Inao, A., et al. (2012) Taguchi’s T Method and Its Improved Method and Performance Comparison of Multiple Regression Analysis. Journal of the Japanese Society for Quality Control, 42, 265-277.
[20]	Hosokawa, T., et al. (2015) A Proposal of Development Methodology Integrating Parameter Design and T-Method. Journal of the Japanese Society for Quality Control, 45, 194-202.
[21]	Taguchi, G. (2002) Evaluation Technology for Optimized Design. Japanese Standards Association, Tokyo, 13-90.
[22]	Taguchi, G. (1988) Parameter Design in New Product Development. Japanese Standards Association, Tokyo, 1-74.
[23]	Taguchi, G. (1999) Quality Engineering for Technological Development. Japanese Standards Association, Tokyo, 1-24.
[24]	Suzuki, M. (2016) Not Difficult Quality Engineering. Nikkan Kogyo Shimbun, 34-48.
[25]	Yano, H. (2002) Introduction to Quality Engineering Calculation Method. Japanese Standards Association, Tokyo, 74-137.
[26]	Tsuruta, H. (2016) Energetic SN Ratio. Union of Japanese Scientists and Engineers, Tokyo, 17-72.
[27]	Taguchi, G. (1993) Design of Experiments and Quality Engineering. Journal of Quality Engineering, 2, 2-8.
[28]	Hirose, K. and Ueda, T. (2015) Introduction to Taguchi Method Analysis Method. Doyukan, 1-24.
[29]	Tatebayashi, K. (2004) Introductory Taguchi Method. Union of Japanese Scientists and Engineers, Tokyo, 9-37.
[30]	Koshimizu, S. and Suzuki, M. (2007) Practical Quality Engineering. Nikkan Kogyo Shimbun, 3-56.
[31]	Taguchi, G. (2012) Mahalanobis-Taguchi System in the 21st Century MTA TS and T Method. Journal of Quality Engineering Forum, 20, 261-268. http://id.ndl.go.jp/bib/000000400206
[32]	Taguchi, G. and Yano, H. (2004) Quality Engineering Application Course: Technology Development of Information Design by Computer-Simulation and MT System. Japanese Standards Association, Tokyo, 310-368.
[33]	Iwasaki, H., et al. (2000) Estimating Optimal Flight Conditions for Water Rockets (Part 1). Proceedings of the 2000 Annual Meeting of the Japan Society of Mechanical Engineers (IV), 517-518. https://doi.org/10.1299/jsmemecjo.2000.4.0_517
[34]	Nakano, R., et al. (2003) Observation of Water Rocket Orbit by Image Processing. Proceedings of the 2003 Annual Meeting of the Japan Society of Mechanical Engineers (V), 387-388. https://doi.org/10.1299/jsmemecjo.2003.5.0_387
[35]	Bairagi, A., et al. (2018) Negative Imaginary Approached High Performance Robust Resonant Controller Design for Single-Phase Islanded Microgrid and Its Voltage Observation on Different Load Condition. Intelligent Control and Automation, 9, 52-63. https://doi.org/10.4236/ica.2018.92004