本发明涉及一种基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法。
背景技术：
近年来，由于计算机水平的迅猛发展，带动潮汐、潮流数值模拟工作也产生了较大的飞跃，使得人们对各海区的潮汐、潮流分布特征的认识正在逐步的深化。水动力模型主要是根据潮波运动方程，由初始条件和边界条件，通过数值方法对原始微分方程进行差分离散求得数值解，从而得到计算海区内的潮汐分布特征。水动力模型一方面能随计算机水平的进步而大大提高模型的分辨率，但是模型计算所需的海底地形数据、底摩擦系数、湍混合参数以及开边界条件等参数很难较准确的通过人为给定，从而给模型的调试带来了很大的麻烦，其结果仍不能达到令人满意的程度。
技术实现要素：
本发明旨在解决上述问题，提供了一种基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法，它基于潮汐潮流数值模拟计算，从而提高多分潮调和常数计算精度，其采用的技术方案如下：
一种基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法，其特征在于，包括:
步骤1)选定分潮，计算边界条件和平衡潮潮高；
步骤2)计算时段及时间步长；
步骤3)各计算点的水位和流速做调和分析，求出各分潮的振幅hl和迟角gl。
在上述技术方案基础上，步骤1)为设所考虑的分潮有l个，其中l1个为主要分潮，后l-l1个为浅水分潮；各分潮的角速率为ωl，(l＝1,…,l)；
其中
边界条件由式(1.1)计算
其中，a0为平均水位；h,g为潮汐调和常数振幅和迟角；l是分潮个数，一般只取l个分潮中的一部分；ω是分潮角速率；t为时间，以小时为单位；
平衡潮潮高由式(1.2)计算
这里，
在上述技术方案基础上，步骤2)为：计算时段分为两部分，前面为启动阶段，记为t1；后面为实际模拟阶段，记为2t；t1的长度取决于模拟结果达到稳定的周期性变化所需的时间；如海域水浅且流速大，则所需时间长度较短，反之较长；2t的长度取决于所模拟潮汐过程的总的会合周期ts，ts可由下式计算：
ts＝2π/min[|ωl-ωm|],l,m＝1,2,…l,m≠l(2.1)
t可取ts/2附近的一个整数；为具体起见，可取
t＝int[ts/2]+1(2.2)
为方便，令
t'＝t+(t1+t)[即t＝t'-(t1+t)](2.3)
并取时间步长
δt＝1/m,m＝整数
这样，以n＝0，1，…代表时间步序数，则n,t,t'有下列关系
即，t'＝n/m；因t＝t'-(t1+t)，故第n步对应的时间t为
t＝n/m-(t1+t)(2.4)
这里，从n＝0至mt1-1步为启动阶段，从n＝mt1开始，每隔m步，即n＝mt1，mt1+m，mt1+2m，…，mt1+2mt计算结果用于调和分析；由公式(2.4)知，这些时间步对应的时间t为-t,-t+1,-t+2,…,t；
对于各时间步，边界条件和平衡潮潮高仍用公式(1.1)和式(1.2)计算，其中t由公式(2.3)计算。
在上述技术方案基础上，步骤3)为对各计算点的水位ζ和流速分量u及v各自作调和分析；水位ζ与(i,j)有关；流速u和v对二维模式与(i,j)有关，对三维模式与(i,j,k)有关；
对每种变量的每一个点，引入f0,f'l,f"l(l＝1,…，l)，它们初始值均取为0；
略去附标(i,j)；
当计算进行至n＝mt1(即t＝-t),作如下计算
当计算进行至n＝mt1+m(即t＝-t+1)，重复公式(3.1)计算(当然此时f0,f'l,f"l,ζn,t等数值已换成新的值)；…直至n＝mt1+2mt,仍用公式(3.1)算得最后的f0,f'l,f"l；
然后计算alm(l,m＝0,1…,l)和blm(l,m＝1,2,…，l)：
这里，{alm}为(l+1)×(l+1)矩阵，{blm}为l×l矩阵；{alm}，{blm}与(i,j)及要素(ζ或u,v)无关，即对各点的ζ,u,v均适用；
再用迭代法计算a0,al,bl(l＝1,2,…,l)：
式中p＝1,2,…为迭代次数序号，初始值
若
有时，由于模式在时间方向采用了交错网格(跳点格式)，(u,v)和ζ的计算值所对应的时间可能相差δt(对二维模式δt一般等于时间步长δt的一半，即0.5δt；对三维模式δt一般是内模态时间步长δt的一半，即0.5δt。要根据具体数值格式给定)；为统一，约定整小时(即t＝-t,-t+1,…,t)对应ζ时间层，(u,v)时间层在ζ层之后δt(即t＝-t+δt,-t+1+δt,…,t+δt),对(u,v)算出al,bl后先用下式进行校正(a0不必校正)：
式(3.5)中右边(al,bl)为校正前值，左边
最后由关系
al＝hlcosgl,bl＝hlsingl(l＝1,…,l)(3.6)
反过来合成求出振幅hl和迟角gl。
一种提高多分潮调和常数计算精度的系统，其特征在于：其计算模块按照上述基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法进行。
本发明的有益效果为：基于潮汐潮流数值模拟的多分潮高精度调和常数推算，从而提高了潮汐潮流分析和预报的准确度。
具体实施方式
下面结合实施例对本发明作进一步说明：
一种基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法，其特征在于，包括:
步骤1)选定分潮，计算边界条件和平衡潮潮高；
步骤2)计算时段及时间步长；
步骤3)各计算点的水位和流速做调和分析，求出各分潮的振幅hl和迟角gl。
优选的，步骤1)为设所考虑的分潮有l个，其中l1个为主要分潮，后l-l1个为浅水分潮；各分潮的角速率为ωl，(l＝1,…,l)；其中
边界条件由式(1.1)计算
其中，a0为平均水位；h,g为潮汐调和常数振幅和迟角；l是分潮个数，一般只取l个分潮中的一部分；ω是分潮角速率；t为时间，以小时为单位；
平衡潮潮高由式(1.2)计算
这里，
优选的，步骤2)为：计算时段分为两部分，前面为启动阶段，记为t1；后面为实际模拟阶段，记为2t；t1的长度取决于模拟结果达到稳定的周期性变化所需的时间；如海域水浅且流速大，则所需时间长度较短，反之较长；2t的长度取决于所模拟潮汐过程的总的会合周期ts，ts可由下式计算：
ts＝2π/min[|ωl-ωm|],l,m＝1,2,…l,m≠l(2.1)
t可取ts/2附近的一个整数；为具体起见，可取
t＝int[ts/2]+1(2.2)
为方便，令
t'＝t+(t1+t)[即t＝t'-(t1+t)](2.3)
并取时间步长
δt＝1/m,m＝整数
这样，以n＝0，1，…代表时间步序数，则n,t,t'有下列关系
即，t'＝n/m；因t＝t'-(t1+t)，故第n步对应的时间t为
t＝n/m-(t1+t)(2.4)
这里，从n＝0至mt1-1步为启动阶段，计算结果一般用不着，从n＝mt1开始，每隔m步，即n＝mt1，mt1+m，mt1+2m，…，mt1+2mt计算结果用于调和分析；由公式(2.4)知，这些时间步对应的时间t为-t,-t+1,-t+2,…,t；
对于各时间步，边界条件和平衡潮潮高仍用公式(1.1)和式(1.2)计算，其中t由公式(2.3)计算。
进一步，步骤3)为对各计算点的水位ζ和流速分量u及v各自作调和分析；水位ζ与(i,j)有关；流速u和v对二维模式与(i,j)有关，对三维模式与(i,j,k)有关；
对每种变量的每一个点，引入f0,f'l,f"l(l＝1,…，l)，它们初始值均取为0；
下面以ζ为例说明，略去附标(i,j)；
当计算进行至n＝mt1(即t＝-t),作如下计算
当计算进行至n＝mt1+m(即t＝-t+1)，重复公式(3.1)计算(当然此时f0,f'l,f"l,ζn,t等数值已换成新的值)；…直至n＝mt1+2mt,仍用公式(3.1)算得最后的f0,f'l,f"l；
然后计算alm(l,m＝0,1…,l)和blm(l,m＝1,2,…，l)：
这里，{alm}为(l+1)×(l+1)矩阵，{blm}为l×l矩阵；{alm}，{blm}与(i,j)及要素(ζ或u,v)无关，即对各点的ζ,u,v均适用；
再用迭代法计算a0,al,bl(l＝1,2,…,l)：
式中p＝1,2,…为迭代次数序号，初始值
若
有时，由于模式在时间方向采用了交错网格(跳点格式)，(u,v)和ζ的计算值所对应的时间可能相差δt(对二维模式δt一般等于时间步长δt的一半，即0.5δt；对三维模式δt一般是内模态时间步长δt的一半，即0.5δt。要根据具体数值格式给定)；为统一，约定整小时(即t＝-t,-t+1,…,t)对应ζ时间层，(u,v)时间层在ζ层之后δt(即t＝-t+δt,-t+1+δt,…,t+δt),对(u,v)算出al,bl后先用下式进行校正(a0不必校正)：
式(3.5)中右边(al,bl)为校正前值，左边
最后由关系
al＝hlcosgl,bl＝hlsingl(l＝1,…,l)(3.6)
反过来合成求出振幅hl和迟角gl。
一种提高多分潮调和常数计算精度的系统，其特征在于：其计算模块按照上述基于潮汐潮流数值模拟提高多分潮调和常数计算精度的方法进行。
此外，应当理解，虽然本说明书按照实施方式加以描述，但并非每个实施方式仅包含一个独立的技术方案，说明书的这种叙述方式仅仅是为清楚起见，本领域技术人员应当将说明书作为一个整体，各实施例中的技术方案也可以经适当组合，形成本领域技术人员可以理解的其他实施方式。