A Research About One-dimensional Chloride Ion Erosion In Concrete Under Drying-wetting Cycles

. Reinforced concrete structures are now widely used in China, and with the exploitation of marine resources, many concrete structures are in a relatively harsh service environment. Concrete in wave and tidal stream areas is in a dry and wet cycle state. When unsaturated concrete is in a dry and wet cycle state, chloride ions will invade the concrete by diffusion and convection, which accelerates the erosion of the reinforcement within the concrete and degrades the concrete performance. Therefore, the study of chloride ion erosion in concrete under dry and wet cycles is particularly important. Fick's law is a good predictor of the diffusion of chloride ions in saturated concrete with stable boundaries, but it is difficult to obtain satisfactory results for concrete under dry and wet cycles. There are many difficult parameters in the microscopic model that make programming the calculations more difficult. In this paper, we propose to use the radial basis function matching point method to solve the problem. It is found that the error is within acceptable limits and can be used to calculate the one-dimensional erosion of chloride ions under dry and wet cycles through solving the validation algorithm.


Introduction
Reinforced concrete structures combine the advantages of steel and concrete, are widely used all over the world. However, due to the concrete has stable chemical properties, which can protect the internal steel bars, people have neglected various factors that affect the durability of concrete structures for a long period of time. As the concrete structure has been used so far, people have found that its service life may not be able to meet expectations. The deterioration of the durability of reinforced concrete structures can be caused by several factors, ranging from inadequate resistance during design to the encounter of extremely adverse and unanticipated load changes. However, the corrosion of steel bars inside concrete is one of the most important factors affecting its durability. The intrusion of chloride ions into concrete, which triggers reinforcement corrosion, is an important cause of deterioration in the durability of concrete structures.
As a large developing country, China is facing unprecedented infrastructure development and cannot afford to tackle the problem of 'major repairs' soon after the 'major construction' that has taken place. The importance of durability is graphically illustrated by the 'fivefold law' of American scholars. When designing a new project, for every $1 saved on reinforcement protection, an additional $5 will be added to the cost of measures taken when reinforcement corrosion is found, an additional $25 will be added to the cost of maintenance when concrete is cracked, and an additional $125 will be added to the cost of maintenance when serious damage occurs [1]. This shows that later maintenance is extremely costly in terms of both money and manpower. Thus, the study of the chloride ion attack process is helpful in improving the durability of concrete.

Background research
With the continuous exploitation of marine resources, concrete structures are gradually being used in wharves, submarine tunnels, etc. As a result, many concrete structures are in the marine environment. They have harsh service conditions due to a variety of factors such as high content of chloride and salt in the water and alternating wet and dry conditions. The highly alkaline pore solution in the concrete passivates the reinforcement and allows it to resist rusting, but the intrusion of chloride ions can partially or destroy this passivation. Even with a thick protective layer, it is only a matter of time before the chloride ions invade the internal reinforcement from the surface layer. These problems caused by chloride ions should not be underestimated. Depending on the location and exposure conditions, such concrete elements can be divided into three zones, which are salt spray zones, wave and tidal zones, and underwater zones. Of these, the wave and tidal zones are subject to a wet and dry cycle. Under such environmental conditions, the penetration rate of chloride ions into the concrete is faster [2]. Under dry and wet cycles, chloride ions invade concrete mainly by diffusion under concentration gradients and capillary absorption under moisture gradients. The main transmission method is different in different situations.
The main absorption is through capillary absorption if the concrete pores are not saturated; if the concrete pores are saturated, diffusion is the main method. For areas such as wave and tidal areas, which are under the wet and dry cycle, the water content on the surface is constantly changing and the internal water content is more stable, it is mainly absorbed by capillary absorption on the surface and diffusion inside.
A lot of work has been done both china and abroad on the transport of chloride ions in concrete. There are generally two ways of trying to establish a computational model, the first one is a macroscopic model using Fick's law and the law of conservation of mass, another one is a microscopic model using the Nernst-Plank equation and the Nernst-Einstein equation. Because the parameters are difficult to determine under microscopic conditions, macroscopic models are used most. For example, Saetta established a control equation for unsaturated concrete that includes water infiltration, heat transfer, and the effect of chloride ions. Costa derived the variation process of chloride ion from Fick's second law in the measured environment. Ababneh et al. also proposed a transport model of chloride ion in terms of unsaturated concrete. Based on this, Suwito proposed a new coupled model of seepage and diffusion. At the same time, there have been some experiments that attempt to use microscopic models. Some scholars suggested that a non-linear diffusion equation can be used to calculate the seepage of water in porous media. Most of above methods use the finite difference method or the finite element method. The finite difference method is poorly adapted to irregular boundaries, while the finite element method, although adapted to irregular boundaries, may require a finer mesh for the calculation, which increases the computational cost. Therefore, a better numerical method needs to be proposed to solve these problems.

Overview of the meshless method
The Finite Element Method is a well-developed numerical method and plays an important role in various research fields and in university teaching. With the development of relevant software in the business, the finite element method has become easy to use and intuitive. The accuracy of the results obtained is relatively good in general. However, as more field of research expand, there are obvious limitations to the finite element method. Because the finite element method relies on meshing, when problems are encountered that do not correspond to the original mesh lines, such as deformation under highspeed impact, the mesh has to be reconstructed, which greatly increases the cost of computing and affects accuracy. Because of the inconvenience of traditional numerical methods, the meshfree method which independent of mesh, has been proposed.
The starting point of the meshless method is to solve without dividing the mesh, by taking discrete points over the physical region. Due to mesh-independent, it is easier to solve for irregular regions than the finite element method and can be better adapted to distortion problems. Furthermore, when using the meshless method for analysis, only node information is required. When there is a need to add nodes to a region, only the node information is added, which greatly facilitates the processing work. The meshless method offers more methods than the finite element method. For example, some methods are moving least squares, compactly radial basis functions, point interpolation, reproducing kernel approximation, unit decomposition, etc. The algebraic equations established by the finite element method are in the use of the variational principle, whereas the solved algebraic equations established by the meshless method are in the use of weighted Residual Method with Compactly Supported Trial Function, etc. The method that meshless used to solve the essential boundary conditions, different with the finite element method, usually using the Lagrange multiplier method and the penalty function method.
Most numerical solution to a partial differential equation first requires that the infinite dimensional problem be approximately transformed into a finite dimensional problem, such that the solution to the problem can be expressed in terms of a finite number of parameters, which are then determined by the partial differential equation and the boundary and initial conditions, or conditions associated with the variational problem. Typically, the true solution u(x) of a problem can be approximated expressed through a finite number of parameters as: Therefore, the shape function in the finite element method is generally an interpolation function. This is one of the advantages of the finite element method, which makes it simple to deal with essential boundary conditions. Earlier shape functions in meshless methods often cannot satisfy this condition, which causes inconvenience to the solution of the meshless method. A great deal of work has been carried out by researchers to overcome this drawback. The most basic and central problem of the meshless method is the construction of the shape function. The construction of the shape function should satisfy the following conditions: Firstly, linear regenerative property, which enables the shape function to satisfy the piecewise test and is an important condition for the convergence of the meshless method. Secondly, it is desirable to satisfy the Kronecker delta function property, which will allow the meshless method to satisfy the essential boundary conditions accurately. In addition, the shape function should be based on a local approximation so that a banded matrix can be generated to improve computational efficiency.

Overview of the radial basis function method
As mentioned above there are many methods of generating shape functions by the meshless method, the radial basis function method will be used in this paper. Radial basis functions (RBFs) are widely used in computational science, such as partial differential equation solving, CAD, computer graphics, etc. Although the radial basis function interpolation method is not constrained by scalar or vector parameters, it is also not easy to use in other interpolation methods.

Domain decomposition method
The original idea can be traced back to the famous Schwarz alternating method proposed by the German mathematician H. A. Schwarz in 1870, but Schwarz's intention was to borrow the alternating method to prove the existence and uniqueness of solutions to non-regular elliptic equations. It was not until the 1950s that Schwarz's method was used in computation, but it failed to attract the attention of mathematicians. In the last decade, traditional algorithms have been challenged by the introduction and increasing popularity of parallel computers, and the classical serial computing landscape is not adapted to parallel computers. How to construct highly parallel algorithms is the key to increasing the speed of computation, and the scientific and engineering problems facing practice are so vast that improvements in computational power depend on developments in both computers and computational methods, and it is in this context that the region decomposition algorithm was produced.
The shape of the sub-domain is as regular as possible, so the solution of the original problem is transformed into a solution on the sub-domain. Region decomposition algorithms are of particular interest because they offer advantages unmatched by other methods:①They reduce the size of the computation by reducing the size of a large problem to many smaller problems; ② They allow the use of familiar fast algorithms such as the Fast Fourier Transform (FFT), spectral methods, etc., if the shape of the subregion is regular (e.g. rectangular), or there is already efficient software available for solving such regular problems;③The use of locally consistent grids is allowed, instead of using the whole consistent grid, and each sub-domain can be computed using different discrete methods;④Different mathematical models are allowed to be chosen for different sub-domains so that the overall model is more suitable for the actual engineering physics; ⑤The algorithms are highly parallel, i.e. the main steps of the computation are performed independently in each subdomain; ⑥ For symmetric regions problems there are simpler region decomposition algorithms. Generally, according to whether the sub-domains overlap or not, region decomposition algorithms can be divided into overlapping region decomposition algorithms and nonoverlapping region decomposition algorithms.
After decades of development the combination of area decomposition methods and finite elements has become more mature, mainly Dirichlet-Neumann methods, Neumann-Neumann methods, Schwarz methods, fictitious area methods and so on. As mentioned before, the meshless method is then not constrained by the shape of the region, but in order to solve it in parallel, on the other hand, it is also hoped that the region decomposition can transform solving large-scale problems into solving multiple small sub-regional problems, thus overcoming the problem of the configuration matrix being an asymmetric full array when solving large-scale problems with the usual full-domain radial basis matching point method.

Overview of the matching point method
In general, the true solution u(x) of a problem can be approximated by a set of discrete nodes (x=1, 2, ...N) and their corresponding functions as The problem we are trying to solve can be considered as a problem given initial and boundary conditions. Then the control equations and boundary conditions   (6) will be brought in, giving  WhenN Ω + N Γ > N, it is a system of super-stationary equations that needs to be solved by least squares. whenN Ω + N Γ = N, then the solution node approximation is a linear system of equations.

Calculation scheme
If the diffusion equation is: and the boundary conditions and initial conditions are S is the number of subdomains to be divided, The problem in the first subdomain can be viewed as A simultaneous equation can be obtained   The basic idea of the matching point method is to make the points in the region satisfy the control equations and the nodes on the boundary satisfy the boundary conditions. The overall time region can first be divided according to ∆t. In each sub-domain of ∆t, the number of points can be determined after determining the time interval and spatial interval of the matching points. When programming the calculation, the information on the mating points is first stored in space and time, and then the radial basis functions are used to calculate the points inside the medium and the points at the boundary, depending on where the mating points are located. The coefficient matrix A and the column vector b of the system of equations are then solved for in the first sub-domain using the PLU decomposition. At this point the coefficient vectors are obtained and can be stored in a separately created matrix. At this point, the result of extracting the last set of matching points within the first sub-domain can be used as the initial condition for the next set of matching points. The loop is thus written so that all subsequent subdomains can be calculated. As the direct result obtained is a vector of coefficients for each group, a result extraction function will also need to be written to obtain the result. Note that if you want to implement a wet and dry loop, the boundary conditions change when calculating the second half of the subdomains, so you need to recalculate part of the coefficient matrix A.
The radial basis function method is used to calculate the validation algorithm, for example: Assume a concrete flat slab with a planar dimension of 1*1, the diffusion mode is one-dimensional diffusion and the humidity diffusion coefficient within the concrete slab has a onedimensional non-homogeneity, i.e.
Boundary conditions and initial conditions: The humidity-valued analytic solution of its onedimensional diffusion takes the form: The Gaussian radial basis function ( ( ) = − 2 , form parameter c=4) was validated using a region decomposition method, in which the time domain was divided into 50 parts, each with time t = 0.02, and N = 33 points (Δ = 0.01, Δ =) were allocated on each time. Using Matlab programming, the non-homogeneous humidity diffusion problem in the region [1,2] *[0,1] was solved.
Plots of the long distribution of humidity with relative errors were calculated.

Scheme for solving one-dimensional wetting processes
For unsaturated concrete in the wave splash zone, the transport of chloride ions during the initial phase of the wetting and drying cycle relies heavily on capillary absorption under a moisture gradient. We therefore simplify the model as follows: during the wetting process, seawater enters the concrete by capillary action; during the drying process, water evaporates outwards through the concrete pores, but the chloride ions remain in the pores, causing the chloride ion concentration to increase. We first obtain a model of the distribution of the moisture field by numerical calculation, and then a model of the distribution of the concentration field of chloride ions in concrete by calculation.
We first consider a one-dimensional model of a single wetting process and, assuming a constant humidity diffusion coefficient, obtain the humidity field distribution by preparing a Matlab program.  Figure 3 plots the calculated approximate humidity field distribution over time and depth. It can be seen from the graph that the relative humidity at the same moment increases with depth; the relative humidity at the same depth increases with time. This is consistent with the law of humidity diffusion.

One-dimensional drying process solution scheme
We then consider a one-dimensional model for a single drying process and, assuming that the humidity diffusion coefficient is constant, obtain the humidity field distribution by programming Matlab.
(1, ) 0; (2, t) e ; ( , 0) 1.  Figure 4 plots the calculated approximate moisture field distribution with time and depth. It can be seen from the graph that the moisture content of the concrete on the drying surface decreases with time, in accordance with the laws of the drying process.

Conclusion
The introduction of radial basis functions into the collocation point method to find the heat conduction equations has the advantages of simplicity of form, convenience of calculation and ease of programming. It overcomes the dependence of the traditional finite element method on the mesh and completely or partially eliminates the mesh, thus reducing the large amount of time and computational cost spent in meshing and reconstruction. Theoretical analysis and numerical experiments show that the global radial basis function configuration method leads to an algebraic system with an excessive number of conditions, and when the number of matching points increases to a certain level, the system of equations appears pathological. The radial basis collocation point method can be easily extended to the high-dimensional case and facilitates the treatment of high-dimensional problems, but the large number of conditions and the computational instability (especially when computing problems with Neumann boundary conditions) limit the method.
The area decomposition algorithm decomposes the computational region into several sub-regions and transforms the large problem into a small problem to be solved. Therefore, the radial basis meshless area decomposition algorithm, which combines the area decomposition algorithm and the radial basis meshless method to form the radial basis meshless area decomposition algorithm, is compatible with the advantages of both. Not only does it have a significant advantage over the finite element area decomposition algorithm in terms of computational time, but it also reduces the condition number of the configuration matrix while maintaining good computational accuracy, which helps to overcome the pathological phenomenon of the system of equations caused by the radial basis matching point method and improves computational efficiency.
There is no perfect method for the selection of the form parameter of the global radial basis function and the radius of the tightly-branched radial basis function, but according to the author's experience, the following methods can be used to find a basic suitable range: roughly estimate the gradient of the function to be solved, the parameters should be selected so that the radial basis function is comparable to the gradient of the function to be solved, such as the gradient is larger, the parameters should be adjusted so that the radial basis function graph is sharper. For example, the use of Gaussian radial basis function can be adjusted to a larger shape, and the use of tightly-branched radial basis function should be tightlybranched radius adjusted to small. In addition, in the choice of radial basis functions, try to choose and solve the form of similar radial basis functions, for example, when the solution function for the exponential type, the use of Gaussian radial basis functions can get better results.
The Hermite type matching point method is an effective way to address the problem of excessive errors in calculating problems with Neumann boundary conditions. the Hermite type approximation function incorporates a derivative term which makes the derivative of the interpolation function at the boundary much more accurate and therefore helps to improve the accuracy at