Abstract
An analytical method based on the wave theory is proposed to calculate the pressure at the interfaces of coated plate subjected to underwater weak shock wave. The method is carried out to give analytical results by summing up the pressure increment, which can be calculated analytically, in time sequence. The results are in very good agreement with the finite element (FE) predictions for the coating case and Taylor’s results for the noncoating case, which validate the method that is suitable for underwater weak shock problem. On the other hand, Taylor’s results for the coating case are invalid, which indicates a potential application field for the method. The extension of the analytical method to q-layer systems and dissipation case is also outlined.
1. Introduction
Coated plates are widely used in the panels of naval vessels in order to make better protection for shock input such as underwater explosion (UNDEX). Therefore, systematic understanding of wave propagation mechanism in multilayered mediums is of great significance in guiding us to a better design for the coatings.
A lot of researches have been performed to study the response of the plane plates that interacted with underwater shock wave. Taylor [1] studied the response of air-backed plate that interacted with weak shock and obtained a one-dimensional analytical solution. Snay and Christian [2] improved Taylor’s solution by considering the strong shock that interacted with the air-backed plates and the results were obtained by method of characteristics. Schechter and Bort [3] considered plates subjected to weak shock with air-backed or water-backed boundary. Huang [4] studied the initial transient response of large elastic plates under spherical shock wave by using the Laplace and Hankel transformations. Rajendran [5] studied the reloading effects on plane plates subjected to noncontact UNDEX. However, the above-mentioned analytical method cannot be applied to coated plates. Chen et al. [6] built an analytical model to investigate the water blast response of one-dimensional marine structures coated with elastic foam. In their studies, the foam coating was modeled by a group of concentrated masses separated by parallel massless nonlinear spring and damper, and the first-order double asymptotic approximation (DAA) method was used to consider the fluid-structure interaction effects. But the results which were only the profile of the real results did not reveal characteristics of wave propagation.
Analytical solutions of transient wave propagation in multilayered mediums have a series of results in the literature. For one-dimensional problem of plane wave propagation in the direction normal to the layered medium, Sun et al. [7, 8] presented continuum theory instead of “effective modulus theory” for determining dispersion relation. Harmonic waves in composites with isotropic layers were studied by Stern et al. [9] and Hegemier and Nayfeh [10]. Transient plane waves propagating in a periodically layered elastic medium were examined by Ting and Mukunoki [11, 12] and Tang and Ting [13]. Chen et al. [14] developed an analytical solution based on Floquet’s theory to the problem of plate impact in layered heterogeneous material systems, and the analytical results and experiment data were in good accordance. Lin and Ma [15, 16] studied the dynamic response in a multilayered medium by the Laplace transform technique. The analytical solutions were presented in the transform domain and the numerical Laplace inversion was performed to obtain the transient response in time domain. The effective material concepts were investigated in their studies. Nevertheless, the methods presented in the above-mentioned paper are either limited to solve late-time solution [14] or given the boundary condition only a function of time, making them not suitable for obtaining the response of the coated plates interacted with underwater shock wave. The main difficulty to solve the problem is that the varying boundary conditions are not only a function of total propagation time in mediums, but also a function of impedance mismatch and the time of the stress wave travelling in each medium for one time since the multiple reflections in the mediums.
Because of the difficulties mentioned above, the problems were always solved numerically in time domain. Desceliers et al. [17] developed a numerical hybrid method to simulate the transient wave propagation in a multilayered semi-infinite medium which could be fluid or solid, subjected to given transient loads. Brasek [18, 19] obtained the dynamic response of one-dimensional aluminum structures that laid coatings under unit step pressure input and studied the response characteristics of metal structures under an exponentially decaying shock wave by FE analysis. Rajendran [20] simulated the response of plane plates subjected to uniform primary shock loading of noncontact UNDEX numerically. The simulation was in good agreement with the experiment for the elastic response. Kwon et al. [21] obtained dynamic response of a metallic cylinder, coated with a rubber material, subjected to an UNDEX, numerically. Gong and Lam [22] studied transient response of layered composite beams subjected to underwater shock by using coupled finite-element and DAA-boundary element codes.
From the above introduction, no analytical method can consider coated plates subjected to underwater shock wave. The primary objective of this paper is to propose an analytical method to calculate the pressure at the interfaces of multilayered mediums, which is suitable to solve the UNDEX problems. FE analysis is applied to validate the method.
The remainder of this paper is organized as follows. Section 2 presents the derivation of basic formulas used throughout this paper. Then, calculation examples and verification of the formulas derived before are treated in Section 3. Comparisons between analytical results and approximate predictions are presented in Section 4. The response of coated plates subjected to the underwater plane wave is solved by using the proposed method in Section 5. Finally, summary and conclusions are given in Section 6.
2. Solution Method
The underwater pressure wave will generate a series of reflection and transmission waves when the wave propagates in the multilayered mediums. As the number of reflection and transmission waves increases, the problem becomes very complicated. Rather than considering the specific reflection and transmission propagation process in time sequence, the method determines the pressure history by accumulating the pressure increments at the interfaces at a specific time, which can be calculated analytically, in time sequence.
2.1. Basic Assumptions
In this paper, wave propagation process in the water, coating, and steel plate 3-medium system is assumed to be a one-dimensional problem as shown in Figure 1 under the condition that the amplitudes of the pressure waves are small. As seen in Figure 1, medium 1 is the semi-infinite water, medium 2 is the coating with finite thick, and medium 3 is the steel plate. Boundary III at the right of medium 3 could be the varying boundary conditions, such as air-backed, water-backed, or rigid boundary. The main assumptions are listed as follows:(1)pressure wave propagates as longitudinal wave in the medium;(2)each medium is homogeneous with fixed linear impedance;(3)rigid body motion is negligible;(4)the varying of the wave propagation time in each medium due to the medium deformation is negligible.
2.2. Solution Procedure without Dissipation
The formulas derived below are all based on the simplified model shown in Figure 1. According to the acoustic theory, the one-dimensional wave function having small amplitudes in the ideal elastic homogeneous medium is presented as where , , , and are pressure, time, spatial dimension, and wave velocity, respectively. Before the reflection wave coming from boundary II reaches boundary I, the traveling wave solution in medium 1 can be obtained as where and are arbitrary functions, and the subscript refers to the different mediums. The transmission wave in medium 2 is given as Since the particle velocity satisfies the Euler Equation where is density of the medium, by substituting (2) and (3) into (4), one can obtain the particle velocity of medium 1 and medium 2, which is
This implies that for any function or , the pressure behaves the same with respect to the particle velocity . Therefore, function could be taken as any arbitrary function, such as harmonic wave, exponentially decaying wave, triangular wave, and step wave. Similar results could also be found in the book [23].
2.2.1. Pressure Increment at Boundary I at a Specific Time
Define that the propagation time of wave propagation in the coating for one time is and that in the steel plate it is and is the thickness of the medium: here one time means the process that the wave passes one boundary, comes to the other boundary, and then reflects back to the initial boundary as shown in Figure 2(a).
(a)
(b)
(c)
(d)
At a specific time, there must be a pressure wave transmitting from the coating into the water. Assuming that the time that the pressure wave has propagated, times in the coating and times in the steel plate, is the specific time , for example, the wave has arrived at boundary I for times from medium 2 and arrived at boundary II for times from medium 3, can be expressed as where and can be taken from 0 to infinity. Assuming that the pressure increment at boundary I at the specific time is and that the initial input pressure is , then the initial reflection wave can be obtained by Hereafter, the subscript 12 indicates the quantity from medium 1 to medium 2 and 21 indicates the quantity from medium 2 to medium 1, and so do 23 and 32. For example, is the transmission coefficient from medium 1 to medium 2 and is the reflection coefficient from medium 3 to medium 2.
The pressure increment at a specific time is related to how many times the wave has transmitted from medium 2 into medium 3. Therefore, the wave that transmits from medium 2 into medium 3 for different times must be calculated.
After the initial time (), if , which means that the wave does not transmit into medium 3, the pressure increment at boundary I denoted by can be written as The propagation process for is displayed in Figure 2(b), where the input wave transmits from medium 1 into medium 2 and propagates times only in medium 2 and finally transmits back into medium 1.
If , then the wave transmission times from medium 2 into medium 3 can be from 1 to . The number can be obtained by the equation .
The wave that transmits one time from medium 2 into medium 3 is presented as where the superscript 1 on indicates one time. The propagation process is shown in Figure 2(c), where the circle indicates the occurrence of the transmission wave which may be possible for any wave that arrives at boundary II from medium 2. Since the wave propagates times in medium 2 and will arrive at boundary II for times, the total times of possible propagation approach from medium 2 to medium 3 are . Therefore, the combinatorial number implies that one kind of propagation approaches is taken from kinds of propagation approaches that the wave transmits from medium 2 into medium 3.
The wave that transmits two times from medium 2 into medium 3 is presented as where the superscript 2 on indicates two times. The propagation process is shown in Figure 2(d), where the circles indicate the occurrence of the transmissions waves which may be possible for any two waves that arrive at boundary II from medium 2. Hence, the combinatorial number implies that two kinds of propagation approaches are taken from kinds of propagation approaches that the wave transmits from medium 2 into medium 3. Since the wave which propagates times in medium 3 is separated into two parts, the combinatorial number accounts for the different propagation approaches in medium 3.
More generally, the wave that transmits times from medium 2 into medium 3 is presented as where the superscript on , indicates times. The combinatorial number implies that the wave transmits times into medium 3 from times of propagation approaches in medium 2. The combinatorial number implies that the wave that propagates times in medium 3 is separated into parts.
To sum up, the pressure increment at boundary I at the specific time can be presented as where , and can be obtain by (12).
2.2.2. Pressure Increment at Boundary II at a Specific Time
The basic idea to calculate the pressure increment at boundary II is the same as that at boundary I. It becomes more complicated here, since the direction of the wave coming to boundary II for the last time can be different. That wave comes to boundary II can be from medium 2 or medium 3 so that the case should be separated into two parts.
The specific time of the pressure increment at boundary II is a little different from that at boundary I which is expressed as where can be taken from 1 to infinity and can be taken from 0 to infinity. Here and mean that the wave arrives at boundary II times from medium 1 and times from medium 2.
Assuming that is the wave that arrives at boundary II from medium 2 and that is the wave that arrives at boundary II from medium 3 at the specific time , the total wave that arrives at boundary II at the specific time , denoted by , is the summation of the wave coming from medium 2 and medium 3.
First, the wave that arrives at boundary II from medium 2 at the specific time will be calculated.
If , that is the wave does not propagate in medium 3, the wave is presented as
If and , the wave that transmits times from medium 2 into medium 3 is presented as where the number can be obtained through the equation .
Then, the wave that arrives at boundary II from medium 3 at the specific time can be calculated.
If and , or and , the wave that transmits 0 time from medium 3 into medium 2 is presented as
If and , the wave that transmits times from medium 3 into medium 2 is presented as where the number can be obtained through the equation .
In conclusion, the wave that arrives at boundary II from medium 2 at the specific time can be presented as where and can be obtained by (16); the wave that arrives at boundary II from medium 3 at the specific time can be described as where and can be obtained by (18). Then the total wave which is the pressure increment at boundary II at the specific time can be obtained by
2.2.3. Pressure Increment at Boundary III at a Specific Time
The procedure to calculate the pressure increment at boundary III is similar to that has been introduced to calculate the pressure increment at boundary I.
The specific time is defined by where and can be taken from 1 to infinity. Here and mean that the wave arrives at boundary II times from medium 2 and arrives at boundary III times from medium 3.
Assuming that is the wave that arrives at boundary III at the specific time , the wave that transmits times from medium 2 into medium 3 is presented as where the number can be obtained through the equation .
To sum up, the pressure increment at boundary III at the specific time is presented as where and can be obtained by (23).
2.2.4. Pressure at the Specific Time
The method to calculate the pressure increment at the specific time has been introduced above. After the specific time and the corresponding pressure increment , for , II, and III, for all different and are obtained, to calculate the pressure, sort the specific time in time sequence, and according to that, sort corresponding . At that time, all the pressure increments are in time sequence. Hence, the pressure at the specific time is the superposition of the pressure increment at the specific time and all the specific time ahead of .
If the input is the exponentially decaying wave, then . The exponentially decaying wave is dealt with in this paper.
Assuming that and are the ordered time series and the corresponding pressure increments for , the pressure at boundary , for , is given as
2.3. Derivation of the Formulas with Dissipation
In this section, the solutions with dissipation are derived.
Considering the damping term, the wave equation is modified as where , , and are the particle displacement, the damping coefficient, and the bulk modulus of the medium.
Assuming that the equation has the harmonic wave solution with the form where and are the imaginary unit and the angular frequency, by substituting (27) into (26), one finds that where Since is a complex number, can be written as where is the acoustic absorption coefficient. The travelling wave solution is presented as Substituting (31) into (32), one finds that Define in (30). By squaring (29) at both sides and substituting using (31), after setting the real parts and the imaginary parts are equal respectively, the wave velocity and the acoustic absorption coefficient are obtained as Here the negative roots of or are eliminated, by considering the positive -direction and the elimination of divergent results. Equations (34) and (35) imply that the wave velocity and the acoustic absorption coefficient are frequency dependent.
It can be proved that the wave equation (26) has the same form for , , and . Therefore, similar results can be obtained for the acoustic pressure and the particle velocity .
The difference between the solution of the wave equation with and without dissipation lies that the decaying term which is frequency dependent is added and the wave velocity turns to be frequency dependent in the solution with dissipation. It implies that despite the varying wave shapes for an arbitrary wave input, the wave shape keeps unchanged for a harmonic input except for the decay of the wave crest. In order to obtain the solutions of the wave propagation in the multilayered mediums, the relation between the acoustic pressure and the particle velocity must be given. Hence, the kinematic equation of the particles is presented as Then, integrating equation (36) as and substituting the solution of , the pressure amplitude associated with the velocity amplitude is obtained as Consequently, the specific acoustic impedance is obtained as
Therefore, the solutions of the wave propagation in the multilayered mediums with dissipation under harmonic input are easily obtained if the wave velocity, the decay of wave crest, and the specific acoustic impedance are modified by (34), (35), and (39), respectively. In order to obtain the solutions of an arbitrary input, Fourier transform is applied. In this paper, the exponentially decaying wave input is under consideration.
For the exponentially decaying wave, the cosine Fourier expansion is presented as where is the time length of the exponentially decaying wave and is the Fourier coefficient which can be given by where .
In the complex exponential form, the Fourier expansion of the exponential decay wave becomes where the Fourier coefficient is obtained as where .
Then, the analytical solutions of the wave propagation in the mediums with dissipation will be derived in the way similar to those without dissipation, except that the solutions with dissipation are the superposition of the Fourier components.
2.4. Extension of the Analytical Method to the -Layer System
The extension of the analytical method to the -layer system is still applicable, though it is very complicated. In this subsection, the main idea to solve the problem is presented, rather than giving the complex form of the formulas. The key point to solve the problem is to calculate the pressure increment at the specific time as it is illustrated before. The specific time of the -layer system is , where and are the times that the pressure wave propagates in the th layer and the time that the pressure wave propagates in the th layer for one time. At the specific time, the pressure increment is the superposition of many wave fragments. For example, in the 3-layer system discussed before, is a wave fragment that transmits one time from medium 2 into medium 3. In the -layer system, the wave fragment that considers the times of the wave transmitting in each interface should be calculated. Then the pressure increment is obtained by accumulating all the wave fragments. Once the pressure increments are obtained, the following procedures are the same as the procedures discussed above in the 3-layer system. The analytical method to the -layer system is obtained.
3. Calculation and Verification
In this section, calculation examples is presented, and comparison between analytical results and FE predictions is conducted to validate the analytical method. The analytical results are achieved by using MATLAB, and FE analysis is performed with FE analysis software ABAQUAS. Medium properties are presented in Table 1.
The model is subjected to an exponentially decaying wave input with a decay time of 0.6 ms and a peak pressure of 1 Pa. (Here 1 Pa is used for simplicity. However, it is enough for weak underwater shock wave. The reasons will be given in Section 5.) The wave propagates through the water, the coating, and the steel plate. Boundary III could be air-backed or water-backed or rigid. FE model is comprised of 4-node bilinear plane strain quadrilateral elements for solid and 4-node 2-D linear acoustic quadrilateral elements for fluid. The major axis is in the -direction. The geometric model established in ABAQUS is shown in Figure 3. The coating and the steel plate consists of 740 bilinear plane strain quadrilateral elements. The water part and the air part are 1m in length and discretized with 10,000 elements, respectively.
3.1. Analytical Calculation Examples
Analytical calculation examples are presented in this subsection to give a general idea of analytical results when the plate is coated.
Using the medium properties presented in Table 1, the analytical results with varying boundaries at different interfaces are obtained in Figure 4. As seen in the figure, pressure jumps appear at the specific time, which indicates that a wave comes to one interface. The peak pressures appear at the very beginning of the wave propagation. As Figure 4(b) shows, the pressure at boundary III is always 0, since it is the free boundary; as the pressure at water-coating surface drops below 0 where local cavitation is formed, analytical method is no longer valid.
(a)
(b)
(c)
3.2. Comparisons between the Analytical Results and FE Predictions
To validate the analytical method, comparisons between the analytical results and FE predictions with the varying boundaries at the different interfaces are conducted. The medium parameters are taken from Table 1. The results are presented in Figures 5–7. To show the reliability of the method, thinner thickness of the steel plates and larger ratio of the acoustic impedance conditions are also investigated. The results are given in Figures 8 and 9.
(a)
(b)
(c)
(a)
(b)
(c)
(a)
(b)
(c)
As seen in Figures 5–7, the analytical results agree considerably well with FE predictions at all three interfaces with all boundaries except for Figure 6(c). For the bad predictions of Figure 6(c), the explanations are as follows: though the nodes at boundary III are selected to see the pressure history, the numerical algorithm in ABAQUS gives the predictions of the boundary element that has a thickness, which cannot present the exact boundary. It is the thickness that triggers the large differences. It can be predicted that as the number of elements increases, the gaps in Figure 6(c) become smaller and finally close up. Then the predictions become coincide with the analytical results.
Figure 8 investigates the comparison results for 1 mm thick steel plate when the other property keeps the same. The FE model is modified as follows: the thickness of the steel plate is modified to 1 mm. The mesh density of the steel plate is 3 times of original one, and there are 30 elements now. The numbers of the coating elements and the water elements are 1,000 and 20,000, respectively, which are twice of the original ones. The FE predictions are the profile of the analytical results, since the numerical diffusions exist. To produce better predictions, more elements are given. However, in this case, numerical diffusions are still large so that the small pressure jumps cannot be distinguished. The FE algorithm for the thinner steel plate accounts for the inconsistent FE predictions that have large numerical diffusions. Hence, we can still say that the analytical results are validated in this case.
Figure 9 displays the comparison results where the wave velocity is 570 m/s. The acoustic impedance which has decreased to a level that is almost 1/3 of original one is less than that of the water. In this case, wave velocity is modified to 570 m/s, while the other medium properties keep unchanged as listed in Table 1. The analytical results are also consistent with FE predictions
To sum up, one can conclude that the comparisons validate the analytical method, since the analytical results agree reasonably well with the FE predictions for all boundaries and interfaces for different medium properties. FE analysis can hardly give perfect predictions under some extreme conditions, while analytical method has no such problem.
4. Approximate Method
The analytical method takes multiple reflections and transmissions into account. However, someone may ask whether the approximate method that only considers several reflections and transmissions at the beginning is enough to obtain reasonable results. In this section, comparison between the analytical results and the approximate predictions at the water-coating surface is presented.
Comparisons at the water-coating surface with rigid, air-backed, and water-backed boundary are investigated in Figure 10. The medium properties are taken from Table 1. The comparisons indicate the following.(1)The approximate predictions that only consider the reflection at boundary I underestimate the results in comparison to the analytical results.(2)The approximate predictions that consider the reflection at boundary II only one time agree well with the analytical results with rigid boundary, but the approximate predictions cannot reveal the oscillation effect of the pressures caused by multiple reflections and transmissions. The approximate predictions that consider the reflection at boundary II only one time are in conflict with the analytical results with air-backed boundary and water-backed boundary. The decaying rate of the case that considers multiple reflections and transmissions is much larger than the case that considers the reflection at boundary II only one time.
(a)
(b)
(c)
The specific examples show that the approximate method that only considers several reflections and transmissions with neglecting the effect of subsequent ones is incorrect. To obtain the correct result, multiple reflections and transmissions should be considered.
5. Application
The analytical method is applied to solve the underwater weak shock problem as Taylor [1] solved it before. The difference of this work is that the coated plate can be dealt with.
To apply the method to underwater weak shock problem, the restrictions of the method must be checked first. The shock wave applied to the coating and the steel plate must be weak enough to guarantee that the medium property is still linear, and the deformation of the mediums is negligible. Since the time span of the incident exponentially decaying wave is always short, the rigid motion can also be ignored.
In this section, the comparisons between Taylor’s results and the analytical results are presented. Two cases will be investigated here. Firstly, the noncoating case where the material of medium 2 and medium 3 that is steel for both of them is conducted. Secondly, the coating case, where the material of medium 2 and medium 3 used is listed in Table 1, is carried out.
Since the analytical results and Taylor’s results both behave linearly with respect to the incident peak pressure, the results can be normalized by dividing the incident peak pressure. Therefore, if the analytical results agree well with Taylor’s results for the noncoating case, one can say the method can be applied to weak underwater shock problem, which is independent of the peak pressure employed. This explains the choice of 1 Pa in Section 3 is without loss of generality to validate the analytical method. In this section, we use the normalized pressure to show the comparison results.
Figure 11 shows the comparisons between the analytical results and Taylor’s results at the water-coating surface where the material of medium 2 and medium 3 is steel. It shows that the analytical results presented in this paper are consistent with Taylor’s results. The oscillations of the analytical results indicate the reflections and transmissions of the pressure wave. The analytical results are more accurate than Taylor’s results with air-backed boundary.
Though Taylor’s results are no longer valid for the coating case, revised results based on Taylor’s method are still obtained by only modifying the mass when the plate is coated. Using the material of medium 2 and medium 3 listed in Table 1, comparisons between analytical results and revised results at the water-coating surface are shown in Figure 12. Though the revised results seem in accordance with the analytical results after a certain period of time, they are in conflict with each other at the beginning. The differences between the analytical results and Taylor’s results may be larger when the ratio of acoustic impedance is larger. It indicates the invalidation of Taylor’s method for the coating case and the potential application field of the analytical method.
6. Summary and Conclusion
An analytical method to calculate the pressure at the interfaces of coated plate subjected to underwater shock wave is obtained in this paper based on the wave theory. The analytical results are in good agreement with the FE predictions for the coating case and Taylor’s results for the noncoating case, which indicates that the method can be used for solving underwater weak shock problem. The failure predictions of Taylor’s method for the coating case show a potential application field of the analytical method. The approximate predictions that only consider several reflections and transmissions are also compared with analytical results. The comparisons show that the results cannot be obtained by simple approximation and multiple reflections and transmissions should be taken into account.
The analytical method is not restricted to two layers and could be executed for -layer systems and dissipation case which have been outlined. The primary drawback to this method is that the shock wave should be weak enough to guarantee the medium property is still linear, and the deformation of the mediums should be negligible. Furthermore, the analytical method solves the response of the coated plate before the formation of the local cavitation and the large rigid motion.
Nomenclature
| : | Wave velocity (m s−1) |
| : | Bulk modulus of medium (N m−2) |
| : | Pressure increment (Pa) |
| : | Arbitrary functions (Pa) |
| : | Thickness of medium (m) |
| : | Impulse (Pa s) |
| : | Pressure (Pa) |
| : | Reflection coefficient |
| : | Specific acoustic impedance (kg s m−2) |
| : | Transmission coefficient |
| : | Time (s) |
| : | Time length (s) |
| : | Energy flow per unit area (J m−2) |
| : | Acoustic absorption coefficient (m−1) |
| : | Damping of medium (N s m−4) |
| : | Decay time (s) |
| : | Density of medium (kg m−3) |
| : | Particle displacement (m). |
| : | Quantity related to the times of transmission. |
| : | Amplitude |
| : | Normalized |
| : | Input |
| : | Specific |
| : | Quantity related to damping |
| I, II, and II: | Quantity related to boundary |
| 1, 2, and 3: | Quantity related to medium |
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
The authors gratefully acknowledge the support for this work by the National Natural Science Foundation of China (NSFC) under Grant no. 11172173 and no. 11272215.