Review of Non-Newtonian Mathematical Models for Rheological Characteristics of Viscoelastic Composites ()
1. Introduction
The rheological properties and extrudate behavior of polymer melts are of central importance in the processing and fabrication of polymer products. The melt-rheological behavior of short sisal, coir, and pineapple fibre-re- inforced polymer systems has been studied in various works [1] [2] . However, there is limited research on the effects of parameters such as the fibre length, the fibre-matrix interaction, and the aspect ratio of green-fibre-re- inforced thermoplastic biopolymer composites on nonlinear rheological behavior [3] [4] .
The mathematical explanation of a viscoelastic fluid is much more complex than its Newtonian counterparts. In addition to the conservation equations of mass and momentum, the constitutive equation or rheological equation of state is required; this relates stress to deformation. For a viscoelastic liquid this relationship is nonlinear and it has no standard form universally valid for each fluid in every flow situation. This reality is one of the reasons why the subject of viscoelasticity is so challenging.
The constitutive equation should not only describe the rheological characteristic of the polymer melt but also give the final fibre orientation of the composite. For this reason it is fundamental to evaluate the role of the biocomposite’s rheology and the natural fibre-polymer interaction. It has been observed that the total stress of the composite increases as fibres are added; consequently a satisfactory constitutive equation could be achieved by adding an extra stress term to an already existing constitutive equation, which then adequately describes the polymer melt [5] -[7] .
Accordingly, constitutive equations found in the literature that adequately describe polymer melt will be explored for their application in biocomposites processing. Particular focus will be given to nonlinear rheological characteristics of viscoelastic materials. The power-law, Cross WFL, Casson, Bird-Carreau and Hershel Bulkley models are among the most preferred rheological models due to their ability to predict velocity and pressure distributions in uniform flows in addition to their simple representation of shear thinning behaviour [7] - [12] . However, in the case of high shear stress of viscoelastic polymer melt, the predictive power of these models is considerably reduced [10] .
In this study, a review of nonlinear rheological models for viscoelastic materials from natural-fibre-reinforced thermoplastic polymers will be presented by a special review of the Upper convected Maxwell, Phan-Thien- Tanner (PTT), K-BKZ, Oldroyd-B, Giesekus and Whhite-Metzner constitutive models.
2. Viscoelastic Characteristics of Materials
Viscoelasticity is the property of a material to demonstrate both viscous and elastic properties under the same conditions when it undergoes deformation. Viscous materials present resistance to shear flow and strain linearly with time when a stress is applied. The shear stress of these materials depends on strain: when strain is applied and then released, they return to their initial configuration. Some common and well-known viscoelastic materials include paint, blood, ketchup, honey, mayonnaise, polymer melt, polymer solution and suspension, shampoo, and corn starch.
At constant temperature, water, air, ethanol, and benzene are represented as Newtonian fluids. This means that the rapport between the shear stress versus shear rate is a straight line with a constant slope for a fixed temperature that is independent of the shear rate. Also, the plot passes through the origin: that is, the shear rate is zero when the shear stress is zero [13] .
A fluid that does not behave in a Newtonian fashion between shear stress and shear rate when it undergoes deformation is commonly termed non-Newtonian. This means that the relation between shear stress and shear is not a straight line but is non-linear. High-molecular-weight liquids, which include polymer melts and solutions of polymers, as well as liquids in which fine particles are suspended, are usually non-Newtonian. In this case, the slope of the shear stress versus shear rate plot will not be constant as we change the shear rate. When viscosity decreases with increasing shear rate, the fluid is called shear-thinning. In the opposite case, where the viscosity increases as the fluid is subjected to a higher shear rate, the fluid is called shear-thickening. Shear-thinning fluids also are called pseudoplastic fluids and shear-thickening fluids are also called dilatants. Shear-thinning behavior is more common than shear-thickening [13] - [16] .
Another type of non-Newtonian fluid is viscoplastic or “yield stress” fluid. This is a fluid that will not flow when only a small shear stress is applied. The shear stress must exceed a critical value known as the yield stress for the fluid to start flowing. Hence, viscoplastic fluids behave like solids when the applied shear stress is less than the yield stress. Once it exceeds the yield stress, the viscoplastic fluid will flow just like an ordinary fluid [14] [15] .
On the other hand, some classes of fluids exhibit time-dependent behavior. This means that even ata given constant shear rate, the viscosity may vary with time. This category of material comprises both thixotropic and rheopectic fluids, whose behaviors in this respect are opposites. The viscosity of a thixotropic liquid will decrease with time under a constant applied shear stress. However, when the stress is removed, the viscosity will gradually recover with time as well. By contrast, rheopectic fluid behavior can be observed when the fluid increases in viscosity with time when a constant shear stress is applied [15] [16] .
3. Rheological Modelling of Viscoelastic Composites
Mathematical models of viscoelasticity are mostly based on a differential or integral representation. From a mathematical point of view, the differential representation is easier to handle than the integral one. However, the integral representation is capable of predicting the time dependence more generally.
The characteristic feature of linear viscoelastic materials is that the stress is linearly proportional to the strain history. Linear viscoelasticity is usually applicable only for small deformations, low rate, low stress, and linear materials. However, in reality about 90% of fluids are nonlinear, with large deformations, and with nonlinear response in the presence of such deformations. Nonlinear viscoelastic behavior is usually exhibited when the deformation is large and most of the time when the material changes its properties under deformations [17] - [24] . Consequently, nonlinear viscoelastic mathematical models are needed [25] -[30] .
Existing nonlinear mathematical rheological models are often constructed through modifications and extensions to higher-order stress or strain terms of the linear theory [31] . As noted earlier, from a mathematical point of view, the integral representation of a viscoelastic constitutive equation is more difficult to perform than the differential form. Thus, several models characterized by elastic, viscous, and inertial nonlinear contributions with various complexities have been developed for describing the nonlinear behavior of these materials [32] . However, in these models, due to the mathematical complications, only the elastic or viscous nonlinearity is often taken into account and the inertial contribution is ignored. Moreover, there are only a few theoretical models formulated with constant-value rheological material parameters.
For these reasons, nonlinear models with constant rheological coefficients are required. The elastic and viscous nonlinearities are taken simultaneously into consideration through a simple nonlinear generalized Maxwell fluid model consisting of a nonlinear spring connected in series with a nonlinear dashpot obeying a power law with constant material coefficients [33] -[36] . According to a previous study by Bauer (1984), suitable constitutive equations for viscoelastic materials must relate stress, strain, and their higher time derivatives: or better said, they must take into consideration the elastic, viscous, and inertial nonlinearities simultaneously [35] . Moreover, various researchers have explored how polymer viscoelasticity affects the diameter distribution of polymer melt- extrudate fibres and have demonstrated that increasing the melt viscosity leads to an increase in fibre diameter but has little effect on the diameter distribution [36] . The commonly used Phan-Thien-Tanner (PTT) and Upper- Convected Maxwell (UCM) constitutive models assume constant shear stress acting on the fibre surface and neglect the effects of heat transfer [37] - [46] . The K-BKZ type of constitutive equation has been widely used in various studies on predicting the rheological behavior of viscoelastic materials. For example, Galante used the constitutive equation to describe viscoelastic effects in an integral equation of the K-BKZ type, suitable for polymer solutions and melts [40] . The problem with the constitutive equation of K-BKZ is that it is not fully applicable to predicting the nonlinear rheological behavior of viscoelastic materials [39] - [41] .
The behavior and properties of a non-Newtonian fluid within compressible flow are provided by the conservation of mass and momentum equations, respectively.
4. Governing Equations
The governing equations for the annulus flows are presented as follows.
Continuity equation for incompressible fluids: 
where
is the velocity vector.
Equation of momentum: 
where
is the density,
is the pressure, and
is the stress tensor.
5. Constitutive Equations
There are many numerical representations for viscoelastic models. The most common models are Upper-con- vected Maxwell (UCM), Phan-Thien-Tanner (PTT), Oldroyd-B, Giesekus, K-BKZ and White-Metzner.
5.1. K-BKZ Model
The K-BKZ integral constitutive equation with multiple relaxation times describes and predicts the stresses within a fibrous suspension, solution, or molten polymer. Also, an extra term is added to the constitutive equation to account for the extra stresses due to the presence of fibres and to predict the orientation of a given fibre undergoing stresses within the suspension or molten [40] . The motivation for developing such a constitutive equation with these two considerations is to present an equation that can describe the rheological behavior of polymeric fibrous solutions and moltens while also to have a model, which can be used in numerical simulations.

The fibre equation is:

where
is the viscosity of the polymer,
is the shear rate,
is the fibre volume fraction,
and
are the diameter and length of the fibre,
is the number of the suspension, and
is the average distance from a given fibre to its nearest neighbor [41] . Dinh and Armstrong proposed the following expression to calculate the distance between two fibres:
for aligned systems
for random systems
Orientation tensor:
![]()
The polymer equation is:
![]()
where
is the shear stress for the polymer,
and
are the relaxation times and relaxation moduli,
is the number of relaxation modes,
is the finger strain tensor, and
is its first and second invariants.
is the strain memory function, and the following formula is used, proposed by Papanastasiou et al.:
![]()
For simple shear flow, the strain-memory function is given as:
![]()
![]()
where
and
are nonlinear model constants to be determined from shear and elongation flow data, respectively.
is the shear strain the Tstrain-memory function in simple shear flow is dependent on
butnoton
. This is expected since
is viewed as a shear-thinning parameter, while
is viewed as an elongational-thinning parameter [40] [41] .
5.2. Upper-Convected Maxwell Model (UCM)
UCM model is a differential generalization of the Maxwell model for the case of large deformations based on the upper-convected time derivative. The model can be written as: ![]()
![]()
is the tensor of the deformation rate
![]()
where
is the relaxation time,
is zero shear viscosity, and
is tis he upper-convected time derivative of the stress tensor, which is expressed as:
![]()
The UCM constitutive model incorporates memory effects of materials, but its viscosity is constant at various shear rates [42] -[45] .
5.3. White-Metzner Model
The White-Metzner model is derived from the network theory of polymers (White and Metzner 1963). Modification of the viscosity and relaxation parameters as a function of the shear rate,
leads to the White-Metzner model. This model exhibits shear thinning, not because of nonaffine motion, but because the relaxation is accelerated at high strain rates, where the relaxation is faster than any deformation [46] . The viscoelastic differential constitutive model takes the form:
![]()
where
can be obtained from the experimental shear viscosity curve and the function
can be obtained from the experimental first normal stress difference curve.
5.4. Phan-Thien-Tanner Model (PTT)
The PTT model refers to a quasi-linear viscoelastic constitutive equation, which is widely used in simulation of polymer solution flows. The original Phan-Thien Tanner equation was written using both of the following modifications simultaneously: the Gordon Schowalter derivative and the segment kinetics. It employs specific forms for the creation and destruction rates of the network junctions in the network theory of Lodge and Yamamoto [47] [48] . Although the Phan-Thien Tanner model over predicts both the shear viscosity at higher shear rates and the transient and extensional properties, it accurately predicts the zero shear viscosity and seems suitable for numerical simulations of polymer melts. It is worth noting that, compared to integral models such as the Bird- Carreau and Wagner models, differential models such as the PTT model provide robust numerical algorithms and exhibit good behavior in FEM simulations [47] . The interested reader is referred or more accurate comparison between the K-BKZ (as integral model) and the PTT model.
In the PTT model, the extra stress tensor is considered as the sum of the viscoelastic component
, and the purely Newtonian component
.
![]()
in which
is given by:
![]()
where
is the strain rate tensor.
The complete form of the PTT constitutive equation for the viscoelastic component
is:
![]()
or
![]()
where the Oldroyd’s upper convective derivative
isdefined by:
![]()
where
represents the velocity matrix,
is the transpose of the velocity matrix, and
is the material derivative of the velocity matrix.
Analyzing the expression above, the first term represents the stress tensor transport and the transient part of the flow. In the second term, the slipping between the fluid polymeric chains is computed. The third term includes the elastic effects. Finally, the term on the right side of the equality represents the diffusive effects.
is the stress tensor,
the deformation-rate tensor,
the fluid relaxation time, and
the relaxation module. The parameter
controls the amount of movement between the fluid polymeric chains [49] . For
the model is named PTT Affine and the slipping between the polymeric chains is neglected. The function
depends on the rate of creation and destruction of the links between the chains.Moreover, the PTT model incorporates the memory effect of materials and its viscosity can vary with the change of the shear rate. When
is zero, PTT constitutive equation reduces to its simplified form (SPTT):
![]()
Phan-Thien and Tanner assumed specific forms for the creation and destruction rates of the network junctions and derived a constitutive equation containing two free parameters,
and
[48] . The exponential constitutive model takes the following form:
![]()
and
are the adjustable parameters of the model.
The parameters
and
are the viscosity and relaxation time respectively, measured from the equilibrium relaxation spectrum of the fluid. They are not considered as adjustable parameters of the model. The PTT model can be solved using a single relaxation time or multiple relaxation times, similar to the Giesekus model. The linear form of the PTT model predicts shear thickening at high elongational rates, after which a plateau is reached [47] -[49] .
5.5. Giesekus-Leonov Model
Giesekus proposed a constitutive model based on a concept of configuration-dependent molecular mobility. In this model, the viscoelastic component of the extra stress tensor is represented with the following parameters
,
and
. Owing to the highly nonlinear nature of the model equations, all of the properties need to be obtained numerically [50] -[55] .
![]()
The
parameter is the dimensionless Giesekus-model mobility factor and controls the extensional viscosity and the ratio of second normal stress difference to the first one. The dimensionless The Giesekus-model mobility factor used to evaluate the anisotropic drag is represented by
. For
the model becomes the isotropic UCM model, while for
the model is merely an anisotropic drag, and for
the model re- presents shear-thinning behavior. The Giesekus model predicts the tension-thickening region for elongational flow, after which a plateau is reached; but it shows the existence of a tension-thinning region at high strain rates [56] -[59] .
5.6. Oldroyd-B Model
The Oldroyd-B model is a principal form of viscoelastic model:
![]()
where
is the rate of deformation tensor,
is the shear viscosity,
is the relaxation time, and
is the second relaxation parameter [60] -[62] . The Oldroyd-B model is mostly used to describe the rheological characteristics of polymer liquids composed at low concentration and moderate shear rates for high-viscosity Newtonian molecular weight polymer.
1) At
the model simplifies to a second-order fluid with a vanishing second normal stress coefficient.
2) At
the model reduces to the convected Maxwell model.
3) At
the model reduces to a Newtonian fluid with viscosity
.
6. Conclusions
Understanding the complex behaviour of polymer materials and interpreting it as one general equation requires a vast knowledge of the characteristics and formation of this complex type of material. Since rheometer does not provide the necessary information for all important rheological properties, constitutive equations are the best available tools for effective process control.
The development of valuable models for composite behavior and the exploration of appropriate constitutive equations to describe this complex behavior have been a high priority for many researchers. However, understanding the rheological behavior of viscoelastic composites is a redoubtable challenge.
The K-BKZ, PTT, Oldroyd-B, and Giesekus models will be widely studied and extended in future works. These models will be applied to the prediction and determination of the shear viscosity of viscoelastic composites as a function of shear stress and shear rate during extrusion and injection moulding processes.
Acknowledgements
This work was carried out with the support of Mitacs funding; the University of Toronto and the Ford Motor Company financially supported this study.
NOTES
*Corresponding author.