Belytschko T. - Introduction (779635), страница 83
Текст из файла (страница 83)
spatial coordinateFig 7.16.1 Problem statement and computational parameters for wave propagationTable 7.16.1 Time step sizes and numbers of time steps for elliptic-plastic wavepropagation example:W.K.Liu, Chapter 7Time Step (∆t )0.040 × 10−20.056 × 10 −20.072 × 10−266∆t / Cra0.50.70.9Number of time steps400286222aCr = ∆x / E / ρ +|c|7.15.2 Breaking of a damThis example is an attempt to model the breaking of a dam or more generally a flowwith large free surface motion by the ALE formulations described in the previous sections.This problem, which has an approximate solution for an inviscid fluid flowing over aperfect frictionless bed, presents a formidable challenge when this solution is applied tomine tailings embankments.
A detailed description of this problem can be found in Huerta& Liu(1988).The problem is solved without the restraints imposed by shallow water theory andonly the case of flow over a still fluid (FSF) is considered. Study on another case of flowover a dry bed (FDB) can be found in the paper of Huerta & Liu(1988). The accuracy ofthe ALE finite element approach is checked by solving the inviscid case, which has ananalytical solution in shallow water theory; then, other viscous cases are studied anddiscussed.Figure 7.8 shows a schematic representation of the flow over a still fluid Thedimensionless problem is defined by employing the following characteristic dimensions:the length scale is the height of the dam, H, over the surface of the downstream still fluid;the characteristic velocity, gH ,is chosen to scale velocities; and ρgH is the pressurescale.
The characteristic time is arbitrarily taken as the length scale over the velocity scale,i.e. H / g . Consequently, if the fluid is Newtonian, the only dimensionless parameterassociated with the field equations is the Reynolds number, Re = H gH / ν , where ν isthe kinematics viscosity. A complete parametric analysis may be found in Huerta (1987).Since the problem is studied in its dimensionless form, H is always set equal to one.Along the upstream and downstream boundaries a frictionless condition isassumed, whereas on the bed perfect sliding is only imposed in the inviscid case (forviscous flows the velocities are set equal to zero).
In the horizontal direction 41 elementsof unit length are usually employed, while in the vertical direction one, three, five, or sevenlayers are taken. depending on the particular case (see Figure 7.9). For the inviscidanalysis, ∆H = H = 1, as in Lohner et al (1984). In this problem both the Lagrange-Eulermatrix method and the mixed formulation are equivalent because an Eulerian description istaken in the horizontal direction; in the vertical direction a Lagrangian description is usedalong the free surface while an Eulerian description is employed everywhere else.Figure 7.10 compares the shallow water solution with the numerical resultsobtained by the one and three layers of elements meshes. Notice how the full integration ofthe Navier-Stokes equations smoothes the surface wave and slows down the initial motionof the flooding wave (recall that the Saint Venant equations predict a constant wave celerity,gH , from t = 0).
No important differences exists between the two discretizations (i.e.one or three elements in depth); both present a smooth downstream surface and a clearlyseparate peak at the tip of the wave. It is believed that this peak is produced in large part bythe sudden change in the vertical component of the particle velocity between still conditionsand the arrival of the wave, instead of numerical oscillations only. Figure 7.11 shows theW.K.Liu, Chapter 767difference between a Galerkin formulation of the rezoning equation, where numerical nodeto node oscillations are clear, and a Petrov-Galerkin integration of the free surface equation(i.e.
the previous 4lx3 element solution). The temporal criterion (Hughes and Tezduyar,1984) is selected for the perturbation of the weighting functions, and, as expected (Hughesand Tezduyar, 1984; Hughes and Mallet, 1986), better results are obtained if the Courantnumber is equal to one. In the inviscid dam-break problem over a still fluid, the secondorder accurate Newmark scheme (Hughes and Liu, 1978) is used (i.e. γ = 0.5 and β =0.25), while in all of the following cases numerical damping is necessary (i.e.
γ > 0.5)because of the small values of ∆H ; this numerical instability is discussed laterThe computed free surfaces for different times and the previous GeneralizedNewtonian fluids are shown in Figures 7.14 and 7.15. It is important to point out that theresults obtained with the Carreau-A model and n = 0.2 are very similar to those of theNewtonian case with Re = 300 , whereas for the Bingham material with µ p = 1× 10 2 Pa • sthe free surface shapes resemble more closely those associated with Re = 3000; this isexpected because the range of shear rate for this problem is from 0 up to 20-30 s −1 Itshould also be noticed that both Bingham cases present larger oscillations at the free surfaceand that even for the µ p = 1× 10 3 Pa • s case the flooding wave moves faster than that forthe Carreau models.
Two main reasons can explain such behavior; first, unlessuneconomical time-steps are chosen, oscillations appear in the areas where the fluid is atrest because of the extremely high initial viscosity (1000 µ p ); second, the larger shearrates occur at the tip of the wave, and it is in this area that the viscosity suddenly drops atleast two orders of magnitude, creating numerical oscillations.Exercise 7.1Observe that if the Jacobian described in Eq. (7.4.3a) is:J = det ∂x ∂x ∂x j ∂x k= ε ijk i ∂X ∂X1 ∂X2 ∂X3(7.4.11a)where ε ijk is the permutation symbol, then J˙ becomes:∂(v1 , x2 ,x 3 )∂(x1,v 2 , x3 )∂(x1 , x2 ,v3 )J˙ =++∂(X1, X 2 , X3 ) ∂(X1 , X2 , X3 ) ∂(X1, X 2 , X3 )where(7.4.11b)W.K.Liu, Chapter 768∂a∂X1∂(a,b,c)∂b=∂( X1 , X2 , X3 ) ∂X1∂c∂X1∂a∂X2∂b∂X2∂c∂X2∂a∂X3∂b∂X3∂c∂X3(7.4.11c)for arbitrary scalars a, b, and c, and vi = x˙ i .∂vUsing the chain rule on 1 , show that:∂X j∂(v1, v2 ,v3 )=∂( X1 , X2 , X3 )3∑m =1∂v1 ∂(x m , x2 ,v3 ) ∂v1=J∂xm ∂( X1 , X2 , X3 ) ∂x1(7.4.12a)Similarly, show that:J˙ = J vk ,k(7.4.12b)Exercise 7.2 Updated ALE Conservation of Angular MomentumThe principle of conservation of angular momentum states that the time rate of changeof the angular momentum of a given mass with respect to a given point, say the origin ofthe reference frame, is equal to the applied torque referred to the same point.
That is:Dx × ρ(x,t)v(x,t)dΩ = ∫ x × b(x,t)dΩ + ∫ x × t(x,t)dΓΩΓDt ∫Ω(7.4.13a)˙.It should be noticed that the left hand side of Eq. (7.4.13a) is simply H(2a) Show that:˙ =H∫Ω0=∫ΩDD(x) × (ρ v)JdΩ0 + ∫ x ×(ρ v J)dΩ0Ω0DtDt(7.4.13b)Dx × (ρ v) + (ρ v)div(v) dΩ DtHint, in deriving the above equation, the following pieces of information have been used of(1) x ,t[X] = v ; (2) v × (ρ v) = 0 and (3) Eq.
(7.3.6b).Now, show that substituting Eqs. (7.4.13b) into Eq. (7.4.13a) yields:W.K.Liu, Chapter 769D∫Ω x × Dt (ρ v) + (ρ v)div(v)dΩ=∫Ω(7.4.14)x × b(x,t)dΩ + ∫ x × (n ⋅ σ)dΓΓ(2b) Show that by employing the divergence theorem and the momentum equations givenin Eq. (7.4.6), the component form of Eq. (7.4.14) is:∫Ω ε ijkσ jkdΩ = 0(7.4.15)(2c) If the Cauchy stress tensor, σ , is smooth within Ω , then the conservation of angularmomentum leads to the symmetry condition of the Cauchy (true) stress via Eq. (7.4.15)and is given as:σ ij = σ ji(7.4.16)Exercise 7.3 Updated ALE Conservation of EnergyEnergy conservation is expressed as (see chapter 3):DρEdΩ = ∫ σ jin j vi dΓ + ∫ ρbi vi dΩ − ∫ qi ni dΓ + ∫ ρsdΩΓΩΓΩDt ∫Ω(7.4.17)where qi is the heat flux leaving the boundary ∂Ω x .
Recall that E is the specific totalenergy density and is related to the specific internal energy e, by:E =e+V22(7.4.18a)where e = e(θ ,ρ ) with θ being the thermodynamic temperature and ρs is the specific heatsource, i.e. the heat source per unit spatial volume and V 2 = vi vi . The Fourier law of heatconduction is:qi = −kijθ , j(7.4.18b)(3a) Show that the energy equation is (hint, use integration by parts and the divergencetheorem):(ρE),t[ χ] + (ρEc j ),j + ρEvˆ j, j = (σ ijvi ), j + bj v j + (kijθ ,j ),i + ρs(7.4.19a)(3b) If there is sufficient smoothness, time differentiate Eq. (7.4.19a) via the chain ruleand make use of the continuity equation to show that Eq.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
