## 2 General Framework of Geomechanical Modeling  
The governing equations for heat and fluid flows, and mechanics are equations for mass and energy conservation, and momentum balance. The mass and energy conservation equations are of a similar form [22]:  
\[\frac{dm_{j}}{dt}+\
abla\cdot\omega_{j}=q_{j} \tag{1}\]  
\[\frac{dm_{\theta}}{dt}+\
abla\cdot f_{\theta}=q_{\theta} \tag{2}\]  
In the mass-conservation Equation (1), the subscript \(j\) denotes a particular fluid phase, \(m_{j}\) (kg/m\({}^{3}\)) is the unit fluid mass, \(w_{j}\) (kg/m\({}^{2}\)/s) is the unit mass-flux, and \(q_{j}\) (kg/m\({}^{3}\)/s) is a unit source term of fluid phase \(j\). In the energy-conservation Equation (2), the subscript \(\theta\) denotes the heat component, \(m_{\theta}\) (/m\({}^{3}\)/s), \(f_{0}\) (/m\({}^{3}\)/s), and \(q_{\theta}\) (/m\({}^{3}\)/s) are unit heat, heat flux, and heat source, respectively. The volumetric flux of fluid phase \(j\), \(v_{j}(m/s)=\omega_{j}/\rho_{j}\), is assumed to follow the multiphase [PERSON]'s law:  
\[\mathbf{v}_{j}=-k\frac{k_{rj}}{\mu_{j}}\big{(}\
abla p_{j}-\rho_{j}\mathbf{g}\big{)} \tag{3}\]  
where \(g\) is the vector of gravitational acceleration, \(k\) (m\({}^{2}\)) is the absolute permeability, \(k_{i}\), \(\mu_{i}\) (Pa\(\cdot\)s), and \(\rho_{j}\) (kg/m\({}^{3}\)) are the relative permeability, viscosity, and density of fluid phase \(j\), respectively. The heat and heat flux are given by  
\[m_{\theta}=-(1-\phi)\rho_{R}C_{R}T+\phi S_{j}\rho_{j}h_{j} \tag{4}\]  
\[f_{\theta}=-\left((1-\phi)K_{R}+\phi S_{j}K_{j}\right)\
abla T+h_{j}\omega_{j} \tag{5}\]  
where \(\phi\), \(\rho_{R}\), \(C_{R}\) (J/kg/K), and \(K_{R}\) (J/m\({}^{2}\)/s/K) are the porosity, density, heat capacity, and thermal conductivity of the porous media, respectively. \(S_{j}\), and \(h_{j}\) (J/kg), are the saturation and specific enthalpy for fluid phase \(j\), respectively. The linear momentum balance equation is given by  
\[\text{Div}\,\sigma+\rho_{b}g=0 \tag{6}\]  
\[\rho_{b}=(1-\phi)\rho_{R}+\phi S_{j}\rho_{j} \tag{7}\]  
where, \(\sigma\) is the Cauchy total-stress tensor, and \(\rho_{b}\) is the bulk density. Based on linear poroelasticity theory [23], the coupled governing equations for pore pressure (here, the pore is saturated with a single \"effective fluid\") and solid mechanics can be written as [24]  
\[(\lambda+2\mu)\
abla(\
abla\cdot u)-\mu\
abla\times(\
abla\times u)-\alpha \
abla p=f \tag{8}\]  
\[Q^{-1}\frac{\partial p}{\partial t}+\alpha\frac{\partial}{\partial t}\
abla \cdot u-D\
abla^{2}p=q \tag{9}\]  
where \(u\) (m) and \(p\) (Pa) are the displacement vector and excess pore pressure, respectively. \(\lambda\) and \(\mu\) are Lame coefficients, \(\alpha\) is the dimensionless coefficient representing the change in pore pressure per unit change in bulk volume under drained conditions, \(Q^{-1}\) (Pa\({}^{-1}\)) is the bulk compressibility, \(D\) is Darcy conductivity. The \(f\) (Nf/m\({}^{2}\)) and \(q\) (m\({}^{3}\)/s) are body force and fluid source, respectively. By assuming the infinitesimal transformation, the strain tensor (\(\varepsilon\)) is given by the symmetric gradient of the displacement vector (\(u\)):  
\[\varepsilon=\frac{1}{2}(\
abla u+\
abla^{T}u) \tag{10}\]  
As seen from these equations, fluid flow, heat transfer, and the mechanic responses of the reservoir systems are coupled with each other. Its solution usually requires a coupled modeling approach, thus, coupled THM analysis is important in geomechanical modeling [25]. Figure 1 shows a schematic of the geomechanical modeling framework, which contains a coupled THM simulation and history matching for calibration of the numerical model. Although full coupling schemes are more accurate, they are also more complex. On the other hand, sequential coupling schemes, in which fluid and heat flows and mechanics are explicitly coupled, may be sufficiently accurate for problems in which coupling is relatively weak.