mathilde.chevreuil@univ- nantes.fr The paper deals with approximations of periodic functions that play a signiﬁcant role in GeM UMR CNRS 6183, harmonic analysis. The approach revisits the trigonometric polynomials, seen as Université de Nantes, Nantes, France combinations of functions, and proposes to extend the class of models of the Full list of author information is combined functions to a wider class of functions. The key here is to use structured available at the end of the article functions, that have low complexity, with suitable functional representation and adapted parametrizations for the approximation. Such representation enables to approximate multivariate functions with few eventually random samples. The new parametrization is determined automatically with a greedy procedure, and a low rank format is used for the approximation associated with each new parametrization. A supervised learning algorithm is used for the approximation of a function of multiple random variables in tree-based tensor format, here the particular Tensor Train format. Adaptive strategies using statistical error estimates are proposed for the selection of the underlying tensor bases and the ranks for the Tensor-Train format. The method is applied for the estimation of the wall pressure for a ﬂow over a cylinder for a range of low to medium Reynolds numbers for which we observe two ﬂow regimes: a laminar ﬂow with periodic vortex shedding and a laminar boundary layer with a turbulent wake (sub-critic regime). The automatic re-parametrization enables here to take into account the speciﬁc periodic feature of the pressure. Keywords: Adapted parametrization, Statistical learning, High-dimensional approximation, Tree tensor networks, Hierarchical tensor format Tensor train format, Adaptive algorithms Introduction The approximation of periodic functions plays a signiﬁcant role in harmonic analysis. In the case of the dynamical response of structures, these responses can notably be highly perturbed by low variability on the model and it then becomes necessary to develop reliable and eﬃcient tools for the prediction of the dynamical random response. We are here interested in constructing an approximation of a multivariate function with periodicity in one or more dimensions based on observations. This is of special interest for instance in uncertainty quantiﬁcation for vibroacoustic problems where the structure is excited with a harmonic wall pressure ﬁeld. The wall pressure ﬁeld is a multivariate function which © The Author(s) 2021. This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. 0123456789().,–: volV Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 2 of 16 depends on the time and on a set of variables, such as the Reynolds number. In practice, the wall pressure is computed for diﬀerent instances of the variables and consequently on diﬀerent discrete time grids depending on the instances. When ﬁne discrete models are involved, the evaluations of the model are costly and we may have access only to sparse information in terms of instances of the variables and of observation time interval. The set of trigonometric polynomials is well adapted for representing periodic functions. Indeed real trigonometric functions of degree m are written in the following form: v(t) = a + (a cos(nt) + b sin(nt)) 0 n n n=1 where the terms are periodic functions with period 2π. This class of functions has nice properties for approximation use, especially if a function f is continuous and given a tolerance ε, there exists a trigonometric polynomial v such that |f (t) − v(t)| <ε for all t. As mentioned above, in many applications we have access to samples of the polyno- mial from which we want to determine its coeﬃcients. The trigonometric polynomials are therefore linked to discrete-time signal processing, e.g. the Discrete-Time Fourier Transform (DTFT) converts a sequence of length N on an equally spaced time grid into a trigonometric polynomial of degree N − 1. The DTFT is extended to the d-dimensional case in the same manner. Given a sample of a multivariate function, the construction of an approximation of the function in the class of trigonometric functions has been widely addressed and the methods for constructing such representation generally depend on the discretization (see [1–3] and the references herein). In the present paper an alternative approach is proposed in order to tackle such prob- lems using a sample that is not necessarily structured. It is based on statistical learning methods [4] for multi-dimensional problems with s variables where the multivariate out- put function u(x , ··· ,x ) of the model, identiﬁed with an order-s tensor, is approximated 1 s in a parametrized subset of functions M = v = (a); a ∈ A (1) where the parameter a belongs to some set of parameters A and is a multilinear function with respect to the parameters a. The key idea is to propose an adapted parametrization with m new variables z = g (x , ... ,x ), i = 1, ... ,m, for the computation of the response i i 1 s so to obtain structured approximations with low complexity by exploiting the periodicity of the function in some dimensions. In the last decade, active subspace [6] and basis adap- tation methods [5] have been proposed to ﬁnd low dimensional structures using adapted parametrizations with reduced dimension. The ﬁrst class of methods consist in detecting the directions of strongest variability of a function using gradient evaluations, and then construct an approximation of the function exploiting the deﬁned low-dimension sub- space. In [7], active subspaces have been advantageously used for quantifying uncertainty of hypersonic ﬂows around a cylinder. The second class of methods, namely basis adap- tation methods, identiﬁes the dominant eﬀects in the form of linear combinations of the input variables and the adapted reduction of the representation is performed through a projection technique. In the current work the change of variables is extended to a wider class of functions g and is done with a method inspired from the projection pursuit method [8] which deﬁnes automatically and sequentially the new variables to add. The approxi- mation with this possibly high-dimensional new set of variables is created by exploiting Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 3 of 16 speciﬁc low-dimensional structures of the function such as sparsity [9]and lowrank[10] structures of the function to be approximated that enable the construction of an approx- imation using few samples as introduced in [11,12]. In the latter reference, the output is approximated in suitable low-rank tensor subsets of the form 1 L l n M = v = (a , ··· , a ); a ∈ R (2) where is a multilinear map with parameters a , l = 1, ... ,L. This is a special case of (1) 1 L n n 1 L with a = (a , ... , a )and A = R × ... × R . The dimension of the parametrization n grows linearly with m and thus makes possible the approximation with high l=1 dimension m. The ﬁrst part of the paper presents an interpretation of the trigonometric functions as a composition of functions h ◦g with speciﬁc structured representations for h, and proposes a generalization of the representation for h. The second part is dedicated to the algorithm for constructing the rank structured approximation combined with the change of variable g(x) used to handle periodic functions. Finally the last part illustrates the method on the pressure of a ﬂow around a cylinder. Periodic functions Let u ∈ H be a multivariate function, with H a Hilbert space, which depends on a set of independent variables X = (X , ... ,X ). In the present paper, we consider the speciﬁc 1 s case where the function u is periodic with respect to one variable denoted τ so that we have X = (, τ), with = (X , ... ,X )and τ = X ,and s = d + 1. d d+1 A variable X has values in X and an associated measure dp ,1 ≤ ν ≤ s. The variable ν ν ν X has values in X = X × ··· × X and an associated measure dp = dp × ··· × 1 s 1 dp . The variables X , ν = 1, ... ,d can be random variables, dp being in that case s ν ν a probability measure on X .Let H be the Hilbert space deﬁned on X,itisatensor space H = H ⊗ ··· ⊗ H with H a Hilbert space deﬁned on X . We consider that 1 s ν ν H ⊂ L (X ) is a ﬁnite dimensional subspace of square-integrable functions equipped ν ν 2 2 2 with the norm u = u dp and Hisasubspaceof L (X ) equipped with the canonical ν p 2 2 ν ν 1 norm u = u dp.Let {ψ } be an orthonormal basis of H and {ψ ⊗ ··· ⊗ X i i=1 i ψ } be an orthonormal basis of H. (i ,...,i )∈[1,...,P ]×···×[1,...,P ] i 1 s 1 d A natural representation of the periodic function u can be obtained using the Fourier series. In the following, x = (ξ,t) will denote an observation of X with ξ an observation of = (X , ... ,X )and t an observation of τ, i.e. a point in the periodic dimension. Trigonometric functions as a composition of functions Let us consider a T-periodic and continuous real valued function u : R → R.Itcan be represented by its truncated Fourier series which is a sum of harmonic functions: m−1 u(t) ≈ a cos(ω t) + b sin(ω t)(3) n n n n n=0 Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 4 of 16 where the circular frequencies are such that ω = nω (multiplicative constraint on ω ), n f n 2π with ω = , and the coeﬃcients a and b deﬁned as follows: f n n π π ω ω ω ω f f f f a = u(t)cos(ω t)dt and b = u(t)sin(ω t)dt, (4) n n n n π π 2π 2π − − ω ω f f for n = 0, ... ,m − 1. The truncated Fourier series can be seen as a composition of functions of the form: v(t) = h(g (t), ... ,g (t)) 1 m where h is an additive model h(z , ... ,z ) = h (z )with h (z ) = a cos(z ) + 1 m n n n n n n n=1 b sin(z ) ∈ V , where V = span 1, cos(·), sin(·) ,and g (t) = ω t for n = 1, ... ,m. n n n n−1 We propose here to extend the Fourier series to a more generalized framework where the function v is a multivariate function depending both on t and ξ, where the circular frequencies ω are chosen adaptively with no multiplicative constraint and where h is chosen in a wider set of functions than additive models. Generalizing Fourier series for the representation of a multivariate function in tensor format Let us focus on a function u(x) = u(ξ,t), periodic with respect to t. A natural representa- tion of u is thus by means of a relevant change of variables: u(x) ≈ h(g (x), ... ,g (x), ξ)(5) 1 m d+1 where g : R → R, n = 1, ... ,m are new variables chosen under the form g (ξ,t) = ω (ξ)t (6) n n ⊗m [0,T] × ×···× 1 d and where h ∈ V ⊗ H ⊗ ··· ⊗ H : R → R. Representing the periodic 1 d function under the form (5) can lead to the deﬁnition of a high number m of new variables so that we will consider subsets of low rank tensor formats for the high dimensional function h of dimension M = m + d. A M-dimensional function v in a subset of tensors can be written: v(z) = (z)(a , ... , a)(7) 1 L where (x) is a multilinear map with parameters (a , ... , a ). We consider here a model 1 L class of rank-structured functions associated with a notion of rank. A well known rank is the canonical rank associated to the sum of multiplicative models. The canonical rank of a function h is the minimal integer rank (v) = r such that 1 M v(z) = v (z ) ... v (z ) 1 M k k k=1 and we deﬁne the subset of canonical tensors: T ={v ∈ H :rank (v) ≤ r}. (8) r Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 5 of 16 r×n It can be associated to the parametrized representation (7)with L = M, where a ∈ R with n the dimension of the functional basis on which the functions v are represented, k = 1, ... ,r and l = 1, ... ,M. We can consider other notions of ranks which provide diﬀerent models with lower complexity. The α-rank of v, denoted by rank (v), is the minimal integer r such that α α v(z) = v (z )v (z ) α α k k k=1 with α ⊂{1, ... ,M},and z and z complementary groups of variables. The T-rank of v, α α denoted by rank (v) ={rank (v): α ∈ T }, is the tuple r ={r } such that T α α α∈T α α v(z) = v (z )v (z ), ∀α ∈ T α α k k k=1 where T is a collection of subsets of {1, ... ,M}. We deﬁne the subset of rank-structured functions: T ={v ∈ H :rank (v) ≤ r}. (9) The complexity of the associated parametrized representations of tensors is linear with the dimension M and polynomial with the ranks. Statistical learning method for approximating a function in tensor format with a change of variables Supervised statistical learning We consider a model that returns a real-valued variable Y = u(X). An approximation v of the function u, also referred to as metamodel, can be obtained by minimizing the risk R(v) = (u(x),v(x))dp(x) over a model class M. The loss function measures a distance between the observation u(x) and the prediction v(x). In the case of the least squares method, it is chosen as (y, v(x)) = (y − v(x)) . k k Let S ={(x ,y ):1 ≤ k ≤ N } be a sample of n realizations of (X, Y ). In practice, the approximation is constructed by minimizing the empirical risk k k R (v) = (y ,v(x )), (10) k=1 taken as an estimator of the risk. A regularization term R can be used for stability reasons when the training sample is small. An approximation u˜ of u is then solution of k k min y − v(x ) + λR(v), (11) v∈M k=1 with the regularization parameter λ ≥ 0, chosen or computed. The accuracy of the metamodel u˜ is estimated using an estimator of the L error. In practice, the number of numerical experiments is too small to sacriﬁce part of it for the error estimation. The Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 6 of 16 error is thus estimated using a k-fold cross validation estimator and more speciﬁcally the leave-one-out cross validation estimator [4] which can be easily evaluated by constructing one single metamodel [13]. Cross validation estimators can be used for model selection. In the following, we present a method to determine u˜(x) in a sequence of model classes M(m) = v(x) = h(g (x), ... ,g (x), ξ); h ∈ M ,g ∈ M for i = 1, ... ,m (12) 1 m h i g with M a linear or multilinear model class and g ∈ M a linear model class. We consider h i g here multilinear models and more speciﬁcally the tensor subset T in order to handle the possibly high dimensional m + d problem. We brieﬂy recall the learning algorithm in a tensor subset [12] and then present the automatic computation of the new variables z = g (x), i = 1, ... ,m. i i Learning with tensor formats Let z ∈ R , an approximation of u in a tensor subset (2) can be obtained by minimizing the empirical least squares risk: N L k k min y − (z )(a , ... , a ) + λ R (a ) (13) 1 L i i i 1 L a ,...,a i=1 k=1 where λ R (a ) are regularization functions. Problem (13) is solved using an alternating i i i minimization algorithm which consists in successively solving an optimization problem on a k k min y − (z )(... , a , ...) + λ R (a)(14) j j j j j N k=1 for ﬁxed parameters a , i = j. Introducing the linear map (z)(a ) = (z)(a , ... , a ), i j 1 L problem (14) yields the following learning problem with a linear model k j k min y − (z )(a ) + λ R (a ). (15) j j j j j N k=1 #v If R (a ) = a ,with v = |v | the -norm, problem (15) is a convex optimiza- j j j i 1 1 i=1 tion problem known as Lasso [14] or basis pursuit [15]. The -norm is a sparsity inducing regularization function that yields a solution a of (15) that may have coeﬃcients equal to zero. The Lasso is solved using the modiﬁed least angle regression algorithm (LARS in [16]). The algorithm to solve Problem (13) is described in [11] for the canonical tensor format, and in [12] for the tree-based tensor format, which is a special case of rank-structured tensor where T is a dimension partition tree. Adaptive algorithms are proposed to auto- matically select the tuple yielding a good convergence of the approximation with respect to its complexity. For tree-based tensor formats, at an iteration i, given an approximation i i v of u with T-rank (r ) , the strategy consists in estimating and studying the truncation α∈T error min R(v) − R(u) for diﬀerent α in T, and choosing to increase the ranks r rank (v)≤r α associated with the indices α yielding the highest errors. The algorithm and more details in the tree-based tensor format case can be found in [12,17]. Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 7 of 16 Learning method with automatic deﬁnition of new variables We now present the method used to automatically search an adapted parametrization of the problem by looking for favored directions in the space of the d + 1 input variables. It consists in writing the approximation under the form (5) where g (x) can be represented with a parametrized linear map g (x) = w ϕ(x) = w ϕ (x) (16) n n,j j j=1 with w = (w , ... ,w ) ∈ R the vector of parameters of the representation of g n n,1 n,p n on an orthonormal functional basis {ϕ } of H,and h a M-dimensional function in j=1 C T the model class of rank structured formats T or T that can be represented with a r r parametrized multilinear map 1 L h(z) = (z)(a , ... , a ) (17) with parameters a , l = 1, ... ,L. The new set of variables z = (z , ... ,z , ξ) is such that 1 m z = g (x), n = 1, ... ,m. n n The method is a Projection Pursuit like method [8] that is generalized to a larger class of models for h than just additive models. It is shown in [18] that under reasonable conditions, T 2 for h ∈ T and g ∈ H we have h ◦ g ∈ L (X ). The approximation u˜ of the form (5) is thus r p parametrized as follows: 1 L u˜(x) = w ϕ(x), ... , w ϕ(x), ξ (a , ... , a ). (18) Let (z , ... ,z ) be an initial set of variables. A new variable z = g (x) is introduced 1 m−1 m m using Algorithm 1. The parameters a , l = 1, ... ,L,and w , n = 1, ... ,m, solve the minimization problem (11) over the model class M(m): N L k k k k 1 L min y − w ϕ(x ), ... , w ϕ(x ), ξ (a , ... , a ) + λ R (a ), i i i L m {a } ,{w } l=1 n=1 k=1 i=1 (19) k k k with x = (ξ ,t ). The solution of this problem is found by alternatively solving N L k k k k 1 L min y − w ϕ(x ), ... , w ϕ(x ), ξ (a , ... , a ) + λ R (a )(20) i i i 1 m a ,...,a 1 L N k=1 i=1 Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 8 of 16 for ﬁxed (w , ... , w ) using a learning algorithm with rank adaptation [11,12]and 1 m k k k k 1 L min y − w ϕ(x ), ... , w ϕ(x ), ξ (a , ... , a ) (21) 1 m p p w ∈R ,...,w ∈R N 1 m k=1 1 L for ﬁxed (a , ... , a ). The optimization problem (21) is a nonlinear least squares problem that is solved with a Gauss-Newton algorithm. The overall algorithm is presented in Algorithm 2. In step 5 of algorithm 2, the parametrization of the model is to be selected in a collection of parametrizations. In this paper, we consider the rank structured function T where T is a dimension partition tree over {1, ... ,m}, which corresponds to the model class of functions in tree-based tensor format [19], a particular class of tensor networks [20]. The new node associated with the new variable z = g (x) is added at the top of the i i tree. The representation of a function v in T (H) requires the storage of a number of parameters that depends both on the collection T and on the associated T-rank r,aswell as on the dimensions of the functional spaces in each dimension. To reduce the number of coeﬃcients that need to be computed during the learning process of a tensor v,one can then represent v ∈ T (H) for diﬀerent collections T and associated T-ranks r, choosing the ones yielding the smallest storage complexity. Furthermore, this adaptation can prove useful when dealing with changes of variables, as introduced in the previous subsection, because it can remove the diﬃculty of how to add a variable: no matter what the initial ordering of the variables is, an adaptation procedure may be able to ﬁnd an optimal one yielding a smaller storage complexity. When T is a dimension partition tree, a stochastic algorithm is presented in [12] for trying to reduce the storage complexity of a tree-based tensor format at a given accuracy. This adaptation is not considered in the paper. Change of variables for periodic functions In this section we present a generalization of the Fourier series where one does not need a structured sample (e.g. a grid) in the variable t. Indeed when the function to approximate is known to be periodic with respect to t, the periodicity of the approximation can be forced. It is done on the one hand by introducing a functional basis ϕ such that g = w ϕ n n Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 9 of 16 in (16) can be identiﬁed to (6) by choosing (ξ) ϕ (x) = ϕ (ξ)t (22) (ξ) where {ϕ } is a d-dimensional tensorized orthogonal basis of H ⊗···⊗ H . The circular 1 d frequencies in (6) are then expressed as: (ξ) (ξ) ω (ξ) = w ϕ (ξ) = w ϕ (ξ). (23) n n,j j=1 n n On the other hand, we choose bases of trigonometric functions {ψ } for the represen- i i=1 tation of h in the dimensions associated with the new variables z , n = 1, ... ,m. Let T be the maximal width of the observation interval in the dimension t. Supposing max this interval is large enough to include the largest periods of the periodic functions that can be learned, the guarantee for the approximation not to have larger periods is to constrain the circular frequencies in (6) such that: 2π ω ≥ ,n = 1, ... ,m. (24) max This constraint is imposed for all values taken by ξ in S. Using expression (23), it is recast under the form −Aw ≤−B (25) (ξ) p×m N ×p where w = (w , ... , w ) ∈ R , A ∈ R is the array of evaluations of ϕ (ξ) for the 1 m k k N values (x , ... ,x ) of ξ in the training set S: d k=1 ⎡ ⎤ (ξ) (ξ) 1 1 1 1 ϕ (x , ... ,x ) ... ϕ (x , ... ,x ) 1 1 d 1 d ⎢ ⎥ . . ⎢ ⎥ . . . A = , (26) . . ⎣ ⎦ (ξ) (ξ) N N N N ϕ (x , ... ,x ) ... ϕ (x , ... ,x ) 1 1 d 1 d N ×m and B ∈ R is a full array with values −2π/T . The optimization problem (21) for the max computation of the parameters w is replaced with the constrained optimization problem k k k k 1 L min y − w ϕ(x ), ... , w ϕ(x ), ξ (a , ... , a ) (27) 1 m p p w ∈R ,...,w ∈R N 1 m k=1 −Aw≤−B which is solved with a NonLinear Programming (NLP) method. Application The method is applied to the prediction of the wall pressure p for a ﬂow over a cylinder for two ranges (low and medium) of Reynolds numbers for which we observe two ﬂow regimes: a laminar ﬂow with periodic vortex shedding and a laminar boundary layer with a turbulent wake (sub-critic regime). The automatic re-parametrization enables here to take into account the speciﬁc periodic feature of the pressure p with the time. The variables of the problem are X = (Re, , τ)with τ the time, Re the Reynolds number and the angular steps, that take values in X = X × X × X and we have d = 2. The pressure Re τ is evaluated on a tensor grid with Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 10 of 16 Fig. 1 Example of a linear tree with ﬁve variables yielding the tensor train format [21]: T ={{1}, {1, 2}, ... , {1, 2, 3, 4}} • at low Reynolds numbers: 300 time steps in X = [27 s, 29.4 s], 50 simulations of Re chosen uniformly on X = [70, 200] and 128 angular steps in X = [0, 2π[, Re • at higher Reynolds numbers: 1500 time steps, 11 simulations of Re chosen on X = Re [7000, 13000] and 320 angular steps in X = [0, 2π[. A representation (5) with new variables g (ξ,t), where ξ = (Re, θ), is computed with Algorithm 2 where the optimization problem (21) of step 6 is replaced with the constrained optimization problem (27). The new variables g are represented on a basis of functions (ξ) as in (22) where {ϕ } is a polynomial basis with maximal partial degree 2 (P = 9). j j=1 We choose for the model class of h the tensor train (TT) format associated to the linear tree T ={{1}, {1, 2}, ... , {1, ... ,M − 1}} (see Fig. 1) where M = m + d, the associated parametrized subset is 1 M l r ×n ×r l−1 l l v = (a , ... , a ); a ∈ R ,l = 1, ... ,M with r = (r ,r , ... ,r ) the TT-rank where r = r = 1, and n is the dimension 0 1 M 0 M l of the functional basis {ψ } in dimension l for the representation of h.Hereweuse l=1 trigonometric bases {1, cos(x), sin(x)} in the dimensions of z = g (t, θ,Re) and polynomial i i bases in the dimensions of θ and Re. The sequential quadratic programming (SQL) method is used to solve (27). Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 11 of 16 The accuracy of the metamodel u˜ is estimated using the unbiased test error estimator based on a test set S of N realizations of (X, Y ) independent of the training set S: test test R (˜u) 2 test e = (28) test R (0) test where R (v) is the empirical risk deﬁned in (10). Considering the least squares estimator, test test k k 2 (y −u˜(x )) N k=1 test we have e = . test N 1 test k 2 (y ) N k=1 test Low Reynolds numbers We ﬁrst consider the approximation of the wall pressure for low Reynolds numbers for which the data is given for 50 values of the Reynolds number in X = [70, 200]. The Re approximation is constructed with the following setting: • for g : polynomial bases with a maximal degree of 1 for t and 2 for θ and Re, • for h: polynomial bases of degree 14 for θ and 3 for Re, • training sample S using 20 simulations of Re and considering only the 181 ﬁrst time steps. The algorithm provided an approximation with dimension M = 6: X = (g (ξ), ... ,g (ξ), ξ) 1 4 and TT-ranks r = [1 3 4442 1]. The model error was estimated using two diﬀerent sam- ples: • estimation of the approximation error on the sample S with 30 simulations of Re test on the train time range (consisting of the 181 ﬁrst time steps): e = 1.14%, test • estimation of the extrapolation error on the sample S with 30 simulations of Re extra on the whole time range (consisting of the 300 available time steps): e = 1.19%. extra On Fig. 2 are plotted the predictions with blue crosses and the observations with circles, the red ones were used for learning the approximation and the green ones for estimating the model error. We observe a very good match of the predictions with the observations even beyond the train time range. The approximation constructed using the proposed change of variable is able to extrap- olate the wall pressure beyond the time range used for training the approximation. This extrapolation was made possible by introducing the constraint (24). As an illustration, Fig. 3 shows the observations (in blue) versus the predictions (in red) on a longer time interval obtained without constraint. The approximation obviously looses its periodicity. Table 1 summarizes the results obtained with the proposed approach and those obtained using the change of variable without using the tensor train format. One can observe that tensor formats enable to break the curse of dimensionality and thus to ease the learning of the approximation based on observations. Indeed the storage complexity of tensor train formats is O(MnR ) where n is the order of the dimension of the representation space in each dimension and R is the order of the rank. That is, it grows only linearly with the dimension M and quadratically with the rank whereas the storage complexity without using tensor formats, i.e. on the polynomial chaos, grows factorially or exponentially with the dimension M. Exploiting low complexity representations as low rank structures is necessary to address the problem when the dimension M increases with the deﬁnition of new variables. Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 12 of 16 Fig. 2 Low Reynolds numbers. Observations and predictions with respect to time for ﬁxed values of the other parameters. The predictions are obtained using tensor formats combined with constrained changes of variables Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 13 of 16 Fig. 3 Low Reynolds numbers. Observations and predictions with respect to time for ﬁxed values of the other parameters. The predictions are obtained using tensor formats combined with changes of variables (unconstrained) Table 1 Storage complexity and test error of the approximations obtained with the change of variable and with and without using Tensor Train (TT) format Storage complexity Test error e (%) test TT format + change of variables 594 1.14 Polynomial chaos + change of variables 8482 5.42 High Reynolds numbers We now consider the approximation of the wall pressure for higher Reynolds numbers for which the data is given for 11 values of the Reynolds number in X = [7000, 13000]. Re The approximation is constructed with the following setting: • for g : polynomial bases with a maximal degree of 1 for t and 2 for θ and Re, • for h: polynomial bases of degree 20 for θ and 6 for Re, • training sample S using 8 simulations of Re and considering only the 1000 ﬁrst time steps. The algorithm provided an approximation with dimension M = 7: X = (g (ξ), ... ,g (ξ), ξ) 1 5 and TT-ranks r = [13555511]. Themodel errorwas estimatedusing twodiﬀerent samples: • estimation of the approximation error on the sample S with 3 simulations of Re on test the train time range (consisting of the 1000 ﬁrst time steps): e = 2.97%, test • estimation of the extrapolation error on the sample S with 3 simulations of Re on extra the whole time range (consisting of the 1500 available time steps): e = 3.02%. extra On Fig. 4 are plotted the predictions with blue crosses and the observations with circles, the red ones were used for learning the approximation and the green ones for estimating the Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 14 of 16 Fig. 4 High Reynolds numbers. Observations and predictions with respect to time for ﬁxed values of the other parameters. The predictions are obtained using tensor formats combined with constrained changes of variables Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 15 of 16 model error. Again we observe a very good match of the predictions with the observations even beyond the train time range. The TT-rank of the approximation is low making possible the approximation using few samples of the Reynolds number R . Conclusion This paper presents a new strategy to approximate multivariate functions with peri- odicity. It gives the principles of the method based on the combination of functions h(g (ξ,t), ... ,g (ξ,t), ξ) chosen in appropriate classes of functions. The functions g (ξ,t) 1 m i deﬁne new variables of the multivariate functions h which is here represented in the class of rank structured functions. Algorithms are proposed for constructing the approxima- tion based on observations of the function, a constraint is added for the deﬁnition of the new parameters to promote periodicity of the representation. The numerical simulations yield good results. An analysis on the convergence of the approximation is to be studied. Acknowledgements Authors knowledge the support of the CNRS GdR 3587 AMORE. Authors’ contributions MC suggested the methodology, MS performed its implementation and run the numerical examples during her post-doctorate position at GeM. All authors contributed to the data analysis and discussed the content of the article. Both authors read and approved the ﬁnal manuscript. Funding Mathilde Chevreuil and Myriam Slama are grateful for the ﬁnancial support provided by the French National Research Agency (ANR) and the Direction Générale de l’Armement (DGA) (MODUL’O π project, Grant ASTRID Number ANR-16-ASTR-0018). Declarations Competing interests The authors declare that they have no competing interests. Author details 1 2 GeM UMR CNRS 6183, Université de Nantes, Nantes, France, LOMC UMR CNRS 6294, Université Le Havre Normandie, Normandie Université, Le Havre, France. Received: 29 December 2020 Accepted: 13 July 2021 References 1. Bass RF, Grochenig K. Random sampling of multivariate trigonometric polynomials. SIAM J. Math. Anal. 2006;36(3):773– 95. https://doi.org/10.1137/S0036141003432316. 2. Kammerer L, Potts D, Volkmer T. Approximation of multivariate periodic functions by trigonometric polynomials based on rank-1 lattice sampling. Journal of Complexity. 2015;31(4):543–76. https://doi.org/10.1016/j.jco.2015.02. 3. Briand T. Trigonometric polynomial interpolation of images. Image Processing On Line. 2019. https://doi.org/10.5201/ ipol.2019.273. 4. Hastie T, Tibshirani R, Friedman J. The Elements of Statistical Learning: Data Mining, Inference, and Prediction. : Springer; 2009. 5. Tipireddy R, Ghanem R. Basis adaptation in homogeneous chaos spaces. Journal of Computational Physics. 2014;259:304–17. 6. Constantine PG, Dow E, Wang Q. Active subspace methods in theory and practice: applications to kriging surfaces. SIAM Journal on Scientiﬁc Computing. 2014;36(4):1500–24. 7. Cortesi A, Constantine P, Magin T, Congedo PM. Forward and backward uncertainty quantiﬁcation with active subspaces: application to hypersonic ﬂows around a cylinder. INRIA: Technical report; 2017. 8. Friedman JH, Stuetzle W. Projection pursuit regression. Journal of the American statistical Association. 1981;76(376):817–23. 9. Hastie T, Tibshirani R, Wainwright M. Statistical Learning with Sparsity: the Lasso and Generalizations. : CRC Press; 2015. 10. Kolda TG, Bader BW. Tensor decompositions and applications. SIAM Review. 2009;51(3):455–500. https://doi.org/10. 1137/07070111X. 11. Chevreuil M, Lebrun R, Nouy A, Rai P. A least-squares method for sparse low rank approximation of multivariate functions. SIAM/ASA Journal on Uncertainty Quantiﬁcation. 2015;3(1):897–921. https://doi.org/10.1137/13091899X. 12. Grelier E, Nouy A, Chevreuil M. Learning with tree-based tensor formats; 2019. arXiv:1811.04455. Chevreuil and Slama Adv. Model. and Simul. in Eng. Sci. (2021) 8:17 Page 16 of 16 13. Cawley GC, Talbot NLC. Fast exact leave-one-out cross-validation of sparse least-squares support vector machines. Neural Networks. 2004;17:1467–75. 14. Tibshirani R. Regression shrinkage and selection via the lasso. Journal of the Royal Statistical Society. 1996;58(1):267–88. 15. Chen SS, Donoho DL, Saunders MA. Atomic decomposition by basis pursuit. SIAM Journal on Scientiﬁc Computing. 1999;20:33–61. 16. Efron B, Hastie T, Johnstone I, Tibshirani R. Least angle regression. The Annals of Statistics. 2004;32(2):407–99. 17. Nouy A, Grelier E, Giraldi L. Approximationtoolbox. Zenodo. 2020. https://doi.org/10.5281/zenodo.3653970. 18. Grelier E Learning with tree-based tensor formats : Application to uncertainty quantiﬁcation in vibroacoustics. PhD thesis, Centrale Nantes; 2019. 19. Falcó A, Hackbusch W, Nouy A. Tree-based tensor formats. SeMA Journal; 2018. 20. Cichocki A, Lee N, Oseledets I, Phan A-H, Zhao Q, Mandic D. Tensor networks for dimensionality reduction and large-scale optimization: Part 1 low-rank tensor decompositions. Foundations and Trends in ® Machine Learning. 2016;9(4–5):249–429. 21. Oseledets I, Tyrtyshnikov E. Recursive decomposition of multidimensional tensors. Doklady Math; 2009. Publisher’s Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional aﬃliations.
"Advanced Modeling and Simulation in Engineering Sciences" – Springer Journals
Published: Aug 4, 2021
Keywords: Adapted parametrization; Statistical learning; High-dimensional approximation; Tree tensor networks; Hierarchical tensor format Tensor train format; Adaptive algorithms
You can share this free article with as many people as you like with the url below! We hope you enjoy this feature!
Read and print from thousands of top scholarly journals.
Already have an account? Log in
Bookmark this article. You can see your Bookmarks on your DeepDyve Library.
To save an article, log in first, or sign up for a DeepDyve account if you don’t already have one.
Copy and paste the desired citation format or use the link below to download a file formatted for EndNote
All DeepDyve websites use cookies to improve your online experience. They were placed on your computer when you launched this website. You can change your cookie settings through your browser.