81 82

Chơng 3
Giải số trị phơng trình cân bằng
năng lợng sóng
3.1. Những vấn đề giải số trị phơng trình tiến triển
năng lợng sóng
Việc giải quyết thnh công bi toán về tính v dự báo sóng
gió phụ thuộc vo chất lợng mô hình vật lý, sự hiện thực hóa
số trị đối với phơng trình cân bằng năng lợng sóng v độ
chính xác của trờng gió cho trớc. Nói chung, phải nhận xét
rằng phần lớn những mô hình sóng gió hiện hnh không phải l
những mô hình đạt trình độ cao về phơng diện sơ đồ giải số trị.
Nhiều công trình [113, 170, 172, 173, 372] đã cho thấy rằng việc
giải số phơng trình cân bằng năng lợng sóng có ý nghĩa
nguyên tắc quyết định độ chính xác v hiệu quả của mô hình.
Những sai số trong khi giải số có thể tơng đơng với sai số liên
quan tới tính thiếu tin cậy của thông tin ban đầu về gió cũng
nh sự cha hon thiện trong quan niệm vật lý về quá trình.
Phơng pháp các đặc trng m t
a đã sử dụng để giải
phơng trình (2.1) cho phép nhận đợc nghiệm phản ánh khá
chính xác những đặc điểm truyền sóng. Ta nhớ lại rằng để xác
định phổ tần số
góc ở một điểm phải "thu lợm" tất cả các
thnh phần phổ từ khắp thủy vực đi tới.
Tuy nhiên, nếu sử dụng phơng pháp đặc trng dới dạng
nh trớc
đây để tính toán sóng gió trên những thủy vực đại
dơng rộng lớn thì sẽ không hợp lý vì hai lý do. Thứ nhất, nó
không có khả năng tính tới sự tơng tác phi tuyến yếu giữa các

365], ở đây đã cố gắng giải phơng trình ny trên mặt cầu, có
tính tới khúc xạ sóng trên nớc nông v ảnh hởng của dòng
chảy. Hiện nay mô hình
WAM có lẽ l mô hình duy nhất trong
đó phơng trình đợc thực hiện đầy đủ nhất. Để xấp xỉ thnh
phần đặc trng cho sự biến thiên mật độ phổ nh một hm
hớng

, đã sử dụng sơ đồ sai phân trung tâm bậc hai.
Tuy nhiên sơ đồ số giải phơng trình (2.1) đã chọn trong mô
hình
WAM cha thật tối u. Vấn đề l ở chỗ khi tính lan truyền
các sóng phổ hớng hẹp thì sai số tính toán thnh phần
83 84

dtdS //

trong phơng trình trở nên đáng kể do khuếch tán
số lớn. Nh vậy, so với phơng trình cân bằng viết cho mặt
phẳng, đã nảy sinh thêm nguồn các sai số, để giảm thiểu chúng
phải tăng khối lợng tính toán lên rất nhiều.
Cách giải số phơng trình cân bằng
năng lợng sóng trên
mặt cầu do V. V. Rvkin đề xuất rất đáng quan tâm. Nếu thông
thờng thì ngời ta chấp nhận rằng tại mỗi nút lới tính trên
mặt cầu số góc hớng
l

nh nhau v thể hiện dới dạng phân
đều khoảng ];[

hạn chế phân giải góc của mô hình sẽ tạo ra sự bất đẳng hớng
nhân tạo trong phân bố không gian độ cao sóng. Hậu quả của
hiện tợng ny l sự lan truyền sóng lừng từ bão xa đợc dự báo
một cách không đạt.
Về sự tồn tại của vấn đề ny đã nhiều lần nhắc
tới trong các
công trình [172, 173, 217, 331, 372, 381], song cha tìm ra một
giải pháp đủ đơn giản. Thiết tởng cách thức giải quyết vấn đề
tự nhiên nhất l tăng thêm số thnh phần phổ. Theo đánh giá
của [217] độ rộng điển hình của khoảng tần )(


v khoảng góc
)( để tính sóng trên thủy vực Bắc Đại Tây Dơng phải bằng

030, v

51, . Sử dụng độ phân giải chi tiết nh vậy
trong thực tế liệu có hợp lý không?
Giải pháp thứ hai cho vấn đề khắc phục "hiệu ứng xé lẻ"
đợc
đề xuất trong công trình của các tác giả [217] cho trờng
hợp truyền sóng trên mặt phẳng. Họ đề xuất đa hai thnh
phần bổ sung vo phơng trình cân bằng năng lợng sóng cho
phép hiệu chỉnh những hiệu ứng liên quan tới sự hữu hạn của
độ rộng khoảng tần v góc khi rời rạc hóa phổ. Nhợc điểm của
các thnh phần hiệu chỉnh ny l ở chỗ phải giải một phơng
trình phức tạp hơn so với phơng trình cân bằng năng lợng
sóng truyền thống. Nó có thêm những thnh phần bổ sung với
các đạo hm riêng bậc hai v ngoi ra còn một phơng trình bổ

Courrant m không mất độ chính xác.
3.2. dẫn lập bi toán về truyền sóng để giải bằng
các phơng pháp khác nhau
Phơng trình xuất phát
. Sự tiến triển của phổ hai chiều
sóng biển ),,,,( tS



một hm số của tần số

, hớng

(ở
đây đợc đo ngợc chiều kim đồng hồ kể từ vĩ tuyến), vĩ độ ,
kinh độ

v thời gian t , đợc mô tả bằng một phơng trình, ở
mô hình
WAM [303, 381] dùng dới dạng









G



),,,,( tGG hm nguồn.
Trên cơ sở hệ (1.86)(1.90) ta viết các phơng trình chuyển động
chùm sóng dọc cung vòng tròn lớn dới dạng:
R
C
g


sin

; (3.2)



cos
cos
R
C
g

; (3.3)
R
C
g


cos
tg

Lr


, trong
đó
r
khoảng cách từ tâm nhiễu động ban đầu đến một điểm
cụ thể,

max
L bán kính tơng quan của nhiễu động ban đầu.
87 88

Trên các thủy vực đại dơng tham số ny đợc ớc lợng bằng
khoảng 1500 km [305]; để khai thác các số liệu vệ tính trên
thủy vực Bắc Hải có thể chấp nhận rằng 150
max

L km [229].
Xem rằng tại thời điểm đầu sóng truyền xuống phía nam.
Hớng truyền sóng tổng quát bằng

90
0
, độ cao sóng lớn ở
tâm bằng 10 m, chu kỳ trung bình bằng 15 s.
Phổ nhiễu động ban đầu xấp xỉ bằng công thức







),(
0
S phổ tần góc
ban đầu của sóng. Hm phân bố không gian chấp nhận tuân
theo mối phụ thuộc hm mũ







]cos[exp,
2
00
2
2
0
F
, (3.6)
trong đó
)/( LR2


v
150L
km những hằng số mô tả mức






















n
n
n
n
n
mnS
maxmax
exp
1

(3.4) vo phơng trình (3.1), phơng trình ny sau
đó có thể dẫn tới dạng bình lu (2.1), tức dới dạng đạo hm
ton phần theo thời gian. Nh đã chứng minh, trong trờng hợp
không có tác động nguồn
0

G , mật độ phổ năng lợng đợc
giữ bảo ton dọc đờng đặc trng. Nh vậy nghiệm giải tích của
phơng trình (2.1) chúng tôi đã nhận đợc ở mục trớc dới
dạng (2.14).
3.3. Khắc phục hiệu ứng "xé lẻ" nghiệm
Những hiệu chỉnh cho phơng
trình động học liên
quan tới sự gián đoạn tần
góc của biểu diễn phổ
Giải số phơng trình cân bằng năng lợng sóng (3.1) đòi hỏi
biểu diễn phổ liên tục dới dạng rời rạc hóa tần
góc. Ta ký hiệu
),(
lk
n
S l thnh phần phổ ứng với tần số
k

v góc
l
tại
bớc thời gian
n .
Năng lợng trung bình

đã cho
89 90











2
2
2
2

11
l
l
k
k
ddSS
lk
,, . (3.9)
Cách truyền thống để ớc lợng tích phân

1
ji
jlikijlk
SbaS
,
,, , (3.10)
trong đó

jj
ba , những hệ số nội suy: 24/1
1111





bbaa ;
12/11
00

ba . Năng lợng trung bình tập trung trong khoảng
ny không chỉ chứa thnh phần
),(
lk
S


, m cả những giá trị
phổ của các thnh phần kế cận.
Bằng cách tơng tự phơng trình cân bằng năng lợng sóng


O
lll
/
cos/coscos


(3.11 a)










,/
sin/sinsin
4
l
2
2
1
cos61
21

























432
11
1
1 O
ii
, (3.12)
thì phơng trình (3.1) có thể biểu diễn dới dạng


































S
g
g
tg
2121
211
2
2
cos
cos
cos
sin
cos
sin
cos
//
cos/
cos








0
22
, (3.13)
trong đó 112 112

)(
2
.
Những số hạng hiệu chỉnh của phơng trình (3.13) có thể
so
sánh với một biểu thức tơng tự đã nhận đợc trong công trình
[217] cho trờng hợp truyền sóng trên mặt phẳng. Số hạng hiệu
chỉnh nhận đợc trong công trình đó phụ thuộc vo "tuổi sóng",
một thứ không xác định địa phơng, v để xác định nó phải giải
một phơng trình bổ sung. Ưu việt chính của số hạng hiệu
chỉnh (3.13) nhận đợc trong công trình ny l ở chỗ nó xác
định địa phơng đối với lan truyền sóng trên mặt cầu. Trên mặt
phẳng phơng trình động học có thể biểu diễn tơng tự


















22















y
S
x
S
C
g
sincos
, (3.14)
trong đó


cos/
ggx
CC 21

. Vì lý do đó chúng tôi giới hạn chỉ xét hiệu ứng gián
đoạn góc, nhng những lập luận của chúng tôi cũng hon ton
có thể áp dụng cho cả trờng hợp khi cần tính tới gián đoạn tần.
Giả sử rằng
L
l quy mô không gian điển hình của lan
truyền sóng ở một vùng no đó. Với đại dơng
L
có bậc một vi
ngn km, với các biển thuộc thèm lục địa nh Bắc Hải
L
có bậc
một vi trăm km. Trên cơ sở kết quả của công trình [217] có thể
chứng minh rằng trong (3.13)













SRSS
L
1 sin















S
R
C

SS
S
t
S
g


~
~
cos
~
cos

21 /
~



;

21 /
~



;

21 /
~



.
Nh đã thấy từ phơng
trình (3.16), vế trái của nó biểu
diễn dới dạng toán tử khuếch tán thông thờng mô tả sự "trao
đổi năng lợng" yếu giữa các thnh phần góc kế cận nhau.
Tham số
A phụ thuộc vo vĩ độ

v góc hớng truyền sóng

.

trong đó

RCA
g
12
2
/
. (3.18)
Nghiệm đơn của phơng trình (3.17),
tuần hon theo biến

, có thể viết dới dạng







2
mimmBS
m
exp, . (3.19)
ở đây

)(mB hm của tham số m , đợc xác định bằng
những điều kiện đầu của bi toán. Nghiệm phơng trình (3.19)
cùng với thời gian sẽ tiến tới phổ đẳng hớng



max
t


, (3.22)
trong đó

max
t thời gian truyền sóng cực đại. Từ đây có thể ớc
lợng giới hạn trên của đại lợng
A


max
)(/ gTtmRA
22
48 . (3.23)
Thí dụ, với chu kỳ sóng trung bình 10 s, phân giải góc
12
/ ,
5

m
v thời gian phát triển sóng cực đại 36 giờ, ta
tìm đợc
2
10A
. Trị số chính xác hơn của tham số ny có thể
tìm theo kết quả thí nghiệm số.
Thuật giải số. Hiện thực hóa số đối với việc hiệu chỉnh các

n chỉ số bớc thời gian; RtAC
g
12/



v /LRA
95 96

l những tham số không thứ nguyên phụ thuộc vo các quy mô
đặc trng mô tả sóng. Nhờ tơng quan (3.21) có thể chứng minh
rằng
không phụ thuộc hiện vo độ phân giải góc


. Ta nhận
thấy rằng phơng trình (3.16) mô tả sự l trơn góc của phổ sóng
hay sự "trao đổi năng lợng" yếu giữa các thnh phần góc diễn
ra trên từng bớc truyền sóng.
Một trong những phơng pháp kiểm tra thuật giải
l kiểm
tra sự bảo ton năng lợng ton phần. Để đánh giá biến thiên
năng lợng cho hai hớng, ta lấy tích phân số đối với phổ theo
tất cả các góc







lk
n
vSSv (3.25)
Vì phổ
),(
lk
S


l hm tuần hon của các hớng
l

, nên có
thể chứng minh rằng toán tử (3.25) bảo tồn năng lợng ton
phần, do







L
l
lk
n
L
l
lk
n

),(
ji


v chùm sóng với tần
số
k
v hớng truyền
l

. Tọa độ điểm đầu ),(
00
ji
m chùm
sóng nằm tại đó ở thời điểm trớc, có thể nhận đợc nếu sử
dụng phơng trình (3.2)
(3.4). Điểm ny sẽ không trung với nút
của mạng lới đều. Vì cần phải xác định giá trị đầu của phổ
0
S
ở điểm
),(
00
ji
, nên ta sẽ dùng nội suy đa thức






),(
qp


ở bớc thời gian trớc
1

n
tt . Trong trờng
hợp nội suy tuyến tính kép, nội suy ny l tối u trong lớp nội
suy đa thức đối với nghiệm bi toán đang xét [170], ta sử dụng
2

L
M
.
Có một sự phức tạp bổ sung trong vấn đề nội suy đó l do
tình huống hớng truyền sóng

không phải l hằng số dọc theo
đờng đặc trng khi truyền sóng trên mặt cầu. Từ tơng quan
(2.8) biến thiên góc

dọc đặc trng có thể viết dới dạng
00
ilil
coscoscoscos . (3.28)
V một lần nữa góc đầu
0
l

, v giá
trị tơng ứng của phổ cũng đợc xác định. Đối với những nút
nằm trên vĩ độ khác, thí dụ
1

j
, góc phân bố tơng đơng có
thể nhận đợc từ tơng quan (2.8).
Để xác định giá trị phổ ứng với góc
0
l
có thể sử dụng nội suy




1
1
11
0
11
m
mljimlji
SaS )()(
,,
(3.29)
trong đó

m
a những hệ số nội suy.

miền lới bằng khoảng 55 km. Miền ny có thể xem nh dạng
đơn giản hóa của thủy vực hai biển Na Uy v Bắc Hải. Hình
dạng bờ đơn giản hóa cho phép có đợc nghiệm giải tích chính
xác của bi toán.
Bi
toán truyền sóng từ nhiễu động ban đầu (3.5)(3.8)
đợc giải tuần tự bằng giải tích v sử dụng sơ đồ số của mô hình
WAM v phơng pháp INTERPOL. Ta sẽ so sánh những kết
quả tính của các phơng pháp bằng cách lập những đại lợng
tích phân liệt kê dới đây:
Trị số trung bình độ cao sóng ),,( th


xác định bằng






ddtSth 2
2
,,,,,, . (3.30)
Trị số quy chuẩn tổng năng lợng trờng sóng lan truyền
xác định bằng
)0(/)()(



tEtEt , (3.31)


)0(/)( TtTt , (3.34)
99 100

trong đó


d d
2
cos,,)( RtFtT , (3.35)
với


),,( tF hm Hevisai








ththFtF
max
,,),,(
3
1
. (3.36)
Sai số bình phơng trung bình (
RMS) tính toán độ cáo sóng




, (3.38)
trong đó
)(th
max
anal
độ cao sóng cực đại trên ton thuỷ vực tại
thời điểm đang xét, tính bằng giải tích.
Các kết quả số trị. Phân bố không gian ban đầu của độ
cao sóng trình by trên hình 3.1 a. Kết quả tính độ cao sóng sau
24 giờ lan truyền bằng nghiệm giải tích, bằng sơ đồ số của mô
hình
WAM
*
v phơng pháp INTERPOL tuần tự thể hiện trên
các hình 3.1 b, c, d. Trong các tính toán số đã sử dụng 12 hớng,
bớc thời gian bằng 120 ph. Phân bố không gian của sai số quy
chuẩn (bản đồ sai số)
ERR (theo (3.38)) tại thời điểm 24

t giờ
thể hiện trên hình 3.2. Các hình 3.1c, d cho ta khái niệm về đặc
điểm diễn biến của nghiệm số trị. Các kết quả tính theo sơ đồ số
thứ nhất v các kết quả của phơng pháp
INTERPOL, ở mức độ

*
Tính toán theo mô hình WAM do Janette Onblee cộng tác viên Viện Khí

tổng quát. Nếu hớng ny chính xác nằm giữa hai hớng cơ sở
của mô hình, thì phân bố không gian độ cao sóng phân rã thnh
cấu trúc hai bớu. Tuy nhiên, mức độ bất đẳng hớng v trị số
sai số tính độ cao sóng trong trờng hợp ny không khác mấy so
với trờng hợp đã biểu diễn trên hình 3.1.
Khác biệt tơng đối giữa các kết quả
số của mô hình WAM,
phơng pháp
INTERPOL v nghiệm giải tích đợc thể hiện
định lợng trên hình 3.2. Mức bất đảng hớng cao gây bởi sự
phân giải góc hạn chế (xem các hình 3.1 b
d) thể hiện rõ nét
trong phân bố sai số quy chuẩn. Điều ny đặc biệt rõ nét trong
trờng hợp sơ đồ số của mô hình
WAM 12 hớng đối với thời
điểm
24t giờ, khi đó giá trị số trị bị cao lên 40% dọc các hớng
cơ sở v bị thấp đi 35% so với các hớng khác. Những sai số ny
tăng dần theo thời gian, đặc biệt ở hớng truyền sóng tổng quát:
kết quả số bị tăng lên +25% tại
12

t giờ, v tại một số điểm sai
số đạt tới 90% sau
48t
giờ. ở các hớng bên giá trị số trị bị
giảm đi
(1520%) tại
12t
giờ v tới 40% tại

Những tham số tích phân của nghiệm: tiến triển thời gian
của tổng năng lợng trờng sóng, vị trí tâm của nó, mức độ
103 104

khuếch tán năng lợng v sai số bình phơng trung bình của
nghiệm số trị (các biểu thức (3.32)
(3.34), (3.37)) đợc biểu diễn
trên hình 3.3. Các kết quả chứng tỏ rằng cả hai sơ đồ số
WAM
v
INTERPOL có khả năng tái hiện đặc điểm diễn biến của hai
tham số đầu tiên trong số các tham số đã nêu một cách khá tốt
(xem hình 3.3 a, b). Các tính toán dựa theo phơng pháp
INTERPOL có phơng sai ít nhiều lớn hơn so với kết quả của
mô hình
WAM. Mức năng lợng sóng tổng cộng nhận đợc theo
mô hình thứ nhất trở thnh nhỏ hơn so với mô hình thứ hai.
Phân bố năng lợng trên miền lới khá phù hợp với nghiệm
chính xác, điều ny có thể thấy rõ theo diễn biến của hm
khuếch tán năng lợng
)(t (xem hình 3.3 c).
Đáng chú ý nhất l đặc điểm diễn biến của sai số
bình
phơng trung bình
RMS (xem hình 3.3d). Tại những giai đoạn
đầu truyền sóng (trớc 24 giờ) xảy ra sự tăng đơn điệu sai số
tính toán theo tất cả các phơng pháp số. Tuy nhiên sau đó sai
số bắt đầu giảm, liên quan tới sự di chuyển nhiễu động ra ngoi
miền lới số trị. Sai số bình phơng trung bình (theo ton miền)
của mô hình

dẫn tới sai số 10% tại
40t
giờ. Nh vậy sai số tính toán giảm
nếu tăng bớc thời gian.

Hình 3.3. Biến trình thời gian các tham số tích phân của nghiệm số trị
a) tổng năng lợng quy chuẩn (t); b) tọa độ vĩ độ trên dịch chuyển tâm miền
nhiễu <
(t)>; c) mức khuếch tán năng lợng theo không gian (t); d) sai số
bình phơng trung bình độ cao sóng tính toán RMS(t)
1
nghiệm giải tích; 2 mô hình WAM 12 hớng, bớc thời gian 20 ph; 3 mô
hình WAM 24 hớng, bớc thời gian 20 ph; 4
mo hình WAM 12 hớng, bớc
thời gian 20 ph, xê dịch các hớng cơ sở


/2; 5 phơng pháp INTERPOL 12
hớng, bớc thòi gian 20 ph; 6
phơng pháp INTERPOL 12 hớng, bớc
thòi gian 3 giờ; 7
phơng pháp INTERPOL 12 hớng, bớc thòi gian 6 giờ.
Để nghiên cứu sự phụ thuộc của kết quả tính vo hình dáng
không gian của nhiễu động ban đầu đã lặp lại tính toán cho một
loạt hm v mức thuyên giảm xa dần tâm nhiễu. Thấy rằng với
miền lới đã cho, hình dáng của nhiễu ban đầu không ảnh
hởng đáng kể tới sự phân bố v tiến triển độ cao sóng. Về định
tính, những chi tiết phân bố không gian độ cao sóng v diễn
biến của các tham số tích phân l tơng tự nh nhau đối với tất
cả các dạng nhiễu ban đầu khác nhau. Tuy nhiên, về định

bậc hai [170, 331], kể cả những sơ đồ Runge-
Kutta. Ta lu ý rằng sử dụng những sơ đồ số bậc cao sẽ trở
thnh không hiệu quả trong các mô hình sóng gió, vì tính cồng
kềnh thực hiện tính toán sóng theo các trờng gió.
Vấn đề còn ở chỗ các vùng khác nhau của phổ sóng gió tiến
triển với tốc độ khác nhau, tốc độ ny về phía mình lại phụ
thuộc vo
độ lớn tốc độ gió v tần số sóng. Điều đó dẫn đến chỗ ở
vùng phổ cao tần các quá trình hình thnh phổ diễn ra khá
nhanh. Để tính toán số trị với chúng đòi hỏi sử dụng bớc thời
gian khá bé, lm tăng số bớc v khối lợng tính toán.
Trong các mô hình số khác nhau thì vấn đề ny đợc giải
quyết một cách khác nhau. Đại đa số trờng hợp, về
sự tăng
trởng phổ ở vùng tần cao ngời ta đề ra những giới hạn nhất
định đối với độ lớn của mật độ phổ, đòi hỏi lm sao để nó không
vợt trội những giá trị của khoảng cân bằng [331].
Trong sơ đồ số của mô hình
WAM [303, 365] phổ đợc chia
ra thnh hai vùng: vùng dự báo, bao gồm vùng cực đại phổ v
vùng tần thấp của nó, v vùng chẩn đoán mô tả phần đuôi tần
cao của phổ. Ngời ta đề ra hai điều kiện đối với độ lớn phổ ở
vùng hai. Điều kiện thứ nhất l: Bắt đầu từ một tần số no đó
mật độ phổ năng lợng đợc cho bằng mối phụ thuộc
)(S
5
.
Ngoi ra, đòi hỏi tốc độ tiến triển mật độ phổ năng lợng sóng
(hay hm nguồn) không vợt trội một giá trị no đó, giá trị ny
trong khi hiện thực số trị mô hình

trở thnh nguồn bổ sung hay nguồn tiêu năng lợng không?
Ta lu ý một điều quan trọng nữa, liên
quan tới tính toán
sóng gió thực tế. Khi thực hiện những tính toán chẩn đoán v
dự báo sóng gió thì thông tin ban đầu về tốc độ gió đợc đa vo
mô hình với một hạn thời gian no đó, thờng trùng với hạn
synop. Hạn ny có thể khác nhau tuỳ thuộc vo một loạt hon
cảnh, thí dụ, tuỳ thuộc vo khoảng thời gian m thông tin đợc
truyền từ những trung tâm thời tiết hạn vừa hay khoảng thời
gian lập bản đồ synop. Thờng khoảng thời gian ny bằng 12
giờ hay 6 giờ, trờng hợp thuận lợi nhất
3 giờ.
Vì vậy, một sơ đồ tối u nhất l
sơ đồ số giải phơng trình
cân bằng năng lợng sóng trong đó với số lần lặp nhỏ nhất nhận
đợc nghiệm số trị chính xác nhất của bi toán cho thời điểm
trùng với hạn synop, tức vo thời điểm kết xuất kết quả tính các
yếu tố trờng sóng. Ta hình dung rằng một sơ đồ số "lý tởng"
phải đa ra giá trị tính toán chính xác nhất vo hạn synop
tơng ứng sau một lần lặp số trị. Tuy nhiên điều ny chắc gì đã
đạt đợc, hoặc do nguyên nhân bất ổn định số, hoặc do độ chính
xác tính toán cha đủ. Có lẽ l ở đây cần sự thoả hiệp giữa bớc
thời gian nhập thông tin ban đầu, số luợng bớc lặp số trị v trị
số sai số tính toán.
Những sơ đồ số tích phân hm
nguồn của phơng trình cân
bằng năng lợng sóng. Nhằm mục đích khảo sát độ chính xác v
tính tối u của các sơ đồ số khác nhau, chúng tôi sẽ dẫn ra
những sơ đồ số thờng dùng để tích phân hm nguồn của
phơng trình cân bằng năng lợng trong các mô hình toán sóng

n
U

.
Sơ đồ hiện hai bớc của Adams (phơng pháp dự báo
sửa
sai) đợc tiến hnh thnh hai bớc tuần tự:
),(
)1(
1
nnn
n
UStGSS



; (3.41a)


),(),(
2
1
1
)1(
1
)2(
1




nửa ẩn chúng tôi dẫn sơ đồ số đợc đề xuất trong mô hình
WAM
[303, 365]. Sơ đồ ny dựa trên sử dụng công thức hình thang ẩn

tUSGUSGSS
nnnnnn


),(),(
2
1
111
. (3.42)
Để có đại lợng
1n
S dới dạng hiện phải giải phơng trình




tUSGtStUSGS
nnnnnn



2
1
2
1
111





S
S
G
GG
n
nn
(3.44)
Đạo hm phiếm hm trong (3.44) đợc biểu diễn dới dạng
ma trận đờng chéo
n

v không đờng chéo
n
N

nn
n
N
S
G



. (3.45)
Thế (3.45) vo (3.43) sẽ đa đến biểu thức sau đây:


đóng góp của các số hạng không đờng chéo trong (3.46) tỏ ra
khá nhỏ để có thể bỏ qua chúng.
Nh vậy biến thiên của mật độ phổ trên bớc thời gian
sẽ
bằng

1
11
)(
2
1
1),(),(
2
1









nnnnnn
UttUSGUSGS . (3.47)
Khi sử dụng sơ đồ nửa ẩn (3.47), ngoi hm nguồn, còn phải
tính đạo hm của nó (3.43).
Lu ý rằng các sơ đồ hiện (3.40)
(3.41) v nửa ẩn (3.47) có
thể đợc dùng với những dạng khá tổng quát của hm nguồn

gió
U

.
Trong tình huống đơn giản nhất, khi có thể giả thiết rằng
phần đóng góp chủ yếu vo hm nguồn l số hạng tuyến
tính
thep phổ
SAG
i

, (trong đó

i
A hệ số tổng quát ), sử dụng
phơng pháp hình thang (3.42) sẽ đa tới biểu thức sau đây để
xấp xỉ đại lợng phổ
1
1
2
1
1
2
1
1






có thể
cũng phụ thuộc vo thời gian t . Thí dụ, với sóng gió có thể chấp
nhận rằng
B
l gia số tăng trởng năng lợng sóng nhờ nạp
năng lợng từ gió cho sóng (theo lý thuyết Miles). Tham số
c
hạn chế sự gia tăng năng lợng do quá trình tiêu tán sóng bằng
một giá trị tới hạn no đó. Nếu ta chấp nhận khoảng cân bằng
),(

S
lm giá trị tới hạn, thì



Sc . Chẳng hạn, ta thấy
trong cơ chế tiêu tán sóng do O. Phillips [336] đề xuất v liên
quan tới quá trình đổ nho đỉnh sóng, thì hm tiêu tán phụ
thuộc vo phổ dới dạng hm lập phơng.
Trong trờng hợp ny
việc xác định công thức của sơ đồ số
ẩn giải phơng trình cân bằng năng lợng sóng quy về việc giải
phơng trình đại số

tUSGStcSBS
nnnnnn




1
4
2
Aabb
A
S
n



, (3.52)
trong đó
1
2
1
12


tBb
n
;
2
12
)(
2
1
tcBa
n



nghiệm các sơ đồ số đã dẫn trên đây, so sánh nghiệm số trị với
giá trị giải tích chính xác cho trờng hợp giá trị ny có thể nhận
đợc dới dạng hiện. Chẳng hạn, có thể nghiệm chính xác cho
trờng hợp gió không đổi v hm nguồn viết dới dạng (3.50).
Vậy ta giả thiết rằng gió tốc độ không đổi 15 m/s thổi trên mặt
biển rộng vô hạn. Theo phơng trình cân bằng năng lợng sóng
(3.1) sẽ không có sự bình lu năng lợng sóng v sự phát triển
sóng diễn ra chỉ với thời gian. Để nhận nghiệm có thể sử dụng
những công thức đã dẫn ở trên.
Nghiệm giải tích chính xác của phơng trình động học với
hm nguồn (3.50) đợc viết dới dạng





/
)(
1
000
1
Bt
ecScSStS , (3.53)
trong đó lấy dấu (+) khi

cS1 lớn hơn không, lấy dấu () trong
trờng hợp ngợc lai.
Trên hình 3.4 dẫn những giá trị mật độ phổ năng lợng
tính
theo một số thuật giải số cho thời điểm 18 phút sau khi bắt

Những tính toán tơng tự đợc lặp lại với các bớc tích
phân số bằng 1, 3, 6 v thậm chí 12 giờ. Những kết quả tính
theo
tất cả ácc sơ đồ tỏ ra ổn định. Tuy nhiên với những bớc
tích phân số lớn thì những kết quả của sơ đồ ẩn (3.52) tỏ ra
đáng chấp nhận hơn cả.
Nghiệm giải tích (3.53) cho phép ớc lợng định lợng sai
số tơng đối của nghiệm số trị. Nó có
thể đợc xác định trê từng
bớc thời gian nh l hiệu số tổng cộng giữa nghiệm số trị v
nghiệm giải tích đợc quy chuẩn với trị số cực đại của nghiệm
giải tích:














MN
kj
kj
kjkji


t , bớc tích phân ph 3


t
115 116 Hình 3.5. Các giá trị mật độ phổ tính theo các thuật toán số trị
với bớc tích phân ph 20t cho thời điểm 3 giờ
1 nghiệm giải tích; 2 phơng pháp nửa ẩn (3.47); 3 phơng pháp bậc tối u
(3.49); 4
sơ đồ ẩn (3.52). Đờng gạch nối chỉ các giá trị khoảng cân bằng
Lu ý rằng quy chuẩn theo trị số cực đại mật độ phổ của
nghiệm giải tích sẽ cho "ớc lợng từ dới". Nếu nh quy chuẩn
không theo trị số cực đại, m theo trị số hiện tại của mật độ
phổ, thì sai số sẽ lớn hơn rất nhiều.
Trên hình 3.6. biểu diễn biến thiên theo thời gian của sai số
tính toán với các phơng pháp khác nhau v bớc tích phân
theo thời gian 1 giờ.
Xu thế chung giảm sai số tính toán liên quan tới sự tiến
dần
các giá trị số trị tới nghiệm giải tích v tăng trị số cực đại
của mật độ phổ dùng để quy chuẩn sai số tính toán trong (3.54).

Hình 3.6. Những giá trị sai số RMS của các phơng pháp,
bớc tích phân theo thời gian 1


t giờ

WAM (3.47), tỏ ra
ổn định đối với những bớc thời gian khá lớn. Tuy nhiên lu ý
rằng khi tăng bớc tích phân theo thời gian độ chính xác giảm.
Sơ đồ ẩn (3.51) tỏ ra ổn định nhất. Nó biểu hiện kết quả ổn
định không những với bớc thời gian
60


t
ph, m cả với
3t
giờ, thậm chí
12t
giờ. Về độ chính xác: Mặc dù trên
những bớc tích phân đầu tiên sai số tính toán có thể khá lớn,
nhng sau đó, sau khi qua một "ngỡng" no đó, sai số giảm
dần v chấp nhận giá trị hon ton thoả mãn. Để có kết quả với
độ chính xác thoả mãn chỉ cần thực hiện 4
6 bớc lặp số.
Phơng pháp phân rã. Từ trớc đến bây giờ đã xét những
nghiệm số trị của phơng trình cân bằng năng lợng sóng trong
ddó không tính tới hm tơng tác phi tuyến yếu
nl
G . Về hm
ny sẽ bn luận ở chơng sau.
ở đây chỉ lu ý rằng tính toán
hm tơng tác phi tuyến yếu
nl
G
l một công việc rất phức tạp

SGtS



/ , (3.56)
v viết dạng triển khai số trị của nó nh sau
nl
nn
nn
GSBcBS
t
SS





)(
1
1
. (3.57)
Sơ đồ ny có thể đợc sử dụng với mọi 10



. Nhớ lại
rằng khi 0


ta có sơ đồ hiện Euler, khi


nn
nn
SBcBS
t
SS
)(
1
1
. (3.59)
Sơ đồ với các bớ
c phân số sẽ trở về sơ đồ (3.57) nếu nh từ
phơng trình thứ nhất (3.58) rút ra

n
S
v thế vo vế trái của
phơng trình thứ hai (3.59) hay lấy tổng từng số hạng hai
phơng trình ny. Sơ đồ sai phân (3.59) thực tế l một sơ đồ
Euler đòi hỏi bớc thời gian nhỏ so với những sơ đồ khác đã biết.
Mặt khác, ta có thể chuyển đổi ngợc lại từ phơng trình
sai
phân sang phơng trình vi phân nếu giải nó cho bớc
t )(1 . Ta viết lại phơng trình (3.59) dới dạng







với điều kiện đầu
n
S tìm theo (3.58). Phơng trình (3.61) có
nghiệm giải tích (3.53) trong đó đại lợng
B
đợc thay thế bằng
)/( 1B .
Sơ đồ với bớc phân số (3.58)
v (3.59), trong đó phơng
trình thứ hai đợc thay thế bởi phơng trình vi phân (3.61) với
0
, từ nay ta sẽ gọi l phơng pháp phân rã hay phơng
pháp bán giải tích, vì nghiệm của phơng trình (3.61) tồn tại
dới dạng giải tích (3.53). Có thể sử dụng nghiệm ny trên từng
bớc thời gian, khi phải giải phơng trình chung (3.57).
Những kết quả thử nghiệm các sơ đồ số giải phơng
trình cân bằng năng lợng sóng có tính tới hm vận
chuyển năng lợng phi tuyến yếu
. Ta sẽ thử nghiệm
phơng pháp phân rã, so sánh các nghiệm số trị với những kết
quả nhận đợc theo một sơ đồ chính xác nhất trong số những sơ
đồ hiện đơn giản
phơng pháp dự báo sửa sai (3.41). Khác
với những thử nghiệm trớc (xem mục 4.1), bây giờ ta không thể
có đợc nghiệm giải tích chính xác của phơng trình cân bằng
năng lợng sóng chứa hm vận chuyển năng lợng phi tuyến
yếu. Vì vậy, thay vì nghiệm giải tích ta sẽ sử dụng một nghiệm
số trị chính xác nhất.
Trong tính toán, hm cung ứng năng lợng từ gió cho
sóng

bằng 1, 3 v thậm chí 6 giờ. Vì vậy không cần thiết phải đa ra
những hạn chế no đó đối với nghiệm hay đối với hm nguồn
nói chung.
Trên hình 3.7 biểu diễn những kết quả tính toán số trị với
phổ tần số cho thời
điểm
30

t giờ. Tính toán đợc thực hiện
theo phơng pháp dự báo
sửa sai v phơng pháp phân rã với
121 122

những bớc tích phân thời gian khác nhau. Nh đã thấy từ
phép so sánh, các kết quả của phơng pháp dự báo
sửa sai
(với bớc tích phân 3 ph) thực tế chính xác trùng với các kết quả
tính của phơng pháp phân rã với bớc tích phân bằng 1 giờ.

Hình 3.7. Những nghiệm số trị của mật độ phổ tại thời điểm 30

t giờ
tình theo các phơng pháp khác nhau
1 nghiệm không tính tới
nl
G ; 2 phơng pháp dự báo sửa sai, bớc tích phân 3 ph;
3
nghiệm với bớc tích phân 1 giờ; 4 nghiệm với bớc tích phân 3 giờ; 5 nghiệm với
bớc tích phân 6 giờ
Ngoi ra, khi tăng bớc tích phân đến 3 v 6 giờ thì độ lớn

để tích phân phơng trình cân bằng năng lợng sóng trong
trờng hợp vế phải có mặt hm tái phân bố phi tuyến yếu trong
phổ sóng gió. Đơng nhiên ở đây nảy sinh những vấn đề bổ sung
về tính căn cứ của hm nguồn khi nó chứa phần nạp năng lợng
v phần tiêu tán, vì dạng giải tích của các phần ấy sẽ cho phép
ta có nghiệm giải tích (3.53). Về vấn đề luận cứ lựa chọn hm
nguồn chúng ta sẽ trở lại trong chơng sau, ở đây chỉ lu ý rằng
hm tiêu tán m ta đã chấp nhận có chứa một loạt tham số tự
do, lựa chọn chúng đúng đắn sẽ lm cho mô hình thích hợp với
nhiều dạng hm mô tả cơ chế phi tuyến tiêu tán năng lợng
sóng gió.

Hòa hợp tích phân hm nguồn với sơ đồ tính truyền

năng lợng sóng
. Những tính toán tiến triển sóng nêu ra ở
trên đã đợc thực hiện cho một điểm không gian v một trờng
gió đồng nhất. Khi đó nảy sinh một câu hỏi tự nhiên liệu có
đúng không nếu sử dụng cách tiếp cận ny vo trờng hợp phức
tạp hơn: tính sóng tại thủy vực cụ thể, khi trờng gió không
đồng nhất v không dừng. Vấn đề còn ở chỗ: khi dùng bớc tích
phân lớn theo thời gian, sóng có thể đi qua những khoảng cách
m trên đó trờng gió không thể đợc xem l đồng nhất v dừng
đợc nữa.
Sử dụng phơng pháp nội suy
tia kết hợp với một phơng
pháp số trị hữu hiệu để tích phân hm nguồn cho phép giảm
nhẹ việc giải bi toán v tăng đáng kể bớc tích phân số. Dĩ
nhiên trong khi đó bớc thời gian có thể bị giới hạn bởi điều
kiện độ chính xác của nghiệm, cũng nh quy mô biến đổi trờng

thì mật độ phổ không giữ nguyên không đổi dọc theo quỹ đạo
truyền chùm sóng. Để lấy tích phân số hm nguồn dọc đờng
đặc trng, có thể dùng một trong những phơng pháp đã mô tả
ở trên (phơng pháp dự báo
sửa sai (3.41), phơng pháp nửa
ẩn (3.47) hay phơng pháp phân rã (3.58), (3.59).
ở đây điều quan trọng phải thấy rằng: tại thời điểm đầu
của mỗi bớc thời gian, ở điểm
),(
00
ji
mật độ phổ năng lợng
0
ij
S , hm nguồn
0
ij
G v tốc độ gió U

đợc xác định bằng phép nội
suy kiểu (3.27) với
n
tt

v tại thời điểm cuối của bớc những
đại lợng ny đợc lấy ở nút
),(
ji



với những phơng pháp nội suy v những sơ đồ sai phân hữu
hạn đơn giản nhất.
Những thí dụ tính toán số về lan truyền sóng lừng thực
hiện trong công trình ny
đã cho thấy rằng sơ đồ số giải phơng
trình cân bằng năng lợng sóng m trong
WAM đã dùng có thể
mắc những sai số đáng kể khi tính lan truyền năng lợng sóng.
Phơng pháp nội suy
tia đợc đề xuất ở đây l một phơng
pháp số thay thế.
Sự khác biệt căn bản giữa sơ đồ số của mô hình
WAM v
phơng pháp nội suy
tia l ở chỗ: phơng pháp thứ nhất dựa
trên việc sử dụng phép xấp xỉ sai phân hữu hạn phơng trình
cân bằng năng lợng sóng v độ ổn định của phơng pháp ny
bị hạn chế bởi điều kiện Courant
Levi. Nó hạn chế bớc tích
phân trong mô hình bởi điều kiện


)/();/cos(
gg
CRCRt




,

ny, ton bộ số bớc trên đó sai số đợc tích luỹ, sẽ nhỏ hơn.
Phơng pháp nội suy
tia u việt không chỉ so với sơ đồ số
bậc một ngợc dòng đợc dùng trong mô hình
WAM. Nó có thể
đợc xem nh một phơng pháp song đề tổng quát để tính sóng
gió cả khi có mặt hm nguồn. Kết hợp nó với phơng pháp tích
phân số ẩn hay nửa ẩn đối với vế phải phơng trình cân bằng
năng lợng sóng cho phép lập ra một sơ đồ tối u nhất để giải
bi toán trong trờng hợp chung nhất. Sơ đồ ny l sơ đồ khá ổn
định đối với những bớc thời gian lớn, cho phép tính đến hm
nguồn dới dạng đầy đủ, bao gồm tơng tác phi tuyến yếu trong
phổ sóng gió, v tạo ra khả năng: bằng một số bớc lặp tối thiểu
nhận đợc nghiệm chính xác nhất cho thời điểm trùng với hạn
synop, tức thời điểm kết xuất kết quả tính các yếu tố sóng.