## 2 Mathematical Model
### Circumferential Grout Pressure Distribution  
A coordinate system is established where the horizontal, vertical, and longitudinal directions are set as the \(x\)-, \(y\)-, and \(z\)-axis, respectively (Figure 2a). The grout pressure distribution of circular tunnels is usually symmetrical by the vertical \(y\)-axis. However, the impact of grout flow direction (top-down or bottom-up) on the grout pressure distribution needs to be discussed separately, owing to the gravity. Force analysis of the grout element during circumferential grout filling when the grout flow direction is upward is shown in Figure 2b. Figure 1: Schematic diagrams of the shield tail void and synchronous grouting ((**a**) three-dimensional annular void, (**b**) circumferential grout filling, (**c**) the longitudinal grout filling, (**d**) radial grout penetration). The rheological equation of Bingham fluid grout considering the temporal variation in grout viscosity can be expressed as [34]:  
\[\tau=\tau_{0}+\mu_{(t)}\gamma \tag{1}\]  
where \(\tau\) is the shear stress (Pa), \(\tau_{0}\) is the yield shear stress (Pa), \(\mu_{(t)}\) is the grout viscosity as a function of time (Pa-s), \(\gamma\) is the shear rate (s\({}^{-1}\)), which is the negative ratio of \(dv\) to \(dr\), \(dv\) is the circumferential grout filling velocity variation vector (m/s), \(dr\) is the grout element height along the tunnel radial direction (m). According to the existing research, the grout viscosity considering the temporal variation can be expressed as [20; 35]:  
\[\mu_{(t)}=\mu_{0}e^{Qt} \tag{2}\]  
where \(\mu_{0}\) is the initial viscosity (Pa-s), \(\xi\) is the time-varying coefficient of viscosity (s\({}^{-1}\)), and \(t\) is the grouting time (s). According to the conservation of mass, the volume of injected grout is equal to the volume of the filled annular void in a certain period of time, described as follows:  
\[qt_{z}=\frac{\Delta\alpha}{2}\Big{(}{R_{0}}^{2}-{R_{1}}^{2}\Big{)}\delta \tag{3}\]  
where \(q\) is the grouting rate (m\({}^{3}\)/s), \(t_{z}\) is the grout viscosity growth time, calculated from the injection of grout into the shield tail void (s), \(\Delta\alpha\) is the diffusion radian (rad) increment which is greater than 0, \(R_{0}\) is the external radius of the shield shell (m), \(R_{1}\) is the external radius of the segment (m), \(\delta\) is the thickness of the filled annular void along the longitudinal direction of the tunnel (m). During synchronous grouting in shield tunnels, the longitudinal grout filling is mainly conducted by the longitudinal movement of the shield machine. Therefore, \(\delta\) can be expressed as follows:  
\[\delta=v_{d}\cdot t_{y} \tag{4}\]  
where \(v_{d}\) is the shield advance velocity (m/s), \(t_{y}\) is the required time for the longitudinal grout filling in the shield tail void (s). The grout viscosity growth time \(t_{z}\) is calculated according to Equation (3):  
\[t_{z}=\frac{\Delta\alpha}{2q}\Big{(}{R_{0}}^{2}-{R_{1}}^{2}\Big{)}\delta \tag{5}\]  
The relationship between the grout viscosity and spatial distribution is established combined with Equation (2):  
\[\mu_{(t_{z})}=\mu_{\big{[}\frac{\Delta\alpha}{2q}({R_{0}}^{2}-{R_{1}}^{2}) \delta\big{]}}=\mu_{0}e^{\xi\big{[}\frac{\Delta q\beta}{2q}({R_{0}}^{2}-{R_{1} }^{2})\big{]}} \tag{6}\]  
Figure 2: Force analysis of grout filling upward ((a) force analysis under the coordinate system, (b) grout element force analysis). The stress equilibrium relationship of the grout element along the streamline is expressed as follows:  
\[ \begin{split}&\rho gdr(Rda)\cos\left(\alpha+\frac{da}{2}\right)- Pdr\cos\left(\frac{da}{2}\right)+(P+dP)dr\cos\left(\frac{da}{2}\right)-\tau\left(R- \frac{dr}{2}\right)da\\ &+(\tau+d\tau)\cdot\left(R+\frac{dr}{2}\right)da=0\end{split}