We describe a method for discretizing planar C2-regular domains immersed in
non-conforming triangulations. The method consists in constructing mappings
from triangles in a background mesh to curvilinear ones that conform exactly to
the immersed domain. Constructing such a map relies on a novel way of
parameterizing the immersed boundary over a collection of nearby edges with its
closest point projection. By interpolating the mappings to curvilinear
triangles at select points, we recover isoparametric mappings for the immersed
domain defined over the background mesh.

In this paper we apply the recently developed mimetic discretization method
[44] to the mixed formulation of the Stokes problem in terms of the
vorticity-velocity-pressure formulation. The mimetic discretization presented
in this paper and in [44] is a higher-order method for curvilinear
quadrilaterals. It relies on the language of differential $k$-forms,
$k$-cochains as its discrete counterpart, and the relations between them in
terms of the mimetic operators: reduction, reconstruction and projection. The
reconstruction consists of a mimetic spectral element method.

The aim of this paper is to discuss some applications of general topology in
computer algorithms including modeling and simulation, and also in computer
graphics and image processing. While the progress in these areas heavily
depends on advances in computing hardware, the major intellectual achievements
are the algorithms. The applications of general topology in other branches of
mathematics are not discussed, since they are not applications of mathematics
outside of mathematics.

We present how entropy estimates and logarithmic Sobolev inequalities on the
one hand, and the notion of quasi-stationary distribution on the other hand,
are useful tools to analyze metastable overdamped Langevin dynamics, in
particular to quantify the degree of metastability. We discuss the interest of
these approaches to estimate the efficiency of some classical algorithms used
to speed up the sampling, and to evaluate the error introduced by some
coarse-graining procedures. This paper is a summary of a plenary talk given by
the author at the ENUMATH 2011 conference.

This paper introduces an interpolation framework for the weighted-H2 model
reduction problem. We obtain a new representation of the weighted-H2 norm of
SISO systems that provides new interpolatory first order necessary conditions
for an optimal reduced-order model. The H2 norm representation also provides an
error expression that motivates a new weighted-H2 model reduction algorithm.
Several numerical examples illustrate the effectiveness of the proposed
approach.

A well-known problem in computing some matrix functions iteratively is the
lack of a clear, commonly accepted residual notion. An important matrix
function for which this is the case is the matrix exponential. Suppose the
matrix exponential of a given matrix times a given vector has to be computed.
We develop the approach of Druskin, Greenbaum and Knizhnerman (1998) and
interpret the sought-after vector as the value of a vector function satisfying
the linear system of ordinary differential equations (ODE) whose coefficients
form the given matrix.

We present a posteriori error estimates for a recently developed
atomistic/continuum coupling method, the Consistent Energy-Based QC Coupling
method. The error estimate of the deformation gradient combines a residual
estimate and an a posteriori stability analysis. The residual is decomposed
into the residual due to the approximation of the stored energy and that due to
the approximation of the external force, and are bounded in negative Sobolev
norms. In addition, the error estimate of the total energy using the error
estimate of the deformation gradient is also presented.

We prove stability estimates for the ENO reconstruction and ENO interpolation
procedures. In particular, we show that the jump of the reconstructed ENO
pointvalues at each cell interface has the same sign as the jump of the
underlying cell averages across that interface. We also prove that the jump of
the reconstructed values can be upper-bounded in terms of the jump of the
underlying cell averages. Similar sign properties hold for the ENO
interpolation procedure.

In this paper we study the phenomenon of nonlinear supratransmission in a
semi-infinite discrete chain of coupled oscillators described by modified
sine-Gordon equations with constant external and internal damping, and subject
to harmonic external driving at the end.

Development of scientific software involves tradeoffs between ease of use,
generality, and performance. We describe the design of a general hyperbolic PDE
solver that can be operated with the convenience of MATLAB yet achieves
efficiency near that of hand-coded Fortran and scales to the largest
supercomputers. This is achieved by using Python for most of the code while
employing automatically-wrapped Fortran kernels for computationally intensive
routines, and using Python bindings to interface with a parallel computing
library and other numerical packages.

We constraint on computer the best linear unbiased generalized statistics of
random field for the best linear unbiased generalized statistics of an unknown
constant mean of random field and derive the numerical generalized
least-squares estimator of an unknown constant mean of random field. We derive
the third constraint of spatial statistics and show that the classic
generalized least-squares estimator of an unknown constant mean of the field is
only an asymptotic disjunction of the numerical one.

In this paper we present a multigrid approach to solve the Poisson equation
in arbitrary domain (identified by a level set function) and mixed boundary
conditions. The discretization is based on finite difference scheme and
ghost-cell method. This multigrid strategy can be applied also to more general
problems where a non-eliminated boundary condition approach is used. Arbitrary
domain make the definition of the restriction operator for boundary conditions
hard to find.

In this article I present a fast and direct method for solving several types
of linear finite difference equations (FDE) with constant coefficients. The
method is based on a polynomial form of the translation operator and its
inverse, and can be used to find the particular solution of the FDE. This work
raises the possibility of developing new ways to expand the scope of the
operational methods.

This paper presents theory for the Lagrange co-rotational (CR) formulation of
finite elements in the geometrically nonlinear analysis of 3D structures. In
this paper strains are assumed to be small while the magnitude of rotations
from the reference configuration is not restricted. A new best fit rotator and
consistent spin filter are derived. Lagrange CR formulation is applied with
Hybrid Trefftz Stress elements, although presented methodology can be applied
to arbitrary problem formulation and discretization technique, f.e.

We consider Implicit-Explicit (IMEX) Runge-Kutta schemes for hyperbolic and
kinetic equations in the diffusion limit. In such regime the system relaxes
towards a parabolic convection-diffusion equation and it is desirable to have a
method that is able to capture the asymptotic behavior with an implicit
treatment of the limiting diffusive terms. To this goal we reformulate the
problem by properly combining the limiting diffusion flux with the convective
flux. This, however, introduces new difficulties due to the dependence of the
stiff source term on the gradient.

RiemCirc is a C++ program which allocates points inside the unit circle for
numerical quadrature on the circle, aiming at homogeneous equidistant
distribution. The weights of the quadrature rule are computed by the area of
the tiles that surround these nodes. The shapes of the areas are polygonal,
defined by Voronoi tessellation.

The condition number of a diagonally scaled matrix, for appropriately chosen
scaling matrices, is often less than that of the original. Approximate
equilibration scales a matrix so that the scaled matrix's row and column norms
are approximately equal. We develop approximate equilibration algorithms for
both nonsymmetric and symmetric matrices that access a matrix only by
matrix-vector products, and we show that approximate equilibration is possible
for all structurally nonsingular matrices.

We seek to approximate a composite function h(x) = g(f(x)) with a global
polynomial. The standard approach chooses points in the domain of f and
computes h(x) at each point, which requires an evaluation of f and an
evaluation of g. We present a Lanczos-based procedure that implicitly
approximates g with a polynomial of f. By constructing a quadrature rule for
the density function of f, we can approximate h(x) using many fewer evaluations
of g. The savings is particularly dramatic when g is much more expensive than f
or the dimension of x is large.

Differential optical absorption spectroscopy (DOAS) is a powerful tool for
detecting and quantifying trace gases in atmospheric chemistry
\cite{Platt_Stutz08}. DOAS spectra consist of a linear combination of complex
multi-peak multi-scale structures. Most DOAS analysis routines in use today are
based on least squares techniques, for example, the approach developed in the
1970s uses polynomial fits to remove a slowly varying background, and known
reference spectra to retrieve the identity and concentrations of reference
gases.

Pseudospectral collocation methods and finite difference methods have been
used for approximating an important family of soliton like solutions of the
mKdV equation. These solutions present a structural instability which make
difficult to approximate their evolution in long time intervals with enough
accuracy. The standard numerical methods do not guarantee the convergence to
the proper solution of the initial value problem and often fail by approaching
solutions associated to different initial conditions.

Finite difference approximations to multi-asset American put option price are
considered. The assets are modelled as a multi-dimensional diffusion process
with variable drift and volatility. Approximation error of order one quarter
with respect to the time discretisation parameter and one half with respect to
the space discretisation parameter is proved by reformulating the corresponding
optimal stopping problem as a solution of a degenerate Hamilton-Jacobi-Bellman
equation.

The aim of this work is to study the asymptotic behavior of a structure made
of plates of thickness $2\delta$ when $\delta\to 0$. This study is carried on
within the frame of linear elasticity by using the unfolding method. It is
based on several decompositions of the structure displacements and on the
passing to the limit in fixed domains. We begin with studying the displacements
of a plate.

In this paper we study the asymptotic behavior of a structure made of curved
rods of thickness 2\delta when \delta rightarrow 0. This study is carried on
within the frame of linear elasticity by using the unfolding method. It is
based on several decompositions of the structure displacements and on the
passing to the limit in fixed domains. We show that any displacement of a
structure is the sum of an elementary rods-structure displacement (e.r.s.d.)
concerning the rods cross sections and a residual one related to the
deformation of the cross-section. The e.r.s.d.

Much work has gone into the construction of quasicontinuum energies that
reduce the coupling error along the interface between atomistic and continuum
regions. The largest consistency errors are typically pointwise
$O(\frac{1}{\eps})$ errors, and in some cases this has been reduced to
pointwise O(1) errors. In this paper we show that one cannot create a coupling
method using a finite-range coupling interface that has o(1)-consistency in the
interface, and we use this to give an upper bound on the order of convergence
in discrete $w^{1,p}$-norms in 1D.

We propose a Multiscale BDDC for a class of saddle-point problems. The method
solves for both flux and pressure variables. The fluxes are resolved in
three-steps: the coarse solve is followed by subdomain solves, and last we look
for a divergence-free flux correction and pressure variables using conjugate
gradients with a Multilevel BDDC preconditioner.

With the current hybridization of treecodes and FMMs, combined with
auto-tuning capabilities on heterogeneous architectures, the flexibility of
fast N -body methods has been greatly enhanced. These features are a
requirement to developing a black- box software library for fast N-body
algorithms on heterogeneous systems, which is our immediate goal.

We investigate the long tim behavior of the following efficient second order
in time scheme for the 2D Navier-Stokes equation in a periodic box: $$
\frac{3\omega^{n+1}-4\omega^n+\omega^{n-1}}{2k} +
\nabla^\perp(2\psi^n-\psi^{n-1})\cdot\nabla(2\omega^n-\omega^{n-1}) -
\nu\Delta\omega^{n+1} = f^{n+1}, \quad -\Delta \psi^n = \om^n. $$ The scheme is
a combination of a 2nd order in time backward-differentiation (BDF) and a
special explicit Adams-Bashforth treatment of the advection term. Therefore
only a linear constant coefficient Poisson type problem needs to be solved at
each time step.

A well known method for convergence acceleration of continued fraction
$\K(a_n/b_n)$ is to use the modified approximants $S_n(\omega_n)$ in place of
the classical approximants $S_n(0)$, where $\omega_n$ are close to tails
$f^{(n)}$ of continued fraction. Recently, author proposed a method of
iterative character producing tail approximations whose asymptotic expansion's
accuracy is improving in each step. This method can be applied to continued
fractions $\K(a_n/b_n)$, where $a_n$, $b_n$ are polynomials in $n$ ($\deg
a_n=2$, $\deg b_n\leq 1$) for sufficiently large $n$.

We consider instability of a two layered solid body of a denser material on
top of a lighter one. This problem is widely known to geoscientist in
sediment-salt migration as salt diapirism. In the literature, this problem has
often been treated as Raleigh-Taylor instability in viscous fluids instead of
solid bodies. In this paper, we propose a successive linear approximation
method for large deformation in viscoelastic solids as a model for salt
migration.

Piezoelectric devices, such as piezoelectric traveling wave rotary ultrasonic
motors, have composite piezoelectric structures. A composite piezoelectric
structure consists of a combination of two or more bonded materials, where at
least one of them is a piezoelectric transducer. Numerical modeling of
piezoelectric structures has been done in the past mainly with the finite
element method.

The main aim of this paper is to document the performance of $p$-refinement
with respect to maximum principles and the non-negative constraint. The model
problem is (steady-state) anisotropic diffusion with decay (which is a
second-order elliptic partial differential equation). We considered the
standard single-field formulation (which is based on the Galerkin formalism)
and two least-squares-based mixed formulations. We have employed non-uniform
Lagrange polynomials for altering the polynomial order in each element, and we
have used $p = 1, ..., 10$.

Partition of unity methods, such as the extended finite element method (XFEM)
allow discontinuities to be simulated independently of the mesh [1]. This
eliminates the need for the mesh to be aligned with the discontinuity or
cumbersome re-meshing, as the discontinuity evolves. However, to compute the
stiffness matrix of the elements intersected by the discontinuity, a
subdivision of the elements into quadrature subcells aligned with the
discontinuity is commonly adopted.

In this paper, the linear free flexural vibration of cracked functionally
graded material plates is studied using the extended finite element method. A
4-noded quadrilateral plate bending element based on field and edge consistency
requirement with 20 degrees of freedom per element is used for this study. The
natural frequencies and mode shapes of simply supported and clamped square and
rectangular plates are computed as a function of gradient index, crack length,
crack orientation and crack location. The effect of thickness and influence of
multiple cracks is also studied.

In this paper we propose some very promissing results in interval arithmetics
which permit to build well-defined arithmetics including distributivity of
multiplication and division according addition and substraction. Thus, it
allows to build all algebraic operations and functions on intervals. This will
avoid completely the wrapping effects and data dependance. Some simple
applications for matrix eigenvalues calculations, inversion of symmetric
matrices and finally optimization are exhibited in the object-oriented
programming language python.

It is described how the Hermite-Pad\'e polynomials corresponding to an
algebraic approximant for a power series may be used to predict coefficients of
the power series that have not been used to compute the Hermite-Pad\'e
polynomials. A recursive algorithm is derived and some numerical examples are
given.

A fuzzy linear system with crisp coefficient matrix and with an arbitrary
fuzzy number in parametric form on the right-hand side is investigated. It is
shown that the well-known existence and uniqueness theorem of a strong fuzzy
solution is equivalent to the following: The coefficient matrix is the product
of a permutation matrix and a diagonal matrix. This means that this theorem can
be applicable only for a special form of linear systems, namely, only when the
system consists of equations, each of which has exactly one variable.

Plotting solution sets for particular equations may be complicated by the
existence of turning points. Here we describe an algorithm which not only
overcomes such problematic points, but does so in the most general of settings.
Applications of the algorithm are highlighted through two examples: the first
provides verification, while the second demonstrates a non-trivial application.
The latter is followed by a thorough run-time analysis. While both examples
deal with bivariate equations, it is discussed how the algorithm may be
generalized for space curves in $\R^{3}$.

Consider being given a mapping \phi from the unit sphere S^{d-1}, d>2, to the
smooth boundary of a simply-connected region \Omega in R^d. We consider the
problem of constructing an extension \Phi from the unit ball B_d to \Omega. The
mapping is required to be 1-1 and continuously differentiable with a
nonsingular Jacobian matrix. We discuss ways of obtaining initial guesses for
such a mapping \Phi and of then improving it by an iteration method.

Achieving computing at the exascale means accelerating today's applications
by one thousand times. Clearly, this cannot be accomplished by hardware alone,
at least not in the short time frame expected for reaching this performance
milestone. Thus, a lively discussion has begun in the last couple of years
about programming models, software components and tools, and algorithms that
will facilitate exascale computing. Among the algorithms that are likely to
play a preeminent role in the new world of computing, the fast multipole method
(F MM) appears as a rising star.

We present a novel way of constructing reduced models for systems of ordinary
differential equations. The reduced models we construct depend on coefficients
which measure the importance of the different terms appearing in the model and
need to be estimated. The proposed approach allows the estimation of these
coefficients on the fly by enforcing the equality of integral quantities of the
solution as computed from the original system and the reduced model.

The classical EMD algorithm has been used extensively in the literature to
decompose signals that contain nonlinear waves. However when a signal contain
two or more frequencies that are close to one another the decomposition might
fail. In this paper we propose a new formulation of this algorithm which is
based on the zero crossings of the signal and show that it performs well even
when the classical algorithm fail. We address also the filtering properties and
convergence rate of the new algorithm versus the classical EMD algorithm.

Total-variation (TV)-based Computed Tomography (CT) image reconstruction has
shown experimentally to be capable of producing accurate reconstructions from
sparse-view data. In particular TV-based reconstruction is very well suited for
images with piecewise nearly constant regions. Computationally, however,
TV-based reconstruction is much more demanding, especially for 3D imaging, and
the reconstruction from clinical data sets is far from being close to
real-time.

In this paper, we develop a new integral representation for the solution of
the time harmonic Maxwell equations in media with piecewise constant dielectric
permittivity and magnetic permeability in R^3. This representation leads to a
coupled system of Fredholm integral equations of the second kind for four
scalar densities supported on the material interface. Like the classical Muller
equation, it has no spurious resonances. Unlike the classical approach,
however, the representation does not suffer from low frequency breakdown.

Approximation of the marginal distribution of the solution of the stochastic
Navier-Stokes equations on the two-dimensional torus by high order numerical
methods is considered. The corresponding rates of convergence are obtained for
a splitting scheme and the method of cubature on Wiener space applied to a
spectral Galerkin discretisation of degree $N$. While the estimates exhibit a
strong $N$ dependence, convergence is obtained for appropriately chosen time
step sizes. Results of numerical simulations are provided, and confirm the
applicability of the methods.

The new concept of numerical smoothness is applied to RKDG methods on the
scalar nonlinear conservation laws. The main result is an a posteriori error
estimate for the RKDG methods of arbitrary order in space and time, with
optimal convergence rate. In this paper, the case of smooth solutions is the
focus point. However, the error analysis framework is prepared to deal with
discontinuous solutions in the future.

We study the half-plane problem for the elastic wave equation subject to a
free surface boundary condition, with particular emphasis on almost
incompressible materials. A normal mode analysis is developed to estimate the
solution in terms of the boundary data, showing that the problem is boundary
stable. The dependence on the material properties, which is difficult to
analyze by the energy method, is made transparent by our estimates. The normal
mode technique is used to analyze the influence of truncation errors in a
finite difference approximation.

Many structured data-fitting applications require the solution of an
optimization problem involving a sum over a potentially large number of
measurements. Incremental gradient algorithms (both deterministic and
randomized) offer inexpensive iterations by sampling only subsets of the terms
in the sum. These methods can make great progress initially, but often slow as
they approach a solution. In contrast, full gradient methods achieve steady
convergence at the expense of evaluating the full objective and gradient on
each iteration.

Given a C^1 path of systems of homogeneous polynomial equations f_t, t in
[a,b] and an approximation x_a to a zero zeta_a of the initial system f_a, we
show how to adaptively choose the step size for a Newton based homotopy method
so that we approximate the lifted path (f_t,zeta_t) in the space of (problems,
solutions) pairs.

The total number of Newton iterations is bounded in terms of the length of
the lifted path in the condition metric.

Low-rank matrix recovery addresses the problem of recovering an unknown
low-rank matrix from few linear measurements. Nuclear-norm minimization is a
tractible approach with a recent surge of strong theoretical backing. Analagous
to the theory of compressed sensing, these results have required random
measurements. For example, m >= Cnr Gaussian measurements are sufficient to
recover any rank-r n x n matrix with high probability. In this paper we address
the theoretical question of how many measurements are needed via any method
whatsoever --- tractible or not.

There are very few results on mixed finite element methods on surfaces. A
theory for the study of such methods was given recently by Holst and Stern,
using a variational crimes framework in the context of finite element exterior
calculus. However, we are not aware of any numerical experiments where mixed
finite elements derived from discretizations of exterior calculus are used for
a surface domain. This short note shows results of our preliminary experiments
using mixed methods for Darcy flow (hence scalar Poisson's equation in mixed
form) on surfaces.

In this paper we introduce a common problem in electronic measurements and
electrical engineering: finding the first root from the left of an equation in
the presence of some initial conditions. We present examples of
electrotechnical devices (analog signal filtering), where it is necessary to
solve it. Two new methods for solving this problem, based on global
optimization ideas, are introduced. The first uses the exact a priori given
global Lipschitz constant for the first derivative. The second method
adaptively estimates local Lipschitz constants during the search.

Much of uncertainty quantification to date has focused on determining the
effect of variables modeled probabilistically, and with a known distribution,
on some physical or engineering system. We develop methods to obtain
information on the system when the distributions of some variables are known
exactly, others are known only approximately, and perhaps others are not
modeled as random variables at all.

In this paper, we focus on two classes of D-invariant polynomial subspaces.
The first is a classical type, while the second is a new class. With matrix
computation, we prove that every ideal projector with each D-invariant subspace
belonging to either the first class or the second is the pointwise limit of
Lagrange projectors. This verifies a particular case of a C. de Boor's
conjecture asserting that every complex ideal projector is the pointwise limit
of Lagrange projectors. Specifically, we provide the concrete perturbation
procedure for ideal projectors of this type.

We present a novel, cell-local shock detector for use with discontinuous
Galerkin (DG) methods. The output of this detector is a reliably scaled,
element-wise smoothness estimate which is suited as a control input to a shock
capture mechanism. Using an artificial viscosity in the latter role, we obtain
a DG scheme for the numerical solution of nonlinear systems of conservation
laws. Building on work by Persson and Peraire, we thoroughly justify the
detector's design and analyze its performance on a number of benchmark
problems.

A self-learning algebraic multigrid method for dominant and minimal singular
triplets and eigenpairs is described. The method consists of two multilevel
phases. In the first, multiplicative phase (setup phase), tentative singular
triplets are calculated along with a multigrid hierarchy of interpolation
operators that approximately fit the tentative singular vectors in a collective
and self-learning manner, using multiplicative update formulas.

We give a complete characterization of the behavior of the Anderson
acceleration (with arbitrary nonzero mixing parameters) on linear problems. Let
n be the grade of the residual at the starting point with respect to the matrix
defining the linear problem. We show that if Anderson acceleration does not
stagnate (that is, produces different iterates) up to n, then the sequence of
its iterates converges to the exact solution of the linear problem. Otherwise,
the Anderson acceleration converges to the wrong solution.

The problem of optimal mass transport arises in numerous applications
including image registration, mesh generation, reflector design, and
astrophysics. One approach to solving this problem is via the Monge-Amp\`ere
equation. While recent years have seen much work in the development of
numerical methods for solving this equation, very little has been done on the
implementation of the transport boundary conditions. In this paper, we propose
a method for solving the transport problem by iteratively solving a
Monge-Amp\`ere equation with Neumann boundary conditions.

A discretization scheme for variable coefficient elliptic PDEs in the plane
is presented. The scheme is based on high-order Gaussian quadratures and is
designed for problems with smooth solutions, such as scattering problems
involving soft scatterers. The resulting system of linear equations is very
well suited to efficient direct solvers such as nested dissection and the more
recently proposed accelerated nested dissection schemes with O(N) complexity.

In this paper we present a numerical scheme for the resolution of matrix
Riccati equation, usualy used in control problems. The scheme is
unconditionnaly stable and the solution is definite positive at each time step
of the resolution. We prove the convergence in the scalar case and present
several numerical experiments for classical test cases.

We introduce a novel technique for constructing higher-order variational
integrators for Hamiltonian systems of ODEs. In particular, we are concerned
with generating globally smooth approximations to solutions of a Hamiltonian
system. Our construction of the discrete Lagrangian adopts Hermite
interpolation polynomials and the Euler-Maclaurin quadrature formula, and
involves applying collocation to the Euler-Lagrange equation and its
prolongation.

We give a new, simple, dimension-independent definition of the serendipity
finite element family. The shape functions are the span of all monomials which
are linear in at least s-r of the variables where s is the degree of the
monomial or, equivalently, whose superlinear degree (total degree with respect
to variables entering at least quadratically) is at most r. The degrees of
freedom are given by moments of degree at most r-2d on each face of dimension
d. We establish unisolvence and a geometric decomposition of the space.

Mesh adaptation for finite element approximation is a procedure used in
numerous applications. The use of thin and long anisotropic triangles improves
the efficiency of the procedure. When piecewise linear finite elements are
used, the aspect ratio for mesh adaptation is generally dictated by the
absolute value of the (estimated) hessian matrix of the approximated function.
We give in this paper the corresponding aspect ratio for piecewise quadratic
finite elements.

The Computation of discrete Contractive semigroups becomes necessary when we
deal with several types of evolution equations in Discretizable Hilbert spaces,
in this work we study some properties of the discrete forms of the contractive
semigroups induced by an approximation scheme in a prescribed Hilbert space, we
also deal with the implementation of computational methods in this Hilbert
Space and apply some of the results presented here in the Heisenberg
representation of quantum dynamical semigroups.

For an elliptic problem with two space dimensions, we propose to formulate
the finite volume method with the help of Petrov-Galerkin mixed finite
elementsthat are based on the building of a dual Raviart-Thomas basis.

We introduce a generalized framework for sampling and reconstruction in
separable Hilbert spaces. Specifically, we establish that it is always possible
to stably reconstruct a vector in an arbitrary Riesz basis from sufficiently
many of its samples in any other Riesz basis. This framework can be viewed as
an extension of that of Eldar et al. However, whilst the latter imposes
stringent assumptions on the reconstruction basis, and may in practice be
unstable, our framework allows for recovery in any (Riesz) basis in a manner
that is completely stable.

In this paper, a new efficient computational algorithm is presented for
solving cyclic heptadiagonal linear systems based on using of heptadiagonal
linear solver and Sherman-Morrison-Woodbury formula. The implementation of the
algorithm using computer algebra systems (CAS) such as MAPLE and MATLAB is
straightforward. Numerical example is presented for the sake of illustration.

Algorithms have two costs: arithmetic and communication. The latter
represents the cost of moving data, either between levels of a memory
hierarchy, or between processors over a network. Communication often dominates
arithmetic and represents a rapidly increasing proportion of the total cost, so
we seek algorithms that minimize communication. In \cite{BDHS10} lower bounds
were presented on the amount of communication required for essentially all
$O(n^3)$-like algorithms for linear algebra, including eigenvalue problems and
the SVD.

We carry out a systematic numerical stability analysis of ZND detonations of
Majda's model with Arrhenius-type ignition function, a simplified model for
reacting flow, as heat release and activation energy are varied. Our purpose
is, first, to answer a question of Majda whether oscillatory instabilities can
occur for high activation energies as in the full reacting Euler equations,
and, second, to test the efficiency of various versions of a numerical
eigenvalue-finding scheme suggested by Humpherys and Zumbrun against the
standard method of Lee and Stewart.

As described in the classic works of Lee--Stewart and Short--Stewart, the
numerical evaluation of linear stability of planar detonation waves is a
computationally intensive problem of considerable interest in applications.
Reexamining this problem from a modern numerical Evans function point of view,
we derive a new algorithm for their stability analysis, related to a much older
method of Erpenbeck, that, while equally simple and easy to implement as the
standard method introduced by Lee--Stewart, appears to be potentially faster
and more stable.

We continue our study on variable order Arnold-Falk-Winther elements for 2D
elasticity in context of both affine and parametric curvilinear elements. We
present an $h$-stability result for affine elements, and an asymptotic
stability result for curvilinear elements. Both theoretical results are
confirmed with numerical experiments.

In this paper we present a steepest descent method with Armijo's rule for
multicriteria optimization in the Riemannian context. The well definedness of
the sequence generated by the method is guaranteed. Under mild assumptions on
the multicriteria function, we prove that each accumulation point (if they
exist) satisfies first-order necessary conditions for Pareto optimality.
Moreover, assuming quasi-convexity of the multicriteria function and
non-negative curvature of the Riemannian manifold, we prove full convergence of
the sequence to a Pareto critical.

In the random case setting, scrambled polynomial lattice rules as discussed
in \cite{BD10} enjoy more favourable strong tractablility properties than
scrambled digital nets. This short note discusses the application of scrambled
polynomial lattice rules to infinite-dimensional integration.

On the one hand the explicit Euler scheme fails to converge strongly to the
exact solution of a stochastic differential equation (SDE) with a superlinearly
growing and globally one-sided Lipschitz continuous drift coefficient. On the
other hand the implicit Euler scheme is known to converge strongly to the exact
solution of such an SDE. Implementations of the implicit Euler scheme, however,
require additional computational effort.

A technique for computing an ILU preconditioner based on the FAPINV algorithm
is presented. We show that this algorithm is well-defined for H-matrices.
Moreover, when used in conjunction with Krylov-subspace-based iterative solvers
such as the GMRES algorithm, results in reliable solvers. Numerical experiments
on some test matrices are given to show the efficiency of the new ILU
preconditioner.

In this paper modified variants of the sparse Fourier transform algorithms
from [14] are presented which improve on the approximation error bounds of the
original algorithms. In addition, simple methods for extending the improved
sparse Fourier transforms to higher dimensional settings are developed. As a
consequence, approximate Fourier transforms are obtained which will identify a
near-optimal k-term Fourier series for any given input function, $f : [0, 2 pi]
-> C, in O(k^2 \cdot D^4)$ time (neglecting logarithmic factors).

This paper develops a first-order system least-squares (FOSLS) formulation
for equations of two-phase flow. The main goal is to show that this
discretization, along with numerical techniques such as nested iteration,
algebraic multigrid, and adaptive local refinement, can be used to solve these
types of complex fluid flow problems. In addition, from an energetic
variational approach, it can be shown that an important quantity to preserve in
a given simulation is the energy law.

The following analog of Bernstein inequality for monotone rational functions
is established: if $R$ is an increasing on $[-1,1]$ rational function of degree
$n$, then $$ R'(x)<\frac{9^n}{1-x^2}\|R\|,\quad x\in (-1,1). $$ The exponential
dependence of constant factor on $n$ is shown, with sharp estimates for odd
rational functions.

This paper proposes a new preconditioning scheme for a linear system with a
saddle-point structure arising from a hybrid approximation scheme on the
sphere, an approximation scheme that combines (local) spherical radial basis
functions and (global) spherical polynomials. Making use of a recently derived
inf-sup condition [13] and the Brezzi stability and convergence theorem for
this approximation scheme, we show that the linear system can be optimally
preconditioned with a suitable block-diagonal preconditioner.

Recently a new class of Monte Carlo methods, called Time Relaxed Monte Carlo
(TRMC), designed for the simulation of the Boltzmann equation close to fluid
regimes have been introduced. A generalized Wild sum expansion of the solution
is at the basis of the simulation schemes.

In recent work (arXiv:1006.3100v1), we have presented a novel approach for
improving particle filters for multi-target tracking. The suggested approach
was based on Girsanov's change of measure theorem for stochastic differential
equations. Girsanov's theorem was used to design a Markov Chain Monte Carlo
step which is appended to the particle filter and aims to bring the particle
filter samples closer to the observations.

The famous Perona-Malik (P-M) equation which was at first introduced for
image restoration has been solved via various numerical methods. In this paper
we will solve it for the first time via applying a new numerical method called
the Variational Iteration Method (VIM) and the correspondent approximated
solutions will be obtained for the P-M equation with regards to relevant error
analysis. Through implementation of our algorithm we will access some effective
results which are deserved to be considered as worthy as the other solutions
issued by the other methods.

The Kaczmarz method is an algorithm for finding the solution to an
overdetermined system of linear equations Ax=b by iteratively projecting onto
the solution spaces. The randomized version put forth by Strohmer and Vershynin
yields provably exponential convergence in expectation, which for highly
overdetermined systems even outperforms the conjugate gradient method. In this
article we present a modified version of the randomized Kaczmarz method which
at each iteration selects the optimal projection from a randomly chosen set,
which in most cases significantly improves the convergence rate.

The present paper introduces an efficient and accurate numerical scheme for
the solution of a highly anisotropic elliptic equation, the anisotropy
direction being given by a variable vector field. This scheme is based on an
asymptotic preserving reformulation of the original system, permitting an
accurate resolution independently of the anisotropy strength and without the
need of a mesh adapted to this anisotropy. The counterpart of this original
procedure is the larger system size, enlarged by adding auxiliary variables and
Lagrange multipliers.

Bidimensional spiking models currently gather a lot of attention for their
simplicity and their ability to reproduce various spiking patterns of cortical
neurons, and are particularly used for large network simulations. These models
describe the dynamics of the membrane potential by a nonlinear differential
equation that blows up in finite time, coupled to a second equation for
adaptation. Spikes are emitted when the membrane potential blows up or reaches
a threshold value $\theta$.

A novel procedure is described for accelerating the convergence of Markov
chain Monte Carlo computations. The algorithm uses an adaptive bootstrap
technique to generate candidate steps in the Markov Chain. It is efficient for
symmetric, convex probability distributions, similar to multivariate Gaussians,
and it can be used for Bayesian estimation or for obtaining maximum likelihood
solutions with confidence limits. As a test case, the Law of Categorical
Judgment (Corrected) was fitted with the algorithm to data sets from simulated
rating scale experiments.

We study a variant of the well known Maxwell model for viscoelastic fluids,
namely we consider the Maxwell fluid with viscosity and relaxation time
depending on the pressure. Such a model is relevant for example in modelling
behaviour of some polymers and geomaterials. Although it is experimentally
known that the material moduli of some viscoelastic fluids can depend on the
pressure, most of the studies concerning the motion of viscoelastic fluids do
not take such effects into account despite their possible practical
significance in technological applications.

We prove quadratic eigenvalue perturbation bounds for generalized Hermitian
eigenvalue problems. The bounds are proportional to the square of the norm of
the perturbation matrices divided by the gap between the spectrums. Using the
results we provide a simple derivation of the first-order perturbation
expansion of a multiple eigenvalue, whose trailing term is tighter than known
results. We also present quadratic bounds for the non-Hermitian case.

We present a symbolic-numeric method to refine an approximate isolated
singular solution $\hat{\mathbf{x}}=(\hat{x}_{1}, \ldots, \hat{x}_{n})$ of a
polynomial system $F=\{f_1, \ldots, f_n\}$ when the Jacobian matrix of $F$
evaluated at $\hat{\mathbf{x}}$ has corank one approximately. Our new approach
is based on the regularized Newton iteration and the computation of approximate
Max Noether conditions satisfied at the singular solution.

We use the work of Milton, Seppecher, and Bouchitt\'{e} on variational
principles for waves in lossy media to formulate a finite element method for
solving the complex Helmholtz equation that is based entirely on minimization.
In particular, this method results in a finite element matrix that is symmetric
positive-definite and therefore simple iterative descent methods and
preconditioning can be used to solve the resulting system of equations. We also
derive an error bound for the method and illustrate the method with numerical
experiments.

The paper introduces the sweeping preconditioner, which is highly efficient
for iterative solutions of the variable coefficient Helmholtz equation
including very high frequency problems. The first central idea of this novel
approach is to construct an approximate factorization of the discretized
Helmholtz equation by sweeping the domain layer by layer, starting from an
absorbing layer or boundary condition. Given this specific order of
factorization, the second central idea of this approach is to represent the
intermediate matrices in the hierarchical matrix framework.

Most three dimensional constitutive relations that have been developed to
describe the behavior of bodies are correlated against one dimensional and two
dimensional experiments. What is usually lost sight of is the fact that
infinity of such three dimensional models may be able to explain these
experiments that are lower dimensional.

The well known algorithm of VOLKER STRASSEN for matrix multiplication can
only be used for $(m2^k \times m2^k)$ matrices. For arbitrary $(n \times n)$
matrices one has to add zero rows and columns to the given matrices to use
STRASSEN's algorithm. STRASSEN gave a strategy of how to set $m$ and $k$ for
arbitrary $n$ to ensure $n\leq m2^k$. In this paper we study the number $d$ of
additional zero rows and columns and the influence on the number of flops used
by the algorithm in the worst case ($d=n/16$), best case ($d=1$) and in the
average case ($d\approx n/48$).

This paper is concerned with a strongly degenerate convection-diffusion
equation in one space dimension whose convective flux involves a non-linear
function of the total mass to one side of the given position. This equation can
be understood as a model of aggregation of the individuals of a population with
the solution representing their local density. The aggregation mechanism is
balanced by a degenerate diffusion term accounting for dispersal.

The elliptic Monge-Amp\`ere equation is a fully nonlinear Partial
Differential Equation which originated in geometric surface theory, and has
been applied in dynamic meteorology, elasticity, geometric optics, image
processing and image registration.

Solutions can be singular, in which case standard numerical approaches fail.
Novel solution methods are required for stability and convergence to weak
solutions.

In this article we build a monotone finite difference solver for the \MA
equation, which we prove converges to the weak (viscosity) solution.

We consider the application of the polyharmonic subdivision wavelets (of
Daubechies type) to Image Processing, in particular to Astronomical Images. The
results show an essential advantage over some standard multivariate wavelets
and a potential for better compression.

Impulse methods are generalized to a family of integrators for Langevin
systems with quadratic stiff potentials and arbitrary soft potentials. Uniform
error bounds (independent from stiff parameters) are obtained on integrated
positions allowing for coarse integration steps. The resulting integrators are
explicit and structure preserving (quasi-symplectic for Langevin systems).

The data in many disciplines such as social networks, web analysis, etc. is
link-based, and the link structure can be exploited for many different data
mining tasks. In this paper, we consider the problem of temporal link
prediction: Given link data for times 1 through T, can we predict the links at
time T+1? If our data has underlying periodic structure, can we predict out
even further in time, i.e., links at time T+2, T+3, etc.? In this paper, we
consider bipartite graphs that evolve over time and consider matrix- and
tensor-based methods for predicting future links.

One can recover sparse multivariate trigonometric polynomials from few
randomly taken samples with high probability (as shown by Kunis and Rauhut). We
give a deterministic sampling of multivariate trigonometric polynomials
inspired by Weil's exponential sum. Our sampling can produce a deterministic
matrix satisfying the statistical restricted isometry property, and also nearly
optimal Grassmannian frames.

We propose a novel numerical method for solving a quadratic vector equation
arising in Markovian Binary Trees. The numerical method consists in a fixed
point iteration, expressed by means of the Perron vectors of a sequence of
nonnegative matrices. A theoretical convergence analysis is performed. The
proposed method outperforms the existing methods for close-to-critical
problems.

Various classes of stable finite difference schemes can be constructed to
obtain a numerical solution. It is important to select among all stable schemes
such a scheme that is optimal in terms of certain additional criteria. In this
study, we use a simple boundary value problem for a one-dimensional parabolic
equation to discuss the selection of an approximation with respect to time. We
consider the pure diffusion equation, the pure convective transport equation
and combined convection-diffusion phenomena.

A mathematical smooth function means that the function has continuous
derivatives to a certain degree C(k). We call it a k-smooth function or a
smooth function if k can grow infinitively. Based on quantum physics, there is
no such smooth surface in the real world on a very small scale. However, we do
have a concept of smooth surfaces in practice since we always compare whether
one surface is smoother than another one. This paper deals with the possible
definitions of "natural" smoothness and their relationship to the original
mathematical definition of smooth functions.

This paper presents a counterexample to the conjecture that the semi-explicit
Lie-Newmark algorithm is variational. As a consequence the Lie-Newmark method
is not well-suited for long-time simulation of rigid body-type mechanical
systems. The counterexample consists of a rigid body in a static potential
field.

Integration of the form $\int_a^\infty {f(x)w(x)dx} $, where $w(x)$ is either
$\sin (\omega {\kern 1pt} x)$ or $\cos (\omega {\kern 1pt} x)$, is widely
encountered in many engineering and scientific applications, such as those
involving Fourier or Laplace transforms. Often such integrals are approximated
by a numerical integration over a finite domain $(a,\,b)$, leaving a truncation
error equal to the tail integration $\int_b^\infty {f(x)w(x)dx} $ in addition
to the discretization error.

An exact, one-to-one transform is presented that not only allows digital
circular convolutions, but is free from multiplications and quantisation errors
for transform lengths of arbitrary powers of two. The transform is analogous to
the Discrete Fourier Transform, with the canonical harmonics replaced by a set
of cyclic integers computed using only bit-shifts and additions modulo a prime
number. The prime number may be selected to occupy contemporary word sizes or
to be very large for cryptographic or data hiding applications.

We describe von Neumann's elegant idea for sampling from the exponential
distribution, Forsythe's generalization for sampling from a probability
distribution whose density has the form exp(-G(x)), where G(x) is easy to
compute (e.g. a polynomial), and my refinement of these ideas to give an
efficient algorithm for generating pseudo-random numbers with a normal
distribution. Later developments are also mentioned.

Factorial series played a major role in Stirling's classic book "Methodus
Differentialis" (1730), but now only a few specialists still use them. This
article wants to show that this neglect is unjustified, and that factorial
series are useful numerical tools for the summation of divergent (inverse)
power series. This is documented by summing the divergent asymptotic expansion
for the exponential integral $E_{1} (z)$ and the factorially divergent
Rayleigh-Schr\"{o}dinger perturbation expansion for the quartic anharmonic
oscillator.

An error bound for the indefinite least squares problem is shown to be an
overestimate, contrary to conjecture. Consequently, the forward error analysis
of the hyperbolic factorization algorithm for those problems does not show the
algorithm to be forward stable.

Eigenvectors of tensors, as studied recently in numerical multilinear
algebra, correspond to fixed points of self-maps of a projective space. We
determine the number of eigenvectors and eigenvalues of a generic tensor, and
we show that the number of normalized eigenvalues of a symmetric tensor is
always finite. We also examine the characteristic polynomial and how its
coefficients are related to discriminants and resultants.

Finding the sparsest solution $\alpha$ for an under-determined linear system
of equations $D\alpha=s$ is of interest in many applications. This problem is
known to be NP-hard. Recent work studied conditions on the support size of
$\alpha$ that allow its recovery using L1-minimization, via the Basis Pursuit
algorithm. These conditions are often relying on a scalar property of $D$
called the mutual-coherence. In this work we introduce an alternative set of
features of an arbitrarily given $D$, called the "capacity sets".

We propose an algorithm for evaluation of the cumulative bivariate normal
distribution, building upon Marsaglia's ideas for evaluation of the cumulative
univariate normal distribution. The algorithm is mathematically transparent,
delivers competitive performance and can easily be extended to arbitrary
precision.

We consider methods for finding high-precision approximations to simple zeros
of smooth functions. As an application, we give fast methods for evaluating the
elementary functions log(x), exp(x), sin(x) etc. to high precision. For
example, if x is a positive floating-point number with an n-bit fraction, then
(under rather weak assumptions) an n-bit approximation to log(x) or exp(x) may
be computed in time asymptotically equal to 13M(n)lg(n), where M(n) is the time
required to multiply floating-point numbers with n-bit fractions. Similar
results are given for the other elementary functions.

Numerical algorithms have two kinds of costs: arithmetic and communication,
by which we mean either moving data between levels of a memory hierarchy (in
the sequential case) or over a network connecting processors (in the parallel
case). Communication costs often dominate arithmetic costs, so it is of
interest to design algorithms minimizing communication.

In this paper, we have given an idea of area specification and its
corresponding sensing of nodes in a dynamic network. We have applied the
concept of Monte Carlo methods in this respect. We have cited certain
statistical as well as artificial intelligence based techniques for realizing
the position of a node. We have also applied curve fitting concept for node
detection and relative verification.

We study in an unified fashion several quadratic vector and matrix equations
with nonnegativity hypotheses. Specific cases of such problems (QBD equations,
nonsymmetric algebraic Riccati equations, Lu's simple equation, Markovian
binary trees equations) have been studied extensively in the past by several
authors. Many of the results appearing here have already been proved for one or
more of the single instances of the problems, resorting to specific
characteristics of the problem.

In this paper approximations of three classes of fractional derivatives (FD)
using modified Gauss integration (MGI) and Gauss-Laguerre integration (GLI) are
considered. The main solutions of these fractional derivatives depend on
inverse of Laplace transforms, which are handled by these procedures. In the
modified form of integration the weights and nodes are obtained by means of a
difference equation, which gives a proper approximation form for the inverse of
Laplace transform and hence the fractional derivatives.

In this paper we present a set of criteria for the choice of the shape
parameter c contained in multiquadrics.

We report here on the recent application of a now classical general reduction
technique, the Reduced-Basis approach initiated in [C. Prud'homme, D. Rovas, K.
Veroy, Y. Maday, A. T. Patera, and G. Turinici. Reliable real-time solution of
parametrized partial differential equations: Reduced-basis output bounds
methods. Journal of Fluids Engineering, 124(1):7080, 2002.], to the specific
context of differential equations with random coefficients. After an elementary
presentation of the approach, we review two contributions of the authors: [S.
Boyaval, C. Le Bris, Y. Maday, N.C.

The exact computation of orbits of discrete dynamical systems on the interval
is considered. Therefore, a multiple-precision floating point approach based on
error analysis is chosen and a general algorithm is presented. The correctness
of the algorithm is shown and the computational complexity is analyzed. As a
main result, the computational complexity measure considered here is related to
the Ljapunow exponent of the dynamical system under consideration.

Quasi-Monte Carlo rules are equal weight quadrature rules defined over the
domain $[0,1]^s$. Here we introduce quasi-Monte Carlo type rules for numerical
integration of functions defined on $\mathbb{R}^s$. These rules are obtained by
way of some transformation of digital nets such that locally one obtains qMC
rules, but at the same time, globally one also has the required distribution.
We prove that these rules are optimal for numerical integration in fractional
Besov type spaces. The analysis is based on certain tilings of the Walsh phase
plane.

Heterogeneous anisotropic diffusion problems arise in the various areas of
science and engineering including plasma physics, petroleum engineering, and
image processing. Standard numerical methods can produce spurious oscillations
when they are used to solve those problems. A common approach to avoid this
difficulty is to design a proper numerical scheme and/or a proper mesh so that
the numerical solution validates the discrete counterpart (DMP) of the maximum
principle satisfied by the continuous solution.

The most critical component of any adaptive numerical quadrature routine is
the estimation of the integration error. Since the publication of the first
algorithms in the 1960s, many error estimation schemes have been presented,
evaluated and discussed. This paper presents a review of existing error
estimation techniques and discusses their differences and their common
features. Some common shortcomings of these algorithms are discussed and a new
general error estimation technique is presented.