Multi-physics coupling simulation method and system for reverse pouring of tunnel arch wall

Abstract

The invention provides a multi-physics coupling simulation method and system for reverse pouring of a tunnel arch wall. The method comprises the steps that a two-dimensional computational domain is constructed and discretized; a Herschel-Bulkley model containing a thixotropic recovery mechanism and a first-order differential equation are adopted to describe the shear stress and thixotropic evolution of the concrete; based on an Arrhenius model and a heat conduction equation, establishing a coupling relationship between the temperature, the hydration degree and the rheological parameters; solving a flow field by an explicit projection method; a VOF-PLIC algorithm is adopted to track a concrete-air interface; calculating template side pressure by combining a flow field and stress, and outputting a distribution curve and an extreme value; dynamically controlling the step length according to a CFL criterion and circulating to a set height or terminating; the system comprises a modeling module, a material model module, a thermal-chemical coupling module, a flow field solving module, an interface tracking module, a side pressure evaluation module and a time control module. According to the method, the multi-field coupling effect is integrated, and the accuracy of prediction and evaluation of the whole reverse pouring process is improved.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of tunnel construction simulation, in particular to a multi-physical field coupling simulation method and system for reverse pouring of tunnel arch wall. BACKGROUND

[0002] In the construction of tunnel lining concrete, the traditional arch wall pouring method is mostly from top to bottom or segmented pouring, and the related process and supporting design method are relatively mature, and a perfect construction specification has been formed in engineering practice. However, with the increasingly complex construction environment, under the conditions of limited terrain, narrow site or special working conditions, a reverse rapid pouring method advancing from arch foot to arch top has been gradually popularized in recent years. This method can significantly reduce high-altitude operation and shorten the construction period, but since the pouring direction is opposite to the traditional method, the concrete flow path, interface advancing mechanism and temperature distribution law are essentially different from the conventional working conditions, which brings new challenges to construction control and structural safety.

[0003] At the same time, the existing reverse pouring construction method lacks professional simulation means that can truly reflect the multi-physical field coupling effect of concrete in the whole process of reverse pouring. This leads to the following problems in engineering practice: 1. The compactness of concrete in the reverse pouring process is difficult to effectively control, and voids and cold joint defects are prone to occur; 2. The formwork side pressure cannot be accurately predicted, resulting in overly conservative or insufficient formwork structure design; 3. The early age temperature field prediction accuracy is low, and the temperature crack risk is difficult to evaluate and control in advance; 4. The free liquid surface advancing form is difficult to predict, and the construction integrity and continuity cannot be effectively evaluated.

[0004] The existing simulation technical means have obvious deficiencies in multi-physical field synchronous coupling, construction dynamic process restoration and key construction parameter prediction accuracy, and cannot directly serve the construction optimization and risk prevention and control under the condition of reverse pouring. For example, the traditional free liquid surface tracking method, such as VOF or Level-set, is mostly applied to general flow conditions and does not consider the time-varying rheological parameters, so the accuracy of the simulation cannot be guaranteed; the conventional CAE fluid simulation process lacks integrated expression of construction dynamic process and multi-physical field coupling effect, so the simulation prediction accuracy cannot be guaranteed.

[0005] Therefore, there is an urgent need for a multi-physical field coupling simulation method and system for reverse pouring of tunnel arch wall, which can comprehensively consider the multi-field coupling effect of concrete and improve the accuracy of prediction and evaluation of the whole process of reverse pouring. SUMMARY

[0006] The application aims to provide a multi-physical field coupling simulation method and system for reverse pouring of a tunnel arch wall, aiming to solve the technical problems of lack of multi-physical field coupling simulation in reverse pouring, inaccurate prediction of density, lateral pressure, temperature control and interface.

[0007] To achieve the above-mentioned purpose, in a first aspect, the application provides a multi-physical field coupling simulation method for reverse pouring of a tunnel arch wall, comprising the following steps: S1, obtaining geometric parameters of the tunnel arch wall, constructing a two-dimensional arch wall elevation calculation domain based on the geometric parameters, discretizing the calculation domain and setting boundary conditions; preferably, MAC staggered grid is used for discretization; S2, defining a concrete rheological constitutive and thixotropic evolution model: using a Herschel-Bulkley model with a thixotropic recovery mechanism to describe the shear stress behavior, and using a first-order differential equation to represent the recovery and destruction process of the thixotropic structure parameter; S3, establishing a coupled calculation model of temperature field and hydration kinetics: using an Arrhenius type exponential model to describe the hydration degree evolution, and combining a heat conduction equation to consider the hydration heat release effect, to establish the coupling relationship among temperature, hydration degree and rheological parameters; S4, solving the flow field by using an explicit projection method; S5, using a geometric reconstruction method (Piecewise Linear Interface Calculation, PLIC) interface tracking algorithm based on a Volume of Fluid (VOF) method to solve the space-time evolution process of the concrete-air free interface; S6, calculating the template lateral pressure according to the flow field and stress distribution, and outputting the change curve and maximum value information of the lateral pressure with height; S7, dynamically controlling the time step based on the CFL (Courant-Friedrichs-Lewy) stability criterion, and cyclically executing steps S2 to S7 until the pouring simulation reaches the set height or the termination time.

[0008] As a further improvement of the above scheme, in step S1, the geometric parameters of the tunnel arch wall include the arch foot position and spacing, the arch wall height H, the analysis width L along the wall direction, the lining thickness b, the pouring port center coordinates and their size and number, the designed pouring height, the segmented pouring sequence and the single pouring amount.

[0009] As a further improvement of the above scheme, in step S1, in the MAC grid structure, the pressure p, the temperature T, the hydration degree , the thixotropic structure parameter and the volume fraction are arranged at the center of the grid element, and the horizontal velocity The component is arranged at the center of the vertical plane of the unit, and the vertical velocity The component is arranged at the center of the horizontal plane of the unit.

[0010] As a further improvement of the above scheme, in step S2, the shear stress model is as follows: ; Wherein: is the shear stress, unit Pa; is the shear rate, calculated from the velocity gradient tensor; is the yield stress related to the thixotropic structure parameter, unit Pa; is the consistency coefficient (viscosity coefficient) related to the thixotropic structure parameter, unit Pa·S n ; is the flow index is a shear thinning fluid, is a Bingham type, is a shear thickening fluid); is the thixotropic structure parameter, indicates that the material structure is completely restored, and the viscosity is high, indicates that the structure is completely destroyed, and the fluidity is the highest.

[0011] As a further improvement of the above scheme, the yield stress and the consistency coefficient are respectively represented as follows: ; ; Wherein are all material constants obtained by test calibration.

[0012] As a further improvement of the above scheme, in step S2, the evolution of the thixotropic structure parameter obeys the following equation: ; Wherein: is the structure recovery coefficient, unit , indicating the speed of structure recovery to under no shear action; is the shear damage coefficient, unit s, indicating the sensitivity of shear rate to structure damage;

[0013] ​The above model can reflect the characteristics of the arch foot area under the reverse pouring working condition, i.e., strong shear, structure destruction, and strong flowability, and the upper part of the arch wall has weak shear, structure recovery, and decreased flowability.

[0014] As a further improvement of the above scheme, in step S3, the hydration kinetics equation is as follows: ; Wherein: is the current hydration degree, dimensionless; is the reference hydration rate constant before temperature correction, unit ; is the apparent activation energy of hydration reaction, unit J / mol; is the gas constant, unit J / (mol·K); is the absolute temperature, unit K; is the final value of hydration, such as 0.8-0.9; The temperature field satisfies the heat conduction and hydration heat coupling control equation, which is specifically as follows: ; Wherein: is the density of concrete, unit ; is the specific heat capacity of concrete, unit J / (kg·K); is the effective thermal conductivity, unit W / (m·K); is the hydration heat term, unit , which can be written as , wherein is the hydration heat potential per unit mass, unit J / kg.

[0015] As a further improvement of the above scheme, in step S3, the coupling relationship of the temperature, hydration degree and rheological parameter includes: Temperature rise leads to the decrease of concrete viscosity, which can be corrected by the following formula: ; Wherein is the temperature sensitivity coefficient, is the reference temperature; Hydration degree rise makes the structure harden and the yield stress increase, which can be corrected by the following formula: ; Wherein Hydration enhancement coefficient; Through the above modification, a strong coupling link of “temperature-hydration-rheology” is formed, and comprehensive reflection on the performance evolution of early age is realized.

[0016] As a further improvement of the above scheme, in step S4, the step method for solving the flow field by using the explicit projection method is as follows: S41, predict the velocity field: Under the premise of knowing the velocity of the previous time step , apparent viscosity , and shear stress , the predicted velocity without pressure term is calculated: ; Wherein: is the time step, in seconds; is the spatial convection term, which is preferably discretized by using the first-order upwind difference or the second-order MUSCL format; is the viscous stress tensor, which is given by the Herschel-Bulkley model, in Pa; is the gravity acceleration vector, in ; S42, solve the pressure Poisson equation: In order to meet the incompressible condition , the pressure Poisson equation is constructed as follows: ; Wherein is the Laplace operator, which is preferably discretized by using the five-point difference format on the MAC grid, and the pressure unknown is solved by Jacobi or Gauss-Seidel iteration.

[0017] S42, velocity correction: After obtaining , the predicted velocity is projected and corrected according to the following formula: ; Wherein, the central difference is used for calculation, so that the corrected velocity field satisfies the mass conservation condition.

[0018] As a further improvement of the above scheme, in step S5, the step method for solving the time-space evolution process of the concrete-air free interface is as follows: S51, the interface transport equation is as follows: ; wherein is the volume fraction of the concrete phase, is the velocity field of the current time step; S52, at each time step, the volume fraction flux is calculated along each grid face, the volume fraction of the cell is updated, and the PLIC method is adopted for the cell of , the linear interface is reconstructed according to the local gradient , and the geometric position and shape of the free surface are obtained; S53, the average height of the free surface and the local rising velocity are obtained from the distribution of , the time when the left and right arch foot jets converge in the middle and the interface shape are obtained, and the volume ratio of the concrete in the wall is obtained.

[0019] As a further improvement of the above scheme, in step S6, the template side pressure is calculated according to the following formula: ; wherein, is the static pressure obtained in the cell, with the unit of Pa; is the viscous stress tensor; is the outward normal vector of the template; is the stress tensor component; is the equivalent side pressure acting on the template; The distribution of is extracted along the height direction of the arch wall, the curve of the template side pressure with the height is drawn, and the maximum side pressure value and its corresponding height are output, which provides quantitative basis for template design, support system configuration and construction safety evaluation.

[0020] As a further improvement of the above scheme, in step S7, in order to ensure the numerical stability of the explicit time marching process, the time step is dynamically controlled according to the CFL condition, which is specifically shown as follows: ; wherein: is the Courant number safety factor, usually taken as 0.2-0.5; is the grid size; is the maximum absolute value of the velocity component of the current time step; is the maximum value of the global apparent viscosity.

[0021] ​In a second aspect, the present application also provides a multi-physical field coupling simulation system for reverse pouring of a tunnel arch wall, which adopts the multi-physical field coupling simulation method for reverse pouring of a tunnel arch wall provided in the first aspect, and comprises: a modeling module for establishing a two-dimensional calculation domain according to tunnel geometric parameters and discretizing the calculation domain and applying boundary conditions; a material model module for implementing a Herschel-Bulkley rheological model with thixotropic recovery and a thixotropic evolution equation; a thermal-chemical coupling module for calculating the interaction of hydration reaction and temperature field and the influence of the interaction on rheological properties; a flow field solving module for performing explicit solving of a velocity field and a pressure field based on a projection method; an interface tracking module for simulating the change of a concrete-air interface by using a VOF-PLIC method; a side pressure evaluation module for calculating the stress distribution of a formwork and outputting a pressure curve and an extreme value;

[0022] a time control module for dynamically adjusting a simulation step length according to a CFL condition and driving the modules to cooperatively operate until the pouring process simulation is completed.

[0023] Due to the above technical solutions, the present application has the following beneficial effects: The present application provides a multi-physical field coupling simulation method for reverse pouring of a tunnel arch wall, which integrates a Herschel-Bulkley rheological model with a thixotropic recovery mechanism, a VOF-PLIC interface tracking algorithm, a temperature-hydration dynamics coupling model and formwork side pressure calculation in the same explicit time advancing framework, and constructs a CAE calculation process specially used for arch wall reverse rapid pouring conditions. This CAE calculation process realizes closed-loop simulation from geometric modeling, material constitutive definition, field coupling solving to construction parameter output, can systematically reflect the interaction of flow field, temperature field and stress field in the reverse pouring process, and thus makes up for the lack of multi-physical field coupling simulation in the prior art.

[0024] Specifically, first, the present application can improve the prediction accuracy of fresh concrete fluidity and compactness. Since the Herschel-Bulkley model with a thixotropic recovery mechanism is adopted, and a first-order differential equation describing the recovery and destruction of structural parameters is introduced into the model, the spatiotemporal reconstruction behavior of thixotropic structure can be depicted according to different conditions of high shear at the arch foot and low shear on the upper segment of the arch wall. In combination with the evolution information of the concrete-air free interface obtained by the VOF-PLIC interface tracking (step S5), the flow front, filling state and possible void defects of concrete in the reverse pouring process can be accurately predicted, so as to improve the prediction reliability of compactness.

[0025] Secondly, the application can realize temperature-hydration-rheology three-field strong coupling, improve the prediction accuracy of lateral pressure and temperature control. Through the Arrhenius-type hydration kinetics equation and the heat conduction equation, and the temperature rise caused by hydration exothermic is fed back to the rheological parameter correction, a closed-loop coupling mechanism of "temperature rise-hydration acceleration-yield stress and viscosity change" is formed. This mechanism can reflect the temperature gradient formed by hydration exothermic in early age concrete and its influence on rheological property and formwork lateral pressure, so that the peak value of lateral pressure at different pouring stages and its distribution along the height can be more accurately predicted, which helps to optimize the formwork design and construction safety control.

[0026] Furthermore, the application can accurately capture the free interface evolution and confluence characteristics of double-arch-foot reverse pouring. The VOF equation combined with PLIC geometric reconstruction can track the concrete-air interface of synchronous injection of left and right arch feet in real time under the explicit solution framework, and obtain the free surface shape, rising speed and spatiotemporal information of the confluence of the two streams in the middle. This technical feature solves the problem that the traditional method cannot predict the intersection position and time of the reverse pouring interface, and can provide quantitative basis for preventing cold joints and ensuring continuous pouring quality.

[0027] In addition, the application can also improve the numerical stability and calculation efficiency. By dynamically controlling the time step based on the CFL stability criterion, and discretizing the calculation domain on the MAC staggered grid, the numerical divergence caused by improper step length can be avoided while ensuring the accuracy of the explicit projection method for solving the flow field, so that the simulation of long-time pouring process is stable and efficient, which is suitable for rapid iteration and scheme comparison in engineering application.

[0028] Through the organic combination of the above technical features, the application forms a multi-physical field coupling simulation method that can simultaneously consider rheology, temperature, hydration and interface behavior, effectively improves the accuracy of density, lateral pressure, temperature control and interface prediction in the reverse pouring process of tunnel arch wall, and provides a reliable analysis tool for construction parameter optimization and safety control. BRIEF DESCRIPTION OF DRAWINGS

[0029] In order to more clearly illustrate the technical solutions in the embodiments of the application or the prior art, the following will briefly introduce the drawings needed to be used in the embodiments or prior art description. Obviously, the drawings in the following description only some embodiments of the application, and for those skilled in the art, other drawings can also be obtained without creative labor based on the structures shown in these drawings.

[0030] Figure 1 A flowchart of a multi-physical field coupling simulation method for reverse pouring of tunnel arch wall is disclosed in the application; Figure 2 A MAC grid and variable arrangement schematic diagram is disclosed in the application; Figure 3 A VOF-PLIC free surface tracking schematic diagram disclosed by the present application; Figure 4 A temperature-hydration-rheology coupling relationship schematic diagram disclosed by the present application; Figure 5 A pouring geometry model and boundary condition schematic diagram disclosed by the present application; Figure 6 A Guangdong Meizhou Changtian tunnel modeling area schematic diagram disclosed by the present application, wherein Figure 6 (a) is a tunnel cross section schematic diagram, Figure 6 (b) is a unilateral arch wall analysis area schematic diagram; Figure 7 A curve schematic diagram of free surface height changing with time in the process of reverse pouring simulation of the Guangdong Meizhou Changtian tunnel disclosed by the present application; Figure 8 A velocity field schematic diagram of the Guangdong Meizhou Changtian tunnel reverse pouring simulation at t = 1.0h disclosed by the present application; Figure 9 A temperature field schematic diagram of the Guangdong Meizhou Changtian tunnel reverse pouring simulation at t = 1.0h disclosed by the present application; Figure 10 A hydration degree distribution schematic diagram of the Guangdong Meizhou Changtian tunnel reverse pouring simulation at t = 1.0h disclosed by the present application; Figure 11 A thixotropic structure parameter distribution schematic diagram of the Guangdong Meizhou Changtian tunnel reverse pouring simulation at t = 1.0h disclosed by the present application; Figure 12 A formwork side pressure distribution schematic diagram of the Guangdong Meizhou Changtian tunnel reverse pouring simulation at t = 1.0h disclosed by the present application.

[0031] The purposes, functional features and advantages of the present application will be further described with reference to the accompanying drawings in conjunction with the embodiments. DETAILED DESCRIPTION

[0032] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the accompanying drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor fall within the protection scope of the present application.

[0033] It should be noted that the technical solutions of the various embodiments of the present invention can be combined with each other, but only if they are based on the ability of those skilled in the art to implement them. When the combination of technical solutions is contradictory or cannot be implemented, it should be considered that such combination of technical solutions does not exist and is not within the scope of protection claimed by the present invention.

[0034] Example 1

[0035] See Figures 1-5 The invention provides a multiphysics coupling simulation method for reverse casting of tunnel arch walls, the steps of which include: S1. Construct and discretize the computational domain: First, the geometric parameters of the tunnel arch wall are obtained. Based on these parameters, a computational domain for the arch wall elevation (height H, width L along the wall, and bottom edge as the arch foot line) is established in a two-dimensional plane, and corresponding boundary conditions are set according to the actual construction boundary. To improve numerical stability and accuracy, a MAC (Marker-and-Cell) staggered mesh is preferred to discretize the computational domain, allowing the velocity and pressure components to be staggered at the mesh nodes, facilitating mass conservation control during subsequent explicit projection method solutions. This method reduces numerical oscillations and ensures stability and computational efficiency in long-term casting simulations.

[0036] In this embodiment, the geometric parameters of the tunnel arch wall are obtained, including the arch foot position and spacing, arch wall height H, analysis width L along the wall direction, lining thickness b, center coordinates of the pouring opening and its size and quantity, design pouring height, segmented pouring sequence and single pouring volume.

[0037] In the MAC grid structure, pressure p, temperature T, and degree of hydration... Thixotropic structural parameters and volume fraction Arranged at the center of the grid cell, horizontal velocity The components are arranged at the center of the vertical plane of the unit, and the vertical velocity is... The components are arranged at the center of the unit's horizontal plane. Wherein: The unit is Pa; Absolute temperature, unit K; The degree of hydration is dimensionless, with a value of 0 to... ; is a dimensionless thixotropic structure parameter, with a value of 0 to 1; This represents the volume fraction of concrete within the element, ranging from 0 (air) to 1 (completely filled with concrete). The staggered arrangement of MACs effectively suppresses numerical oscillations caused by pressure-velocity coupling, improving the stability of the solutions for the free surface and pressure field. The boundary conditions are set as follows: at the arch foot pouring inlet, a given vertical inlet velocity... and mold temperature VOF value ; no-slip condition at the template and surrounding rock boundary , temperature is treated as convective heat transfer condition; upper opening, atmospheric pressure boundary is adopted.

[0038] S2, define the rheological constitutive and thixotropic evolution model of concrete: In view of the characteristics of concrete experiencing high shear (arch foot) and low shear (upper segment of arch wall) alternately in the reverse pouring process, the Herschel-Bulkley model with thixotropic recovery mechanism is used to describe the shear stress behavior, and the yield stress and apparent viscosity of which evolve with the thixotropic structure parameter, hydration degree and temperature. Further, the "recovery-damage" process of thixotropic structure is described by a first-order differential equation, so as to reflect the reconstruction and disintegration of the structure network of fresh concrete in time and space. This model can more truly predict the change of concrete fluidity at different pouring stages, help to improve the prediction accuracy of compactness, and reduce the risk of voids or cold joints caused by the deviation of rheological property estimation.

[0039] S3, establish a coupled calculation model of temperature field and hydration kinetics: According to the hydration reaction law of concrete, the evolution of hydration degree with time is calculated by using the Arrhenius type exponential model, and the temperature rise effect caused by hydration heat release is taken into account by combining with the heat conduction equation, so as to establish the coupling relationship between temperature field, hydration degree and rheological parameters, see Figure 4 . Specifically, the hydration heat release is taken as a heat source into the heat conduction equation, and after the temperature distribution is obtained, the yield stress and viscosity parameters of the Herschel-Bulkley model are feedback corrected, so as to form a closed-loop coupling mechanism of "temperature rise→hydration acceleration→rheological property change". This way can reflect the temperature gradient of early age concrete and its influence on hardening process and template side pressure, so as to improve the accuracy of temperature control and side pressure prediction.

[0040] S4, solve the flow field by using explicit projection method: After defining the coupling relationship between the constitutive and field, the Navier-Stokes equation is solved by using explicit projection method for time advancement, and the velocity field and pressure field distribution at each time are obtained. The explicit format is easy to synchronize with VOF-PLIC interface tracking and thixotropic evolution equation, and can well maintain mass conservation with MAC grid, which is suitable for numerical simulation of concrete reverse pouring with large deformation and free surface flow.

[0041] S5, solve the free interface by using VOF-PLIC interface tracking algorithm: For the working condition of synchronous reverse pouring of left and right arch feet, the volume fraction method (VOF) is combined with the piecewise linear interface calculation (PLIC) geometric reconstruction method to track the space-time of the free interface between concrete and air. This method obtains the volume fraction distribution by solving the VOF transport equation, and then uses PLIC to geometrically reconstruct the interface, so as to obtain the position, shape, rising speed of the free surface, and the time and position of the convergence of the two streams in the middle, see Figure 3 Such a setting can accurately predict the advancing process of the pouring interface, provide a basis for judging whether there is a cold joint and evaluating the filling uniformity, and further improve the engineering reliability of the density prediction.

[0042] S6, calculate the template side pressure and output the result: Based on the flow field and corresponding stress distribution obtained in step S4, the side pressure borne by the template is calculated according to the superposition principle of fluid statics and dynamics, and the curve of the side pressure changing with the pouring height is drawn, and the maximum value information is extracted. Since the side pressure is affected by the coupling of rheological properties, temperature field and interface shape, this method can accurately reflect the maximum side pressure and its distribution trend under different working conditions, thereby providing quantitative reference for template strength design and construction safety control.

[0043] S7, dynamically control the time step and loop calculation: To ensure the stability of numerical calculation, the time step is dynamically determined according to the CFL (Courant-Friedrichs-Lewy) stability criterion, so that the explicit solution process can meet the accuracy requirements while avoiding divergence. Then steps S2 to S7 are looped to continuously update the hydration degree, temperature, rheological parameters, flow field and interface position until the simulated pouring height reaches the set value or the calculation time is reached. This way can maintain stable and efficient operation in long-time pouring simulation, facilitating rapid iteration and comparative analysis of engineering schemes.

[0044] Through the organic combination of the above steps, the embodiment realizes the multi-physical field coupling simulation for tunnel arch wall reverse pouring: the Herschel-Bulkley rheological model with thixotropy recovery, VOF-PLIC interface tracking, temperature-hydration dynamics coupling and template side pressure calculation are unified in the explicit time marching framework, forming a complete CAE calculation process. This process can simultaneously reflect the interaction of flow field, temperature field, stress field and free interface, effectively overcoming the problem of lack of multi-physical field coupling in the prior art, significantly improving the prediction accuracy of density, side pressure, temperature control and interface position, and can guide the template design and construction parameter optimization.

[0045] As a preferred embodiment, the present embodiment defines and calculates the shear constitutive and thixotropic evolution process of fresh concrete in the following manner when implemented, so as to solve the problem of inaccurate prediction of fluidity and compactness due to differences in shear environment during reverse pouring.

[0046] S21, Herschel-Bulkley model with thixotropic recovery mechanism: In step S2, the Herschel-Bulkley (HB) model with a thixotropic recovery mechanism is selected to describe the shear stress behavior of fresh concrete, and its expression is as follows: ; wherein: is the shear stress, unit Pa; is the shear rate, calculated from the velocity gradient tensor, unit s -1 ; is the yield stress related to the thixotropic structure parameter , unit Pa; is the consistency coefficient (viscosity coefficient) related to the thixotropic structure parameter , unit Pa· ; is the flow index is a shear thinning fluid, is a Bingham type, is a shear thickening fluid; is the thixotropic structure parameter, indicates that the material structure is completely restored, and the viscosity is high, indicates that the structure is completely destroyed, and the fluidity is the highest.

[0047] S22, parameterization expression of yield stress and consistency coefficient: In order to introduce the structure parameter λ into the constitutive relation, the present embodiment further expresses the yield stress and the consistency coefficient as follows: ; ; wherein are all material constants calibrated by indoor rheological test. This expression makes the yield stress and viscosity change continuously with the change of thixotropic structure state, and can reflect the response of concrete under different shear histories in numerical calculation.

[0048] S23, evolution equation of thixotropic structure parameter:​ To depict the evolution of the structure parameter λ in time and space, the embodiment introduces the following first-order differential equation in step S2: ; Wherein: is the structure recovery coefficient, unit , representing the speed of the structure recovering to under the action of no shear; is the shear damage coefficient, unit s, representing the sensitivity of shear rate to structure damage; From the above: when the shear rate is large (such as the arch foot area), the second term dominates, λ decreases rapidly, the structure damage is intensified, and the concrete shows low yield stress and low viscosity, i.e. the flowability is enhanced; when the shear rate is small (such as the upper part of the arch wall), the first term A(1-λ) acts significantly, λ gradually increases, the structure tends to recover, the viscosity rises, and the flowability decreases.

[0049] Through the introduction of the above specific expression and evolution equation, the rheological model of the embodiment can accurately reflect the spatial shear difference of the reverse pouring condition, so as to realize the fine prediction of the change of the flowability of the concrete with time and position in the numerical simulation. Combined with the subsequent VOF-PLIC interface tracking and side pressure calculation, the judgment accuracy of the dense degree during pouring can be effectively improved, and the risk of cavities or cold joints caused by flowability estimation deviation can be reduced. In addition, after coupling the model with temperature-hydration dynamics, the temperature rise caused by hydration heat can be fed back to and correction, further improving the reliability of side pressure and temperature control analysis.

[0050] As a preferred embodiment, the coupling model is established in the following manner to overcome the problem that the mutual influence of the temperature gradient, hydration process and rheological properties of the early-age concrete during reverse pouring is difficult to accurately predict, thereby improving the reliability of temperature control, side pressure and dense degree analysis.

[0051] S31, hydration kinetics equation: the evolution of hydration degree with time is described by using an Arrhenius-type exponential model, which has the form: ; Wherein: is the current hydration degree, dimensionless; is the reference hydration rate constant before temperature correction, unit ; is the apparent activation energy of hydration reaction, unit J / mol; R is the gas constant, with the unit of J / (mol·K); T is the absolute temperature, with the unit of K; a is the final hydration value, with the range of 0.8-0.9; The above equation reflects the significant influence of temperature on the hydration rate, and the increase of temperature accelerates the hydration process, thus forming a direct coupling between the temperature field and the hydration degree in numerical value.

[0052] S32, early-age temperature field control equation: the temperature field is determined by heat conduction and hydration heat release, and satisfies the following coupling equation: ; wherein: ρ is the density of concrete, with the unit of kg / m3; ; Cp is the specific heat capacity of concrete, with the unit of J / (kg·K); k is the effective thermal conductivity, with the unit of W / (m·K); Q is the hydration heat release term, with the unit of J / (m3·s); which can be written as wherein λ is the hydration heat release potential per unit mass, with the unit of J / kg; The above equation converts the hydration kinetics calculation result of into a heat source term introduced into the heat conduction process, so as to reflect the temperature gradient formed inside the pouring body due to hydration heat release and its development over time.

[0053] S33, coupling relationship between temperature, hydration and rheological parameters: in order to form a strong coupling link of “temperature-hydration-rheology”, the present embodiment further establishes the following correction relationship: Influence of temperature on rheological parameters: the increase of temperature will reduce the viscosity of concrete, and the following correction is adopted: ; wherein α is the temperature sensitivity coefficient, Tref is the reference temperature; this correction enables the temperature rise caused by hydration heat release to be fed back to the consistency coefficient of the Herschel-Bulkley model, thereby affecting the flow resistance and lateral pressure distribution.

[0054] Influence of hydration degree on rheological parameters: with the progress of hydration, the structure of concrete gradually hardens, and the yield stress increases, which is expressed as: ; wherein β is the hydration strengthening coefficient; this relationship introduces the hydration degree into the yield stress calculation, so that the hardening process is reflected in numerical simulation.

[0055] Through the above modification, a strong coupling link of "temperature-hydration-rheology" is formed, and comprehensive reflection on the performance evolution in early age is realized.

[0056] Through the modeling of S31-S33, the embodiment realizes strong coupling of temperature field, hydration kinetics and rheological properties: hydration heat is included in the heat conduction calculation, which can accurately predict the temperature distribution and gradient change in the early age concrete; temperature change directly affects viscosity through exponential correction, and the increase of hydration degree increases the yield stress, which together change the flow and filling performance of concrete; the closed-loop coupling mechanism formed can reflect the chain reaction of "temperature rise → hydration acceleration → yield stress and viscosity evolution", so as to more realistically reproduce the side pressure change and interface advancing characteristics caused by hydration heat and structure hardening in the simulation of reverse casting. Combined with the thixotropic HB model in step S2 and the side pressure calculation in step S6, this coupling method significantly improves the accuracy of temperature control and side pressure prediction, and helps to predict the influence of local flow difference caused by temperature gradient on the dense degree, providing a feasible numerical way to solve the problems of missing multi-physical field coupling simulation and inaccurate key parameter prediction in reverse casting.

[0057] As a preferred embodiment, to further clarify the implementation details of the explicit projection method in step S4, the embodiment completes the flow field solution in the following order according to the numerical characteristics of the reverse casting condition, so as to ensure that the calculation stability and mass conservation are considered on the premise of not introducing excessive complexity, thereby improving the prediction accuracy of the formwork side pressure, interface advancing and dense degree.

[0058] S41, predict the velocity field: On the basis of the known velocity , apparent viscosity and shear stress at the last time step, the predicted velocity without pressure term is calculated first, and the discrete momentum equation is: ; Wherein: is the time step, unit s; is the spatial convection term, and in the embodiment, first-order upwind difference or second-order MUSCL format is preferably used for discretization to suppress numerical oscillation caused by convection dominance; is the viscous stress tensor given by the Herschel-Bulkley model, unit Pa; is the gravity acceleration vector, unit ; This step realizes the explicit promotion of inertia, convection, viscous shear and gravity, providing a basic field for subsequent pressure correction.

[0059] S42, solve the pressure Poisson equation: To meet the incompressible condition , a pressure Poisson equation needs to be constructed and solved, as shown below: ; Where is the Laplace operator, which is preferably discretized using a five-point difference format on the MAC staggered grid, which can improve the discrete precision and facilitate boundary condition processing; Pressure unknown is solved by Jacobi or Gauss-Seidel iteration until the residual meets the set convergence criteria.

[0060] The role of this step is to adjust the predicted velocity through the pressure field to meet the continuity equation, thereby avoiding numerical divergence and accurately reflecting the change in lateral pressure due to flow restriction during the reverse pouring process.

[0061] S42, velocity correction: After obtaining , the predicted velocity is projected and corrected according to the following formula: ; Where, Central difference is used for calculation to ensure that the velocity field after correction has second-order accuracy in space and strictly satisfies mass conservation.

[0062] Through the above three-step explicit projection method, the present embodiment realizes stable solution of the Navier-Stokes equation on the MAC staggered grid. The prediction step clearly shows the contribution of inertia, viscosity and gravity, which can accurately capture the difference between high flow velocity at the arch foot and slow flow at the upper section of the arch wall during reverse pouring; the solution of the pressure Poisson equation ensures the mass conservation of the velocity field, making the template lateral pressure calculated from the flow field (step S6) more reliable; the combination of explicit format and CFL dynamic time step control (step S7) can maintain numerical stability under conditions of large deformation and free surface flow, avoiding non-physical oscillation.

[0063] The above embodiment can effectively improve the accuracy and stability of flow field calculation, and cooperates with the thixotropic HB model of step S2 and the temperature-hydration-rheological coupling model of step S3, so that the overall simulation method can more accurately predict the compactness, lateral pressure distribution and interface advancement during the reverse pouring process, thereby directly providing a numerical solution to the technical problems of "lack of multi-physical field coupling simulation and inaccurate prediction of related parameters".

[0064] As a preferred embodiment, to further clarify the calculation method of the spatio-temporal evolution of the concrete-air free interface in step S5, this embodiment targets the working condition of reverse pouring with synchronous feeding at the left and right springings, and implements interface tracking in the following order to improve the prediction accuracy of interface advancement, material flow confluence, and filling state, thereby improving the judgment of compactness and cold joint risk.

[0065] S51. Interface transport equation: The volume fraction method is used to describe the distribution of the concrete phase in the grid cells, and its interface transport equation is: ; where is the volume fraction of the concrete phase. F = 1 indicates that the cell is filled with concrete, F = 0 indicates that it is all air, and 0 < F < 1 indicates that there is an interface in the cell; u is the velocity field at the current time step, obtained by the explicit projection method in step S4. The above equation tracks the overall position of the concrete-air interface through convective transport and can adapt to the severe deformation of the free surface during reverse pouring.

[0066] S52. Volume fraction update and PLIC interface reconstruction: Within each time step, first calculate the volume fraction fluxes on each grid surface and update the distribution accordingly. For the mixed cells with 0 < F < 1, the Piecewise Linear Interface Calculation (PLIC) method is used for geometric reconstruction: Based on the local volume fraction gradient fit the plane equation of the interface within the cell to obtain the precise geometric position and shape of the interface.

[0067] Compared with the approximation that only depends on the volume fraction, PLIC can retain the interface curvature and normal information with higher accuracy, reduce the interface diffusion error; it can clearly present the confluence process of the left and right springing jets in the middle, providing an intuitive basis for judging whether there is a cold joint.

[0068] S53. Interface and filling state parameter extraction; Based on the updated distribution, the following can be quantitatively calculated: Average height and local rising velocity of the free liquid surface: Statistically obtain the top position of the cells occupied by concrete to get the average elevation of the liquid surface along the height direction of the arch wall and its change rate, reflecting the overall trend of the pouring front advancement; Time and interface shape of the confluence of the left and right springing jets in the middle: Identify the moment when the leading edges of the volume fractions of the two side material flows first meet in the middle of the wall and the corresponding interface shape, which can be used to evaluate the continuity of the double-stream lap; Volume ratio of concrete in the wall: Integrate The value is the proportion of the total volume of the poured concrete, which helps to determine the pouring progress and the degree of dense filling.

[0069] In this embodiment, high-precision tracking of the concrete-air free interface is realized in the explicit time-marching framework. Specifically, the VOF equation ensures the conservation of the overall convective transport of the interface, and the PLIC reconstruction preserves the geometric details of the interface. The combination of the two can accurately predict the free surface morphology, the rising speed, and the convergence characteristics of the left and right arch foot flow in the middle. The obtained interface and filling parameters can be directly used to evaluate the dense degree in the reverse pouring process, and the cold joints or void defects that may be formed due to the delayed flow lapping can be found in advance. In conjunction with the thixotropic HB model in step S2, the temperature-hydration-rheological coupling model in step S3, and the calculation of the side pressure in step S6, the simulation method can simultaneously reflect the mutual influence of rheology, thermodynamics, and interface behavior in the same calculation process, thereby effectively solving the problems of missing multi-physical field coupling simulation and inaccurate prediction of key parameters in the prior art, and providing reliable numerical basis for form design, construction parameter optimization, and quality control.

[0070] As a preferred embodiment, for the calculation of the form side pressure in step S6, this embodiment further provides a specific mechanical model and implementation method to accurately obtain the distribution of the side pressure along the arch wall height in the reverse pouring simulation, thereby overcoming the problem of inaccurate side pressure prediction in the prior art and providing a quantitative basis for form design and construction safety control.

[0071] S61, form side pressure calculation: On the basis of the completed flow field and stress field solution, the equivalent side pressure acting on the unit form area is calculated according to the static pressure and the viscous stress tensor in the element by the following formula: wherein, is the static pressure obtained in the element, with the unit of Pa; is the viscous stress tensor, which is calculated by the Herschel-Bulkley constitutive model in step S2; is the outward normal unit vector of the form; is the component of the stress tensor in the local coordinate system of the form; is the equivalent side pressure acting on the form; The form side pressure is composed of the combined effect of the fluid static pressure and the viscous shear stress in the normal direction, which can reflect the comprehensive effect of the dynamic pressure and viscous resistance generated by the concrete flow on the side form during the reverse pouring process. ​​

[0072] S62, lateral pressure distribution extraction and output: extracting in the direction of the height of the arch wall (y direction) per unit distribution, forming a curve of the change of lateral pressure with height, and determining the maximum lateral pressure value and its corresponding height position from it. This data can be directly used for: template structure strength and rigidity design, to ensure the safety of the template under the maximum lateral pressure; support system arrangement and reinforcement scheme development, to reasonably distribute the load; safety assessment and monitoring of the construction process, to prevent the risk of structural instability or form explosion caused by lateral pressure overrun.

[0073] Since the calculation combines static pressure and viscous shear stress, and considers the directionality of the outward normal of the template, it can more truly reflect the difference in the effect of the flow state of concrete at different heights on the side mold during reverse pouring.

[0074] As a preferred embodiment, for the numerical stability control of the explicit time marching process in step S7, this embodiment uses the CFL stability criterion to dynamically determine the time step Δt to avoid numerical divergence caused by too large a step size, ensuring that the calculation is stable and the accuracy is controllable in long-time reverse pouring simulation, thereby supporting the reliability of the results of multi-physical field coupling.

[0075] In the calculation steps that need to be advanced synchronously, such as explicit projection method and thixotropic evolution solution, the time step Δt must satisfy the following CFL condition: : where: is the Courant number safety factor, usually taken as 0.2-0.5, to leave enough numerical stability margin; is the size of the calculation domain grid in the x and y directions; is the maximum absolute value of the velocity component at the current time step, reflecting the fastest scale of convection transport; is the maximum value of the apparent viscosity of the whole domain, derived from the calculation results of the Herschel-Bulkley model, is the density of the concrete.

[0076] The three terms of the above formula correspond to the limiting conditions of convection stability, cross-grid propagation time and viscous diffusion stability respectively, and the minimum value is multiplied by the safety factor to obtain the maximum Δt allowed for this time step.

[0077] At the beginning of each cycle, the , and , and the step length is used to promote the solution of steps S2-S6. In this way, a larger step length can be used to improve the calculation efficiency while ensuring the stability of the explicit format.

[0078] To further illustrate the inventive concept of the present application, see Figures 6-12 , which is described in conjunction with a specific application scenario.

[0079] I. Geometric modeling and boundary condition setting: Taking the Changtian tunnel project in Meizhou City, Guangdong Province as the application background, see Figure 6 , a two-dimensional arch wall facade calculation domain with a height of 5.0 m and a width of 6.0 m is constructed, the left wall is selected as the analysis area, the wall thickness is processed according to the two-dimensional unit thickness, the bottom edge = arch foot position (as the pouring inlet), the left and right sides are template or initial support side wall, and the top edge is the arch top free surface.

[0080] See Figure 7 , 121x101 MAC staggered grids are used for discretization, in which the pressure p, temperature T, hydration degree a, thixotropic parameter l and volume fraction VOF(F) are arranged at the center of the element, and the velocity components u, v are arranged at the vertical edge and horizontal edge of the element, respectively, to avoid pressure oscillation and improve the stability of the coupled field.

[0081] The boundary conditions are set as shown in the following table: ; II. Material model and physical property parameter setting: 1. Basic physical property parameters of concrete: density p = 2350 kg / m³, specific heat c = 900 J / kg·K, thermal conductivity k = 1.7 W / m·K.

[0082] 2. Herschel-Bulkley thixotropic constitutive model: the shear stress is represented as , wherein τ y0 = 100 Pa, K0 = 18 Pa·sⁿ, n = 0.45; the evolution equation of thixotropic structure parameter l is , the recovery coefficient = 0.20 s -1 , and the damage coefficient = 1.8.

[0083] 3. Hydration kinetics model: the Arrhenius type equation is used , wherein = 3x10 4 , E = 40 kJ / mol, = 0.85.

[0084] 4. Temperature field control equation: , hydration heat source , hydration latent heat H h = 280 kJ / kg.

[0085] III. Explicit time-marching calculation procedure: Step S1, domain initialization: Generate 2D geometry model and MAC mesh, initialize velocity field, temperature field, hydration degree and thixotropic parameters.

[0086] Step S2, shear rate calculation: Calculate shear rate based on velocity gradient tensor , provide input for thixotropic model update.

[0087] Step S3, thixotropic structure parameter λ n+1 update: Update λ distribution according to evolution equation . See Figure 11 , simulation results show that λ in high shear zone at arch foot drops to 0.1-0.3 (structure failure), and λ in low shear zone at upper part recovers to 0.95-1.0, reflecting the spatio-temporal differentiation characteristics of "shear failure-structure recovery".

[0088] Step S4, hydration degree α evolution calculation: Update hydration degree through Arrhenius model , used for temperature field and rheological parameter coupling calculation. See Figure 10 .

[0089] Step S5, rheological parameter correction: Update yield stress and consistency coefficient according to λ and α, where C1=1.2, C2=2.6.

[0090] Step S6, flow field solution (explicit projection method): Predict velocity: ; Pressure Poisson equation: , solved by Jacobi iteration; Velocity correction: .

[0091] Step S7, VOF-PLIC interface tracking: Solve volume fraction transport equation , and use PLIC geometric reconstruction interface to extract free liquid surface height h(t). See Figure 7 , simulation results show that the liquid surface rises nearly linearly in the early stage (0-1h), and the speed slows down in the later stage due to the increase of thixotropic recovery, and finally stabilizes at the arch top height of 5m.

[0092] Step S8, temperature field update: By heat conduction and hydration heat release coupling equation Calculate temperature distribution, see Figure 9 , simulation presents a "lower hot upper cold" stratification trend, the bottom temperature rises to 28-29℃, and the hydration degree reaches 0.25-0.30.

[0093] Step S9, lateral pressure extraction and output: According to the formula Calculate the template lateral pressure to generate a distribution curve along the height. See Figure 12 , in the simulation of Changtian tunnel, the maximum lateral pressure is 60-65kPa, and it rises approximately linearly in the range of 0-2m depth, and the thixotropic recovery leads to a significant decrease in the upper lateral pressure.

[0094] Four, time step control and loop termination: Based on the CFL stability criterion to dynamically control the time step, to ensure the numerical stability of the explicit solution process. Loop execution steps S2 to S9 until the pouring height reaches 5m or the simulation time is terminated.

[0095] Example 2

[0096] The application also provides a multi-physical field coupling simulation system for tunnel arch wall reverse pouring, which adopts the simulation method described in example 1 and is composed of several functional modules. Each module works cooperatively according to the predetermined logic to realize the whole process calculation from geometric modeling to pouring process simulation, so as to solve the problems of lack of multi-physical field coupling simulation, density, lateral pressure, temperature control and inaccurate interface prediction in reverse pouring. The composition and implementation of the system are described below.

[0097] Specifically, it includes: Modeling module, for establishing a two-dimensional arch wall elevation calculation domain according to the obtained tunnel arch wall geometric parameters, and applying boundary conditions consistent with actual construction in the calculation domain. This module preferably uses MAC staggered grid to discretize the calculation domain, to obtain the numerical advantage of staggered arrangement of velocity field and pressure field at grid nodes, which is beneficial to the mass conservation control in subsequent explicit projection method solution. This processing method can reduce numerical oscillation, and ensure the stability and calculation efficiency of long-time simulation.

[0098] Material model module, for realizing Herschel-Bulkley rheological model with thixotropic recovery mechanism and thixotropic evolution equation. The module is built-in yield stress And consistency coefficient The parameterized form of the Bingham plastic model is adopted, and the recovery and damage processes of the thixotropic structure parameter λ are described by first-order differential equations. Through this module, the simulation can reflect the working condition characteristics of arch foot high shear leading to structural damage and increased fluidity, and low shear in the upper section of the arch wall promoting structural recovery and decreased fluidity, thereby improving the prediction accuracy of concrete fluidity and density.

[0099] A thermo-hydraulic coupling module is used to calculate the interaction of hydration reaction and temperature field and its influence on rheological properties. The module solves the evolution of hydration degree α based on the Arrhenius-type exponential model, and incorporates the hydration heat release effect through the heat conduction equation to form a coupled calculation of temperature field and hydration degree. Further, temperature and α are fed back to the rheological parameter (viscosity and yield stress) correction formula to establish a closed-loop coupling mechanism of "temperature rise → accelerated hydration → change of rheological properties", thereby improving the accuracy of the prediction of lateral pressure and temperature control under the influence of early age temperature gradient and hardening process.

[0100] A flow field solving module is used to solve the velocity field and pressure field based on the explicit projection method. The module completes the predicted velocity calculation without pressure term, pressure field solving to meet the incompressible condition, and velocity field projection correction in sequence according to the steps of prediction-pressure Poisson-velocity correction, ensuring mass conservation and accurately capturing the difference between high-speed flow at the arch foot and slow flow in the upper section during reverse pouring. In cooperation with the material model module, the Herschel-Bulkley stress tensor can be introduced into the momentum equation to improve the consistency of the flow field with the actual pouring process.

[0101] An interface tracking module is used to simulate the space-time evolution of the concrete-air free interface using the VOF-PLIC method. The module solves the volume fraction transport equation at each time step, updates the concrete phase distribution, and obtains the interface position and morphology through PLIC geometric reconstruction for mixed cells, which can accurately output the average height of the free liquid surface, the local rising velocity, the time of left and right arch foot jets converging in the middle, and the interface morphology, etc., providing a basis for judging the filling uniformity and cold joint risk, thereby improving the reliability of the density prediction.

[0102] A lateral pressure evaluation module is used to calculate the stress distribution of the formwork and output the pressure curve and extreme value. The module calculates the equivalent lateral pressure based on the flow field and stress distribution according to the formula for the synthesis of static pressure and viscous shear stress in the outer normal direction of the formwork, extracts the distribution data along the height direction of the arch wall, draws the change curve and outputs the maximum lateral pressure value and the corresponding height, providing quantitative basis for formwork structure design, support system configuration and construction safety evaluation.

[0103] The time control module is used for dynamically adjusting the simulation step length according to the CFL stability criterion, and driving the cooperative operation of each module until the pouring process simulation is completed. The module calculates the allowed maximum time step length according to the current maximum velocity component and the apparent viscosity according to the CFL condition at each cycle step, so as to ensure the numerical stability of the explicit solution process, and at the same time, as large as possible step length is used to improve the calculation efficiency. The stability control of the module makes the long-time multi-physical field coupling simulation feasible, and avoids the influence of numerical divergence on the side pressure and interface prediction accuracy.

[0104] The system combines the functions of geometric modeling, material constitutive, thermal-hydraulic coupling, flow field solving, interface tracking, side pressure evaluation and time control through modular design, and forms a complete simulation platform that can run in the same explicit time marching framework. The parameter transmission and coupling calculation between the modules make the rheology, temperature, hydration, interface and stress field update synchronously, which makes up for the problem of the lack of multi-physical field coupling simulation in the prior art; the material model and the thermal-hydraulic coupling module are linked, which can truly reproduce the changes of fluidity and side pressure of early-age concrete caused by hydration heat release and structure hardening; the interface tracking module and the side pressure evaluation module are combined, which can provide quantitative data of filling state and formwork stress, and significantly improve the prediction accuracy of compactness, side pressure, temperature control and interface position; the time control module guarantees the stability and efficiency of long-time simulation, so that the system can be practically applied to engineering scheme comparison and construction safety evaluation. The system corresponds to the simulation method of embodiment 1 in structure, and realizes the multi-physical field coupling simulation for tunnel arch wall reverse pouring through the cooperative work of each functional module, which can effectively solve the technical problem of inaccurate prediction of key parameters in reverse pouring, and provides a reliable analysis tool for formwork design and construction control.

[0105] The above is only the preferred embodiment of the present application, and does not limit the patent scope of the present application, and any equivalent structural transformation made under the inventive concept of the present application, or direct / indirect application in other related technical fields is included in the patent protection scope of the present application.

Claims

1. A multi-physics coupling simulation method for reverse pouring of a tunnel arch wall, characterized in that, The steps include: S1, acquiring geometric parameters of the tunnel arch wall, constructing a two-dimensional arch wall facade calculation domain based on the geometric parameters, discretizing the calculation domain and setting boundary conditions; S2, using a Herschel-Bulkley model containing a thixotropic recovery mechanism to describe shear stress behavior, and using a first-order differential equation to represent the recovery and destruction process of the thixotropic structure parameter; S3, using an Arrhenius-type exponential model to describe the hydration degree evolution, and considering the hydration heat release effect by combining a heat conduction equation to establish the coupling relationship of temperature, hydration degree and rheological parameters; S4, solving the flow field by using an explicit projection method; S5, solving the space-time evolution process of the concrete-air free interface by using a PLIC interface tracking algorithm based on the VOF method; S6, calculating the template side pressure according to the flow field and stress distribution, and outputting the side pressure curve and maximum value information with height; S7, dynamically controlling the time step based on the CFL stability criterion, and cyclically executing steps S2 to S7 until the pouring simulation reaches the set height or the termination time.

2. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1, characterized in that, In step S1, the geometric parameters of the tunnel arch wall are acquired, including the arch foot position and spacing, the arch wall height H, the analysis width L along the wall direction, the lining thickness b, the pouring port center coordinates and their size and number, the designed pouring height, the segmented pouring sequence and the single pouring amount.

3. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1 or 2, characterized in that, In step S2, the shear stress model is as follows: ； wherein: is the shear stress, is the shear rate, is the yield stress related to the thixotropic structure parameter, is the consistency coefficient related to the thixotropic structure parameter, is the flow index, is the thixotropic structure parameter.

4. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 3, characterized in that, In step S2, the thixotropic structure parameters evolve according to the following equations: ； wherein, is a structural recovery coefficient, is a shear failure coefficient.

5. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1 or 2, characterized in that, In step S3, the hydration kinetics equation uses an Arrhenius-type exponential model, as shown below: ； wherein, is the current degree of hydration, is the reference rate constant for hydration before temperature correction, is the apparent activation energy for hydration reaction, is the gas constant, is the absolute temperature, is the final value of hydration. The temperature field satisfies the heat conduction and hydration heat release coupling control equation, which is as follows: ； wherein, is the density of the concrete, is the specific heat capacity of the concrete, is the effective thermal conductivity, is the hydration heat release term.

6. The multi-physics coupling simulation method for reverse pouring of a tunnel vault wall according to claim 3, characterized in that, In step S3, the temperature correction relationship of the rheological parameter is introduced: The concrete consistency coefficient is corrected by the following formula: ； The yield stress is corrected by the following formula: ； wherein, is the temperature sensitivity coefficient, is the reference temperature; is the hydration enhancement coefficient, is the absolute temperature.

7. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1 or 2, characterized in that, In step S4, the step method for solving the flow field by using an explicit projection method is as follows: S41, predicting velocity field: under the premise of knowing the velocity of the last time step , apparent viscosity and shear stress , calculate the predicted velocity without pressure term : ； where: is the time step, is the spatial advection term, is the viscous stress tensor, is the gravitational acceleration vector; S42, construct and solve the pressure Poisson equation: To satisfy the incompressibility condition ; construct the pressure Poisson equation as follows: ； wherein is the Laplace operator; Pressure unknown Solved by Jacobi or Gauss-Seidel iteration; S43, velocity correction: after obtaining the projected correction of the predicted velocity is made as follows: ； wherein, The central difference is used to calculate the corrected velocity field, which satisfies the mass conservation condition.

8. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1 or 2, characterized in that, In step S5, the step method for solving the space-time evolution process of the concrete-air free interface is as follows: S51, the interface transport equation is as follows: ； wherein is the volume fraction of the concrete phase, is the velocity field of the current time step; S52、In each time step, the volume fraction flux is calculated along each grid face, and the is updated, and the unit of adopts the PLIC method to reconstruct the linear interface according to the local gradient , and the geometric position and shape of the free surface are obtained. S53、 by The distribution calculation obtains: the average height of the free surface and the local rising velocity; the time and interface shape of the left and right arch foot jets converging in the middle; the volume ratio of the concrete in the wall.

9. The multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to claim 1 or 2, characterized in that, In step S6, the template side pressure The calculation formula is as follows: ； wherein, P is the static pressure found within the cell, is the viscous stress tensor; is the outward normal vector to the template; is the stress tensor component; is the equivalent side pressure acting on the template; Extracting along the direction of the height of the arch wall Distributing, drawing the curve of the template side pressure with the height, and outputting the maximum side pressure value and the corresponding height.

10. A multi-physics coupling simulation system for reverse pouring of a tunnel arch wall, adopting the multi-physics coupling simulation method for reverse pouring of a tunnel arch wall according to any one of claims 1-9, characterized in that, It includes: A modeling module for establishing a two-dimensional calculation domain according to the tunnel geometric parameters and discretizing the calculation domain and applying boundary conditions; A material model module for implementing a Herschel-Bulkley rheological model containing thixotropic recovery and a thixotropic evolution equation; A thermal-hydration coupling module for calculating the interaction of hydration reaction and temperature field and its influence on rheological properties; A flow field solving module for explicit solving of velocity field and pressure field based on projection method; An interface tracking module for simulating the change of concrete-air interface by using VOF-PLIC method; A side pressure evaluation module for calculating the stress distribution of the template and outputting the pressure curve and extreme value; A time control module for dynamically adjusting the simulation step according to the CFL condition, and driving the collaborative operation of each module until the pouring process simulation is completed.

Images