# Fractional K-BKZ Numerical Model of the Start-Up Flow for a Viscoelastic Shock Absorber

Authors: Baoquan Mao; Zhiqian Wang; Yuying Yang; Rui Zhu; Qijin Zhao; Bowen Zheng
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

This paper has proposed a fractional K-BKZ numerical model by adopting the framework of the classical K-BKZ model and the relaxation modulus of the fractional Maxwell model with quasiproperties to study the start-up flow of a viscoelastic shock absorber. The start-up flows in both the orifice and the gap of a shock absorber were simplified to unidirectional accelerated flows in a pipe and between two parallel plates where one plate is accelerating and the other is at rest. The fractional K-BKZ numerical model was then developed using the finite difference method with real-world initial and boundary conditions. Numerical simulation was then performed, and the results were validated through laboratory testing, based on a comparison of the maximum fluid level and the contact angle. The proposed fractional K-BKZ numerical model successfully simulated the characteristics of the viscoelastic material passing through the orifice or the gap of a shock absorber, as demonstrated by accurately capturing the change of the shape of the flow. This fractional K-BKZ numerical model provided better accuracy for the fluid’s viscoelasticity and can be used for shock absorber design.


---

## Body

## 1. Introduction

Viscoelastic materials belong to the class of non-Newtonian pseudoplastic fluids. According to the mechanical structure of a shock absorber, the start-up flow forms of the viscoelastic material in the shock absorber can be divided into two types: the gap flow and the orifice flow. Due to the unique nonlinearity of viscoelastic fluid, many researchers have reviewed different research studies and viscoelastic materials and chose different mathematical models and constitutive relations to study the flow characteristics of viscoelastic fluids in a shock absorber. The research can be divided into two main categories: one involves empirical or semiempirical formula and the other involves classical models and their derived models.Jia et al. [1, 2] used a linear fluid model and a power function to study the flow characteristics of a viscoelastic fluid between two parallel plates. The semiempirical formula suitable for the impact environment was derived experimentally, which can provide guidance for the damper design. Zhang et al. [3] combined the fractional Maxwell model with the semiempirical formula reported in the literature [1] to design a damper for heavy vehicles and validated the model’s accuracy experimentally. Constantinou and Symans [4] studied a fluid damper that is used to make buildings earthquake-resistant and derived the empirical formula for experimental vibration environments. The model’s results agreed well with the experimental results. However, due to the different mechanical structures and the operating environments of the damper and absorber in the above research, as well as the formula that was based on a large amount of experimental data, the derived empirical formulas had different forms as well as more restrictions in their use and limited scopes of application. Furthermore, the quantitative measurement for the start-up flow process of a viscoelastic fluid is difficult, so it is usually obtained through numerical simulation.In the current research, the classical Maxwell model and its derivative models have attracted more and more attention. Especially with the introduction of the fractional derivative, the fractional Maxwell model has been validated and shown to have good applicability and can better simulate the flow process of a viscoelastic fluid [5, 6]. Yang and Zhu [7] used the fractional Maxwell model to study the unidirectional start-up flow of a viscoelastic fluid from a state of rest and under a constant pressure gradient in an infinitely long pipe and obtained an exact solution. Furthermore, it has been proposed that when the fractional element exhibits solid-like behavior, the fluid that is simulated using the fractional Maxwell model also exhibits the characteristics of a solid. Yan et al. [8] adopted the fractional Maxwell method to study the oscillating flow in a pipe in rolling motion to simulate the conditions of a ship’s pipelines subjected to seawater fluctuation; the velocity gradient and its variation of the oscillating flow of the pipe were then obtained. Tan et al. [9, 10] used the Laplace inverse transform with a continuous fractional derivative to discuss the fractional Maxwell model of a non-Newtonian viscoelastic fluid between two parallel plates when one plate undergoes random motion, and the analytical solutions for both the velocity and stress were then obtained. Carrera et al. [11] studied the fractional Maxwell model containing the relaxation process presenting non-Newtonian viscous behavior, and they analyzed and calculated the storage and loss modulus corresponding to the frequency and the amplitude function. Furthermore, they analyzed the stability of the function and found that the model’s results agreed well with the experimental results. On the basis of the converse ladder Maxwell model, Yao [12] introduced the tensor and different fractional derivative definitions to apply the fractional dashpot model to the simulation for the flow process of a nonlinear viscoelastic fluid under large deformation. Through this method, the transient viscosity and modulus could be simulated by varying the parameters. The aforementioned studies can be divided into two categories: a velocity category and a modulus category. The former is used to obtain the relationship between velocity, time, stress, etc. Then the fluid’s flow characteristics are characterized by the variation in velocity [1, 2, 7–10]; the latter is used to obtain the relationship between the modulus, viscosity, etc. as well as the frequency using the Fourier transform and other methods in order to characterize the rheological properties of the fluids [11, 12]. The two categories are different since the former is more intuitive for the comprehension of the motion of a viscoelastic fluid, while the latter can better explain the rheological characteristics of a viscoelastic fluid through the perspective of its mechanical properties. However, when the fluid’s velocity is higher, the fractional Maxwell model will gradually develop into a Newtonian fluid and lose viscoelasticity [13]. Similarly, when a viscoelastic fluid has a large deformation, the fractional Maxwell model is less nonlinear, which can lead to incorrect results [14]. A viscoelastic material belongs to the class of pseudoplastic fluids and the unique working characteristics of a shock absorber causes the viscoelastic fluid to have a higher velocity, so there are certain limitations to the fractional Maxwell model when simulating the flow characteristics of a viscoelastic material in a shock absorber.The K-BKZ constitutive equation is a relatively new model [15]. Many researchers have found that the K-BKZ model can describe nonlinear viscoelastic behavior very well [16–19]. Some researchers [20, 21] modified the K-BKZ model and derived the Wagner integral constitutive equation by ignoring the Cauchy tensor in the constitutive equation. The accuracy and consistency of the Wagner integral model prediction were demonstrated experimentally. Yamaguchi [22] used both the Oldroyd-B model and the K-BKZ model to simulate the flow conditions of an incompressible viscoelastic fluid within the inlet region of a circular tube. Bower et al. [23] used the K-BKZ model to simulate a fluid’s shear thinning effect for a specific fluid type, which can provide substantial insights about the combined effect of inertia and shear thinning in an orthogonal rheometer and provide predictive guidance for similar fluid models. Ebrahimi [24] et al. used the K-BKZ constitutive models to predict the transient shear response of high-density polyethylene melt. The comparison between the simulated results and the experimental data, obtained by a rotational rheometer, showed good agreement. By using the finite difference method, Tomé [25] et al. obtained the numerical solution to the K-BKZ integral constitutive equation with 2D free-surface motion, and the convergence of the mesh refinement in the numerical solution was validated and evaluated. The K-BKZ model with the relaxation kernel is described by power law function has the accuracy of the viscoelasticity model only at short times [14] and [26]. Jaishankar and McKinley [14] proposed a fractional K-BKZ model, which provided an accurate method to complete the transformation between the linear and nonlinear rheological study of multiscale materials. The model’s accuracy has been validated experimentally, which proved that the fractional K-BKZ model can effectively capture the shear thinning phenomenon of a viscoelastic fluid and maintain the good nonlinearity of the viscoelastic fluid.A viscoelastic material belongs to the class of pseudoplastic fluids and its flow inside a shock absorber has a larger flow velocity and a strong viscoelasticity, which can be studied using a K-BKZ model, especially the fractional K-BKZ model. Most of the related studies were focused on obtaining the analytical solutions of the K-BKZ model, except for two studies [13, 25], where numerical solutions were derived to obtain the models. However, the solution of the fractional K-BKZ model involves fractional calculus equations, and the process for solving the analytical solution is difficult and more complicated, so computational simulation is not easy. Comparatively, simulation on the basis of a numerical solution is easy to implement. The use of the fractional K-BKZ model to simulate the unidirectional accelerated start-up flows of viscoelastic fluid in a pipe has rarely been reported. Meanwhile, in the engineering field, it has been rarely reported that the fractional K-BKZ model was applied to a study on the start-up flow of a shock absorber.Therefore, aimed at the start-up flow of a viscoelastic material inside the orifice and gap of a shock absorber, the present paper has simplified the flow processes in both the orifice and the gap to accelerated motion in a pipe and the flow between a moving plate and a static plate. Subsequently, the fractional K-BKZ numerical model was constructed by adopting the framework of the classical K-BKZ model and the relaxation modulus of the fractional Maxwell model with quasiproperties. Comparative analysis of classical K-BKZ numerical model, the fractional Maxwell numerical model, and the Newtonian fluid model have been performed, and the unidirectional accelerated start-up flows of viscoelastic fluid in the orifice and the gap have been studied. A start-up flow test system was designed, and it captured the flow process of a viscoelastic material in the orifice of a shock absorber when the flow starts up, which validated the accuracy of the fractional K-BKZ numerical model.


## 2. Development of the Fractional K-BKZ Start-Up Flow Numerical Model

### 2.1. The Fractional K-BKZ Model

The classical K-BKZ constitutive equation is expressed as [26, 27](1)σt=∫−∞t2∂Wt−s,I1,I2∂I1C−1t,s−2∂Wt−s,I1,I2∂I2Ct,sds,where C−1=F−1T is the relative Finger strain tensor [26]. C is the relative Cauchy strain tensor. I1 and I2 are the first and second invariants, respectively, of C−1. Wt−s,I1,I2 is the potential function related to the strain energy function of the material can be factored into a time-dependent memory function mt−s and a scalar potential wI1,I2, i.e., Wt−s,I1,I2=mt−swI1,I2 [20, 28]. Wagner proposed [20, 26–28] that the second normal shear stress difference is not considered; that is, the relative Cauchy strain tensor is ignored, i.e., ∂wI1,I2/∂I2=0 [20, 28], the Wagner model can be written as [20, 26–28]:(2)σt=∫−∞tmt−shI1,I2C−1t,sds,where hI1,I2=2∂wI1,I2/∂I1 [20, 28]. In the unidirectional flow, the relative Finger strain tensor can be written as [27, 29, 30](3)C−1=1+ε2t,sεt,s0εt,s10001,then the expression of the shear stress can be simplified to [14, 20, 21, 27, 30](4)σt=∫−∞tmt−shεεt,sds,where mt−s is related to the relaxation modulus Gt−s and can be written as [14, 28](5)mt−s=dGt−sds.Using the relaxation modulus of the fractional Maxwell constitutive equation gives [14](6)Gt=ψt−βEα−β,1−β−ψκtα−β,where κ,α and ψ,β denote two mechanical elements of the spring-pots, α and β are fractional orders, and 0≤β<α≤1; κ and ψ are quasiproperties, and their units are Pa⋅sα and Pa⋅sβ, respectively; Eb,cz is the generalized Mittag-Leffler function, and b=α−β>0, c=1−β>0, z∈C, which is defined as [31, 32](7)Eb,cz=∑i=0∞ziΓbi+c.Through simultaneous solutions (6) and (7), the memory function (5) can be written as [14](8)mt−s=−ψt−s−1−βEα−β,−β−ψκt−sα−β,and here the form of the damping function is [28](9)hε=11+aε2t,s,where εt,s is the relative shear strain between any two moments, which can be defined as [14](10)εt,s=εt−εs=ε˙t−s.Substituting (8) and (9) into (4), the fractional K-BKZ model can be written as [14](11)σt=−∫−∞tψt−s−1−βEα−β,−β−ψκt−sα−βεt,s1+aε2t,sds.


### 2.2. Structure and Simplification

The geometrical schematic of a viscoelastic shock absorber has been shown in Figure1. The shock absorber was installed as a stressed structure by the shock absorber’s stand (Figure 1(h)). The viscoelastic material was injected into the shock absorber’s chamber (Figure 1(b)) through the viscoelastic material inlet (Figure 1(a)). The volume of the viscoelastic material filling was adjusted as needed to give the shock absorber a certain prestress, which was called the initial preload force. When the external impact force was greater than the initial preload force, the piston rod (Figure 1(e)) was pushed and drove the piston (Figure 1(d)) to compress the viscoelastic material in the shock absorber’s chamber (Figure 1(b)). The viscoelastic material flowed through the gap (Figure 1(c)) and the orifice (Figure 1(f)) of the piston after compression that generated the viscous damping force which hindered the piston from moving forward, and this process converted part of the external force into heat energy which was lost and potential energy which was stored. When the external force was withdrawn, that is, when the impact force was less than the initial preload force, the viscoelastic material expanded by itself and released the potential energy stored from the compression and some part of the potential energy converted to heat energy, while the rest potential energy pushed the piston back to its original position. The viscoelastic material flowed back to its initial state before compression through the orifice and the gap, and was ready for the next impact.Figure 1
The geometrical schematic of a viscoelastic shock absorber. (a) The viscoelastic material inlet. (b) The shock absorber’s chamber. (c) The gap. (d) The piston. (e) The piston rod. (f) The orifice. (g) The shock absorber’s chamber wall. (h) The shock absorber’s stand.In order to meet the actual working conditions of the start-up flow of a viscoelastic material in a shock absorber, the flow relations of the start-up flow of a viscoelastic material in the orifice and gap can be simplified into the two-dimensional space. The former (Figure1(f)) can be considered as the accelerated flow of a fluid in a pipe (see Figure 2(a)) while the latter (Figure 1(c)) can be considered as the flow of a fluid between two parallel plates while one plate is at rest and the other is accelerating (see Figure 2(b)).Figure 2
Simplifying the relationship of the motion between a shock absorber’s structure and the viscoelastic material. (a) The motion relationship between an orifice and the viscoelastic fluid and (b) a gap and the viscoelastic fluid.
(a)
(b)


### 2.3. Governing Equation

The time rate of change of the total momentum in the fluid’s motion is equal to the addition of its total volume force and its interactive surface force. The acceleration can be considered as the time rate of the momentum of the fluid having unit mass. According to Newton’s law of motion and the momentum balance theorem, the following can be obtained:(12)MA=∑F.The vector expression of the motion equation that is irrelevant to the coordinate system can be expressed as [33](13)ρ∂∂tU=ρf+∇⋅σ,where U is the flow velocity, f represents the unit volume force, ρ is the density, and σ is the stress tensor.Thus, the volume force is negligible and can be ignored. Then, the following can be obtained:(14)ρ∂∂tU=∇⋅σ.The velocity is only distributed in they direction (see Figure 2), so(15)ρ∂ut,y∂t=∂σt,y∂y.Whent≤0, εt=0, (11) can be written as(16)σt=−ψ∫0tθ−1−βEα−β,1−β−ψκθα−βεt−εt−θ1+aεt−θ−εt2dθ+ψ∫t∞θ−1−βEα−β,1−β−ψκθα−βεt1+aε2tdθ.Since the definition of the velocity gradient isdu/dy=ε˙, the strain definition is εt,y=∫0t∂εs,y/∂sds, and the definition of the Mittag-Leffler function derivatives is d/dsθ−βEα−β,1−β−κ/ψθα−β=θ−1−βEα−β,−β−κ/ψθα−βds [32]; the fractional K-BKZ governing equation can be written as(17)σt,y=−ψ∫0tθ−1−βEα−β,−β−ψκθα−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ψt−βEα−β,1−β−ψκtα−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,(18)∂σt,y∂y=ρ∂ut,y∂t=−ψ∫0tθ−1−βEα−β,−β−ψκθα−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫tt−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ψt−βEα−β,1−β−ψκtα−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.


### 2.4. The Initial Conditions and the Boundary Conditions

Whent≤0, ut=0. In the orifice, the flow velocity of the start-up flow at the center gradually increases with acceleration A, while the velocity near the pipe wall is always zero. In the gap, the flow velocity near the side of the moving plate increases with acceleration A, while the velocity near the stationary plate is always zero. So the initial conditions and the boundary conditions are as follows:(1)
The orifice:(19)ut,0=0,t>0,ut,L=0,t>0,ut,L2=At,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0.(2)
The gap:(20)ut,0=At,t>0,ut,L=0,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0,where L is the distance of the gap and the orifice and A is the acceleration.


### 2.5. Dimensionless Analysis

According to theΠ theorem, the dimensionless variables of the fractional K-BKZ model can be obtained as(21)u∗=ψκ1/α−βuA,y∗=L−1y,t∗=ψκ1/α−βt,σ∗=σLρA,a∗=aAL2κψ4/α−β,ϑ=ψL2ρκψ2−β/α−β.Substituting the dimensionless variables in (21) in (17) and (18), the fractional K-BKZ start-up flow dimensionless equation can be obtained (for brevity, the dimensionless mark “∗” is omitted):(22)σt,y=−ϑ∫0tθ−1−βEα−β,−β−θα−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ϑt−βEα−β,1−β−tα−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,(23)∂ut,y∂t=−ϑ∫0tθ−1−βEα−β,−β−θα−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫t−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ϑt−βEα−β,1−β−tα−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.The corresponding dimensionless initial conditions and boundary conditions of the orifice and gap are given as(24)ut,0=0,t>0,ut,1=0,t>0,ut,12=t,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0,(25)ut,0=t,t>0,ut,1=0,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0.


### 2.6. Numerical Solution

Numerical discretization was then applied to the fractional derivative using the finite difference method. The following were defined:tm=mΔt, m=0,1,2,…,K, 0≤k≤K−1; yn=nΔy, n=0,1,2,…,P. Δy=1/P is the space step, Δt=T/K is the time step. Applying the finite difference method to discretize the integrodifferential operators of (22) and (23), these can be approximated to(26)∫tatb∂uw,y∂ydw≈∑i=a+1buti,yn+1−uti,ynΔyΔt,∫tatb∂2uw,y∂y2dw≈∑i=a+1buti,yn+1−2uti,yn+uti,yn−1Δy2Δt.Letuji and σji denoteuti,yj and σti,yj. As Ftk=tk−1−βEα−β,−β−tkα−β and Htk=tk−βEα−β,1−β−tkα−β, the fractional K-BKZ start-up flow numerical model can be written as(27)σnm=−ϑ∑i=1mFtiΔt2/Δy∑j=m−i+1mun+1j−unj1+aΔt/Δy2∑j=m−i+1mun+1j−unj2+ϑHtmΔt/Δy∑j=1mun+1j−unj1+aΔt/Δy2∑j=1mun+1j−unj2,um+1,n−um,nΔt=−ϑ∑i=1mFtiΔt/Δy2∑j=m−i+1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21−2aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n2+ϑHtmΔt/Δy2∑j=1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=1muj,n+1−uj,n21−2aΔt/Δy2∑j=1muj,n+1−uj,n21+aΔt/Δy2∑j=1muj,n+1−uj,n2.


## 2.1. The Fractional K-BKZ Model

The classical K-BKZ constitutive equation is expressed as [26, 27](1)σt=∫−∞t2∂Wt−s,I1,I2∂I1C−1t,s−2∂Wt−s,I1,I2∂I2Ct,sds,where C−1=F−1T is the relative Finger strain tensor [26]. C is the relative Cauchy strain tensor. I1 and I2 are the first and second invariants, respectively, of C−1. Wt−s,I1,I2 is the potential function related to the strain energy function of the material can be factored into a time-dependent memory function mt−s and a scalar potential wI1,I2, i.e., Wt−s,I1,I2=mt−swI1,I2 [20, 28]. Wagner proposed [20, 26–28] that the second normal shear stress difference is not considered; that is, the relative Cauchy strain tensor is ignored, i.e., ∂wI1,I2/∂I2=0 [20, 28], the Wagner model can be written as [20, 26–28]:(2)σt=∫−∞tmt−shI1,I2C−1t,sds,where hI1,I2=2∂wI1,I2/∂I1 [20, 28]. In the unidirectional flow, the relative Finger strain tensor can be written as [27, 29, 30](3)C−1=1+ε2t,sεt,s0εt,s10001,then the expression of the shear stress can be simplified to [14, 20, 21, 27, 30](4)σt=∫−∞tmt−shεεt,sds,where mt−s is related to the relaxation modulus Gt−s and can be written as [14, 28](5)mt−s=dGt−sds.Using the relaxation modulus of the fractional Maxwell constitutive equation gives [14](6)Gt=ψt−βEα−β,1−β−ψκtα−β,where κ,α and ψ,β denote two mechanical elements of the spring-pots, α and β are fractional orders, and 0≤β<α≤1; κ and ψ are quasiproperties, and their units are Pa⋅sα and Pa⋅sβ, respectively; Eb,cz is the generalized Mittag-Leffler function, and b=α−β>0, c=1−β>0, z∈C, which is defined as [31, 32](7)Eb,cz=∑i=0∞ziΓbi+c.Through simultaneous solutions (6) and (7), the memory function (5) can be written as [14](8)mt−s=−ψt−s−1−βEα−β,−β−ψκt−sα−β,and here the form of the damping function is [28](9)hε=11+aε2t,s,where εt,s is the relative shear strain between any two moments, which can be defined as [14](10)εt,s=εt−εs=ε˙t−s.Substituting (8) and (9) into (4), the fractional K-BKZ model can be written as [14](11)σt=−∫−∞tψt−s−1−βEα−β,−β−ψκt−sα−βεt,s1+aε2t,sds.


## 2.2. Structure and Simplification

The geometrical schematic of a viscoelastic shock absorber has been shown in Figure1. The shock absorber was installed as a stressed structure by the shock absorber’s stand (Figure 1(h)). The viscoelastic material was injected into the shock absorber’s chamber (Figure 1(b)) through the viscoelastic material inlet (Figure 1(a)). The volume of the viscoelastic material filling was adjusted as needed to give the shock absorber a certain prestress, which was called the initial preload force. When the external impact force was greater than the initial preload force, the piston rod (Figure 1(e)) was pushed and drove the piston (Figure 1(d)) to compress the viscoelastic material in the shock absorber’s chamber (Figure 1(b)). The viscoelastic material flowed through the gap (Figure 1(c)) and the orifice (Figure 1(f)) of the piston after compression that generated the viscous damping force which hindered the piston from moving forward, and this process converted part of the external force into heat energy which was lost and potential energy which was stored. When the external force was withdrawn, that is, when the impact force was less than the initial preload force, the viscoelastic material expanded by itself and released the potential energy stored from the compression and some part of the potential energy converted to heat energy, while the rest potential energy pushed the piston back to its original position. The viscoelastic material flowed back to its initial state before compression through the orifice and the gap, and was ready for the next impact.Figure 1
The geometrical schematic of a viscoelastic shock absorber. (a) The viscoelastic material inlet. (b) The shock absorber’s chamber. (c) The gap. (d) The piston. (e) The piston rod. (f) The orifice. (g) The shock absorber’s chamber wall. (h) The shock absorber’s stand.In order to meet the actual working conditions of the start-up flow of a viscoelastic material in a shock absorber, the flow relations of the start-up flow of a viscoelastic material in the orifice and gap can be simplified into the two-dimensional space. The former (Figure1(f)) can be considered as the accelerated flow of a fluid in a pipe (see Figure 2(a)) while the latter (Figure 1(c)) can be considered as the flow of a fluid between two parallel plates while one plate is at rest and the other is accelerating (see Figure 2(b)).Figure 2
Simplifying the relationship of the motion between a shock absorber’s structure and the viscoelastic material. (a) The motion relationship between an orifice and the viscoelastic fluid and (b) a gap and the viscoelastic fluid.
(a)
(b)


## 2.3. Governing Equation

The time rate of change of the total momentum in the fluid’s motion is equal to the addition of its total volume force and its interactive surface force. The acceleration can be considered as the time rate of the momentum of the fluid having unit mass. According to Newton’s law of motion and the momentum balance theorem, the following can be obtained:(12)MA=∑F.The vector expression of the motion equation that is irrelevant to the coordinate system can be expressed as [33](13)ρ∂∂tU=ρf+∇⋅σ,where U is the flow velocity, f represents the unit volume force, ρ is the density, and σ is the stress tensor.Thus, the volume force is negligible and can be ignored. Then, the following can be obtained:(14)ρ∂∂tU=∇⋅σ.The velocity is only distributed in they direction (see Figure 2), so(15)ρ∂ut,y∂t=∂σt,y∂y.Whent≤0, εt=0, (11) can be written as(16)σt=−ψ∫0tθ−1−βEα−β,1−β−ψκθα−βεt−εt−θ1+aεt−θ−εt2dθ+ψ∫t∞θ−1−βEα−β,1−β−ψκθα−βεt1+aε2tdθ.Since the definition of the velocity gradient isdu/dy=ε˙, the strain definition is εt,y=∫0t∂εs,y/∂sds, and the definition of the Mittag-Leffler function derivatives is d/dsθ−βEα−β,1−β−κ/ψθα−β=θ−1−βEα−β,−β−κ/ψθα−βds [32]; the fractional K-BKZ governing equation can be written as(17)σt,y=−ψ∫0tθ−1−βEα−β,−β−ψκθα−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ψt−βEα−β,1−β−ψκtα−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,(18)∂σt,y∂y=ρ∂ut,y∂t=−ψ∫0tθ−1−βEα−β,−β−ψκθα−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫tt−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ψt−βEα−β,1−β−ψκtα−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.


## 2.4. The Initial Conditions and the Boundary Conditions

Whent≤0, ut=0. In the orifice, the flow velocity of the start-up flow at the center gradually increases with acceleration A, while the velocity near the pipe wall is always zero. In the gap, the flow velocity near the side of the moving plate increases with acceleration A, while the velocity near the stationary plate is always zero. So the initial conditions and the boundary conditions are as follows:(1)
The orifice:(19)ut,0=0,t>0,ut,L=0,t>0,ut,L2=At,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0.(2)
The gap:(20)ut,0=At,t>0,ut,L=0,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0,where L is the distance of the gap and the orifice and A is the acceleration.


## 2.5. Dimensionless Analysis

According to theΠ theorem, the dimensionless variables of the fractional K-BKZ model can be obtained as(21)u∗=ψκ1/α−βuA,y∗=L−1y,t∗=ψκ1/α−βt,σ∗=σLρA,a∗=aAL2κψ4/α−β,ϑ=ψL2ρκψ2−β/α−β.Substituting the dimensionless variables in (21) in (17) and (18), the fractional K-BKZ start-up flow dimensionless equation can be obtained (for brevity, the dimensionless mark “∗” is omitted):(22)σt,y=−ϑ∫0tθ−1−βEα−β,−β−θα−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ϑt−βEα−β,1−β−tα−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,(23)∂ut,y∂t=−ϑ∫0tθ−1−βEα−β,−β−θα−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫t−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ϑt−βEα−β,1−β−tα−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.The corresponding dimensionless initial conditions and boundary conditions of the orifice and gap are given as(24)ut,0=0,t>0,ut,1=0,t>0,ut,12=t,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0,(25)ut,0=t,t>0,ut,1=0,t>0,u0,y=0,y>0,∂u0,y∂t=0,y>0.


## 2.6. Numerical Solution

Numerical discretization was then applied to the fractional derivative using the finite difference method. The following were defined:tm=mΔt, m=0,1,2,…,K, 0≤k≤K−1; yn=nΔy, n=0,1,2,…,P. Δy=1/P is the space step, Δt=T/K is the time step. Applying the finite difference method to discretize the integrodifferential operators of (22) and (23), these can be approximated to(26)∫tatb∂uw,y∂ydw≈∑i=a+1buti,yn+1−uti,ynΔyΔt,∫tatb∂2uw,y∂y2dw≈∑i=a+1buti,yn+1−2uti,yn+uti,yn−1Δy2Δt.Letuji and σji denoteuti,yj and σti,yj. As Ftk=tk−1−βEα−β,−β−tkα−β and Htk=tk−βEα−β,1−β−tkα−β, the fractional K-BKZ start-up flow numerical model can be written as(27)σnm=−ϑ∑i=1mFtiΔt2/Δy∑j=m−i+1mun+1j−unj1+aΔt/Δy2∑j=m−i+1mun+1j−unj2+ϑHtmΔt/Δy∑j=1mun+1j−unj1+aΔt/Δy2∑j=1mun+1j−unj2,um+1,n−um,nΔt=−ϑ∑i=1mFtiΔt/Δy2∑j=m−i+1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21−2aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n2+ϑHtmΔt/Δy2∑j=1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=1muj,n+1−uj,n21−2aΔt/Δy2∑j=1muj,n+1−uj,n21+aΔt/Δy2∑j=1muj,n+1−uj,n2.


## 3. Development of the Classical K‐BKZ Start‐Up Flow Numerical Model

Larson expressed the relaxation modulus describing polymers in the K-BKZ model as [34](28)Gt=ct−g,where c and g are constants obtained through experiments, respectively.The determination of the values ofc and g is usually based on a grate quantity of experimental data. Therefore, it would be difficult to use equation (28) for accurate numerical calculation in this research due to the lack of experimental data.Consequently, Jaishankar et al. [14, 30, 35] introduced the fractional element model into the K-BKZ model and used fractional power law function to represent relaxation modulus as follows:(29)Gt=ψt−βΓ1−β,The memory function (5) can be written as(30)mt−s=ψβΓ1−βt−s−1−β.Thus, the K-BKZ model can be written as(31)σt=∫−∞tψβΓ1−βt−s−1−βεt,s1+aε2t,sds.Using the methods in Section2, the classical K-BKZ governing equation can be written as(32)σt,y=ψβΓ1−β∫0tθ−1−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ψΓ1−βt−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,∂σt,y∂y=ρ∂ut,y∂t=ψβΓ1−β∫0tθ−1−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫t−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ψΓ1−βt−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.and the classical K-BKZ start-up flow dimensionless equation can be obtained (for brevity, the dimensionless mark “∗” is omitted):(33)σt,y=ϑβΓ1−β∫0tθ−1−β∫t−θt∂uw,y/∂ydw1+a∫t−θt∂uw,y/∂ydw2dθ+ϑΓ1−βt−β∫0t∂uw,y/∂ydw1+a∫0t∂uw,y/∂ydw2,∂ut,y∂t=ϑ∫0tθ−1−βEα−β,−β−θα−β∫t−st∂2uw,y/∂y2dw1+a∫t−st∂uw,y/∂ydw21−2a∫t−st∂uw,y/∂ydw21+a∫t−st∂uw,y/∂ydw2ds+ϑt−β∫0t∂2uw,y/∂y2dw1+a∫0t∂uw,y/∂ydw21−2a∫0t∂uw,y/∂ydw21+a∫0t∂uw,y/∂ydw2.Then the classical K-BKZ start-up flow numerical model can be written as(34)σnm=ϑβΓ1−β∑i=1mti−1−βΔt2/Δy∑j=m−i+1mun+1j−unj1+aΔt/Δy2∑j=m−i+1mun+1j−unj2+ϑtm−βΓ1−βΔt/Δy∑j=1mun+1j−unj1+aΔt/Δy2∑j=1mun+1j−unj2,um+1,n−um,nΔt=ϑβΓ1−β∑i=1mti−1−βΔt/Δy2∑j=m−i+1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21−2aΔt/Δy2∑j=m−i+1muj,n+1−uj,n21+aΔt/Δy2∑j=m−i+1muj,n+1−uj,n2+ϑtm−βΓ1−βΔt/Δy2∑j=1muj,n+1−2uj,n+uj,n−11+aΔt/Δy2∑j=1muj,n+1−uj,n21−2aΔt/Δy2∑j=1muj,n+1−uj,n21+aΔt/Δy2∑j=1muj,n+1−uj,n2.The corresponding dimensionless initial conditions and boundary conditions of the orifice and gap are given as (24) and (25).


## 4. Development of the Fractional Maxwell Start-Up Flow Numerical Model

The fractional Maxwell model is proposed in [18, 30, 36, 37]:(35)σ+κψDtα−βσ=κDtαε.According to the methods in Section2, the fractional Maxwell governing equation can be written as(36)σt,y+κψDtα−βσt,y=κJt1−α∂ut,y∂y,ρ∂ut,y∂t+ρκψDtα−β+1ut,y=κJt1−α∂2ut,y∂2y.and the fractional Maxwell start-up flow dimensionless equation can be obtained (for brevity, the dimensionless mark “∗” is omitted):(37)σ+Dtα−βσt,y=ϑJt1−α∂ut,y∂y,(38)∂ut,y∂t+Dtα−β+1ut,y=ϑJt1−α∂2ut,y∂2y.Applying the finite difference method in the literature [38] to discretize and approximate (37) and (38), the fractional Maxwell start-up flow numerical model can be written as(39)σnk+1+Δtβ−αΓ2−α+βσnk+1+∑j=1kj+11−α+β−2j1−α+β+j−11−α+βσnk−j+1=ϑΔt1−αΓ2−α⋅Δy∑j=0kk−j+11−α−k−j1−αun+1j−unj,(40)unk+1−unkΔt+Δtβ−α−1Γ2−α+βunk+1−unk+∑j=1kj+11−α+β−2j1−α+β+j−11−α+βunk−j+1−unk−j=ϑΔt1−αΓ2−α⋅Δy2∑j=0kk−j+11−α−k−j1−αun+1j+un−1j−2unj,where (39) and (40) should satisfy 0<β<α<1. The corresponding dimensionless initial conditions and boundary conditions of the orifice and gap are given as (24) and (25).When there is only one spring-pot element andα=1, then the difference approximation to the (37) and (38) can be written as(41)σnm+1=ϑΔyun+1m−unm,(42)unm+1−unmΔt=ϑΔy2un+1m+un−1m−2unm,(41) and (42) are numerical expressions of the Newtonian fluid model.


## 5. Numerical Results Analysis

The time step was set as Δt = 0.001 and the space step was set as Δy = 0.025. The flow characteristics of the fractional K-BKZ numerical model, the classical K-BKZ numerical model, the fractional Maxwell numerical model, and the Newtonian fluid model were then compared.


### 5.1. Numerical Simulation of the Orifice

According to the initial conditions and the boundary conditions of (24), the numerical results of the flow velocity distribution of the start-up flow of the fractional K-BKZ numerical model, the classical K-BKZ numerical model, the fractional Maxwell numerical model, and the Newtonian fluid model in the orifice as the time increased have been shown in Figure 3. The basic parameters were α = 0.8, β = 0.3, ϑ = 0.2, a = 0.2. It can be seen from Figure 3 that the flow velocities of the four models in the orifice had a maximum value at y = 0.5, and the velocity distribution curves were symmetric along the line where y = 0.5, which conformed to the characteristics of the fluid flow in the pipe. By comparing the three-dimensional diagrams of the flow velocity distribution (Figures 3(a), 3(c), 3(e), and 3(g)), it can be seen that when the time and the maximum velocity growth rate had the same value, the flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model had been gradually developing parabolic shapes. These were more rounded than the fractional Maxwell numerical model and the Newtonian fluid model, and the curvatures of the curves changed greatly. However, the variation trends of the fractional Maxwell numerical model and the Newtonian fluid model were not obvious. By comparing the three-dimensional flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model (Figures 3(a) and 3(c)), it can be found that the curvature of the classical K-BKZ numerical model was larger than the fractional K-BKZ numerical model, and the former was more rounded than the latter. Figures 3(b), 3(d), 3(f), and 3(h) show the two-dimensional flow velocity distributions of the four models with respect to time. It can be seen from Figures 3(b) and 3(d) that, as the time increased, the flow velocity distribution curves of the fractional K-BKZ numerical model and the classical K-BKZ numerical model gradually changed from a flat shape to a parabolic shape which fitted the flow characteristics of the actual viscoelastic fluid. Moreover, regardless of how large the maximum flow velocity was, the distribution curve always maintained good nonlinear characteristics and the curve’s curvature was always more obvious than that of the fractional Maxwell numerical model and the Newtonian fluid model. This meant that the fluid that was simulated using the K-BKZ model (either the fractional K-BKZ numerical model or the classical K-BKZ numerical model) always maintained stronger viscoelasticity. However, the flow velocity distribution curve of the fractional Maxwell numerical model showed a flat shape at t = 0.5. With the increase of time and the maximum velocity, the distribution curves gradually resembled that of the Newtonian fluid model, which meant that in the interval where 0 < y < 0.5 or 0.5 < y < 1, the flow velocity distribution curve of the Newtonian fluid model had been developing a straight shape and the curvature of the curve was basically unchanged. Besides, when time was large, the curve’s curvature of the classical K-BKZ numerical model was larger than the fractional K-BKZ numerical model at the same time, which was consistent with the conclusions of the three-dimensional velocity distribution comparison.Figure 3
Flow velocity distributions of the four models in the orifice for different time values. (a) The three-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (b) The two-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (c) The three-dimensional flow velocity distribution of the classical K-BKZ numerical model. (d) The two-dimensional flow velocity distribution of the classical K-BKZ numerical model. (e) The three-dimensional flow velocity distribution of the fractional Maxwell numerical model. (f) The two-dimensional flow velocity distribution of the fractional Maxwell numerical model. (g) The three-dimensional flow velocity distribution of the Newtonian fluid model. (h) The two-dimensional flow velocity distribution of the Newtonian fluid model.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure4 shows the two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the orifice at low velocity (that is, at short time) as the time increased. As can be seen from Figure 4, with the increase of time, the flow velocity distribution of the fractional Maxwell numerical model gradually changed from a flat shape to a parabolic shape, while the curvature of the curve declined compared with the previous time at t = 0.8, and the fractional Maxwell numerical model gradually changed to the Newtonian fluid model. According to Figures 3 and 4, it can be found that the fractional Maxwell numerical model had good nonlinearity at low velocity or short time, while when the velocity was higher or the time was longer, it lost nonlinearity and viscoelasticity, which was consistent with the conclusions of studies [13] and [14].Figure 4
The two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the orifice for different time values.Figure5 shows the stresses variation of the fluid of four models at the center of the orifice. With the increase of the time, the stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly, but the stresses of the fractional K-BKZ numerical model and the classical K-BKZ numerical model increased slowly and tended to reach a steady state. The reason is that the K-BKZ model is able to effectively capture the phenomenon of shear thinning, while the Maxwell model and the Newtonian fluid model are not [14]. Since the classical K-BKZ numerical model has a damping function that converges to zero with time increasing, which can effectively describe the shear thinning phenomenon of viscoelastic fluid; therefore, the fractional K-BKZ numerical model derived from the model can also effectively capture the shear thinning phenomenon. Compared with the classical K-BKZ numerical model, the fractional K-BKZ numerical model employs the Mittag-Leffler function as the relaxation kernel instead of a power-law function. The Mittag-Leffler function is a generalized exponential function, which can improve the accuracy of the model [38], and it is more suitable to describe the relaxation process of polymer materials. The main component of the viscoelastic materials that were injected into the shock absorber was either methyl silicon rubber or methyl siloxane. These two materials belong to the class of high polymer materials and pseudoplastic fluids, which possesses a shear thinning phenomenon. Therefore, the fractional K-BKZ numerical model is able to capture the nonlinear viscoelastic properties of this viscoelastic material more accurately.Figure 5
Stress distributions at the center point of the orifice.


### 5.2. Numerical Simulation of the Gap

According to the initial conditions and the boundary conditions from (25), the flow velocity distributions of the start-up flow of four models in the gap of the shock absorber as the time increased have been shown in Figure 6. As can be seen from Figure 6, the flow velocity of the four models in the gap gradually decreased to zero from a certain value, which was consistent with the characteristics of the fluid flow between the plates. By comparing the three-dimensional diagrams of the flow velocity distribution of four models (Figures 6(a), 6(c), 6(e), and 6(g)), it can be seen that when the time and the maximum velocity growth rate had the same value, the fractional K-BKZ numerical model and the classical K-BKZ numerical model showed a fuller velocity distribution and greater curvature than the fractional Maxwell numerical model and the Newtonian fluid model. Figures 6(b), 6(d), 6(f), and 6(h) show the two-dimensional flow velocity distributions of the four models varying with time. It can be seen from these graphs that the velocity distribution of the Newtonian fluid model was linear while the velocity distribution of the fractional K-BKZ numerical model maintained good nonlinearity as the time increased. When the time was short or the velocity was small, the classical K-BKZ numerical model and the fractional Maxwell numerical model were both similar to the fractional K-BKZ numerical model; all can accurately simulate the viscoelasticity of fluid. However, when the time was long or the velocity was large, the fractional Maxwell numerical model gradually developed to the Newtonian fluid model, while the classical K-BKZ numerical model lost the accuracy of viscoelasticity model. It can be seen in Figure 6(d) that the curves of t = 2.5 and t = 3.5 overlapped each other when y > 0.5. The K-BKZ model with the relaxation kernel is described by power law function has the accuracy of the viscoelasticity model only at short times [14] and [26]. This phenomenon was proved in Figure 6, which is why the curve’s curvature and shape of the classical K-BKZ numerical model were larger and more rounded than the fractional K-BKZ numerical model at the same time, when time was large.Figure 6
Flow velocity distributions of the four models in the gap for different time values. (a) The three-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (b) The two-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (c) The three-dimensional flow velocity distribution of the classical K-BKZ numerical model. (d) The two-dimensional flow velocity distribution of the classical K-BKZ numerical model. (e) The three-dimensional flow velocity distribution of the fractional Maxwell numerical model. (f) The two-dimensional flow velocity distribution of the fractional Maxwell numerical model. (g) The three-dimensional flow velocity distribution of the Newtonian fluid model. (h) The two-dimensional flow velocity distribution of the Newtonian fluid model.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure7 shows the two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the gap at low velocity or at short time as the time increased. As can be seen from Figure 7, with the increase of time, the flow velocity distribution of the Newtonian fluid model had a trend to develop to a linear shape, while the fractional Maxwell numerical model had a good nonlinearity. By combining Figures 6 and 7, it can be seen that at low velocity or short time, the fractional Maxwell numerical model had good nonlinearity, while the velocity is higher or the time is longer, it lost the nonlinearity and the viscoelasticity. The results further verified the simulation results of the start-up flow in the orifice: the fractional K-BKZ numerical model always keeps good nonlinearity and viscoelasticity as the time or the speed increased. However, the classical K-BKZ numerical model and the fractional Maxwell numerical model will respectively lose the viscoelastic accuracy and nonlinearity at long times. Besides, the Newtonian fluid model has no viscoelastic characteristics.Figure 7
The two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the gap for different time values.For the start-up flow in the gap, the stress distribution curves of four models near the side of the moving plate (y = 0) have been shown in Figure 8. The stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly with the increasing time while the stresses of the fractional K-BKZ numerical model and the classical K-BKZ numerical model gradually became steady, which was consistent with the increasing trend of the orifice. When t = 10, the stress values of the start-up flow of the fractional Maxwell numerical model and the Newtonian fluid in the gap were 2.37 and 2.32, respectively, while the values of those two models in the orifice were 4.24 and 4.15 respectively. As can be seen from the values, the stresses of the start-up flow of the fractional Maxwell numerical model and the Newtonian fluid in the orifice all were 1.78 times that in the gap. The maximum velocity of the start-up flow in the orifice occurred at the center point (y = 0.5) while the maximum velocity of the start-up flow in the gap occurred near the side of the moving plate (y = 0), and there were also different boundary conditions. When the time step and the space step were the same, the calculated data of the velocity of the start-up flow in the orifice was half that of the data in the gap. As the stresses were calculated using the difference in velocity, and the stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly, the stresses of the start-up flow in the orifice were 1.78 times that in the gap. Although the velocity data could be kept consistent by reducing the number of steps, the adjustment of the steps would also cause a variation in the difference coefficient, so the result of the calculation would remain unchanged. The stresses of the start-up flow of the fractional K-BKZ numerical model and the classical K-BKZ numerical model in the orifice and the gap both were 0.50 and 0.48, respectively. Due to the decrease in the difference of the data, the stresses were increased although the variations were not large. The reason was that the fractional K-BKZ numerical model and the classical K-BKZ numerical model both had a damping function that converged to zero with the increase in the deformation, which made the models’ stress tend to be a steady state. As a result, the variation in the stress of the fractional K-BKZ numerical model was not large and the result was more stable; the fractional K-BKZ numerical model had a same characteristic with the K-BKZ model that it was able to effectively capture the shear thinning phenomenon.Figure 8
Stress distribution of four models near the side of the moving plate in the gap.


### 5.3. Influence Analyses of Parameters Change

Figures9–11 show the flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model in the gap for different parameters a, β, and ϑ, respectively. By combining Figures 9–11, we can find that when the time was short, the flow velocity distributions of the classical K-BKZ numerical model had a good consistency with the fractional K-BKZ numerical model, while when the time was large, the former was different from the latter. It can be seen from Figure 9 that as the time increased, no matter what the value of a was, the fractional K-BKZ numerical model always maintained good viscoelasticity, while when t = 2.5 and t = 3.5, the curves of the classical K-BKZ numerical model overlapped each other. The viscoelastic accuracy of the classical K-BKZ numerical model was low at long times, and the accuracy did not change with the values of parameter a increased. We can found from Figures 10 and 11 that the viscoelastic accuracies of the classical K-BKZ numerical model were strengthened with the increase of values of parameters β and ϑ, and the effect of parameter β was greater than that of parameter ϑ. Parameters β and ϑ were the greater, the accuracy of the classical K-BKZ numerical model was the better. Comparatively, the fractional K-BKZ numerical model had good accuracy.Figure 9
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parametera: (a) a = 0.1; (b) a = 0.3; (c) a = 0.5; (d) a = 0.7.
(a)
(b)
(c)
(d)Figure 10
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parameterβ: (a) β = 0.1; (b) β = 0.2; (c) β = 0.3; (d) β = 0.4.
(a)
(b)
(c)
(d)Figure 11
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parameterϑ: (a) ϑ = 0.1; (b) ϑ = 0.3; (c) ϑ = 0.5; (d) ϑ = 0.7.
(a)
(b)
(c)
(d)The flow velocity distributions of the fractional K-BKZ numerical model with change of parameters are shown in Figure12. With the increase of the parameter a, the nonlinearity and the flow velocity of the fractional K-BKZ numerical model gradually increased (as shown in Figure 7). The parameter a is related with the damping function [28]. The less the parameter a is, the more the damping function tends to zero and the larger the viscosity decreases. Because of the damping function tends to zero, the molecular chain entanglement and network structure of the viscoelastic material are mostly destroyed [28]. The parameters α and β were the fractional orders. It can be seen from Figure 7 that as α and β increased, the trends of the flow velocity decreased and the curvature of the flow velocity distribution curve gradually decreased, so the viscosity decreased. As shown in Figure 7 that the dependence of the parameter ζ on the flow velocity field had the same trend as that of α and β. The parameter ζ was smaller, the velocity curve curvature was larger, and the nonlinearity was stronger.Figure 12
Flow velocity distributions of the orifice of the fractional K-BKZ numerical model for different parameters. (a) Parameter a. (b) Parameterα. (c) Parameter β. (d) Parameter ϑ.
(a)
(b)
(c)
(d)


## 5.1. Numerical Simulation of the Orifice

According to the initial conditions and the boundary conditions of (24), the numerical results of the flow velocity distribution of the start-up flow of the fractional K-BKZ numerical model, the classical K-BKZ numerical model, the fractional Maxwell numerical model, and the Newtonian fluid model in the orifice as the time increased have been shown in Figure 3. The basic parameters were α = 0.8, β = 0.3, ϑ = 0.2, a = 0.2. It can be seen from Figure 3 that the flow velocities of the four models in the orifice had a maximum value at y = 0.5, and the velocity distribution curves were symmetric along the line where y = 0.5, which conformed to the characteristics of the fluid flow in the pipe. By comparing the three-dimensional diagrams of the flow velocity distribution (Figures 3(a), 3(c), 3(e), and 3(g)), it can be seen that when the time and the maximum velocity growth rate had the same value, the flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model had been gradually developing parabolic shapes. These were more rounded than the fractional Maxwell numerical model and the Newtonian fluid model, and the curvatures of the curves changed greatly. However, the variation trends of the fractional Maxwell numerical model and the Newtonian fluid model were not obvious. By comparing the three-dimensional flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model (Figures 3(a) and 3(c)), it can be found that the curvature of the classical K-BKZ numerical model was larger than the fractional K-BKZ numerical model, and the former was more rounded than the latter. Figures 3(b), 3(d), 3(f), and 3(h) show the two-dimensional flow velocity distributions of the four models with respect to time. It can be seen from Figures 3(b) and 3(d) that, as the time increased, the flow velocity distribution curves of the fractional K-BKZ numerical model and the classical K-BKZ numerical model gradually changed from a flat shape to a parabolic shape which fitted the flow characteristics of the actual viscoelastic fluid. Moreover, regardless of how large the maximum flow velocity was, the distribution curve always maintained good nonlinear characteristics and the curve’s curvature was always more obvious than that of the fractional Maxwell numerical model and the Newtonian fluid model. This meant that the fluid that was simulated using the K-BKZ model (either the fractional K-BKZ numerical model or the classical K-BKZ numerical model) always maintained stronger viscoelasticity. However, the flow velocity distribution curve of the fractional Maxwell numerical model showed a flat shape at t = 0.5. With the increase of time and the maximum velocity, the distribution curves gradually resembled that of the Newtonian fluid model, which meant that in the interval where 0 < y < 0.5 or 0.5 < y < 1, the flow velocity distribution curve of the Newtonian fluid model had been developing a straight shape and the curvature of the curve was basically unchanged. Besides, when time was large, the curve’s curvature of the classical K-BKZ numerical model was larger than the fractional K-BKZ numerical model at the same time, which was consistent with the conclusions of the three-dimensional velocity distribution comparison.Figure 3
Flow velocity distributions of the four models in the orifice for different time values. (a) The three-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (b) The two-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (c) The three-dimensional flow velocity distribution of the classical K-BKZ numerical model. (d) The two-dimensional flow velocity distribution of the classical K-BKZ numerical model. (e) The three-dimensional flow velocity distribution of the fractional Maxwell numerical model. (f) The two-dimensional flow velocity distribution of the fractional Maxwell numerical model. (g) The three-dimensional flow velocity distribution of the Newtonian fluid model. (h) The two-dimensional flow velocity distribution of the Newtonian fluid model.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure4 shows the two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the orifice at low velocity (that is, at short time) as the time increased. As can be seen from Figure 4, with the increase of time, the flow velocity distribution of the fractional Maxwell numerical model gradually changed from a flat shape to a parabolic shape, while the curvature of the curve declined compared with the previous time at t = 0.8, and the fractional Maxwell numerical model gradually changed to the Newtonian fluid model. According to Figures 3 and 4, it can be found that the fractional Maxwell numerical model had good nonlinearity at low velocity or short time, while when the velocity was higher or the time was longer, it lost nonlinearity and viscoelasticity, which was consistent with the conclusions of studies [13] and [14].Figure 4
The two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the orifice for different time values.Figure5 shows the stresses variation of the fluid of four models at the center of the orifice. With the increase of the time, the stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly, but the stresses of the fractional K-BKZ numerical model and the classical K-BKZ numerical model increased slowly and tended to reach a steady state. The reason is that the K-BKZ model is able to effectively capture the phenomenon of shear thinning, while the Maxwell model and the Newtonian fluid model are not [14]. Since the classical K-BKZ numerical model has a damping function that converges to zero with time increasing, which can effectively describe the shear thinning phenomenon of viscoelastic fluid; therefore, the fractional K-BKZ numerical model derived from the model can also effectively capture the shear thinning phenomenon. Compared with the classical K-BKZ numerical model, the fractional K-BKZ numerical model employs the Mittag-Leffler function as the relaxation kernel instead of a power-law function. The Mittag-Leffler function is a generalized exponential function, which can improve the accuracy of the model [38], and it is more suitable to describe the relaxation process of polymer materials. The main component of the viscoelastic materials that were injected into the shock absorber was either methyl silicon rubber or methyl siloxane. These two materials belong to the class of high polymer materials and pseudoplastic fluids, which possesses a shear thinning phenomenon. Therefore, the fractional K-BKZ numerical model is able to capture the nonlinear viscoelastic properties of this viscoelastic material more accurately.Figure 5
Stress distributions at the center point of the orifice.


## 5.2. Numerical Simulation of the Gap

According to the initial conditions and the boundary conditions from (25), the flow velocity distributions of the start-up flow of four models in the gap of the shock absorber as the time increased have been shown in Figure 6. As can be seen from Figure 6, the flow velocity of the four models in the gap gradually decreased to zero from a certain value, which was consistent with the characteristics of the fluid flow between the plates. By comparing the three-dimensional diagrams of the flow velocity distribution of four models (Figures 6(a), 6(c), 6(e), and 6(g)), it can be seen that when the time and the maximum velocity growth rate had the same value, the fractional K-BKZ numerical model and the classical K-BKZ numerical model showed a fuller velocity distribution and greater curvature than the fractional Maxwell numerical model and the Newtonian fluid model. Figures 6(b), 6(d), 6(f), and 6(h) show the two-dimensional flow velocity distributions of the four models varying with time. It can be seen from these graphs that the velocity distribution of the Newtonian fluid model was linear while the velocity distribution of the fractional K-BKZ numerical model maintained good nonlinearity as the time increased. When the time was short or the velocity was small, the classical K-BKZ numerical model and the fractional Maxwell numerical model were both similar to the fractional K-BKZ numerical model; all can accurately simulate the viscoelasticity of fluid. However, when the time was long or the velocity was large, the fractional Maxwell numerical model gradually developed to the Newtonian fluid model, while the classical K-BKZ numerical model lost the accuracy of viscoelasticity model. It can be seen in Figure 6(d) that the curves of t = 2.5 and t = 3.5 overlapped each other when y > 0.5. The K-BKZ model with the relaxation kernel is described by power law function has the accuracy of the viscoelasticity model only at short times [14] and [26]. This phenomenon was proved in Figure 6, which is why the curve’s curvature and shape of the classical K-BKZ numerical model were larger and more rounded than the fractional K-BKZ numerical model at the same time, when time was large.Figure 6
Flow velocity distributions of the four models in the gap for different time values. (a) The three-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (b) The two-dimensional flow velocity distribution of the fractional K-BKZ numerical model. (c) The three-dimensional flow velocity distribution of the classical K-BKZ numerical model. (d) The two-dimensional flow velocity distribution of the classical K-BKZ numerical model. (e) The three-dimensional flow velocity distribution of the fractional Maxwell numerical model. (f) The two-dimensional flow velocity distribution of the fractional Maxwell numerical model. (g) The three-dimensional flow velocity distribution of the Newtonian fluid model. (h) The two-dimensional flow velocity distribution of the Newtonian fluid model.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)Figure7 shows the two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the gap at low velocity or at short time as the time increased. As can be seen from Figure 7, with the increase of time, the flow velocity distribution of the Newtonian fluid model had a trend to develop to a linear shape, while the fractional Maxwell numerical model had a good nonlinearity. By combining Figures 6 and 7, it can be seen that at low velocity or short time, the fractional Maxwell numerical model had good nonlinearity, while the velocity is higher or the time is longer, it lost the nonlinearity and the viscoelasticity. The results further verified the simulation results of the start-up flow in the orifice: the fractional K-BKZ numerical model always keeps good nonlinearity and viscoelasticity as the time or the speed increased. However, the classical K-BKZ numerical model and the fractional Maxwell numerical model will respectively lose the viscoelastic accuracy and nonlinearity at long times. Besides, the Newtonian fluid model has no viscoelastic characteristics.Figure 7
The two-dimensional flow velocity distribution comparison of the fractional Maxwell numerical model and the Newtonian fluid model in the gap for different time values.For the start-up flow in the gap, the stress distribution curves of four models near the side of the moving plate (y = 0) have been shown in Figure 8. The stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly with the increasing time while the stresses of the fractional K-BKZ numerical model and the classical K-BKZ numerical model gradually became steady, which was consistent with the increasing trend of the orifice. When t = 10, the stress values of the start-up flow of the fractional Maxwell numerical model and the Newtonian fluid in the gap were 2.37 and 2.32, respectively, while the values of those two models in the orifice were 4.24 and 4.15 respectively. As can be seen from the values, the stresses of the start-up flow of the fractional Maxwell numerical model and the Newtonian fluid in the orifice all were 1.78 times that in the gap. The maximum velocity of the start-up flow in the orifice occurred at the center point (y = 0.5) while the maximum velocity of the start-up flow in the gap occurred near the side of the moving plate (y = 0), and there were also different boundary conditions. When the time step and the space step were the same, the calculated data of the velocity of the start-up flow in the orifice was half that of the data in the gap. As the stresses were calculated using the difference in velocity, and the stresses of the fractional Maxwell numerical model and the Newtonian fluid model increased linearly, the stresses of the start-up flow in the orifice were 1.78 times that in the gap. Although the velocity data could be kept consistent by reducing the number of steps, the adjustment of the steps would also cause a variation in the difference coefficient, so the result of the calculation would remain unchanged. The stresses of the start-up flow of the fractional K-BKZ numerical model and the classical K-BKZ numerical model in the orifice and the gap both were 0.50 and 0.48, respectively. Due to the decrease in the difference of the data, the stresses were increased although the variations were not large. The reason was that the fractional K-BKZ numerical model and the classical K-BKZ numerical model both had a damping function that converged to zero with the increase in the deformation, which made the models’ stress tend to be a steady state. As a result, the variation in the stress of the fractional K-BKZ numerical model was not large and the result was more stable; the fractional K-BKZ numerical model had a same characteristic with the K-BKZ model that it was able to effectively capture the shear thinning phenomenon.Figure 8
Stress distribution of four models near the side of the moving plate in the gap.


## 5.3. Influence Analyses of Parameters Change

Figures9–11 show the flow velocity distributions of the fractional K-BKZ numerical model and the classical K-BKZ numerical model in the gap for different parameters a, β, and ϑ, respectively. By combining Figures 9–11, we can find that when the time was short, the flow velocity distributions of the classical K-BKZ numerical model had a good consistency with the fractional K-BKZ numerical model, while when the time was large, the former was different from the latter. It can be seen from Figure 9 that as the time increased, no matter what the value of a was, the fractional K-BKZ numerical model always maintained good viscoelasticity, while when t = 2.5 and t = 3.5, the curves of the classical K-BKZ numerical model overlapped each other. The viscoelastic accuracy of the classical K-BKZ numerical model was low at long times, and the accuracy did not change with the values of parameter a increased. We can found from Figures 10 and 11 that the viscoelastic accuracies of the classical K-BKZ numerical model were strengthened with the increase of values of parameters β and ϑ, and the effect of parameter β was greater than that of parameter ϑ. Parameters β and ϑ were the greater, the accuracy of the classical K-BKZ numerical model was the better. Comparatively, the fractional K-BKZ numerical model had good accuracy.Figure 9
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parametera: (a) a = 0.1; (b) a = 0.3; (c) a = 0.5; (d) a = 0.7.
(a)
(b)
(c)
(d)Figure 10
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parameterβ: (a) β = 0.1; (b) β = 0.2; (c) β = 0.3; (d) β = 0.4.
(a)
(b)
(c)
(d)Figure 11
Flow velocity distributions of the fractional K-BKZ numerical model and classical K-BKZ numerical model in the gap for parameterϑ: (a) ϑ = 0.1; (b) ϑ = 0.3; (c) ϑ = 0.5; (d) ϑ = 0.7.
(a)
(b)
(c)
(d)The flow velocity distributions of the fractional K-BKZ numerical model with change of parameters are shown in Figure12. With the increase of the parameter a, the nonlinearity and the flow velocity of the fractional K-BKZ numerical model gradually increased (as shown in Figure 7). The parameter a is related with the damping function [28]. The less the parameter a is, the more the damping function tends to zero and the larger the viscosity decreases. Because of the damping function tends to zero, the molecular chain entanglement and network structure of the viscoelastic material are mostly destroyed [28]. The parameters α and β were the fractional orders. It can be seen from Figure 7 that as α and β increased, the trends of the flow velocity decreased and the curvature of the flow velocity distribution curve gradually decreased, so the viscosity decreased. As shown in Figure 7 that the dependence of the parameter ζ on the flow velocity field had the same trend as that of α and β. The parameter ζ was smaller, the velocity curve curvature was larger, and the nonlinearity was stronger.Figure 12
Flow velocity distributions of the orifice of the fractional K-BKZ numerical model for different parameters. (a) Parameter a. (b) Parameterα. (c) Parameter β. (d) Parameter ϑ.
(a)
(b)
(c)
(d)


## 6. Test Verification

The test system consisted of an AC servo motor, a machine controller, a shock absorber made from a glass tube, a high-speed camera, image acquisition software, and a spotlight (as shown in Figure13(a)). To facilitate the capture of the start-up flow process of the viscoelastic fluid, the plastic tube was used as the piston and the acrylic glass tube was used as the shock absorber’s chamber (as shown in Figure 13(b)). The test device was designed according to the actual sizes of the shock absorber’s piston and orifice. The chamber’s inner diameter was 20 mm, the thickness of the piston was 20 mm, and the diameter of orifice was 1.5 mm. In order to facilitate the process of capturing the flow of the orifice, the gap between the piston and the chamber was removed and a sealing ring was installed on the head of the plastic piston tube in order to prevent the viscoelastic fluid from squeezing into the gap during motion.Figure 13
The test system for the start-up flow measure. (a) The test system layout for the start-up flow measurement. (b) The partial enlarged picture of the glass tube shock absorber.
(a)
(b)A Mitsubishi AC servo motor was used as the motion source to drive the piston. The motor controller was used to adjust the rotational speed to keep the piston at a constant speed or acceleration. A Fastec HiSpec 5 high-speed camera and the accompanying image acquisition software were used to capture the flow condition of the viscoelastic fluid in the orifice at the start-up of the piston. In the shooting process, the spotlight was used as an auxiliary light source to enhance the camera’s image capture.A viscoelastic material, with methylsiloxane as the main component with a dynamic viscosity of 5 × 105 cSt and a density of 1.1 g/cm3, was injected into the chamber of the shock absorber in the glass tube. The motor controller was adjusted to start the piston from 0 mm/s, and the acceleration was maintained at 1.0 mm/s2. The experiment was carried out at ambient temperature. The acquisition frequency of the high-speed camera was 500 frames per second, so it could record the flow conditions of the viscoelastic fluid at different moments, as shown in Figure 14. The black region in the figure shows the viscoelastic material while the white region is the piston of the plastic tube. During the experiment, the auxiliary light source was on the left side of the shock absorber in the glass tube so that the reflective point was easily visible at the left half when the high-speed camera was shooting. As a result, some areas on the left side of the test results were missing and became white due to overexposure. Due to the small acceleration, the fluid velocity was also small in the first 0.6 s, and the change in the viscoelastic material was not obvious.Figure 14
The test results of the start-up flow in the orifice: (a)t = 0 s; (b) t = 0.6 s; (c) t = 0.8 s; (d) t = 1.0 s; (e) t = 1.2 s; (f) t = 1.4 s.
(a)
(b)
(c)
(d)
(e)
(f)According to the test results, the flow velocity was the largest at the orifice’s center and was smaller around the orifice’s edge. As the time increased and the piston’s speed increased, the shape of the viscoelastic material gradually changed from a flat shape to a parabolic shape. The start-up flow in the orifice always maintained stronger viscoelasticity during the flow process.The velocity-time variation of the start-up flow in the orifice was obtained by substituting the model’s parameters listed in Table1 into the fractional K-BKZ numerical model. The fluid level distribution curves of the fractional K-BKZ numerical model with increasing time were obtained by integration of the velocity, as shown in Figure 15. The shape of the fluid level distribution, simulated by the fractional K-BKZ numerical model, gradually changed from a flat shape to a parabolic shape as the time increased, which was basically the same as the trend of the test results.Table 1
Parameters of the fractional K-BKZ numerical model.
Parameter
Physical description
Numerical value
A
Acceleration
1.0 mm/s2
Ρ
Density
1.1 kg/m3
L
Orifice distance
1.5 mm
Λ
Relaxation time [39]
1.94 s
G
t
Relaxation modulus [39]
132.83 Pa
Α
Fractional order [39]
0.35
κ
Quasiproperty
1036.2 Pa·sα
Β
Fractional order [39]
0.8
ψ
Quasiproperty
231.9 Pa·sβ
A
A parameter related with the damping function [28]
0.3
The quasipropertiesκ and ψ can be approximated by the relaxation modulus Gt and the fractional orders α and β, which can be written as [14](43)
κ
≈
R
t
λ
−
α
Γ
1
−
α
,
ψ
≈
R
t
λ
−
β
Γ
1
−
β
.Figure 15
The fluid level distribution of the fractional K-BKZ numerical model in the orifice.The measurement method of the contact angle used in the literature [40] was applied in this paper, and the test results were combined for modification. As an example from the test results of the viscoelastic material’s fluid level at t = 1.0 s, two parameters of the maximum fluid level H and the contact angle γ were proposed, as shown in Figure 16. The maximum fluid level was taken as the vertical distance from the horizontal line of the lowest fluid level (i.e., the initial position of the viscoelastic fluid) to the highest point at the chosen time. This could verify the consistency of the numerical model results and the test results at the same time. An auxiliary tangent was drawn from the lowest point of the fluid level. The contact angle was the included angle between the auxiliary tangent and the horizontal line of the lowest fluid level, which was mainly used to investigate the viscoelasticity of the fluid. Some areas on the left side of the test results become white due to overexposure, so the contact angles were subjected to measurement on the right. According to the scale length, the maximum fluid level and the contact angle of the simulation results (as shown in Figure 15) and the test results (as shown in Figure 14) were then measured; the results have been shown in Table 2, and the trends have been shown in Figures 17 and 18.Figure 16
Viscoelastic fluid level measurement method (L is the orifice distance, H is the maximum fluid level, and γ is the contact angle).Table 2
Comparison of the maximum fluid level and the contact angle.
Test results
H (mm)
γ (°)
The results of the fractional K-BKZ numerical model
H (mm)
Γ (°)
The results of the Newtonian numerical model
H (mm)
γ (°)
t = 0.6 s
0.1889
44.4
t = 0.6 s
0.1800
45.6
t = 0.6 s
0.1800
29.5
t = 0.8 s
0.3306
55.8
t = 0.8 s
0.3173
58.2
t = 0.8 s
0.3173
41.6
t = 0.1 s
0.4665
62.3
t = 0.1 s
0.4889
65.2
t = 0.1 s
0.4889
52.8
t = 1.2 s
0.6995
70.8
t = 1.2 s
0.6971
72.4
t = 1.2 s
0.6971
61.3
t = 1.4 s
0.9075
74.4
t = 1.4 s
0.9402
76.7
t = 1.4 s
0.9402
67.7Figure 17
Comparison of the maximum fluid level.Figure 18
Comparison of the contact angle.As has been shown in Figure17 and Table 2, the maximum liquid level of the two models was basically consistent with the tests, indicating that the models and the tests were well matched at the same time. It can be seen from Figure 18 and Table 2 that, as the time increased, the contact angles of the fractional K-BKZ numerical model were in good agreement with the tests, while those of the Newtonian fluid model were greatly different (the average error of the contact angle was 19.3%), which meant that the fractional K-BKZ numerical model showed strong viscoelasticity, while the Newtonian fluid model showed no viscoelasticity. The maximum error of the maximum fluid level was 4.79%, and the maximum error of the contact angle was 4.65% (as shown in Figure 19). In summary, the fractional K-BKZ numerical model was able to accurately describe the viscoelastic characteristics of the viscoelastic material and accurately simulate the actual flow of the start-up flow in the shock absorber.Figure 19
The error bars of the fractional K-BKZ numerical model.


## 7. Conclusions

In this paper, the fractional K-BKZ numerical model was proposed in order to study the start-up flows of the viscoelastic material in a shock absorber. According to the actual flow conditions of a viscoelastic material in the orifice and the gap of a shock absorber, the flow behavior was simplified to two forms: the accelerated flow of the viscoelastic fluid in a pipe and the flow of the viscoelastic fluid between two parallel plates with one plate accelerating and the other at rest.The fractional K-BKZ numerical model has been established by adopting the framework of the classical K-BKZ model and the relaxation modulus of the fractional Maxwell model with quasiproperties. According to the actual working environment, the initial conditions and the boundary conditions were set. Compared with the results given by the classical K-BKZ numerical model, the fractional Maxwell numerical model, and the Newtonian fluid model, the variations in the flow velocity and stress were analyzed. Finally, the test system for the start-up flow was designed to capture the flow process of the viscoelastic material in the orifice during the start-up flow. The maximum fluid level and the contact angle were defined and measured. The accuracy of the fractional K-BKZ numerical model was verified by comparing the results from the model and the tests.The numerical results and the test results proved that(1)
When the viscoelastic fluid had an accelerated flow in the pipe or between the two parallel plates, the flow velocity distribution simulated by the fractional K-BKZ numerical model always maintained good viscoelasticity. The classical K-BKZ numerical model and the fractional Maxwell numerical model lost the viscoelasticity and the accuracy of viscoelasticity model, respectively, when the time was long or the velocity was large. The Newtonian fluid model displayed no viscoelasticity.(2)
The fluid stresses increase for the fractional K-BKZ numerical model and the classical K-BKZ numerical model tended to be a steady state, which could effectively capture the shear-thinning phenomenon.(3)
Although the model accuracy did not change with the increase of the value of parametera, the higher the values of parameter β and ϑ were, the better the accuracy of the classical K-BKZ numerical model would be. With the increase of the parameters α, β, and ϑ, the velocity curve curvatures and the nonlinearity of the fractional K-BKZ numerical model gradually decreased, while the influence of the parameter a was opposite to that of other parameters.(4)
The shape of the start-up flow of the viscoelastic material in the orifice developed from a flat shape to a parabolic shape. The comparison between the tests and the fractional K-BKZ numerical model proved that they had basically identical variation trends.(5)
The results of the maximum fluid level and the contact angle of the fractional K-BKZ numerical model were in good agreement with the test results. The maximum errors were 4.79% and 4.65%, respectively, which proved the accuracy of the model. Thus, the fractional K-BKZ numerical model was able to accurately simulate the start-up flow of the viscoelastic material in a shock absorber, and the fractional K-BKZ numerical model always maintained strong viscoelasticity.



---

*Source: 5412361-2020-02-22.xml*