Computation of the analysis error covariance in variational data assimilation problems with nonlinear dynamics

The problem of variational data assimilation for a nonlinear evolution model is formulated as an optimal control problem to find the initial condition function. The data contain errors (observation and background errors), hence there will be errors in the optimal solution. For mildly nonlinear dynamics, the covariance matrix of the optimal solution error can often be approximated by the inverse Hessian of the cost functional. Here we focus on highly nonlinear dynamics, in which case this approximation may not be valid. The equation relating the optimal solution error and the errors of the input data is used to construct an approximation of the optimal solution error covariance. Two new methods for computing this covariance are presented: the fully nonlinear ensemble method with sampling error compensation and the ‘effective inverse Hessian’ method. The second method relies on the efficient computation of the inverse Hessian by the quasi-Newton BFGS method with preconditioning. Numerical examples are presented for the model governed by Burgers equation with a nonlinear viscous term.     

托卡马克中由线积分数据重建分布的神经网络方法

用多层前馈神经网络,处理托卡马克中由物理量沿观察弦的线积分值重建其空间分布的反演问题。用BFGS拟牛顿法大大加快了标准误差反向传播算法(BP)的收敛速度     

On the accuracy of covariance matrix: Hessian versus Gauss-Newton methods in atmospheric remote sensing with infrared spectroscopy

Most retrieval schemes use a linear approximation of the radiative transfer function within each iteration as well as for error analysis. Like most standard methods, the improved Hessian method relies on a quadratic form of the cost function and linear approximation in the error analysis. Often, there is no robust criterion in determining step size that can be used to calculate covariance matrix by discrete perturbation of the cost function in the Hessian approach. The Hessian method improved recently, however, overcomes this problem by employing adaptive algorithm which uses small step sizes in steep directions and large step sizes in flat directions of the cost function. The results of retrievals of atmospheric trace gases from simulated limb emission spectra show that Gauss-Newton algorithm and the improved Hessian generally give nearly identical volume mixing ratios and error covariance matrices in the original state vector space. Due to interlevel correlations, however, the agreement in the uncertainities in the original state vector coordinate system is partly lost in a space in which the elements of state vector are independent after orthogonal coordinate transformation. The significant discrepancies between the estimated uncertainities by the two methods are found to be related with elements of state vector that are dominantly controlled by flattest eigenvector directions of the inverse covariance matrix. The improved Hessian method determines the uncertainities in those shallowest directions with better accuracy than Gauss-Newton approach. The performance of the Hessian method is also found to be better in resolving structures related to the shallowest eigenvector directions as revealed by better vertical resolutions in the retrieved profiles of the trace species.     

Analysis of the accuracy of the inverse problem solution for a differential heterodyne microscope as applied to rectangular plasmonic waveguides

The existence, uniqueness, and stability of the inverse problem solution for a scanning differential heterodyne microscope as applied to rectangular plasmonic waveguides have been analyzed. The consideration is based on an algorithm using a trial-and-error method that we proposed previously to characterize plasmonic waveguides with a triangular profile. The error of the inverse problem (IP) solution is calculated as dependent on the initial data and with allowance for their errors. Instability domains are found for the IP solution, where the solution error sharply increases. It is shown that the instability domains can be eliminated and the accuracy of the IP solution can be significantly improved in the entire range of initial data by taking initial data in the form of two phase responses of the microscope at different wavelengths.     

基于无导数优化方法的数值模式误差估计

初始场误差和模式误差是制约数值预报准确率的两个关键因素,本文主要考虑利用历史观测资料实现时空演变的模式误差的估计问题.通过把模式误差综合考虑成为准确模式中的未知项,把历史资料看作是带有未知项的准确模式的特解,构造了求解时空演变的模式误差项的反问题及其最优控制问题.给出了一个解决最优控制问题的无导数优化方法,该方法的优点是不需要建立原数值模式的切线性模式与伴随模式,它只需在增加一个外强迫项的基础上运行原数值模式即可实现模式误差项的最优估计.关于Burgers方程的算例表明,无论模式的初始状态是否准确已知,无导数优化方法都能有效解决时空演变的模式误差的最优估计问题,它为实际业务模式利用历史数据提取模式误差信息并显著地改进预报效果提供了一种方便可行的数值方法与理论依据.     

Inverse radiation problem of temperature distribution in one-dimensional isotropically scattering participating slab with variable refractive index

In this paper, an inverse analysis is performed for estimation of source term distribution from the measured exit radiation intensities at the boundary surfaces in a one-dimensional absorbing, emitting and isotropically scattering medium between two parallel plates with variable refractive index. The variation of refractive index is assumed to be linear. The radiative transfer equation is solved by the constant quadrature discrete ordinate method. The inverse problem is formulated as an optimization problem for minimizing an objective function which is expressed as the sum of square deviations between measured and estimated exit radiation intensities at boundary surfaces. The conjugate gradient method is used to solve the inverse problem through an iterative procedure. The effects of various variables on source estimation are investigated such as type of source function, errors in the measured data and system parameters, gradient of refractive index across the medium, optical thickness, single scattering albedo and boundary emissivities. The results show that in the case of noisy input data, variation of system parameters may affect the inverse solution, especially at high error values in the measured data. The error in measured data plays more important role than the error in radiative system parameters except the refractive index distribution; however the accuracy of source estimation is very sensitive toward error in refractive index distribution. Therefore, refractive index distribution and measured exit intensities should be measured accurately with a limited error bound, in order to have an accurate estimation of source term in a graded index medium.     

Design of an optimal wave-vector filter for enhancing the resolution of reconstructed source field by near-field acoustical holography (NAH)

In near-field acoustical holography using the boundary element method, the reconstructed field often diverges due to the presence of small measurement errors. In order to handle this instability in the inverse problem, the reconstruction process should include some form of regularization for enhancing the resolution of source images. The usual method of regularization has been the truncation of wave vectors associated with small singular values, although the determination of an optimal truncation order is difficult. In this article, an iterative inverse solution technique is suggested in which the mean-square error prediction is used. A statistical estimation of the minimum mean-square error between measured pressures and the model solution is required for yielding the optimal number of iterations. The continuous curve of an optimal wave-vector filter is designed, for suppressing the high-order modes that can produce large reconstruction errors. Experimental results from a baffled radiator reveal that the reconstruction errors can be reduced by this form of regularization, by at least 48% compared to those without any regularization. In comparison to results using the optimal truncation method of regularization, the new scheme is shown to give further reductions of truncation error of between 7% and 39%, for the example in this article.     

An approach to estimating and extrapolating model error based on inverse problem methods:towards accurate numerical weather prediction

Model error is one of the key factors restricting the accuracy of numerical weather prediction(NWP). Considering the continuous evolution of the atmosphere, the observed data(ignoring the measurement error) can be viewed as a series of solutions of an accurate model governing the actual atmosphere. Model error is represented as an unknown term in the accurate model, thus NWP can be considered as an inverse problem to uncover the unknown error term. The inverse problem models can absorb long periods of observed data to generate model error correction procedures. They thus resolve the deficiency and faultiness of the NWP schemes employing only the initial-time data. In this study we construct two inverse problem models to estimate and extrapolate the time-varying and spatial-varying model errors in both the historical and forecast periods by using recent observations and analogue phenomena of the atmosphere. Numerical experiment on Burgers’ equation has illustrated the substantial forecast improvement using inverse problem algorithms. The proposed inverse problem methods of suppressing NWP errors will be useful in future high accuracy applications of NWP.     

遥感数据误差对地表参数定量反演可靠性的影响——以太湖叶绿素a反演为例

地表参数的遥感反演误差的大小由遥感数据误差和反演模型误差共同构成,数据误差不是简单地“加减”到反演误差中,而是经过反演模型改造后融入到反演误差中。因此,在地表参数定量反演过程中,用回归系数最大或均方根误差最佳代价函数来描述地表参数与遥感反射光谱之间的关系将不太可靠。从理论上指出了最高回归系数或最小均方根误差评价方法失效的根源在于反演模型的形式,以2003年10月27日太湖实测数据为例进行了论证。研究表明,虽然TM2/TM3算法比TM2/TM1算法的回归系数高,但其对数据误差的放大效果是TM2/TM算法的2.28倍,这导致了反演结果的均方根误差比TM2/TM1算法大了7.938 5 μg·L-1;另外从定量反演结果来看,基于TM2/TM3算法和基于TM2/TM1算法的反演结果完全相反,与以往研究成果对比可知,基于TM2/TM1算法的反演结果更符合实际。因此,数据误差应该作为一个约束条件,加入到代价函数的求解过程中,才能增加反演结果的可靠性。     

透平叶栅三维形状反问题研究

随着CFD技术的发展,基于伴随方法的求解Euler和NS方程的气动优化设计已成为流体力学形状反问题研究中的热门领域.本文应用该方法对透平叶栅进行三维气动优化设计,详细推导了Euler方程伴随系统的偏微分方程组及其各类边界条件,首次给出了透平内流伴随方程边界条件的具体形式,并给出伴随变量的物理意义.结合拟牛顿算法发展了三维透平叶栅形状反问题气动优化算法,并给出了算法的流程.     

Second order gradient ascent pulse engineering

We report some improvements to the gradient ascent pulse engineering (GRAPE) algorithm for optimal control of spin ensembles and other quantum systems. These include more accurate gradients, convergence acceleration using the Broyden–Fletcher–Goldfarb–Shanno (BFGS) quasi-Newton algorithm as well as faster control derivative calculation algorithms. In all test systems, the wall clock time and the convergence rates show a considerable improvement over the approximate gradient ascent.     

基于非对称模型的变电所数据融合误差优化研究

牵引变电所作为电力系统中实现变换电能的场所,给电力机车提供稳定可靠的电能。由于机车具有非线性运载特点,在吸收电网电流的同时也在输入谐波电流,引起电网干扰波动和噪声误差、电压畸变等现象,对供电系统安全构成威胁。本文通过建立变电所三相匹配与失配模型,利用自适应数据融合算法仿真检测控制误差的效果和优化误差的最佳方案,合理控制参数的状态和优化参数的规律。结果表明数据融合可以解决局部电网误差波动电流超标和电压畸变,在保障供电质量和电网安全方面效果显著。     

Flux-corrected transport. III. Minimal-error FCT algorithms

This paper presents an error analysis of numerical algorithms for solving the convective continuity equation using flux-corrected transport (FCT) techniques. The nature of numerical errors in Eulerian finite-difference solutions to the continuity equation is analyzed. The properties and intrinsic errors of an “optimal” algorithm are discussed and a flux-corrected form of such an algorithm is demonstrated for a restricted class of problems. This optimal FCT algorithm is applied to a model test problem and the error is monitored for comparison with more generally applicable algorithms. Several improved FCT algorithms are developed and judged against both standard flux-uncorrected transport algorithms and the optimal algorithm. These improved FCT algorithms are found to be four to eight times more accurate than standard non-FCT algorithms, nearly twice as accurate as the original SHASTA FCT algorithm, and approach the accuracy of the optimal algorithm.     

A Meshless Regularization Method for a Two-Dimensional Two-Phase Linear Inverse Stefan Problem

In this paper, a meshless regularization method of fundamental solutions is proposed for a two-dimensional, two-phase linear inverse Stefan problem. The numerical implementation and analysis are challenging since one needs to handle composite materials in higher dimensions. Furthermore, the inverse Stefan problem is ill-posed since small errors in the input data cause large errors in the desired output solution. Therefore, regularization is necessary in order to obtain a stable solution. Numerical results for several benchmark test examples are presented and discussed.     

基于拟牛顿法的连铸结晶器热传导反问题求解

在实际的工业过程中,热传导反问题广泛的存在于各个应用领域。由于连铸的生产环境恶劣,给连铸钢坯的测温造成了很大的困难,导致连铸结晶器热传导反问题的求解是很困难的。在研究了热传导反演算法的基础上,基于连铸结晶器的有限的测量点,通过热传导反演算法可以求出结晶器的热流密度,给连铸的优化和控制提供依据。在建立了连铸热传导反问题模型的基础上,提出了基于拟牛顿法的连铸结晶器热传导反问题求解的改进算法,该算法可以有效的利用当前搜索过程中的信息,通过近似目标函数的二阶导数来加快算法的搜索效率,采用了Wolf步长选择规则克服了算法收敛速率慢的缺点,并且该算法无需计算Hessen矩阵,能够有效降低计算的复杂度。实验结果表明该算法能够快速有效的求解该问题。该方法用于连铸结晶器热传导反问题是可行的、有效的。     

Neuronal event detection in fMRI time series using iterative deconvolution techniques

An iterative estimation algorithm for deconvolution of neuronal activity from Blood Oxygen Level Dependent (BOLD) time series data is presented. The algorithm requires knowledge of the hemodynamic impulse response function but does not require knowledge of the stimulation function. The method uses majorization-minimization of a cost function to find an optimal solution to the inverse problem. The cost function includes penalties for the l₁ norm, total variation and negativity. The algorithm is able to identify the occurrence of neuronal activity bursts from BOLD time series accurately. The accuracy of the algorithm was tested in simulations and experimental fMRI data using blocked and event-related designs. The simulations revealed that the algorithm is most sensitive to contrast-to-noise ratio levels and to errors in the assumed hemodynamic model and least sensitive to autocorrelation in the noise. Within normal fMRI conditions, the method is effective for event detection.     

变分同化中水平误差函数的正交小波模拟新方法

背景误差协方差是变分资料同化系统中的一个重要组成部分,能将观测信息从观测点传播到周围的模式格点和垂直层上.为了模拟背景误差协方差中水平误差函数的非均匀性和各向异性,提出了一种用正交小波变换模拟水平误差函数的新方法.试验结果表明：新方法能模拟出水平误差函数中固有的非均匀性和各向异性,极好地表示了背景误差协方差中固有的结构和特征. 关键词： 变分资料同化 背景误差协方差 二维正交小波 水平误差函数     

Numerically accurate computational techniques for optimal estimator analyses of multi-parameter models

Modelling unclosed terms in partial differential equations typically involves two steps: First, a set of known quantities needs to be specified as input parameters for a model, and second, a specific functional form needs to be defined to model the unclosed terms by the input parameters. Both steps involve a certain modelling error, with the former known as the irreducible error and the latter referred to as the functional error. Typically, only the total modelling error, which is the sum of functional and irreducible error, is assessed, but the concept of the optimal estimator enables the separate analysis of the total and the irreducible errors, yielding a systematic modelling error decomposition. In this work, attention is paid to the techniques themselves required for the practical computation of irreducible errors. Typically, histograms are used for optimal estimator analyses, but this technique is found to add a non-negligible spurious contribution to the irreducible error if models with multiple input parameters are assessed. Thus, the error decomposition of an optimal estimator analysis becomes inaccurate, and misleading conclusions concerning modelling errors may be drawn. In this work, numerically accurate techniques for optimal estimator analyses are identified and a suitable evaluation of irreducible errors is presented. Four different computational techniques are considered: a histogram technique, artificial neural networks, multivariate adaptive regression splines, and an additive model based on a kernel method. For multiple input parameter models, only artificial neural networks and multivariate adaptive regression splines are found to yield satisfactorily accurate results. Beyond a certain number of input parameters, the assessment of models in an optimal estimator analysis even becomes practically infeasible if histograms are used. The optimal estimator analysis in this paper is applied to modelling the filtered soot intermittency in large eddy simulations using a dataset of a direct numerical simulation of a non-premixed sooting turbulent flame.     

Analytic Hessian derivation for the quasi-one-dimensional Euler equations

The Hessian for the quasi-one-dimensional Euler equations is derived. A pressure minimization problem and a pressure matching inverse problem are considered. The flow sensitivity, adjoint sensitivity, gradient and Hessian are calculated analytically using a direct approach that is specific to the model problems. For the pressure minimization problem we find that the Hessian exists and it contains elements with significantly larger values around the shock location. For the pressure matching inverse problem we find at least one case for which the gradient as well as the Hessian do not exist. In addition, two formulations for calculating the Hessian are proposed and implemented for the given problems. Both methods can be implemented in industrial applications such as large scale aerodynamic optimization.     

Trans-dimensional geoacoustic inversion

This paper develops a general trans-dimensional Bayesian methodology for geoacoustic inversion. Trans-dimensional inverse problems are a generalization of fixed-dimensional inversion that includes the number and type of model parameters as unknowns in the problem. By extending the inversion state space to multiple subspaces of different dimensions, the posterior probability density quantifies the state of knowledge regarding inversion parameters, including effects due to limited knowledge about appropriate parametrization of the environment and error processes. The inversion is implemented here using a reversible-jump Markov chain Monte Carlo algorithm and the seabed is parametrized using a partition model. Unknown data errors are addressed by including a data-error model. Jumps between dimensions are implemented with a birth-death methodology that allows transitions between dimensions by adding or removing interfaces while maintaining detailed balance in the Markov chain. Trans-dimensional inversion results in an inherently parsimonious solution while partition modeling provides a naturally self-regularizing algorithm based on data information content, not on subjective regularization functions. Together, this results in environmental estimates that quantify appropriate seabed structure as supported by the data, allowing sharp discontinuities while approximating smooth transitions where needed. This approach applies generally to geoacoustic inversion and is illustrated here with seabed reflection-coefficient data.