A Chebyshev collocation based sequential matrix exponential method for the generalized density evolution equation

Abstract

In perspective of global approximation, this study presents a numerical method for the generalized density evolution equation (GDEE) based on spectral collocation. A sequential matrix exponential solution based on the Chebyshev collocation points is derived in consideration of the coefficient or velocity term of GDEE being constant in each time step, then the numerical procedure could be successively implemented without implicit or explicit difference schemes. The results of three numerical examples indicate that the solutions in terms of the sequential matrix exponential method for GDEE have good agreement with the analytical results or Monte Carlo simulations. For sufficiently smooth cases, there need no more than one hundred representative points to achieve a satisfied solution by the proposed method, whereas for the case in presence of severe discontinuity a few more sampling points are required to keep numerical stability and accuracy.

Introduction

Many scientific and engineering problems are inevitably related to uncertain factors. Currently, uncertainty quantification has attracted more and more attention. In the framework of probability theory, besides Monte Carlo simulations, there are still a few popular alternative methods, e.g. FPK equation [1], [2], path integral methods [3], [4], stochastic finite element [5], and so on [6]. In recent decades, starting from the principle of preservation of probability [7], Li and Chen developed probability density evolution method (PDEM) for analyzing structural responses with uncertain parameters and random excitations. At first, dynamic responses and reliability analysis of structures subjected to earthquake were performed in [8], [9], [10], [11]. Thereafter, PDEM has been successively applied to other research fields, structural buckling analysis [12], [13], fatigue life and reliability evaluation [14], [15], collapsing and failure analysis of structures [16], [17], etc. The key of PDEM lies in the generalized density evolution equation (GDEE), which reflects the probability preservation law of a reduced-dimensional stochastic dynamic system in view of the random event description [18]. PDEM is also combined with other methods [19], [20], or solved in an integral form [21]. Theoretically, the GDEE holds for any arbitrary dimension [22], whereas one-dimensional GDEE is often of interest in most cases. Hence, in order to carry out PDEM, one has to solve the one-dimensional GDEE numerically.

GDEE is a first order hyperbolic partial differential equation with time-varying coefficients. The first numerical method for PDEE proposed by Li and Chen [23] was a total variation diminishing (TVD) based finite difference method, in which a flux limiter was imposed on a hybrid strategy in combination of the upwind scheme with Lax–Wendroff scheme to eliminate numerical oscillation and keep the non-negative nature of probability density function. Thereafter, Papadopoulos and Kalogeris [24] discussed the conversation properties of the GDEE, and proposed a time-marching discontinuous Galerkin scheme and the Streamline Upwind/Petrov Galerkin (SUPG) scheme for GDEE. It was found that the SUPG method not restricted by the CFL condition improved the computational efficiency. To overcome the mesh sensitivity of finite difference for GDEE, a difference-nonlinear wavelet density estimation method was proposed in [25]. Taking advantages of the multi-resolution property of wavelet functions and choosing the optimal scale at each time step, a better numerical stability and results could be achieved. In addition, Tao and Li [26] improved the original finite difference method by replacing the point velocity with the ensemble velocity incorporating a first correction term. In this way, the influence of the additionally generated samples in each probability subdomain could be taken into consideration to improve the numerical results. Recently, a new enrichment strategy based on reproducing kernel particle method (RKPM) was proposed for solving the generalized probability density evolution equation [27]. By taking the RKPM as a surrogate model instead of performing time-consuming deterministic dynamic analysis, more representative points were efficiently produced, so that one can get a more accurate solution to PDF. The extra computation cost is low, because there needs just one RKPM based surrogate model for the refinement process.

So far, all numerical methods for PDEM are in the framework of finite difference, which are local interpolation approximation for GDEE. Whereas, since 1960s spectral methods [28] as a global approximation have been successively applied to solve differential equation in many subjects, particularly, fluid dynamics [29], [30]. For this issue, spectral methods for hyperbolic problems were in details reviewed in [31]. Theoretically, for smooth problems spectral methods have a more efficient convergence rate than finite element method and finite difference method. The most popular spectral methods include Galerkin method, Tau method, collocation method or pseudospectral method [32]. Fundamental Galerkin method and Tau method are implemented with the expansion coefficients as unknowns, while collocation method is implemented with the nodal values as unknowns. Among these spectral methods, the collocation method is more suitable for problems with variable coefficient and convenient to deal with non-linear terms. The probability density evolution equation is a linearly hyperbolic partial differential equation, generally with time-dependent coefficients, which brings most difficulties for solving it. In this study, taking the coefficient or velocity term of GDEE as a constant in each time step, we derive a sequential matrix exponential solution to GDEE based on the Chebyshev collocation method, which could be successively implemented in replace of implicit or explicit difference schemes.

The reset of this paper is organized as follows. In Section 2, a brief introduction to probability density evolution method is presented. In Section 3, the Chebyshev polynomials and theoretical solution form of linearly partial differential equation are introduced. Based on these results, collocation point formula for GDEE are derived and corresponding numerical algorithm is presented. As numerical examples, the proposed method for an undamped oscillator with random frequency and Riccati equation is validated by analytical solutions, and a clamped bar with random imperfections and impulse load is verified by Monte Carlo simulations in Section 4. Finally, some conclusions and discussions are drawn in Section 5.