On the expansion of solutions of Laplace-like equations into traces of separable higher dimensional functions

This paper deals with the equation $-\Delta u+\mu u=f$, $\mu$ a positive constant, on high-dimensional spaces $\mathbb{R}^m$. If the right-hand side $f$ is a rapidly converging series of separable functions, the solution $u$ can be represented in the same way. These constructions are based on the approximation of the function $1/r$ by sums of exponential functions. We derive results of similar kind for more general right-hand sides $f(x)=F(Tx)$ that are composed of a separable function $F$ on a space of a dimension $n$ greater than $m$ and a linear mapping given by a matrix $T$ of full rank. These results are based on the observation that in the high-dimensional case the euclidian norm of $T^t\omega$ behaves for $\omega$ in most of the $\mathbb{R}^n$ like the euclidian norm of $\omega$.


Introduction
The approximation of high-dimensional functions, whether they be given explicitly or implicitly as solutions of differential equations, represents a grand challenge for applied mathematics. High-dimensional problems arise in many fields of application such as data analysis and statistics, but first of all in the sciences. The Schrödinger equation, which links chemistry to physics and describes a system of electrons and nuclei that interact by Coulomb attraction and repulsion forces, forms one of the most important examples. The present work is partly motivated by applications in the context of the Schrödinger equation and is devoted to the equation with µ > 0 a given constant, on the R m for high dimensions m. Provided the righthand side f of this equation possesses an integrable Fourier transform, is a solution of this equation, and the only solution that tends uniformly to zero as x goes to infinity. If the right-hand side f of the equation is a tensor product of univariate functions or a rapidly converging series of such tensor products, the same holds for the Fourier transform of f . Using expansions of the function 1/r into sums of exponential functions [1], [2], one can show that the solution (1.2) can then, independent of the space dimension, be well approximated by a finite sum of such tensor products. This effect is meanwhile well-known and widely used and forms the basis of a considerable number of studies; [4] and [5] are predecessors of many of them and [3] is a more recent example. The aim of the present paper is to generalize such results and to represent the solution for right-hand sides like that are traces of higher dimensional separable functions in form of a series of functions of the same type, plus possibly a very small smooth remainder. Some basic facts about functions with Fourier transforms in L 1 and their lower dimensional traces are collected in the first section of this paper. In the second section we show that the solution (1.2) can for right-hand sides f (x) = F(T x), with T a given matrix, be written in the form u(x) = U(T x), where U is the solution of a degenerate elliptic equation L U = F in the higher dimensional space. Approximations to U are then iteratively generated. For corresponding right-hand sides F, this leads to expansions of the solution u into a series of traces of separable functions.

Functions with integrable Fourier transform and their traces
We consider in this paper functions u : R d → R, d a varying and potentially high dimension, that possess a then unique representation in terms of a function u ∈ L 1 (R d ), their Fourier transform. The space W 0 (R d ) of these functions (side remark, one of the Wiener algebras) becomes under the norm a Banach space. The functions in W 0 (R d ) are uniformly continuous and tend uniformly to zero as x goes to infinity, the Riemann-Lebesgue theorem. The norm (2.2) dominates their maximum norm. If the functions (iω) β u(ω), β ≤ α, in multi-index notation, are integrable as well, the partial derivative D α u of u exists, is given by and is as the function u itself uniformly continuous and tends uniformly to zero. For partial derivatives of first order, this follows from the Fourier representation of the corresponding difference quotients and the dominated convergence theorem, and for derivatives of higher order then by induction. Let T be a from now on fixed (n × m)-matrix of full rank m < n and let be the trace of a function U : R n → R with an integrable Fourier transform. We first calculate the Fourier transform of such trace functions. Let q m+1 , . . . , q n be an orthonormal basis of the orthogonal complement of the range of the matrix T . We add the columns q m+1 , . . . , q n to the m linearly independent columns of T , yielding an invertible matrix of dimension n that we denote by T . To avoid technical problems, we restrict ourselves here at first to rapidly decreasing functions.
Proof That the trace of a rapidly decreasing function is a rapidly decreasing function is rather obvious. To derive the representation (2.5) of its Fourier transform, first we introduce the rapidly decreasing function Its Fourier transform possesses after transformation of the variables the representation The trace function u(x) = V (x, 0) can be recovered from the Fourier transform of the full function V (x, y) via the inverse Fourier transformation as That is, the Fourier transform of u is the integrable function which is, as Fourier transform of a rapidly decreasing function, itself rapidly decreasing. Inserting the given representation of V , the proposition follows.

Lemma 2.2
The norm of the trace function (2.4) of the function U from Lemma 2.1 can be estimated by the norm of the function U itself as follows: Proof This is a rather immediate consequence of the representation (2.5) of the Fourier transform of u and of Fubini's theorem. Because of the estimate (2.6) for the scaled L 1 -norms (2.2) of u and U, respectively, follows again from the transformation theorem for multivariate integrals.
Remarkably the matrix T does not enter into the constant. As the rapidly decreasing functions are dense in W 0 (R n ), the trace functions (2.4) of functions U in W 0 (R n ) thus are functions in W 0 (R m ) of a norm not greater than that of U. Another consequence of Lemma 2.1 is that the traces of functions with Fourier transform vanishing outside of a strip around the kernel of T t are bandlimited. Proof It suffices to show the result for functions U with infinitely differentiable Fourier transform with bounded support inside the set under consideration as these are dense in the given subspace of W 0 (R n ). We split the vectors in R n again into the parts ω ∈ R m and η ∈ R n−m . Because of our assumptions on the support of U, holds for the ω and η for which the integrand in the representation (2.5) takes a value different from zero, which means that u(ω) must vanish for ω > Ω .

Shifted Laplace equations with trace functions as right-hand sides
We now return to the equation −∆ u + µu = f from Sect. 1. We show that its solution can, for a right-hand side f (x) = F(T x) that is trace of a function F : R n → R with Fourier transform in L 1 , be written as trace u(x) = U(T x) of a function U : R n → R that solves a degenerate elliptic equation.
with ω the euclidian norm of ω, is then the only solution of the equation on the R m that tends uniformly to zero as x goes to infinity.
Proof That u is a twice continuously differentiable function that solves the equation and tends uniformly to zero as x goes to infinity follows from the remarks made in the last section about functions with integrable Fourier transform. The maximum principle states that it is the only solution of the equation with this property.
The solution (3.1) of the equation (3.2) can also be characterized in a different way. We call a function u ∈ W 0 (R m ) a weak solution of this equation if holds for all rapidly decreasing functions ϕ.
and Fubini's theorem. If u 1 and u 2 are weak solutions of the equation, we have for all rapidly decreasing functions ϕ. As the equation −∆ ϕ + µϕ = χ possesses for all rapidly decreasing functions χ a rapidly decreasing solution ϕ, with a Fourier representation as above, for all rapidly decreasing functions χ then holds. The difference of u 1 and u 2 must therefore vanish.
Let the right-hand side now be the trace f (x) = F(T x) of a function F in W 0 (R n ). As such, it is by the results of the previous section a function in W 0 (R m ). The crucial observation is that we can lift the equation from R m into R n .
mapping the higher dimensional space to the real numbers.
Proof Using Fubini's theorem and observing that for rapidly decreasing ϕ holds, one recognizes that the trace function u is a weak solution of equation (3.2). As such, it coincides by Lemma 3.1 with the solution (3.1) of this equation.
The function (3.4) is in the domain of the operator L given by By definition, it solves the equation Because the expression µ + T t ω 2 is a second order polynomial in the components of ω, L can be considered as a second order differential operator and this equation therefore as a degenerate elliptic equation. If the Fourier transform of F and then also that of the solution U have a bounded support, U is infinitely differentiable and its derivatives can be obtained by differentiation under the integral sign. In this case, the function (3.4) is a classical solution of the equation (3.6). This is, however, irrelevant for the following considerations. The question under which conditions this generally holds is therefore not the subject of the present study. Instead of attacking the original equation (3.2) directly, we will approximate the solution of the higher dimensional equation (3.6). We are here primarily interested in right-hand sides F that are products of lower-dimensional functions or sums or rapidly converging series of such functions. The problem is now that the function is more difficult to expand into a series of separable functions than its counterpart in the representation (3.1) of the solution of the original equation, a problem that will find in the truly high-dimensional case a surprisingly simple solution.
4 The iterative solution of the degenerate elliptic equation To start with, let α : R n → R be a measurable, bounded function for which holds for all ω ∈ R n and assign to it the operator α : W 0 (R n ) → W 0 (R n ) given by The function U = αF can then serve as a first approximation of the solution (3.4). In general, this approximation will not be very precise but can be iteratively improved.  Proof The errors possess the representation The L 1 -norm of the Fourier transform of the errors is therefore The integrands are by the assumption made above bounded by the absolute value of the integrable function U and tend almost everywhere to zero as k goes to infinity. The dominated convergence theorem yields therefore which proves the proposition.
The same kind of result obviously also holds in the norm which dominates the derivatives up to second order of the trace of U, and in other weighted norms, corresponding properties of the right-hand side provided. Provided that on the support of F and thus also of the solution and the iterates 1 − α(ω) µ + T t ω 2 ≤ q < 1 (4.5) holds, one obtains from the Fourier representation of the errors the estimate in the norm (2.2) or (4.4). In this case, one can use Chebyshev acceleration to speed up the convergence of the iteration or, in other words, to improve the quality of the approximations of the solution. That means, one replaces the k-th iterate by a weighted average of the iterate itself and all previous ones. The errors of the recombined iterates U k possess then again a Fourier representation but now not with the polynomials P k (λ ) = (1 − λ ) k but polynomials The polynomial P k of degree k that satisfies the normalization condition P k (0) = 1 and attains the minimum maximal absolute value on an interval 0 < a ≤ λ ≤ b is, up to a linear transformation of the variable and the multiplication by a constant factor, the Chebyshev polynomial of degree k, a fact that is widely used in the analysis of Krylov-space methods. In our case, 1 − q ≤ λ ≤ 1 + q. Choosing the averaging coefficients α k accordingly, one gets the error estimate for the recombined iterates, where This is in comparison to the original iteration with the convergence rate a potentially big and very substantial improvement.

A particular class of iterative methods
Let ρ : R n → R be a measurable, locally bounded function for which holds for all ω ∈ R n . We assign to this ρ the iteration (4.3) based on the function Our primary example is the function ρ(ω) = T ω . The function (5.2) can then in a second step again be approximated by a sum of Gauss functions. A construction that keeps the relative error under control is analyzed in [6,Sect. 5]. In any case, That is, the condition (4.1) is because of (5.1) everywhere satisfied and the error tends in the norms (2.2) and (4.4) as well to zero. The key feature for understanding the convergence of the iteration comes from the analysis of the sets If ρ(ω) is a norm or seminorm, S(δ ) is a cone, that is, with ω every scalar multiple of ω is contained in S(δ ). Independent of the size of ρ(ω), on the set S(δ ) holds. If the Fourier transform of the right-hand side and with that also the Fourier transform of the solution vanish outside of S(δ ), this implies the error estimate The crucial point is that the regions S(δ ) fill in case of high dimensions m almost the complete space once δ falls below a certain bound, without any sophisticated adaption of ρ. We demonstrate this here for the function ρ(ω) = T ω .
Theorem 5.1 Let ρ(ω) = T ω and let κ be the condition of the matrix T , the ratio of its maximum and minimum singular value. If κδ < 1, for all R > 0 then holds, where λ is the volume measure on the R n and Equality holds if and only if κ = 1, that is, if all singular values of T coincide.
Proof We split the vectors in R n again into parts ω ∈ R m and η ∈ R n−m and start from a singular value decomposition T t = UΣV t of the matrix T t . As the multiplication with the orthogonal matrices U and V t , respectively, does not change the euclidian length of a vector, the set whose volume has to be estimated consists of the points with components ω and η in the ball of radius R around the origin for which holds. Because the volume is invariant under orthogonal transformations, the volume of this set coincides with the volume of the set of all points in the ball for which holds. Let 0 < σ 1 ≤ . . . ≤ σ m be the diagonal elements of the diagonal matrix Σ , the singular values of T t as well as of T , and let κ = σ m /σ 1 . Because of the given set is a subset of the set of all those points in the ball for which ω < κδ ω η holds. If κ = 1, that is, if all singular values of the matrix T are equal, The two sets then coincide and nothing is lost up to here. If κ > 1, there exists a vector with components ω 0 and η 0 inside of the ball under consideration for which holds. For this vector and thus also for all vectors sufficiently close to it we have That means that the two sets then differ and that the second one has a greater volume.
In what follows, we will calculate the volume of the latter one and compare it with the volume of the ball. We can restrict ourselves here to the radius R = 1. Let The set consists then of the points in the n-dimensional unit ball for which ω < ε η holds. Its volume can be expressed as double integral where H(t) = 0 for t ≤ 0, H(t) = 1 for t > 0, χ(t) = 1 for t ≤ 1, and χ(t) = 0 for arguments t > 1. It tends to the volume of the unit ball as ε goes to infinity. In terms of polar coordinates, this double integral reads Dividing this by the volume ν n of the unit ball itself and remembering that this completes the proof of the estimate (5.7) and shows that equality holds if and only if κ = 1, that is, if all singular values of T coincide. Moreover, we have shown that the function (5.8) tends to one as ε goes to infinity, and the bound for the ratio of the two volumes itself therefore to one as κδ goes to one. The ratio of the two volumes thus tends with the m-th power of δ to zero. It can also for rather large δ still attain extremely small values, which means that in most of the frequency space T t ω behaves like the norm of ω. For m = 128 and n = 256, for example, the ratio is for κδ ≤ 1/4 less than 1.90 · 10 −42 , and for κδ ≤ 1/2 still less than 6.95 · 10 −10 . Figure 5.1 shows the bound (5.7) as a function of κδ < 1 for n = 2m, with m = 2, 4, 8, . . . , 512. As long as the Fourier transform of the solution is not strongly concentrated around the kernel of T t , the regions of worse convergence will in such cases hardly affect the iterates and can basically be ignored.
Polynomial acceleration as described in Sect. 4 does not change this picture. If one inserts for the polynomial P k in (4.7) and (4.8) the polynomial that attains the minimum maximal absolute value on the interval δ 2 ≤ λ ≤ 1, the error estimate holds for the recombined iterates as long as the Fourier transform of the solution U vanishes outside of S(δ ). In the general case, the corresponding part of the error is thus reduced by a much larger factor than with the basic iteration. As the polynomial with the described properties satisfies for 0 ≤ λ ≤ 1 the estimate |P k (λ )| ≤ 1, the remaining part of the error does not blow up. We conclude that the iteration error can also in the described cases be reduced very substantially replacing the original iterates by linear combinations of all iterates up to the given one.
6 The limit behavior for fixed dimension ratios pointwise to zero and for arguments ε right of it pointwise to one. The maximum distance of the function values ψ(ε) to zero and one, respectively, tends for all ε outside any given interval around the jump discontinuity exponentially to zero.
Proof The proof is based on Stirling's formula in its logarithmic version If the ratio m/n is kept fixed and m goes to infinity, it leads to the representation ln 2Γ (n/2) where the constants β and β 0 depend only on the ratio m/n and are given by We keep the integers m and n in the following fixed, assume that they are relatively prime, and study with help of this representation the limit behavior of the functions as k goes to infinity. The with the prefactor multiplied integrands can be written as where C(k) remains bounded and tends to the limit lim k→∞ C(k) = m π e β m−β 0 .
The term in the brackets attains its global maximum at the point t = ε 0 specified above. As it increases strictly for t < ε 0 , takes the value one at t = ε 0 , and decreases strictly for t > ε 0 , there exists for every open interval around ε 0 a q < 1 with 0 ≤ e β m t m (1 + t 2 ) n/2 ≤ q for all t ≥ 0 outside of it. For arguments ε right of the interval we have therefore and for arguments ε > 0 left of the interval As the integrals are uniformly bounded in ε, this proves the proposition.
If the ratio m/n of the dimensions m and n is kept fixed and n and with that also the dimension m tend to infinity, the bound (5.7) from Theorem 5.1 for the ratio of the two volumes thus tends for values of κδ left of the jump discontinuity to zero and for values of κδ right of it to one, uniformly and exponentially outside every given interval around ξ 0 and the faster, the larger this interval is. The more the ratio of the dimensions approaches one, the closer the jump discontinuity comes to the endpoint of the admissible interval. We remark that one can, with help of the representation of the prefactor from the proof, also derive the estimate of the function ψ by the first term of its Taylor expansion at ε = 0, which describes the behavior of the bound when κδ approaches zero. The scaling factor ε 1 = e −β depends only on the ratio m/n and c(m) remains bounded as m goes to infinity.
7 Back to the motivating example The components of T t ω can be calculated via the relation T t ω · x = ω · T x. The square of the norm of T x can in matrix notation, with the vector e = (1, . . . , 1) t , be written as T x 2 = x · T t T x, T t T = 2 I − 1 m ee t .  The motivation behind this example is the electronic Schrödinger equation. In applications in conjunction with the electronic Schrödinger equation, the x i and ω i j are no longer real numbers but vectors in R 3 . That is, the matrix T can be reordered to a block matrix that is composed of three identical copies of the matrix above. λ ω T t ω < δ ω , ω ≤ R λ ω ω ≤ R ≤ ψ δ √ 1 − δ 2 , (7.7) as a function of δ < 1 for N = 2, 3, . . . , 18.