Lijian Jiang

Xuehan Zhang, Lijian Jiang, Eric T. Chung

Deep learning-based surrogate models have been extensively developed for efficiently approximating multiscale systems with random input fields. However, most existing approaches require retraining neural networks from scratch when source terms, boundary conditions, or differential operators change, resulting in significant computational costs and limited adaptability. To address this challenge, we integrate our previous CNN-based reduced-order model (ROM) framework with the multiscale finite element method (MsFEM) and propose an MsFEM-inspired transfer learning strategy, termed MITL. The CNN-based ROM consists of two components: Basis CNNs, which learn reduced basis functions, and Coef CNNs, which predict the corresponding linear combination coefficients. To enhance the transferability of learned multiscale representations, global MsFEM basis problems are employed as source tasks during pretraining. For new target problems, MITL requires training only lightweight adaptation networks to construct task-specific reduced bases and coefficients, thereby substantially reducing the computational burden. Numerical experiments demonstrate that MITL achieves accurate and efficient predictions across a range of target tasks, with particularly significant advantages in data-scarce scenarios.

Mengnan Li, Eric Chung, Lijian Jiang

In this paper, we present a Constraint Energy Minimizing Generalized Multiscale Finite Element Method (CEM-GMsFEM) for parabolic equations with multiscale coefficients, arising from applications in porous media. We will present the construction of CEM-GMsFEM and rigorously analyze its convergence for the parabolic equations. The convergence rate is characterized by the coarse grid size and the eigenvalue decay of local spectral problems, but is independent of the scale length and contrast of the media. The analysis shows that the method has a first order convergence rate with respect to coarse grid size in the energy norm and second order convergence rate with respect to coarse grid size in $L^2$ norm under some appropriate assumptions. For the temporal discretization, finite difference techniques are used and the convergence analysis of full discrete scheme is given. Moreover, a posteriori error estimator is derived and analyzed. A few numerical results for porous media applications are presented to confirm the theoretical findings and demonstrate the performance of the approach.

Lijian Jiang, Qiuqi Li

In this paper, we propose a model's sparse representation based on reduced mixed generalized multiscale finite element (GMsFE) basis methods for elliptic PDEs with random inputs. Mixed generalized multiscale finite element method (GMsFEM) is one of the accurate and efficient approaches to solve multiscale problem in a coarse grid with local mass conservation. When the inputs of the PDEs are parameterized by the random variables, the GMsFE basis functions usually depend on the random parameters. This leads to a large number degree of freedoms for the mixed GMsFEM and substantially impacts on the computation efficiency. In order to overcome the difficulty, we develop reduced mixed GMsFE basis methods such that the multiscale basis functions are independent of the random parameters and span a low-dimensional space. To this end, a greedy algorithm is used to find a set of optimal samples from a training set scattered in the parameter space. Reduced mixed GMsFE basis functions are constructed based on the optimal samples using two optimal sampling strategies: basis-oriented cross-validation and proper orthogonal decomposition. Although the dimension of the space spanned by the reduced mixed GMsFE basis functions is much smaller than the dimension of the original full order model, the online computation still depends on the number of coarse degree of freedoms. To significantly improve the online computation, we integrate the reduced mixed GMsFE basis methods with sparse tensor approximation and obtain a sparse representation for the model's outputs. The sparse representation is very efficient for evaluating the model's outputs for many instances of parameters. To illustrate the efficacy of the proposed methods, we present a few numerical examples for multsicale problems with random inputs.

Lijian Jiang, Na Ou

This work presents a model reduction approach to the inverse problem in the application of subsurface flows. For the Bayesian inverse problem, the forward model needs to be repeatedly computed for a large number of samples to get a stationary chain. This requires large computational efforts. To significantly improve the computation efficiency, we use generalized multiscale finite element method and least-squares stochastic collocation method to construct a reduced computational model. To avoid the difficulty of choosing regularization parameter, hyperparameters are introduced to build a hierarchical model. We use truncated Karhunen-Loeve expansion (KLE) to reduce the dimension of the parameter spaces and decrease the mixed time of Markov chains. The techniques of hyperparameter and KLE are incorporated into the model reduction method. The reduced model is constructed offline. Then it is computed very efficiently in the online sampling stage. This strategy can significantly accelerate the evaluation of the Markov chain and the resultant posterior distribution converges fast. We analyze the convergence for the approximation between the posterior distribution by the reduced model and the reference posterior distribution by the full-order model. A few numerical examples in subsurface flows are carried out to demonstrate the performance of the presented model reduction method with application of the Bayesian inverse problem.

Yuming Ba, Lijian Jiang, Na Ou

Ensemble Kalman filter (EnKF) has been widely used in state estimation and parameter estimation for the dynamic system where observational data is obtained sequentially in time. To reduce uncertainty and accelerate posterior inference, a two-stage ensemble Kalman filter is presented to improve the sequential analysis of EnKF. It is known that the final posterior ensemble may be concentrated in a small portion of the entire support of the initial prior ensemble. It will be much more efficient if we first build a new prior by some partial observations, and construct a surrogate only over the significant region of the new prior. To this end, we construct a very coarse model using generalized multiscale finite element method (GMsFEM) and generate a new prior ensemble in the first stage. GMsFEM provides a set of hierarchical multiscale basis functions supported in coarse blocks. This gives flexibility and adaptivity to choosing degree of freedoms to construct a reduce model. In the second stage, we build an initial surrogate model based on the new prior by using GMsFEM and sparse generalized polynomial chaos (gPC)-based stochastic collocation methods. To improve the initial surrogate model, we dynamically update the surrogate model, which is adapted to the sequential availability of data and the updated analysis. The two-stage EnKF can achieve a better estimation than standard EnKF, and significantly improve the efficiency to update the ensemble analysis (posterior exploration). To enhance the applicability and flexibility in Bayesian inverse problems, we extend the two-stage EnKF to non-Gaussian models and hierarchical models. In the paper, we focus on the time fractional diffusion-wave models in porous media and investigate their Bayesian inverse problems using the proposed two-stage EnKF.

Lijian Jiang, Ilya D. Mishev

We develop a framework for constructing mixed multiscale finite volume methods for elliptic equations with multiple scales arising from flows in porous media. Some of the methods developed using the framework are already known \cite{jennylt03}; others are new. New insight is gained for the known methods and extra flexibility is provided by the new methods. We give as an example a mixed MsFV on uniform mesh in 2-D. This method uses novel multiscale velocity basis functions that are suited for using global information, which is often needed to improve the accuracy of the multiscale simulations in the case of continuum scales with strong non-local features. The method efficiently captures the small effects on a coarse grid. We analyze the new mixed MsFV and apply it to solve two-phase flow equations in heterogeneous porous media. Numerical examples demonstrate the accuracy and efficiency of the proposed method for modeling the flows in porous media with non-separable and separable scales.

Lijian Jiang, J. David Moulton, Jia Wei

Stochastic modeling has become a popular approach to quantify uncertainty in flows through heterogeneous porous media. The uncertainty in heterogeneous structure properties is often parameterized by a high-dimensional random variable. This leads to a deterministic problem in a high-dimensional parameter space and the numerical computation becomes very challengeable as the dimension of the parameter space increases. To efficiently tackle the high-dimensionality, we propose a hybrid high dimensional model representation (HDMR) technique, through which the high-dimensional stochastic model is decomposed into a moderate-dimensional stochastic model in a most active random space and a few one-dimensional stochastic models. The derived low-dimensional stochastic models are solved by incorporating sparse grid stochastic collocation method into the proposed hybrid HDMR. The porous media properties such as permeability are often heterogeneous. To treat the heterogeneity, we use a mixed multiscale finite element method (MMsFEM) to simulate each of derived stochastic models. To capture the non-local spatial features of the porous media and the important effects of random variables, we can hierarchically incorporate the global information individually from each of random parameters. This significantly enhances the accuracy of the multiscale simulation. The synergy of the hybrid HDMR and the MMsFEM reduces the stochastic model of flows in both stochastic space and physical space, and significantly decreases the computation complexity. We carefully analyze the proposed HDMR technique and the derived stochastic MMsFEM. A few numerical experiments are carried out for two-phase flows in random porous media and support the efficiency and accuracy of the MMsFEM based on the hybrid HDMR.

Lijian Jiang, Dylan Copeland, J. David Moulton

We develop a family of expanded mixed Multiscale Finite Element Methods (MsFEMs) and their hybridizations for second-order elliptic equations. This formulation expands the standard mixed Multiscale Finite Element formulation in the sense that four unknowns (hybrid formulation) are solved simultaneously: pressure, gradient of pressure, velocity and Lagrange multipliers. We use multiscale basis functions for the both velocity and gradient of pressure. In the expanded mixed MsFEM framework, we consider both cases of separable-scale and non-separable spatial scales. We specifically analyze the methods in three categories: periodic separable scales, $G$- convergence separable scales, and continuum scales. When there is no scale separation, using some global information can improve accuracy for the expanded mixed MsFEMs. We present rigorous convergence analysis for expanded mixed MsFEMs. The analysis includes both conforming and nonconforming expanded mixed MsFEM. Numerical results are presented for various multiscale models and flows in porous media with shales to illustrate the efficiency of the expanded mixed MsFEMs.

Xiaoyan Song, Guanghui Zheng, Lijian Jiang

In this paper, we present an inverse problem of identifying the reaction coefficient for time fractional diffusion equations in two dimensional spaces by using boundary Neumann data. It is proved that the forward operator is continuous with respect to the unknown parameter. Because the inverse problem is often ill-posed, regularization strategies are imposed on the least fit-to-data functional to overcome the stability issue. There may exist various kinds of functions to reconstruct. It is crucial to choose a suitable regularization method. We present a multi-parameter regularization $L^{2}+BV$ method for the inverse problem. This can extend the applicability for reconstructing the unknown functions. Rigorous analysis is carried out for the inverse problem. In particular, we analyze the existence and stability of regularized variational problem and the convergence. To reduce the dimension in the inversion for numerical simulation, the unknown coefficient is represented by a suitable set of basis functions based on a priori information. A few numerical examples are presented for the inverse problem in time fractional diffusion equations to confirm the theoretic analysis and the efficacy of the different regularization methods.

Lingling Ma, Qiuqi Li, Lijian Jiang

In this paper, a local-global model reduction method is presented to solve stochastic optimal control problems governed by partial differential equations (PDEs). If the optimal control problems involve uncertainty, we need to use a few random variables to parameterize the uncertainty. The stochastic optimal control problems require solving coupled optimality system for a large number of samples in the stochastic space to quantify the statistics of the system response and explore the uncertainty quantification. Thus the computation is prohibitively expensive. To overcome the difficulty, model reduction is necessary to significantly reduce the computation complexity. We exploit the advantages from both reduced basis method and Generalized Multiscale Finite Element Method (GMsFEM) and develop the local-global model reduction method for stochastic optimal control problems with PDE constraints. This local-global model reduction can achieve much more computation efficiency than using only local model reduction approach and only global model reduction approach. We recast the stochastic optimal problems in the framework of saddle-point problems and analyze the existence and uniqueness of the optimal solutions of the reduced model. In the local-global approach, most of computation steps are independent of each other. This is very desirable for scientific computation. Moreover, the online computation for each random sample is very fast via the proposed model reduction method. This allows us to compute the optimality system for a large number of samples. To demonstrate the performance of the local-global model reduction method, a few numerical examples are provided for different stochastic optimal control problems.

Lijian Jiang, Michael Presho

In this paper we use a splitting technique to develop new multiscale basis functions for the multiscale finite element method (MsFEM). The multiscale basis functions are iteratively generated using a Green's kernel. The Green's kernel is based on the first differential operator of the splitting. The proposed MsFEM is applied to deterministic elliptic equations and stochastic elliptic equations, and we show that the proposed MsFEM can considerably reduce the dimension of the random parameter space for stochastic problems. By combining the method with sparse grid collocation methods, the need for a prohibitive number of deterministic solves is alleviated. We rigorously analyze the convergence of the proposed method for both the deterministic and stochastic elliptic equations. Computational complexity discussions are also offered to supplement the convergence analysis. A number of numerical results are presented to confirm the theoretical findings.

Na Ou, Guang Lin, Lijian Jiang

This work presents a multiscale model reduction approach to discontinuous fields identification problems in the framework of Bayesian inference. An ensemble-based variable separation (VS) method is proposed to approximate multiscale basis functions used to build a coarse model. The variable-separation expression is constructed for stochastic multiscale basis functions based on the random field, which is treated Gauss process as prior information. To this end, multiple local inhomogeneous Dirichlet boundary condition problems are required to be solved, and the ensemble-based method is used to obtain variable separation forms for the corresponding local functions. The local functions share the same interpolate rule for different physical basis functions in each coarse block. This approach significantly improves the efficiency of computation. We obtain the variable separation expression of multiscale basis functions, which can be used to the models with different boundary conditions and source terms, once the expression constructed. The proposed method is applied to discontinuous field identification problems where the hybrid of total variation and Gaussian (TG) densities are imposed as the penalty. We give a convergence analysis of the approximate posterior to the reference one with respect to the Kullback-Leibler (KL) divergence under the hybrid prior. The proposed method is applied to identify discontinuous structures in permeability fields. Two patterns of discontinuous structures are considered in numerical examples: separated blocks and nested blocks.

Xiaoyan Song, Lijian Jiang, Guanghui Zheng

This paper presents an improved implicit sampling method for hierarchical Bayesian inverse problems. A widely used approach for sampling posterior distribution is based on Markov chain Monte Carlo (MCMC). However, the samples generated by MCMC are usually strongly correlated. This may lead to a small size of effective samples from a long Markov chain and the resultant posterior estimate may be inaccurate. An implicit sampling method proposed in [11] can generate independent samples and capture some inherent non-Gaussian features of the posterior based on the weights of samples. In the implicit sampling method, the posterior samples are generated by constructing a map and distribute around the MAP point. However, the weights of implicit sampling in previous works may cause excessive concentration of samples and lead to ensemble collapse. To overcome this issue, we propose a new weight formulation and make resampling based on the new weights. In practice, some parameters in prior density are often unknown and a hierarchical Bayesian inference is necessary for posterior exploration. To this end, the hierarchical Bayesian formulation is used to estimate the MAP point and integrated in the implicit sampling framework. Compared to conventional implicit sampling, the proposed implicit sampling method can significantly improve the posterior estimator and the applicability for high dimensional inverse problems. The improved implicit sampling method is applied to the Bayesian inverse problems of multi-term time fractional diffusion models in heterogeneous media. To effectively capture the heterogeneity effect, we present a mixed generalized multiscale finite element method (mixed GMsFEM) to solve the time fractional diffusion models in a coarse grid, which can substantially speed up the Bayesian inversion.

Qiuqi Li, Lijian Jiang

In this paper, we propose a novel variable-separation (NVS) method for generic multivariate functions. The idea of NVS is extended to to obtain the solution in tensor product structure for stochastic partial differential equations (SPDEs). Compared with many widely used variation-separation methods, NVS shares their merits but has less computation complexity and better efficiency. NVS can be used to get the separated representation of the solution for SPDE in a systematic enrichment manner. No iteration is performed at each enrichment step. This is a significant improvement compared with proper generalized decomposition. Because the stochastic functions of the separated representations obtained by NVS depend on the previous terms, this impacts on the computation efficiency and brings great challenge for numerical simulation for the problems in high stochastic dimensional spaces. In order to overcome the difficulty, we propose an improved least angle regression algorithm (ILARS) and a hierarchical sparse low rank tensor approximation (HSLRTA) method based on sparse regularization. For ILARS, we explicitly give the selection of the optimal regularization parameters at each step based on least angle regression algorithm (LARS) for lasso problems such that ILARS is much more efficient. HSLRTA hierarchically decomposes a high dimensional problem into some low dimensional problems and brings an accurate approximation for the solution to SPDEs in high dimensional stochastic spaces using limited computer resource. A few numerical examples are presented to illustrate the efficacy of the proposed methods.

Yuming Ba, Lijian Jiang

In the paper, we develop an ensemble-based implicit sampling method for Bayesian inverse problems. For Bayesian inference, the iterative ensemble smoother (IES) and implicit sampling are integrated to obtain importance ensemble samples, which build an importance density. The proposed method shares a similar idea to importance sampling. IES is used to approximate mean and covariance of a posterior distribution. This provides the MAP point and the inverse of Hessian matrix, which are necessary to construct the implicit map in implicit sampling. The importance samples are generated by the implicit map and the corresponding weights are the ratio between the importance density and posterior density. In the proposed method, we use the ensemble samples of IES to find the optimization solution of likelihood function and the inverse of Hessian matrix. This approach avoids the explicit computation for Jacobian matrix and Hessian matrix, which are very computationally expensive in high dimension spaces. To treat non-Gaussian models, discrete cosine transform and Gaussian mixture model are used to characterize the non-Gaussian priors. The ensemble-based implicit sampling method is extended to the non-Gaussian priors for exploring the posterior of unknowns in inverse problems. The proposed method is used for each individual Gaussian model in the Gaussian mixture model. The proposed approach substantially improves the applicability of implicit sampling method. A few numerical examples are presented to demonstrate the efficacy of the proposed method with applications of inverse problems for subsurface flow problems and anomalous diffusion models in porous media.

Fuchen Chen, Lijian Jiang, Guanghui Zheng

In this paper, we present a dimension reduction method to reduce the dimension of parameter space and state space and efficiently solve inverse problems. To this end, proper orthogonal decomposition (POD) and radial basis function (RBF) are combined to represent the solution of forward model with a form of variable separation. This POD-RBF method can be used to efficiently evaluate the model's output. A gradient regularization method is presented to solve the inverse problem with fast convergence. A generalized cross validation method is suggested to select the regularization parameter and differential step size for the gradient computation. Because the regularization method needs many model's evaluations. This is desirable for POD-RBF method. Thus, the POD-RBF method is integrated with the gradient regularization method to provide an efficient approach to solve inverse problems. We focus on the coefficient inversion of diffusion equations using the proposed approach. Based on different types of measurement data and different basis functions for coefficients, we present a few numerical examples for the coefficient inversion. The numerical results show that accurate reconstruction for the coefficient can be achieved efficiently.

Dengfei Zeng, Lijian Jiang, Shuyu Sun et al.

A likelihood-free transport filtering method is proposed based on the couplings between state and observation variables. By exploiting a block-triangular structure in the transport map, the analysis step of filtering is reformulated as the minimization of the maximum mean discrepancy (MMD) between the true joint measure and its transport-based approximation. To circumvent the non-convexity in the MMD optimization, we introduce a training-free transport filter method via gradient flows, which leads to an analytic computation for the transport map that implies the steepest descent direction of the MMD. The proposed approach accurately approximates non-Gaussian filtering posteriors and avoids particle collapse. We provide a convergence analysis for the expectation of the MMD between the approximated posterior and the truth posterior. Finally, we extend the method to high-dimensional problems through domain localization. Numerical examples demonstrate the superior performance of our approach over conventional filtering methods in nonlinear, non-Gaussian scenarios.

Eric T. Chung, Lijian Jiang, Mengnan Li et al.

Simulating diffusion in heterogeneous media presents a significant computational challenge, as resolving microscopic physical scales traditionally demands excessively fine computational grids. To overcome this barrier, we extend the Localized Subspace Iteration (LSI) framework to multiscale parabolic equations. The proposed method constructs optimal, low-dimensional trial spaces by iteratively approximating the dominant eigenspaces of local inverse operators via Localized Standard Subspace Iteration (LSSI) or Localized Krylov Subspace Iteration (LKSI). Because these LSI basis functions are inherently tailored to capture the slow-decaying, low-frequency modes of the parabolic solution, they naturally suppress error accumulation over long-term integration. To further improve computational efficiency, we decouple the basis construction into an offline phase and implement a contrast-independent, partially explicit temporal splitting scheme for online time-stepping. By explicitly advancing the dominant macroscopic modes while implicitly treating high-frequency microscopic corrections, this scheme guarantees stability without imposing restrictive time-step constraints. We establish rigorous a priori error estimates in both the energy and $L^2$ norms. Numerical experiments illustrate the accuracy and efficiency of the LSI framework, particularly highlighting the LKSI method's advantages in handling high-contrast, complex multiscale media.

Xiaofei Guan, Lijian Jiang, Yajun Wang et al.

This paper proposes localized subspace iteration (LSI) methods to construct generalized finite element basis functions for elliptic problems with multiscale coefficients. The key components of the proposed method consist of the localization of the original differential operator and the subspace iteration of the corresponding local spectral problems, where the localization is conducted by enforcing the local homogeneous Dirichlet condition and the partition of the unity functions. From a novel perspective, some multiscale methods can be regarded as one iteration step under approximating the eigenspace of the corresponding local spectral problems. Vice versa, new multiscale methods can be designed through subspaces of spectral problem algorithms. Then, we propose the efficient localized standard subspace iteration (LSSI) method and the localized Krylov subspace iteration (LKSI) method based on the standard subspace and Krylov subspace, respectively. Convergence analysis is carried out for the proposed method. Various numerical examples demonstrate the effectiveness of our methods. In addition, the proposed methods show significant superiority in treating long-channel cases over other well-known multiscale methods.

Dengfei Zeng, Lijian Jiang

In this paper, we present a new ensemble-based filter method by reconstructing the analysis step of the particle filter through a transport map, which directly transports prior particles to posterior particles. The transport map is constructed through an optimization problem described by the Maximum Mean Discrepancy loss function, which matches the expectation information of the approximated posterior and reference posterior. The proposed method inherits the accurate estimation of the posterior distribution from particle filtering while gives an extension to high dimensional assimilation problems. To improve the robustness of Maximum Mean Discrepancy, a variance penalty term is used to guide the optimization. It prioritizes minimizing the discrepancy between the expectations of highly informative statistics for the reference posteriors. The penalty term significantly enhances the robustness of the proposed method and leads to a better approximation of the posterior. A few numerical examples are presented to illustrate the advantage of the proposed method over ensemble Kalman filter.

Ziyi Wang, Lijian Jiang

Data assimilation (DA) has increasingly emerged as a critical tool for state estimation across a wide range of applications. It is signiffcantly challenging when the governing equations of the underlying dynamics are unknown. To this end, various machine learning approaches have been employed to construct a surrogate state transition model in a supervised learning framework, which relies on pre-computed training datasets. However, it is often infeasible to obtain noise-free ground-truth state sequences in practice. To address this challenge, we propose a novel method that integrates reinforcement learning with ensemble-based Bayesian ffltering methods, enabling the learning of surrogate state transition model for unknown dynamics directly from noisy observations, without using true state trajectories. Speciffcally, we treat the process for computing maximum likelihood estimation of surrogate model parameters as a sequential decision-making problem, which can be formulated as a discretetime Markov decision process (MDP). Under this formulation, learning the surrogate transition model is equivalent to ffnding an optimal policy of the MDP, which can be effectively addressed using reinforcement learning techniques. Once the model is trained offfine, state estimation can be performed in the online stage using ffltering methods based on the learned dynamics. The proposed framework accommodates a wide range of observation scenarios, including nonlinear and partially observed measurement models. A few numerical examples demonstrate that the proposed method achieves superior accuracy and robustness in high-dimensional settings.

Xuehan Zhang, Lijian Jiang

In this article, we present a data-driven method for parametric models with noisy observation data. Gaussian process regression based reduced order modeling (GPR-based ROM) can realize fast online predictions without using equations in the offline stage. However, GPR-based ROM does not perform well for complex systems since POD projection are naturally linear. Conditional variational autoencoder (CVAE) can address this issue via nonlinear neural networks but it has more model complexity, which poses challenges for training and tuning hyperparameters. To this end, we propose a framework of CVAE with Gaussian process regression recognition (CVAE-GPRR). The proposed method consists of a recognition model and a likelihood model. In the recognition model, we first extract low-dimensional features from data by POD to filter the redundant information with high frequency. And then a non-parametric model GPR is used to learn the map from parameters to POD latent variables, which can also alleviate the impact of noise. CVAE-GPRR can achieve the similar accuracy to CVAE but with fewer parameters. In the likelihood model, neural networks are used to reconstruct data. Besides the samples of POD latent variables and input parameters, physical variables are also added as the inputs to make predictions in the whole physical space. This can not be achieved by either GPR-based ROM or CVAE. Moreover, the numerical results show that CVAE-GPRR may alleviate the overfitting issue in CVAE.

Lijian Jiang, Na Ou

In the paper, we present a strategy for accelerating posterior inference for unknown inputs in time fractional diffusion models. In many inference problems, the posterior may be concentrated in a small portion of the entire prior support. It will be much more efficient if we build and simulate a surrogate only over the significant region of the posterior. To this end, we construct a coarse model using Generalized Multiscale Finite Element Method (GMsFEM), and solve a least-squares problem for the coarse model with a regularizing Levenberg-Marquart algorithm. An intermediate distribution is built based on the approximate sampling distribution. For Bayesian inference, we use GMsFEM and least-squares stochastic collocation method to obtain a reduced coarse model based on the intermediate distribution. To increase the sampling speed of Markov chain Monte Carlo, the DREAM$_\text{ZS}$ algorithm is used to explore the surrogate posterior density, which is based on the surrogate likelihood and the intermediate distribution. The proposed method with lower gPC order gives the approximate posterior as accurate as the the surrogate model directly based on the original prior. A few numerical examples for time fractional diffusion equations are carried out to demonstrate the performance of the proposed method with applications of the Bayesian inversion.