Finite element simulations have been used to solve various partial differential equations (PDEs) that model physical, chemical, and biological phenomena. The resulting discretized solutions to PDEs often do not satisfy requisite physical properties, such as positivity or monotonicity. Such invalid solutions pose both modeling challenges, since the physical interpretation of simulation results is not possible, and computational challenges, since such properties may be required to advance the scheme. We, therefore, consider the problem of computing solutions that preserve these structural solution properties, which we enforce as additional constraints on the solution. We consider in particular the class of convex constraints, which includes positivity and monotonicity. By embedding such constraints as a postprocessing convex optimization procedure, we can compute solutions that satisfy general types of convex constraints. For certain types of constraints (including positivity and monotonicity), the optimization is a filter, i.e., a norm-decreasing operation. We provide a variety of tests on one-dimensional time-dependent PDEs that demonstrate the method's efficacy, and we empirically show that rates of convergence are unaffected by the inclusion of the constraints.
Various iterative reconstruction algorithms for inverse problems can be unfolded as neural networks. Empirically, this approach has often led to improved results, but theoretical guarantees are still scarce. While some progress on generalization properties of neural networks have been made, great challenges remain. In this chapter, we discuss and combine these topics to present a generalization error analysis for a class of neural networks suitable for sparse reconstruction from few linear measurements. The hypothesis class considered is inspired by the classical iterative soft-thresholding algorithm (ISTA). The neural networks in this class are obtained by unfolding iterations of ISTA and learning some of the weights. Based on training samples, we aim at learning the optimal network parameters via empirical risk minimization and thereby the optimal network that reconstructs signals from their compressive linear measurements. In particular, we may learn a sparsity basis that is shared by all of the iterations/layers and thereby obtain a new approach for dictionary learning. For this class of networks, we present a generalization bound, which is based on bounding the Rademacher complexity of hypothesis classes consisting of such deep networks via Dudley's integral. Remarkably, under realistic conditions, the generalization error scales only logarithmically in the number of layers, and at most linear in number of measurements.
This paper proposes a level set-based topology optimization method for designing acoustic structures with viscous and thermal boundary layers in perspective. It is known that acoustic waves propagating in a narrow channel are damped by viscous and thermal boundary layers. To estimate these viscothermal effects, we first introduce a sequential linearized Navier-Stokes model based on three weakly coupled Helmholtz equations for viscous, thermal, and acoustic pressure fields. Then, the optimization problem is formulated, where a sound-absorbing structure comprising air and an isothermal rigid medium is targeted, and its sound absorption coefficient is set as an objective function. The adjoint variable method and the concept of the topological derivative are used to obtain design sensitivity. A level set-based topology optimization method is used to solve the optimization problem. Two-dimensional numerical examples are provided to support the validity of the proposed method. In addition, the mechanisms that lead to the high absorption coefficient of the optimized design are discussed.
We provide a method to compute the entropy-satisfying weak solution to the eikonal equation in an arbitrary-order polynomial space. The method uses an artificial viscosity approach and is demonstrated for the signed distance function, where exact solutions are available. The method is designed specifically for an existing high-order discontinuous-Galerkin framework, which uses standard convection, diffusion, and source terms. We show design order of accuracy and good behavior for both shocks and rarefaction type solutions. Finally the distance function around a complex multi-element airfoil is computed using a high-order-accurate representation.
We present a discontinuous Galerkin internal-penalty scheme that is applicable to a large class of linear and non-linear elliptic partial differential equations. The scheme constitutes the foundation of the elliptic solver for the SpECTRE numerical relativity code. As such it can accommodate (but is not limited to) elliptic problems in linear elasticity, general relativity and hydrodynamics, including problems formulated on a curved manifold. We provide practical instructions that make the scheme functional in a production code, such as instructions for imposing a range of boundary conditions, for implementing the scheme on curved and non-conforming meshes and for ensuring the scheme is compact and symmetric so it may be solved more efficiently. We report on the accuracy of the scheme for a suite of numerical test problems.
In this paper, we propose and analyze a temporally second-order accurate, fully discrete finite element method for the magnetohydrodynamic (MHD) equations. A modified Crank--Nicolson method is used to discretize the model and appropriate semi-implicit treatments are applied to the fluid convection term and two coupling terms. These semi-implicit approximations result in a linear system with variable coefficients for which the unique solvability can be proved theoretically. In addition, we use a decoupling projection method of the Van Kan type \cite{vankan1986} in the Stokes solver, which computes the intermediate velocity field based on the gradient of the pressure from the previous time level, and enforces the incompressibility constraint via the Helmholtz decomposition of the intermediate velocity field. The energy stability of the scheme is theoretically proved, in which the decoupled Stokes solver needs to be analyzed in details. Optimal-order convergence of $\mathcal{O} (\tau^2+h^{r+1})$ in the discrete $L^\infty(0,T;L^2)$ norm is proved for the proposed decoupled projection finite element scheme, where $\tau$ and $h$ are the time stepsize and spatial mesh size, respectively, and $r$ is the degree of the finite elements. Existing error estimates of second-order projection methods of the Van Kan type \cite{vankan1986} were only established in the discrete $L^2(0,T;L^2)$ norm for the Navier--Stokes equations. Numerical examples are provided to illustrate the theoretical results.
We consider the problem of approximating a function in general nonlinear subsets of $L^2$ when only a weighted Monte Carlo estimate of the $L^2$-norm can be computed. Of particular interest in this setting is the concept of sample complexity, the number of samples that are necessary to recover the best approximation. Bounds for this quantity have been derived in a previous work and depend primarily on the model class and are not influenced positively by the regularity of the sought function. This result however is only a worst-case bound and is not able to explain the remarkable performance of iterative hard thresholding algorithms that is observed in practice. We reexamine the results of the previous paper and derive a new bound that is able to utilize the regularity of the sought function. A critical analysis of our results allows us to derive a sample efficient algorithm for the model set of low-rank tensors. The viability of this algorithm is demonstrated by recovering quantities of interest for a classical high-dimensional random partial differential equation.
Combining three themes: port-Hamiltonian energy-based modelling, structural analysis as used in the circuit world, and structural analysis of general differential-algebraic equations, we form a new model for electrical circuits, the compact port-Hamiltonian equations. They have remarkable simplicity and symmetry, and always have index at most 1 and other good numerical properties. The method has been implemented in Matlab. We give proofs and numerical results.
We study constrained reinforcement learning (CRL) from a novel perspective by setting constraints directly on state density functions, rather than the value functions considered by previous works. State density has a clear physical and mathematical interpretation, and is able to express a wide variety of constraints such as resource limits and safety requirements. Density constraints can also avoid the time-consuming process of designing and tuning cost functions required by value function-based constraints to encode system specifications. We leverage the duality between density functions and Q functions to develop an effective algorithm to solve the density constrained RL problem optimally and the constrains are guaranteed to be satisfied. We prove that the proposed algorithm converges to a near-optimal solution with a bounded error even when the policy update is imperfect. We use a set of comprehensive experiments to demonstrate the advantages of our approach over state-of-the-art CRL methods, with a wide range of density constrained tasks as well as standard CRL benchmarks such as Safety-Gym.
Graph signals are signals with an irregular structure that can be described by a graph. Graph neural networks (GNNs) are information processing architectures tailored to these graph signals and made of stacked layers that compose graph convolutional filters with nonlinear activation functions. Graph convolutions endow GNNs with invariance to permutations of the graph nodes' labels. In this paper, we consider the design of trainable nonlinear activation functions that take into consideration the structure of the graph. This is accomplished by using graph median filters and graph max filters, which mimic linear graph convolutions and are shown to retain the permutation invariance of GNNs. We also discuss modifications to the backpropagation algorithm necessary to train local activation functions. The advantages of localized activation function architectures are demonstrated in four numerical experiments: source localization on synthetic graphs, authorship attribution of 19th century novels, movie recommender systems and scientific article classification. In all cases, localized activation functions are shown to improve model capacity.
Graph neural networks (GNNs) are a popular class of machine learning models whose major advantage is their ability to incorporate a sparse and discrete dependency structure between data points. Unfortunately, GNNs can only be used when such a graph-structure is available. In practice, however, real-world graphs are often noisy and incomplete or might not be available at all. With this work, we propose to jointly learn the graph structure and the parameters of graph convolutional networks (GCNs) by approximately solving a bilevel program that learns a discrete probability distribution on the edges of the graph. This allows one to apply GCNs not only in scenarios where the given graph is incomplete or corrupted but also in those where a graph is not available. We conduct a series of experiments that analyze the behavior of the proposed method and demonstrate that it outperforms related methods by a significant margin.