The computation of strain rate tensor in multiple-relaxation-time lattice Boltzmann model

TL;DR

This paper compares two formulas for computing the strain rate tensor in MRT-LB models, analyzing their relationship, advantages, and applicability to different lattice patterns, to guide better choice in complex flow simulations.

Contribution

It clarifies the relationship between Yu and Chai formulas and discusses their respective advantages, aiding in selecting the appropriate method for MRT-LB simulations.

Findings

Yu formula is more computationally efficient.

Chai formula applies to more lattice patterns.

Yu formula can be derived from Chai formula in specific cases.

Abstract

Multiple-relaxation-time (MRT) lattice Boltzmann (LB) model is an important class of LB model with lots of advantages over traditional single-relaxation-time (SRT) LB model. In addition, the computation of strain rate tensor is crucial in MRT-LB simulations of some complex flows. Up to now, there are only two formulas to compute the strain rate tensor in the MRT LB model. One is to compute the strain rate tensor by using non-equilibrium parts of macroscopic moments (Yu formula). The other is to compute the strain rate tensor by using non-equilibrium parts of density distribution functions (Chai formula). The mathematical expressions of these two formulas are so different that we do not know which formula to choose for computing the strain rate tensor in MRT LB model. In this paper, we study the relationship of these two formulas. It is found that Yu formula can be deduced from Chai formula in a particular procedure. However, these two formulas have their own advantages and disadvantages. Yu formula is more efficient in the computation aspect while Chai formula can be applied to more lattice patterns of MRT LB models. It is also found that deducing Yu formula for a particular lattice pattern from Chai formula is more convenient than the way proposed by Yu et al.

Equations

$$f_{\alpha }(\mathbf{x}+{\mathbf{c}}_{\alpha }{\delta }_t,t+{\delta }_t)-f_{\alpha }(\mathbf{x},t)=-\sum _{i=0}^{18}{\Lambda }_{\alpha i}(f_i(\mathbf{x},t)-f_i^{(eq)}(\mathbf{x},t)),\alpha =0-18,$$

$$f_i^{(eq)}(\mathbf{x},t)={\omega }_i\left[\rho +{\rho }_0\left(\frac{{\mathbf{c}}_i\cdot \mathbf{u}}{c_s^2}+\frac{{({\mathbf{c}}_i\cdot \mathbf{u})}^2}{2c_s^4}-\frac{|\mathbf{u}{|}^2}{2c_s^2}\right)\right],$$

$${\omega_{i}}=\left\{{\begin{array}[]{*{20}{l}}{1/3,\;\;\;\;\;\;i=0,}\\ {1/18,\;\;\;\;\;\;i=1-6,}\\ {1/36,\;\;\;\;i=7-18,}\end{array}}\right.$$

$$\footnotesize\begin{array}[]{c}\{\emph{{c}}_{0},\emph{{c}}_{1},\ldots,\emph{{c}}_{18}\}=\\ \left\{\begin{array}[]{ccccccccccccccccccc}0&1&-1&0&0&0&0&1&-1&1&-1&1&-1&1&-1&0&0&0&0\\ 0&0&0&1&-1&0&0&1&1&-1&-1&0&0&0&0&1&-1&1&-1\\ 0&0&0&0&0&1&-1&0&0&0&0&1&1&-1&-1&1&1&-1&-1\end{array}\right\}c,\end{array}$$

$$f_i^{(eq)}(\mathbf{x},t)={\omega }_i\left[\delta \rho +{\rho }_0\left(\frac{{\mathbf{c}}_i\cdot \mathbf{u}}{c_s^2}+\frac{{({\mathbf{c}}_i\cdot \mathbf{u})}^2}{2c_s^4}-\frac{|\mathbf{u}{|}^2}{2c_s^2}\right)\right].$$

$$|f(\mathbf{x}+{\mathbf{c}}_{\alpha }{\delta }_t,t+{\delta }_t)⟩-|f(\mathbf{x},t)⟩=-\mathbf{\Lambda }(|f(\mathbf{x},t)⟩-|f^{(eq)}(\mathbf{x},t))⟩,$$

$$\mathbf{T}=\left({\begin{array}[]{*{20}c}1&1&1&1&1&1&1&1&1&1&1&1&1&1&1&1&1&1&1\\ -30&{-11}&{-11}&{-11}&{-11}&{-11}&{-11}&8&8&8&8&8&8&8&8&8&8&8&8\\ 12&{-4}&{-4}&{-4}&{-4}&{-4}&{-4}&1&1&1&1&1&1&1&1&1&1&1&1\\ 0&1&{-1}&0&0&0&0&1&{-1}&1&{-1}&1&{-1}&1&{-1}&0&0&0&0\\ 0&{-4}&4&0&0&0&0&1&{-1}&1&{-1}&1&{-1}&1&{-1}&0&0&0&0\\ 0&0&0&1&{-1}&0&0&1&1&{-1}&{-1}&0&0&0&0&1&{-1}&1&{-1}\\ 0&0&0&{-4}&4&0&0&1&1&{-1}&{-1}&0&0&0&0&1&{-1}&1&{-1}\\ 0&0&0&0&0&1&{-1}&0&0&0&0&1&1&{-1}&{-1}&1&1&{-1}&{-1}\\ 0&0&0&0&0&{-4}&4&0&0&0&0&1&1&{-1}&{-1}&1&1&{-1}&{-1}\\ 0&2&2&{-1}&{-1}&{-1}&{-1}&1&1&1&1&1&1&1&1&{-2}&{-2}&{-2}&{-2}\\ 0&{-4}&{-4}&2&2&2&2&1&1&1&1&1&1&1&1&{-2}&{-2}&{-2}&{-2}\\ 0&0&0&1&1&{-1}&{-1}&1&1&1&1&{-1}&{-1}&{-1}&{-1}&0&0&0&0\\ 0&0&0&{-2}&{-2}&2&2&1&1&1&1&{-1}&{-1}&{-1}&{-1}&0&0&0&0\\ 0&0&0&0&0&0&0&1&{-1}&{-1}&1&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&{-1}&{-1}&1\\ 0&0&0&0&0&0&0&0&0&0&0&1&{-1}&{-1}&1&0&0&0&0\\ 0&0&0&0&0&0&0&1&{-1}&1&{-1}&{-1}&1&{-1}&1&0&0&0&0\\ 0&0&0&0&0&0&0&{-1}&{-1}&1&1&0&0&0&0&1&{-1}&1&{-1}\\ 0&0&0&0&0&0&0&0&0&0&0&1&1&{-1}&{-1}&{-1}&{-1}&1&1\\ \end{array}}\right),$$

$$|f(\mathbf{x}+{\mathbf{c}}_{\alpha }{\delta }_t,t+{\delta }_t)⟩-|f(\mathbf{x},t)⟩=-{\mathbf{T}}^{-1}\hat{\mathbf{\Lambda }}(\mathbf{m}(\mathbf{x},t)-{\mathbf{m}}^{(eq)}(\mathbf{x},t)),$$

$$\mathbf{m}=(\delta \rho ,e,ϵ,{\rho }_0{u}_x,q_x,{\rho }_0{u}_y,q_y,{\rho }_0{u}_z,q_z,3p_{xx},3{\pi }_{xx},p_{\omega \omega },{\pi }_{\omega \omega },p_{xy},p_{yz},p_{xz},t_x,t_y,t_z{)}^{\prime },$$

$${\mathbf{m}}^{(eq)}=(\delta \rho ,e^{(eq)},{ϵ}^{(eq)},{\rho }_0{u}_x,q_x^{(eq)},{\rho }_0{u}_y,q_y^{(eq)},{\rho }_0{u}_z,q_z^{(eq)},3p_{xx}^{(eq)},3{\pi }_{xx}^{(eq)},p_{\omega \omega }^{(eq)},{\pi }_{\omega \omega }^{(eq)},p_{xy}^{(eq)},p_{yz}^{(eq)},p_{xz}^{(eq)},t_x^{(eq)},t_y^{(eq)},t_z^{(eq)}{)}^{\prime }.$$

$$\left.\begin{array}[]{c}e^{(eq)}=-11\delta\rho+\frac{19}{\rho_{0}}\bm{j}\cdot\bm{j},\\ \epsilon^{(eq)}=\omega_{\epsilon}\delta\rho+\frac{\omega_{\epsilon j}}{\rho_{0}}\bm{j}\cdot\bm{j},\end{array}\right\}$$

$$\left.\begin{array}[]{c}q_{x}^{(eq)}=-2j_{x}/3,\\ q_{y}^{(eq)}=-2j_{y}/3,\\ q_{z}^{(eq)}=-2j_{z}/3,\\ \end{array}\right\}$$

$$\left.\begin{array}[]{c}3p_{xx}^{(eq)}=\frac{1}{\rho_{0}}(3j_{x}^{2}-\bm{j}\cdot\bm{j}),\\ p_{\omega\omega}^{(eq)}=3\omega_{xx}p_{xx}^{(eq)},\\ \end{array}\right\}$$

$$\left.\begin{array}[]{c}3\pi_{xx}^{(eq)}=\frac{1}{\rho_{0}}(j_{y}^{2}-j_{z}^{2}),\\ \pi_{\omega\omega}^{(eq)}=3\omega_{xx}\pi_{xx}^{(eq)},\\ \end{array}\right\}$$

$$\left.\begin{array}[]{c}p_{xy}^{(eq)}=\frac{1}{\rho_{0}}j_{x}j_{y},\\ p_{yz}^{(eq)}=\frac{1}{\rho_{0}}j_{y}j_{z},\\ p_{xz}^{(eq)}=\frac{1}{\rho_{0}}j_{x}j_{z},\\ \end{array}\right\}$$

$$\left.\begin{array}[]{c}t_{x}^{(eq)}=0,\\ t_{y}^{(eq)}=0,\\ t_{z}^{(eq)}=0,\\ \end{array}\right\}$$

$$\hat{\mathbf{\Lambda }}\equiv diag(s_0,s_1,s_2,s_3,s_4,s_5,s_6,s_7,s_8,s_9,s_{10},s_{11},s_{12},s_{13},s_{14},s_{15},s_{16},s_{17},s_{18}),$$

$$\hat{\mathbf{\Lambda }}\equiv diag(s_c,s_e,s_{ϵ},s_c,s_q,s_c,s_q,s_c,s_q,s_{\nu },s_{\pi },s_{\nu },s_{\pi },s_{\nu },s_{\nu },s_{\nu },s_t,s_t,s_t),$$

$$S_{xx}\approx -\frac{\epsilon }{38{\rho }_0{\delta }_t}(s_1{m}_1^{(1)}+19s_9{m}_9^{(1)}),$$

$$S_{yy,zz}\approx -\frac{\epsilon }{76{\rho }_0{\delta }_t}\left[2s_1{m}_1^{(1)}-19s_9(m_9^{(1)}\mp 3m_{11}^{(1)})\right],$$

$$S_{xy,yz,zx}\approx -\frac{3\epsilon s_9}{2{\rho }_0{\delta }_t}{m}_{13,14,15}^{(1)},$$

$$\mathbf{S}=-\frac{{\sum }_i{\mathbf{c}}_i{\mathbf{c}}_i{\left({\mathbf{T}}^{-1}\hat{\mathbf{\Lambda }}\mathbf{T}\right)}_{ij}{f}_j^{(neq)}}{2{\rho }_0{c}_s^2{\delta }_t},$$

$$\mathbf{S}=-\frac{{\sum }_i{\mathbf{c}}_i{\mathbf{c}}_i\left(f_i-f_i^{(eq)}\right)}{2{\rho }_0{c}_s^2{\delta }_t\tau },$$

$$f_i=f_i^{(0)}+\epsilon f_i^{(1)}+{\epsilon }^2{f}_i^{(2)}+\cdots$$

$${\partial }_t=\epsilon {\partial }_{t_1}+{\epsilon }^2{\partial }_{t_2}$$

$${\partial }_{\alpha }=\epsilon {\partial }_{\alpha 1}.$$

$${\epsilon }^0:f_i^{(0)}=f_i^{(eq)},$$

$${\epsilon }^1:D_{1i}{f}_i^{(0)}=-\frac{1}{{\delta }_t}{\Lambda }_{ij}{f}_j^{(1)},$$

$${\epsilon }^2:{\partial }_{t_2}{f}_i^{(0)}+D_{1i}{f}_i^{(1)}+\frac{{\delta }_t}{2}{D}_{1i}^2{f}_i^{(0)}=-\frac{1}{{\delta }_t}{\Lambda }_{ij}{f}_j^{(2)},$$

$${\epsilon }^0:{\mathbf{m}}^{(0)}={\mathbf{m}}^{(eq)},$$

$${\epsilon }^1:{\tilde{\mathbf{D}}}_1{\mathbf{m}}^{(0)}=-{\hat{\mathbf{\Lambda }}}^{\prime }{\mathbf{m}}^{(1)},$$

$${\epsilon }^2:{\partial }_{t_2}{\mathbf{m}}^{(0)}+{\tilde{\mathbf{D}}}_1(\mathbf{I}-\frac{{\hat{\mathbf{\Lambda }}}^{\prime }}{2}){\mathbf{m}}^{(1)}=-{\hat{\mathbf{\Lambda }}}^{\prime }{\mathbf{m}}^{(2)},$$

$$\mathbf{T}{\mathbf{D}}_1{\mathbf{T}}^{-1}=\mathbf{T}({\partial }_{t_1}\mathbf{I}+{\mathbf{C}}_x{\partial }_{x1}+{\mathbf{C}}_y{\partial }_{y1}+{\mathbf{C}}_z{\partial }_{z1}){\mathbf{T}}^{-1}={\partial }_{t_1}\mathbf{I}+{\hat{\mathbf{C}}}_k{\partial }_{k1},$$

$$\mathbf{\hat{C}}_{x}=\left({\begin{array}[]{*{20}c}0&0&0&1&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&\frac{21}{5}&\frac{19}{5}&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ \frac{10}{19}&\frac{1}{57}&0&0&0&0&0&0&0&\frac{1}{3}&0&0&0&0&0&0&0&0&0\\ 0&\frac{4}{63}&\frac{10}{63}&0&0&0&0&0&0&-\frac{2}{9}&\frac{5}{9}&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\ 0&0&0&\frac{6}{5}&-\frac{1}{5}&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0\\ 0&0&0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&-\frac{1}{2}&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&\frac{1}{2}\\ 0&0&0&0&0&0&0&0&0&0&0&\frac{2}{3}&\frac{1}{3}&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&-1&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\ \end{array}}\right).$$

$$\mathbf{\hat{C}}_{y}=\left({\begin{array}[]{*{20}c}0&0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&\frac{21}{5}&\frac{19}{5}&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0&0\\ \frac{10}{19}&\frac{1}{57}&0&0&0&0&0&0&0&-\frac{1}{6}&0&\frac{1}{2}&0&0&0&0&0&0&0\\ 0&\frac{4}{63}&\frac{10}{63}&0&0&0&0&0&0&\frac{1}{9}&-\frac{5}{18}&-\frac{1}{3}&\frac{5}{6}&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0\\ 0&0&0&0&0&-\frac{3}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&-\frac{3}{2}&0\\ 0&0&0&0&0&0&-\frac{1}{2}&0&0&0&0&0&0&0&0&0&0&-\frac{3}{2}&0\\ 0&0&0&0&0&\frac{3}{5}&-\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&-\frac{1}{2}&0\\ 0&0&0&0&0&0&\frac{1}{2}&0&0&0&0&0&0&0&0&0&0&-\frac{1}{2}&0\\ 0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&0&\frac{1}{2}&0&0\\ 0&0&0&0&0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&-\frac{1}{2}\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&-\frac{1}{3}&-\frac{1}{6}&-\frac{1}{3}&-\frac{1}{6}&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&-1&0&0&0&0\\ \end{array}}\right).$$

$$\mathbf{\hat{C}}_{z}=\left({\begin{array}[]{*{20}c}0&0&0&0&0&0&0&1&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&\frac{21}{5}&\frac{19}{5}&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&1&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0\\ \frac{10}{19}&\frac{1}{57}&0&0&0&0&0&0&0&-\frac{1}{6}&0&-\frac{1}{2}&0&0&0&0&0&0&0\\ 0&\frac{4}{63}&\frac{10}{63}&0&0&0&0&0&0&\frac{1}{9}&-\frac{5}{18}&\frac{1}{3}&-\frac{5}{6}&0&0&0&0&0&0\\ 0&0&0&0&0&0&0&-\frac{3}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&\frac{3}{2}\\ 0&0&0&0&0&0&0&0&-\frac{1}{2}&0&0&0&0&0&0&0&0&0&\frac{3}{2}\\ 0&0&0&0&0&0&0&-\frac{3}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&-\frac{1}{2}\\ 0&0&0&0&0&0&0&0&-\frac{1}{2}&0&0&0&0&0&0&0&0&0&-\frac{1}{2}\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0&0\\ 0&0&0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&\frac{1}{2}&0\\ 0&0&0&\frac{2}{5}&\frac{1}{10}&0&0&0&0&0&0&0&0&0&0&0&-\frac{1}{2}&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&-1&0&0&0\\ 0&0&0&0&0&0&0&0&0&0&0&0&0&0&1&0&0&0&0\\ 0&0&0&0&0&0&0&0&0&\frac{1}{3}&\frac{1}{6}&-\frac{1}{3}&-\frac{1}{6}&0&0&0&0&0&0\\ \end{array}}\right).$$

$$\tiny{\partial_{{t_{1}}}}\left[{\begin{array}[]{*{20}{c}}{\delta\rho}\\ {-11\delta\rho+\frac{19\textbf{\emph{j}}\cdot\textbf{\emph{j}}}{\rho_{0}}}\\ {\omega_{\epsilon}\delta\rho+\frac{\omega_{\epsilon j}\textbf{\emph{j}}\cdot\textbf{\emph{j}}}{\rho_{0}}}\\ {j_{x}}\\ {-2j_{x}/3}\\ {j_{y}}\\ {-2j_{y}/3}\\ {j_{z}}\\ {-2j_{z}/3}\\ {(2j_{x}^{2}-j_{y}^{2}-j_{z}^{2})/\rho_{0}}\\ {\omega_{xx}(2j_{x}^{2}-j_{y}^{2}-j_{z}^{2})/\rho_{0}}\\ {(j_{y}^{2}-j_{z}^{2})/\rho_{0}}\\ {\omega_{xx}(j_{y}^{2}-j_{z}^{2})/\rho_{0}}\\ {j_{x}j_{y}/\rho_{0}}\\ {j_{y}j_{z}/\rho_{0}}\\ {j_{x}j_{z}/\rho_{0}}\\ {0}\\ {0}\\ {0}\end{array}}\right]+{\partial_{x1}}\left[{\begin{array}[]{*{20}{c}}{{j_{x}}}\\ {5{j_{x}}/3}\\ {-2{j_{x}}/3}\\ {j_{x}^{2}/\rho_{0}+\delta\rho/3}\\ A\\ {{j_{x}}{j_{y}}/\rho_{0}}\\ {{j_{x}}{j_{y}}/\rho_{0}}\\ {{j_{x}}{j_{z}}/\rho_{0}}\\ {{j_{x}}{j_{z}}/\rho_{0}}\\ {4j_{x}/3}\\ {-2j_{x}/3}\\ 0\\ {0}\\ {j_{y}/3}\\ {0}\\ {j_{z}/3}\\ {\frac{(j_{y}^{2}-j_{z}^{2})(2+\omega_{xx})}{3\rho_{0}}}\\ {-j_{x}j_{y}/\rho_{0}}\\ {j_{x}j_{z}/\rho_{0}}\end{array}}\right]+{\partial_{y1}}\left[{\begin{array}[]{*{20}{c}}{{j_{y}}}\\ {5{j_{y}}/3}\\ {-2{j_{y}}/3}\\ {{j_{x}}{j_{y}}/\rho_{0}}\\ {{j_{x}}{j_{y}}/\rho_{0}}\\ {j_{y}^{2}/\rho_{0}+\delta\rho/3}\\ {B}\\ {j_{y}j_{z}/\rho_{0}}\\ {j_{y}j_{z}/\rho_{0}}\\ {-2j_{y}/3}\\ {j_{y}/3}\\ {2j_{y}/3}\\ {-j_{y}/3}\\ {j_{x}/3}\\ {j_{z}/3}\\ {0}\\ {j_{x}j_{y}/\rho_{0}}\\ {\frac{(j_{z}^{2}-j_{x}^{2})(2+\omega_{xx})}{3\rho_{0}}}\\ {-j_{y}j_{z}/\rho_{0}}\end{array}}\right]+{\partial_{z1}}\left[{\begin{array}[]{*{20}{c}}{{j_{z}}}\\ {5{j_{z}}/3}\\ {-2{j_{z}}/3}\\ {{j_{x}}{j_{z}}/\rho_{0}}\\ {{j_{x}}{j_{z}}/\rho_{0}}\\ {{j_{y}}{j_{z}}/\rho_{0}}\\ {{j_{y}}{j_{z}}/\rho_{0}}\\ {j_{z}^{2}/\rho_{0}+\delta\rho/3}\\ {C}\\ {-2{j_{z}}/3}\\ {j_{z}/3}\\ {-2{j_{z}}/3}\\ {j_{z}/3}\\ {0}\\ {j_{y}/3}\\ {j_{x}/3}\\ {-j_{x}j_{z}/\rho_{0}}\\ {j_{y}j_{z}/\rho_{0}}\\ {\frac{(j_{x}^{2}-j_{y}^{2})(2+\omega_{xx})}{3\rho_{0}}}\end{array}}\right]=-\left[{\begin{array}[]{*{20}{c}}0\\ {{s_{1}^{{}^{\prime}}}{m^{(1)}_{1}}}\\ {{s_{2}^{{}^{\prime}}}{m^{(1)}_{2}}}\\ {0}\\ {{s_{4}^{{}^{\prime}}}{m^{(1)}_{4}}}\\ {0}\\ {{s_{6}^{{}^{\prime}}}{m^{(1)}_{6}}}\\ {0}\\ {{s_{8}^{{}^{\prime}}}{m^{(1)}_{8}}}\\ {{s_{9}^{{}^{\prime}}}{m^{(1)}_{9}}}\\ {{s_{10}^{{}^{\prime}}}{m^{(1)}_{10}}}\\ {{s_{11}^{{}^{\prime}}}{m^{(1)}_{11}}}\\ {{s_{12}^{{}^{\prime}}}{m^{(1)}_{12}}}\\ {{s_{13}^{{}^{\prime}}}{m^{(1)}_{13}}}\\ {{s_{14}^{{}^{\prime}}}{m^{(1)}_{14}}}\\ {{s_{15}^{{}^{\prime}}}{m^{(1)}_{15}}}\\ {{s_{16}^{{}^{\prime}}}{m^{(1)}_{16}}}\\ {{s_{17}^{{}^{\prime}}}{m^{(1)}_{17}}}\\ {{s_{18}^{{}^{\prime}}}{m^{(1)}_{18}}}\\ \end{array}}\right]$$

Full text

††journal:

1 Introduction

The lattice Boltzmann equation (LBE), as a mesoscopic numerical method, has been widely used to simulate various complex fluid flows and has gained significant success because of its distinct advantages such as the natural parallelism of algorithm, simplicity of programming and ease of dealing with complex boundary conditions [1, 2, 3].

Recently, the LBE method has also been used to the studies of blood flows[4, 5, 6, 7], non-Newtonian fluid flows[8, 9, 10, 11, 12, 13, 14, 15], large eddy simulation of turbulent flows[16, 17, 18], multiphase flows[19, 20] and so on. In all the above studies, the computation of strain rate tensor is of key importance for simulations, therefore the computation of strain rate tensor has received increasing attention over the past decade [16, 12, 21, 22].

For SRT LB model, it compute the strain rate tensor by using non-equilibrium parts of density distribution functions [4, 6, 19, 21]. This particular computation method is local and has second-order accuracy in space, so it is very suitable for the studies of flows in complex geometries and has great advantage over the finite difference method. Although the computation of strain rate tensor for SRT LB (also called LBGK) model has been studied a lot, the computation of strain rate tensor for MRT LB model is little investigated.

Yu et al. proposed to compute the strain rate tensor with non-equilibrium parts of macroscopic moments (Yu formula) for large eddy simulation (LES) of turbulent jets by using MRT LB model with nineteen discrete velocities in three dimensions (D3Q19) [16]. Based on Yu formula, Premnath et al. developed the formula to compute the strain rate tensor with the external force effect for LES of turbulent flows [23]. Pattison et al. used the developed formula by Premnath et al. to compute the strain rate tensor for LES of turbulent flow in a straight square duct driven by a pressure gradient [24]. Some other works also used the formula developed by Yu et al. and Premnath et al. to compute the strain rate tensor for the LES of turbulent flows [25, 26].

More recently, Chai et al. proposed the other formula to compute the strain rate tensor with non-equilibrium parts of density distribution functions (Chai formula) for simulation of non-Newtonian fluid flows by using MRT LB model with nine discrete velocities in two dimensions (D2Q9) [12]. After this work, Chai formula was widely used to compute the strain rate tensor in many studies on non-Newtonian fluid flows by MRT-LB model [27, 28, 29].

To the knowledge of the authors, Yu formula and Chai formula are the only two formulas to compute the strain rate tensor in the MRT LB model up to now. The computation formulas are so different that we are not sure which formula to choose for current simulations. In this paper, we firstly analyze the derivation of these two formulas. Then we study the characteristics of these two formulas. Finally, we give our recommendation to the formula, which should be chosen for computing the strain rate tensor in the MRT LB model.

In the following, we firstly present He-Luo type MRT LB model with D3Q19 lattice [30] as the starting point. Based on the Chapman-Enskog analysis, we then deduce two formulas in computing the strain rate tensor in D3Q19 MRT LB model, i.e., Yu formula and Chai formula. Thirdly, we study the relationship of above two formulas and analyze the advantages and disadvantages of these two formulas. The suggestions to the choice of formulas when computing the strain rate tensor in MRT LB model are also given. Finally, we give our conclusions.

2 He-Luo D3Q19 MRT LB model

In this paper, we compare the only two formulas (Yu formula and Chai formula) in computing the strain rate tensor in MRT LB model. The macroscopic equilibrium moments of the MRT LB model are chosen to be derived from the equilibrium distribution functions of SRT LB (also called LBGK) model proposed by He and Luo [31], so we call this type of MRT LB model as He-Luo MRT LB model. In the following, we take He-Luo MRT LB model with D3Q19 lattice as an example to carry out the comparison.

The evolution equation of He-Luo D3Q19 MRT LB model is

[TABLE]

where $f_{l}(\bm{x},t)$ and $f_{l}^{(eq)}(\bm{x},t)$ ($l=\alpha,i$) are the distribution function and equilibrium distribution function of particles with velocity $\bm{c}_{l}$ at node $\bm{x}$ and time $t$, $\Lambda_{\alpha i}$ is the element located in $\alpha$ row and $i$ column of $19\times 19$ collision matrix $\bm{\Lambda}$. For He-Luo MRT model, the equilibrium distribution function is chosen as

[TABLE]

where

[TABLE]

$\bm{c}_{i}$ is defined as

[TABLE]

where $c=\delta_{x}/\delta_{t}$ is the particle velocity and $\delta_{x}$ and $\delta_{t}$ are the lattice spacing and time step, respectively. $c_{s}=c/\sqrt{3}$ is the sound speed. $\rho$ and $\bm{u}$ are the density and velocity of fluid. $\rho_{0}$ is the mean density. The density fluctuation $\delta\rho=\rho-\rho_{0}$ is usually used instead of $\rho$ in the equilibrium distribution function to reduce the numerical effects due to the round-off error [30, 32], i.e., the following equilibrium distribution function is commonly used,

[TABLE]

In the following, we suppose $c=1$ such that the relevant quantities are dimensionless.

The above equation can also be written in a vector form:

[TABLE]

where $|f(\bm{x},t)\rangle=(f_{0}(\bm{x},t),f_{1}(\bm{x},t),\cdots,f_{18}(\bm{x},t))^{{}^{\prime}}$ is a column vector, $|f^{(eq)}(\bm{x},t)\rangle$ and $|f(\bm{x}+{\bm{c}_{\alpha}}{\delta_{t}},t+{\delta_{t}})\rangle$ have similar definitions and the superscript ′ represents the transpose operator. For D3Q19 MRT model, a $19\times 19$ transformation matrix is defined as [30]:

[TABLE]

which can transform the distribution function into the moment with the linear mapping $\bm{m}=\mathbf{T}{|f\rangle}$ and $\bm{m}^{(eq)}=\mathbf{T}{|f^{(eq)}\rangle}$, and simultaneously convert the collision matrix into a diagonal one by $\bm{\hat{\Lambda}}=\mathbf{T}\bm{\Lambda}\mathbf{T}^{-1}$. Thus, we can further write Eq. (5) as

[TABLE]

where $\bm{m}$ is defined as

[TABLE]

and $\bm{m}^{(eq)}$ is written as

[TABLE]

where the equilibrium moments are

[TABLE]

where $\bm{j}=(\rho_{0}u_{x},\rho_{0}u_{y},\rho_{0}u_{z})^{\prime}$, $j_{i}=\rho_{0}u_{i}\;(i=x,y,z)$, $\omega_{\epsilon}=3$, $\omega_{\epsilon j}=-11/2$ and $\omega_{xx}=-1/2$. It should be noted that $\omega_{\epsilon}$, $\omega_{\epsilon j}$ and $\omega_{xx}$ do not have much effect on the recovered Navier-Stokes equation. To attain an optimized stability of the model, we can adjust these parameters through linear analysis [30]. In the following, we can see from Yu formula that these parameters do not have effect on the computation of strain rate tensor either.

The diagonal collision matrix $\bm{\hat{\Lambda}}$ is

[TABLE]

where $s_{0}$, $s_{3}$, $s_{5}$ and $s_{7}$ are the relaxation parameters corresponding to the conserved moments. The values of these parameters do not affect the recovered macroscopic N-S equations, which are always set to be zeros. In addition, Eq. (11) can also be written as

[TABLE]

where $s_{c}$, $s_{e}$, $s_{\epsilon}$, $s_{q}$, $s_{\nu}$, $s_{\pi}$ and $s_{t}$ are the parameters corresponding to the conserved moments, the moments related to kinetic energy, energy square, energy flux, strain rate tensor, et al. He-Luo D3Q19 MRT LB model can recovered to the incompressible Navier-Stokes equation by the Chapman-Enskog analysis. In the following, the discussion is based on this model.

3 Two formulas in computing the strain rate tensor in MRT LB model

Yu et al. proposed one formula to compute the strain rate tensor $S_{\alpha\beta}=(\partial_{\alpha}u_{\beta}+\partial_{\beta}u_{\alpha})/2$ for He-Luo D3Q19 MRT LB model [16]. In this formula, the strain rate tensor is computed from non-equilibrium parts of moments. The computation formula is

[TABLE]

where $\varepsilon m_{i}^{(1)}\approx m_{i}-m_{i}^{(eq)}$ is the non-equilibrium part of moment $m_{i}$, $m_{i}$ and $m_{i}^{(eq)}$ are the $i$-th element of vector $\bm{m}$ and $\bm{m}^{(eq)}$.

Chai et al. proposed the other formula to compute the strain rate tensor for D2Q9 MRT LB model [12, 22]. This formula can also be used for D3Q19 MRT LB model. In this formula, the strain rate tensor is computed from non-equilibrium parts of density distribution functions. The computation formula is

[TABLE]

where $\bm{S}$ is strain rate tensor, $f_{j}^{(neq)}={{f}_{j}}-f_{j}^{(eq)}$. If $\bm{{\hat{\Lambda}}}=\frac{1}{\tau}\mathbf{I}$, where $\tau$ is the relaxation time in LBGK model and $\mathbf{I}$ is the unit matrix, the above expression becomes

[TABLE]

which is exactly the strain rate tensor computation formula for LBGK model. Therefore, Chai formula is an extension of strain rate tensor computational formula for LBGK model.

4 Derivation of two formulas in computing the strain rate tensor in MRT LB model

4.1 Derivation of Yu formula in computing the strain rate tensor in MRT LB model

We first expand the density distribution function and the derivatives of time and space as

[TABLE]

Substituting the above expansions into Eq. (1), we obtain the zero-, first-, and second-order equations in $\varepsilon$,

[TABLE]

where ${{D}_{1i}}={{\partial}_{{{t}_{1}}}}+{{c}_{i\gamma}}{{\partial}_{\gamma 1}}$.

If we rewrite Eq. (17) in the vector form and multiply the matrix $\mathbf{T}$ on both sides of them, the corresponding equations in the moment space can be derived

[TABLE]

where ${\hat{\mathbf{\Lambda}}}^{\prime}=\hat{\mathbf{\Lambda}}/{{\delta}_{t}}$, $\hat{\mathbf{\Lambda}}={\mathbf{T}}{\mathbf{\Lambda}}\mathbf{T}^{-1}$, ${{{\tilde{\mathbf{D}}}}_{1}}=\mathbf{T}{{\mathbf{D}}_{1}}{{\mathbf{T}}^{-1}}$ and ${\mathbf{D}_{1}}=diag(\partial_{t_{1}},\partial_{t_{1}}+{{c}_{1\gamma}}{{\partial}_{\gamma 1}},...,\partial_{t_{1}}+{{c}_{18\gamma}}{{\partial}_{\gamma 1}}$). ${{{\tilde{\mathbf{D}}}}_{1}}$ can be computed as

[TABLE]

where ${\mathbf{\hat{C}}_{k}}=\mathbf{T}{\mathbf{C}_{k}}{\mathbf{T}^{-1}},\ k=x,y,z$, $\mathbf{C}_{k}$ is a diagonal matrix with the $k$ component of every discrete velocity $\bm{c}_{i}$ as the diagonal element. ${\mathbf{\hat{C}}_{k}}$ are computed as

[TABLE]

[TABLE]

[TABLE]

After calculating ${\mathbf{\hat{C}}_{k}}\bm{m}^{(eq)},k=x,y,z$ and substituting them into Eq. (18b), we have

[TABLE]

where $s^{\prime}_{i}=s_{i}/\delta t$ and

[TABLE]

From Eq. (23), we have

[TABLE]

From the rows corresponding to the conserved moments in Eq. (23), we have

[TABLE]

where $c_{s}^{2}=1/3$. In fact, above equation is satisfied for both LBGK and MRT LB models, and can be applied to most lattices, such as D2Q9, D3Q15 and D3Q19 lattices. From above equation, we have

[TABLE]

indicating that ${\partial_{{t_{1}}}}{\rho_{0}}{u_{\alpha}}{u_{\beta}}\ (\alpha,\beta=x,y,z)$ is of order $O(Ma^{3})$. Therefore, the terms $\frac{1}{{{\rho}_{0}}}{{\partial}_{t_{1}}}\left({{j}_{\alpha}}{{j}_{\beta}}\right),\alpha,\beta=x,y,z$ in Eqs. (25) can be omitted for the incompressible flows. After multiplying $\varepsilon$ to both sides of Eq. (25), we obtain

[TABLE]

From the equations of $\varepsilon m_{1}^{(1)}$, $\varepsilon m_{9}^{(1)}$ and $\varepsilon m_{11}^{(1)}$, we can compute the diagonal elements of strain rate tensor. And from the equations of $\varepsilon m_{13}^{(1)}$, $\varepsilon m_{14}^{(1)}$ and $\varepsilon m_{15}^{(1)}$, we can calculate the off-diagonal elements of strain rate tensor. The final formula is Eq. (13).

4.2 Derivation of Chai formula in computing the strain rate tensor in MRT LB model

Equilibrium distribution functions for He-Luo MRT LB model satisfy

[TABLE]

From Eqs. (17a), (17b), (26a), (29) and noting ${\mathbf{\Lambda}}={\mathbf{T}^{-1}}\hat{\mathbf{\Lambda}}\mathbf{T}$, we obtain

[TABLE]

Omitting the terms $O(Ma^{3})$ (see Eq. (27)) and multiplying $\varepsilon$ to the both sides of Eq. (30), we obtain

[TABLE]

From above equation, we can easily deduce the Chai formula (Eq. (14)) used for computing the strain rate tensor in MRT LB model.

5 Discussion on two formulas in computing the strain rate tensor in MRT LB model

5.1 The equivalence of two formulas in computing the strain rate tensor in MRT LB model

In view of the form of formulas, Chai formula (Eq. (14)) is very different from Yu formula (Eq. (13)). Yu formula is based on non-equilibrium moments, while Chai formula is based on non-equilibrium density distribution functions. To see if these two formulas are equivalent, we denote Chai formula using non-equilibrium moments. The proof is as follows.

Obviously, ${\bm{m}^{(neq)}}=\mathbf{T}\left|{{f}^{(neq)}}\right\rangle$ and $\varepsilon{\bm{m}^{(1)}}={\bm{m}^{(neq)}}$ are commonly used. Then we calculate $\mathbf{M}={\mathbf{T}^{-1}}\bm{\hat{\Lambda}}{\bm{m}^{(1)}}$, where $\mathbf{M}$ is column vector. After that, we compute ${{N}_{\alpha\beta}}=\sum\nolimits_{i}{{{c}_{i\alpha}}{{c}_{i\beta}}{{M}_{i}}}$, where ${{M}_{i}}$ is the $i$-th element of column vector $\mathbf{M}$. The computation results are

[TABLE]

Finally, we calculate ${{S}_{\alpha\beta}}=-\varepsilon{{N}_{\alpha\beta}}/\left(2{{\rho}_{0}}c_{s}^{2}{{\delta}_{t}}\right)$, where ${{S}_{\alpha\beta}}$ is the element of strain rate tensor. Surprisingly, the final formula obtained from above computation is the same with Yu formula (Eq. (13)).

In Ref. [22], the accuracy of Chai formula is proved to be second-order accurate in space. However, the accuracy of Yu formula has not been studied in previous literature. Because the equivalence of these two formulas is proved, it can be concluded that Yu formula is also second-order accurate in space.

5.2 The computational efficiency of two formulas in computing the strain rate tensor in MRT LB model

To calculate the strain rate tensor with Chai formula, we must firstly calculate all the elements of column vector ${\bf{T}}\left|{{f^{(neq)}}}\right\rangle$, where $\bf{T}$ is a $19\times 19$ transformation matrix, $\left|{{f^{(neq)}}}\right\rangle$ is a column vector and is composed of non-equilibrium distribution functions in nineteen discrete velocity directions. We then calculate ${{\bf{T}}^{-1}}\left({{\bf{\hat{\Lambda}T}}\left|{{f^{({{neq}})}}}\right\rangle}\right)$, where $\bm{\hat{\Lambda}}$ is a diagonal matrix and ${{\bf{T}}^{-1}}$ is a $19\times 19$ matrix. Next, we have to calculate $\sum\limits_{i=1}^{18}{{c_{i\alpha}}{c_{i\beta}}{{\left[{{{\bf{T}}^{-1}}\left({\bm{\hat{\Lambda}}{\bf{T}}\left|{{f^{(neq)}}}\right\rangle}\right)}\right]}_{i}}}$, where $\alpha,\beta=x,y,z$. Finally, the computation results are divided by $-2\rho_{0}c_{s}^{2}{{\delta}_{t}}$ to obtain the strain rate tensor as Eq. (14).

Compared with calculating the strain rate tensor by Chai formula, calculating the strain rate tensor by Yu formula requires much less computation quantity. Through Yu formula, we firstly calculate only 6 elements of column vector ${\bf{T}}\left|{{f^{(neq)}}}\right\rangle$, which are 2nd, 10th, 12th, 14th, 15th and 16th elements of column vector ${\bf{T}}\left|{{f^{(neq)}}}\right\rangle$. The above elements are approximately equal to $\varepsilon m_{1}^{(1)}$, $\varepsilon m_{9}^{(1)}$, $\varepsilon m_{11}^{(1)}$, $\varepsilon m_{13}^{(1)}$, $\varepsilon m_{14}^{(1)}$, $\varepsilon m_{15}^{(1)}$. Then we do some arithmetic operations to obtain the strain rate tensor as Eq. (13). From quantitative analysis, it is known that 957 times of multiplication and 792 times of addition are needed for Chai formula to work out the strain rate tensor at one grid point while only 155 times of multiplication/division and 113 times of addition/subtraction are required for Yu formula. Therefore, to calculate the strain rate tensor in MRT LB model more efficiently, Yu formula is a better choice.

5.3 The generality of two formulas in computing the strain rate tensor in MRT LB model

Yu formula is based on Eq. (28), which is originated from Eq. (23) or Eq. (25). Obviously, Eq. (23) and Eq. (25) are only satisfied for D3Q19 lattice. Therefore, Yu formula can only be applied to D3Q19 MRT LB model to compute the strain rate tensor.

Compared with Yu formula, Chai formula is more general to compute the strain rate tensor for MRT LB model. Chai formula is derived mainly from Eq. (30), which is stemmed from Eqs. (17), (26) and (29). Interestingly, Eqs. (17), (26) and (29) are satisfied for most lattices in the MRT LB model. Thus, Chai formula can be used to calculate the strain rate tensor not only for D3Q19 lattice, but also for other lattices, such as D2Q9 and D3Q15 lattices in the MRT LB model. From above analysis, it can be inferred that, if we want to compute the strain rate tensor for MRT LB model in a unified framework, Chai formula is a better choice.

5.4 Two ways to obtain Yu formula for other lattices in MRT LB model

Because Yu formula is more efficient in computing the strain rate tensor in MRT LB model, and we do not always use D3Q19 lattice, it is worthy to deduce the Yu formula for other lattices, such as the commonly used D2Q9 and D3Q15 lattices. For D$d$Q$q$ lattice, the derivation of Yu formula requires the computation of $d$ matrices ${\mathbf{\hat{C}}_{k}}=\mathbf{T}{\mathbf{C}_{k}}{\mathbf{T}^{-1}}$ and $d$ vectors ${\mathbf{\hat{C}}_{k}}\bm{m}^{(eq)}$, where $\mathbf{T}$, ${\mathbf{C}_{k}}$ and ${\mathbf{T}^{-1}}$ are $q\times q$ matrices. For example, the derivation of Yu formula for D3Q15 lattice needs to compute ${\mathbf{\hat{C}}_{k}}=\mathbf{T}{\mathbf{C}_{k}}{\mathbf{T}^{-1}}$ and $\ {\mathbf{\hat{C}}_{k}}\bm{m}^{(eq)}\ (k=x,y,z)$, where $\mathbf{T}$, ${\mathbf{C}_{k}}$ and ${\mathbf{T}^{-1}}$ are $15\times 15$ matrices. After above calculation, we can obtain similar equations to Eq. (23), Eq. (25) and Eq. (28). Finally, we do some tedious arithmetic operations to deduce the final Yu formula. It should be noted that above derivation can not be conducted only by symbol calculation but needs lots of manual work, which make the derivation very time-consuming.

Compared with above procedure, the derivation of Yu formula from Chai formula is more convenient. The derivation steps are the same with that in section 5.1. Firstly, $\varepsilon{\bm{m}^{(1)}}$ is used to instead of $\mathbf{T}\left|{{f}^{(neq)}}\right\rangle$ in Eq.(14). Then we compute $\mathbf{M}={\mathbf{T}^{-1}}\bm{\hat{\Lambda}}{\bm{m}^{(1)}}$ and ${{N}_{\alpha\beta}}=\sum\nolimits_{i}{{{c}_{i\alpha}}{{c}_{i\beta}}{{M}_{i}}}$, where ${{M}_{i}}$ is the $i$-th element of column vector $\mathbf{M}$. Finally, we compute ${{S}_{\alpha\beta}}=-\varepsilon{{N}_{\alpha\beta}}/\left(2{{\rho}_{0}}c_{s}^{2}{{\delta}_{t}}\right)$ and thus obtain the Yu type of computational formula for stain rate tensor. Fortunately, above derivation way can be conducted completely by symbol calculation, which makes this way more convenient to deduce Yu formula.

6 Conclusion

In this paper, we take He-Luo D3Q19 MRT LB model for an example to study the only two existing computational formulas for strain rate tensor in MRT LB model. These two formulas are named Yu formula and Chai formula in this paper. In view of the expression forms of these two formulas, it seems difficult to establish their correlation. However, through the theoretical analysis and symbol computation, it is found that these two formulas are actually equal to each other. In addition, the advantages and disadvantages of these two formulas are found. Yu formula is deduced through much more matrix computation and arithmetic operations, thus this formula is more specific to one type of MRT LB model and is more efficient in computing the strain rate tensor compared with Chai formula. In the other hand, Chai formula is deduced from some basic equations for most MRT LB models, so Chai formula can be applied to more lattice patterns of MRT LB models than Yu formula. Finally, it is found that to deduce the Yu type of formulas for other lattices, such as D2Q9 lattice and D3Q15 lattice, Chai formula can be used and this way is more convenient than the way proposed by Yu et al. All in all, in this paper it is found that the only two formulas in computing the strain rate tensor for MRT LB model are equal to each other but have their own advantages. To compute the strain rate tensor more efficiently, Yu formula is recommended. To compute the strain rate tensor for different lattice patterns of MRT LB models in a unified framework, Chai formula is suggested.

Acknowledgment

This work is supported by the National Natural Science Foundation of China (Grant No.11502124, Grant No.11302073), the Natural Science Foundation of Zhejiang Province (Grant No.LQ16A020001), the Scientific Research Fund of Zhejiang Provincial Education Department (Grant No.Y201533808), Natural Science Foundation of Ningbo (Grant No. 2016A610075) and is sponsored by K.C. Wong Magna Fund in Ningbo University.