A water-salt simulation method considering multiple chemical reaction processes

Abstract

The present invention relates to a method for simulating soil water and salt, including: collecting soil data, meteorological data, irrigation data, and groundwater data for driving the model; determining soil hydraulic parameters, evaporation stress parameters, solute transport parameters, ion adsorption parameters, precipitation / dissolution parameters, and simulation time; constructing a water and salt transport model considering multi-component salt competitive adsorption and precipitation / dissolution reactions; setting the initial conditions and upper and lower boundary conditions for model operation, inputting soil information, meteorological data, groundwater data, and required parameters, and solving the above-constructed model using numerical algorithms to obtain the water and salt distribution characteristics on the soil profile at a specified time. Based on the traditional water and salt transport model, the present invention considers the physical and chemical reaction processes of multi-component salts and the influence of salt crystallization on soil evaporation during the water and salt transport process, and can reflect the water movement and salt reaction and transport processes in the soil, making it more suitable for simulating the water and salt transport process in saline soil.

Description

Technical Field

[0001] The present invention relates to the technical field of agricultural production, and particularly to a water-salt simulation method considering multiple chemical reaction processes. Background Art

[0002] Studying the process of water-salt migration in soil is of great significance for understanding the distribution law of soil water-salt in saline-alkali areas and ensuring the sustainable development of agriculture. At present, a large number of experimental and simulation studies have been carried out in this regard. Among them, experimental studies can observe the law of soil water-salt migration under a certain specific treatment. However, there are many influencing factors in experimental studies, the experimental conditions are not easy to control, and long-term experimental observations are time-consuming and laborious. In the research of simulation models, the simulation input variables are easy to control, which has a good supplementary effect on long-term experimental studies. Soil saturated-unsaturated water-salt migration models based on kinetic processes, such as RZWQM, HYDRUS, and SWAP, etc., deal with the state transformation in the process of multi-component salt migration relatively simply, and do not consider the complex physical and chemical reaction processes enough, and it is difficult to describe the physical and chemical reaction processes such as ion adsorption, precipitation / dissolution, etc. in water-salt movement. And the models involving chemical reaction processes in the saturated zone are mainly hydrogeochemical models. For example, PHREEQC and MINTEQA2 are generally applicable to uniform flow and are not suitable for simulating water-salt reaction migration under complex unsaturated water flow movement. In addition, UNSATCHEM, HP1 (HYDRUS1D-PHREEQC), and LEACHM are currently commonly used hydrochemical models that can simulate salt reaction migration in unsaturated soil. These models consider the water stress on soil evaporation, but do not consider the influence of salt crystallization on soil evaporation enough. Existing water-salt migration models and hydrogeochemical models each have their own drawbacks. For example, currently, the water-salt migration models in unsaturated soil rarely consider the complex physical and chemical reaction processes of multi-component salts and their influence on soil evaporation, and it is difficult to describe the chemical processes such as precipitation, dissolution, ion exchange, adsorption, crystallization, etc. involved in salt migration. While hydrogeochemical models can simulate the chemical reaction processes of multi-component salts, but are mostly applicable to the saturated zone and are not suitable for simulating water-salt reaction migration in the unsaturated zone.

[0003] Therefore, how to provide a water-salt simulation method considering multiple chemical reaction processes and more comprehensively and accurately consider the water-salt migration process in the saturated-unsaturated state is a technical problem that needs to be solved urgently by those skilled in the art. Summary of the Invention

[0004] The present invention aims to solve at least one of the technical problems in the related art to a certain extent, and proposes a water-salt simulation method considering multiple chemical reaction processes. On the basis of the traditional water-salt migration model, multiple chemical reaction processes are considered, such as ion competitive adsorption, precipitation / dissolution, etc.; and on the basis of the reaction and migration of multi-component salts in the soil, the influence of salt crystallization on soil evaporation is considered, which can reflect the water movement and salt reaction and migration processes in the soil, making it more suitable for the simulation of water-salt migration processes in saline soils.

[0005] In view of this, the present invention proposes a soil water-salt simulation method considering the chemical reaction process of multi-component salts, including

[0006] Obtaining the specific values of the parameters associated with the soil water-salt simulation process;

[0007] Constructing a water-salt migration model considering multi-component salt competitive adsorption, precipitation / dissolution reactions, setting the initial conditions and upper and lower boundary conditions for model operation, and inputting the specific values of all parameters into the water-salt migration model for solution to obtain the water-salt distribution characteristics on the soil profile at a specified time; wherein, the water-salt migration model is used to describe the functional relationship between the parameters and the soil water-salt content.

[0008] In the soil water-salt simulation method of the present invention, the specific values of the parameters associated with the soil water-salt simulation process are soil data, meteorological data, irrigation data, and groundwater data for driving the model; determining soil hydraulic parameters, evaporation stress parameters, solute transport parameters, ion adsorption parameters, precipitation / dissolution parameters, and simulation time.

[0009] In an embodiment of the present invention, the above steps further include discretizing the soil profile into a number of spatial grid nodes, setting the simulation time and step size, and solving and calculating the moisture and salt of the discrete nodes by the implicit finite difference method.

[0010] In an embodiment of the present invention, the water-salt migration model includes: simulating water movement by the one-dimensional vertical mixed Richards equation; simulating the convection and mechanical dispersion processes of salts in the soil by the one-dimensional vertical convection-dispersion equation.

[0011] The expression of the one-dimensional vertical mixed Richards equation is

[0012]

[0013] In formula (1), t is time (d); z is the spatial coordinate in the vertical direction, positive downward (cm); θ is the soil volumetric water content (cm 3 / cm 3 ); h is the soil matrix potential or pressure head (cm); K(h) is the soil hydraulic conductivity (cm / d); S is the source-sink term.

[0014] The van Genuchten-Mualem model (abbreviated as VG-M) is adopted to describe the soil hydraulic properties. The relationships between the soil water content θ and the hydraulic conductivity K(h) with respect to the matric potential h are expressed as follows:

[0015]

[0016]

[0017] where S e is the effective saturation of the soil (-); θ r is the residual water content of the soil (cm 3 / cm 3 ); θ s is the saturated water content of the soil (cm 3 / cm 3 ); α (cm -1 ), m (-), and n (-) are soil shape coefficients, and m = 1 - 1 / n; K s is the saturated hydraulic conductivity of the soil (cm / d); l is the soil pore connectivity parameter (-), generally taking a value of 0.5.

[0018] The expression form of the one-dimensional vertical convection-dispersion equation is:

[0019]

[0020] where c k is the solution concentration of salt component k (mg / cm 3 ); S k is the adsorption concentration of salt component k in the soil (mg / g); is the precipitation concentration of salt component k in the soil (mg / g); S sk is the absorption amount of salt component k by the roots (mg / (cm 3 ·d)); ρ is the dry bulk density of the soil (g / cm 3 ); D sh (v, θ) is the hydrodynamic dispersion coefficient of the solute (cm 2 / d); v is the average flow velocity of the soil solution in the soil pores (cm / d); N is the total number of salt components; q is the unit-width flow rate of water in the vertical direction (cm / d), calculated according to Darcy's law for unsaturated water movement:

[0021]

[0022] The hydrodynamic dispersion coefficient D sh (v, θ) reflects the molecular diffusion and mechanical dispersion processes of the solute in the soil solution, and its expression is:

[0023]

[0024] In the formula, D L is the longitudinal mechanical dispersion coefficient (cm); D0 is the diffusion coefficient of the solute in free water (cm 2 / d).

[0025] The amount of solute absorbed by the root system is related to the root water uptake rate S and the solution concentration c of the solute k as follows:

[0026] S sk = S r ·S·c k (7)

[0027] In the formula, S r is the empirical parameter for the root system to absorb solute (-).

[0028] When considering the competitive adsorption among salt ions (Na + , K + , Mg 2+ and Ca 2+ ), assuming that the cation exchange capacity CEC (meq / 100 g soil) in the soil remains unchanged, the cation exchange is a reversible process, and the charges carried by exchangeable sodium (NaX), potassium (KX), calcium (Ca 1 / 2X), and magnesium (Mg 1 / 2 X) are equal to CEC, that is:

[0029] CEC = NaX + KX + Ca 1 / 2 X + Mg 1 / 2 X (8)

[0030] Then, according to the modified Gapon equation:

[0031]

[0032]

[0033]

[0034]

[0035]

[0036]

[0037] The charge number (meq / 100 g soil) of the exchangeable cations adsorbed is obtained by solving:

[0038]

[0039]

[0040]

[0041] [KX] = CEC - [NaX] - [Ca 1 / 2 X] - [Mg 1 / 2 X] (18)

[0042] Wherein, K1, K2, K3, K4, K5 and K6 are ion exchange selectivity coefficients (-); (·) represents ion activity (mol / L).

[0043] According to the number of adsorbed charges, the exchangeable cation concentration can be further calculated, and the obtained exchangeable cation concentration is substituted into the adsorption term S in the convective-dispersion equation (4) k , thus obtaining a multi-component solute transport model considering ion competitive adsorption.

[0044] In the present invention, the precipitation-dissolution processes of two salts, CaSO4 and NaCl, are also considered, and they satisfy the solubility product principle in the soil solution. Taking CaSO4 as an example, the dissolution equilibrium reaction equation can be expressed as:

[0045]

[0046] Its reaction rate R (mol / (m 2 ·min)) satisfies the second-order kinetic process:

[0047]

[0048] Wherein, k d and k p respectively represent the dissolution constant and the precipitation constant (dm 6 / (mol·m 2 ·min)); K sp is the solubility product constant (-); γ is the ion activity coefficient (-); is the unit jump function:

[0049]

[0050] Therefore, it is obtained that:

[0051]

[0052] Wherein, S a is the precipitation surface area (m 2 / dm 3 ); ξ is the unit conversion coefficient (-).

[0053] The salt precipitation on the soil surface will inhibit the evaporation rate of the soil, and its reduction ratio coefficient is γ k The empirical formula is adopted:

[0054]

[0055] In the formula, is the empirical parameter (mg / g).

[0056] In this study, the calculation of the activity coefficient γ adopts the modified Debye-Hückel formula, and its calculation results in high-concentration solutions still have high reliability.

[0057]

[0058] In formula (24), Z k is the valence of the k-th salt ion; a s and b s are the activity coefficient parameters of a certain ion in the formula (-), and the values are shown in Table 1; I is the ionic strength (-), and its definition is:

[0059]

[0060] Table 1 Activity coefficient parameters of ions

[0061]

[0062] In the present invention, the initial conditions and boundary conditions of soil moisture and salt need to be set. Among them, the initial moisture conditions can be expressed by water content or matrix potential:

[0063] θ(z) = θ i (z) 0 ≤ z ≤ L, t = 0 (26)

[0064] h(z) = h i (z) 0 ≤ z ≤ L, t = 0 (27)

[0065] In the formula, θ i (z) and h i (z) are the water contents (cm 3 / cm 3 ) and matrix potentials (cm) at different depths of the soil profile respectively; L is the maximum depth of the soil profile (cm).

[0066] The boundary conditions of moisture can be expressed as constant head or flux boundary:

[0067] h(0, t) = h0(t) z = 0, t ≥ 0 (28)

[0068]

[0069] where \(h_0(t)\) and \(q_0(t)\) are the pressure head (cm) and flow rate (cm / d) at the upper boundary, respectively.

[0070] The lower boundary can be expressed as a constant head boundary (Equation 29), a flow rate boundary (Equation 30), or a free drainage boundary (Equation 31):

[0071] \(h(L, t)=h\) L (t) \(z = L, t\geq0\) (30)

[0072]

[0073]

[0074] where \(h\) L (t) and \(q\) L (t) are the pressure head (cm) and flow rate (cm / d) at the lower boundary, respectively.

[0075] The initial conditions for soil salt simulation are the distributions of the concentrations of each component salt in the soil profile at the initial simulation time:

[0076]

[0077] where \(c\) ki (z), \(S\) ki (z) and are the solution concentration (mg / cm 3 )), adsorption concentration (mg / g), and precipitation concentration (mg / g) of salt at different depths in the soil profile, respectively.

[0078] The upper boundary of salt simulation can be expressed as a concentration boundary (Equation 34) or a concentration flow rate boundary (Equation 35):

[0079] \(c\) k (0, t)=c k0 (t) \(z = 0, t\geq0\) (34)

[0080]

[0081] where \(c\) ki (z), \(S\) ki (z) and are the solution concentration (mg / cm 3 )), adsorption concentration (mg / g), and precipitation concentration (mg / g) of salt at different depths in the soil profile, respectively.

[0082] where \(c\) k0 (t) is the soil solution concentration (g / L) at the upper boundary or the solution concentration in the water flux at the upper boundary (g / L).

[0083] The lower boundary of the salt simulation can be expressed as:

[0084] c k (L, t) = c kL (t) z = L, t ≥ 0 (36)

[0085]

[0086]

[0087] In the formula, c kL (t) is the solute concentration of the soil solution at the lower boundary (g / L) or the solution concentration in the water flux at the lower boundary (g / L). The soil profile is discretized into several spatial grid nodes, and the simulation time and step size are set.

[0088] In the present invention, the simulation of soil water movement adopts the one-dimensional vertical mixed Richards equation because this equation can simulate the water movement process in saturated-unsaturated soils and can reduce the mass conservation error. The salt simulation uses the one-dimensional vertical convection-dispersion equation to describe the convection and mechanical dispersion processes of salts in the soil, and at the same time considers the competitive adsorption, precipitation / dissolution and other chemical reaction processes of multi-component salt ions (Na + , K + , Mg 2+ , Ca 2+ , Cl- and SO4 2 -) in the soil and their influence processes on evaporation.

[0089] In an embodiment of the present invention, the process of solving the water movement using the one-dimensional vertical mixed Richards equation at least includes:

[0090] (1) Input the solution variables hydraulic conductivity, specific water storage capacity and upper and lower boundary conditions. Among them, for the upper boundary evaporation, the inhibition effect of soil surface salt crystallization is considered, and the chasing method is used to program and calculate the discrete linear tridiagonal equations to obtain the matrix potential or soil profile water content at the end of the time period;

[0091] (2) Judge whether the matrix potential meets the iterative accuracy requirement. "Yes", use the obtained matrix potential as the initial value for the next moment, and repeat the above loop for the calculation of the next moment; "No", go to step (3);

[0092] (3) Use the matrix potential to solve the hydraulic conductivity and specific water storage capacity again, and repeat the above calculation method until the accuracy requirement is met.

[0093] In one embodiment of the present invention, step (4) is further included. When the number of iteration cycles in step (3) exceeds 20 and fails to meet the iteration accuracy requirement, the step size of the calculation period is reduced to one third, and the above calculation method is repeated again for calculation.

[0094] In one embodiment of the present invention, the hydraulic conductivity and the specific water density are calculated by directly obtaining the parameters associated with the soil water and salt simulation process.

[0095] The process of using the one-dimensional vertical hybrid Richards equation to solve water movement in the present invention is to calculate potential evaporation, hydraulic conductivity K and specific water density C based on input soil information, meteorological data, irrigation data, etc., and determine upper and lower boundary conditions. The variables and boundary conditions solved above are substituted into the discrete Richards equation, and the linear tridiagonal equation system is solved to obtain the matrix potential h1 at the end of the time period. If the iteration accuracy requirement is met, the h1 obtained at the end of this time period is used as the initial value of the next moment, and the above cycle is repeated to perform the calculation at the next moment. If the accuracy requirement is not met, the hydraulic conductivity K1 and specific water density C1 are solved again based on the calculated h1, and the above calculation method is repeated until the accuracy requirement is met. If the number of iterations exceeds 20 times and still does not meet the iteration accuracy requirement, the time step DT of the calculation period is reduced to one-third, and the above calculation method is repeated again to calculate.

[0096] In one embodiment of the present invention, using a one-dimensional vertical convection-diffusion equation to simulate and describe the convection and mechanical diffusion process of salt in soil at least includes:

[0097] (1) First, determine whether to consider the chemical reaction process of salt; if "yes", proceed to step (2); if "no", proceed to step (3);

[0098] (2) The adsorption concentration and precipitation concentration of each salt component during the chemical reaction and the soil profile moisture content are input, and the tridiagonal coefficient matrix is solved by the chase method through simultaneous equations to obtain the salt solution concentration value at the end of the period; the salt adsorption and precipitation amount of the soil profile are calculated based on the salt solution concentration value to obtain the updated salt solution concentration value after the reaction;

[0099] (3) Solve the corresponding variables of the parameters directly obtained from the soil water-salt simulation process and input them, discretize the equation, and solve the tridiagonal equation to obtain the salt solution concentration value at the end of the time period. The salt solution concentration value at the end of this time period is used as the initial salt solution concentration at the next moment, and the above cycle is repeated to perform calculations at the next moment.

[0100] In one embodiment of the present invention, the adsorption concentration and the precipitation concentration are respectively calculated by directly obtaining the associated parameters of the soil water-salt simulation process using the competitive adsorption equation and the precipitation-dissolution equation.

[0101] In one embodiment of the present invention, step (2) further includes that when the updated salt solution concentration value meets the iterative accuracy requirement, the calculation for the next moment is carried out; otherwise, the adsorption concentration and precipitation concentration are recalculated according to the updated salt solution concentration value, and the salt solution concentration value is updated again until the iterative accuracy requirement is met.

[0102] In the present invention, the one-dimensional vertical convection-dispersion equation is used to simulate and describe the convection and mechanical dispersion processes of salts in the soil: in this simulation method, the chemical reaction process can be optionally considered or not considered. First, it is judged whether to consider the chemical reaction process of salts. If chemical reactions such as adsorption and precipitation are considered, the adsorption concentration and precipitation concentration of each salt component in the soil during the chemical reaction in this period are first calculated according to the competitive adsorption equation and the precipitation dissolution equation. The above calculation results of chemical reactions and the calculated soil water content results are substituted into the one-dimensional vertical convection-dispersion equation, and then the tridiagonal equation set is solved to obtain the salt solution concentration value c1 at the end of this period. According to the value of c1, the adsorption amount and precipitation amount of the soil profile are calculated to obtain the updated salt solution concentration value c2 after the reaction. If the iterative accuracy requirement is met, the calculation for the next moment is carried out; otherwise, it is necessary to recalculate the adsorption concentration and precipitation concentration of each salt component in the soil during the chemical reaction according to c2, and update the salt solution concentration value until the iterative accuracy requirement is met. If chemical reactions such as adsorption and precipitation are not considered, the second term (adsorption term) and the third term (precipitation phase) on the left side of equation (4) in the convection-dispersion equation are both zero. The set soil parameter values and solute transport parameters are substituted into the one-dimensional vertical convection-dispersion equation, the equation is discretized, and the tridiagonal equation is solved to obtain the salt solution concentration value at the end of the period. Then, the salt solution concentration value at the end of this period is used as the initial salt solution concentration for the next moment, and the above cycle is repeated to carry out the calculation for the next moment.

[0103] The above steps are cyclically calculated until the set last simulation time point is reached and stopped, and finally the spatio-temporal distributions of water and salts can be obtained.

[0104] Through the above technical solutions, the present invention proposes a water-salt simulation method considering multiple chemical reaction processes, having the following technical effects:

[0105] 1) By numerical simulation, the distributions of water and salts in the soil, especially the distributions of various salt ions, are calculated, which is convenient for analyzing the distribution laws of water and salts in the soil.

[0106] 2) On the basis of the traditional water-salt transport model, multiple chemical reaction processes such as ion adsorption and precipitation / dissolution are considered. The description of water-salt transport in the soil is more comprehensive and accurate, which is conducive to truly reflecting the water-salt transport process in farmland saline soil.

[0107] 3) On the basis of the water stress on soil evaporation, the influence of soil salt crystallization on the evaporation rate is considered, making the calculation of soil evaporation more accurate. BRIEF DESCRIPTION OF THE DRAWINGS

[0108] The above and / or additional aspects and advantages of the present invention will become apparent and be readily understood from the following description of embodiments in conjunction with the drawings, in which:

[0109] Figure 1 It is a flowchart of water simulation calculation provided in Embodiments 1 and 2.

[0110] Figure 2 It is a salt simulation calculation process provided in Embodiments 1 and 2.

[0111] Figure 3 It is a schematic diagram of the soil column test device provided in Embodiment 1.

[0112] Figure 4 It is the simulated and measured values of the salt solution in Example 1.

[0113] Figure 5 It is the simulated and measured values of the salt solution in Embodiment 2.

[0114] Figure 6 It is the simulated and measured values of the actual evaporation rate and cumulative evaporation amount of the soil in Embodiment 3.

[0115] Figure 7 For Na in the soil in Embodiment 3 + , Ca 2+ , Cl - and SO4 2- It is a diagram showing the profile distribution of ions after 260 h of evaporation. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0116] In order to more clearly understand the above objects, features and advantages of the present invention, the present invention will be further described in detail below in conjunction with the drawings and specific embodiments. It should be noted that, without conflict, the embodiments of the present invention and the features in the embodiments can be combined with each other.

[0117] In the following description, many specific details are set forth in order to fully understand the present invention. However, the present invention can also be implemented in other ways different from those described herein. Therefore, the protection scope of the present invention is not limited by the specific embodiments disclosed below.

[0118] Embodiment 1

[0119] As Figure 1 and 2, in view of this, this embodiment proposes to provide a water-salt simulation method considering multiple chemical reaction processes, including the following steps:

[0120] Step 1: Collect the initial condition information of the soil profile in the soil column experiment (water content, salt content, basic physical parameters of the soil, solute parameters, etc.).

[0121] In this embodiment, the initial water content of the soil is the saturated water content of 0.409 cm 3 / cm 3 , the initial salt ion solution concentration is 0.008 g / L (Na + ), 0.031 g / L (Ca 2+ ), and 0.01 g / L (Cl - ), and the initial adsorbed salt ion content of the soil is 0.025 g / kg soil (Na + ), 2.0 g / kg soil (Ca 2+ ).

[0122] The test content is the study of the competitive adsorption of multi-component salt ions in the soil, mainly studying the penetration process of different salt ions (Na + , Ca 2+ , and Cl - ) in the soil under the action of competitive adsorption. As Figure 3 , the test device of this embodiment mainly includes a peristaltic pump, a soil column, and a collection device. Among them, the peristaltic pump is a water supply device with a constant flow rate. The soil column is 11.6 cm high and 2.6 cm in inner diameter. The collection device is used to collect the salt solution flowing through the soil column. The infiltrating salt solution is a mixed solution of 2 g / L NaCl and 5 g / L CaCl2. The test soil in the experiment is sandy loam, taken from the surface soil of the field in Aksu area, Xinjiang (the soil sampling depth is 0-30 cm). First, air-dry the soil in the laboratory and pass it through a sieve with a pore size of 2 mm. Measure the water content of the air-dried soil after sieving, fill the soil column according to the designed soil bulk density. First, use the peristaltic pump to transport deionized water to the soil column, collect the outflow liquid at the outlet end, wait until the outflow liquid flow rate is constant, then replace the deionized water with the salt solution, start recording the time, and continue to use the peristaltic pump to transport the salt solution to the soil column at a constant flow rate (0.02 cm / min). At the same time, collect the outflow liquid every 0.5 h at the outlet end and continuously collect the outflow liquid for 7 h. The salt ion concentration in the outflow liquid is measured using an ion chromatograph of model ECO-IC-940. The main salt ions measured include: Na + , Ca 2+ , and Cl - . Table 2 shows the basic soil parameters used in the soil column experiment. Table 3 shows the solute parameters.

[0123] Table 2 Basic soil parameters used in the soil column experiment

[0124]

[0125] Table 3 Solute Parameters

[0126]

[0127] Step 2: Set the boundary conditions for the simulation.

[0128] Since the upper outlet of the soil column is connected to the atmosphere, a constant head boundary (matrix potential is 0) is adopted for the upper boundary in the simulation, and a flux boundary is adopted for the lower boundary. The infiltration velocity is kept consistent with the actual experimental velocity (0.02 cm / min). The simulated infiltration time is 450 minutes, and an observation point is set at the upper boundary outlet to output the concentrations of Na + , Ca 2+ and Cl - in the outflow liquid.

[0129] Step 3: Calculate the hydraulic conductivity K and the specific water capacity C according to the soil parameter information in Step 1. Combine the initial soil moisture content condition in Step 1 and the soil moisture boundary condition in Step 2, and substitute the above results into the one-dimensional vertical mixed Richards equation. Discretize it according to the implicit difference format to obtain:[[]]

[0130]

[0131] where Δz is the spatial step (cm); Δt is the time step (d); the subscript i is the spatial node number; the superscript j is the time node number.

[0132] Both the soil moisture content and the matrix potential are contained in the one-dimensional vertical mixed Richards equation, and certain transformation treatments are required during the solution process. The modified Picard iteration can achieve this purpose. That is, assume that the soil moisture content and the matrix potential in the M-th iteration process are θ M and h M , respectively. Use the formula:[[]]

[0133]

[0134] where C is the specific water capacity (cm -1 ). The difference format only about the matrix potential h can be obtained:[[]]

[0135]

[0136] where

[0137]

[0138]

[0139]

[0140] When the upper boundary is a constant head, a1 = 0, b1 = 1, c1 = 0, d1 = h0(t); when the upper boundary is a flux boundary, Darcy's law can be discretized as:

[0141] a1 = 0, However, when the flux at the upper boundary changes violently, Darcy's law is prone to cause non - convergence of the numerical solution. Therefore, the water balance formula is often used to replace the solution:

[0142]

[0143] After discretizing the formula at the upper boundary, we get:

[0144]

[0145]

[0146]

[0147] In this embodiment, the upper boundary is a constant - head boundary. Therefore, a1 = 0, b1 = 1, c1 = 0, d1 = h0(t).

[0148] When the lower boundary is a constant head, a N+1 = 0, b N+1 = 1, c N+1 = 0, d N+1 = h L (t); when the lower boundary is a flux boundary, Darcy's law can be discretized as: c1 = 0, When the lower boundary is free - drainage, a N+1 = - 1, b N+1 = 1, c N+1 = 0, d N+1 = 0.

[0149] In this embodiment, the lower boundary is a flux boundary. Therefore c1 = 0,

[0150] After the above - mentioned discretization process of the one - dimensional vertical mixed - type Richards equation, N + 1 linear equations are obtained, and all the equations form a tridiagonal matrix:

[0151]

[0152] The tridiagonal matrix is solved by the chasing method for programming calculation.

[0153] Step 4: Calculate the exchangeable cation concentration S according to the solute transport parameters in Step 1 kSubstitute it into the one-dimensional vertical convection-dispersion equation, and then the multi-component solute transport model considering ion competitive adsorption can be obtained.

[0154] The convection-dispersion equation is discretized according to the implicit difference scheme to obtain:

[0155]

[0156] After organizing the formula, a tridiagonal equation is obtained:

[0157]

[0158] Among them,

[0159]

[0160]

[0161]

[0162]

[0163] When the upper boundary is a concentration boundary, a′1 = 0, b′1 = 1, c′1 = 0, When the upper boundary is a concentration flux boundary, formula (35) can be discretized as:

[0164] a′1 = 0, However, formula (35) is prone to numerical instability, and better results can be obtained through the salt balance formula of the upper boundary:

[0165]

[0166] After discretization, it is obtained:

[0167]

[0168]

[0169] In this embodiment, the upper boundary is a concentration boundary, so a′1 = 0, b′1 = 1, c′1 = 0,

[0170] When the lower boundary is a concentration boundary, a′ N+1 = 0, b′ N+1 = 1, c′ N+1 = 0, When the lower boundary is a concentration flux boundary, formula (37) can be discretized as:

[0171] c′ N+1 = 0, When the lower boundary is free of salt discharge, a′N+1 =-1, b′ N+1 =1, c′ N+1 =0,d′ N+1 =0.

[0172] After the above discretization process, N+1 linear equations are obtained. The tridiagonal coefficient matrix is solved by the simultaneous equations using the pursuit method.

[0173] Step 5: Obtain the soil profile moisture content and salt (Na) at the end of the first period through steps 3 and 4. + , Ca 2+ and Cl - ) information, and the soil profile moisture content, salt (Na + , Ca 2+ and Cl - ) information as the initial value of the next stage, and then enter the calculation of the next period until the last simulation time point set is stopped.

[0174] This example simulates the competitive adsorption of Na + , Ca 2+ and Cl - The penetration process in the soil. The simulated value and measured value of the salt solution are as follows Figure 4 It can be seen that the water-salt simulation method of the present invention, which takes multiple chemical reaction processes into consideration, can better describe and simulate the migration process of multi-component salt ions in soil under the influence of competitive adsorption.

[0175] Example 2

[0176] The initial soil moisture content in this example is 0.409 cm 3 / cm 3 , the initial salt ion concentration is 0.008g / L (Na + )、0.031g / L(Ca 2+ ) and 0.01 g / L (Cl - ), the initial adsorption salt ion content of the soil is 0.025g / kg soil (Na + ), 2.0 g / kg soil (Ca 2+ ).

[0177] The experimental device of this embodiment mainly includes a peristaltic pump, a soil column and a collection device. The peristaltic pump is a water supply device with a constant flow rate. The soil column is 11.6 cm high and has an inner diameter of 2.6 cm. The collection device is used to collect the salt solution after it flows through the soil column. The infiltration salt solution is a mixture of 2g / L NaCl and 2g / L CaCl2. The simulated value and the measured value of the salt solution are shown in Figure 2. Figure 5 .like Figure 4-5 According to the statistical analysis results, Na +, Ca 2+ and Cl - The determination coefficient R of the simulated and measured values of the ion penetration curve 2 is above 0.92, and the model efficiency NSE is greater than 0.9.

[0178] Other technical features in this embodiment are the same as those in Embodiment 1 and will not be elaborated here one by one.

[0179] Embodiment 3

[0180] In this embodiment, the initial soil moisture content is the saturated moisture content of 0.409 cm 3 / cm 3 , Na + , Ca 2+ , Cl - and the initial values of the liquid phase and adsorbed phase concentrations of SO4 2- are shown in Table 4.

[0181] Table 4 The initial concentrations of the liquid phase and adsorbed phase of Na + , Ca 2+ , Cl - and SO4 2- in the model

[0182]

[0183] In this embodiment, the inner diameter of the soil column is 19 cm and the height is 35 cm. An inlet chamber and a gravel layer are provided at the lower part of the soil column. Three salt treatments S0, S1, and S2 are set, representing soil salt contents of 5.4 (background value), 10, and 15 g / kg respectively, and 6 replicates are set, with a total of 3×6 = 18 soil columns. The test soil is air-dried and sieved through a 2 mm sieve, and the moisture content of the air-dried soil is measured by the drying method to be 0.01 cm 3 / cm 3 , and the soil column is filled according to the designed bulk density (1.55 g / cm 3 ), and the filling height of the soil column is 30 cm. Then, a certain amount of Xinjiang salt crust is mixed with deionized water to prepare infiltration solutions with different concentrations. The infiltration solution for the S0 treatment is deionized water without salt, and the soil columns are saturated from the bottom with a Mariotte bottle respectively. After saturation, the soil columns are placed under an infrared lamp (200 W), and the height of the infrared lamp from the soil surface is 50 cm to simulate the evaporation process in saline soil considering precipitation / dissolution reactions. Before evaporation, the initial contents of salt ions in the saturated soil with different concentrations of salt solutions are shown in Table 5.

[0184] Table 5 Initial contents of salt ions in saturated soil under different test treatments (mmol / kg soil)

[0185]

[0186] Step 2: Set the simulation boundary conditions.

[0187] In the simulation of soil moisture, the upper boundary adopts a flux boundary, and the flux is the evaporation rate after the measured maximum soil evaporation rate is affected by soil water-salt stress. The lower boundary adopts the first-type boundary, and its boundary water content is the measured soil water content. In the simulation of multi-component salts, the upper boundary adopts a solute flux boundary, and the ion concentration in the water flux at the upper boundary during the evaporation process is 0. The lower boundary adopts the third-type boundary. The simulation evaporation time is 260 h, and the liquid phase, adsorbed phase, and precipitation phase concentrations of salt ions and the actual soil evaporation rate are output.

[0188] Step 3

[0189] Calculate the hydraulic conductivity K and the specific water capacity C according to the soil parameter information in Step 1, and combine the initial soil water content condition in Step 1 and the soil moisture boundary condition in Step 2. Substitute the above results into the one-dimensional vertical mixed Richards equation and discretize it according to the implicit difference format to obtain:

[0190]

[0191] In the formula, Δz is the spatial step (cm); Δt is the time step (d); the subscript i is the spatial node number; the superscript j is the time node number.

[0192] The one-dimensional vertical mixed Richards equation contains both soil water content and matrix potential, and certain transformation treatments need to be carried out during the solution process. The modified Picard iteration can achieve this purpose. That is, assume that the soil water content and matrix potential in the M-th iteration process are θ M and h M , respectively, and use the formula:

[0193]

[0194] where

[0195]

[0196]

[0197]

[0198] When the upper boundary is a constant head, a1 = 0, b1 = 1, c1 = 0, d1 = h0(t); when the upper boundary is a flux boundary, Darcy's law can be discretized as:

[0199] a1 = 0, However, when the upper boundary flow rate changes drastically, Darcy's law is prone to causing non-convergence of numerical solutions. Therefore, the water balance formula is often used to replace the solution:

[0200]

[0201] After discretizing the upper boundary, the formula becomes:

[0202]

[0203]

[0204]

[0205] In this embodiment, the upper boundary is a flux boundary. Therefore, a1 = 0.

[0206] When the lower boundary is a constant head, a N+1 = 0, b N+1 = 1, c N+1 = 0, d N+1 = h L (t); when the lower boundary is a flux boundary, Darcy's law can be discretized as: c1 = 0. When the lower boundary is free drainage, a N+1 = -1, b N+1 = 1, c N+1 = 0, d N+1 = 0.

[0207] The lower boundary adopts the first type of boundary, and its boundary water content is the measured soil water content.

[0208] After the above discretization process of the one-dimensional vertical mixed Richards equation, N + 1 linear equations are obtained, and all the equations form a tridiagonal matrix:

[0209]

[0210] The tridiagonal matrix is solved by the chasing method for programming calculation.

[0211] Step 4: Calculate the exchangeable cation concentration S k and the ion precipitation concentration according to the solute transport parameters in Step 1 and substitute them into the one-dimensional vertical convection-dispersion equation, then the multi-component solute transport model considering ion competitive adsorption is obtained.

[0212] Discretize according to the implicit difference format to obtain:

[0213]

[0214] The formula is rearranged to obtain a tridiagonal equation:

[0215]

[0216] Wherein,

[0217]

[0218]

[0219]

[0220]

[0221] When the upper boundary is a concentration boundary, a′1 = 0, b′1 = 1, c′1 = 0, When the upper boundary is a concentration flux boundary, formula (35) can be discretized as:

[0222] a′1 = 0, However, formula (35) is prone to numerical instability, and better results can be obtained through the salt balance formula of the upper boundary:

[0223]

[0224] After discretization, we get:

[0225]

[0226]

[0227] The upper boundary adopts a solute flux boundary, and the ion concentration in the water flux of the upper boundary during the evaporation process is 0.

[0228] When the lower boundary is a concentration boundary, a′ N+1 = 0, b′ N+1 = 1, c′ N+1 = 0, When the lower boundary is a concentration flux boundary, formula (37) can be discretized as:

[0229] c′ N+1 = 0, When the lower boundary is free salt drainage, a′ N+1 = -1, b′ N+1 = 1, c′ N+1 = 0, d′ N+1 = 0.

[0230] The lower boundary adopts the third type of boundary, that is, the concentration flux boundary. Therefore

[0231] c′ N+1 =0,

[0232] After the above discretization process, N+1 linear equations are obtained. The tridiagonal coefficient matrix is solved by the simultaneous equations using the pursuit method.

[0233] This example simulates the evaporation process in saline soil taking into account the precipitation / dissolution reaction of CaSO4 and NaCl. The simulated values of the actual evaporation rate and cumulative evaporation of the soil are compared with the measured values. Figure 6 The simulated values are basically consistent with the measured values, and the three stages of evaporation can be simulated. + , Ca 2+ , Cl- and SO4 2- The profile distribution of ions after 260 hours of evaporation is as follows Figure 7 The sum of the simulated concentrations of each ion in the liquid phase and the precipitation phase is basically consistent with the measured ion concentration in the soil.

[0234] The foregoing description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Those skilled in the art will readily appreciate that various modifications and variations of the present invention are possible. Any modifications, equivalent substitutions, or improvements made within the spirit and principles of the present invention are intended to be within the scope of protection of the present invention.

Claims

1. A soil water and salt simulation method considering the chemical reaction process of multi-component salts, characterized in that include Obtain the specific values of parameters associated with the soil water and salt simulation process; A water-salt transport model is constructed that considers competitive adsorption and precipitation / dissolution reactions of multi-component salts. Initial conditions and upper and lower boundary conditions are set for the model operation. Specific values of all the parameters are input into the water-salt transport model for solution to obtain the water-salt distribution characteristics on the soil profile at a specified time. The water-salt transport model is used to describe the functional relationship between the parameters and the soil water-salt content. The water-salt transport model includes a one-dimensional vertical mixed Richards equation to simulate water movement and a one-dimensional vertical convection-diffusion equation to simulate the convection and mechanical diffusion process of salt in the soil. The process of solving the water movement using the one-dimensional vertical mixed Richards equation includes at least: (1) Input the variables to be solved, namely, hydraulic conductivity, specific water density, and upper and lower boundary conditions. The upper boundary evaporation takes into account the inhibitory effect of salt crystallization on the soil surface. The chase method is used to program and calculate the discretized linear tridiagonal equation system to obtain the matrix potential or soil profile moisture content at the end of the time period. (2) Determine whether the matrix potential meets the iteration accuracy requirement. If yes, use the obtained matrix potential as the initial value at the next moment and repeat the above cycle to perform the calculation at the next moment. If no, proceed to step (3). (3) using the matrix potential to solve the hydraulic conductivity and specific water density again, repeating the above calculation method until the accuracy requirements are met; The one-dimensional vertical convection-diffusion equation is used to simulate the convection and mechanical diffusion process of salt in the soil, which at least includes: (1) First, determine whether to consider the chemical reaction process of salt; if "yes", proceed to step (2); if "no", proceed to step (3); (2) The adsorption concentration and precipitation concentration of each salt component during the chemical reaction and the soil profile moisture content are input, and the tridiagonal coefficient matrix is solved by the chase method through simultaneous equations to obtain the salt solution concentration value at the end of the period; the salt adsorption amount and precipitation amount of the soil profile are calculated based on the salt solution concentration value to obtain the updated salt solution concentration value after the reaction; (3) Solve the corresponding variables of the parameters directly obtained from the soil water-salt simulation process and input them, discretize the equation, and solve the tridiagonal equation to obtain the salt solution concentration value at the end of the time period. The salt solution concentration value at the end of this time period is used as the initial salt solution concentration at the next moment, and the above cycle is repeated to perform calculations at the next moment.

2. The soil water and salt simulation method according to claim 1, wherein It includes discretizing the soil profile into several spatial grid nodes, setting the simulation time and step size, and using the implicit finite difference method to solve and calculate the moisture and salt of the discrete nodes.

3. The soil water and salt simulation method according to claim 2, characterized in that The process of solving the water movement using the one-dimensional vertical mixed Richards equation also includes step (4). When the number of iteration cycles in step (3) exceeds 20 and fails to meet the iteration accuracy requirement, the step size of the calculation period is reduced to one third, and the above calculation method is repeated again for calculation.

4. The soil water and salt simulation method according to claim 1 or 3, characterized in that In the process of solving the water movement using the one-dimensional vertical mixed Richards equation, the hydraulic conductivity and the specific water density are calculated by directly obtaining the associated parameters of the soil water-salt simulation process.

5. The soil water and salt simulation method according to claim 1, characterized in that A one-dimensional vertical convection-diffusion equation is used to simulate the convection and mechanical diffusion of salt in the soil. The adsorption concentration and the precipitation concentration are respectively calculated by directly obtaining the associated parameters of the soil water-salt simulation process using the competitive adsorption equation and the precipitation dissolution equation.

6. The soil water and salt simulation method according to claim 1 or 5, characterized in that, The convection and mechanical diffusion of salt in the soil are simulated using a one-dimensional vertical convection-diffusion equation. Step (2) further includes: updating the concentration value of the salt solution to meet the iteration accuracy requirement and performing calculation at the next moment; Otherwise, the adsorption concentration and the precipitation concentration are recalculated according to the updated salt solution concentration value, and the salt solution concentration value is updated again until the iteration accuracy requirement is met.

Images