阅读说明：本技术 一种基于松鼠搜索算法的电阻率溶洞识别方法 (Resistivity karst cave identification method based on squirrel search algorithm ) 是由 陈建华 罗盈洲 韩庆忠 邓亮 于 2019-11-07 设计创作，主要内容包括：本发明公开了一种基于松鼠搜索算法的电阻率溶洞识别方法,该方法主要步骤如下：S1.构建溶洞有限元模型,利用电阻抗法正问题相关公式求得结构测量点的电势分布数据；S2.构建溶洞结构的目标函数,即待优化的目标函数；S3.利用松鼠搜索算法不断的迭代优化目标函数,满足停止准则后,保存最优解,便可识别出溶洞位置。该方法对比传统的灵敏度或梯度矩阵方法可以准确地检测出各种溶洞的位置,无需初值和梯度信息,并且对噪声不敏感。本发明通过单个溶洞和多个溶洞两个实施例说明了该方法的有效性和鲁棒性,与粒子群算法相比,本发明即使在一定的测量噪声下也具有较好的探测精度,具有良好的工程应用能力。(The invention discloses a resistivity karst cave identification method based on a squirrel search algorithm, which mainly comprises the following steps: s1, constructing a karst cave finite element model, and solving potential distribution data of a structural measurement point by using an electric impedance method positive problem correlation formula; s2, constructing a target function of the karst cave structure, namely the target function to be optimized; and S3, continuously and iteratively optimizing the objective function by using a squirrel search algorithm, and storing the optimal solution after the stopping criterion is met so as to identify the karst cave position. Compared with the traditional sensitivity or gradient matrix method, the method can accurately detect the positions of various karst caves without initial values and gradient information and is insensitive to noise. The effectiveness and robustness of the method are illustrated through two embodiments of a single karst cave and a plurality of karst caves, and compared with a particle swarm algorithm, the method has better detection precision even under certain measurement noise and good engineering application capability.)

1. A resistivity karst cave identification method based on a squirrel search algorithm is characterized by comprising the following steps:

the method comprises the following steps: establishing a finite element model of a karst cave soil layer, determining an electrode scheme, and obtaining potential distribution data of a measuring point by using an electrical impedance method;

step two: constructing an objective function of the karst cave structure, namely an optimized objective function, wherein the objective function is as follows:

where g (c) is the objective function, | () | represents the 2 norm of the vector, i.e.c is the apparent resistivity of each cell, Su represents the data on the set of measured potential points,

step three: and continuously optimizing the objective function by using a squirrel search algorithm, and finally obtaining the identification position of the karst cave after meeting the stop criterion.

2. The resistivity cavern identification method based on the squirrel search algorithm as claimed in claim 1, wherein: the concrete process of optimizing the objective function by utilizing the squirrel search algorithm in the third step is as follows:

s1, random initialization: assuming n squirrels (FS), the position of the ith squirrel can be represented by a vector, and the positions of all squirrels can be represented by the following matrix:

wherein, FS_i.jRepresenting the j-th dimension of the ith squirrel, the initial position of each squirrel in the forest was assigned with a uniform distribution:

FS_i＝FS_L+U(0，1)×(FS_U-FS_L) (3)

wherein, FS_LAnd FS_URespectively the lower bound and the upper bound of the ith squirrel in the j dimension, and U (0,1) is a value range [0,1 ]]Uniformly distributed random numbers inside;

s2, fitness evaluation: the values of the decision variables (solution vectors) are put into a user-defined fitness function, the position fitness of each squirrel is calculated, and the corresponding values are stored in the following array:

s3, sorting, declaring, and randomly selecting: after storing the fitness value of the position of each squirrel, the arrays are sorted in ascending order;

s4, generating a new location: the mathematical model of the new position has the following three cases:

wherein d is_gFor random glide distance, R₁、R₂And R₃Is [0,1 ]]Random number in the range, FS_nt、FS_atAnd FS_htThe squirrel positions of the common tree, the oak tree and the pecan tree are reached respectively, t is the current iteration number, G_cHas a value of 1.9, P_dpThe value of (A) is 0.1;

s5, seasonal monitoring conditions: seasonal monitoring conditions are introduced into the algorithm, so that the algorithm is prevented from falling into a local optimal solution, and the modeling behavior of the algorithm involves the following steps:

(i) first, a seasonal constant S is calculated by equation (13)_c

Wherein t is 1, 2, 3;

(ii) checking seasonal monitoring conditions, i.e.

where t and t_mRespectively a current iteration value and a maximum iteration value;

(iii) if the seasonal monitoring condition is true, randomly relocating the squirrels;

s6, randomly moving at the end of winter; the migration of squirrel is expressed by the following equation:

the mathematical method of the Levy distribution is as follows:

where μ, γ > 0.γ is a scale parameter, μ is a transfer parameter, and the Levy flight mode is calculated by:

Lévy(x)＝0.01×r_a×σ/|r_b|^1/β(16)

wherein r is_aAnd r_bIs [0,1 ]]Two normally distributed random numbers in (1), consider β constant to be 1.5 in this work, and σ is calculated as:

wherein Γ (x) ═ x-1! (ii) a

S7, stop criterion: defining an allowed small threshold between the last two continuous results, considering the maximum iteration number as a stopping criterion, storing the current optimal solution, judging whether the stopping criterion is met, outputting the optimal solution if the stopping criterion is met, otherwise, repeating the steps S2-S6 until the algorithm is finished.

Technical Field

The invention relates to a resistivity karst cave identification method based on a squirrel search algorithm, and belongs to the technical field of underground karst cave detection.

Background

At present, in areas where limestone underground karst caves develop, the underground karst caves of construction engineering fields are often located below underground water level, each construction party only needs to invest a large amount of manpower, financial resources and time to carry out exploration work of the underground karst caves of the construction engineering fields of the limestone areas in order to ensure safety of building foundations, traditional engineering geological drilling is only one-hole observation, or traditional electrical exploration is adopted, and due to the fact that the underground karst caves are large in buried depth or small in scale, the resistivity of the underground karst caves and the resistivity of surrounding rocks are small in difference and cannot be accurately explored, and serious potential safety hazards are left for buildings. The underground karst cave can cause certain harm to building construction, subway construction and the like, the underground karst cave is detected and processed in time, and the method has important significance to safety in construction, use and the like.

The resistivity method for detecting the karst cave problem is an optimization problem in the identification research of inverse problems, and the basic idea is as follows: the detection area, the electrode scheme and the karst cave are determined, the measured electric potential distribution data are changed, and the karst cave position can be positioned by using the changes in response to the change of the apparent resistivity of the soil layer. Namely, the method is realized by defining an objective function related to the karst cave structure and then utilizing various optimization methods to detect the position of the karst cave. The traditional resistivity cavern detection method introduces a regularization method, adopts sensitivity or a gradient matrix to identify the position of the cavern, and has the following defects: the method is sensitive to initial values and noise, and is easy to fall into a local optimal solution, so that the identification effect is poor.

Disclosure of Invention

The invention overcomes the defects of the prior art and provides a resistivity karst cave identification method based on a squirrel search algorithm, which has practicability, effectiveness and accuracy.

In order to achieve the purpose, the technical scheme of the invention is as follows: a resistivity karst cave identification method based on a squirrel search algorithm comprises the following steps:

the method comprises the following steps: establishing a finite element model of a karst cave soil layer, determining an electrode scheme, and obtaining potential distribution data of a measuring point by using an electrical impedance method;

step two: constructing an objective function of the karst cave structure, namely an optimized objective function, wherein the objective function is as follows:

where g (c) is the objective function, | () | represents the 2 norm of the vector, i.e.

c is the apparent resistivity of each cell, Su represents the data on the set of measured potential points,

is the potential distribution data measured by the ith set of electrode schemes,

is the potential distribution data identified by the ith group of electrode schemes,

is a matrix formed by s sets of measurement data, and R (c) is a matrix formed by s sets of identification data;

step three: continuously optimizing a target function by using a squirrel search algorithm, and finally obtaining the identification position of the karst cave after meeting a stopping criterion;

the concrete process of optimizing the objective function by utilizing the squirrel search algorithm in the third step of the method is as follows:

s1, random initialization: assuming n squirrels (FS), the position of the ith squirrel can be represented by a vector, and the positions of all squirrels can be represented by the following matrix:

wherein, FS_i.jRepresenting the j-th dimension of the ith squirrel, the initial position of each squirrel in the forest was assigned with a uniform distribution:

FS_i＝FS_L+U(0，1)×(FS_U-FS_L) (3)

wherein, FS_LAnd FS_URespectively the lower bound and the upper bound of the ith squirrel in the j dimension, and U (0,1) is a value range [0,1 ]]Uniformly distributed random numbers inside;

s2, fitness evaluation: the values of the decision variables (solution vectors) are put into a user-defined fitness function, the position fitness of each squirrel is calculated, and the corresponding values are stored in the following array:

s3, sorting, declaring, and randomly selecting: after storing the fitness value of the position of each squirrel, the arrays are sorted in ascending order;

s4, generating a new location: the mathematical model of the new position has the following three cases:

wherein d is_gFor random glide distance, R₁、R₂And R₃Is [0,1 ]]Random within rangeNumber, FS_nt、FS_atAnd FS_htThe squirrel positions of the common tree, the oak tree and the pecan tree are reached respectively, t is the current iteration number, G_cHas a value of 1.9, P_dpThe value of (A) is 0.1;

s5, seasonal monitoring conditions: seasonal monitoring conditions are introduced into the algorithm, so that the algorithm is prevented from falling into a local optimal solution, and the modeling behavior of the algorithm involves the following steps:

(i) first, a seasonal constant S is calculated by equation (13)_c

Wherein t is 1, 2, 3;

(ii) checking seasonal monitoring conditions, i.e.

Wherein S_minFor the minimum of the seasonal constants, the calculation formula is:

where t and t_mRespectively a current iteration value and a maximum iteration value;

(iii) if the seasonal monitoring condition is true, randomly relocating the squirrels;

s6, randomly moving at the end of winter; the migration of squirrel is expressed by the following equation:

the mathematical method of the Levy distribution is as follows:

where μ, γ > 0.γ is a scale parameter, μ is a transfer parameter, and the Levy flight mode is calculated by:

Lévy(x)＝0.01×r_a×σ/|r_b|^1/β(16)

wherein r is_aAnd r_bIs [0,1 ]]Two normally distributed random numbers in (1), considered in this work

β constant is 1.5, σ is calculated as:

wherein Γ (x) ═ x-1! (ii) a

S7, stop criterion: defining an allowed small threshold between the last two continuous results, considering the maximum iteration number as a stopping criterion, storing the current optimal solution, judging whether the stopping criterion is met, outputting the optimal solution if the stopping criterion is met, otherwise, repeating the steps S2-S6 until the algorithm is finished. Compared with the prior art, the technical scheme of the invention has the beneficial effects that: according to the method, the target function is constructed through the potential distribution data of the karst cave structure measuring points, the karst cave position is identified by using the metaheuristic algorithm, and the method is insensitive to the initial value and the noise, is not easy to fall into the local optimal solution, and has better efficiency and precision.

Drawings

FIG. 1 is a flow chart of a solution cavity identification problem normalized to an optimization problem;

FIG. 2 is a schematic diagram of an approximate model of the sliding behavior of the squirrel search algorithm;

FIG. 3 is a schematic diagram of the implementation flow of the squirrel search algorithm;

FIG. 4 is a finite element model of a single karst cave in example 1 of the present invention;

FIG. 5 is a finite element model of a plurality of karsts in embodiment 2 of the present invention;

FIG. 6 shows the result of the detection in the absence of noise by the method of the present invention in example 1 of the present invention;

FIG. 7 shows the results of the detection in the presence of 0.1% noise in example 1 of the present invention using the method of the present invention;

FIG. 8 is a graph showing the detection result of the particle swarm optimization under the noise-free condition in embodiment 1 of the present invention;

FIG. 9 shows the detection result of the particle swarm optimization in the case of 0.1% noise in embodiment 1 of the present invention;

FIG. 10 is a comparison graph of the convergence curves of the algorithm for a single cavern with 0.1% noise in example 1 of the present invention;

FIG. 11 shows the result of the detection in the absence of noise by the method of the present invention in example 2 of the present invention;

FIG. 12 shows the result of the detection in the presence of 0.1% noise in example 2 of the present invention;

FIG. 13 is a graph showing the detection result of the particle swarm optimization in the noise-free case in embodiment 2 of the present invention;

FIG. 14 shows the detection result of the particle swarm optimization in the case of 0.1% noise in embodiment 2 of the present invention;

FIG. 15 is a comparison graph of the convergence curves of the algorithm for multiple cavities with 0.1% noise in example 2 of the present invention.

Detailed Description

In order to make the objects, technical solutions and advantages of the present invention more apparent, the present invention will be further described in detail with reference to the accompanying drawings and examples.

Referring to fig. 1 and 2, the resistivity cavern identification method based on the squirrel search algorithm of the invention comprises the following specific steps:

(1) objective function

The karst cave model is simplified into a rectangular detection area. The current points are not identical to the potential points, but are equal in number. The positive problem solution is carried out by using finite elements, and the following equation is provided

Wherein u ═ u₁；u₂；…；u_n]Including all measurement node potentials, f ═ f₁；f₂；…；f_n]The current vector is characterized, different electrode schemes, and f is different. K is a stiffness matrix and the apparent resistivity c of each unit is ═ c₁；c₂；…；c_m]Is linearly related to where K_jIs a matrix of cell stiffness at unit resistivity.

When no karst cave exists in the soil layer, the apparent resistivity of all the units is taken as c₀When the solution cavity occurs, the apparent resistivity of the unit is reduced, and the apparent resistivity of the unit with the solution cavity is 3, namely c_jThat is, 3( j 1, 2, 3.., m) indicates that the jth unit has a cavity. Through the formula of the finite element, c has a functional relation with u, and the electric potential distribution data calculated by c deduction can be identified

The objective function constructed based on the potential distribution data is as follows:

where c is the apparent resistivity of each cell, Su represents the data on the set of measured potential points,is the potential distribution data measured by the ith set of electrode schemes,

is the potential distribution data identified by the ith group of electrode schemes,

is a matrix of s sets of measured data, and R (c) is a matrix of s sets of calculated data. When the measured potential distribution data and the identified potential distribution data are completely matched, the numerical value is minimum, the karst cave identification problem is equivalent to an optimization problem, a certain individual position of the squirrel is an apparent resistivity distribution condition, and when the objective function reaches a minimum value, the position of the karst cave can be reflected through the identified apparent resistivity c, namely an optimal solution.

(2) And optimizing the target function by using a squirrel search algorithm to obtain an identification result.

Referring to fig. 3, the squirrel search algorithm can be divided into the following 7 stages:

stage 1: random initialization

The parameter setting comprises the initial population number n of the algorithm, the maximum iteration number t and the upper and lower bounds of all dimensions of the squirrel. The position of the ith squirrel can be represented by a vector. The positions of all squirrels can be represented by the following matrix:

wherein FS_i.jRepresents the j-th dimension of the ith squirrel. The initial position of each squirrel in the forest was assigned with a uniform distribution (formula (2)).

FS_i＝FS_L+U(0，1)×(FS_U-FS_L) (4)

Wherein, FS_LAnd FS_URespectively the lower bound and the upper bound of the ith squirrel in the j dimension, and U (0,1) is a value range [0,1 ]]Uniformly distributed random numbers within.

And (2) stage: the degree of fitness is evaluated and the degree of fitness,

the values of the decision variables (solution vectors) are put into a user-defined fitness function, the position fitness of each squirrel is calculated, and the corresponding values are stored in the following array:

the fitness value of the position of each squirrel describes the quality of the food sources it searches, i.e. the optimal food source (pecan tree), the normal food source (oak tree) and no food source (squirrels on normal trees), and therefore also their survival probability.

And (3) stage: ordering, declaring and randomly selecting

After storing the fitness value for each squirrel's position, the arrays are sorted in ascending order. Squirrels with the smallest fitness value are considered to be on walnut trees. The next three best squirrels were considered to be on the oak tree, which were considered to beMoving towards the pecan tree. The remaining squirrels should be on normal trees. Further to random selection, some squirrels were considered to walk towards the pecan tree, assuming they had met their daily energy requirements. The remaining squirrels will go to acorn nut trees (to meet their daily energy requirements). This foraging behavior of squirrels is always affected by predators. Using predator presence probability (P)_dp) The location update mechanism of (a) models this natural behavior.

And (4) stage: generating new locations

As mentioned previously, three conditions may occur during the dynamic foraging of squirrels. In each case, it is assumed that in the absence of predators, the squirrel will glide and effectively find its favorite food throughout the forest, and the presence of predators will make it discreet and be forced to search nearby hidden places using small random walks. The mathematical model of dynamic foraging behavior is as follows:

1) squirrels on acorn nut trees may move towards the walnut tree. In this case, the new position of the squirrel can be obtained as follows:

wherein d is_gFor random glide distance, R₁Is [0,1 ]]Random number in the range, FS_htTo reach the squirrel position of the hickory tree, t is the current iteration number. In the mathematical model, a sliding constant G is used_cThe balance of exploration and development is realized. Its value has a significant impact on the performance of the proposed algorithm. In the present work, G_cThe value of (a) was regarded as 1.9 and was obtained by rigorous analysis.

2) Squirrel in normal tree FS_ntMay move towards acorn nut trees to meet their daily energy requirements. In this case, the new position of the squirrel can be obtained as follows:

wherein R is₂Is [0,1 ]]Random numbers within a range.

3) Some squirrels are on normal trees and have eaten acorn nuts, which may fly towards the walnut tree to store walnuts that may be eaten when food is scarce. In this case, the new position of the squirrel can be obtained as follows:

wherein R is₃Is [0,1 ]]Random numbers within a range. In the present work, the probability P of predator presence_dpIn all cases considered to be 0.1.

4) The gliding mechanism of squirrels is described by equilibrium gliding, in which the sum of the lift (L) and drag (D) produces a resultant force (R) equal and opposite in magnitude to the squirrel's body weight (Mg). The invention uses an approximate model of the sliding behavior in the design of the optimization algorithm (fig. 1). The steady glide speed of the squirrel always falls at an angular level and lift-drag ratio or glide rate, defined as follows:

squirrels can increase their glide landing length through a small glide angle (phi) and thus increase the lift drag ratio. Lift is generated by the downward deflection of air as it passes over the wing, defined herein as:

where ρ (═ 1.204 kgm)^-3) Is the density of air, C_LIs the coefficient of lift, V (═ 5.25 ms)^-1) Is the velocity, S (═ 154 cm)²) Is the limb surface area, and the frictional resistance is defined as:

wherein C is_DIs the coefficient of frictional resistance. At low speeds, this drag coefficient is very high, while at high speeds it becomes smaller. The steady state glide angle is:

approximate glide distance d_gCalculated as follows:

wherein h is_g(═ 8m) is the loss in height after glide. Calculating d_gAll parameter values required, including C_LAnd C_DAre considered from the actual data. Therefore, the squirrel can change the glide path length or d simply by changing the lift-drag ratio according to the landing position_g. In the simulation process, C is added_LAt 0.675 ≦ C_LRandom variation in the range of ≦ 1.5, consider C_DFixed at 0.60.

The horizontal glide distance of the squirrel is usually 5-25 m in one glide. The sliding distance of the model is 9-20 m, d_gLarger values may introduce larger perturbations in generating new positions, resulting in less than ideal algorithm performance. Will d_gIs scaled down to achieve acceptable performance of the algorithm. d_gDivided by a suitable non-zero value, called the scale factor (sf), is obtained by rigorous experimentation with the basis function. In the present work sf 18 provides sufficient d_gPerturbation range, helps achieve the desired balance between exploration and development phases.

And (5) stage: seasonal monitoring conditions

Seasonal variations significantly affect the foraging activity of squirrels. They suffer a great deal of heat loss at low temperatures because they are very hot and small, which makes foraging very costly and risky due to the presence of active predators. Compared to autumn, the climatic conditions force them to be less active in winter. Therefore, the squirrel's movements are affected by weather changes, and the inclusion of such behavior may provide a more realistic approach to optimization. Therefore, seasonal monitoring conditions are introduced into the algorithm, and the algorithm is prevented from falling into a local optimal solution. The modeling behavior involves the following steps:

(i) first, a seasonal constant S is calculated_c

Wherein t is 1, 2, 3.

(ii) Checking seasonal monitoring conditions, i.e.

Wherein S_minFor the minimum of the seasonal constants, the calculation formula is:

where t and t_mRespectively, a current iteration value and a maximum iteration value. S_minAffecting the exploration and development capabilities of the method. S_minThe larger the algorithm, the higher the exploration capacity; s_minThe smaller the algorithm, the more powerful the algorithm can be developed. For any metaheuristic to be effective, there must be a proper balance between these two phases. Albeit by the slip constant G_cTo maintain this balance, but by adaptively changing S during the iteration process_minThe value of (b) can be improved.

(iii) If the seasonal monitoring condition is true, squirrels that cannot find the best winter food source in the forest are randomly relocated.

And 6: winter end random relocation

As previously mentioned, the termination of the winter season causes squirrels to become active due to the lower cost of foraging. Squirrels that survive without finding the best food source in the forest in winter may forage in a new direction. Incorporating this behavior into the model may improve the exploratory capabilities of the algorithm. It is speculated that only those squirrels that survive without being able to find a source of hickory will move in different directions to find a better source of food. This squirrel migration is illustrated by the following equation:

the levy distribution encourages better and more efficient search space exploration. Levy flight patterns are a powerful mathematical tool used by researchers to improve the global search capabilities of various metaheuristic algorithms. The levy flight mode helps to find new candidate solutions, far from the best solution at present. It is a random walk in which the step size is derived from the levy distribution. The distribution expression is often expressed by a power law formula L(s) to | s^-1-βWhere 0 < β < 2 is an indicator the mathematical method for the Levy distribution is as follows:

where μ, γ > 0.γ is the scale parameter and μ is the transfer parameter. The levy flight pattern is calculated by:

Lévy(x)＝0.01×r_a×σ/|r_b|^1/β(16)

wherein r is_aAnd r_bIs [0,1 ]]Two normally distributed random numbers in (1), consider β constant to be 1.5 in this work, and σ is calculated as:

wherein Γ (x) ═ x-1! .

And (7) stage: stopping criterion

Functional tolerance is a common convergence criterion in which an allowable small threshold is defined between the last two consecutive results. The maximum execution time is sometimes used as a stop condition. In this study, the maximum number of iterations is considered the stopping criterion.

As the squirrel search algorithm belongs to a swarm intelligence algorithm, the squirrel search algorithm can be used for solving a nonlinear problem, and the resistivity method for detecting the karst cave problem is a typical nonlinear problem. When a squirrel search algorithm is combined, a finite element model of a karst cave is identified by a resistivity method according to a detection region, the apparent resistivity of the unit is used for calibrating the corresponding karst cave position, the set input and output currents (electrode scheme) are used for solving the potential distribution of the detection area according to the positive problem formula, in practical engineering, we can only obtain the situation of potential distribution, and the problem of non-linearity of the resistivity distribution to be determined by the potential distribution, therefore, the apparent resistivity distribution condition of the detection area is solved by combining a squirrel search algorithm, so that the position of the karst cave is obtained (an objective function is established according to the potential distribution condition obtained by the positive problem, a certain body position of the squirrel is the apparent resistivity distribution condition, parameters such as the initial position, the number of individuals and the like of the squirrel are set, and then the optimal squirrel position finally obtained is the identified apparent resistivity distribution condition after the stop criterion is met through continuous iteration of the algorithm).