65-XX NUMERICAL ANALYSIS
Refine
Year of publication
Document Type
- ZIB-Report (35)
- Master's Thesis (5)
- Doctoral Thesis (2)
- Software (2)
Is part of the Bibliography
- no (44)
Keywords
- finite element method (5)
- optimal control (4)
- reduced basis method (4)
- rigorous optical modeling (3)
- Algebraic multigrid (2)
- Automatic Differentiation (2)
- Lipschitz Continuity (2)
- Monte Carlo (2)
- Nonsmooth (2)
- PDE (2)
Institute
- Numerical Mathematics (23)
- Computational Medicine (7)
- Mathematical Optimization (6)
- Modeling and Simulation of Complex Processes (6)
- Visual and Data-centric Computing (6)
- Computational Nano Optics (5)
- Visual Data Analysis (5)
- Computational Systems Biology (4)
- Energy Network Optimization (3)
- Mathematical Optimization Methods (3)
Quantitative photoacoustic tomography aims to recover maps of the local concentrations of tissue chromophores from multispectral images. While model-based inversion schemes are promising approaches, major challenges to their practical implementation include the unknown fluence distribution and the scale of the inverse problem. This paper describes an inversion scheme based on a radiance Monte Carlo model and an adjoint-assisted gradient optimization that incorporates fluence-dependent step sizes and adaptive moment estimation. The inversion is shown to recover absolute chromophore concentrations, blood oxygen saturation and the Grüneisen parameter from in silico 3D phantom images for different radiance approximations. The scattering coefficient was assumed to be homogeneous and known a priori.
In many applications, geodesic hierarchical models are adequate for the study of temporal observations. We employ such a model derived for manifold-valued data to Kendall's shape space. In particular, instead of the Sasaki metric, we adapt a functional-based metric, which increases the computational efficiency and does not require the implementation of the curvature tensor. We propose the corresponding variational time discretization of geodesics and apply the approach for the estimation of group trends and statistical testing of 3D shapes derived from an open access longitudinal imaging study on osteoarthritis.
In many applications, geodesic hierarchical models are adequate for the study of temporal observations. We employ such a model derived for manifold-valued data to Kendall's shape space. In particular, instead of the Sasaki metric, we adapt a functional-based metric, which increases the computational efficiency and does not require the implementation of the curvature tensor. We propose the corresponding variational time discretization of geodesics and employ the approach for longitudinal analysis of 2D rat skulls shapes as well as 3D shapes derived from an imaging study on osteoarthritis. Particularly, we perform hypothesis test and estimate the mean trends.
This thesis firstly presents a nonlinear extended deterministic model for the transmission dynamics of tuberculosis, based on realistic assumptions and data collected from the WHO. This model enables a comprehensive qualitative analysis of various aspects in the outbreak and control of tuberculosis in Sub-Saharan Africa countries and successfully reproduces the epidemiology of tuberculosis in Cameroon for the past (from 1994-2010). Some particular properties of the model and its solution have been presented using the comparison theorem applied to the theory of differential equations. The existence and the stability of a disease free equilibrium has been discussed using the Perron-Frobenius theorem and Metzler stable matrices.
Furthermore, we computed the basic reproduction number, i.e. the number of cases that one case generates on average over the course of its infectious period. Rigorous qualitative analysis of the model reveals that, in contrast to the model without reinfections, the full model with reinfection exhibits the phenomenon of backward bifurcation, where a stable disease-free equilibrium coexists with a stable endemic equilibrium when a certain threshold quantity, known as the basic reproduction ratio (R0), is less than unity. The global stability of the disease-free equilibrium has been discussed using the concepts of Lyapunov stability and bifurcation theory.
With the help of a sensitivity analysis using data of Cameroon, we identified the relevant parameters which play a key role for the transmission and the control of the disease. This was possible applying sophisticated numerical methods (POEM) developed at ZIB. Using advanced approaches for optimal control considering the costs for chemoprophylaxis, treatment and educational campaigns should provide a framework for designing realistic cost effective strategies with different intervention methods. The forward-backward sweep method has been used to solve the numerical optimal control problem. The numerical result of the optimal control problem reveals that combined effort in education and chemoprophylaxis may lead to a reduction of 80\% in the number of infected people in 10 years.
The mathematical and numerical approaches developed in this thesis could be similarly applied in many other Sub-Saharan countries where TB is a public health problem.
In this thesis, adaptive algorithms in optimization under PDE constraints have been inves-
tigated. In its application, the aim of optimization is to increase the longevity of implants,
namely the hip joint implant, and in doing so to minimize stress shielding and simultaneously minimize the influence of locally high stresses, that, above a threshold value, are malign to the bone structure. Under the constraint of the equilibrium of forces, describing an elastodynamic setup, coupled with a contact inequality condition, a computationally expensive problem formulation is given.
The first step to make the solution of the given problem possible and efficient was to change over to the spatial equilibrium equation, thus rendering an elastostatic setup. Subsequently the intrinsically dynamic motions – trajectories in the load domain – were converted to the static setup. Thus, the trajectories are marginalized to the load domain and characterized with probability distributions. Therefore the solving of the PDE constraint, the contact problem, is simplified.
Yet in the whole optimization process, the solving of the PDE, the spatial equilibrium equation together with the contact condition has the most expensive contribution still and hence needed further reduction. This was achieved by application of Kriging interpolation to the load responses of the integrated distribution of stress difference and the maximum stresses. The interpolation of the two response surfaces only needs comparatively few PDE solves to set up the models. Moreover, the Kriging models can be adaptively extended by sequentially adding sample-response pairs. For this the Kriging inherent variance is used to estimate ideal new sample locations with maximum variance values. In doing so, the overall interpolation variance and therefore the interpolation error is reduced.
For the integration of the integrated stress differences and penalty values on the relative high dimensional load domain Monte Carlo integration was implemented, averting the curse of dimension. Here, the motion’s probability distribution combined with patient specific data of motion frequencies is taken advantage of, making obsolete the use of the otherwise necessary importance sampling.
Throughout the optimization, the FE-discretization error and the subsequently attached errors entering the solution process via PDE discretization and approximative
solving of the PDE, Kriging interpolation and Monte Carlo integration need to decrease. While the FE-discretization error and the solution of the elastostatic contact problem were assumed precise enough, numerics showed, that the interpolation and integration errors can be controlled by adaptive refinement of the respective methods. For this purpose comparable error quantities for the particular algorithms were introduced and effectively put to use.
For the implant position’s optimization, the derivative of the objective function was derived using the implicit function theorem. As the FE-discretization changes with implant position modifications big enough, a special line search had to be used to deal with the discontinuities in the objective function.
The interplay and performance of the subalgorithms was demonstrated numerically on a reduced 2D setup of a hip joint with and without the implant. Consequently the load domain and the control variable were also limited to the 2D case.
Adaptive sampling strategies for efficient parameter scans in nano-photonic device simulations
(2014)
Rigorous optical simulations are an important tool in optimizing scattering properties of nano-photonic devices and are used, for example, in solar cell optimization. The finite element method (FEM) yields rigorous, time-harmonic, high accuracy solutions of the full 3D vectorial Maxwell's equations [1] and furthermore allows for great flexibility and accuracy in the geometrical modeling of these often complex shaped 3D nano-structures. A major drawback of frequency domain methods is the limitation of single frequency evaluations. For example the accurate computation of the short circuit current density of an amorphous silicon / micro-crystalline multi-junction thin film solar cell may require the solution of Maxwell's equations for over a hundred different wavelengths if an equidistant sampling strategy is employed. Also in optical metrology, wavelength scans are frequently used to reconstruct unknown geometrical and material properties of optical systems numerically from measured
scatterometric data.
In our contribution we present several adaptive numerical integration and sampling routines and study their efficiency in the context of the determination of generation rate profiles of solar cells. We show that these strategies lead to a reduction in the computational effort without loss of accuracy. We discuss the employment of tangential information in a Hermite interpolation scheme to achieve similar accuracy on coarser grids. We explore the usability of these strategies for scatterometry and solar cell simulations.
Our focus is on Maxwell's equations in the low frequency range; two specific applications we aim at are time-stepping schemes for eddy current computations and the stationary double-curl equation for time-harmonic fields. We assume that the computational domain is discretized by triangles or tetrahedrons; for the finite element approximation we choose N\'{e}d\'{e}lec's $H(curl)$-conforming edge elements of the lowest order. For the solution of the arising linear equation systems we devise an algebraic multigrid preconditioner based on a spatial component splitting of the field. Mesh coarsening takes place in an auxiliary subspace, which is constructed with the aid of a nodal vector basis. Within this subspace coarse grids are created by exploiting the matrix graphs. Additionally, we have to cope with the kernel of the $curl$-operator, which comprises a considerable part of the spectral modes on the grid. Fortunately, the kernel modes are accessible via a discrete Helmholtz decomposition of the fields; they are smoothed by additional algebraic multigrid cycles. Numerical experiments are included in order to assess the efficacy of the proposed algorithms.
An automatic adaptive importance sampling algorithm for molecular dynamics in reaction coordinates
(2017)
In this article we propose an adaptive importance sampling scheme for dynamical quantities of high dimensional complex systems which are metastable. The main idea of this article is to combine a method coming from Molecular Dynamics Simulation, Metadynamics, with a theorem from stochastic analysis, Girsanov's theorem. The proposed algorithm has two advantages compared to a standard estimator of dynamic quantities: firstly, it is possible to produce estimators with a lower variance and, secondly, we can speed up the sampling. One of the main problems for building importance sampling schemes for metastable systems is to find the metastable region in order to manipulate the potential accordingly. Our method circumvents this problem by using an assimilated version of the Metadynamics algorithm and thus creates a non-equilibrium dynamics which is used to sample the equilibrium quantities.
In this article we analyze a generalized trapezoidal rule for initial value problems with piecewise smooth right hand side F:IR^n -> IR^n.
When applied to such a problem, the classical trapezoidal rule suffers from a loss of accuracy if the solution trajectory intersects a non-differentiability of F. In such a situation the investigated generalized trapezoidal rule achieves a higher convergence order than the classical method. While the asymptotic behavior of the generalized method was investigated in a previous work, in the present article we develop the algorithmic structure for efficient implementation strategies
and estimate the actual computational cost of the latter.
Moreover, energy preservation of the generalized trapezoidal rule is proved for Hamiltonian systems with piecewise linear right hand side.
Convergence Analysis of Smoothing Methods for Optimal Control of Stationary Variational Inequalities
(2011)
In the article an optimal control problem subject to a stationary variational inequality
is investigated. The optimal control problem is complemented with pointwise control constraints.
The convergence of a smoothing scheme is analyzed. There, the variational inequality
is replaced by a semilinear elliptic equation. It is shown that solutions of the regularized optimal
control problem converge to solutions of the original one. Passing to the limit in the
optimality system of the regularized problem allows to prove C-stationarity of local solutions of the original problem.
Moreover, convergence rates with respect to the regularization parameter for the error in the control are obtained.
These rates coincide with rates obtained by numerical experiments, which are included in the paper.