The subdiffusion equation with a Caputo fractional derivative of order alpha is an element of(0,1) in time arises in a wide variety of practical applications, and it is often adopted to model anomalous subdiffusion processes in heterogeneous media. The L1 scheme is one of the most popular and successful numerical methods for discretizing the Caputo fractional derivative in time. The scheme was analysed earlier independently by Lin and Xu (2007, Finite difference/spectral approximations for the time-fractional diffusion equation. J. Comput. Phys., 225, 1533-1552) and Sun and Wu (2006, A fully discrete scheme for a diffusion wave system. Appl. Numer. Math., 56, 193-209), and an O(tau(2) (alpha)) convergence rate was established, under the assumption that the solution is twice continuously differentiable in time. However, in view of the smoothing property of the subdiffusion equation, this regularity condition is restrictive, since it does not hold even for the homogeneous problem with a smooth initial data. In this work, we revisit the error analysis of the scheme, and establish an O(tau) convergence rate for both smooth and nonsmooth initial data. The analysis is valid for more general sectorial operators. In particular, the L1 scheme is applied to one-dimensional space-time fractional diffusion equations, which involves also a Riemann-Liouville derivative of order beta is an element of (3/2, 2) in space, and error estimates are provided for the fully discrete scheme. Numerical experiments are provided to verify the sharpness of the error estimates, and robustness of the scheme with respect to data regularity.

A symmetric discretization scheme for heterogeneous anisotropic diffusion problems on general meshes is developed and studied. The unknowns of this scheme are the values at the centre of the control volumes and at some internal interfaces that may, for instance, be chosen at the diffusion tensor discontinuities. The scheme is therefore completely cell centred if no edge unknown is kept. It is shown to be accurate for several numerical examples. The convergence of the approximate solution to the continuous solution is proved for general (possibly discontinuous) tensors and general (possibly nonconforming) meshes and with no regularity assumption on the solution. An error estimate is then deduced under suitable regularity assumptions on the solution.

We construct a preconditioned modified Hermitian and skew-Hermitian splitting (PMHSS) iteration scheme for solving and preconditioning a class of block two-by-two linear systems arising from the Galerkin finite element discretizations of a class of distributed control problems. The convergence theory of this class of PMHSS iteration methods is established and the spectral properties of the PMHSS-preconditioned matrix are analysed. Numerical experiments show that the PMHSS preconditioners can be quite competitive when used to precondition Krylov subspace iteration methods such as GMRES.

As long as a square non-negative matrix A has total support then it can be balanced, that is, we can find a diagonal scaling of A that has row and column sums equal to one. A number of algorithms have been proposed to achieve the balancing, the most well known of these being Sinkhorn-Knopp. In this paper, we derive new algorithms based on inner-outer iteration schemes. We show that Sinkhorn-Knopp belongs to this family, but other members can converge much more quickly. In particular, we show that while stationary iterative methods offer little or no improvement in many cases, a scheme using a preconditioned conjugate gradient method as the inner iteration converges at much lower cost (in terms of matrix-vector products) for a broad range of matrices; and succeeds in cases where the Sinkhorn-Knopp algorithm fails.

We present, in a unified framework, new conforming and nonconforming virtual element methods for general second-order elliptic problems in two and three dimensions. The differential operator is split into its symmetric and nonsymmetric parts and conditions for stability and accuracy on their discrete counterparts are established. These conditions are shown to lead to optimal H-1- and L-2-error estimates, confirmed by numerical experiments on a set of polygonal meshes. The accuracy of the numerical approximation provided by the two methods is shown to be comparable.

The matrix completion problem is to complete an unknown matrix from a small number of entries, and it captures many applications in diversified areas. Recently, it was shown that completing a low-rank matrix can be successfully accomplished by solving its convex relaxation problem using the nuclear norm. This paper shows that the alternating direction method (ADM) is applicable for completing a low-rank matrix including the noiseless case, the noisy case and the positive semidefinite case. The ADM approach for the matrix completion problem is easily implementable and very efficient. Numerical comparisons of the ADM approach with some state-of-the-art methods are reported.

The moving least squares (MLS) method provides an approximation u of a function u based solely on values u(x(j)) of u on scattered 'meshless' nodes x(j). Derivatives of u are usually approximated by derivatives of u. In contrast to this, we directly estimate derivatives of u from the data, without any detour via derivatives of u. This is a generalized MLS technique, and we prove that it produces diffuse derivatives as introduced by Nyroles et al. (1992, Generalizing the finite element method: diffuse approximation and diffuse elements. Comput. Mech., 10, 307-318). Consequently, these turn out to be efficient direct estimates of the true derivatives, without anything 'diffuse' about them, and we prove optimal rates of convergence towards the true derivatives. Numerical examples confirm this, and we finally show how the use of shifted and scaled polynomials as basis functions in the generalized and standard MLS approximation stabilizes the algorithm.

We consider numerical methods for thermodynamic sampling, i.e., computing sequences of points distributed according to the Gibbs-Boltzmann distribution, using Langevin dynamics and overdamped Langevin dynamics (Brownian dynamics). A wide variety of numerical methods for Langevin dynamics may be constructed based on splitting the stochastic differential equations into various component parts, each of which may be propagated exactly in the sense of distributions. Each such method may be viewed as generating samples according to an associated invariant measure that differs from the exact canonical invariant measure by a stepsize-dependent perturbation. We provide error estimates A la Talay-Tubaro on the invariant distribution for small stepsize, and compare the sampling bias obtained for various choices of the splitting method. We further investigate the overdamped limit and apply the methods in the context of driven systems where the goal is sampling with respect to a nonequilibrium steady state. Our analyses are illustrated by numerical experiments.

In this paper, we consider the initial-boundary value problem for the three-dimensional incompressible magnetohydrodynamic equations. We first establish certain regularity results for the solution (u, p, H) under the assumption that u, H is an element of L-4(0, T;H-1 (Omega)(3)). Next, we consider the fully discrete first-order Euler semi-implicit scheme based on the mixed finite element method. Then we obtain the L-2-unconditional convergence of this scheme by using the negative norm technique.

We consider the initial-boundary value problem for an inhomogeneous time-fractional diffusion equation with a homogeneous Dirichlet boundary condition, a vanishing initial data and a nonsmooth right-hand side in a bounded convex polyhedral domain. We analyse two semidiscrete schemes based on the standard Galerkin and lumped mass finite element methods. Almost optimal error estimates are obtained for right-hand side data f(x,t) is an element of L-infinity(0,T;(H)over bar(q)(Omega)), -1 < q <= 1, for both semidiscrete schemes. For the lumped mass method, the optimal L-2(Omega)-norm error estimate requires symmetric meshes. Finally, two-dimensional numerical experiments are presented to verify our theoretical results.

Abstract In this work, we analyse a Crank-Nicolson type time-stepping scheme for the subdiffusion equation, which involves a Caputo fractional derivative of order $\alpha\in (0,1)$ in time. It hybridizes the backward Euler convolution quadrature with a $\theta$-type method, with the parameter $\theta$ dependent on the fractional order $\alpha$ by $\theta=\alpha/2$ and naturally generalizes the classical Crank–Nicolson method. We develop essential initial corrections at the starting two steps for the Crank–Nicolson scheme and, together with the Galerkin finite element method in space, obtain a fully discrete scheme. The overall scheme is easy to implement and robust with respect to data regularity. A complete error analysis of the fully discrete scheme is provided, and a second-order accuracy in time is established for both smooth and nonsmooth problem data. Extensive numerical experiments are provided to illustrate its accuracy, efficiency and robustness, and a comparative study also indicates its competitive with existing schemes.

We develop and analyse a new family of virtual element methods on unstructured polygonal meshes for the diffusion problem in primal form, which uses arbitrarily regular discrete spaces V-h subset of C-alpha, alpha is an element of N. The degrees of freedom are (a) solution and derivative values of various degrees at suitable nodes and (b) solution moments inside polygons. The convergence of the method is proved theoretically and an optimal error estimate is derived. Numerical experiments confirm the convergence rate that is expected from the theory.

We establish a class of accelerated Hermitian and skew-Hermitian splitting (AHSS) iteration methods for large sparse saddle-point problems by making use of the Hermitian and skew-Hermitian splitting (HSS) iteration technique. These methods involve two iteration parameters whose special choices can recover the known preconditioned HSS iteration methods, as well as yield new ones. Theoretical analyses show that the new methods converge unconditionally to the unique solution of the saddle-point problem. Moreover, the optimal choices of the iteration parameters involved and the corresponding asymptotic convergence rates of the new methods are computed exactly. In addition, theoretical properties of the preconditioned Krylov subspace methods such as GMRES are investigated in detail when the AHSS iterations are employed as their preconditioners. Numerical experiments confirm the correctness of the theory and the effectiveness of the methods.

In this paper we introduce and analyze a virtual element method (VEM) for a mixed variational formulation of the Stokes problem in which the pseudostress and the velocity are the only unknowns, whereas the pressure is computed via a postprocessing formula. We first recall the corresponding continuous variational formulation, and then, following the basic principles for mixed-VEM, define the virtual finite element subspaces to be employed, introduce the associated interpolation operators, and provide the respective approximation properties. In particular, the latter includes the estimation of the interpolation error for the pseudostress variable measured in the H(div)-norm. We remark that a Bramble-Hilbert type theorem for averaged Taylor polynomials plays a key role in the respective analysis. Next, and in order to define calculable discrete bilinear forms, we propose a new local projector onto a suitable space of polynomials, which takes into account the main features of the continuous solution and allows the explicit integration of the terms involving the deviatoric tensors. The uniform boundedness of the resulting family of local projectors is established and, using the aforementioned compactness theorem, its approximation properties are also derived. In addition, we show that the global discrete bilinear forms satisfy all the hypotheses required by the Babu. ska-Brezzi theory. In this way, we conclude the well-posedness of the actual Galerkin scheme and derive the associated a priori error estimates for the virtual solution as well as for the fully computable projection of it. Finally, several numerical results illustrating the good performance of the method and confirming the theoretical rates of convergence are presented.

We investigate the discretization of Darcy flow through fractured porous media on general polyhedral meshes. We consider a hybrid-dimensional model, invoking a complex network of planar fractures. The model accounts for matrix-fracture interactions and fractures acting either as drains or as barriers, i.e., we have to deal with pressure discontinuities at matrix-fracture interfaces. The numerical analysis is performed in the general framework of gradient discretizations which is extended to the model under consideration. Two families of schemes, namely the vertex approximate gradient scheme and the hybrid finite volume scheme, are detailed and shown to fit the gradient scheme framework, which yields, in particular, convergence. Numerical tests confirm the theoretical results.

Abstract In this work, we present a mass conservative numerical scheme for two-phase flow in porous media. The model for flow consists of two fully coupled, nonlinear equations: a degenerate parabolic equation and an elliptic one. The proposed numerical scheme is based on backward Euler for the temporal discretization and mixed finite element method for the spatial one. A priori stability and error estimates are presented to prove the convergence of the scheme. A monotone increasing, Hölder continuous saturation is considered. The convergence of the scheme is naturally dependant on the Hölder exponent. The nonlinear systems within each time step are solved by a robust linearization method, called the $L$-scheme. This iterative method does not involve any regularization step. The convergence of the $L$-scheme is rigorously proved under the assumption of a Lipschitz continuous saturation. For the Hölder continuous case, a numerical convergence is established. Numerical results (two-dimensional and three-dimensional) are presented to sustain the theoretical findings.

In this paper, we devise and analyse an unconditionally stable, second-order-in-time numerical scheme for the Cahn-Hilliard equation in two and three space dimensions. We prove that our two-step scheme is unconditionally energy stable and unconditionally uniquely solvable. Furthermore, we show that the discrete phase variable is bounded in L-infinity (0,T;L-infinity) and the discrete chemical potential is bounded in L-infinity (0,T;L-2), for any time and space step sizes, in two and three dimensions, and for any finite final time T. We subsequently prove that these variables converge with optimal rates in the appropriate energy norms in both two and three dimensions. We include in this work a detailed analysis of the initialization of the two-step scheme.

In this paper, we will present a generalized convolution quadrature for solving linear parabolic and hyperbolic evolution equations. The original convolution quadrature method by Lubich works very nicely for equidistant time steps while the generalization of the method and its analysis to nonuniform time stepping is by no means obvious. We will introduce the generalized convolution quadrature allowing for variable time steps and develop a theory for its error analysis. This method opens the door for further development towards adaptive time stepping for evolution equations. As the main application of our new theory, we will consider the wave equation in exterior domains which is formulated as a retarded boundary integral equation.

In this paper, we propose a modified Polak-Ribiere-Polyak (PRP) conjugate gradient method. An attractive property of the proposed method is that the direction generated by the method is always a descent direction for the objective function. This property is independent of the line search used. Moreover, if exact line search is used, the method reduces to the ordinary PRP method. Under appropriate conditions, we show that the modified PRP method with Armijo-type line search is globally convergent. We also present extensive preliminary numerical experiments to show the efficiency of the proposed method.

A two-point boundary value problem whose highest order term is a Caputo fractional derivative of order delta is an element of (1, 2) is considered. Al-Refai's comparison principle is improved and modified to fit our problem. Sharp a priori bounds on derivatives of the solution u of the boundary value problem are established, showing that u '' (x) may be unbounded at the interval endpoint x = 0. These bounds and a discrete comparison principle are used to prove pointwise convergence of a finite difference method for the problem, where the convective term is discretized using simple upwinding to yield stability on coarse meshes for all values of d. Numerical results are presented to illustrate the performance of the method.