技術レポート

A Porosity Formation and Flattening Model of an Impinging Molten Particle

H.Fukanuma

Abstract

Thermal spray coatings have porosity;however,reasons for the production of porosity during the coating process are not known.This paper proposes a physical and mathematical model for the production of porosity by considering deformation of a molten particle during thermal spray coating processes.The theoretical model shows that the impinging velocity,the ambient gas pressure,the particle diameter,and the molten material viscosity contribute to producing porosity.This paper also proposes that there is a porosity distribution along the splat radius and that most of the porosity exists in the periphery of the splat.Also,a flattening model proposed in this work agrees well with the results of Engel(Ref 1).

1. Introduction

THE boundary structure between splats and the substrate or previously coated layers greatly affects the physical characteristics of thermal spray coatings.In particular,porosity in the thermal spray coating is an important factor in determining the characteristics of the layers.At present,the structures of splat interfaces are being clarified by various studies(Ref 2-6),but there is not yet an established theory of the mechanism of porosity formation in the layers.In this paper,the porosity formation mechanism of thermal spray coatings is mathematically analyzed under several assumptions.

Although there have been some studies concerning the flattening process,they have not theoretically described the flattening process as a function of time(Ref 7,8).In the present work,a particle flattening model is proposed from a new physical viewpoint.The model includes the radius expansion as a function of time,so that the flattening time can be calculated.

2. Theory

2.1 Definition of Porosity

The structure of a thermal spray coating can be pictured as lamellar(Fig.1),where a first layer of splats is contact with the substrate.It is known that gaps exist between the substrate and the splats,and between the splats themselves.The schematic of an enlarged interface is shown in Fig.2.It is assumed that the interface consists of two parts,the gap and the contact region.The true contact region is defined as the interface where chemical bonding exists.It may be a very small proportion of the interface(Ref 3).The contact distance should be at the atomic level such that gas molecules or atoms do not pass through the interfaces.Arata et al.have shown that these gaps are continuous(Ref 5,6).

null

These studies suggest that there is little probability that closed pores exist in the interface.

The porosity ƒp of the thermal spray coating is defined as:

null

where Vg is the total gap volume in the intersplat and splat-substrate interface in the unit volume of thermal spray coating and Vs is the total solid volume in the unit volume of the thermal spray coating.The total gap volume,Vg,is defined as follows:

null

where z is a vector parallel to the z-axis with a length |z|=z(the z-axis is perpendicular to the substrate surface); dS is a vector perpendicular to the small area dS on the splat surface,positive to the outer direction of the splat and with length equal to the surface area ofdS;and S is the total surface area of all splats in the unit volume of the coating.Vg is the volume sum of every dV in the splat interface of the unit volume.Since DV = zdS , Vg can be obtained by the integration of zdS.However,when integration is carried out for the surface area of every splat,integration is performed twice(i.e.,both the lower and upper surfaces of each splat are integrated)and therefore the resultant value is divided by two to calculate the total gap volume.

null

The “z-distribution” of the interface can describe thermal spray coatings better than porosity.In other words,to know the “z-distribution” is to know the interface structure.Although it is considered to be a more important factor than porosity,it is difficult to measure and analyze theoretically.

It is reasonable to assume that Vs and Vg can be replaced by a splat volume and the gap volume under the splat,respectively,as shown in Fig.3,because the cylindrical surface area of the splat(i.e.,the surface area of the splat edges)is negligible compared to its flat surface area when it is a thin disk.When the diameter and the thickness of the splat are Dm andh,respectively,the two flat surface areas of the top and bottom disks are 2π(Dm/2)2 and the cylindrical surface area is πDmh.The ratio of πDmh to 2π(Dm/2)2 is expressed as:

null

where d0 is the diameter of the molten particle at the moment of impingement.Actually,because Eq 3 is a negligibly small value compared to Eq 1,when Dm / d0 is larger than 3,the cylindrical area is negligible compared to the area of the flat surface.Hereafter,the porosity in the case of a single splat is discussed instead of the unit volume of the thermal sprayed coating.

2.2 A Mathematical Model of Porosity Formation

2.2.1 Gas Compression in a Hole by a Rigid Body Impingement

Initially,consider a gas compressed by a moving rigid body.The rigid body plunges into a hole on a substrate at initial velocity v0 and compresses the gas,as shown in Fig.4.When the rigid body velocity is reduced to v,the gas pressure is expressed by the polytropic equation

null

where P is gas pressure,k is the gas constant,V is gas volume,and n is the polytropic exponent.It is assumed that the gas is ideal,there is no gas leakage through the boundary between the rigid body and the hole wall,and there is no friction between the body and the wall.In this case,n = 1 for isothermal compression,and n is γ = Cp / CV for a reversible adiabatic compression,where Cp and CV are the specific heat of gas at constant pressure and constant volume,respectively.When the gas volume is compressed to V from V0 and the velocity of the rigid body is reduced to v from 0,under the reversible adiabatic condition,the energy transfer to the gas from the rigid body equals the energy that the rigid body loses.So the next equation holds:

null

where m is the mass of the rigid body.

null

2.2.2 Gas Compression in a Hole by a Molten Disk

Next,consider the gas compression in the hole by a molten disk impinging on the substrate instead of the rigid body,as shown in Fig.5.If there is no friction in the molten material and at the interface between the molten material and the wall,and there is no flow parallel to the substrate in the liquid and no gas leakage through the interface,the molten material ABCD compresses the gas in the hole as in the case of the rigid body.Since the mass of ABCD can be replaced by k1ρ(V – V0),where k1 is a constant and ρ is the density of the molten material,the following equation holds:

null

There are other conditions that may affect the gas volume V and must be considered in the case of a molten material compressing the gas:

・Gas leakage through the boundary at the hole wall during the compression process ・Friction between the molten material and the hole wall ・Viscous work in the molten material ・Heat transfer from the liquid to the gas and from the gas to the substrate ・Flow of the molten material parallel to the substrate surface ・The contact angle between the molten material and the substrate or previously coated layers ・The vapor pressure of the molten material(unless this is negligible)

null

The parameter n can be used in place of γ and represents the difficulty of compressibility.The parameter n is influenced by the above factors and is larger than 1.The molten material stops in the hole when v = 0,and the next equation holds:

null

Integrating Eq 7 and substituting kV0-n = P0 yields the equation:

null

where P0 is the gas pressure at gas volume V0.

When V/V0 << 1,because 1 - V/V0 ≅ 1 and (V/V0) 1-n - 1 ≅ (V/V0) 1-n,Eq 8 can be simplified as: null

The ratio of the compressed volume V to the initial volume V0 decreases as the impinging velocity increases and the ambient pressure decreases;that is,the porosity decreases under these conditions.Figure 6 shows the relationship between the gas compression ratio and the impinging velocity in Eq 8 and 9 when ρ = 800 kgm/3 , P = 105Pa , k1=10 and n = 1.5 or 2.It also shows that V/V0/sub> in Eq 9 agrees with V/V0 in Eq 8.The expression k1=10 means that when the molten disk of 10μm thickness impinges on the substrate,the molten material enters into the hole to 1μm depth.

null

2.2.3 Solidification Time of Molten Material in a Hole

The material is pressed away from the developing gas interface if the molten material does not solidify and the gas pressure is high.If the molten material does solidify in the hole,then the solidification time of the material must be considered.Figure 7 is a schematic of the solidification of molten material.When the temperature of the substrate and the molten material are TS0 and T0,respectively,at time 0,and TS0 is held constant and the molten material thickness is semi-infinite,the solidified thickness X at time t is expressed as:

null

where λ is a numerical constant(detailed in Eq 11),and kS and t are the thermal diffusivity of the solidified material and time,respectively(Ref 9).The parameter λ is found from:

null

where KS and CS are the thermal conductivity and the specific heat of the solidified material,and KL and kL,Tm and L are the thermal conductivity.the thermal diffusivity,the melting temperature,and the latent heat of the solidified material,respectively.

null

null

Figure 8 shows the results calculated for a copper substrate sprayed with copper particles at 2000 and 2500K using the physical property values in Table 1.Now,consider a cylindrical hole with a diameter and depth of 1μm,and assume that the solidification thickness is half of the hole diameter.Because the solidification thickness is 0.5μm in a 1μm-diam hole,the solidification time of the molten copper in the hole is about 3 ×10-9s when the initial temperature of the molten copper is 2000

null

K,and about 5×10-9 at 2500K in Fig.7.When the solidification time is 5×10-9s,the molten copper enters the hole to a depth of 0.5μm or less at an impinging velocity of 100m/s,and to a depth of 1μm or less at an impinging velocity of 200m/s.The true solidification time is greater than the calculated results because the heat transfer coefficient at the interface between the substrate and the solidified material is far smaller than the assumption used for a thermal spray coating.

Consider a 10μm thick disk of molten copper compressing gas in the hole(1μm deep)when the impinging velocity is 200m/s.Substituting v0=200 m/s,ρ=8780 kg/m3,P0=1,atm=1.013×105Pa, k1 = 10 ,and n=2,V/V0 = 5.8×10-4 is obtained.Because V/V0= 5.8×10-5,the molten copper almost reaches to the hole bottom.It takes 5×10-9s for molten copper to reach the bottom at speeds of at least 200m/s.However,because the velocity of the particle is decelerated by the compressed gas in the hole,it takes longer than 5×10-9s for the molten material to reach the bottom.Because solidification begins from the hole wall during the compression process,when the top of the molten material reaches to the hole bottom most material in the hole is solidified.When the hole diameter and depth are of the order of 1μm or less,it is reasonable to assume that the solidification of the top of the molten material in the hole occurs during the compression process.The remaining highpressure gas leaks through the boundary,and the pressure equalizes at ambient pressure.

Only one simply shaped hole has been considered for porosity formation with respect to the trapped gas phenomenon.It is reasonable to assume that the hole volume V0 represents the sum of hole volumes on a unit surface area of the substrate or previously coated layers.The argument also holds when the hole shapes are not cylindrical.

2.3 A Mathematical Model of Flattening

Flattening models of liquid droplets that impinge on a rigid flat surface have been proposed(Ref 1,7,8),and this paper develops a prior investigation(Ref 10).Because molten particles impinge on a substrate,flatten,and solidify,the particle-gas interaction during flattenig must be considered to understand porosity formation.

2.3.1 Radial Expansion Shortly After Impingement

It is assumed that immediately after a molten particle impinges on a flat substrate,a thin radial disk begins to spread,and that the splat thickness is a constant during the flattening,as shown in Fig.9.The disk spreads under the kinetic energy that is mainly dissipated by viscous work.The disk thickness is probably a function of the impinging velocity,the particle diameter,and the viscosity of the molten material.The assumption agrees with the results of high-speed photographs of a water drop spreading in Engel’s investigation(Ref 1).However,in his investigation Engel assumed that the disk thickness at center is half of the peripheral thickness.The following are also assumed:

null

・The flattening process is isothermal(i.e.,the viscosity and density of the molten material are constant).
・The kinetic energy that the particle possesses at impingement is dissipated only by viscous friction in the molten material.
・The friction work in the disk only contributes to the energy dissipation,but the friction work in the part of the sphere on the disk can be negligibly small.
・The surface tension of molten material,the interface tension between the material and the substrate,and gravity can also be negligible,as they are generally far less than the kinetic energy experienced during thermal spraying.

Additional assumptions from the present work are:

・The impinging particle is completely molten
・The particles impinge perpendicularly on the substrate or the previously coated layers,which are smooth.
・Flattening is completed before solidification.
The particle does not rotate at impingement.
・There is no flow in the particle before impingement.

Consider a molten particle of diameter d0 that impinges at velocity v0.At time t=0,a bottom portion of the particle is deformed to a flat disk whose thickness and radius are h andR,respectively.When cylindrical coordinates are used,the volume of 0AB,now identified as ω,is expressed as:

null

where γ is the cord radius made by the sphere of AB length,and H is the height of the bottom portion of the sphere,as shown in Fig.9.The parameter γ can be expressed as:

null

Substituting Eq 13 in Eq 12,the expression becomes:

null

When the portion 0AB in the sphere is deformed to the flat disk whose radius and thickness are R and h,respectively,the volume 0AB equals the volume of the disk.Then:

null

When H<<d0 then:

null

Differentiation of Eq 15 by t leads to an expression for the expansion rate of the disk.Because H = v0t,the expansion rate can be expressed as:

null

When H/d0<<1,Eq 17 is simplified as:

null

null

Consequently,the radial expanding rate of the disk is obtained shortly after impingement.

2.3.2 The Flow Field in the Disk

Consider a jet impinging on a flat surface as shown in Fig.10,the flow field consists of four distinct zones.The maximum radial velocity Um along the γ axis is(Ref 11-13):

null

Flow zones III and IV are divided at r = R0.Um increases proportionally with r and becomes maximum at r = R0,then decreases in inverse proportion to r.

When the flow field in the disk of the impinged particle is defined as u = u(ur,uθ,uz),thenur,uθ,and uz are the r,θ,and z components of the flow field,respectively,in cylindrical coordinates.From the analogy of impinging jets,it is assumed that the radial velocity ur is proportional to the radius when r≤R0 and is inversely proportional to the radius when r > R0.The parameter R0 is a function of the impinging velocity,the particle diameter,and viscosity,and it is assumed that R0 is close to the radius of the impinging particle(d0/2).Therefore the flow field u = u(ur,uθ,uz) is:

null

where h is the disk thickness and C1I(t) and C2(t) are coefficients that are a function of time.

null

The flow field must satisfy the equation of continuity,namely div u = 0,and meet the following requirements:(a)no shear stress on the free surface;(b)no flow on the substrate surface;and(c)when r → ∞,ur → 0.The flow field u satisfies all of these conditions.The flow field is symmetrical around the z-axis(because no particle rotation at impingement is assumed),so uθ equals 0.

The flow field could be assumed to be laminar in the disk when the impinging velocity is less than 250m/s(Fig.11),because the Reynolds number is substantially below 5000.For example,Raynolds number 5000 occurs when density ρ,the mean radial flow velocity ur, representative length h,and viscosity μ are 104kg/m3,250m/s,2μm,and 10-3Pa・s,respectively.

Because it is expected that the flow field decreases rapidly with time,then:

null

where C10,C20,and α are constants.There is a volume balance of material about R0 such that the material that flows from the inside region of R0 is the same volume that flows into the outside region of R0.Therefore:

null

Equation 23 is simplified as:

null

The flow field can be determined when C10,C20,and α are known under certain boundary conditions.For example,when the disk radius becomes R0 and R at times t0 andt,respectively,the fluid volume outside R0 equals the volume that flows out of R0 for the time from t0 to t.Hence at the time difference τ = t – t0:

null

Equation 25 is simplified as:

null

When Eq 26 is differentiated and τ=0 is substituted,then the expression becomes:

null

When τ=0,t = t0,and R = R0,substituting t = t0 in Eq 17,the expression becomes:

null

When R0≤d0,v0t0d0 is negligible in Eq 28,then:

null

Then,C20 can be determined as follows:

null

Substituting Eq 30 in Eq 24,the expression becomes:

null

In Eq 26,letting R = Rm when τ=∞,the expression becomes:

null

Because C10,C20 and α are found above,the flow field is determined.

2.3.3 A Flattening Model

From an energy balance,we find:

null

where

null

From the flow field u,because the disk thickness h is thin(i.e.,h is negligibly small compared to disk radius R),∂uz∂r is negligibly small compared to ∂ur∂z and Φ is written as:

null

where μ is the viscosity of molten material and V is the entire volume of molten material.Substituting Eq 34 and then substituting the result in Eq 33,assuming that the surface and interfacial energy are negligible,the expression becomes:

null

Substituting Eq 30,Eq 31,and Eq 32,and h = d30/6R2m in Eq 35,the expression becomes:

null

null

Letting Dm = 2Rm,R0 = εd0/2,where ε is in the neighborhood of 1,and substituting them in Eq 36,the expression becomes:

null

When ε=1,Eq 37 is approximated as follows:

null

When Reynolds number(Re)is ρd0v0/μ,the expression becomes:

null

The disk diameter at t = ∞ is determined.

Madejski’s theory,Jones’theory,and this theory are shown in Fig.12,where the flattening ratios as a function of Reynolds number are shown.Madejski’s and Jones’theorise are respectively(Ref 8,9):

null

null

In actual flattening process,it is considered that R0/d0 must change with time and that therefore ε probably changes from 1 to 2.When the impinging velocity is faster,R0/d0 must be larger.However,Dm/d0 changes little with ε between 1 and 2,as shown in Fig.13,so ε can be considered a constant.

Substituting Eq 30 and Eq 32 in Eq 26,when R =D/2,R0=d0/2,and Rm = Dm/2,Eq 26 is simplified as:

null

The spreading disk diameter is determined as a function of time.When 1/(Dm/d0)2<<1,Eq 42 is simplified as: null

Substituting Eq 39 in Eq 43,the expression becomes:

null

null

From Eq 43,the spreading time is determined as:

null

When R = R0,t = t0.Substituting H = v0t,h = 2d30/3D2m,and R0 = d0/2 in Eq 15,because H/d0is negligible,t0 is determined as follows:

null

Because t = t0+τ,t is expressed as:

null

When the first term in Eq 47 is negligible,t becomes:

null

Substituting Eq 39 in Eq 48,the time for the disk to spread to D = 0.9Dm is:

null

and the time for the disk to spread to D = 0.99Dm is:

null

null

where t0.9 and t0.99 are the times when the disk diameter becomes 0.9Dm and0.99Dm,respectively.Equations 49 and 50 show that determining the real spreading time is difficult,because the spreading rate is slower in the neighborhood of the maximum diameter.Figure 14 shows that the theory agrees well with the experimental results of Engel(Ref 1)and shows Madejski’s and Jones’flattening ratios.

2.4 Porosity Formation Model with Flattening Process

The simple theory of porosity formation was described in section 2.2 for the case where a molten disk impinges on a flat substrate having small holes.This section describes a porosity formation model that considers the flattening process.When V(R) is the gap volume per unit area at radius R of the splat,the entire gap volume Vg between the splat and the substrate surface or previously sprayed layers is expressed as:

null

It is assumed that the impinging velocity v0 affects the gas between the splat and substrate surface in the region inside the radius R0 and that the radial spreading rate dR/dt = vRaffects the gas in the region outside R0,as shown in Fig.15.Because Eq 9 is the approximation of Eq 8,when v0→0,V →V0 in Eq 8.However,V →∞ in Eq 9.When R →Rm,vR→0,so that V(R)→∞.In order to avoid that,Eq 9 diverges at R = Rm,the equation V(R) is substituted for Eq 9,using a parameter ζ that is a numerical value close to and smaller than 1.So V(R) is:

null

null

Substituting Eq 52 in Eq 51,the expression becomes:

null

The radial spreading rate vR can be determined by substituting Eq 30 and Eq 32 after differentiating Eq 26 with respect to time as follows:

null

The gap volume VgR from R = 0 to R is defined as:

null

The accumulative distribution of gap volume VgR/Vg is determined by integrating Eq 53 and Eq 55,after substituting Eq 54 in Eq 53 and Eq 52 in Eq 55 because the values in 0≤ R ≤ R0of Eq 53 and Eq 55 are negligible compared with those in R0 ≤ R ≤ Rm.Also,when R20/R2mis negligible,the expressions become,respectively,

null

null

null

Figures 16 and 17 show VgR/Vg as a function of the normalized radius R/Rm,when n = 1.5and n = 2.They show that the gap volume increases with increasing radius and becomes extremely large in the periphery of the splat.Figure 16 shows that when n = 1.5,most gap volume exists outside 0.9Rm(i.e.,in the peripheral region of the splat).When n=2,the gap volume begins to increase from about 0.4Rm;it becomes about 10% up to approximately0.7Rm.It is shown that most gap volume exists in the periphery of the splat.

Porosity is defined by Eq 1.When same-size molten particles impinge at different velocities and all other conditions are the same,then the ratio at impinging velocity v0 to the porosity atv00 is expressed as:

null

where Vg and Vg0 are the gap volumes at impinging velocities v0 and v00,respectively;ƒpand ƒp0, are the porosities at v0 and v00,respectively.Vs is a splat volume;and Vg for0≤R<ζRm is simplified in Eq 59 and 60.It is assumed that R20/R2m and the first term on the right-hand side of Eq 53 and are negligible,therefore,integrating the second term,the expression becomes: null

null

null

The porosity ratio is shown as a function of impinging velocity in Fig.18 as ƒp0 = 0.1 at 100m/s and n = 1 and 2.Figure 19 shows the porosity ratio as a function of ambient pressure as ƒp0=0.1 at 1 atm.Figures 18 and 19 show that porosity decreases with increasing impinging velocity and increases with increasing ambient pressure.Equations 59 and 60 show that porosity decreases as the density or viscosity of the molten material increases.These results seem to be confirmed by empirical knowledge.

3. Conclusion

According to the theory developed in this article,the gap volume between the splat and the substrate or previously coated layers is minimal in the central portion and relatively larger in the periphery of the splat.It is believed that the interaction force between the splat and the substrate within the radius R0 is very different form that outside it,because the flow field uoutside R0 shows that pressure is zero on the substrate surface.The pressure at the impingement point is aρC0v0 at the moment of impingement(Ref 1),where a is a constant(about 0.5)and C0 is the velocity of a compression wave in the molten particle.This gradually becomes smaller as the flattening particle expands to R0.The bonding force insideR0 is probably stronger than that outside R0.Interaction between the molten particle and the substrate is strongly affected by the impinging velocity inside R0,but not that outsideR0.Although wettability of coating material onto the substrate or of coating material onto prior solidified material is not directly considered in this theory,it might influence the interface interaction very much outside R0.This would be particularly true in the periphery of the splat,when the molten material wets the substrate.The radial spreading rate is as fast as the impinging velocity in the neighborhood of R0,so it is questionable whether the rate at which molten material west the sub-

strate is faster than the spreading rate and whether the wetting action contributes to the interaction when the spreading rate is fast.However,because the spreading rate is slow in the outer part of the splat,wettability becomes important for porosity and bonding.It is believed that the velocity affects porosity and bonding when r < R0 and that the wettability is of prime concern when r > R0.

This theory illustrates that gaps must be continuous.It is believed that closed pores are not produced in the flattening process,because it is difficult for molten material to seal the hole completely against the high pressure.Gas probably escapes from the hole opening somewhere in the interface between the molten material and the substrate.Many paths for gas to escape through the interface are made during flattening.With the same reasoning,there are few true contact areas in the interface on the atomic level,and thus the bonding force is produced by an anchor effect or mechanical bond in thermal spray coating.If the molten particle temperature is sufficient to melt the substrate material,or if the substrate temperature is so high that it is molten,closed pores could conceivably be produced as molten particles impinge,and metal or chemical bonding may be achieved.

REFERENCES

1. O.G.Engel,Waterdrop Collisions With Solid Surfaces,J.Res.Nat.Bur.Stand.,Vol 54(No.5),May 1995,p281-298
2. R.McPherson,The Relationship Between the Mechanism of Formation,Microstructure,and Properties of Plasma Sprayed Coatings,Thin Solid Films,Vol 83,1981,p297-310
3. R.McPherson and B.V.Shafer,Interlamellar Contact Within Plasma Sprayed Coatings,Thin Solid Films,Vol 97,1982,p201-204
4. S.Safai and H.Herman,Microstructural Investigation of Plasma Sprayed Coatings,Thin Solid Films,Vol 45,1977,p295-307
5. Y.Arata,A.Ohmori,and C.Li,Electrochemical Method to Evaluate the Connected Porosity in Ceramic Coatings,Thin Solid Films,Vol 156,1988,p315-325
6. Y.Arata,A.Ohmori,and C.Li,Study on the Structure of Plasma Sprayed Ceramic Coating by Using Copper Electroplating,Proc.ATTC(Osaka,Japan),May 1988,p205-210
7. H.Jones,Cooling,Freezing and Substrate Impact of Droplets Formed by Rotary Atomization,J.Appl.Phys.,Vol 4,1971,p1657-1660
8. J.Madejski,Solidification Droplets on a Cold Surface,Int.J.Heat Transfer,Vol 19,1976,p1009-1013
9. H.S.Carlow and J.C.Jaeger,Conduction of Heat in Solids,2nd ed.,Oxford University Press,Oxford,England,1959,p285
10. H.Fukanuma,An Analysis of the Porosity Producing Mechanism,Thermal Spray:International Advances in Coatings Technology,C.C.Berndt,Ed.,ASM International,1992,p767-772
11. N.Rajaratnam,Chapter 11,Turbulent Jets,Elsevier Scientific Publishing Co.,1981(translated to Japanese)
12. M.poreh and J.E.Cermark,Flow Characteristics of a Circular Submerged Jet Impinging Normally on a Smooth Boundary,Proc.Sixth Mid-Western Conf.on Fluid Mechanics,1967,p198-212
13. S.Beltaos and N.Rajaratnam,Impinging Circular Turbulent Jet,Proc.A.S.C.E.J.Hydraulics Division,1974,p1313-1328
14. L.D.Landau and E.M.Lifshitz,Fluid Mechanics,2nd ed.,Pergamon Press,Oxford,England,1987,p50-51