{'fraction_non_alphanumeric': 0.038171115999762174, 'fraction_numerical': 0.017509958974968785, 'mean_word_length': 4.4169082125603865, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An adaptive reduced basis collocation method based on PCM ANOVA decomposition for anisotropic stochastic PDEs

The combination of reduced basis and collocation methods enables efficient and accurate evaluation of the solutions to parameterized PDEs. In this paper, we study the stochastic collocation methods that can be combined with reduced basis methods to solve high-dimensional parameterized stochastic PDEs. We also propose an adaptive algorithm using a probabilistic collocation method (PCM) and ANOVA decomposition. This procedure involves two stages. First, the method employs an ANOVA decomposition to identify the effective dimensions, i.e., subspaces of the parameter space in which the contributions to the solution are larger, and sort the reduced basis solution in a descending order of error. Then, the adaptive search refines the parametric space by increasing the order of polynomials until the algorithm is terminated by a saturation constraint. We demonstrate the effectiveness of the proposed algorithm for solving a stationary stochastic convection-diffusion equation, a benchmark problem chosen because solutions contain steep boundary layers and anisotropic features. We show that two stages of adaptivity are critical in a benchmark problem with anisotropic stochasticity.


INTRODUCTION
Numerical methods for high-dimensional anisotropic parameterized partial differential equations (PDEs) have received increasing interest in various areas of science and engineering. These problems are costly to solve because of requirements of high fidelity in space or high dimensionality in the parameterized random space. This study concerns some new approaches for reducing costs, by combining reduced basis methods for the spatial component of the problem with specialized quadrature rules for the parameterized component.
An approach to reduce spatial complexity is reduced basis methods, developed to cope with parameterized systems that require repeated use of a high-fidelity solver. Costly computations for high-resolution solutions for multiple parameter values are replaced by an approach with lower complexity, which computes parameterized solutions with a cost that is independent of the dimension of the original problem [2,19,21,[27][28][29][30]32,33]. Instead of solving the full discrete PDE system, a cheaper system is obtained by a projection into a lower-dimensional parametric space determined by the full-system solution at a distinguished set of parameters. Convergence analysis is available for various types of PDEs [27,31,33] using so-called certified reduced basis methods. Extension to parameterized stochastic PDEs is straightforward [7,8], and applications in stochastic control [6] have been studied as well.
The two components of reduced basis methods are referred to as an offline step (construction of the reduced space by computation of required full-system solutions) and an online step (computations on the reduced space) designed to be cheap. In addition to requirements of spatial resolution, the other contribution to cost comes from the high dimensionality of the stochastic space. In order to parametrize noisy input data for a given PDE, random fields can often be expanded as an infinite combination of random variables by, for instance, Karhunen-Loève [15,26] or polynomial chaos [40,41] expansions. Such random fields can be well approximated using a finite number of variables [34], but the dimensionality may be high. For such high dimensions in the stochastic space, high-dimensional integration techniques [3,16], based for example on Smolyak sparse-grid [22,24] or ANOVA decomposition [5,13,36], can be employed to compute the statistical moments of the solution. Although these methods are often efficient compared to Monte-Carlo methods, the curse of dimensionality remains a challenge. Thus, new ideas are necessary to increase efficiency and facilitate rapid evaluation of quantities of interest.
The idea of reduced basis methods combined with high-dimensional collocation schemes [10] brings the aforementioned efforts together. Instead of searching arbitrary samples, the candidate reduced basis points are chosen to be the sparse-grid probabilistic collocation points. If the solution vector is not of full rank in the parametric space, the computational cost can be reduced with essentially no loss of accuracy. However, the online stage involves a cost that is proportional to the number of collocation points. Thus, the algorithm suffers from the curse of dimensionality with costs that are less than for using the full collocation scheme, but still depending on the stochastic dimensionality. This cost can be further reduced by using dimension adaptation approaches, for instance, anisotropic sparse-grid [7] and ANOVA decomposition [18,20].
In this paper, we study the reduced basis collocation method for stochastic PDEs, and propose an adaptive procedure both in the dimensions of the stochastic space and the degrees of freedom in each parametric direction. An outline of the paper is as follows. In section 2, we review the reduced basis collocation method based on ANOVA collocation methods and examine the quality of different interpolation points. Then, we introduce the adaptive ANOVA procedure employed in the reduced basis algorithm, emphasizing its advantage for reducing computational cost and automatic sorting. Then, in section 3, we propose an adaptive reduced basis algorithm based on a probabilistic collocation method (PCM). In addition to the adaptive dimension procedure described in section 2, resolution in the parametric space is managed by adapting the polynomial order p discussed in section 3. We test these ideas on the stationary stochastic convection-diffusion equation chosen as a benchmark problem because of its interesting features such as solutions with steep boundary layers and discontinuities. The reduced basis ANOVA approach is suitable to study the interaction of these characteristics with the underlying noise that yields anisotropy in the stochastic space.

REDUCED BASIS COLLOCATION METHODS BASED ON ANOVA DECOMPOSITION
In this section, we review the reduced basis approximation for the solution computed using collocation in the stochastic space. Consider a stochastic problem consisting of a partial differential operator Ł and a boundary operator b, on D × Ω, where D ⊂ R d is the spacial domain, assumed to be a bounded and connected domain with polygonal boundary ∂D, and Ω is a complete probability space of (Ω, Σ, P) with σ-algebra Σ and probability measure P. Assume that the stochastic excitation can be represented by a finite number of random variables ξ = [ξ 1 , ..., . This could come from a variety of sources, for example, a truncated Karhunen-Loève (KL) expansion [15,26], a partitioning of D into subdomains, or uncertain boundary conditions.
In this paper, we consider the steady-state convection-diffusion equation, where a(x, ξ) is a parametrized random diffusion coefficient, w is a constant convective velocity, and f (x) is the forcing term. The boundary conditions are imposed as Dirichlet conditions with g D (x). We will consider choices of a(x, ξ) and g D (x) that involve discontinuity. The stochasticity in the diffusion coefficient a(x, ξ) is taken as a piecewise constant function on partitions, i.e., D is divided into N D equal-sized subdomains of horizontal strips, FIG. 1: Examples of subdomains for test problem (3) of horizontal strips, squares, and vertical strips with N D = 16.
with a constant ν and random variables {ξ m } N D m=1 . The dimensionality of the random space is identical to the number of subdomains, M = N D . The boundary conditions are imposed as Dirichlet conditions with The conditions yield jump discontinuities at the points (0, −1) and (−1, 1).
Examples of solutions to this system with a constant convective velocity w = (sin(π/6), cos(π/6)) and a constant forcing f = 1 are shown in Figure 2. The mean and standard deviation are plotted for the case of N D = 4 with a 2 × 2 partition and uniform random variables {ξ m } N D m=1 on Γ m = [0.01, 1] with diffusion parameter ν = 1/2 and 1/20. There is an exponential boundary layer near the top boundary x 2 = 1. The discontinuity from the boundary condition is smeared by the presence of diffusion, producing an internal layer. In addition, the standard deviation is less regular on the edges of the subdomains associated with a(·, ·). The diffusion coefficient lies in the range [0.01ν, ν], and the images on the right show results for a smaller diffusion parameter ν = 1/20. For this case, the overall variance on the domain tends to be smaller than for ν = 1/2, but the magnitude is larger in a narrow region near the top boundary layer. We will further study the impact of randomness on each subdomain on the irregular features of the solution using the reduced basis method combined with ANOVA decomposition. This approach will reveal the effects of the anisotropic nature of the problem and show the effectiveness of the method for anisotropic systems while obtaining good spatial resolution with a reasonable amount of computational cost.

Collocation methods
Collocation methods approximate the solution by interpolation at a set of Lagrange polynomials at collocation points, where Θ ⊂ Γ is a finite sample set (the collocation points), {L ξ (k) (ξ)} are the Lagrange interpolation polynomials on Γ, and the coefficient function u c x, ξ (k) is the solution of a deterministic problem corresponding to a given realization ξ (k) of the random vector ξ. We denote the weak form of the system (2) at a given realization ξ as where h is a grid parameter of the spatial domain. We consider a quadrilateral finite element space with the streamline diffusion method [12] to enhance stability when boundary layers are present. Nonhomogeneous boundary conditions are interpolated using addition basis functions at the boundary. Then the finite element solution u h (·, ξ) ∈ X h can be computed as where B (sd) ξ incorporates streamline-diffusion terms (see [12], section 6.3 for details). For fixed ξ, we refer to this sample solution as a snapshot. Given a finite set of parameters Θ, we can gather all the snapshots computed on Θ as This can be associated with a matrix S Θ ∈ R N h ×|Θ| , where each column of S Θ is a vector of coefficients of basis functions of the finite element solution. Here, we assume that the finite element approximation is sufficiently fine so that the spatial discretization error is acceptable. The full collocation solution (8) can be computed by applying a deterministic solver |Θ| times. If Θ consists of a full tensor product of one-dimensional parameter values, i.e., Θ = Θ i 1 1 ⊗ · · · ⊗ Θ i M M , where Θ im m is the set of collocation points in the m-th direction of order i m , costs grow exponentially with respect to the dimension M . The sparse-grid collocation method [3,16,22,24] reduces this cost using a subset of the full tensor grid collocation points where |i| = i 1 + · · · + i M and the level is the parameter that truncates the level of interaction. The choice of sparsegrid collocation points are typically taken to be nodes used for quadrature in M -dimensional space. The most common ones are Clenshaw-Curtis and Gauss abscissae. The convergence behavior of the reduced basis collocation method using Clenshaw-Curtis quadrature is studied in [10]. Here, we compare Clenshaw-Curtis and Gauss-Legendre points used in the reduced basis collocation method that will be described later. Figure 3 shows the error in the moments with respect to the size of the reduced basis N r . It is evident that the choice of the interpolation point is critical and Gauss-Legendre is the more effective choice for this example. In the sequel, we restrict our attention to Gauss-Legendre quadrature.

ANOVA decomposition for collocation methods
An alternative to sparse-grid collocation uses an ANOVA decomposition of an M -dimensional function. This decomposition has the form This representation decomposes multivariate functions according to the degree of interaction among the variables. For instance, u j1j2 ( ·, ξ j1 , ξ j2 ) describes the second-order or two-body interaction between the variables ξ j1 and ξ j2 . Similar to sparse-grid collocation, by truncating the series at a certain level, an M -dimensional function can be approximated using a function on a lower-dimensional subspace of the M -dimensional parameter space. Such a decomposition can be computed by u K ( ·, ξ K ) .
This construction can be put in a form similar to that of collocation by replacing the Lebesgue measure with a Dirac measure dµ = δ(ξ − c). This is called the anchored-ANOVA method where c is the anchor point [13,36]. The first few terms in the anchored-ANOVA decomposition become u 0 (·) = u( ·, c) = u( ·, By considering the ANOVA approximation in the random space discretized by the probabilistic collocation method (PCM) [1,39] in each parametric space, the solution approximated up to level can be written as where κ M,|K|, is a signed number counting the frequency of the sample solutions at ξ (k) ∈ Ξ p K appearing in the ANOVA expansion up to level , that is, κ M,j, = r=j (−1) [13,44] in direction K up to polynomial order p for all the parameters.
are the PCM points in the m-th parameter of polynomial order p m , then We denote the set of collocation points up to level and its size as respectively. For stochastic problems, it was shown in [23] that an ANOVA decomposition of level two has good accuracy if the anchor point is chosen as the mean of the probability density function considered, and we will use this technique.
An error indicator for the ANOVA truncation can be constructed using the relative size of either the mean or variance in each term [22,42]. In our work, we take the ANOVA indicator for the K-th term using the mean as

Reduced basis collocation method
Reduced basis algorithms aim to reduce the computational cost by approximating the solution using full discrete solutions computed only on a small subset of parameters Θ r ⊂ Ξ . Let N r .
= |Θ r |. The goal is to find a set of finite element solutions Q = {q 1 , · · · , q Nr } computed on Θ r that give an accurate representation of the column space of S Ξ . Such an approximate solution u rb (·, ξ) ∈ span(Q) satisfies The size N r of the reduced problem (14) is typically much smaller than N h from (7), and the utility of reduced basis methods derives from this fact [27,31]. Here we summarize the procedure. Assume that the operator Ł is linear and is affinely dependent on the vector of parameters ξ, that is, Then, the weak form (7) requires solution of a linear system of order N h , where For each q i ∈ Q, let q i denote the vector of coefficients associated with q i . Then Q = [q 1 , · · · , q Nr ] ∈ R N h ×Nr is the matrix representation of Q, and the linear system obtained from a Galerkin condition corresponding to the reduced problem (14) can be written as where Q u ξ is the reduced basis solution approximating u ξ . If the parameter-independent quantities Q T A i Q and Q T f are precomputed, the cost to assemble the reduced system for each ξ becomes O(N 2 r ) timesn a . The reduced basis algorithm computed with collocation points Θ is summarized in Algorithm 1. Either sparse-grid or ANOVA collocation can be employed by updating the reduced basis matrix while increasing the level from 0 to , that is, Q = RBM(Q, Ξ l , ), for 0 ≤ l ≤ (see [10,20] for details). Here, RB is an acceptance tolerance for the reduced basis solution based on an error indicator. We use a residual based error indicator, where [6,38]. An alternative error indicator is a dual-based indicator developed in [27,33], which depends on having an a posteriori error estimate for the associated problem. The total cost of a reduced basis algorithm can be divided into two parts as C(M, ) · (cost of reduced system solves and residual calculation) + (18) N r · (cost of full system solves and construction of offline quantities), where the former and latter terms are on the order of N r and N h , respectively. Here, C(M, ) is the number of sparsegrid nodes or ANOVA collocation points. The efficiency of the reduced basis algorithm comes from avoiding the full system solve on some collocation points. However, the cost becomes large as M and increase. Dimension-adaptive approaches can further improve the efficiency for anisotropic problems. The ANOVA indicator γ K (13) can be used Algorithm 1: (Q, S Θr ) = RBM(Q, Θ, RB ) : Reduced basis algorithm to update the reduced basis matrix Q considering the collocation points Θ and tolerance RB .

3:
if η Q,ξ (k) < RB then 4: Use the reduced solution u R (·, ξ (k) ) = Q u ξ (k) to serve as u c (·, ξ (k) ) in (5). 5: else 6: Compute the full solution vector u ξ (k) by solving (15). 7: Use u h (·, ξ (k) ) derived from u ξ (k) to serve as u c (·, ξ (k) ) in (5). 8: Update S Θr and Q by augmenting u ξ (k) and the orthogonal complement. 9: Reconstruct the offline reduced matrices and vectors, to select only the important directions. After each level l, we choose the ANOVA terms u K (10) such that γ K is greater than a tolerance A , and refer to the component of these terms as the effective dimensions J l , i.e., We use these indices to build the next level ANOVA approximation (see lines 19-20 of Algorithm 3 for details). The number of search points is reduced to C( M , ), where M is the number of the effective dimensions. This is a standard procedure in adaptive ANOVA schemes [22,42] and has been studied in the context of reduced basis methods in [18,20].
Here, we emphasize that we also employ the ANOVA indicator γ K (13) to sort the search direction J and the newly added reduced basis solution S Nr . Both are sorted in each level according to the descending order of γ K . Since the ANOVA indicator implies the priority of the reduced basis, this procedure enhances the quality of the reduced basis and saves the additional cost that occurs from sorting the reduced basis afterwards.

Numerical results
In this section, we use the reduced basis algorithm based on ANOVA to examine the features of the benchmark problem (2) and the accuracy of the method. We consider the spatial domain D = [−1, 1] × [−1, 1], with the system parameters described in the beginning of section 2. We take the random variables ξ m in the piecewise constant diffusion coefficient a(x, ξ) of (3) to be independent and uniformly distributed on Γ m = [0.01, 1]. The magnitude of the diffusion parameter ν is either 1/2 or 1/20 for examples of moderate and small diffusion. We choose the convective velocity w to correspond to a 30-degree angle right of vertical, and the Dirichlet boundary condition g D (x 1 , x 2 ) as in (4). For each sample point, we find a weak discrete solution u h (x, ξ) ∈ X h . For the finite element discretization of D, we use a bilinear (Q 1 ) finite element approximation [4]. The deterministic solver is based on the IFISS software package [11,35] and the spatial discretization is done on a uniform 128 × 128 grid, so that N h = 128 2 .
To explore accuracy with respect to parameterization, we compute a reference solution using a quasi Monte-Carlo (MC) method with 10 5 samples of a Halton sequence, and denote the computed mean and standard deviation as E[u mc ] and σ[u mc ], respectively. We then estimate the relative L 2 -error of reduced basis solutions u rb in the mean and standard deviation as We first present the results where the ANOVA decomposition helps us to understand the impact of anisotropy on features of the solution. We use the reduced basis tolerance RB = 10 −4 (see Algorithm 1) and the full ANOVA decomposition of the solution without truncation by setting A = 10 −15 ( see (19)). Figure 4 plots magnitudes of the means of the first and second-order ANOVA terms, E[u j (x, ξ j )] and E[u jk (x, ξ j , ξ k )] , for a 1 × 16 partition (see Figure 1). The image on the left shows that the ANOVA terms for components of ξ j near the top boundary (large j), where there is a boundary layer, and near the bottom (small j), where there is a discontinuity along the inflow boundary, are large. Among the second order terms, the pairs involving the aforementioned directions and the terms with physically adjacent subdomains, e.g, u j,j+1 (x, ξ j , ξ j+1 ), have large magnitudes as well. These points are further substantiated in Figure 5, where the absolute values of the means |E[u j (x, ξ j )]| are plotted for several j. In particular, there is significant activity near discontinuities, u 16 , u 14 , and u 1 . This is also true in the case of smaller diffusion parameter ν = 1/20 plotted in Figure 6. In addition, the system is more anisotropic compared to ν = 1/2   The dimension-adaptive procedure using ANOVA decomposition has an advantage for anisotropic problems, e.g., 1 × 16.   and the contribution is more localized in each subdomain. We emphasize that this feature will make the reduced basis algorithm based on ANOVA more efficient, since the existence of relatively important parameters indicates a possible truncation in the ANOVA decomposition. In constrast, Figures 7 and 8 show results for a 16×1 partition. The maximum values of the ANOVA terms are all of the same order, in other words, the system is less anisotropic with respect to the parameters than the case of 1 × 16, since all subdomains are attached to the top boundary layer. This is more apparent when ν = 1/20 in Figure 8, where the effect of the noise from each subdomain is more localized. Next, we examine the quality of the reduced basis sorted using the ANOVA indicator γ K (13) as in Algorithm 2. Figure 9 plots the relative L 2 -error in the moments with respect to the size of the reduced basis N r for the problems on 1×16 and 16×1 partitions with ν = 1/2. The errors are computed every two indices up to 100. The pink dashed lines show the error using the sparse-grid points of level two without the dimension-adaptive procedure, typically, generated in increasing order of index, e.g., ξ 1 , ξ 2 , · · · , and ξ 12 , ξ 13 , · · · . The initial saturation of the error in the 1 × 16 case is due to this ordering, and the error drops suddenly after obtaining the reduced basis that has a large contribution to the solution, for instance, when the points using index ξ 16 are considered. However, our algorithm sorts the solution by the magnitude of γ K and the error effectively decays while avoiding an additional computation induced either by adopting the greedy procedure or sorting the reduced basis afterwards. This computation becomes particularly expensive for large M and N r , which makes our algorithm more preferable. This advantage is not apparent in the 16 × 1 case shown in the second column, due to the fact that all directions are influenced in essentially the same way by the top boundary layer. However, our approach will show its effectiveness once the random variables {ξ j } have distinctive distributions.
We close this section by presenting a convergence study of the PCM-ANOVA points in Figure 10 and Table 1, which motivates the necessity of adapting the resolution in each direction. Unlike the sparse-grid points, the ANOVA points have another independent parameter, the polynomial order p, in addition to the level . Figure 10 shows the relative L 2 -error of the mean with respect to these two parameters for a 1 × 4 partition and ν = 1/20. The results are computed by either fixing the level = 3 or the polynomial order p = 9. Although in both cases, the errors decrease with respect to p and , they tend to stagnate after certain maximal values are reached. Thus, the polynomial order and the level should be increased together with a carefully selected tolerance regarding the system [13].

ADAPTIVE ALGORITHM FOR RBM-ANOVA
In this section, we introduce an adaptive reduced basis algorithm by increasing the number of collocation points used to compute each ANOVA term, in other words, the polynomial order p in each effective dimension J (19). As shown in Figure 10 and also later in Figure 11, fixed p may limit accuracy. However, it is not practical to choose large p from the beginning, since that leads to a large number of collocation points. Thus, we propose an adaptive PCM-ANOVA reduced basis algorithm that increases p only among the active dimensions J and with an appropriate criterion that limits the increases in p.
The idea is shown in Algorithm 3. The method initializes the collocation points with polynomial order p 0 , typically p 0 = 3, and level 0 . This choice will make the cost in the initial step to be smaller than that using a sparse-grid with the same level 0 . Then, among the effective dimensions J computed by the adaptive ANOVA criterion γ K (13), we increase p in specific directions until a termination criterion holds. For this, we use a saturation condition that measures the relative difference of ANOVA terms in direction K after increasing the polynomial order, that is, where u p K K is the ANOVA term in direction K computed with the increased polynomial order p K , and up K K is from the previous iteration, withp K < p K . If this quantity is smaller than a prescribed tolerance p = O( RB ), the reduced basis solution added in the current iteration is excluded from S Nr and the direction K is excluded from the set of effective dimensions J . This condition is based on the error estimate of the stochastic spectral solution between large polynomial orders [38]. In addition, if all the error indicators (17) computed from the collocation points in the direction K are less than the reduced basis tolerance RB , we exclude K from J . A stronger terminating criterion can be imposed by restricting the maximum order of polynomial at level K, |K| = l, from the previous level as p K ≤ max |S|=l−1 p S .

Numerical results
We first study the adaptive strategies in Algorithm 3 for benchmark problems on N D = 16 subdomains with 1 × 16, 4 × 4, and 16 × 1 partitions. The initial parameters are selected as p 0 = 3 and 0 = 2, and we successively increase while J = ∅ do 6: StoreJ = J . 7: for K ∈ J do 8: Storeū = u K andQ = Q. 9: Update Q in direction K, (Q, S K ) = RBM(Q, Ξ p K K , RB ) in Algorithm 1. 10: Compute the ANOVA term u K and indicator γ K in (13).

11:
if Q =Q or γ K < A or ρ p K K < p in (21) then 12: Exclude term K, i.e., J = J \ K, and reset u K =ū, Q =Q. p K by 2. We have three different tolerances to investigate, the reduced basis acceptance RB , the truncation of the ANOVA decomposition A , and the polynomial order p .
We first point out in Figure 11 the advantage of adapting the polynomial order p. The tolerance RB varies from 10 −2.5 to 10 −4 , while p = 10 −5 , and A = 10 −15 to prevent ANOVA truncation. The error of the reduced basis solution computed with fixed polynomial orders p = 3 or 5 does not decay after reaching a certain level. The accuracy is not improved while the size of the reduced basis increases. The adaptive procedure for p overcomes this problem.
Next, we study the sensitivity of the method with respect to the tolerance. One of the tolerances among RB , A , and p is varied from 10 −2.5 to 10 −4 , while others are fixed as 10 −4 . We also test the case where all the tolerances are reduced simultaneously as RB = A = p . The results are plotted in Figure 12 which shows that all the errors eventually converge to the level attained by RB = A = p = 10 −4 . Moreover, the selected reduced basis solution is more accurate when p is less than RB . In further simulations, we take p = RB /C, for some constant C, e.g., C = 2. On the other hand, results for A < 10 −3.5 show no difference since the effective dimensions J include all directions (see Figure 4 for the magnitude of the ANOVA terms in the 1 × 16 case). Thus, A should be chosen carefully by taking into consideration the trend of γ K to balance accuracy and computational complexity.
We now investigate the performance of our algorithm in higher dimensional parametric space, considering 6 × 6, 8 × 8, and 10 × 10 subdomain partitions. The reduced basis tolerance RB is varied from 10 −2.5 to 10 −4 , and A and p are taken to be proportional to RB , A = p = RB /2. In Figure 13 shows the solution on a 6 × 6 partition with diffusion parameter ν = 1/2 using tolerance RB = 10 −4 . On each subdomain D j with parameter ξ j , the figures plot the magnitude of the first order ANOVA terms E[u j ] , the adapted maximum polynomial order p j in direction j, and the number of reduced basis N rj selected in direction j. They show that the algorithm automatically increases the polynomial order and selects more reduced basis solutions on the subdomains near the top boundary layer and the bottom inflow discontinuity, where the ANOVA term is relatively large.
Let us compare the efficiency of the reduced basis algorithm using a fixed polynomial order p = 9 and adaptive determination of p. Table 2 shows the reduced basis size and the error in the moments for diffusion parameter ν = 1/20 (top) and ν = 1/2 (bottom). With a slightly smaller value of RB , the adaptive p approach yields a reduced basis set that achieves the same level of accuracy with the fixed p approach. Figure 14 shows accuracy versus computational time for both fixed p and adaptively determined p (These results were obtained in Matlab 2016 on a 3.33GHz Intel core i5 processor). It can be seen that changing p adaptively significantly reduces the computational cost. The computational cost saving of the adaptive p approach is primarily due to its use of fewer search points C(M, ) in (18), as shown in Table 3. We compare the total number of points on which the reduced basis is searched between using the PCM-ANOVA collocation points in the entire M -dimensions with a fixed p = 9, adapting directions with a fixed p = 9, and adapting both directions and p. Compared to considering a tensor product, which yields p M points, the number of PCM-ANOVA points C(M, ) is reduced by the truncation in the interaction level with = 2. However, the number is large due to the high dimensionality of size M , and the adaptive approach in both dimension and p reduces the computational cost. For the case of ν = 1/20, due to its strong anisotropy, the cost saving is significant. When ν = 1/2, although the cost saving from the truncation of directions is less than for ν = 1/20, the procedure of adaptive p reduces the cost further by increasing p only among the effective dimensions. The size of reduced basis N r and the error in the moments e µ and e σ for N D ×N D partitions and diffusion parameter ν = 1/2 (top) and ν = 1/20 (bottom), comparing the reduced basis algorithm with adaptive p and fixed p = 9. Corresponding computational time is plotted in Figure 14. 153 e µ 9.159e-03 5.574e-03 9.333e-04 8.447e-03 6.093e-03 9.052e-04 e σ 7.037e-01 3.208e-01 2.363e-02 6.196e-01 3.237e-01 2.448e-02

CONCLUSION
In this paper, we developed an adaptive reduced basis method for anisotropic parametrized stochastic PDEs based on the PCM and ANOVA approximation. The proposed algorithm enhances the reduced basis collocation methods by first obtaining the effective dimensions and then increasing the resolution in those directions. The dimensional adaptive procedure based on the ANOVA decomposition reduces the computational cost significantly for an anisotropic problem. In addition, the ANOVA indicator rearranges both of the search directions and the reduced basis solution, so that the reduced basis is sorted in a descending order of error with lower computational cost than the cost that typically occurs in greedy type approaches. The second adaptive procedure increases the resolution in each direction among the effective dimensions. This eliminates the stagnation of the error that may occur by fixing the polynomial order in PCM. Moreover, this adaptive procedure can be terminated when the ANOVA term converges. The proposed adaptive reduced basis algorithm is implemented on the stationary stochastic convection-diffusion problem. The algorithm effectively detects the directions in the stochastic space that are associated with the irregularity of the solution, for instance, discontinuous boundary conditions and steep boundary layers. By emphasizing those directions and by increasing the number of collocation points until necessary, the algorithm reduces the computational cost while improving accuracy in the moments. Our future work includes analysis of the proposed adaptive reduced basis algorithm using reduced basis error indicators based on posterior error estimates of parametrized PDEs [27] and combining them with analytical studies of adaptive PCM [14,38] and anchored ANOVA method [44]. In addition, we aim to further investigate the effect of the quadrature rules, e.g., Clenshaw-Curtis and Gauss abscissae, as candidate points to the accuracy of the reduced basis solution. Since different choice for the interpolation and anchor points in PCM-ANOVA yields significantly distinct set of collocation points, strategies to utilize the error estimates involving different interpolation points [25,37] and to choose the anchor point that minimizes the error estimates of the ANOVA approximation [43,44] may be useful. Finally, extending our algorithm to nonlinear and time-dependent problems, for example, stochastic Navier-Stokes equation, is another future direction. Various techniques including linearized iterative solvers such as Newton iterations [12,17] and empirical interpolation methods [2,9] can be employed.