## 2 Methodology
### Improved Mass-Conserved and Interface-Corrected DI Method
#### 2.1.2 Improved Cahn-Hilliard Model with Mass Conservation and Interface Correction  
This section presents a development of the improved Cahn-Hilliard model to effectively capture the phase interface. In order to maintain the equilibrium shape, as well as the thickness of the fluid interface, and guarantee mass conservation, an interface correction term \(q_{1}\) and a mass conservation term \(q_{2}\) are adopted in the C-H model. The modified C-H model with mass conservation and interface correction can be written as:  
\[\frac{\partial C}{\partial t}+\
abla\cdot(\mathbf{u}C)=\Gamma\
abla^{2}\mu_{C}-q_ {1}+q_{2}, \tag{5}\]  
where \(q_{1}\) is introduced to maintain the shape of the phase interface and \(q_{2}\) is used to guarantee mass conservation. The interfacial shape can be obtained analytically in the one-dimensional case by using the free energy density model for immiscible two-phase flows:  
\[C(d)=\frac{1}{2}\Big{[}1+\tanh\left(\frac{2d}{\xi}\right)\Big{]}, \tag{6}\]  
where \(d\) is the signed distance along the gradient of \(C\). It ought to be noticed that the surface tension force is obtained by assuming that the interface shape follows Equation (6) during the evolution of the fluid interfaces. This means that, when the shape of the interface is different from that described by Equation (6), the surface tension cannot be obtained accurately. However, the DI model is not able to guarantee the local shape of the interface in unsteady flow simulations, especially if the density ratio of fluids is huge and the convective effect is strong. It should be noted that the level set approach can quite easily maintain the interfacial shape, since the signed distance function \(d\) is obtained accurately by applying the Hamilton-Jacobian equation, which can be written as:  
\[\frac{\partial d}{\partial\tau}+S(d)(|\
abla d|-1)=0, \tag{7}\]  
\[S(d)=\frac{d}{\sqrt{(d)^{2}+(|\
abla d|\varepsilon)^{2}}}, \tag{8}\]  
where \(\tau\) is the pseudo time, \(S(d)\) is a smoothed signed function, and \(\varepsilon\) is equal to the mesh size \(\Delta x\) in this paper. Inspired by the level set approach, where the Hamilton-Jacobian equation is applied to accurately calculate the signed distance equation, the interface correction term \(q_{1}\) is introduced into the Cahn-Hilliard model, which can be written as:  
\[q_{1}=kS(d)(|\
abla d|-1), \tag{9}\]  
Here, \(k\) is a constant connected with the mesh size \(\Delta x\) and the characteristic length \(L\). It ought to be noticed that \(q_{1}\) only works in the interfacial layer, while it is zero elsewhere. To get the expression of the source term \(q_{1}\) explicitly in terms of the order parameter, the signed distance function \(d\) can be expressed by using the order parameter with the aid of Equation (6):  
\[d=\frac{\varepsilon}{2}\tanh^{-1}[2(C-0.5)], \tag{10}\]  
Equation (10) implies that the signed \(d=0\) and the interface of the order parameter \(C=0.5\) both represent the phase interface. In our research, in order to avoid a numerical error, the interface correction term \(q_{1}\) is applied in the zone of \(0.01\leq C\leq 0.99\). Substituting Equation (10) into Equation (9) gives the final expression of the source term \(q_{1}\). \[q_{1}= \begin{cases}\frac{k\ln\left(\frac{C}{1-C}\right)}{\sqrt{\left(\ln\left( \frac{C}{1-C}\right)\right)^{2}+}\left(\left|\frac{1}{C-C}\right|\left|\
abla C \right|\right)^{2}}\left(\left|\frac{\xi}{4C-4C^{2}}\right|\left|\
abla C \right|-1\right),&\text{if}(0.01\leq C\leq 0.99)\\ 0,&\text{elsewhere}\end{cases}  \tag{11}\]  
In addition to the interface shape correction term \(q_{1}\), the mass conservation term \(q_{2}\) is obtained by assuming that this term exists averagely in the interfacial layer. It can thus be obtained by integrating Equation (5) in the zone \(\Omega_{L}\) filled with lighter fluid.