## 2 Materials and Methods
### Thermo-Mechanical Coupling Governing Equations
#### 2.1.2 Governing Equations of the Stress-Strain Field  
To better describe the orthotropic anisotropic deformation of frozen soil, a local coordinate system was established, with directions 1, 2, and 3 representing the three principal axes. In the local coordinate system, the relationship between the thermal strain component caused by temperature and the freeze-thaw volume deformation is as follows [29]:  
\[\begin{bmatrix}\varepsilon_{11}^{T}&\varepsilon_{12}^{T}&\varepsilon_{13}^{T} \\ \varepsilon_{21}^{T}&\varepsilon_{22}^{T}&\varepsilon_{23}^{T}\\ \varepsilon_{31}^{T}&\varepsilon_{32}^{T}&\varepsilon_{33}^{T}\end{bmatrix}= \begin{bmatrix}\eta&0&0\\ 0&(1-\eta)/2&0\\ 0&0&(1-\eta)/2\end{bmatrix}\varepsilon_{\
u} \tag{4}\]  
Assuming that frozen soil is a linear elastic body, when the soil reaches the freezing or thawing temperature, the total strain of frozen soil is composed of two parts--the elastic strain caused by the load and the volume deformation caused by the freezing-thawing--which can be expressed as  
\[\left\{\begin{array}{l}\varepsilon_{11}=\frac{1}{E}[\sigma_{11}-\mu(\sigma_{22}+ \sigma_{33})]+\eta\varepsilon_{v}\\ \varepsilon_{22}=\frac{1}{E}[\sigma_{22}-\mu(\sigma_{11}+\sigma_{33})]+\frac{1 }{2}(1-\eta)\varepsilon_{v}\\ \varepsilon_{33}=\frac{1}{E}[\sigma_{33}-\mu(\sigma_{11}+\sigma_{22})]+\frac{1 }{2}(1-\eta)\varepsilon_{v}\\ \varepsilon_{12}=\frac{2}{E}(1+\mu)\tau_{12}\\ \varepsilon_{23}=\frac{2}{E}(1+\mu)\tau_{23}\\ \varepsilon_{31}=\frac{2}{E}(1+\mu)\tau_{31}\end{array}\right. \tag{5}\]  
where \(E\) is the elastic modulus of the frozen soil, \(\mu\) is the Poisson's ratio of the frozen soil, and \(\eta\) is the deformation characteristic coefficient, which is a dimensionless quantity. When \(\eta=1/3\), the deformation characteristic of frozen soil is isotropic; \(\varepsilon_{v}\) is the volume deformation of soil caused by frost heave or thawing settlement due to temperature changes.  
For the plane strain problem, given by \(\varepsilon_{33}=\varepsilon_{23}=\varepsilon_{31}=0\), the following equation can be used:  
\[\left\{\begin{array}{l}\varepsilon_{11}=\frac{1-\mu^{2}}{E}\Big{(}\sigma_{11 }-\frac{\mu}{1-\mu}\sigma_{22}\Big{)}+\big{[}\eta+\frac{\mu}{2}(1-\eta)\big{]} \varepsilon_{v}\\ \varepsilon_{22}=\frac{1-\mu^{2}}{E}\Big{(}\sigma_{22}-\frac{\mu}{1-\mu} \sigma_{11}\Big{)}+\frac{1}{2}(1-\eta)(1+\mu)\varepsilon_{v}\\ \varepsilon_{12}=\frac{2}{E}(1+\mu)\tau_{12}\end{array}\right. \tag{6}\]  
Suppose that the angle between heat flow direction 1 and the \(x\)-axis of the global coordinate system is \(\theta\); then, the thermal strain component caused by the temperature in the global coordinate system is  
\[\left\{\begin{array}{l}\varepsilon_{x}^{T}=\Big{\{}\big{[}\eta+\frac{\mu}{2 }(1-\eta)\big{]}\cos^{2}\theta+\Big{[}\frac{1}{2}(1-\eta)(1+\mu)\Big{]}\sin^{2 }\theta\Big{\}}\varepsilon_{v}\\ \varepsilon_{y}^{T}=\Big{\{}\big{[}\eta+\frac{\mu}{2}(1-\eta)\big{]}\sin^{2} \theta+\Big{[}\frac{1}{2}(1-\eta)(1+\mu)\Big{]}\cos^{2}\theta\Big{\}} \varepsilon_{v}\\ \varepsilon_{xy}^{T}=\frac{1}{2}(3\eta-1)\sin\theta\cos\theta\varepsilon_{v} \end{array}\right. \tag{7}\]  
In summary, the governing equation for the thermo-mechanical coupling stress-strain field of the 2D tunnel cross section is as follows:  
\[\left\{\begin{array}{l}\varepsilon_{x}=\frac{1-\mu(T)^{2}}{E(T)}\Big{(} \sigma_{x}-\frac{\mu(T)}{1-\mu(T)}\sigma_{y}\Big{)}+\varepsilon_{x}^{T}\\ \varepsilon_{y}=\frac{1-\mu(T)^{2}}{E(T)}\Big{(}\sigma_{y}-\frac{\mu(T)}{1-\mu (T)}\sigma_{x}\Big{)}+\varepsilon_{y}^{T}\\ \varepsilon_{xy}=\frac{2}{E(T)}[1+\mu(T)]\tau_{xy}+\varepsilon_{xy}^{T}\end{array}\right. \tag{8}\]