本发明涉及一种河床床面冲淤变形的求解方法，具体涉及一种河床调整方法。
背景技术：
当前的科学研究中，河床床面冲淤变形的求解方法一般为：通过计算近床泥沙悬移质和推移质的质量，并以泥沙质量守恒方程为基础求解。悬移质和推移质的数值模拟，则大多基于经验公式，这些经验公式大多是建立泥沙浓度和底部流速的对应关系。在采用直接数值模拟法(dns)求解高雷诺数水流时，流速常常具有随机性。将流速代入经验公式后会导致悬移质和推移质的量在空间和时间上都存在很大的变异性。因此，在dns中用传统泥沙数值模拟方式计算河床底部变化，河床常常会凹凸不平，且这些凹凸不平的床面形态通常不满足泥沙颗粒本身在水中的稳定要求。泥沙在水中的休止角，作为决定底部河床深坑和沙脊的最大坡度的物理量，理应在数值模拟中被充分考虑。
技术实现要素：
发明目的：本发明的目的在于针对现有技术中存在的问题，提供一种基于泥沙水下休止角的河床稳定床面调整方法。
技术方案：本发明提供了一种基于泥沙水下休止角的河床调整方法，包括以下步骤：
(1)计算某节点与相邻节点的坡度，将坡度与泥沙水下休止角的正切值相比较并赋值坡比系数；
(2)根据坡比系数计算分类系数，并利用分类系数对节点进行分类；
(3)计算相邻节点的理想调整高度，将理想调整高度排序；
(4)结合节点类型，分阶段对河床高度进行调整；
(5)重复上述步骤数次后可得到满足泥沙水下休止角要求的河床形态。
进一步，步骤(1)首先计算在结构网格下某i节点与其相邻j节点的坡度kj：
kj＝(yj-yi)/dj(1)
式中，yi是i节点的河床高程，yj是j节点的河床高程，dj为j节点与i节点的水平间距，j＝1，2，3，4；
坡比系数由下式得到：
aj为坡比系数，α为泥沙水下休止角。
进一步，步骤(2)分类系数c＝(cp，cn)采用下式计算：
式中，计算得出的cp，cn实其本质含义分别为坡比系数为正和为负的个数。
进一步，步骤(2)节点的分类包括：s型节点、d型节点和t型节点；其中，
s类节为分类系数中至少一项为0，包括s0、s1、s2、s3、s4，数字表示除0之外的另一分类系数，该类节点的物理特征是相邻节点的不稳定坡度类型仅为正坡或者负坡；
d类节点为分类系数中至少一项为1且另一项不为0，包括d1、d2、d3，数字表示除1之外的另一分类系数；
t型节点为分类系数两项均为2。
进一步，步骤(3)理想调整高度采用下式计算：
式中，
进一步，步骤(4)用角标a、b、c、d表示相邻节点标号1、2、3、4的一个全排列，该排列满足理想调整高度dya0，dyb0，dyc0，dyd0按从小到大排序，节点的调整方法包括以下步骤：：
a.若为s0型节点：则dya0，dyb0，dyc0，dyd0均为0，无需调整；
b.若为s型节点：则dya0，dyb0，dyc0，dyd0不为0的数量与其子类型有关，无论何种子类型，均可用四阶段计算s型节点，各节点调整公式如下：
阶段一以a号节点为基本调整节点，调整计算公式如下：
其中，ya，yb，yc，yd分别表示与i相邻四个节点的河床高度，上标0表示起始高度，即应用调整方法前的高度，上标1表示该第一阶段末更新后的河床高度；
阶段二以b号节点为基本调整节点，调整计算公式如下：
其中，上标2表示该第二阶段末更新后的河床高度，上标1表示该阶段的初始值；
阶段三以c号节点为基本调整节点，调整计算公式如下：
其中，上标3表示第三段末更新后的河床高度，上标2表示该阶段的初始值；
阶段四以d号节点为基本调整节点，调整计算公式如下：
其中，上标4表示第四段末更新后的河床高度，3表示该阶段的初始值；
c.若为d1型节点，dyc0，dyd0不为0，dya0，dyb0为0，此时c，d两节点需要调整，该类型节点调整仅有一个阶段：
阶段一以c，d两节点同时为基本调整节点，调整计算公式如下：
d.若为d2型节点，dyb0，dyc0，dyd0不为0，其调整阶段数与相邻节点特点有关：
依据以下计算式分类：
l＝(ab2)acad
(d1)若l值为1，阶段一以b号节点为基本调整节点，调整的计算公式为：
随后更改对应坡比系数ab为0，依次采用s型节点计算公式(7)、(8)进行计算，采用式(7)计算时用该阶段末的值作为式(7)的初始值，采用式(8)计算时用式(7)阶段末的值作为式(8)的初始值；
(d2)若l值为-1，阶段一以b号节点为基本调整节点，调整的计算公式为：
随后更改对应坡比系数ab为0，采用d1型节点计算公式计算阶段二，当采用式(9)计算时，需要将上标对应为0的地方改成1，1改成2，以表示代入阶段一的计算结果进行进一步调整。
e.若为d3型节点，dya0，dyb0，dyc0，dyd0均不为0，相邻4个节点均需要调整；
(e1)若坡比系数aa与其他三个符号不同
阶段一以a号节点为基本调整节点，调整的计算公式如下：
随后更改对应坡比系数aa为0，并依次采用s型节点的阶段计算方法计算，即式(6)、(7)、(8)；采用式(6)计算时用该阶段末的值作为式(6)的初始值，采用式(7)计算时用式(6)阶段末的值作为式(7)的初始值，采用式(8)计算时用式(7)阶段末的值作为式(8)的初始值；
(e2)若坡比系数aa与其他三个符号存在相同
阶段一以a号节点为基本调整节点，调整的计算公式如下：
随后更改对应坡比系数aa为0，并采用d2型节点的阶段计算方法计算剩余阶段；
f.若为t型节点，阶段一以a号节点为基本调整节点，调整的计算公式如下：
随后更改对应的坡比系数aa为0，并采用d2型节点的计算方法计算剩余阶段，在采用式(11)或者式(12)计算时，需要将上标对应为0的地方改成1，1改成2，以表示代入阶段一的计算结果进行进一步调整。
进一步，重复步骤(1)-(5)直至达到一定循环次数或者最终i节点变为s0型节点。
有益效果：本发明主要应用于矩形网格中，通过计算各个空间点与其相邻点的坡度，利用泥沙水下休止角判断该点与其周围各点的高程关系是否稳定，针对不满足稳定的点，比较该点与其相邻点的坡度和泥沙水下休止角正切值的大小，根据大小关系将其分类，将可能出现的不同的河床形态进行分类；同时，本发明具有分阶段调整河床的特点，即将某点的河床调整过程依据其类别分成不同的阶段，在每个阶段，采用不同的计算公式进行河床的调整，且一阶段只针对相邻的一至二个节点，各个阶段的河床调整均满足土体质量守恒，即具有守恒性，且在单一阶段内，各个调整具有同时性，守恒性和同时性使得整个河床的调整过程更符合物理实际，综合考虑了相邻节点之间的相互影响以及河床土体质量守恒等因素，使得调整过程具有良好的物理仿真性，调整结果更符合物理规律，也有利于提高河床演变模拟的效率。虽然本发明只给出单节点即其临近点的调整方式，但其易于扩展应用至所有空间网格节点，可以与结构化网格的水力学计算程序较好结合。
附图说明
图1是本发明物理概念图；
图2是实施例最大理想高度之和随循环次数的变化图；
图3是河床高程云图对比，图中(a)为原图，(b)为应用本发明调整后的云图。
具体实施方式
下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。
一种基于泥沙水下休止角的河床调整方法，定义i节点和其相邻四节点的水平距离均为单位长度1，泥沙休止角为45度，i点与其各相邻节点的高度为：假设i节点的高度为5，与其周围节点高度分别设为12.8、11.2、5.5、8.9。
调整河床过程如下：
(1)计算坡比系数(a1，a2，a3，a4，)＝(1，1，0，1)。
(2)计算分类系数c＝(3，0)，依据分类依据可知为s3型节点。
(3)计算理想调整高度并将其从小到大排序如下，式(4)的物理概念如图1所示：(dy03，dy04，dy02，dy01)＝(0，1.45，2.6，3.4)，则三个调整阶段的依次调整3，4，2，1四个节点。
(4)给出s3型节点阶段调整公式分阶段调整：
阶段一以3号阶段一以4号节点为基本调整节点并采用式(5)计算，可得各项均为0。
阶段二以4号节点为基本调整节点并采用式(6)计算：
代入后可计算得阶段二4号节点需要调整的高度为0.725，阶段二末各点调整后的高度为(yi2，y12，y22，y32，y42)＝(7.175，12.075，10.475，5.5，8.175)。
阶段三采用以2号节点为基本调整节点并采用式(7)计算：
将阶段二末的结果代入上式后，可计算得阶段三2号节点需要调整的高度为0.767，阶段三末各点调整后的高度为(yi3，y13，y23，y33，y43)＝(8.709，11.308，9.708，5.5，8.175)。
阶段四以1号节点为基本调整节点并采用式(8)计算：
将阶段三末的结果代入上式后，可计算得阶段四1号节点需要调整的高度为0.8，阶段四末各点调整后的高度为(yi4，y14，y24，y34，y44)＝(9.509，10.508，9.708，5.5，8.175)。
(6)将调整后各节点的高程求和为43.4，与起初河床高度求和相等，说明该发明具有质量守恒性。
(7)用更新后的节点值再次重复上述步骤，首先计算坡比系数为(a1，a2，a3，a4，)＝(0，0，-1，-1)；再计算分类系数c＝(0，2)为s2型节点。计算理想调整高度(dy01，dy02dy04，dy03)＝(0，0，0.167，1.505)，因此四个阶段依次调整1，2，4，3节点。
(8)阶段一和阶段二分别以1，2号节点为基本节点进行调整，代入后得到各点高程不改变，阶段三以4号节点为基本调整节点并采用式(7)计算，可以计算得阶段三4号节点需要调整的高度为0.112，阶段三末各点调整后的高度为(yi3，y13，y23，y33，y43)＝(9.285，10.508，9.708，5.612，8.287)。
(9)阶段四以3号节点为基本调整节点并采用式(8)计算，可以计算得阶段四3号节点需要调整的高度为1.337，阶段二末各点调整后的高度为(yi4，y14，y24，y34，y44)＝(7.948，10.508，9.708，6.949，8.287)。同样的，调整后各床面节点求和为43.4，依旧保持质量守恒。
通过分析上述计算结果可知，在本发明调整河床的过程中，河床土体质量始终保持守恒，且随着调整过程的继续，整体河床趋于平缓。此外，可用理想调整高度之和来判定河床是否平滑。如图2为上述算例各个循环初始时刻的理想调整高度之和，从图中可知，经过6次循环后，理想调整高度之和由起初的7.45变为0.058，且逐渐趋于0，说明河床经过调整逐渐趋于平整。此外，本发明在各阶段给出的计算公式，保证了每个调整阶段不会存在过度调整的情况，即不会出现原本的正坡由于调整而变成负坡的现象。最后通过图3展示本发明在计算水力学中的实际应用，图3为应用本发明前后的河床高度云图，图a为某数值模拟结果计算河床变形的原始结果，河床高程点存在突兀和不连续型，图b为在数值计算中应用该发明进行河床调整后的模拟结果，可以发现原本凹凸不平的河床经过本发明方法应用后变得更加平整，说明本发明在实际计算水力学中的应用成果合理性。