The numerical simulation of hydraulic fracturing fracture propagation is the core content of hydraulic fracturing design and construction. Its simulation results directly affect the effect of fracturing, and can effectively guide the fracturing construction plan and reduce the construction risk. At present, two-dimensional or quasi-three-dimensional models are mainly used, but most of them are used to simulate the vertical fracture of hydraulic fracturing. There are errors in the application process. In this paper, a three-dimensional mathematical model, including an elastic rock mechanics equation and a material flow continuity equation, is established to simulate horizontal fracture propagation in shallow reservoirs. The emphasis of this paper is to propose a new method for solving equations. The basic idea of the iteration method has been proposed by previous scholars: Firstly, assuming that the initial pressure of each point in the fracture is uniform, the fracture height of each initial point can be obtained by using Equation (20). Using the initial height values, the pressure values at each point of continuous variation are calculated by Equation (16), and then the new fracture height values are obtained by Equation (20). Because of the equal initial pressure, this method leads to too many iterations in the later stages, which makes the calculation more complicated. In this paper, a new Picca iteration method is proposed. The iteration parameters are changed sequentially. Firstly, the distribution value of fracture height is assumed. Then, the pressure distribution value is calculated according to Equation (16). Then, the new distribution value of fracture height is obtained by bringing the obtained pressure distribution value into Equation (20). Then, the new distribution value of the fracture height is calculated according to Equation (16). The pressure distribution value completes an iteration process until the iteration satisfies the convergence condition. In addition, Sneddon’s model is introduced into the hypothesis of fracture height to obtain the maximum fracture height and assume that the initial fracture profile is a parabola. Finally, the proposed method can rapidly improve the convergence rate. Next, on the basis of investigating the solutions of previous equations, the Galerkin finite element method is used to solve the above equations. The new Picard iteration sequence method is applied to solve the height and pressure at different points in the fracture. By calculating the stress intensity factor, we can judge whether the fracture continues to extend or not, and then simulate the full three-dimensional horizontal fracture of the hydraulic fracturing expansion process. The infiltration process of three types of oil reservoirs in Daqing Changyuan oilfield is simulated. The results show that during the initial fracture stage, the radius and height of fractures increase rapidly, and the rate of increase slows down with the increase of construction time. The height and net pressure of each point in the fracture are unequal. The height and net pressure of the fracture in the wellbore reach the maximum, and gradually decrease to the front of the fracture. Compared with conventional fracturing, the fracturing-flooding percolation process has the characteristics of short fracture-making and large vertical percolation distance, which can greatly increase the swept volume of flooding fluid and thus enhance oil recovery. With the increase in the rock modulus of elasticity, the radius of fractures decreases and the height of fractures increases. With the increase in construction displacement, the radius of fractures hardly changes, the height of fractures increases, and the vertical infiltration distance of the fractures increases. It is suggested that the construction displacement should be 4.0 m3/min. In the range of fracturing fluid viscosity in the studied block, with the change of fracturing fluid viscosity, the change of fracture radius and height is not obvious. In order to further increase sweep volume, the fracturing fluid viscosity should be further reduced.
