Concrete three-dimensional mesoscopic modeling method based on simulated annealing algorithm

Abstract

The invention relates to the technical field of concrete material analogue simulation, in particular to a concrete three-dimensional mesoscopic modeling method based on a simulated annealing algorithm. All aggregate particles are generated according to aggregate assembly, all the aggregate particles are initially put into the model space, and aggregate overlapping is allowed; the overlapping ratio between the aggregates is simulated as extrusion deformation, so that the aggregate putting problem is converted into the optimization problem of the elastic potential energy of the aggregate configuration; and dynamically adjusting the spatial position of the aggregate according to the mutual repulsive force among the aggregate particles, gradually eliminating the overlap among the aggregate through variable step iteration, and finally realizing the generation of the concrete three-dimensional mesoscopic model. The aggregate shape is implicitly represented through support mapping, wide aggregate shapes can be described, unification of collision detection methods of aggregates in different shapes can be realized, and the collision detection efficiency is improved through wide-phase and narrow-phase two-stage collision detection. The problems of single aggregate shape, low volume fraction and long consumed time in concrete three-dimensional mesoscopic modeling are effectively solved.

Description

Technical Field

[0001] The invention relates to the technical field of concrete material simulation technology, and in particular to a concrete three-dimensional microscopic modeling method based on a simulated annealing algorithm. Background Art

[0002] From a microscopic perspective, concrete is a multiphase composite material consisting of aggregate, mortar, and interfacial transition zone (ITZ). The complex mechanical properties of concrete depend largely on its heterogeneous mesostructural characteristics. Therefore, numerical simulation methods based on concrete mesoscopic models have gradually become an important means of concrete material research.

[0003] Concrete mesoscopic modeling mainly involves the generation of random spatial structure of aggregates. Random order placement method is a classic mesoscopic modeling method that can accurately control the gradation and volume fraction of aggregates. This method places aggregates into the target space one by one according to particle size from large to small, and ensures that the current aggregate does not overlap with the previously placed aggregates during the placement process. However, as the space occupied by the placed aggregates increases, the effective placement area of ​​the newly added aggregates will decrease, resulting in very inefficient placement in the later stage. In order to improve the placement efficiency, the Chinese patent "A three-dimensional concrete mesoscopic simulation model and its establishment method (CN102136157A)" proposes a method for establishing a three-dimensional concrete mesoscopic model using a hexahedral background grid, and the space occupied by the placed aggregates is marked by the background grid to improve the success rate of aggregate placement. In the Chinese patent "A method for constructing a concrete mesoscopic model (CN117195382A)", the finite element grid is used as a global background grid for block division, and the complexity and calculation amount of the construction are reduced by parallel placement of blocks and secondary subdivision of blocks. In the Chinese patent "Mesoscopic stochastic modeling method of four-graded concrete based on four-fold discrete point set (CN115798645A)", the four-fold discrete point set structure is used to reduce the number of attempts to place aggregates and increase the upper limit of aggregate content. However, in these improved methods, the positions of the placed aggregates remain unchanged. When the spatial position of the previously placed aggregates is inappropriate, it will seriously affect the placement of the remaining aggregates, making it difficult to generate a mesoscopic model with a high aggregate volume fraction. In addition, in these modeling methods, the collision detection between aggregates is mostly global, which makes the placement process still time-consuming, and the collision detection algorithm for aggregates is often only applicable to a single aggregate shape. Summary of the invention

[0004] In order to overcome the problems of single applicable aggregate shape, low aggregate placement efficiency, low aggregate volume fraction and the like in existing modeling methods, the present invention provides a three-dimensional concrete mesoscopic modeling method based on a simulated annealing algorithm, which dynamically adjusts the spatial position of aggregates according to the mutual repulsive force between aggregate particles, iteratively eliminates the overlap between aggregates, and thus completes the establishment of a concrete mesoscopic model.

[0005] The present invention is realized by the following technical scheme: a concrete three-dimensional microscopic modeling method based on simulated annealing algorithm, comprising:

[0006] Step S1, pre-generate all aggregate particles according to the mesoscopic model size, aggregate gradation and aggregate volume fraction;

[0007] Step S2, generating an initial aggregate configuration S=S0, placing aggregate particles in order from large to small particle sizes, placing the selected aggregate into the model space at random center coordinates and rotation angles, performing collision detection between the current aggregate and the previously placed aggregate and recording collision information, calculating the elastic potential energy caused by the overlap between the current aggregate and the surrounding aggregates based on the collision information of the current aggregate, repeating several times and taking the position with the minimum potential energy as the final placement position of the current aggregate, and updating the collision information of the aggregate and the aggregates overlapping with it;

[0008] Step S3, calculating the elastic potential energy U0 of the initial aggregate configuration S0;

[0009] Step S4, setting the initial values ​​of the simulated annealing algorithm control parameters, including the aggregate translation length L = L0 and its reduction coefficient α, the aggregate rotation angle θ = θ0 and its reduction coefficient β, the sampling number M = M0 and its magnification coefficient γ, the historical minimum potential energy U of the aggregate configuration min =U0, sampling counter n=0;

[0010] Step S5, updating the position of each aggregate in order from large to small particle size to update the aggregate configuration S;

[0011] Step S6, calculating the elastic potential energy U of the current aggregate configuration S;

[0012] Step S7, under the condition of U=0, output the shape and spatial position information of each aggregate in the current aggregate configuration to complete the three-dimensional microscopic modeling of concrete; if U≠0, proceed to step S8;

[0013] Step S8: Compare the elastic potential energy U of the current aggregate configuration with the historical minimum potential energy U min Compare, when U<U min Update U under the condition min = U and reset the sampling counter to n = 0, and when U ≥ U min Under the condition of , the sampling counter is reset to n = n + 1;

[0014] Step S9, based on the comparison between n and M, determine whether to terminate the current sampling, update the parameters L=αL, θ=βθ, M=γM under the condition that n=M and reset the sampling counter n=0 to cool down, and perform step S5 under the condition that n≠M;

[0015] Step S10, based on L and L min The comparison of L<L min Under the condition of L≥L min Under the condition of, it is determined to jump to step S5 to perform the next sampling process, where L min The preset minimum translation length.

[0016] Furthermore, the step S5 comprises:

[0017] Step S501, calculate the aggregate P to be moved currently according to the following formula: i The resultant force of repulsion,

[0018]

[0019] Where V j Aggregate P j The volume, d ij Indicates that from aggregate P j Point to aggregate P i The penetration depth vector, N represents the number of aggregates;

[0020] Step S502, according to the aggregate P i The resultant force F i Calculation of aggregate P i The new position after moving i ′, where si′=R(si)+D, where s i ′ and s i P i The new position and current position of , R(si) represents the current position s i Perform rotation transformation according to the rotation matrix R; D represents the translation vector of the aggregate position, D = λ t Lu, where u·F i ≥0, where u is a random unit vector, λ t is a random number between 0 and 1; the rotation matrix R = I + sin (λ r θ)M+(1-cos(λ r θ))M 2 , where I represents the unit matrix, λ r is a random number in the range of [-0.5,0.5], and M is a random rotation axis k (k x ,k y ,k z ),

[0021]

[0022] Step S503, check the boundary of the aggregate. If the aggregate moves beyond the model space after moving to the new position, keep the original position of the aggregate unchanged and start the next aggregate P. i+1 If the aggregate does not exceed the model space range after moving to the new position, proceed to step S504;

[0023] Step S504: i Perform collision detection, update the collision information of the aggregates that collide, and then start the next aggregate P i+1 Location updates;

[0024] Step S505: After all aggregates complete the updating of aggregate positions in sequence, an updated aggregate configuration S is obtained.

[0025] Furthermore, the collision detection is based on the support map representation of the aggregate shape and traverses the minimum gap δ between aggregates. a and the thickness of the interface transition zone t, where the support map of the aggregate shape

[0026]

[0027] In the formula, s A (v) represents the original support map of aggregate A, v is an arbitrary direction, s′ A (v) is the support map after the minimum gap between aggregates and the thickness of the interface transition zone, δ a is the minimum gap between aggregates, and t is the thickness of the interface transition zone.

[0028] Furthermore, the collision detection includes wide-phase collision detection and narrow-phase collision detection, wherein the wide-phase collision detection uses the axis-aligned bounding box of the aggregate for interference detection, and the narrow-phase collision detection uses the GJK algorithm to detect and determine the penetration depth vectors of the two aggregates when two aggregates intersect, and updates the collision information of the two aggregates respectively.

[0029] Furthermore, the step S6 comprises:

[0030] Step S601: The overlap between aggregates is regarded as extrusion deformation. If the aggregate P i and aggregate P j There is overlap, by P j The effect of P i The elastic potential energy u ij Calculated by the following formula,

[0031]

[0032] Step S602: aggregate P i Elastic potential energy U iThe sum of contributions from all surrounding aggregates, which is calculated by the following formula,

[0033]

[0034] Step S603, the elastic potential energy U of the current entire aggregate configuration is the sum of the elastic potential energies of all aggregates, which is calculated by the following formula:

[0035] The beneficial effects of the present invention are:

[0036] (1) The aggregate position is dynamically adjusted through the simulated annealing algorithm, which improves the aggregate placement efficiency and the upper limit of aggregate volume fraction.

[0037] (2) The efficiency of aggregate collision detection is improved through two-level collision detection: wide phase and narrow phase.

[0038] (3) The support mapping is used to implicitly represent the aggregate shape, which can describe a wide range of aggregate shapes and achieve the unification of collision detection methods for aggregates of different shapes. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Figure 1 It is a flow chart of a concrete three-dimensional mesoscopic modeling method based on a simulated annealing algorithm according to an embodiment of the present invention;

[0040] Figure 2 It is a schematic diagram of the repulsive force caused by overlapping aggregates in an embodiment of the present invention;

[0041] Figure 3A A cross-sectional view of a three-dimensional mesoscopic model of concrete with spherical aggregates according to an embodiment of the present invention;

[0042] Figure 3B A cross-sectional view of a three-dimensional mesoscopic model of concrete with ellipsoidal aggregate according to an embodiment of the present invention;

[0043] Figure 3C A cross-sectional view of a concrete three-dimensional mesoscopic model of polyhedral aggregate according to an embodiment of the present invention;

[0044] Figure 3D This is a cross-sectional view of a three-dimensional mesoscopic model of concrete with mixed-shaped aggregates according to an embodiment of the present invention. DETAILED DESCRIPTION

[0045] The implementation process of the present invention is described in detail below in conjunction with the accompanying drawings and specific embodiments. The following embodiment is only a specific application example of the present invention, but the protection scope of the present invention is not limited to this embodiment.

[0046] like Figure 1 As shown, a concrete three-dimensional microscopic modeling method based on a simulated annealing algorithm comprises:

[0047] Step S1, pre-generate all aggregate particles according to the mesoscopic model size, aggregate gradation and aggregate volume fraction;

[0048] In this embodiment, the size of the mesoscopic model is set as a cube with a side length of 50 mm; spherical aggregate, ellipsoidal aggregate, polyhedral aggregate and these three types of mixed aggregates are used to establish concrete mesoscopic models respectively to verify the modeling method proposed in the present invention; the aggregate volume content is taken as 50%, the aggregate particle size range is taken as 2.36 mm to 19 mm, and the aggregate particle size distribution adopts Fuller grading.

[0049] Step S2, generating an initial aggregate configuration S=S0, placing aggregate particles in order from large to small in particle size, placing the selected aggregate into the model space at random center coordinates and rotation angles, performing collision detection between the current aggregate and the previously placed aggregate and recording collision information, calculating the elastic potential energy caused by the overlap between the current aggregate and the surrounding aggregates according to the collision information of the current aggregate, repeating several times and taking the position with the minimum potential energy as the final placement position of the current aggregate, and updating the collision information of the aggregate and the aggregates overlapping with it;

[0050] In this embodiment, 10 trials were conducted for each aggregate;

[0051] Step S3, calculating the elastic potential energy U0 of the initial aggregate configuration S0;

[0052] Step S4, setting the initial values ​​of the simulated annealing algorithm control parameters, including the aggregate translation length L = L0 and its reduction coefficient α, the aggregate rotation angle θ = θ0 and its reduction coefficient β, the sampling number M = M0 and its magnification coefficient γ, the historical minimum potential energy U of the aggregate configuration min =U0, sampling counter n=0;

[0053] In this embodiment, L0=2.5 ​​mm, θ0=90°, M0=10, α=0.8, β=0.8, and γ=1.2 are set.

[0054] Step S5, updating the position of each aggregate in order from large to small particle size to update the aggregate configuration S;

[0055] Step S6, calculating the elastic potential energy U of the current aggregate configuration S;

[0056] Step S7, under the condition of U=0, output the shape and spatial position information of each aggregate in the current aggregate configuration to complete the three-dimensional microscopic modeling of concrete; if U≠0, proceed to step S8;

[0057] Step S8: Compare the elastic potential energy U of the current aggregate configuration with the historical minimum potential energy U min Compare, when U<U min Update U under the conditionmin = U and reset the sampling counter to n = 0, and when U ≥ U min Under the condition of , the sampling counter is reset to n = n + 1;

[0058] Step S9, based on the comparison between n and M, determine whether to terminate the current sampling, update the parameters L=αL, θ=βθ, M=γM under the condition that n=M and reset the sampling counter n=0 to cool down, and perform step S5 under the condition that n≠M;

[0059] Step S10, based on L and L min The comparison of L<L min Under the condition of L≥L min Under the condition of, it is determined to jump to step S5 to perform the next sampling process, where L min To preset the minimum translation length, set L min =0.01mm.

[0060] Specifically, step S5 includes:

[0061] Step S501, calculate the aggregate P to be moved currently according to the following formula: i The resultant force of repulsion,

[0062]

[0063] Where V j Aggregate P j The volume, d ij Indicates that from aggregate P j Point to aggregate P i The penetration depth vector, N represents the number of aggregates;

[0064] Step S502, according to the aggregate P i The resultant force F i Calculation of aggregate P i The new position after moving i ′, where si′=R(si)+D, where s i ′ and s i P i The new position and current position of , R(si) represents the current position s i Perform rotation transformation according to the rotation matrix R; D represents the translation vector of the aggregate position, D = λ t Lu, where u·F i ≥0, where u is a random unit vector, λ t is a random number between 0 and 1; the rotation matrix R = I + sin (λ rθ)M+(1-cos(λ r θ))M 2 , where I represents the unit matrix, λ r is a random number in the range of [-0.5,0.5], and M is a random rotation axis k (k x ,k y ,k z ),

[0065]

[0066] Step S503, check the boundary of the aggregate. If the aggregate moves beyond the model space after moving to the new position, keep the original position of the aggregate unchanged and start the next aggregate P. i+1 If the aggregate does not exceed the model space range after moving to the new position, proceed to step S504;

[0067] Boundary checking can be implemented using the axis-aligned bounding box (AABB) of the aggregate. Consider the model space of a cube, whose length, width, and height are represented by X, Y, and Z respectively. First, according to the support mapping of the aggregate, a total of 6 support mapping points along the positive and negative directions of the three coordinate axes of the Cartesian coordinate system are calculated. According to the maximum values ​​of the coordinates (x, y, z) of these points, the AABB of the aggregate can be obtained, which is expressed as [x min ,x max ]×[y min ,y max ]×[z min ,z max ]. Considering that there is a minimum distance δ between the aggregate and the model boundary b If the aggregate is in the model space, the following conditions must be met:

[0068]

[0069] Where: δ0 = δb-(t+δa / 2); if not all of them are met, the aggregate exceeds the model boundary.

[0070] Step S504: i Perform collision detection, update the collision information of the aggregates that collide, and then start the next aggregate P i+1 Location updates;

[0071] Step S505: After all aggregates complete the updating of aggregate positions in sequence, an updated aggregate configuration S is obtained.

[0072] Specifically, the collision detection is based on the support map representation of the aggregate shape and traverses the minimum gap δ between aggregates. a and the thickness of the interface transition zone t, where the support map of the aggregate shape

[0073]

[0074] In the formula, s A (v) represents the original support map of aggregate A, v is an arbitrary direction, s′ A (v) is the support map after the minimum gap between aggregates and the thickness of the interface transition zone, δ a is the minimum gap between aggregates, t is the thickness of the interface transition zone, and the present invention uses support mapping to represent the shape of aggregates, and implements aggregate collision detection based on this. In addition, the minimum gap between aggregates δ is also considered in the support mapping. a and the thickness of the interface transition zone t (ITZ), the support mapping diagram of the aggregate shape is shown in Figure 2 As shown, in this embodiment, δ a and t are both set to 0.08 mm.

[0075] Specifically, the collision detection includes wide-phase collision detection and narrow-phase collision detection. The wide-phase collision detection uses the axis-aligned bounding box of the aggregate for interference detection, and the narrow-phase collision detection uses the GJK algorithm to detect and determine the penetration depth vectors of the two aggregates when two aggregates intersect, and updates the collision information of the two aggregates respectively.

[0076] Specifically, step S6 includes:

[0077] Step S601: The overlap between aggregates is regarded as extrusion deformation. If the aggregate P i and aggregate P j There is overlap, by P j The effect of P i The elastic potential energy u ij Calculated by the following formula,

[0078]

[0079] Step S602: aggregate P i Elastic potential energy U i The sum of contributions from all surrounding aggregates, which is calculated by the following formula,

[0080]

[0081] Step S603, the elastic potential energy U of the current entire aggregate configuration is the sum of the elastic potential energies of all aggregates, which is calculated by the following formula:

[0082] Using the concrete three-dimensional microscopic modeling method of the present invention, this embodiment establishes concrete microscopic models for four different aggregate shapes, wherein the spherical three-dimensional view and the cross-sectional view are as follows: Figure 3A, the 3D view and the cross-sectional view of the ellipsoid are as follows Figure 3B , the three-dimensional view and the cross-sectional view of the polyhedron are as follows Figure 3C , the 3D view and cutaway view of the blended shape are as follows Figure 3D The number of aggregate particles and the time required for aggregate placement in the four mesoscopic models are shown in Table 1:

[0083] Table 1 Aggregate quantity and calculation time of four concrete mesoscopic models

[0084]

Claims

1. A three-dimensional concrete micro-modeling method based on simulated annealing algorithm, characterized in that: include: Step S1, pre-generate all aggregate particles according to the mesoscopic model size, aggregate gradation and aggregate volume fraction; Step S2, generating an initial aggregate configuration S=S0, placing aggregate particles in order from large to small in particle size, placing the selected aggregate into the model space at random center coordinates and rotation angles, performing collision detection between the current aggregate and the previously placed aggregate and recording collision information, calculating the elastic potential energy caused by the overlap between the current aggregate and the surrounding aggregates according to the collision information of the current aggregate, repeating several times and taking the position with the minimum potential energy as the final placement position of the current aggregate, and updating the collision information of the aggregate and the aggregates overlapping with it; Step S3, calculating the elastic potential energy U0 of the initial aggregate configuration S0; Step S4, setting the initial values ​​of the simulated annealing algorithm control parameters, including the aggregate translation length L = L0 and its reduction coefficient α, the aggregate rotation angle θ = θ0 and its reduction coefficient β, the sampling number M = M0 and its magnification coefficient γ, the historical minimum potential energy U of the aggregate configuration min =U0, sampling counter n=0; Step S5, updating the position of each aggregate in order from large to small particle size to update the aggregate configuration S; Step S6, calculating the elastic potential energy U of the current aggregate configuration S; Step S7, under the condition of U=0, output the shape and spatial position information of each aggregate in the current aggregate configuration to complete the three-dimensional microscopic modeling of concrete; If U≠0, proceed to step S8; Step S8: Compare the elastic potential energy U of the current aggregate configuration with the historical minimum potential energy U min Compare, when U<U min Update U under the condition min = U and reset the sampling counter to n = 0, and when U ≥ U min Under the condition of , the sampling counter is reset to n = n + 1; Step S9, based on the comparison between n and M, determine whether to terminate the current sampling, update the parameters L=αL, θ=βθ, M=γM under the condition that n=M and reset the sampling counter n=0 to cool down, and perform step S5 under the condition that n≠M; Step S10, based on L and L min The comparison of L<L min Under the condition of L≥L min Under the condition of, it is determined to jump to step S5 to perform the next sampling process, where L min The preset minimum translation length.

2. The method for three-dimensional concrete micro-modeling based on simulated annealing algorithm according to claim 1 is characterized in that: The step S5 comprises: Step S501, calculate the aggregate P to be moved currently according to the following formula: i The resultant force of repulsion, Where V j Aggregate P j The volume, d ij Indicates that from aggregate P j Point to aggregate P i The penetration depth vector, N represents the number of aggregates; Step S502, according to the aggregate P i The resultant force F i Calculation of aggregate P i The new position after moving i ′, where si′=R(si)+D, where s i ′ and s i P i The new position and current position of , R(si) represents the current position s i Perform rotation transformation according to the rotation matrix R; D represents the translation vector of the aggregate position, D = λ t Lu, where u·F i ≥0, where u is a random unit vector, λ t is a random number between 0 and 1; the rotation matrix R = I + sin (λ r θ)M+(1-cos(λ r θ))M 2 , where I represents the unit matrix, λ r is a random number in the range of [-0.5,0.5], and M is a random rotation axis k (k x ,k y ,k z ), Step S503, check the boundary of the aggregate. If the aggregate moves beyond the model space after moving to the new position, keep the original position of the aggregate unchanged and start the next aggregate P. i+1 If the aggregate does not exceed the model space range after moving to the new position, proceed to step S504; Step S504: i Perform collision detection, update the collision information of the aggregates that collide, and then start the next aggregate P i+1 Location updates; Step S505: After all aggregates complete the updating of aggregate positions in sequence, an updated aggregate configuration S is obtained.

3. The method for three-dimensional concrete micro-modeling based on simulated annealing algorithm according to any one of claims 1 or 2, characterized in that: The collision detection is based on the support map representation of the aggregate shape and traverses the minimum gap δ between aggregates. a and the thickness of the interface transition zone t, where the support map of the aggregate shape In the formula, s A (v) represents the original support map of aggregate A, v is an arbitrary direction, s′ A (v) is the support map after the minimum gap between aggregates and the thickness of the interface transition zone, δ a is the minimum gap between aggregates, and t is the thickness of the interface transition zone.

4. The method for three-dimensional concrete micro-modeling based on simulated annealing algorithm according to claim 3 is characterized in that: The collision detection includes wide-phase collision detection and narrow-phase collision detection, wherein the wide-phase collision detection uses the axis-aligned bounding box of the aggregate to perform interference detection, and the narrow-phase collision detection uses the GJK algorithm to detect and determine the penetration depth vectors of the two aggregates when two aggregates intersect, and updates the collision information of the two aggregates respectively.

5. The method for three-dimensional concrete micro-modeling based on simulated annealing algorithm according to claim 1, characterized in that: The step S6 comprises: Step S601: The overlap between aggregates is regarded as extrusion deformation. If the aggregate P i and aggregate P j There is overlap, by P j The effect of P i The elastic potential energy u ij Calculated by the following formula, Step S602: aggregate P i Elastic potential energy U i The sum of contributions from all surrounding aggregates, which is calculated by the following formula, Step S603, the elastic potential energy U of the current entire aggregate configuration is the sum of the elastic potential energies of all aggregates, which is calculated by the following formula:

Images