（水声工程专业论文）耦合简正波快速算法研究.pdf

a b s t r a c t i nt h el i g h to ft h eh i g hp r e c i s i o no ft h ec o u p l e dn o r m a lm o d e ( c n m ) ，i ti s w i d e l ya d o p t e d a st h ec r i t e r i o no ft h ea c o u s t i cf i e l dc o m p u t a t i o ni nc o m p l e xo c e a n e n v i r o n m e n t s h o w e v e r , t h e r ee x i s t s ag r e a tm o u n to fc o m p u t a t i o n s t h e c o m p u t i n gs p e e di sn o ta c c e p t a b l ef o rp r a c t i c e ，e i t h e r t h i sp a p e rf o c u so nt h e p r o b l e mm e n t i o n e da b o v ev i ar e s e a r c h i n gt h ef a s ta l g o r i t h m o fs o l v i n gt h ec n m a m p l i t u d ee q u a t i o n t h ep a p e rs t a r t sf r o mi n t r o d u c i n gt h ec n m ，p a r a b o l i ce q u a t i o nm e t h o d ( p e l v oa n dr a ya c o u s t i c sm e t h o d ( r a m ) t h e nt h ed i f f e r e n t i a le q u a t i o n sf o r s o l v i n gc n ma m p l i t u d ea r eo b t a i n e d t h ec n m r a mm o d e la n dc n m p e m m o d e lf o rh i g hs p e e dn o r m a lm o d es o l u t i o na r ed e s c r i b e di nd e t a i l b e s i d e s ，t h e p r o g r a mf o rc n m p e mn u m e r i cc o m p u t a t i o n i s d e v e l o p e d o nt h ef o r t r a n p l a t f o r m i t sc o m p u t a t i o ne r r o r so f t r a n s i t i o nl o s si nc o n d i t i o n so fr a n g ed e p e n d e n t a n di n d e p e n d e n ta r ev e r i f i e dt h r o u g hc o m p a r i s o nw i t hb e n c h m a r kp r o b l e m r e a s o n sf o re r r o r sa r ea n a l y z e dv i ac o m p a r i n gw i t ht h ec o u p l e 0 7r e s u l t s k e yw o r d s ：k r a k e n ；c n m p e m ；a c o u s t i cp r o p a g a t i o n ；m o d e l i n g 哈尔滨工程大学 学位论文原创性声明 本人郑重声明：本论文的所有工作，是在导师的指导下， 由作者本人独立完成的。有关观点、方法、数据和文献的引用 已在文中指出，并与参考文献相对应。除文中已注明引用的内 容外，本论文不包含任何其他个人或集体已经公开发表的作品 成果。对本文的研究做出重要贡献的个人和集体，均己在文中 以明确方式标明。本人完全意识到本声明的法律结果由本人承 担。 作者( 签字) ： 芟厦 日期： 功d7 年岁月刃日 哈尔滨工程大学 学位论文授权使用声明 本人完全了解学校保护知识产权的有关规定，即研究生在校 攻读学位期间论文工作的知识产权属于哈尔滨工程大学。哈尔滨 工程大学有权保留并向国家有关部门或机构送交论文的复印件。 本人允许哈尔滨工程大学将论文的部分或全部内容编入有关数据 库进行检索，可采用影印、缩印或扫描等复制手段保存和汇编本 学位论文，可以公布论文的全部内容。同时本人保证毕业后结合 学位论文研究课题再撰写的论文一律注明作者第一署名单位为哈 尔滨工程大学。涉密学位论文待解密后适用本声明。 本论文( 口在授予学位后即可口在授予学位1 2 个月后口 解密后) 由哈尔滨工程大学送交有关部门进行保存、汇编等。 作者( 签字) ： 李目绍 导师( 签字) ：吲a 皖 e l 期： 砂，年多月刃日年月 日。 哈尔滨工程大学硕士学位论文 第1 章绪论 1 1 研究意义 海洋蕴含着丰富的资源，海洋环境对全球的气候也有重要的影响。随着 社会经济的发展，人们逐渐重视海洋的研究和开发。而我国拥有辽阔的海域， 其海底有丰富的资源，也是水下目标活动的主要区域，因此对我国海域的传 播特性研究，无论在军用还是民用上都有非常重要的意义。又因为不可能在 所有海域做实际的实验而掌握其传播特性，因此建立符合实际海洋环境的声 场模型，进行声场预报显得尤为重要。 从人们开始研究声传播问题至今，已经开发了许多声场计算模型。依据 数学处理方法的不同，描述海洋中声传播模型主要有射线方法，简正波理论 和抛物近似方法等。 射线法是波动方程的高频近似解，计算速度快，物理意义明确。对低频 以及一些焦散问题处理困难。随着国际上水声学的研究重点转向低频问题， 现在射线法的应用越来越少，基本限制在某些特殊的高频问题。 抛物方程法是波动方程的窄角近似解，近年来发展很快，先后出现了可 以处理宽角、向后散射等问题。 如果将所有的简正波都计算进来，耦合简正波理论的局部剖面解是精确 的，而且在水平方向满足应力和位移连续条件，因此耦合简正波理论的计算 精度最高，它对所有频率都适用。缺点是当需要计算的简正波号数很大时， 计算量大，计算速度慢。当海底倾斜程度加大时，计算步长必须大幅度减小， 从而使计算速度大为降低，计算误差也将增加。因此提高耦合简正波方法的 计算效率得到了研究人员的关注。 1 2 历史背景和研究现状 人们虽然很早以前就发现海洋是声传播的良好载体，但却在二十世纪初 才开始对声波在海洋中传播问题进行真正的研究。第一篇关于水声的论文是 1 9 1 9 年德国人发表的。这篇论文描述了由于海水中弱温度梯度和盐度梯度而 引起声线弯曲，并得到了实验验证。二战后，由于重视的人多起来，因此水 哈尔滨工程大学硕士学位论文 声传播的研究有了突飞猛进的发展。 声场计算模型的研究是从二十世纪六十年代开始起步的，当时只有射线 理论和水平分层的简正波理论。射线理论模型通过声线图可给声场以直观形 象的理解，数学形式简单，至今仍广泛应用于声场计算1 。但若声线在一个 波长内发生弯曲或声强发生变化，这种方法不能给出信赖的声场图像，在浅 海和水文变化剧烈的海洋环境不适用。而水平分层的简正波理论只能处理水 平不变的问题，而且分析的精度比较差。 七十年代开始出现了多种计算模型。虽然抛物型近似方法在波传播的应 用可追溯到4 0 年代，但到了1 9 7 3 年该方法才被t a p p e r t t u 应用于水声传播问 题。抛物型近似方法通过对算子的分解，获得辐射波方程，当赋予距离初值 时，该方程的求解非常有效。当海洋环境与距离无关时，因式分解是准确的。 抛物型近似法要求声场的垂直开角较小( 1 ，得： 州辱忙) ( 2 - 7 1 ) 窘搬。考+ 0 瑟2 _ _ u u ：+ 瑶z ) 一1 k = 。 ( 2 7 2 ) 唆一一) 甜3 ， = 一瓴( 昙瓜一雨f 导卜 。7 非x = 洲嗉等 ( 昙+ 氓一瓜产= o ( 2 7 4 ) 2 1 射线理论是通过射线计算声场的，即对声场的描述是由射线传递声能量。 从声源出发的每一条射线都通过各自的路程到达接收点，因此在接收点的声 场是所有到达射线的叠加。且每一条射线相应的有各自的到达时间和相位， 决定一条射线有两个基本方程：程函方程、强度方程h 1 。 假设： p = a e 却 ( 2 7 5 ) 将( 2 7 5 ) 代入亥姆霍兹方程，分开实部和虚部，得： 去删么一妒) 2 榭= o ( 2 - 7 6 ) 2 ( 跗v o ) + a v 2 妒= 0 ( 2 _ 7 7 ) 方程( 2 7 6 ) 为程函方程确定射线的几何形状。而方程( 2 - 7 7 ) 为强度 方程，确定声波的振幅。这种函数分离是基于一个假定，必须满足近似条件： 在一个波长内，声速不能有太大变化。在近似下，方程( 2 7 6 ) 简化为： ) 2 = k 2 ( 2 7 8 ) 、弋一黻 、： ： 上； r 峪 r波阵面 图2 4 射线和波阵面 等相位面( 9 ：常数) 即为波阵面，波阵面的法线为声线，如图2 4 所 示。考虑边界条件和体积效应，通常本征声线分为四种：直达路径( d p ) 、 折射海面反射( r s r ) 、折射海底反射( r b r ) 、折射一海面反射 哈尔滨工程大学硕+ 学位论文 海底反射( r s r b r ) ，如图2 5 所示。 声源 器 图2 5 本征声线的四种类型 2 5 本章小结 本章在波动理论的基础上，对简正波模型和抛物方程模型进行了研究。 在不同近似条件下，给出距离无关和距离有关简正波解。距离无关的简正波 解介绍了w k b 近似解、w k b z 近似解以及波束位移射线简正波解。详细讨 论了耦合简正波的计算过程。对抛物方程和射线模型做了简单介绍。 哈尔滨工程大学硕士学位论文 第3 章耦合简正波快速算法模型 水声工作者不断追求的目标是快速而精确计算声传播。实际海洋环境非 常复杂，对声传播影响很大。目前计算水平变化声场的模型有射线法、耦合 简正波法、抛物方程方法等。如何提高现有算法的速度和精度，扩展现有算 法的计算能力，寻求新方法具有重要的意义口3 1 。本章在耦合简正波的基础上， 推导了简正波幅度方程的快速求解方法：耦合简正波射线理论和耦合简 正波丢勉物方程方法。 3 1 耦合简正波射线模型 本小节在亥姆霍兹方程的基础上，基于射线声学理论推导出了求解简正 波幅度的微分方程组。 由亥姆霍兹方程： 垮考+ 磐+ 羽0 3 2 p = 。( 3 - 1 0 z ， ，毋a ， 2 c 2 l ，z ) 1 拉普拉斯算子为： v 2 = 吾拿o r ( 舀o r + 专著+ 罢o z ( 3 - 2 ) rr 。o 秽。 。 式( 3 - 1 ) 忽略了方位角方向。 与射线模型类似，用包络函数i ，z ) 表示p 设伫4 1 ： p ：下vp 衍( 3 3 ) 2 。j = p 。 l j j j 则： 考= 眵刖2 1 2 v r - 3 1 2 + i k 井咖 c 3 q 勿i 毋 2 rj 、7 等= 钞坨一2 + 三詈r - 5 1 2 + 2 i ko _ 勿z r _ l ，2 p 5 ， 一i k v r 驯2 一k 2 ；) p 鼢 q r 哈尔滨丁程大学硕士学位论文 可a 2 p = 万i 掣8 船 ( 3 - 6 ) 可2 万丁8 p 哪 将式( 3 - 4 ) ，( 3 5 ) ，( 3 6 ) 代入式( 3 一1 ) ，得： 计算远场，忽略对，的二阶偏导，式( 3 7 ) 的近似描述为： 2 征业+ 堡一2 k 2 u y ：o( 3 8 ) 舯叱班弘南 声线轨迹由h a m i l t o n t 凋方程定义： 鱼：丝丝：一望( 3 9 ) 一= 一一= 一一 、j 一， d ra f d r d g z 、7 其中哈密顿算子为： 日= 譬+ u ( 3 - 1 0 ) f 与射线掠射角z 的关系为： f = d z l d r = t a n z( 3 - 1 1 ) 与距离有关的海洋环境可以分成n 个距离段，在每一个距离段内介质均 匀，此时本征函数、本征值都与距离有关。利用耦合简正波法可设彻： p ，z ) = e 丸o 砂肌，z ) ( 3 1 2 ) 式( 3 1 2 ) 强调本征函数妙与距离，有关。由式( 3 8 ) ，( 3 1 2 ) 得： 詈= ( 扩2 圳心) p 曲；帆+ ，曲( 善等善警九) p 功 m 一r 铲 群 他 v了。 刮 铲一瑟 ， l 屿加 3 4 一 旦万伊砖 哈尔滨工程大学硕士学位论文 标准的本征方程为： 由式( 3 1 0 ) 得： 0 眩2 v = r l 2 e - i k rz 。妒m 磐a z眩 m z 可d 2 u j m + 僻一，；。= 。 由式( 3 1 1 ) 和角度关系得： 日。一u = 了2 后2 f 2 = k 2t a n 2z = k 2 - 7 2 ( 3 1 4 ) ( 3 1 5 ) ( 3 1 6 ) ( 3 1 7 ) 将式( 3 - 1 7 ) 代入式( 3 - 1 5 ) 整理得： 丢争埘呱一u 杪肼= o ( 3 - 1 8 ) 将式( 3 1 8 ) 代入( 3 1 4 ) 得： 笔= r l 2 e - 鼢2 后z 妙一h 。玩少册 ( 3 1 9 ) 知z 厶一 、h “”h 1 、 将( 3 1 3 ) 和( 3 1 9 ) 代入( 3 8 ) 得： ( 万1 一舫) 善屯妙用+ 善警九+ 善笔豫+ k i n m z 埘c m 少讲= 。( 3 - 2 。) 两端乘少。，由正交性，得： 警+ 蝇九+ ( 去一舫弦一军丸謦啪 p 2 ， 其中初始条件为。( o ) = 。( o ，忍) ，即本征函数在声源位置的值。依据初始条 件幅度方程( 3 2 1 ) 可以做数值求解。 3 2 耦合简正波抛物方程模型 依据彭朝晖提出的基于w k b z 理论的耦合简正波抛物方程理论，本 哈尔滨工程大学硕士学位论文 小节在抛物方程的基础上推导出了求解简正波幅度的微分方程组。设如图3 1 所示的声道中，其海面为自由平面，海底为分段变化的不平海底，密度p ， 声速c ，只随深度z 变化，海底界面只有一层，柱坐标系中简谐点源作用下 水平变化声道中声场需满足的波动方程为： 越 赚 0 z v 2 p 一吉窘一吉卯v p - o 和相应的边界条件： f 尸p ，o ) = 0 南掣三裂掣 【 p ，z 咒一= 0 式中：h - 一海水深度 善_ 海底斜面法线方向的偏导数 2 7 ( 3 - 2 2 ) ( 3 - 2 3 ) 哈尔滨工程大学硕士学位论文 旦- o n 高丽 刳 ( 3 - 2 4 ) 一= ；= = = = = l 一一仃一l i j - z 斗i 1 + 7 ) 2l 瑟 务 。 式中吼冽，h = d h ( r ) d r 引入波函数： p = p l p j p = p 4 ； ( 3 2 5 ) 将( 3 2 5 ) 式代入( 3 2 2 ) 得： v p + o 可2 p + f 2 p = o ( 3 - 2 6 ) 式中： 护= 后2 + 万v 2 p 一可3 ( v p ) 2 ( 3 2 7 ) 2 p ap i 、 令： 则得到： y = l o g p ( 3 2 8 ) 排“i 1 矿0 2 yj 4 ( 老) 22 瑟2 k 瑟， ( 3 - 2 9 ) 求解满足边界条件( 3 2 3 ) 的波动方程是非常困难的，通常的方法是祸 合简正波法3 1 1 ，设： p ( z ) = 丸( x ，j ，) 。( z ；x ，y ) ( 3 3 0 ) 式中：丸b ，y ) 简正波幅值函数 婶；x ，川本地简正波本征函数 g ( z ；x ，y ) 满足以下方程： 绰+ 石2 ：砰y 。 ( 3 。3 1 ) a ，三 1“hih 、 和相应的边界条件： l 旨尔浜工程大字坝士字伍论又 fy o ，o ) = 0 1a 雅莎时j 巩) 啪 ) ( 3 - 3 2 ) i 瓦葡t 2 碉t 【y ，z ) l ：= o 本征函数z ；x ，y ) 为完整系，满足正交归一性： p n d z = 6 m n ( 3 - 3 3 ) 式中：地n e c k e rd e l t a 函数。 将( 3 - 3 0 ) 代入( 3 2 6 ) 得： ( 2 v v 丸坝v 2 丸+ p n v 2 丸幅等糟丸= 。 ( 3 - 3 4 ) 随z ，j ，变化缓慢，只保留眠对水平方向的一阶导数。忽略掉第三项， 铲死+ 珐碍2 v v 馋 ( 3 3 5 ) + 碍) 篇警+ 等等) p 3 6 ， 式( 3 3 6 ) 两边乘在深度上积分，利用( 3 - 3 3 ) 得到水平波动方程： + 霸) ( 一誓峨。，挚 ( 3 - 3 7 ) = p 疗爹出( 3 - 3 8 ) b = 札哪- & m - 出 写成矩阵形式，得： 尝+ 磐+ k2 + 2 b j 丝+ 2 b 孚：o ( 3 - 3 9 ) o x 2 a 1 ，2“o x y a 1 ， 上式中，k 都是大写表示向量和矩阵。 哈尔浜i 程大学硕十学位论文 忽略高阶耦合项p 2 3 3 3 4 1 ，得： ( 丢峨) 2 + 等y 詈一= 。 p 4 。， 采用抛物方程方法，做因式分解，得： 睁侈+ 码纠怡壤+ 侈+ 磅册2 = 0 ( 3 州， 只取前向散射，得： 罢= 也+ 侈+ 辑弘卜 p 4 2 ， 在柱坐标系下，拉普拉斯算子为： v 2 = 卫r0 r 8 r + 三r 要0 0 专 ( 3 - 4 3 ) kj 22 瑟2 、7 则在柱坐标系下的波动方程为： c 等弓昙+ 专等聃后2 p = 。( 3 - 4 4 )( 矿+ 7 万+ 7 万+ 萨) p + 后2 0 为了满足能量守恒p 5 3 q ，需要做一修正。定义： 兹= 屈7 旎 ( 3 4 5 ) 则得到能量守恒耦合形式解： 比鲥= 矿乏k 广破他概伍罅o(346)co 1 4 # r l 式中：吒本地本征值 ( z ；秒，) 相应的本征函数 甜。寺求的简正波系数 与在直角坐标系解法相同，得到幅度方程： 啪一倍窘+ 2 砉饬冲= 。 p 4 7 ， 哈尔滨工程大学硕+ 学位论文 式中e 是对角元素为e x p ( 矽 的对角阵，耳，删和岛，。为耦合系数，有以 下形式： b r , n m = ”j 警c 3 如= k 等比 p 4 8 ) 幅度方程( 3 4 7 ) 依据初始条件九( o ) = ( o ，乙) 可做数值计算。 3 3 本章小结 耦合简正波算法中解决耦合系数需要花费很长时间，针对这种问题本章 介绍了两种耦合简正波的解法：简正波射线方法以及简正波抛物方 程方法。这两种方法分别用射线理论和抛物方程理论解决耦合系数问题。而 着重讨论了本论文中要采用的方法：简正波抛物方程方法。 3 l 哈尔滨工程大学硕士学位论文 第4 章耦合简正波抛物方程方法数值计算 4 1 本地简正波解 耦合简正波抛物方程方法数值计算时需要解决本地简正波和简正波 的幅度方程。本小节求解本征值和本征函数，采用有限差分法。 z = 0 z = d 图4 1 单层海洋环境 正如图4 1 所示，我们将区间0 z d 等分成个间隔，构成等间隔 z ，= 办d = 0 , 1 n ) 网格，其中h = d m 为网格宽度。而且为了保证充分抽样 数目n 应选择得足够大，通常每个波长取1 0 个点就足够了。 假定密度是常数，由此简正波满足方程嘲： 1 南“；p - o ( 4 - 1 ) 式中撇号表示对于z 的导数。依照推导有限差分方程的标准方法，利用 泰勒级数展开可得： 。嘲也厅也 等+ 妙j 等+ ( 4 - 2 ) 由此可以得到一阶导数的前向差分近似式： 二二工二 峁 历 勿 缈 卵 哈尔滨丁程大学硕士学位论文 缈，：坚竽嘶吐h ( 4 - 3 ) 2缈，。了_ 二一y _ ，一+ 则一阶导数的d 0 ) 近似式为： 鼍半 利用方程( 4 。1 ) 计算前向差分近似式可以获得： ( 4 4 ) 胖一呙卅荆 ， 由此产生。g 2 ) 近似式为： 罕+ 南斗害 件6 ， 类似地可以获得反向差分近似式。从以下泰勒级数出发 。粥嘶办嵋”等一y ；等+ ( 4 - 7 ) v j d 2 v j v j h 七v j _ 最一v j _ 聂七 u 一| ) 由此产生。o ) 近似式和。g 2 ) 近似式分别为： 。j vj v 1 ， i 一 华一南叫 ， 最后利用式( 4 2 ) 和( 4 7 ) ，可以得n - - 阶导数的中心差分近似式： ( 4 8 ) ( 4 - 9 ) ，”= 堕笋+ 。g 2 ) ( 4 - 1 0 ) 综合这些已有的有限差分近似式，就可以用离散的模拟量代替连续问题 中的导数。重新写出连续i - j 题： 肿高卅刑= 。 厂r 僻( 0 ) + 掣掣= 。 厂7 僻掣掣= 。 利用常微分方程及顶部和底部边界条件鲫，可以得到： 。+ - 2 + h 高卅卜川一o ，一， 知+ 吉降+ 高0 3 2 卅扣 件均 毕+ 高h 一- 扣+ 一扣：。 于是由这些差分方程得到以下形式的代数特征值问题： 4 僻眵= 0 这里妒为向量，其元素为妒0 ，少 ，y 。这些元素是方程( 4 1 9 ) 劓罢高钭+ 者卟阱。 么= d oe 1 p ld le 2 1 e 2d 2e 3 e 一2d ，一2e 一1 e | - l 九一1 e n 3 4 p d n ( 4 1 7 ) ( 4 1 8 ) ( 4 - 1 9 ) ( 4 2 0 ) 一p + 少 厂7 诎 缸程方理整 哈尔滨工程大学硕士学位论文 式中系数嘭和g 的定义为： 巩= 型嘎笋划+ ( 4 - 2 1 ) d ，：二三二! ! 垦竺型，歹：1 ，一1 ( 4 2 2 ) j h p “ 、。 0 = 尘生挚趔+ ( 4 2 3 ) p ，= ，= 1 ，n ( 4 - 2 4 ) 对于压力释放表面，边界条件中出现的比值f g 将变为无限大。在这种 情况下为零，可以从特征值问题矩阵中直接删除第一行第- - n 。解出a 矩 阵的特征值与特征函数就是简正波的本征值与本征函数。 海洋声学问题常常包含声速或密度的不连续性，例如从海水过度到海底 就存在不连续性。处理这类不连续问题时可以把海洋环境分成若干层，使得 在每一层内介质特性平滑变化。在一层内，前面所述的有限差分方程仍是适 用的。然后在界面处导出边界条件把各个层的解结合在一起。 本文中考虑两层界面，一层是海水，另一层是海底。在每一层内构成独 立的有限差分网格，网格间距分别为h 。和h b ，如图4 2 所示。 在海水中，方程( 4 1 ) 的有限差分近似式为： 十棚高一砰，一o ，川，肌 在海底中，有限差分近似式为： 十+ 磁南叫。一o ，剐扎 在界面处，利用边界条件即声压连续和法向振速连续。其中法向振速连 续条件为： 垡笙鱼丛垄：垡丝鱼丛垄( 4 - 2 7 ) p wp 、 3 5 哈尔滨工程大学硕+ 学位论文 知一l z n 1 工吃 z = 0 z = d 1 譬甏辫动= 竿一 南卅妙封成 一 堕+ 二丝丝丝：鱼：) 二笙k 丛篮+ 竺匝鱼：墨缸+ 监：。 二笙型k ：笙鱼：丛k 箧+ 监：o 叫 二二工 五 乞 哈尔滨工程大学硕士字何论文 以得到对称矩阵，与单层的情况完全一样。 在计算声压时需要使用归一化模式。归一化常数是： 虬= 搿出+ 击掣荆+ 砭1 百d ( f g ) 荆 其中积分项可以用梯形法则计算，即： 辫出号( 扣：+ m m + 剐 ， 舯：旷为 4 2 简正波幅度计算 本小节用抛物方程方法求解简正波幅度方程( 3 4 7 ) ，设方位角变化不太 大，简正波沿方位角方向的耦合相对于水平方向弱p 8 3 ，即忽略方位角方向的 耦合，得： 考：聊刊+ 倍窘彬卜 式中矩阵d 的元素可以按下式计纠州1 1 ： 一 九。= 一b m 口： ( 4 3 3 ) f f “一) 毋 墨。一掣 = 讣挑( 扇一罢y 乏) + 誓警( ，一到 ：蓐批勿( 4 - 3 6 ) 式( 4 3 5 ) 中下标1 ，2 分别代表海底界面上下处的参数值。以= 0 孺， 儿：蕊分别为海底界面上下处的垂直波数。 3 7 哈尔滨工程大学硕+ 学位论文 方程( 4 3 2 ) 可用数值方法求解。该方程可以采用分裂步方法求解，这 就意味着需要分别求解以下两个方程： 丝：d 矽 加 警呻川陆鲁搿卜 ( 4 - 3 7 ) ( 4 3 8 ) 方程( 4 3 7 ) 可以采用c r a n k - n i c o l s o n 积分法求解，即： r + a r ( - 一等。) 一( + 等。) 7 方程( 4 3 8 ) 中的矩阵均为对角矩阵，所以可以针对每一号简正波进行 分别求解，即： 挚：弛。- 1 + 订西诊。，朋：1 ，2 ，( 4 - 4 0 ) o r 非肛寿著慨 ，一“二d 一 7 7 。= g ：一三) ： 采用与抛物方程模型f o r 3 d 相同放入隐式有限差分方法对方程( 4 4 0 ) 进行求解，可得： ( 1 + 重x + 廿= ( 1 + p y 为i 二( 4 - 4 1 ) 式中：p = 0 + 8 ) 1 4 ，万= ( 1 8 ) 4 ，万= 啦m a r 。 在本论文中不计及方位角方向将方程( 4 3 7 ) 和( 4 3 8 ) 整理为一个方 程： 掣：p + 一) f 眵 跏 。、“ ( 4 4 2 ) 此方程也可采用c r a n k - n i c o l s o n 积分法求解，即： 矿廿= ( ，一等a ) 。1 ( + 钏矽7 4 3 ， a 是主对角线上的元素为( 后一) f ，其他元素为d 的矩阵。初始场为： 3 8 哈尔滨t 程大学硕士学位论文 4 3 程序流程 ( o ) = y ( z ，0 ) 输入环境参数 j r 本地简正波分岳计算起始场c f 。， l ，= ，+ 血i 0 解方程( 4 4 3 ) 0 计算声压和传播损失 4 ；藩否 ( 4 - 4 4 ) 图4 3 程序流程图 如图4 3 所示，本论文设计的算法程序如下：第一，做距离上的循环。 第二，在循环内输入海水深度、声源深度、接收深度、海水密度、海底介质 密度等环境参数，由于海水深度随着距离变化，因此输入的参数与距离有关。 第三，调用k r a k e n 或k r a k e n c 使用输入的环境参数计算本征值和本征函数。 第四，利用声源位置的本征函数作为初始值，计算幅度方程。第五，结合幅 度值和本征函数计算声场。 4 4 算例与分析 本文编写f o r t r a n 程序实现了耦合简正波抛物方程方法的数值仿真， 3 9 哈尔滨工程大学硕士学位论文 首先需要计算本征值、本征函数，结合计算的幅度，最后计算传播损失。在 计算本征值和本征函数时分别调用f o r t r a n 版k r a k e n 和k r a k e n c 。 4 4 1 水平海底模型 勺= 1 5 0 0 r a l s 见= 1 o g c r ：? 龟= 1 7 0 0 n s p b = 1 5 9 c ? 距离k m 图4 4 具有自由表面和半无限海底的均匀水层 如图4 4 所示，海底为水平时，水声传播建模软件包a c h 腰( a c o u s t i c t o o l b o xu s e ri n t e r f a c ea n dp o s tp r o c e s s o r ) 可以用k r a k e n 和k r a k e n c 计算该海 洋环境下的本征值、本征函数、传播损失。以下为本方法与a c t u p 计算的结 果的比较与分析。 图4 5 是声源频率为2 5 h z ，声源位置为1 0 0 m ，海底吸收为0 5 d b a 时， 本方法中调用k r a k e n 计算的声压传播损失与a c t u p 中k r a k e n 计算的声压 传播损失比较图。从图4 5 中可以看出，两种方法计算的声压传播损失在近 场和远场都具有很高的一致性。 哈尔滨_ t 程大学硕+ 学位论文 仁2 5 h z z s = 1 0 0 m ，z r = 3 0 m r a n g e ( m ) 图4 5k r a k e n 传播损失比较图 仁2 5 h z ，z s = 1 0 0 m ，z r = 3 0 m r a n g e ( m ) 图4 6k r a k e n c 传播损失比较图 图4 6 是频率为2 5 h z ，声源位置为1 0 0 m ，海底吸收为o 5 d b 2 , 时，本方 法中调用k r a k e n c 计算的传播损失与a c t u p 中k r a k e n c 计算的传播损失比较 一p)篮j cois一ec巴卜 一mp)蚜)1 cols一ec坦_l 哈尔滨工程大学硕士学位论文 图。从图中可以看出，这两种方法计算的传播损失在近场有一些差别，在远 场基本一致。在近场有差别的原因是f o r t r a n 版的k r a k e n c 计算本征值和本征 函数与a c t u p 中k r a k e n c 计算的结果比较存在差别。 图4 7 是频率为2 5 h z ，声源位置为1 0 0m ，海底吸收为0 5 d b 2 时，本方 法调用a c t u p 中k r a k e n c 计算出来的本征值和本征函数，计算的传播损失与 a c t u p 中k r a k e n c 计算的传播损失比较图。 图4 7k r a k e n c 传播损失比较图 从图4 7 可以看出，两种方法计算的传播损失无论在近场还是在远场都 具有很高的一致性。 表4 1 是a c t u p 中k r a k e n 计算的本征值与f o r t r a n 中调用k r a k e n 计算 的本征值的比较，表4 2 是a c t u p 中k r a k e n 计算的声源位置本征函数与 f o r t r a n 中调用k r a k e n 计算的声源位置本征函数的比较。表4 3 是a c t u p 中 k r a k e n c 计算的本征值与f o r t r a n 中调用k r a k e n e 计算的本征值的比较，表 4 4 是a c t