A Three-Dimensional Velocity Field Related to a Generalized Third-Grade Fluid Model

: In this work, we propose a new three-dimensional constitutive equation related to a third-grade fluid. This proposal is based on experimental work for which the viscosity term and the terms related to viscoelasticity may depend on the shear rate, in accordance with a power-law type model. The numerical implementation of this fluid model is rather demanding in terms of computational calculation and, in this sense, we use the Cosserat theory related to fluid dynamics, which makes the transition from the three-dimensional fluid model to a one-dimensional fluid model for a specific geometry under study which, in this case, is a straight tube with constant circular cross-section. Based on this approximation theory, the one-dimensional fluid model is solved by assuming an ordinary differential equation involving: an unsteady mean pressure gradient; an unsteady volume flow rate; the Womersley number; and the viscosity and viscoelasticity parameters. Consequently, for specific data, and using the Runge–Kutta method, we can obtain the solution for the unsteady volume flow rate and we can present simulations to the three-dimensional velocity field.


Introduction
This article is based on the works of Caulk and Naghdi [1] and Carapau et al. [2,3] applied to a new three-dimensional model related to a third-grade non-Newtonian fluid.Here, we consider the viscosity and the viscoelasticity parameters to be dependent on the shear rate.In order to develop this new three-dimensional fluid model, we start with the standard constitutive equation for a third-grade non-Newtonian fluid (see Truesdell and Noll [4]), which is given by In Equation (1), the tensors A 1 , A 2 and A 3 are given by (see Rivlin and Ericksen [5]) and In Equations ( 2) and (3), the vector ϑ = ϑ(x, t) is the three-dimensional velocity field (the vector x = (x 1 , x 2 , z = x 3 ) is the vector space coordinates and the parameter t is the time variable).∇ϑ is the velocity gradient, and ∇ϑ T denotes the transpose of ∇ϑ.In Equations (2) 2 and (3), the expression d dt (•) is the material time derivative, defined by: Furthermore, in Equation (1), "tr" denotes the trace operator, µ is the viscosity constant, and −pI is the undetermined part of the stress due to incompressibility, where p is the pressure and α 1 , α 2 , β 1 , β 2 and β 3 are the viscoelasticity constants, also called the normal stress coefficients.
The stability and thermodynamics of constitutive Equation (1) were studied in detail by Fosdick and Rajagopal [6].Based on that work, the solution stability is ensured under the following conditions: Using the condition (5) in Equation (1), we obtain the following constitutive equation: In Equation (6), the different material parameters µ, α 1 , α 2 and β 3 may depend on several factors, for example temperature, pressure, or shear rate (see Truesdell and Noll [4]).In this work, we consider that the viscosity (µ) and the viscoelasticity terms (α 1 , α 2 , β 3 ) depend on the shear rate.Therefore, we present a new constitutive equation: where the positive function (the set R + is the set of positive real numbers), is the shear-dependent viscoelasticity function, γ the traditional scalar measure of the rate of shear given by | γ| = √ 2D : D where ":" denotes the tensor product, and the tensor D is given by Considering the experimental work of Beracea et al. [7], Mall-Gleissle et al. [8], and Tao et al. [9], related to polymers, suspensions, and liquid crystals, respectively, we can conclude that there is a variation in the viscosity and in the viscoelasticity terms associated with the shear rate, with this variation being of the power-law type.Therefore, we consider the positive power-law function (8) in Equation (7), defined as follows: where n > 0 is the flow index.By Equation (9), and for different values of the flow index n, we obtain two distinct and relevant situations in the current study, i.e., shear-thinning (or pseudoplastic) viscoelastic fluid situation and shear-thickening (or dilatant) fluid situation.Therefore, in condition (9) and we obtain the shear-thinning (or pseudoplastic) viscoelastic fluid case, i.e., the viscoelasticity decreases when we increase the shear rate; see Figure 1a.In the same sense, considering n > 1 in condition (9), the result is lim and we obtain the shear-thickening (or dilatant) viscoelastic fluid situation, i.e., the viscoelasticity increases when we increase the shear rate; see Figure 1b.Finally, considering n = 1 in condition (9) applied to Equation (7), we recover constitutive Equation (6) for a third-grade non-Newtonian fluid with the stability condition (5).Furthermore, considering β 3 = 0, Equation (6) becomes the standard constitutive equation for a second-grade non-Newtonian fluid; see Coleman and Noll [10].Also, if α 1 = α 2 = 0 and β 3 = 0, Equation (6) becomes the constitutive equation for a Newtonian fluid.Based on the work of Chien et al. [11], we can consider that Equation (7), with the condition (9) for n < 1, may be relevant for the study of blood flow in small vessels, where the action of the shear rate associated with the phenomena of aggregation and deformability of red blood cells produce variational effects on the viscosity and viscoelasticity.The study in terms of classical numerical implementation (for example, by finite element methods) for a three-dimensional model associated with an incompressible fluid that follows constitutive Equation (7) with the condition (9) in a tube of a circular crosssection with a variable radius is, in most situations, an unfeasible study in computational terms.In this sense, to get around the difficulty related to the space dimensions (x 1 , x 2 , x 3 ) and time t, we will use the Cosserat theory related to fluid dynamics, as developed by Caulk and Naghdi [1].This theory, under very specific conditions, makes the transition from a three-dimensional model (with three space variables and one time variable) to a onedimensional model (one space variable and one time variable), making the computational implementation more accessible.Furthermore, in terms of numerical implementation, we will only consider the simplest geometry, i.e., a straight tube with a circular crosssection of constant radius throughout the flow.Using the Cosserat theory associated with fluid dynamics (see [1,[12][13][14]), in which the three-dimensional velocity field of the fluid is approximated in a very specific way, it is possible, by integration of the linear momentum equation, to transform the three-dimensional model of the fluid into a one-dimensional model.This procedure allows us to obtain a partial (or ordinary) differential equation.In most cases, the difficulty of this differential equation is related to the geometry under study.Therefore, considering [1], let us assume that the unsteady three-dimensional velocity field ϑ(x, t) = ϑ i (x, t)e i , is approximated by where Latin indices take the values 1, 2, 3, Greek indices the values 1, 2, and we use the convention of summing over repeated indices.Moreover where the function v denotes the velocity along the symmetric axis z (axis relative to the flow) at variable time t, and x θ 1 . . .x θ N are the polynomial weighting functions with order k; this number is the number of directors.The vectors W θ 1 ...θ N are the director velocities associated with specific physical characteristics of the fluid (see [1]), and e i are the unit basis vectors.
In Section 2, we consider a three-dimensional model for an incompressible fluid that follows the Equation (7) with conditions (9), where the unsteady three-dimensional velocity field is approximated by the expression (10) with nine-directors, i.e., k = 3 in (10) (see [1]).Therefore, we obtain a one-dimensional fluid model, making it possible to obtain for specific flow regime data, applying the Runge-Kutta method to the solution for the unsteady volume flow rate.Consequently, we can present simulations to the unsteady three-dimensional velocity field, including an analysis on perturbed flows.

Proposed Flow Model
Considering the arguments presented in the previous section, we propose a new three-dimensional model for a non-Newtonian incompressible fluid, where viscosity and elasticity vary depending on the shear rate, with this variation being of the power-law type.Therefore, for our proposed generalized third-grade fluid model, the equations of motion stating the conservation of linear momentum and mass are given in Ω × (0, T) by where ρ is the constant fluid density and B is the external body forces.In our study, we neglect the external forces; in this sense, B ≡ 0. The balance of linear momentum is given by Equation ( 12) 1 , and the incompressibility condition of the fluid is given by ( 12) 2 .The first expression of Equation ( 12) 3 is the proposed constitutive Equation (7) with the given condition (9).Finally, the second expression of Equation ( 12) 3 is the stress vector on the surface, given by t w , where η = η i e i denotes the outward unit normal vector.In this section, we consider the general straight tube geometry Ω contained in R 3 (see Figure 2) of the circular cross-section, where the scalar function ϕ(z, t) is related to the variable radius of the tube by is the general geometry fluid domain Ω; functions p e and τ 1 , τ 2 are the normal and tangential components of the surface traction vector, respectively, where τ 1 is the wall shear stress.Also, Γ 1 is the proximal cross-section, Γ 2 is the distal cross-section, and Γ w is the lateral wall of the tube.
Then, considering the scalar function ϕ(z, t), the components η i of the outward unit normal vector η to the domain Ω are given by where a subscripted variable denotes partial differentiation.Applying condition (4) to Equation (13), then, on the boundary of the domain Ω, we have: Now, using the unsteady three-dimensional velocity approach (10) with k = 3, based on the work [1], we obtain where ξ(z, t), ω(z, t), γ(z, t), σ(z, t) and ψ(z, t) are physical scalar functions with the following meaning: ξ and σ are related to transverse elongation; γ is related to transverse shearing motion, while ω and ψ are related to rotational motion (also called swirling motion).Using velocity Equation (16), the lateral boundary (15) can be rewritten as and the incompressibility condition (12) 2 reduces to Consequently, for Equation (18) to be valid for any point in the fluid domain, we have to impose the following conditions: In order to simplify our study, we neglect the time variable in the scalar function i.e., we consider the fluid flowing in a straight tube with a rigid wall.Also imposing a no-slip condition on the boundary, velocity Equation (16) is identically equal to zero on surface (13).Therefore, we obtain from (16): Furthermore, using condition (20), kinematic Equation (17) is satisfied, and we can rewrite the incompressibility condition (19) 2 as Throughout our work, we are interested in understanding the behavior of the unsteady volume flow rate; in this sense, we present the definition to the volume flow rate Q, i.e., where S(z, t) is an arbitrary cross-section of the tube geometry.Now, using the definition for unsteady volume flow rate (23), the component of the three-dimensional velocity field (16) related to e 3 and conditions (21) 3 , (22), we conclude that the unsteady volume flow rate Q does not depend on the space variable z, and is given by Starting from the balance of linear momentum (12) 1 without external body forces, we impose the following integral assumptions, where N = 1, 2, 3 (see [1]): S(z,t) Now, applying the divergence theorem and integration by parts, the integral conditions ( 25) and ( 26) can be reduced to the equations: The forces terms n, h θ 1 ...θ N and m θ 1 ...θ N are given by and The inertia terms a and b θ 1 ...θ N are given by and the surface traction terms f and l θ 1 ...θ N , are given by Considering work [1], the stress vector t w , given by (12) 3 on the lateral surface related to Equation (32), can be rewritten in terms of the outward unit normal vector η, and by the normal and tangential components p e and τ 1 , τ 2 of the surface traction, i.e., where the scalar function τ 1 is the wall shear stress.Moreover, taking into account a flow without rotation (i.e., ψ ≡ ω ≡ 0 in Equation ( 16)), conditions (19) 1 , (20), ( 21) 1,3 , and the unsteady volume flow rate (24), then the unsteady three-dimensional velocity field (16) reduces to The model (12), with n = 1 on the power-law function given in condition (12) 3 reduces to a three-dimensional model for a non-Newtonian incompressible fluid of the third grade, where the constitutive equation is given by (6).The third-grade fluid model has been studied for different frameworks (see, e.g., [15][16][17][18][19][20]) and, in particular, by the Cosserat theory (see, e.g., [2,3]).Analogously, we can reduce Model (12), with n = 1 and β 3 = 0 in condition (12) 3 , to the three-dimensional model for a non-Newtonian incompressible fluid of the second grade.The second-grade fluid model has been studied for different frameworks (see, e.g., [21][22][23][24][25][26][27][28]) and, in particular, by the Cosserat theory (see, e.g., [29][30][31]).Finally, Model (12) can be reduced to a three-dimensional model for a Newtonian incompressible fluid, with n = 1 and α 1 = α 2 = β 3 = 0 in condition (12) 3 .The Newtonian fluid model has already been studied by the Cosserat theory, considering the viscosity constant and non-constant (see, e.g., [32,33]), respectively.Remark 1.For concrete flow regimes of Newtonian and non-Newtonian fluids in specific geometries, it is possible to find exact solutions in the scientific literature.For these cases, this Cosserat theory proved to be an approximation theory to take into account.More specifically, the Cosserat theory was validated for Newtonian fluids in the work of Caulk and Naghdi [1] and Robertson and Sequeira [32].Also, the theory was validated for non-Newtonian fluids in the work of Carapau et al. [30,31,33].

One-Dimensional Results
Due to the complexity of the model under study, we only present results for the simplest geometry case, i.e., the scalar function ϕ given by ( 20) is constant.Therefore, the unsteady three-dimensional velocity field (34) reduces to ϑ(x, t) = 2Q(t) and the unsteady stress vector (33) reduces to Now, under the conditions of the system (12) with ϕ constant, and applying the velocity field (35) and the stress vector (36) to Equations ( 29)- (32), we can calculate the results related to the quantities n, k α 1 ...α n , m α 1 ...α n , a, b α 1 ...α N and f , l α 1 ...α N .Therefore, using the solutions related to the forces, the inertia and surface traction terms in system Equations ( 27) and ( 28) with the average pressure (the quantity A(z, t) is the area of the section S(z, t)), we obtain the unsteady average pressure differential equation and the unsteady wall shear stress equation Integrating differential Equation (38) over the cross-section of tube interval [z 1 , z 2 ], we obtain the unsteady mean pressure gradient equation, given by Let us consider the following dimensionless variables (the characteristic frequency for unsteady flow is given by the parameter ω 0 ): (42) Now, substituting this new variables (41) and (42) into Equations ( 39) and (40), we obtain the unsteady nondimensional mean pressure gradient equation Ĝ Moreover, we obtain the unsteady nondimensional wall shear stress equation In Equations (43 and (44), the constant n is the power index, α1 , β3 are viscoelastic coefficients, and W 0 is the Womersley number, given by which reflects the pulsatility of the flow frequency in relation to viscous effects, which is an unsteady phenomenon (see [34]).Finally, using the dimensionless variables in Equation ( 35), we obtain the unsteady nondimensional three-dimensional velocity field equation At this stage, it is important to mention that, for flow index n = 1 in Equations ( 43) and ( 44), we recover the results obtained by Carapau and Correia [2].In next section, for specific flow regimes data, we present numerical illustrations to the unsteady volume flow rate and the unsteady three-dimensional velocity field.

Numerical Simulations
For specific flow regimes, we could obtain numerical results for both Equation ( 43) and ( 44), but since we intend to present the results for the unsteady three-dimensional velocity field (46), we will focus our numerical simulations on Equation (43).Therefore, considering Equation (43) for specific flow regime data, with constant and non-constant mean pressure gradient, and applying the Runge-Kutta method, we obtain the unsteady solution for the volume flow rate.Consequently, we can present numerical simulations to solve the unsteady three-dimensional velocity field (46), for both shear-thinning and shear-thickening viscoelastic fluids, considering a straight tube with a constant circular cross-section.

Constant Mean Pressure Gradient
In this case, let us consider a specific flow regime, i.e., Ĝ = 1 in Equation (43), Q(0) = 1, α1 = 1 and β3 = 1.Applying the Runge-Kutta method to solve Equation (43), we can present the behavior of the unsteady volume flow rate Q for various values of the Womersley number and power index n.Figure 3 shows the behavior of the shear-thinning viscoelastic fluid, and Figure 4 shows the behavior of the shear-thickening viscoelastic fluid.
Based on the solutions presented in Figures 3 and 4, we can conclude that after the initial transient phase the unsteady volume flow (43) converges in time to the stationary solution.The convergence is faster for the case of shear-thickening viscoelastic fluid (see Figure 4) than for the shear-thinning case.Finally, we illustrate in Figures 5 and 6 the behavior of the unsteady three-dimensional velocity field (46), related to shear-thinning and shear-thickening viscoelastic fluids.The behavior of the solution presented for three-dimensional velocity field θ is associated with the behavior of the volume flow rate solution Q during and after the initial transient phase.We can remark that the shown behavior of the unsteady volume flow rate solution Q has an influence on the choice of time parameters in the illustrations presented in Figures 5 and 6.

Non-Constant Mean Pressure Gradient
Given the following equation for the non-constant mean pressure gradient Ĝ we can verify that it presents a strong variation during the initial transient phase and smaller variation in time after the initial transient phase; see Figure 7. Let us consider a specific flow regime, i.e., Q(0) = 1, α1 = 1, β3 = 1 and W 0 = (0.2, 0.4, 0.6, 0.8).Using the Runge-Kutta method for solving ordinary differential Equation (43), we can observe the behavior of the unsteady volume flow rate Q for different power index values, related to shear-thinning and shear-thickening viscoelastic fluids; see Figures 8 and 9, respectively.Analyzing the behavior of the volume flow rate illustrated in Figures 8 and 9, we can conclude that, over time, the unsteady solution converges to the steady solution, as in the case where the mean pressure gradient was constant.Also, in these conditions, the convergence is faster for the case of shear-thickening viscoelastic fluid (see Figure 9) than for the shear-thinning case.Finally, considering the data solution for the unsteady volume flow rate on Equation (46), we can illustrate the behavior of the three-dimensional velocity field (46), related to shear-thinning and shear-thickening viscoelastic fluids; see Figures 10 and 11, respectively.As in the case of a constant mean pressure gradient, the behavior of the solution presented for the three-dimensional velocity field θ is associated with the behavior of the volume flow rate solution Q during and after the initial transient phase.Again, the choice of time parameters in Figures 10 and 11

Perturbed Flows
When considering the new constitutive Equation (7), we lose the guarantee of stability solution mentioned in condition (5).Therefore, in this section, we intend to take a first approach to the study of the stability solution related to the unsteady volume flow rate Q, obtained by the proposed flow model (12).In this sense, we only study the stability solution, where Q(0) = 1, W 0 = 0.8, α1 = 1 and β3 = 1 are fixed, for different power index parameters.Now, let us consider the perturbation function, given by Ĝ± where ϵ > 0 is the magnitude perturbation and Q± ϵ is the perturbed volume flow rate related to the perturbation Ĝ± ϵ .From the previous section, we know that, over time, the solution to the unsteady volume flow rate Q for constant and non-constant mean pressure gradient converges to the steady solution.As a result, assuming Qt ( t) = 0 in Equation (43), considering specific perturbation Ĝ± ϵ , then it is not possible to calculate explicitly the exact expression to the perturbed unsteady volume flow rate Q± ϵ .However, we can overcome this difficulty by calculating the time evolution of the perturbed flow for fixed magnitude ϵ > 0. Consequently, Figure 12 shows the perturbation given by (49) with magnitude ϵ = 0.1 for shear-thinning and shear-thickening viscoelastic fluid.More properly, Figure 12a illustrates the evolution in time of the perturbation (49) for different power index values, where the mean pressure gradient is constant.In the same way, Figure 12b shows the evolution in time of the perturbation (49) for different power index values, where the mean pressure gradient is non-constant and given by Function (47).Therefore, for specific flow regimes associated with the solutions, we can conclude that, after the initial transition phase, the unsteady volume flow rate behavior is stable in both situations of shear-thinning and shear-thickening viscoelastic fluids.

Conclusions
In this work, we proposed a new three-dimensional model for a non-Newtonian incompressible fluid, where viscosity and elasticity vary depending on the shear rate being this dependence a power-law type.The study carried out allowed us to obtain relevant results with applications in different areas of scientific development; for example, in biology, physics, industry, engineering, chemistry, hemodynamics and mathematics.Below, we highlight the main conclusions, namely:

•
We obtained the ordinary differential equations involving unsteady mean pressure gradient (unsteady wall shear stress, respectively), unsteady volume flow rate, Womersley number, power index, and material parameters; • We obtained numerical simulations for the unsteady volume flow rate in specific flow regimes; • Based on the unsteady solution obtained for the volume flow rate, we presented numerical simulations to the unsteady three-dimensional velocity field in specific flow regimes; • We presented a study on the perturbed flows.
In addition to the innovation in the proposed model in the field of applications, we highlight the innovative use of Cosserat theory associated with fluid dynamics.This allowed us to move from a three-dimensional model to a one-dimensional model, thus obtaining information about the behavior of the flow volume, which was decisive in the numerical simulations presented for three-dimensional velocity.Finally, for future work, we intend to apply this approximation theory to investigate open problems related to curved tubes and fluid-structure interaction.

Figure 7 .
Figure 7. Mean pressure gradient dependent on time given by (47).
is associated with the shown behavior of the unsteady volume flow rate solution Q during and after the initial transient phase.