XMa N.Zabares/Joumal f Computational Physics()44-471 8453 where the expectation value(can be evaluated by any numerical integration rule. Next.we consider a general non-linear function g(x.).where e is the random porosity.We need to express this function (27) 28) Here.))))ornd)oMThis is beca heterm Hermite poly- nomials are justo()By writing the k-Le jection onto each basis element)we can obtain from Eq.(27)the following: (2）-明 29 Thus,we obtain 哥-gK62 (30) ( This expression may be evaluated usng discusd-quraturebased approach suchas the ased on the Latin-Hypercubestatethiscase wer be The ide. subdivide the space of then calculate the value of the integrand in the numerator of Eq.(30).Summing all the values,the expectation value ()is just the arithmetic mean of these realizations.Also.the value of(can be pre-computed using quadrature rule. 4.2.GPCE-based formulation 院 By no nedia.n the sto cussed in the pre 议可. 1 31 --含 32) e(Xω) (33) (X.0 where Pis the number of expansion terms determined by the stochastic dimension and expansion order as in Eq.(23) hus.the stochastic problem is to fir the stochas ch tha e loe aand (150Da)e(x,t,) +PrV2v(x.t.@)-c(x.t.w)Vp(x.t.@)-c(x.t.@)PrRa0(x.t.w)eg. (34 7vx,t,=0. 65 .vx.t.)-V..)-vo.t.) (36)
where the expectation value h; i can be evaluated by any numerical integration rule. Next, we consider a general non-linear function gðx;Þ, where is the random porosity. We need to express this function as gðx;Þ ¼ XP i¼0 giWi; ð27Þ where gi is the expansion coefficient onto the polynomial basis and the porosity ðx;xÞ is written here based on the K–L expansion as follows: ðx;xÞ ¼ ðxÞ þXM i¼1 ffiffiffiffi ki p fiðxÞniðxÞ ¼ XP i¼0 iðxÞWiðxÞ: ð28Þ Here, 0ðxÞ ¼ ðxÞ, iðxÞ ¼ ffiffiffiffi ki p fiðxÞ, for i ¼ 1; ... ; M and iðxÞ ¼ 0, for i > M. This is because the first M þ 1 term Hermite polynomials are just W0 ¼ 1; W1 ¼ n1ðxÞ; ... ; WM ¼ nMðxÞ. By writing the K–L expansion as Eq. (28) instead of Eq. (16), it is easy to formulate and perform the polynomial chaos calculations. Using the same method as before (i.e., performing a Galerkin projection onto each basis element), we can obtain from Eq. (27) the following: g x; XP i¼0 iWi !; Wj * + ¼ gjhW2 j i: ð29Þ Thus, we obtain gj ¼ hgðx; PP i¼0iWiÞ; Wji hW2 j i : ð30Þ This expression may be evaluated using quadrature rule as discussed in [43] or a non-quadrature-based approach such as the integration, Taylor series and sampling approach discussed in [42]. Here, we employ a Monte-Carlo based sampling approach based on the Latin-Hypercube sampling (LHS) strategy [44]. In this case, we first generate samples of uncorrelated standard normal variables n using LHS. The idea of LHS is to subdivide the stochastic support space of the joint PDF of n into N subintervals along each stochastic dimension and to ensure that one sample of n lies in each subinterval. For each sample, we then calculate the value of the integrand in the numerator of Eq. (30). Summing all the values, the expectation value h; i is just the arithmetic mean of these realizations. Also, the value of hW2 j i can be pre-computed using quadrature rule. 4.2. GPCE-based formulation By now, we have developed all the tools we need to formulate natural convection in random porous media. In the stochastic natural convection problem, the input uncertainties are due to the Gaussian random field of porosity. Note that since there are non-linear functions ofðx; xÞ in the governing equation (1), we need to first express them in the polynomial basis using the method discussed in the previous section: 1 ðx;xÞ ¼ XP i¼0 ^iWi; ð31Þ ð1 ðx;xÞÞ2 ðx;xÞ 2 ¼ XP i¼0 iWi; ð32Þ 1 ðx;xÞ ðx;xÞ 2 ¼ XP i¼0 ~iWi; ð33Þ where P is the number of expansion terms determined by the stochastic dimension and expansion order as in Eq. (23). Thus, the stochastic problem is to find the stochastic functions that describe the velocity field vðx;t;xÞ : D ½0; T X ! Rd , the pressure field pðx;t;xÞ : D ½0; T X ! R and the temperature field hðx;t;xÞ : D ½0; T X ! R, such that the following equations are satisfied: ovðx;t;xÞ ot þ vðx;t;xÞ ðx;t;xÞ rvðx;t;xÞ¼ Pr Da ð1 ðx;t;xÞÞ2 ðx;t;xÞ 2 vðx;t;xÞ 1:75kvðx;t;xÞkð1 ðx;t;xÞÞ ð150DaÞ 1=2 ðx;t;xÞ 2 vðx;t;xÞ þ Prr2 vðx;t;xÞ ðx;t;xÞrpðx;t;xÞ ðx;t;xÞPrRahðx;t;xÞeg; ð34Þ r vðx;t;xÞ ¼ 0; ð35Þ ohðx;t;xÞ ot þ vðx;t;xÞ rhðx;t;xÞ ¼ r2 hðx;t;xÞ: ð36Þ X. Ma, N. Zabaras / Journal of Computational Physics 227 (2008) 8448–8471 8453
8454 XMa,N.Zabaras/Joumal of Computational Physics 7(0)8448-8471 Since the input uncertainty is taken as a Gaussian random field.we use Hermite polynomials to represent: z..-...) xto-xΨ⑤.t,o-4xΨ（⑤ (37) Substitution of Eqs.(28).(31)-(33)and(37)into the stabilized second-order projection method formulation Eqs.(7)and (10)-(13).results in the following: 2n-+会会会42wW-g网m 品2会-w+容92名w-g (38) w+加玄以"男=A三+会+会 (39 "-g (40) 立名w1-9%+玄-名防=0 2w-明+22四-呀m-云g。 (42 (0dtheroonity of thepoly -0,得言若名e听- -高含名-名2a网-+n-名名 -prRaeneu (43) 0-0 (At+t)p"AtV+ 4) 暖=Vp吸l (45) 心1-+暖1-哦=0， (46) -用+客呀时g (47) wheree0.1...P.This results in (P+1)(3d+)decoupled deterministic equations. Remark 2.In Eq.(43)it is time consuming to evaluate the fourth-order product term ()directly using the method discussed in Section 4.1.To simplify this calculation,we introduce an auxiliary random variable as follows:
Since the input uncertainty is taken as a Gaussian random field, we use Hermite polynomials to represent the solution: vðx;t;xÞ ¼ XP i¼0 viðx;tÞWiðnÞ; pðx;t;xÞ ¼ XP i¼0 piðx;tÞWiðnÞ; hðx;t;xÞ ¼ XP i¼0 hiðx;tÞWiðnÞ; pðx;t;xÞ ¼ XP i¼0 piðx;tÞWiðnÞ: ð37Þ Substitution of Eqs. (28), (31)–(33) and (37) into the stabilized second-order projection method formulation Eqs. (7) and (10)–(13), results in the following: 1 Dt XP i¼0 ðvnþ1=2 i vn i ÞWi þXP i¼0 XP j¼0 XP l¼0 ^lð2vn i rvn j vn1 i rvn1 j ÞWiWjWl ¼ Pr Da XP i¼0 XP j¼0 ið2vn j vn1 j ÞWiWj þ PrXP i¼0 r2 ðvnþ1=2 i ÞWi 1:75kvnk ð150DaÞ 1=2 XP i¼0 XP j¼0 ~ið2vn j vn1 j ÞWiWj XP i¼0 XP j¼0 irpn j WiWj PrRaXP i¼0 XP j¼0 ihn j WiWjeg; ð38Þ ðDt þ sÞr2XP i¼0 pnþ1 i Wi ¼ Dt XP i¼0 r2 pn i Wi þXP i¼0 r vnþ1=2 i Wi þ s XP i¼0 r pn i Wi; ð39Þ XP i¼0 pnþ1 i Wi ¼ XP i¼0 rpnþ1 i Wi; ð40Þ 1 Dt XP i¼0 ðvnþ1 i vnþ1=2 i ÞWi þXP i¼0 rpnþ1 i Wi XP i¼0 rpn i Wi ¼ 0; ð41Þ 1 Dt XP i¼0 ðhnþ1 i hn i ÞWi þXP i¼0 XP j¼0 ð2vn i rhn j vn1 i rhn1 j ÞWiWj ¼ XP i¼0 r2 hnþ1 i Wi: ð42Þ Then performing a Galerkin projection of each equation by h; Wki [8,10], and using the orthogonality of the polynomial basis Eq. (19), we obtain 1 Dt ðvnþ1=2 k vn k Þ þ hWiWjWlWki hW2 k i XP i¼0 XP j¼0 XP l¼0 ^lð2vn i rvn j vn1 i rvn1 j Þ ¼ Pr Da XP i¼0 XP j¼0 ið2vn j vn1 j Þeijk 1:75kvnk ð150DaÞ 1=2 XP i¼0 XP j¼0 ~ið2vn j vn1 j Þeijk þ Prr2 vnþ1=2 k XP i¼0 XP j¼0 irpn j eijk PrRaXP i¼0 XP j¼0 ihn j eijkeg; ð43Þ ðDt þ skÞr2 pnþ1 k ¼ Dtr2 pn k þ r vnþ1=2 k þ skr pn k ; ð44Þ pnþ1 k ¼ rpnþ1 k ; ð45Þ 1 Dt ðvnþ1 k vnþ1=2 k Þ þ rpnþ1 k rpn k ¼ 0; ð46Þ 1 Dt ðhnþ1 k hn k Þ þXP i¼0 XP j¼0 ð2vn i rhn j vn1 i rhn1 j Þeijk ¼ r2 hnþ1 k ; ð47Þ where eijk ¼ hWiWjWki hW2 k i , k ¼ 0; 1; ... ; P. This results in ðP þ 1Þð3d þ 2Þ decoupled deterministic equations. Remark 2. In Eq. (43), it is time consuming to evaluate the fourth-order product term 2 hW2 k i PP i¼0 PP j¼0 PP l¼0^lvn i rvn j hWiWjWlWki directly using the method discussed in Section 4.1. To simplify this calculation, we introduce an auxiliary random variable as follows: 8454 X. Ma, N. Zabaras / Journal of Computational Physics 227 (2008) 8448–8471
X Ma N.Zabares/Joumal of Computational Physics 27(0)844-471 8455 d-宫4-名名新m (48) such that 4-22am (49 =0.u==0 0=1 u=v=0 Free fluid Porous Medium Fg.Schematic of atural saturated variable porous medium .Comparison of th
d ¼ XP m¼0 dmWm ¼ XP i¼0 XP l¼0 ^lvn i WiWl; ð48Þ such that dm ¼ XP i¼0 XP l¼0 ^lvn i elim: ð49Þ Free fluid Porous Medium = 1 u = v =0 = 0 u = v =0 θ θ = 0, u = v =0 n ∂ θ ∂ = 0, u = v =0 n ∂ θ ∂ l ε =1 w ε Fig. 1. Schematic of natural convection in a fluid saturated variable porous medium. 0. 0 1 .2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 -1.73 -1.73 -1.73 -0.09 -0.09 -0.09 -0.09 -0.01 -0.01 -0.01 -0.01 -0.01 -5.48 -9 -8 -7 -6 -5 -4 -4 -3 -3 -2 -2 -1 -1 1 1 2 2 3 4 5 Fig. 2. Comparison of the streamlines (left) and isotherms (right). Top row: first-order projection method. Bottom row: second-order projection method. X. Ma, N. Zabaras / Journal of Computational Physics 227 (2008) 8448–8471 8455