一种颗粒流微观力学参数反演方法
阅读说明：本技术 一种颗粒流微观力学参数反演方法 (Particle flow micromechanics parameter inversion method ) 是由 王敏 万文 赵延林 于 2021-06-29 设计创作，主要内容包括：本发明公开了一种颗粒流微观力学参数反演方法,首先通过物理实验获取岩石宏观力学参数,然后基于物理宏观力学参数反演出颗粒流微观力学参数,接着通过数值模拟计算获取数值模拟宏观力学参数,最后采用模拟退火算法调整颗粒流微观力学参数,使得数值模拟宏观力学参数与物理实验宏观力学参数不断接近,当数值模拟宏观力学参数与物理实验宏观力学参数之间的误差小于10%时,则此时对应的颗粒流微观力学参数则是所需要确定的微观力学参数。本发明操作简便,微观力学参数调整过程中无主观性,收敛速度较好,适用于颗粒流数值模拟计算领域。(The invention discloses a particle flow micromechanics parameter inversion method, which comprises the steps of firstly obtaining rock micromechanics parameters through a physical experiment, then performing particle flow micromechanics parameters based on the physical micromechanics parameters, then obtaining numerical simulation macroscopic mechanics parameters through numerical simulation calculation, and finally adjusting the particle flow micromechanics parameters by adopting a simulated annealing algorithm, so that the numerical simulation macroscopic mechanics parameters are continuously close to the physical experiment macroscopic mechanics parameters, and when the error between the numerical simulation macroscopic mechanics parameters and the physical experiment macroscopic mechanics parameters is less than 10%, the corresponding particle flow micromechanics parameters are the micromechanics parameters to be determined. The method is simple and convenient to operate, has no subjective property in the adjustment process of the micromechanics parameters, has good convergence rate, and is suitable for the field of particle flow numerical simulation calculation.)
技术领域
本发明涉及一种颗粒流微观力学参数反演方法。
背景技术
在开展颗粒流数值模拟计算之前需要确定颗粒流数值模拟模型相应的微观力学参数，该微观力学参数不能随意选取，为了使所选取的微观力学参数对应的数学模型与工程实践更加贴切，一般通过物理试验获取岩石的相关力学参数：单轴抗压强度、弹性模量及泊松比，然后调节颗粒流微观力学参数，使其数值模拟模型对应的宏观力学参数：单轴抗压强度、弹性模型、泊松比与物理试验宏观力学参数靠近，当数值模拟实验得到的宏观力学参数与物理试验宏观力学参数之间误差较小时，则此时对应的微观力学参数即是所求的微观力学参数。目前一般采用试错法不断调整颗粒流微观力学参数使数值模拟计算的宏观力学参数与物理实验宏观力学参数不断接近，当其误差量小于某一值时，则停止调试。采用试错法主要依赖于试验者对颗粒流软件掌握的熟练程度，且调试过程中存在一定的盲目性，同时获取好的微观力学参数具有较大的偶然性。
发明内容
为了解决上述技术问题，本发明提供一种算法简单、收敛速度好的颗粒流微观力学参数反演方法。
本发明解决上述问题的技术方案是：一种颗粒流微观力学参数反演方法，包括以下步骤：
步骤一：通过物理实验获取岩石物理实验宏观力学参数；
步骤二：基于物理宏观力学参数反演出颗粒流微观力学参数；
步骤三：通过数值模拟计算获取数值模拟宏观力学参数；
步骤四：采用模拟退火算法调整颗粒流微观力学参数，使得数值模拟计算宏观力学参数与物理实验宏观力学参数不断接近，当数值模拟宏观力学参数与物理实验宏观力学参数之间的误差小于10％时，则此时对应的颗粒流微观力学参数则是所需要确定的微观力学参数。
上述颗粒流微观力学参数反演方法，所述步骤一中，物理实验宏观力学参数包括单轴抗压强度UCSexperimental、弹性模量Eexperimental、泊松比vexperimental。
上述颗粒流微观力学参数反演方法，所述步骤二中，颗粒流微观力学参数包括颗粒密度ρ、颗粒最大半径与最小半径比值Rmax/Rmin、颗粒最小半径Rmin、颗粒接触刚度Ec、颗粒法向刚度与切向刚度比值kn/ks、平行连接刚度平行连接法向刚度与切向刚度比值颗粒摩擦系数μ、平行连接法向刚度平均值σc-mean、平行连接法向刚度方差σc-std、平行连接剪切刚度平均值τc-mean、平行连接剪切刚度方差τc-std。
上述颗粒流微观力学参数反演方法，所述步骤三中，将步骤二所得到的流微观力学参数开展数值模拟计算，并获取该组微观力学参数所对应的数值模拟宏观力学参数，即单轴抗压强度UCSnumerical、弹性模量Enumerical、泊松比vnumerical。
上述颗粒流微观力学参数反演方法，所述步骤四中，采用模拟退火算法调整颗粒流微观力学参数的具体步骤为：
步骤1：初始化模拟退火算法超参数：温度Temperature、降温系数Decay、马可夫链长度Markov及数据跳动因子Stepfactor，马可夫链长度Markov表示每一温度条件下循环次数，迭代次数Iteration＝0；
初始数值模拟计算时，微观力学参数：颗粒密度、颗粒最大半径与最小半径比值、颗粒最小半径、颗粒接触刚度、颗粒法向刚度与切向刚度比值、平行连接刚度、平行连接法向刚度与切向刚度比值、颗粒摩擦系数、平行连接法向刚度平均值、平行连接法向刚度方差、平行连接剪切刚度平均值、平行连接剪切刚度方差的最大值和最小值分别为ρmax、Rmax-max/Rmin-max、Rmin-max、Ec-max、kn-max/ks-max、μmax、σc-mean-max、σc-std-max、τc-mean-max、τc-std-max和ρmin、Rmax-min/Rmin-min、Rmin-min、Ec-min、kn-min/ks-min、μmin、σc-mean-min、σc-std-min、τc-mean-min、τc-std-min，然后在微观力学参数最大值与最小值之间选取随机值ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre作为微观力学参数初始值，采用该组微观力学参数开展数值模拟单轴压缩实验，获得数值模拟宏观力学参数单轴抗压强度、弹性模量和泊松比分别为UCSnumerical-pre、Enumerical-pre、vnumerical-pre；
计算获得的数值模拟宏观力学参数与物理试验获得的宏观力学参数之间的相对误差最大值Judgepre，Judgepre即作为判别值表示为：
同时将初始的微观力学参数视为最优的微观力学参数：
其对应的判别值也是目前最优的相对误差值：
Judgebest＝Judgepre (3)
步骤2：基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)确定下一组微观力学参数的数值，其计算公式如下：
上式中rand为-1到1之间的随机数，通过式(4)基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)，产生当前的微观力学参数集合(ρnext、Rmax-next/Rmin-next、Rmin-next、Ec-next、kn-next/ks-next、μnext、σc-mean-next、σc-std-next、τc-mean-next、τc-std-next)，采用当前的微观力学参数集合开展单轴压缩数值模拟实验获取相应的数值模拟模型宏观力学参数：单轴抗压强度UCSnumerical-next、弹性模量Enumerical-next和泊松比vnumerical-next，并计算判别函数的值Judgenext：
步骤3：如果计算得到的Judgenext<Judgebest，则更新最优的微观力学参数集合及最优的目标值：
Judgebest＝Judgenext (7)
如果Judgenext>Judgebest，则无需更新最优微观力学参数集合；
步骤4：如果计算得到的Judgenext<Judgepre，则更新先前颗粒流微观力学参数及先前判断函数值：
Judgepre＝Judgenext (9)
而当Judgenext≥Judgepre时，先计算概率值p1，p1表示为：
然后产生一个0到1的随机数p2，如果p1>p2则更新先前颗粒流微观力学参数及先前判断函数值，即执行公式(8)和公式(9)，否则不执行任何操作；
步骤5：迭代次数更新，Iteration＝Iteration+1，当Iteration>Markov时：Iteration＝0，当前温度下降Temperature＝Temperature×Decay；执行完步骤5之后，再次执行步骤2，依此循环。
上述颗粒流微观力学参数反演方法，所述步骤四中，数值模拟计算终止条件是数值模拟计算得到的宏观力学参数与物理实验得到的宏观力学参数差值相对误差小于10％；数值模拟计算过程中，每迭代一次就对Judgebest进行判定，当Judgebest≤10％时，则停止计算，其中Judgebest表示为：
式中，UCSnumerical-best，Enumerical-best，vnumerical-best表示计算过程中最佳的单轴抗压强度、弹性模量及泊松比组合。
本发明的有益效果在于：本发明首先通过物理实验获取岩石宏观力学参数，然后基于物理宏观力学参数反演出颗粒流微观力学参数，接着通过数值模拟计算获取数值模拟宏观力学参数，最后采用模拟退火算法调整颗粒流微观力学参数，使得数值模拟宏观力学参数与物理实验宏观力学参数不断接近，当数值模拟宏观力学参数与物理实验宏观力学参数之间的误差小于10％时，则此时对应的颗粒流微观力学参数则是所需要确定的微观力学参数。本发明操作简便，微观力学参数调整过程中无主观性，收敛速度较好，适用于颗粒流数值模拟计算领域。
附图说明
图1为本发明的流程图。
图2为本发明采用模拟退火算法调整颗粒流微观力学参数的流程图。
具体实施方式
下面结合附图和实施例对本发明做进一步的说明。
如图1所示，一种颗粒流微观力学参数反演方法，包括以下步骤：
步骤一：通过物理实验获取岩石物理实验宏观力学参数，物理实验宏观力学参数包括单轴抗压强度UCSexperimental、弹性模量Eexperimental、泊松比vexperimental。
步骤二：基于物理宏观力学参数反演出颗粒流微观力学参数，颗粒流微观力学参数包括颗粒密度ρ、颗粒最大半径与最小半径比值Rmax/Rmin、颗粒最小半径Rmin、颗粒接触刚度Ec、颗粒法向刚度与切向刚度比值kn/ks、平行连接刚度平行连接法向刚度与切向刚度比值颗粒摩擦系数μ、平行连接法向刚度平均值σc-mean、平行连接法向刚度方差σc-std、平行连接剪切刚度平均值τc-mean、平行连接剪切刚度方差τc-std。
步骤三：将步骤二所得到的流微观力学参数开展数值模拟计算，并获取该组微观力学参数所对应的数值模拟宏观力学参数，即单轴抗压强度UCSnumerical、弹性模量Enumerical、泊松比vnumerical。
步骤四：采用模拟退火算法调整颗粒流微观力学参数，使得数值模拟计算宏观力学参数与物理实验宏观力学参数不断接近，当数值模拟宏观力学参数与物理实验宏观力学参数之间的误差小于10％时，则此时对应的颗粒流微观力学参数则是所需要确定的微观力学参数。
采用模拟退火算法调整颗粒流微观力学参数的具体步骤为：
步骤1：初始化模拟退火算法超参数：温度Temperature、降温系数Decay、马可夫链长度Markov及数据跳动因子Stepfactor，马可夫链长度Markov表示每一温度条件下循环次数，迭代次数Iteration＝0；
初始数值模拟计算时，微观力学参数：颗粒密度、颗粒最大半径与最小半径比值、颗粒最小半径、颗粒接触刚度、颗粒法向刚度与切向刚度比值、平行连接刚度、平行连接法向刚度与切向刚度比值、颗粒摩擦系数、平行连接法向刚度平均值、平行连接法向刚度方差、平行连接剪切刚度平均值、平行连接剪切刚度方差的最大值和最小值分别为ρmax、Rmax-max/Rmin-max、Rmin-max、Ec-max、kn-max/ks-max、μmax、σc-mean-max、σc-std-max、τc-mean-max、τc-std-max和ρmin、Rmax-min/Rmin-min、Rmin-min、Ec-min、kn-min/ks-min、μmin、σc-mean-min、σc-std-min、τc-mean-min、τc-std-min，然后在微观力学参数最大值与最小值之间选取随机值ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre作为微观力学参数初始值，采用该组微观力学参数开展数值模拟单轴压缩实验，获得数值模拟宏观力学参数单轴抗压强度、弹性模量和泊松比分别为UCSnumerical-pre、Enumerical-pre、vnumerical-pre；
计算获得的数值模拟宏观力学参数与物理试验获得的宏观力学参数之间的相对误差最大值Judgepre，Judgepre即作为判别值表示为：
同时将初始的微观力学参数视为最优的微观力学参数：
其对应的判别值也是目前最优的相对误差值：
Judgebest＝Judgepre (3)
步骤2：基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)确定下一组微观力学参数的数值，其计算公式如下：
上式中rand为-1到1之间的随机数，通过式(4)基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)，产生当前的微观力学参数集合(ρnext、Rmax-next/Rmin-next、Rmin-next、Ec-next、kn-next/ks-next、μnext、σc-mean-next、σc-std-next、τc-mean-next、τc-std-next)，采用当前的微观力学参数集合开展单轴压缩数值模拟实验获取相应的数值模拟模型宏观力学参数：单轴抗压强度UCSnumerical-next、弹性模量Enumerical-next和泊松比vnumerical-next，并计算判别函数的值Judgenext：
步骤3：如果计算得到的Judgenext<Judgebest，则更新最优的微观力学参数集合及最优的目标值：
Judgebest＝Judgenext (7)
如果Judgenext>Judgebest，则无需更新最优微观力学参数集合；
步骤4：如果计算得到的Judgenext<Judgepre，则更新先前颗粒流微观力学参数及先前判断函数值：
Judgepre＝Judgenext (9)
而当Judgenext≥Judgepre时，先计算概率值p1，p1表示为：
然后产生一个0到1的随机数p2，如果p1>p2则更新先前颗粒流微观力学参数及先前判断函数值，即执行公式(8)和公式(9)，否则不执行任何操作；
步骤5：迭代次数更新，Iteration＝Iteration+1，当Iteration>Markov时：Iteration＝0，当前温度下降Temperature＝Temperature×Decay；执行完步骤5之后，再次执行步骤2，依此循环。
数值模拟计算终止条件是数值模拟计算得到的宏观力学参数与物理实验得到的宏观力学参数差值相对误差小于10％；数值模拟计算过程中，每迭代一次就对Judgebest进行判定，当Judgebest≤5％时，则停止计算，其中Judgebest表示为：
式中，UCSnumerical-best，Enumerical-best，vnumerical-best表示计算过程中最佳的单轴抗压强度、弹性模量及泊松比组合。
实施例
一种颗粒流微观力学参数反演方法，包括基于颗粒流微观力学参数开展单轴压缩实验获取数值模拟模型宏观力学参数、根据数值模拟宏观力学参数与物理实验宏观力学参数差值采用模拟退火算法调整颗粒流微观力学参数、数值模拟计算终止条件；结合具体实例其步骤如下：
步骤1：通过室内实验获取岩石的单轴抗压强度UCSexprimental、弹性模量Eexperimental、泊松比vexperimental分别为：37MPa、23.3GPa和0.17。同时确定模拟退火算法的超参数温度Temperature为100、降温系数Decay为0.99、马可夫链长度Markov(每一温度条件下循环次数)为10000，数据跳动因子Stepfactor为0.02。并确定微观力学参数(颗粒密度ρ、颗粒最大半径与最小半径比值Rmax/Rmin、颗粒最小半径Rmin、颗粒接触刚度Ec、颗粒法向刚度与切向刚度比值kn/ks、平行连接刚度平行连接法向刚度与切向刚度比值颗粒摩擦系数μ、平行连接法向刚度平均值σc-mean、平行连接法向刚度方差σc-std、平行连接剪切刚度平均值τc-mean、平行连接剪切刚度方差τc-std)的最小值和最大值分别为：(1000,1,0.1e-3,1e9,0.1,1e9,0.1,0.1,1e6,1e6,1e6,1e6)和(5000,10,10e-3,1000e9,10,1000e9,10,10,1000e6,100e6,1000e6,100e6)。
步骤2：初始微观力学参数(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)是通过在微观力学参数最小值和最大值直接随机选取，然后得到初始微观力学参数组合。然后采用初始微观力学参数组合开展数值单轴压缩实验获得宏观力学参数：单轴抗压强度、弹性模量和泊松比分别为UCSnumerical-pre、Enumerical-pre、vnumerical-pre。计算出目标函数：
同时将初始的微观力学参数视为最优的微观力学参数：
其对应的判别值也设置为目前最优的判别值：
Judgebest＝Judgepre
步骤3：基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)确定下一组微观力学参数的数值，其计算公式如下：
式中rand为-1到1的随机数。通过上式可基于先前的微观力学参数集合(ρpre、Rmax-pre/Rmin-pre、Rmin-pre、Ec-pre、kn-pre/ks-pre、μpre、σc-mean-pre、σc-std-pre、τc-mean-pre、τc-std-pre)，产生当前的微观力学参数集合(ρnext、Rmax-next/Rmin-next、Rmin-next、Ec-next、kn-next/ks-next、μnext、σc-mean-next、σc-std-next、τc-mean-next、τc-std-next)，采用当前的微观力学参数集合开展单轴压缩数值模拟实验获取相应的数值模拟模型宏观力学参数：单轴抗压强度UCSnumerical-next、弹性模量Enumerical-next和泊松比vnumerical-next。并计算判别函数的值Judgenext：
步骤4：如果计算得到的Judgenext<Judgebest，则更新最优的微观力学参数集合及最优的目标值：
Judgebest＝Judgenext
如果Judgenext>Judgebest，则无需更新最优微观力学参数集合。
步骤5：如果计算得到的Judgenext<Judgepre，则更新先前颗粒流微观力学参数及先前判断函数值：
Judgepre＝Judgenext
而当Judgenext≥Judgepre时，先计算概率值p1，p1可表示为：
然后产生一个0到1的随机数p2，如果p1>p2则需更新先前颗粒流微观力学参数及先前判断函数值，即执行公式(17)和公式(18)，否则不执行任何操作。
步骤6：迭代次数更新，Iteration＝Iteration+1，当Iteration>10000时：Iteration＝0，当前温度下降Temperature＝Temperature×0.99。执行完步骤6之后，再次执行步骤3，依次循环。
所述数值模拟计算终止条件是数值模拟计算得到的宏观力学参数与物理实验得到的宏观力学参数差值相对误差小于10％时，则终止数值模拟计算。数值模拟计算过程中，每迭代依次就对Judgebest进行判定，当Judgebest≤10％时，则停止计算。其中Judgebest可表示为：
式中，UCSnumerical-best，Enumerical-best，vnumerical-best表示计算过程中最佳的单轴抗压强度、弹性模量及泊松比组合。通过计算得到UCSnumerical-best，Enumerical-best，vnumerical-best分别为35.82MPa、23.18GPa和0.16，其中对应的Judgebest为5.8％，小于实验终止条件10％。
同时得到对应最佳的微观力学参数集合(ρbest、Rmax-best/Rmin-best、Rmin-best、Ec-best、kn-best/ks-best、μbest、σc-mean-best、σc-std-best、τc-mean-best、τc-std-best)为(2411，5.54，0.76e-3，1.76e9，0.57，55.96e9，2.94,6.28，77.63e6，77.76e6，64.04e6，33.52e6)。
- 上一篇：一种医用注射器针头装配设备
- 下一篇：态势图生成方法、装置、非易失性存储介质及处理器