近海系統(tǒng)化生態(tài)－沉積－環(huán)境動(dòng)力學(xué)數(shù)值模式——理論基礎(chǔ)＊

蘆 靜，劉學(xué)海，滕 涌，蒲新明，崔迎春，辛 明

（1.國(guó)家海洋局 第一海洋研究所，山東 青島266061；2.海洋環(huán)境科學(xué)與數(shù)值模擬國(guó)家海洋局重點(diǎn)實(shí)驗(yàn)室，山東 青島266061；3.中國(guó)海洋大學(xué) 海洋環(huán)境學(xué)院，山東 青島266003；4.海洋生態(tài)環(huán)境科學(xué)與工程國(guó)家海洋局重點(diǎn)實(shí)驗(yàn)室，山東 青島266061；5.海洋沉積與環(huán)境地質(zhì)國(guó)家海洋局重點(diǎn)實(shí)驗(yàn)室，山東 青島266061）

近幾十年來(lái)，各類海洋開發(fā)、海水增養(yǎng)殖、污染排放等使我國(guó)近海生態(tài)功能退化。大量陸源入海沉積物、通過(guò)再懸浮進(jìn)入水體的海底沉積物及泥沙顆粒吸附物，也直接或者間接影響了生態(tài)環(huán)境；養(yǎng)殖規(guī)模和養(yǎng)殖密度的增加及其不合理布局影響了養(yǎng)殖自身的健康，造成養(yǎng)殖環(huán)境壓力增大［1］；作為當(dāng)前研究熱點(diǎn)的氣候變化和碳循環(huán)同樣與生態(tài)環(huán)境密切相關(guān)。以上科學(xué)問(wèn)題都需要加強(qiáng)對(duì)生態(tài)及與之相關(guān)的沉積過(guò)程的研究，以評(píng)價(jià)人類活動(dòng)對(duì)生境要素的影響及資源利用的合理配置。如何在耦合物理、生物、地質(zhì)、化學(xué)各過(guò)程的基礎(chǔ)上，建立用以有效評(píng)價(jià)和預(yù)測(cè)生態(tài)系統(tǒng)平衡、演變的數(shù)學(xué)模型成為沿海各國(guó)可持續(xù)發(fā)展的科研戰(zhàn)略。

目前常用的水動(dòng)力學(xué)模式如ROMS、ECOMSED、FVCOM等都添加了一些功能模塊。如ECOMSED模式以泥沙模塊為優(yōu)勢(shì)，同時(shí)包含水質(zhì)模塊；FVCOM的三角形網(wǎng)格在模擬近岸的過(guò)程中可以更貼近岸線，其借鑒了ROMS泥沙模塊，并增加UNSW－sed泥沙模塊可供選擇。另外，工程上應(yīng)用非常普遍的商業(yè)模式如Delft3d、MIKE21、EFDC等軟件型模式包括的功能模塊更為全面，但由于代碼不開放，不便于其他用戶和研究人員根據(jù)自己的需求進(jìn)行二次開發(fā)和模型設(shè)計(jì)。以往生態(tài)模式多限于傳統(tǒng)的NPZD（浮游植物－浮游動(dòng)物－營(yíng)養(yǎng)鹽－碎屑）模型和簡(jiǎn)單的水質(zhì)模型，且獨(dú)立性的生態(tài)、沉積物輸運(yùn)和環(huán)境動(dòng)力學(xué)模型在模塊間相互方面存在不足。實(shí)際上，泥沙及其他懸浮物的沉積懸浮過(guò)程直接影響著各生態(tài)要素的變化變遷，如不進(jìn)行耦合顯然是不足的。此外，我國(guó)使用的多是國(guó)外開發(fā)的各種模式，但針對(duì)中國(guó)海域發(fā)展自己的模型卻相當(dāng)不足。另外，海水養(yǎng)殖是我國(guó)近岸海域一個(gè)重要的海洋開發(fā)利用產(chǎn)業(yè)，建立能夠刻畫養(yǎng)殖在生態(tài)過(guò)程中的作用并能對(duì)合理的養(yǎng)殖進(jìn)行評(píng)價(jià)評(píng)估的模式也是滿足于研究中國(guó)獨(dú)特的海洋生態(tài)環(huán)境所急需的。

基于以上考慮，本研究擬發(fā)展一套海洋生態(tài)－沉積－環(huán)境動(dòng)力學(xué)模式，特別是將生態(tài)過(guò)程與懸沙及底質(zhì)的懸浮及沉積過(guò)程、養(yǎng)殖與生態(tài)過(guò)程及水動(dòng)力過(guò)程耦合起來(lái)，且較完備地考慮了海浪作用。模式將代碼開放，便于用戶使用和二次開發(fā)。本文作為系統(tǒng)模式研究的先行工作，先給出建立模式的理論基礎(chǔ)，包括構(gòu)建模型系統(tǒng)的概念框架，模式中典型的物理、生物、化學(xué)、沉積過(guò)程的理論介紹及其主要參數(shù)化表達(dá)方式，并給出初步模擬結(jié)果。

1 模型系統(tǒng)框架

系統(tǒng)化海洋生態(tài)－沉積－環(huán)境動(dòng)力學(xué)耦合數(shù)值模式建立在水動(dòng)力模式基礎(chǔ)上，包含6個(gè)功能模塊（圖1），各模塊通過(guò)相互作用（圖2），協(xié)同構(gòu)成一個(gè)功能全面的系統(tǒng)模式。水動(dòng)力模式：模擬三維潮流、環(huán)流、溫度、鹽度、水位。模式有自帶的波浪模塊，也可將外部海浪模式的輸出結(jié)果作為本模式的輸入文件以備調(diào)用。波浪模塊輸入要素包括波致混合系數(shù)、有效波高、波周期。水動(dòng)力模塊輸出的流速和混合系數(shù)作為計(jì)算生態(tài)、沉積物擴(kuò)散的基礎(chǔ)。

沉積模塊：可模擬懸浮物濃度、底床升降、底切應(yīng)力。波浪模式輸入的有效波高和波周期用于波致底切應(yīng)力的計(jì)算，波流耦合底切應(yīng)力用于計(jì)算物質(zhì)的再懸浮或沉積通量。泥沙濃度模擬結(jié)果參與水動(dòng)力模型中水體密度的計(jì)算，并用于生態(tài)模塊中光衰減系數(shù)的計(jì)算。

生態(tài)模塊：基本變量為NPZD（浮游植物－浮游動(dòng)物－營(yíng)養(yǎng)鹽－碎屑），以氮循環(huán)或碳循環(huán)為物質(zhì)和能量流，考慮了外源輸入及界面交換過(guò)程。變量中，浮游植物能包括多個(gè)藻種，營(yíng)養(yǎng)鹽包括磷酸鹽、硅酸鹽及不同形態(tài)的無(wú)機(jī)氮和溶解有機(jī)物。浮游生物考慮了其生長(zhǎng)、死亡、代謝、捕食等過(guò)程。另外，還設(shè)置了赤潮發(fā)生的閾值，具有赤潮模擬功能；碳循環(huán)模擬中，可給出pCO2、CO2通量和堿度。

精細(xì)化水質(zhì)模塊：在生態(tài)模型的基礎(chǔ)上對(duì)一些變量和生物化學(xué)過(guò)程細(xì)化，增設(shè)了溶解有機(jī)物、DO（溶解氧）、化學(xué)耗氧量（COD）、微生物（主要是細(xì)菌）等要素及硝化、氨化等過(guò)程。考慮了底棲環(huán)境（生物及化學(xué)物質(zhì)）、微量元素如鐵對(duì)浮游植物生長(zhǎng)的限制、各化學(xué)要素不同形態(tài)的轉(zhuǎn)化及與生物之間的物質(zhì)轉(zhuǎn)化，生物化學(xué)要素的海底交換過(guò)程利用泥沙沉積模塊通過(guò)應(yīng)力關(guān)系間接計(jì)算沉降和懸浮率。

養(yǎng)殖模塊：模型中考慮養(yǎng)殖對(duì)水動(dòng)力的影響及養(yǎng)殖生物體參與的生物化學(xué)過(guò)程，即養(yǎng)殖的下行控制作用。本模塊可與其他所有模塊耦合。

擴(kuò)散輸運(yùn)及粒子追蹤模塊：主要針對(duì)不詳細(xì)考慮生物化學(xué)過(guò)程的各種污染物、溫排水和濃鹽水等的擴(kuò)散，可通過(guò)較為簡(jiǎn)單的參數(shù)化（而不是多變量之間轉(zhuǎn)化的模擬）考慮降解、風(fēng)化等必要的非保守過(guò)程，并可給出擴(kuò)散穩(wěn)定后的濃度包絡(luò)。

評(píng)估評(píng)價(jià)模塊：給出養(yǎng)殖容量和環(huán)境容量的估算及海水質(zhì)量的判別，具備水體富營(yíng)養(yǎng)狀況評(píng)價(jià)功能。

圖1 系統(tǒng)化生態(tài)－沉積－環(huán)境動(dòng)力學(xué)耦合模式總體構(gòu)架Fig.1 Framework of the systematic ecology－sediment－environment dynamical numerical coupled model

圖2 模塊間的耦合過(guò)程圖Fig.2 Sketch of coupling relationship of modules

2 各模塊的理論基礎(chǔ)

2.1 水動(dòng)力模式

2.1.1 環(huán)流模式

作為各模塊載體的水動(dòng)力模式采用國(guó)家海洋局第一海洋研究所基于POM（Princeton Ocean Model）發(fā)展的浪－潮－流耦合模式，以三維流速、溫度、鹽度、湍動(dòng)能和混合長(zhǎng)為預(yù)報(bào)變量，通過(guò)嵌入湍封閉模型對(duì)垂向混合過(guò)程進(jìn)行參數(shù)化。該模式為用FORTRAN語(yǔ)言編寫的開源模式，采用了曲線正交網(wǎng)格，接口能力強(qiáng)，便于今后進(jìn)一步的開發(fā)完善和模塊嵌入。各模塊的時(shí)間和空間積分及網(wǎng)格均與水動(dòng)力模式一致。

2.1.2 海浪的耦合

海浪的作用（圖3）包括以下方面：

1）將波致混合Bv［2］引入到泥沙擴(kuò)散方程的湍混合系數(shù)KH之中，合理地增強(qiáng)了層泥沙混合；

2）在波流耦合的底邊界層模型［3－4］的基礎(chǔ)上，加入海底沙波的計(jì)算；

底邊界模型中，波浪邊界層和物理底粗糙度聯(lián)合引起的表觀粗糙度z0c為

式中，u＊cw為波流耦合的摩擦速度；Ab為波浪近底振幅；β為波所占的比例系數(shù)其中u＊c為流摩擦速度；ub為波浪軌道速度；w為波頻率。z0c經(jīng)過(guò)波浪修正，用z0c得到新的拖曳系數(shù)cd為

式中，zr為參考高度；卡曼常數(shù)κ取0.4。cd用以計(jì)算流致底切應(yīng)力，進(jìn)而得到波流耦合的底切應(yīng)力。

在此底邊界模型基礎(chǔ)上考慮海底沙波的存在，波浪剪切摩擦速度穩(wěn)定流的摩擦速度之比，u＊ws／u＊cs，可用來(lái)定義聯(lián)合流作用下不同種類的沙波。當(dāng)u＊ws／u＊cs＞1.25時(shí)，波浪起主導(dǎo)作用。而當(dāng)u＊ws／u＊cs＜1.25時(shí)，應(yīng)用波流聯(lián)合流用來(lái)預(yù)測(cè)沙波的幾何形態(tài)。

（1）波浪起主導(dǎo)作用下海底沙波模型［5］

沙波波高η用近底偏移振幅Ab和遷移數(shù)φ來(lái)計(jì)算：

φ與Shields夾帶參數(shù)有關(guān)，且φ由波浪近底軌道速度決定，當(dāng)φ＜10時(shí)沙波的波高為不規(guī)則波。沙波的波長(zhǎng)λr由波陡計(jì)算。

（2）波流聯(lián)合作用下的海底沙波模型

分為3種狀態(tài)：

①弱傳輸狀態(tài)，u＊cws＜u＊cr

②平衡狀態(tài)，u＊cws≥u＊cr且u＊cwb＜u＊bf

③破碎狀態(tài)，u＊cwb≥u＊bf

其中，u＊cws為僅考慮泥沙粒徑時(shí)波流聯(lián)合作用的摩擦速度；u＊cr為推移質(zhì)運(yùn)動(dòng)臨界摩擦速度；u＊bf為破碎狀態(tài)臨界摩擦速度。分別計(jì)算出沙波波高η，沙波波陡η／λ。沙波增強(qiáng)的Nikuradse粗糙度ks由式（4）計(jì)算。其中，θs為Sheilds參數(shù)，由波流共同作用的底切應(yīng)力決定。ks用于更新波致底切應(yīng)力的波浪摩擦因子fw。

3）引入由波致孔隙壓力建立的再懸浮通量E1，考慮由此導(dǎo)致的底沙液化［6］。

其中，HS為有效波高；A1為只取決于海底細(xì)沙特征的經(jīng)驗(yàn)常數(shù)；F是波頻率ω、水深H的函數(shù)。

圖3 海浪的作用框架圖Fig.3 Sketch of wave＇s effect

2.2 沉積模塊

沉積模塊融合了POMSED和ECOMSED的優(yōu)勢(shì)。參考目前國(guó)際上較先進(jìn)的沉積物輸運(yùn)模式POMSED［7－8］，通過(guò) Richardson數(shù)引入底邊界層化效應(yīng)；通過(guò)嵌入 ECOMSED波流耦合底邊界層模型［3－4］，考慮了波致底沙液化［6］對(duì)底沙的侵蝕作用；改進(jìn)了泥沙濃度輻射邊界條件；考慮了泥沙密度對(duì)水動(dòng)力的作用，可模擬異重流等情況。

2.2.1 泥沙對(duì)水動(dòng)力的反作用

主要是考慮泥沙對(duì)水體密度的影響。根據(jù)霧狀層動(dòng)力學(xué)，考慮懸浮泥沙的貢獻(xiàn)，可模擬異重流等情況。海水的密度用體積關(guān)系進(jìn)行計(jì)算［10］：

式中，ρw是純水的密度；ρs是泥沙的密度；C是懸浮泥沙濃度。

2.2.2 底邊界有機(jī)顆粒物的沉積和再懸浮

海底沉積物里面包含碎屑顆粒物，其濃度擴(kuò)散方程中考慮了和波浪混合和底部的再懸浮及沉積過(guò)程，沉積通量SD為

式中，WS為沉降速度；P為沉降概率；D為碎屑顆粒物濃度。

總再懸浮通量E為

式中，fD為有機(jī)碎屑所占分?jǐn)?shù)；E1為由波致孔隙壓力建立的再懸浮通量；E2為底部剪切引起的再懸浮通量，，其中τb為波流聯(lián)合底切應(yīng)力，τce為啟動(dòng)臨界應(yīng)力。

2.2.3 泥沙對(duì)生態(tài)的光效應(yīng)

將懸浮物濃度模擬結(jié)果耦合到生態(tài)模塊浮游植物光效應(yīng)的計(jì)算中。考慮懸浮物濃度對(duì)可見(jiàn)光的阻隔作用，浮游植物生長(zhǎng)的光限制因子：

式中，I為有效光照強(qiáng)度，采用國(guó)際上最新且通過(guò)觀測(cè)進(jìn)行了驗(yàn)證的R－修正雙指數(shù)參數(shù)化公式［11］：

其中，I0為有效光合作用的海表太陽(yáng)輻射；R′、R是表層和水體分配常數(shù)；k1、k2是表層和水體的垂向衰減系數(shù)，光衰減系數(shù)k（即k1和k2）隨SSC的變化為

其中，kew為水自身衰減系數(shù)；P，D和SSC分別為浮游植物、碎屑和懸浮物質(zhì)量濃度。

2.3 生態(tài)模塊

生態(tài)模塊主要模擬浮游生態(tài)系統(tǒng)及與之相關(guān)的關(guān)鍵環(huán)境因子的時(shí)空分布，大體可將生物化學(xué)要素分為4類：浮游植物（P）、浮游動(dòng)物（Z）、營(yíng)養(yǎng)鹽（N）、有機(jī)碎屑（D）。浮游植物為浮游植物生物量，并可分為多個(gè)藻種?？紤]中國(guó)近海藻類生物種類和數(shù)量情況，主要給出甲藻和硅藻兩類，對(duì)于較少的藻類如藍(lán)藻、金藻、綠藻等暫不添加；浮游動(dòng)物個(gè)體差異很大，考慮目前生態(tài)模型主要用于研究浮游植物和營(yíng)養(yǎng)鹽，暫將浮游動(dòng)物作為一個(gè)功能群考慮，以總生物量表示，若要研究某一個(gè)或幾個(gè)種類，如單細(xì)胞的原生動(dòng)物（鞭毛類、有孔蟲類、纖毛蟲類）、浮游甲殼動(dòng)物、水母類、毛顎動(dòng)物、幼形類及其他種類，可針對(duì)其獨(dú)特的攝食和代謝特點(diǎn)建立模型；營(yíng)養(yǎng)鹽包括總無(wú)機(jī)氮（可分為氨氮、亞硝酸態(tài)氮、硝酸態(tài)氮）、磷酸鹽和硅酸鹽；有機(jī)碎屑主要包括浮游生物尸體和糞便等。用戶可根據(jù)需要，在生態(tài)模塊中添加游泳動(dòng)物。

生態(tài)要素之間關(guān)鍵的生物化學(xué)轉(zhuǎn)化過(guò)程表示為

其中，浮游植物生長(zhǎng)、代謝和死亡分別為

浮游動(dòng)物的攝食、生長(zhǎng)、代謝和死亡分別為

有機(jī)碎屑的分解為

各式中，N，DIP和Si分別為總無(wú)機(jī)氮、磷酸鹽濃度和硅酸鹽濃度；KN，KP和KSi分別為無(wú)機(jī)氮、磷酸鹽和硅酸鹽半飽和濃度常數(shù)；gPm和mPm分別為浮游植物0℃時(shí)的最大生長(zhǎng)率和呼吸率；gZm和mZm分別為浮游動(dòng)物0℃最大生長(zhǎng)率和代謝率；μP和γP分別為浮游植物生長(zhǎng)率和呼吸率隨溫度變化的系數(shù)；μZ和γZ分別為浮游動(dòng)物生長(zhǎng)率和代謝率隨溫度變化的系數(shù)；dP為浮游植物自然死亡率；PK為浮游植物的半飽和死亡率；β和dZ分別為浮游動(dòng)物同化率和自然死亡率；λ為Ivlev常數(shù)；e′為有機(jī)碎屑分解率；Ws為浮游植物或碎屑的沉降系數(shù)；Source為包括陸源輸入、大氣干濕沉降和沉積物礦化溶入的營(yíng)養(yǎng)鹽外源輸入。

模型還考慮了高等動(dòng)物和大型藻對(duì)生態(tài)要素的貢獻(xiàn)和消耗，高等動(dòng)物包括游泳動(dòng)物和養(yǎng)殖動(dòng)物。Egst＿Hz，Eat＿Hz，Exc＿Hz分別為高等級(jí)動(dòng)物對(duì)浮游植物的攝食、對(duì)浮游動(dòng)物的攝食、代謝過(guò)程中對(duì)營(yíng)養(yǎng)鹽的釋放；Grow＿La為大型藻對(duì)營(yíng)養(yǎng)鹽的吸收；Mort＿LaHz表示高等動(dòng)物及大型藻生長(zhǎng)死亡過(guò)程對(duì)碎屑的釋放。針對(duì)養(yǎng)殖動(dòng)物（如貝類）和養(yǎng)殖大型藻（如海帶）涉及的生物化學(xué)過(guò)程設(shè)置在養(yǎng)殖模塊中。

2.4 精細(xì)化水質(zhì)模塊

區(qū)別于生態(tài)模塊，本模塊主要是通過(guò)模擬化學(xué)要素參量來(lái)研究水質(zhì)情況。除生態(tài)模塊中的生物化學(xué)過(guò)程外，還包括以下過(guò)程：

1）微生物（細(xì)菌）生態(tài)過(guò)程

其生長(zhǎng)方程為

其中，B為微生物濃度；Grow＿B，Rsp＿B，Egst＿Z，Rs＿B分別為細(xì)菌的生長(zhǎng)（伴隨著碎屑分解過(guò)程）、代謝、死亡、被動(dòng)物的攝食及底部進(jìn)入。

2）底棲生態(tài)過(guò)程

主要變量包括大型底棲動(dòng)物、小型底棲動(dòng)物、底棲細(xì)菌、有機(jī)碎屑、間隙水營(yíng)養(yǎng)鹽等。有機(jī)碎屑主要來(lái)源于水體中顆粒有機(jī)物（碎屑、浮游植物）的沉降及底棲生物的死亡，間隙水營(yíng)養(yǎng)鹽通過(guò)細(xì)菌分解和礦化得以再生并溶入水體。

3）耗氧或釋氧過(guò)程

溶解氧濃度源匯方程為

其中，O為溶解氧濃度；GrowP＿O為浮游植物光合作用釋放氧氣、Air＿O為海氣界面的曝氣還氧、NitN＿O為營(yíng)養(yǎng)鹽的氨化作用產(chǎn)生的氧氣；RspZ＿O，RspP＿O，RspB＿O分別為浮游動(dòng)物、浮游植物呼吸作用和細(xì)菌分解有機(jī)物所消耗的氧氣；COD＿O為污染排放引起的有機(jī)物氧化所耗的氧；OxN＿O為硝化耗氧；Botm＿O為底泥耗氧（涉及生物呼吸、硝化、沉積物氧化分解）。

4）溶解無(wú)機(jī)物、有機(jī)碎屑與溶解無(wú)機(jī)物之間的轉(zhuǎn)化過(guò)程

有機(jī)碎屑（POM）被細(xì)菌分解轉(zhuǎn)化為溶解有機(jī)物（DOM）和溶解無(wú)機(jī)物（DIM），DOM 礦化為DIM，POM、DOM又來(lái)源于浮游動(dòng)植物和細(xì)菌的代謝（分泌、排泄、死亡）及外源輸入。各要素通過(guò)碳、氮、磷元素的含量可表達(dá)為對(duì)應(yīng)的變量，如碳形式的顆粒有機(jī)碳（POC）、溶解無(wú)機(jī)碳（DOC）、溶解無(wú)機(jī)碳（DIC）。

5）營(yíng)養(yǎng)鹽之間的硝化、氨化過(guò)程

2.5 養(yǎng)殖模塊

本模塊設(shè)定養(yǎng)殖區(qū)和水層、養(yǎng)殖種類后，將能模擬養(yǎng)殖對(duì)水動(dòng)力的作用及養(yǎng)殖對(duì)生物化學(xué)過(guò)程的作用。

1）養(yǎng)殖對(duì)水動(dòng)力的作用

通過(guò)在運(yùn)動(dòng)方程中加入養(yǎng)殖引起的動(dòng)量損失實(shí)現(xiàn)。動(dòng)量損失的參數(shù)化方案主要參考已有工作［12］并結(jié)合近年積累的大量觀測(cè)資料進(jìn)行完善，分別對(duì)不同的養(yǎng)殖種類和密度給出合適的動(dòng)量損失系數(shù)（或摩擦系數(shù)）。模型主要考慮采用懸繩與吊籠養(yǎng)殖的海帶和貝類，及采用網(wǎng)箱養(yǎng)殖的魚類。

2）養(yǎng)殖對(duì)生物化學(xué)過(guò)程的影響

主要考慮養(yǎng)殖生物對(duì)生態(tài)要素的吸收、釋放及對(duì)水質(zhì)的影響等。

海帶的影響：包含養(yǎng)殖期間對(duì)營(yíng)養(yǎng)鹽的凈吸收、成熟收割期海帶腐爛對(duì)營(yíng)養(yǎng)鹽或有機(jī)物的釋放，及對(duì)氧氣的釋放或消耗。涉及的有關(guān)參數(shù)如干濕重比、生長(zhǎng)率、C／N／P的摩爾比、養(yǎng)殖密度、光合作用放氧量等可根據(jù)有關(guān)文獻(xiàn)［11，12－15］或特定海域的測(cè)量來(lái)確定。

貝類的影響：主要包括中上層水體和底播養(yǎng)殖的濾食性雙殼貝類。典型的貝類有吊籠養(yǎng)殖的扇貝，懸繩養(yǎng)殖的牡蠣、貽貝，底播養(yǎng)殖的菲律賓蛤仔、縊蟶等。貝類對(duì)浮游植物和有機(jī)碎屑的攝食率及代謝過(guò)程中對(duì)營(yíng)養(yǎng)鹽、溶解有機(jī)物、有機(jī)碎屑的釋放率和對(duì)氧氣的消耗等參考相關(guān)文獻(xiàn)［12，15－19］。另外，吊籠養(yǎng)殖還要考慮濾食性附著動(dòng)物的攝食和排泄［14，20］。

魚類的影響：主要考慮養(yǎng)殖魚類對(duì)浮游動(dòng)物的攝食、排泄出的碎屑及對(duì)氧氣的消耗等，有關(guān)攝食率和同化率參考有關(guān)文獻(xiàn)［12，19］。

設(shè)定養(yǎng)殖區(qū)的養(yǎng)殖種類和規(guī)模后，本模塊將計(jì)算出對(duì)生態(tài)要素（含關(guān)鍵水質(zhì)變量）的吸收率或釋放率。

2.6 物質(zhì)擴(kuò)散及粒子追蹤模塊

該模塊主要計(jì)算常用的、典型性的物質(zhì)或環(huán)境參數(shù)的擴(kuò)散及運(yùn)移情況。包括對(duì)溫排水、濃鹽水、污染物排放的模擬，可給出濃度或溫度的最大包絡(luò)、特定位置的溫升（或濃度升高）及擴(kuò)散范圍，以及粒子運(yùn)移路徑。

物質(zhì)擴(kuò)散方程為

其中，C代表溫度、鹽度或污染物濃度；t為時(shí)間；u，v，w分別為x，y，z方向的流速分量；Ah和Av分別為水平和垂向湍混合系數(shù)；R為外源輸入；λ為衰減系數(shù)。對(duì)溫排水，考慮蒸發(fā)和風(fēng)引起的散熱，水面綜合散熱系數(shù)采用Gunnerberg經(jīng)驗(yàn)公式。對(duì)污染物擴(kuò)散，考慮其濃度衰減，將對(duì)不同的污染物給定不同的衰減系數(shù)，可以為一個(gè)常數(shù)也可為經(jīng)驗(yàn)公式。

粒子追蹤方程采用了拉格朗日粒子追蹤模型

其中，x為坐標(biāo)，→v（x，t）為隨空間和時(shí)間變化的速度矢量。另外，將生態(tài)模塊中的大型藻生長(zhǎng)模型和粒子追蹤結(jié)合起來(lái)可模擬大型藻的生長(zhǎng)和漂流軌跡。

2.7 評(píng)估評(píng)價(jià)模塊

2.7.1 環(huán)境容量

通過(guò)物質(zhì)擴(kuò)散模塊、生態(tài)模塊或水質(zhì)模塊計(jì)算某一污染物源強(qiáng)下形成的濃度場(chǎng)，確定水質(zhì)控制點(diǎn)和污染源之間的濃度響應(yīng)系數(shù)場(chǎng)和污染分擔(dān)率，按分擔(dān)率法和排海通量最優(yōu)法估算環(huán)境容量。

采用分擔(dān)率法計(jì)算的污染源的入海最大負(fù)荷為

其中，F(xiàn)i為某一個(gè)污染源的排海通量；Ci（x，y，z）為該污染源單獨(dú)作用造成的目標(biāo)海域的污染濃度為該污染源單獨(dú)作用而使目標(biāo)海域濃度處于一定等級(jí)水質(zhì)標(biāo)準(zhǔn)的濃度。各污染源入海最大負(fù)荷之和為研究海域的總環(huán)境容量。

為優(yōu)化環(huán)境容量，使各污染源的允許入海負(fù)荷分配之和達(dá)到最大且各水質(zhì)控制點(diǎn)仍滿足一定等級(jí)水質(zhì)標(biāo)準(zhǔn)要求，水質(zhì)標(biāo)準(zhǔn)約束條件為

式中，為水質(zhì)控制點(diǎn)處污染物背景濃度為控制點(diǎn)處的水質(zhì)標(biāo)準(zhǔn)濃度；m為水質(zhì)控制點(diǎn)數(shù)目；n為污染源數(shù)；αij為第i個(gè)污染源對(duì)第j個(gè)水質(zhì)控制點(diǎn)的分擔(dān)率。求解方程所界定的線性規(guī)劃問(wèn)題，得出污染物最大允許排放量。

2.7.2 養(yǎng)殖容量

由生態(tài)模塊得到浮游植物生物量和初級(jí)生產(chǎn)力，并基于營(yíng)養(yǎng)收支平衡原理或生態(tài)系統(tǒng)中的營(yíng)養(yǎng)級(jí)關(guān)系研究養(yǎng)殖容量。

對(duì)貝類的養(yǎng)殖容量，基于餌料收支平衡關(guān)系，計(jì)算濾食性貝類的單位面積養(yǎng)殖容量：

其中，CC為養(yǎng)殖區(qū)內(nèi)某一養(yǎng)殖貝類的平均養(yǎng)殖容量（ind／m2）；PP為養(yǎng)殖區(qū)內(nèi)平均初級(jí)生產(chǎn)量；PPRZ，PPRaddit，PPRout分別為養(yǎng)殖區(qū)平均單位面積1d內(nèi)的浮游動(dòng)物攝食量、附著生物的濾食量、及餌料的流出（或流入）量，IR為養(yǎng)殖區(qū)內(nèi)某一養(yǎng)殖貝類的濾食率；PPD、PPB分別為養(yǎng)殖區(qū)內(nèi)平均碎屑中的有機(jī)碳的供應(yīng)、底部微藻的初級(jí)生產(chǎn)供應(yīng)。

對(duì)大型藻的養(yǎng)殖容量，其養(yǎng)殖容量的計(jì)算公式

其中，TNK為該海域單位面積可供大型養(yǎng)殖藻生長(zhǎng)的某種營(yíng)養(yǎng)鹽的總量；PCK為大型藻類的該種營(yíng)養(yǎng)鹽的含量。TNK涉及通過(guò)海水交換、陸源輸入、大氣沉降釋放、海洋動(dòng)物排泄物及海底沉積物釋放的營(yíng)養(yǎng)鹽，及浮游植物、野生大型藻、潮灘微型藻類的生長(zhǎng)和海洋動(dòng)物通過(guò)食物鏈傳遞吸收的營(yíng)養(yǎng)鹽［12］。

對(duì)魚類的養(yǎng)殖容量，基于營(yíng)養(yǎng)動(dòng)態(tài)模型［21］估算魚類營(yíng)養(yǎng)階層的生產(chǎn)量

式中，P為一個(gè)區(qū)域內(nèi)增養(yǎng)或放養(yǎng)的魚類容量，以鮮重計(jì)，單位為t／a；B為魚類所需要的初級(jí)生產(chǎn)量，以浮游植物的濕重計(jì)；E為生態(tài)效率；n為魚類的營(yíng)養(yǎng)階層（營(yíng)養(yǎng)級(jí)轉(zhuǎn)換級(jí)數(shù)）。3個(gè)重要參數(shù)的具體計(jì)算方法參考文獻(xiàn)［12］。對(duì)于網(wǎng)箱魚類養(yǎng)殖，可以通過(guò)計(jì)算養(yǎng)殖水域網(wǎng)箱內(nèi)魚類生長(zhǎng)的最大餌料供應(yīng)力（包括人工添加餌料），確定研究海域的魚類增養(yǎng)能力。對(duì)于放養(yǎng)的魚類，可以通過(guò)計(jì)算餌料的總生產(chǎn)量并減去養(yǎng)殖貝類及其他動(dòng)物對(duì)浮游植物的攝食，以得到放養(yǎng)魚類的最大供餌力。

3 模式的初步結(jié)果

目前的耦合模型已在水動(dòng)力模型的基礎(chǔ)上完成了沉積物輸運(yùn)模塊、生態(tài)模塊、養(yǎng)殖模塊、物質(zhì)擴(kuò)散模塊、評(píng)估評(píng)價(jià)模塊幾個(gè)獨(dú)立模塊的構(gòu)建，并已實(shí)現(xiàn)了波浪－環(huán)流的耦合、生態(tài)和養(yǎng)殖模塊的耦合、養(yǎng)殖和評(píng)估評(píng)價(jià)模塊的耦合。以下示例給出幾個(gè)模擬結(jié)果。

3.1 浮游生態(tài)過(guò)程和赤潮生消的模擬

考慮浪－流－潮的物理作用，建立了南黃海浮游動(dòng)力學(xué)模型，模擬了南黃海的浮游生態(tài)過(guò)程，刻畫出了浮游植物生物量的時(shí)、空變化規(guī)律，揭示了春華和次表層最大化等典型生態(tài)特征及其受控的環(huán)境動(dòng)力機(jī)制。圖4為模擬的典型月份葉綠素質(zhì)量濃度的表層和垂直分布。

圖4 南黃海典型月份的葉綠素a質(zhì)量濃度的平面和35°N斷面垂直分布Fig.4 Chl－a concentration in the South Yellow Sea for typical months in horizontal distribution and vertical along 35°N section

在生態(tài)模型的基礎(chǔ)上，建立了赤潮模型，并按照假定的環(huán)境條件模擬的青島奧帆賽場(chǎng)及臨近海域的赤潮生消過(guò)程，如圖5所示，圖中a～i為生消過(guò)程中葉綠素a質(zhì)量濃度的變化。

3.2 養(yǎng)殖海域生態(tài)模型及容量評(píng)估結(jié)果

對(duì)典型養(yǎng)殖海域桑溝灣，考慮養(yǎng)殖對(duì)水動(dòng)力的影響、養(yǎng)殖生物參與的生物化學(xué)過(guò)程，建立了該灣的生態(tài)動(dòng)力學(xué)模型，并估算了主要養(yǎng)殖生物的養(yǎng)殖容量，圖6給出了模型的理論框架、全灣平均的葉綠素a質(zhì)量濃度變化及與實(shí)測(cè)的對(duì)比、全年不同時(shí)間、不同規(guī)格扇貝的養(yǎng)殖容量。

圖5 青島近岸海域赤潮生消過(guò)程控制實(shí)驗(yàn)（改自文獻(xiàn)［22］）Fig.5 Numerical experiment of evolvement of HAB along the cosat of Qingdao（from modifying the Number 22literature）

圖6 養(yǎng)殖海域生態(tài)模型及桑溝灣葉綠素a質(zhì)量濃度與扇貝養(yǎng)殖容量模型Fig.6 Ecosystem model for aquaculture areas and model results for Chl－a mass concentration and carrying capacity of shellfish of Sanggou Bay

3.3 渤黃海沉積動(dòng)力學(xué)模擬結(jié)果

圖7給出2006－08—2009－04期間10個(gè)航次中渤海底層懸浮物調(diào)查結(jié)果的平均值及以氣候態(tài)模擬的渤海底層懸浮物濃度的年平均值。對(duì)比可見(jiàn)，模擬與觀測(cè)在分布上大體一致。懸浮物高質(zhì)量濃度區(qū)域主要在渤海灣南部至黃河口區(qū)域，此外，遼東灣灣頂、遼東半島西側(cè)及山東半島東端也有高質(zhì)量濃度懸浮泥沙分布。這些懸浮物高值區(qū)域與底切應(yīng)力高值區(qū)域?qū)?yīng)，因此這些區(qū)域水動(dòng)力作用較強(qiáng)。

圖8給出了年沉積厚度沖淤的模擬結(jié)果分布。沉積厚度集中的區(qū)域主要集中在萊州灣，其次在北黃海西北部。發(fā)生沖刷的主要區(qū)域在渤海灣、遼東半島西側(cè)及朝鮮半島西部。此外，還模擬了現(xiàn)代黃河的沉積厚度分布，結(jié)果顯示，在山東半島東南端水下三角洲頂積層，現(xiàn)代黃河沉積速率模擬結(jié)果約為0.1cm／a，這與Alexander［23］基于210Pb的調(diào)查結(jié)果（0.1～0.2cm／a）在量級(jí)上一致。

圖7 觀測(cè)與模擬的底層懸浮物質(zhì)量濃度對(duì)比Fig.7 Validation of suspended sediment concentration of bottom layer

圖8 模擬的年沉積速率Fig.8 Simulated deposition rates

4 結(jié)論

本文闡述了系統(tǒng)化的海洋生態(tài)－沉積－環(huán)境動(dòng)力學(xué)耦合數(shù)值模式的理論基礎(chǔ)，給出了生態(tài)模塊、沉積物輸運(yùn)模塊、精細(xì)化水質(zhì)模塊、養(yǎng)殖模塊、物質(zhì)擴(kuò)散及粒子追蹤模塊、評(píng)估評(píng)價(jià)模塊之間的耦合方式，給出了作為各模塊載體的水動(dòng)力模式的浪流耦合方案，以及功能模塊的表達(dá)方程和參數(shù)化方案或?qū)崿F(xiàn)的計(jì)算方式。對(duì)各功能模塊之間的耦合重點(diǎn)考慮了：沉積動(dòng)力過(guò)程對(duì)生物化學(xué)過(guò)程的作用，包括懸沙導(dǎo)致的光衰減、及考慮層化和波流耦合的底邊界層中有機(jī)顆粒的沉積和再懸?。缓＠说鸟詈献饔?，包括上層的波致垂向混合、波流耦合的底邊界層過(guò)程、三維輻射應(yīng)力作用下的波致流沿岸物質(zhì)輸運(yùn)、及波致底沙液化等；泥沙活動(dòng)和養(yǎng)殖設(shè)施對(duì)水動(dòng)力的作用。

本研究已初步建立了主要的功能模塊，實(shí)現(xiàn)了一些模塊的耦合，并用于中國(guó)近海海域得到了較好的模擬結(jié)果。下一步我們將主要建立并完善精細(xì)化水質(zhì)模塊和物質(zhì)擴(kuò)散及粒子追蹤模塊，完成各模塊之間的有效耦合。另外還著重于泥沙過(guò)程和生態(tài)過(guò)程的銜接，并在生態(tài)模塊實(shí)現(xiàn)碳循環(huán)的模擬功能。最后在完成模塊組合和模型驗(yàn)證的基礎(chǔ)上，給出標(biāo)準(zhǔn)化的系統(tǒng)模式，寫出模式系統(tǒng)介紹和使用說(shuō)明書。發(fā)展的模式以期能方便、可靠地用于海洋生態(tài)系統(tǒng)動(dòng)力過(guò)程、赤潮生消過(guò)程、懸浮物濃度、海底沉積物沖淤、海洋碳循環(huán)、水質(zhì)、物質(zhì)擴(kuò)散、環(huán)境容量、養(yǎng)殖容量等方面的研究和評(píng)估評(píng)價(jià)。

（References）：

［1］SHU T F，LUO L，WEN Y M.Effects of mariculture on coastal ecological environment［J］.Marine Environmental Science，2002，21（2）：74－79.舒廷飛，羅琳，溫琰茂.海水養(yǎng)殖對(duì)近岸生態(tài)環(huán)境的影響［J］.海洋環(huán)境科學(xué)，2002，21（2）：74－79.

［2］QIAO F L，YUAN Y L，YANG Y Z，et al.Wave－induced mixing in the upper ocean：Distribution and application to a globalocean circulation model［J］.Geophysical Research Letters，2004，31（11）：doi：10.1029／2004GL019824

［3］GRANT W D，MADSEN O S.Combined wave and current interaction with a rough bottom［J］.Journal of Geophisical Research，1979，84（C4）：1797－808.

［4］SIGNELL R P，BEARDSLEY R C，GRABER H C，et al.Effect of wave－current interaction on steady wind－driven circulation in narrow，shallow embayments［J］.Journal of Geophisical Research，1990，95（C6）：9671－9678.

［5］SMITH G A，BABANIN A V，RIEDEL P，et al.Introduction of a new friction routine into the SWAN model that evaluates roughness due to bedform and sediment size changes［J］.Coastal Engineering，2011，58（4）：317－326.

［6］LI M Z，AMOS C L.Predicting ripple geometry and bed roughness under combined waves and currents in a continental shelf environment［J］.Continental Shelf Research，1998，18（9）：941－970.

［7］LAMBRECHTS J，HUMPHREY C，MCKINNA L，et al.Importance of wave－induced bed liquefaction in the fine sediment budget of Cleveland Bay，Great Barrier Reef［J］.Estuarine，Coastal and Shelf Science，2010，89（2）：154－162.

［8］WANG X H，PINARDI N.Modeling the dynamics of sediment transport and resuspension in the Northern Adriatic Sea［J］.Journal of Geophysical Research，2002，107（C12）：doi：10.1029／2001JC001303.

［9］WANG X H，PINARDI N，MALACIC V.Sediment transport and resuspension due to combined motion of wave and current in the northern Adriatic Sea during a Bora event in January 2001：a numerical modeling study［J］.Continental Shelf Research，2007，27（5）：613－633.

［10］WANG X H.Tide－induced sediment resuspension and the bottom boundary layer in an idealized estuary［J］.Journal of Physical Oceanography，2002，32：3113－3131.

［11］BYUN D S，WANG X H，HARTC D E，et al.Review of parameterizations in ocean ecosystem models［J］.Estuarine，Coastal and Shelf Science，2014，151：318－323.

［12］LIU X H.Ecosystem dynamic models for South Yellow Sea and an aquaculture bay［D］.Qingdao：Ocean University of China，2009.劉學(xué)海.南黃海及養(yǎng)殖功能海域生態(tài)動(dòng)力學(xué)模型研究［D］.青島：中國(guó)海洋大學(xué)，2009.

［13］FANG J G，KUANG S H，SUN H L.Accessing the carrying capacity of Sanggou Bay for culture of kelp Laminaria Japonica［J］.Marine Fisheries Research，1996，17（2）：7－17.方建光，匡世煥，孫慧玲.桑溝灣海帶養(yǎng)殖容量的研究［J］.海洋水產(chǎn)研究，1996，17（2）：7－17.

［14］FANG J G，KUANG S H，SUN H L.Study on the carrying capacity of Sanggou Bay for the culture of scallop［J］.Marine Fisheries Research，1996，17（2）：18－31.方建光，匡世煥，孫慧玲.桑溝灣櫛孔扇貝養(yǎng)殖容量的研究［J］.海洋水產(chǎn)研究，1996，17（2）：18－31.

［15］ZHU M Y，ZHANG X L，TANG Y T.Application of ecological model in studing the sustainable development of coastal shellfish culture［J］.Advances in Marine Science，2002，20（4）：34－42.朱明遠(yuǎn)，張學(xué)雷，湯庭耀.應(yīng)用生態(tài)模型研究近海貝類養(yǎng)殖的可持續(xù)發(fā)展［J］.海洋科學(xué)進(jìn)展，2002，20（4）：34－42.

［16］NUNES J P，F(xiàn)ERREIRA J G，GAZEAU F.A model for sustainable management of shellfish polyculture in coastal bays［J］.Aquaculture，2003，219：257－277.

［17］KUANG S H，SUN H L，LI F.Feeding and growth of scallop Chlanys Farreri begore and after spawning［J］.Marine Fisheries Research，1996，17（2）：80－86.匡世煥，孫慧玲，李鋒.櫛孔扇貝生殖活動(dòng)前后的濾食和生長(zhǎng)［J］.海洋水產(chǎn)研究，1996，17（2）：80－86.

［18］SONG X L，SUN Y，WU J X.Effects of temperature and shell length on the excretion of nitrogen and phoshorus of Chlamys Farreri［J］.Journal of Fisheies of China，2007，31（1）：38－44.宋嫻麗，孫耀，武晉宣.櫛孔扇貝的氮、磷排泄及其與溫度、殼長(zhǎng)的關(guān)系［J］.水產(chǎn)學(xué)報(bào)，2007，31（1）：38－44.

［19］WU J X.Development and applications of dynamic model for dissovel introgen and phosphate budgets in Sungo Bay［D］.Qingdao：Ocean University of China，2005.武晉宣.桑溝灣養(yǎng)殖海域氮、磷收支及環(huán)境容量模型［D］.青島：中國(guó)海洋大學(xué)，2005.

［20］SUN H L，F(xiàn)ANG J G，KUANG S H.Testing of the filtation rate of Ascidiceain different seasons［J］.Marine Fisheries Research，1996，17（2）：103－107.孫慧玲，方建光，匡世煥.不同季節(jié)海鞘濾水率的測(cè)定［J］.海洋水產(chǎn)研究，1996，17（2）：103－107.

［21］PARSONS T R，TAKAHASHI M，HARGRAVE B.Biological Oceanography Processes［J］.Oxford：Pergamon Press，1973.

［22］LIU X H，LI R X，YUAN Y L.Numerical experiments of HAB inducement in Qingdao adjacent coastal area［J］.Acta Oceanologica Sinica，2010，29（5）：109－120.

［23］ALEXANDERCR，DEMASTER D J，NITTROUER C A.Sediment accumulation in a modern epicontinental－shelf setting：the Yellow Sea［J］.Marine Geology，1991，98（1）：51－72.