• 沒有找到結果。

含渦流作用之重力水波研究(1/2)

N/A
N/A
Protected

Academic year: 2021

Share "含渦流作用之重力水波研究(1/2)"

Copied!
11
0
0

加載中.... (立即查看全文)

全文

(1)

行政院國家科學委員會專題研究計畫 期中進度報告

含渦流作用之重力水波研究(1/2)

計畫類別: 個別型計畫 計畫編號: NSC92-2611-E-002-030- 執行期間: 92 年 08 月 01 日至 93 年 07 月 31 日 執行單位: 國立臺灣大學水工試驗所 計畫主持人: 黃良雄 報告類型: 精簡報告 處理方式: 本計畫可公開查詢

中 華 民 國 93 年 5 月 31 日

(2)

•

•

•Å

Å

Åo

o

»

»



I

I.

I

.

.

õ

õ

õº

º

ºÝ

Ý

ÝÞ

Þ

Þ@

@

@~

~

Œ

Œi

i

i









—

—

— 

 

 ×

×

×

â

â

â:

:

ø

ø®

®

®à

à

à



¥

¥æ

æ

æi

i

®

®@

@

@~

~

~

(1/2)

The Study of Gravitational Water Wave

with Vortex Effect (1/2)

Œi_r

: NSC 92-2611-E-002-030

ƕ§

: 92

O

8

`

1

93

O

7

`

31

^

x

x

¹

¹ß

ß

ß

:

?

?

ˆ

ˆ.

.

.>

>

>0

0

0

»

»

Ȗ

ñ

ñ



È

È



.

.

þa

þ

a

a





.

.

.





Abstract

In this study, a grid-free numerical model for the simulation of gravitational water wave with vortex effect is proposed. The computational algorithm is established by combining a potential flow field with a vortex flow field. The potential at the free surface and the free-surface displacement are solved by the time integrations. The vorticity field is the solution of the vorticity equation, which is solved using the vortex method and the vortex sheet method. The potential and velocity fields are solved by the boundary integral method.

Up to the present we have finished the establishment of mathematical formulations, the numerical models about potential flow and linear wave, and the numerical tests of vorticity generation from solid boundary. The unfinished parts are the numerical tests of vorticity generation from free surface, model integration, verifications and applica-tions.

Keywords: gravitational water wave, vortex, boundary integral equation method,

vor-tex method, numerical simulation

1

Introduction

Because of the demand for security of constructions and the complex phenomenon in fluid dynamics, the research in wave interaction with structures, such as bridge piers and submerged dikes, is an important topic for numerical simulation.

Water wave passing a structure causes wave reflection and diffraction as well as the generation and shedding of vortices. Vortex shedding will lead to severe form drag force, and the wave height will be influenced by the energy dissipation and the vorticity near free surface. Owing to the complication of wave-vortex interaction and the consideration for the computational efficiency, it is difficult to simulate the water wave with vortex effect.

In this study, a grid-free numerical model using boundary integral equation method for the simulation of gravitational water wave with vortex effect is proposed. The computational

(3)

algorithm is established by combining a potential flow field with a vortex flow field. In the first year, we develop the vortex separating mechanisms from rigid and free surface boundaries, and set up a two-dimensional vortex model. In the second year, we will improve the boundary conditions; deduce the boundary integral equations and calculating procedures to combine the two-dimensional potential model with vortex model. This hybrid model will be applied to the simulations of viscous water wave passing over a submerged dike, and the propagating of wave groups.

2

Mathematical Formulations

2.1

General Formulations of Two-Dimensional Linear Wave

Prob-lems

In the two-dimensional fluid domain V (t) with boundary S(t), the linearized Navier-Stokes equations for incompressible fluid with constant density ρ and kinematic viscosity ν are expressible in terms of the velocity u and the pressure p as

∂u

∂t = −

1

ρ∇p − gˆj + ν∇

2u in V , (1)

where ˆj is the unit vector which points vertically upwards, and x = (x, y), u = (u, v). Together with the continuity equation

∇ · u = 0 in V , (2)

the Navier-Stokes equations form a set of three scalar differential equations (for two-dimensional problems) sufficient for the determination of u and p, provided that adequate initial and boundary conditions are known.

It is, however, advantageous to introduce the concept of the vorticity ω = (0, 0, ω) defined by

ω = ∇ × u , (3)

and to consider the linearized vorticity equation

∂ω

∂t = ν∇

2ω (4)

obtained from equation (1) by taking the curl of each term in that equation and using equation (3). The set of equations (2) to (4), with u and ω as dependent variables, replaces the set of equations (1) and (2) in which u and p are dependent variables. There are several reasons for using ω in the formulation of the problem. The principal one is that the set of equations in terms of u and ω separate conveniently into two aspects: a dynamic aspect which deals with the change of the vorticity field ω with time described by equation (4), and a kinematic aspect which deals with the continuity of velocity field u governed by (2) and relates u at any instant to the vorticity field ω at that instant.

(4)

1. Kinematic free surface boundary condition:

∂ζ

∂t = v ; (5)

2. Dynamic free surface boundary conditions, normal stress:

−p + 2µ∂v

∂y = −p

s; (6)

3. Dynamic free surface boundary conditions, tangential stress:

∂u ∂y + ∂v ∂x = τs µ , (7)

where ζ(x, y, t) is the vertical displacement of the free surface from its unperturbed config-uration consisting of the xy-plane, ps is the atmospheric pressure, τs is the wind stress, and

µ is the dynamic viscosity. Because the vorticity ω = ∂v/∂x − ∂u/∂y, (7) can be rearranged

to give the vorticity on the surface,

ω = 2∂v

∂x

τs

µ . (8)

2.2

The Application of the Helmholtz Decomposition Theorem

The Helmholtz decomposition theorem states that, an arbitrary differentiable velocity field u can be decomposed into its irrotational part ue and rotational part uv as

u = ue+ uv, (9)

and there exists a scalar function φ, called the scalar potential, and a vector function, Bv = (0, 0, Ψ) (for two-dimensional problems), called the vector potential, which satisfies

ue = ∇φ , (10)

and

uv = ∇ × Bv. (11)

From the continuity equation in the incompressible flow, the scalar potential φ satisfies the Laplace equation

2φ = 0 . (12)

From the definition of vorticity, we get that Ψ satisfies the Poisson’s equation

2Ψ = −ω . (13) Substituting the decomposition (9) to (11) into the linearized Navier-Stokes equations, we get

(5)

p ρ = − ∂φ ∂t − gy , (14) and ∂Ψ ∂t = ν∇ 2Ψ . (15)

Using (13), (15) can be rewritten as

∂Ψ

∂t = −νω . (16)

Substituting for the surface vorticity from the shear stress boundary condition (8), at y = 0 the above equation becomes

∂Ψ ∂t = −2ν ∂v ∂x + τs ρ at y = 0 . (17)

Also, we can substitute (5) into (17) and have

∂Ψ ∂t = −2ν 2ζ ∂t∂x + τs ρ at y = 0 . (18)

These are equivalent to the stress boundary condition, and can be regarded as the boundary condition of Ψ at the free surface.

Applying the decomposition (9) to (11), the KSFBC (5) can be expressed as

∂ζ ∂t = ∂φ ∂y ∂Ψ ∂x at y = 0 . (19)

and the normal stress boundary condition becomes

∂φ ∂t + gζ − 2ν 2φ ∂x2 = − ps ρ + 2ν 2Ψ ∂x∂y at y = 0 . (20)

The results are summarized below. The problem requires to solve

2φ = 0 (21)

and the linearized vorticity equation

∂ω

∂t = ν∇

2ω (22)

in the flow field, and at y = 0 the associated boundary conditions are

∂Ψ ∂t = −2ν 2ζ ∂t∂x+ τs ρ , (23) ∂ζ ∂t = ∂φ ∂y ∂Ψ ∂x , (24) ∂φ ∂t + gζ − 2ν 2φ ∂x2 = − ps ρ + 2ν 2Ψ ∂x∂y. (25)

(6)

Now we focus on the special case where the wind stress τs is zero. From (8) the vorticity

at the free surface is

ω = 2 2φ

∂x∂y − 2

2Ψ

∂x2 , (26)

and (23) can be integrated with respect to time Ψ = −2ν∂ζ

∂x + Ψ0 at y = 0 . (27)

If the fluid was initially at rest, the constant of integration Ψ0 can be taken to be zero;

therefore

Ψ = −2ν∂ζ

∂x at y = 0 . (28)

When these are substituted into (24), it becomes

∂ζ ∂t = ∂φ ∂y + 2ν 2ζ ∂x2 at y = 0 . (29)

2.3

Integral Representations

2.3.1 Irrotational Flow Field

The governing equation of irrotational velocity field ue is the Laplace equation (12).

Applying the Green’s second identity, the integral representation of Laplace equation (12) can be obtained as

αφ(x) =

Z

S

00G − G∇0φ0) n0 dS0, (30)

where the subscript “0” indicates a variable, or a differentiation, or an integration evaluated in the x0 space. α is the solid angle at x, n is the unit outward vector normal to S, φ and

G are finite and continuous and possess continuous first and second partial derivative in V . G is the fundamental solution of Laplace equation,

G =

(

1

ln1r for two - dimensional problems

1

4πr for three - dimensional problems

, (31)

where r is the magnitude of the vector (x − x0).

2.3.2 Rotational Flow Field

The second Green’s identity for vector states that, if vectors P and Q are single-valued, finite, and have continuous second derivatives, then

Z V [P · (∇ × ∇ × Q) − Q · (∇ × ∇ × P)] dV = Z S n · (Q × ∇ × P − P × ∇ × Q) dS , (32)

(7)

where S is the boundary of the fluid region V , n is the unit normal vector on S directed outward from V . If we let

Q = Bv and P = ∇ ׳a

r

´

,

where a is a constant unit vector, we can get the integral presentation of the rotational part of velocity vector from equation (32) as (Wu and Thompson, 1973)

uv = Z V ω0× ∇0G dV0 Z S [(uv0· n0)∇0G − (uv0× n0) × ∇0G] dS0, (33)

here the subscript “0” indicates that a variable or an integration is in the r0 space, G is the

fundamental solution of Laplace equation. The gradient of G is

0G = −

(x − x0)

2(d − 1)πrd, (34)

where d = 3 for three-dimensional problems and d = 2 for two-dimensional problems. Thus, the general expression for the rotational part of velocity is of the form

uv(x) = 1 2(d − 1)π ·Z V ω0× (x − x0) |x − x0|d dV0 Z S (uv0· n0)(x − x0) |x − x0|d dS0 + Z S (uv0× n0) × (x − x0) |x − x0|d dS0 ¸ . (35)

The category of flows which were studied in Wu and Thompson (1973) is fluid in infinite domain bounded internally by a moving solid surface. In other words, besides the free stream velocity, the flow field is only driven by vorticity, and the boundary condition is simply the no-slip condition on the solid surface. In the viscous wave problem, there exists an important potential effect – pressure gradient due to gravity waves, and the boundary conditions are more complicated. Owing to this, we have to formulate the relationship between the integral representation of potential flow and that of rotational flow.

In the present research, we find that the irrotational part of velocity can be expressed as ue(x) = − 1 4π∇ × Z S µ(x0)∇0 µ 1 r× n0dS0, (36)

where µ(x0) is the strength of dipole located at x0. In other words,

ue(x) ≡ ∇ × Be, (37) where Be = − 1 Z S µ(x0)∇0 µ 1 r× n0dS0. (38)

The above corollary indicates the result: the irrotational velocity field can also be expressed as a solenoidal field

(8)

i.e., the velocity field can be expressed as

u = ∇ × (Be+ Bv) = ∇ × B

without loss of generality. Therefore, we can let Q = B instead of Q = Bv in (32), and get

the integral representation of velocity field u as u = Z V ω0× ∇0G dV0 Z S [(u0· n0)∇0G − (u0× n0) × ∇0G] dS0. (39)

This indicates that if the vorticity field in the domain and the velocity on the boundaries are known, the velocity field in the domain can be evaluated from the integral equation (39).

3

Numerical Methods

The potential φ at the free surface and the free-surface displacement ζ are solved by the time integrations of (25) and (29) respectively (ignoring the effects of wind stress and atmospheric pressure). The vorticity field ω is the solution of the vorticity equation (4), which is solved using the vortex method and the vortex sheet method. The potential field

φ and the velocity field u are the solution of the integral representations (30) and (39)

respectively, solved by the boundary integral method.

The flow chart of computational scheme for two-dimensional linear wave problems is shown in Fig. 1, and the brief introductions of numerical methods for vortex field simulation are shown in the following subsections.

3.1

The Vortex Method

The vortex method is a particle method in which fluid particles carrying concentrations of vorticity are followed as their positions and concentrations evolve with the motion of the fluid. Let

ω(x, t) =X

j

fσ(x − xj(t))Γj (40)

be an approximation to the vorticity at some time t = k∆t. The jth term on the right hand side of (40) is referred to as the jth vortex element or jth vortex blob, xj(t) is its position at

time t and Γj its strength. fσ is the vorticity distribution function. Here we use the Gaussian

distribution function fσ(x − xj) = 1 πσ2 j exp µ −|x − xj| 2 σ2 j, (41)

where the quantity σj is clearly a measure of the core size of the vortex element.

3.2

The Vortex Sheet Method

Near solid and free-surface boundaries, vortex blobs do not provide good approximations of the vorticity because the flow has a large velocity gradient in the normal direction that is not well-captured by the radially symmetric blob representation. A better approximation of

(9)

the vorticity in these regions consists of a collection of smoothed one-dimensional segments parallel to the boundary called vortex sheets.

In the vortex sheet method the vorticity at time t = k∆t is approximated by a sum of linear concentrations of vorticity

ω(x, t) =X

j

ωjbl(x − xj(t))δ(yj(t) − y) , (42)

where (x, y) denote coordinates which are parallel and perpendicular to the boundary re-spectively. Each term of the sum in (42) is referred to as a vortex sheet. The jth sheet has center (xj, yj) and strength ωj. Here δ is the Dirac delta function, and bl = b(x/l) is what

we refer to as the smoothing or cutoff function in analogy with the vortex method, and l is the length of sheet. The most commonly used cutoff is the hat or tent function originally proposed by Chorin (1978)

b(x) =

(

1 − |x| |x| ≤ 1,

0 otherwise. (43)

3.3

Vorticity Creation and Shedding

Besides the diffusion of vorticity in the flow field, another important viscous effect is the vorticity creation at the boundary. At the free surface, amount of vorticity creation is determined by (26). At solid boundary, the creation of vorticity is to maintain the no-slip condition at the surface. While the vorticity field in the domain V (t) (except the boundary) is solved from the step 2 shown in Fig. 1, and the vorticity distribution at the free surface is determined from (26), the vorticity distribution at solid boundary can be solved from (39). Fig. 2 is the velocity field of vortex sheets, which are just shedding from the flat plate. Fig. 3 is the numerical result of a uniform flow passing a rectangular dike. Compared to the experimental result (Van Dyke, 1982), the flow pattern seems to be reasonable.

References

Batchelor, G. K. (1967), An Introduction to Fluid Dynamics, Cambridge University Press, Cambridge.

Chorin, A. J. (1978), “Vortex Sheet Approximation of Boundary Layers”, J. Comp. Phys., 27, 428-442.

Van Dyke, M. (1982), An Album of Fluid Motion, The Parabolic Press, California.

Wu, J. C. and Thompson, J. F. (1973), “Numerical Solutions of Time-Dependent Incom-pressible Navier-Stokes Equations Using An Integral-Differential Formulation”, Computers

(10)

Input data and data preparation.

?

1. Time Integrations of φ,t and ζ,t.

Solve φ and ζ on the free surface.

?

2. Vorticity Transportation. Perform Euler convection and viscous diffusion of core radius

of all discrete vortices.

?

3. Vorticity Creation and Shedding. Evaluate the vorticity on boundaries. Shed the discrete vortices into flow field.

?

4. Potential Flow Computation. Solve the Laplace equations:

2φ = 0 and ∇2φ

,t = 0.

?

Advance time by ∆t.

¾

(11)

x (m) y (m ) -1 -0.5 0 0.5 1 0 0.001 0.002 0.003 0.004

Figure 2: Vortex sheets created from the solid boundary

U = 1.0, B = 3.0, H = 1.4, Re = 1000000 Time = 40 sec.

Figure 3: Numerical simulation of a uniform flow passing a rectangular dike

數據

Figure 1: Flow chart for numerical methods
Figure 4: Experimental result of a uniform flow passing a rectangular dike (Van Dyke, 1982)

參考文獻

相關文件

• One technique for determining empirical formulas in the laboratory is combustion analysis, commonly used for compounds containing principally carbon and

For ASTROD-GW arm length of 260 Gm (1.73 AU) the weak-light phase locking requirement is for 100 fW laser light to lock with an onboard laser oscillator. • Weak-light phase

Elsewhere the difference between and this plain wave is, in virtue of equation (A13), of order of .Generally the best choice for x 1 ,x 2 are the points where V(x) has

For a 4-connected plane triangulation G with at least four exterior vertices, the size of the grid can be reduced to (n/2 − 1) × (n/2) [13], [24], which is optimal in the sense

It is concluded that the proposed computer aided text mining method for patent function model analysis is able improve the efficiency and consistency of the result with

To solve this problem, this study proposed a novel neural network model, Ecological Succession Neural Network (ESNN), which is inspired by the concept of ecological succession

The objective of this study is to establish a monthly water quality predicting model using a grammatical evolution (GE) programming system for Feitsui Reservoir in Northern

Therefore, a study of the material (EPI) re-issued MO model for an insufficient output of the LED chip manufacturing plant is proposed in this paper.. Three material