## Abstract

In this chapter, we present a brief description of existing viscoelastic models, starting with the classical differential and integral models, and then focusing our attention on new models that take advantage of the enhanced properties of the Mittag-Leffler function (a generalization of the exponential function). The generalized models considered in this work are the fractional Kaye-Bernstein, Kearsley, Zapas (K-BKZ) integral model and the differential generalized exponential Phan-Thien and Tanner (PTT) model recently proposed by our research group. The integral model makes use of the relaxation function obtained from a step-strain applied to the fractional Maxwell model, and the differential model generalizes the familiar exponential Phan-Thien and Tanner constitutive equation by substituting the exponential function of the trace of the stress tensor by the Mittag-Leffler function. Since the differential model is based on local operators, it reduces the computational time needed to predict the flow behavior, and, it also allows a simpler description of complex fluids. Therefore, we explore the rheometric properties of this model and its ability (or limitations) in describing complex flows.

### Keywords

- Mittag-Leffler
- viscoelastic
- memory function
- fractional calculus
- rheology

## 1. Introduction

Since viscoelastic materials are abundant in nature and present in our daily lives (examples are paints, blood, polymers, biomaterials, etc.), it is important to study and understand viscoelastic behavior. Therefore, in this chapter, we further develop the modeling of viscoelasticity making use of fractional calculus tools.

We start this section with some basic concepts that are needed to derive and understand classical and fractional viscoelastic models. These are trivial concepts such as force, stress, viscosity, Hooke’s law of elasticity and also Newton’s law of viscosity. Later, we evolve to more complex concepts of viscoelasticity that involve the knowledge of fractional calculus, integral and differential models.

It is well known that a *force* is any interaction that when unopposed will change the motion of an object/body. *Stress* is an internal resistance provided by the body itself whenever it is under deformation. Stress is defined as the intensity of internal forces developed in the material. The intensity of any quantity is defined as the ratio of the quantity to the area on which it is acting, leading to: Average Stress = Force/Area. If we want to know the stress in one material point, then we must take the limit of the area to zero. A good example on how stress works is given by imagining a person lying on top of thin layer of ice. When the person is lying down on the ice, the force (weight) divided by the area of the surface of the person in contact with the ice is smaller, when compared to the case when someone is standing up (the weight is the same, but the area in contact with the ice is smaller). Therefore, eventually, the ice will break due to the high internal stresses when the person is standing. Finally, we refer to *elasticity* as the ability of a body to resist a distorting influence and to return to its original size and shape when that influence or force is removed. See for example Figure 1 where three springs are stretched. If we remove the weights attached to the springs, the spring would ideally return to its initial/natural position.

Figure 1(b) also shows an experiment where we observe that the force (mass times gravity) applied to the spring (increasing weight) is proportional to the displacement. This is known as *Hooke’s law* (the force

We will now explore the concept of *viscosity* in fluids. The viscosity of a fluid is a measure of the internal resistance to the rate of deformation:

As an example, imagine that we have a thin film of fluid in between two parallel plates, as shown in Figure 2. The fluid is at rest, and suddenly the upper plate starts moving with constant velocity *U*. This velocity will be felt at the bottom layer due to diffusion of momentum, and to keep the bottom wall fixed, we must exert a restraining force, that is measured with a force gage or dynamometer attached to that wall. Note that if we take the view of this portion of fluid as infinitesimally thin layers, we observe that each layer will drag the underlying layer due to the action of viscosity (internal resistance). The higher the viscosity, the more force will be required to deform the fluid at a given speed *U*.

Since the velocity of the thin layer adjacent to the top wall is *U* and the velocity of the bottom layer is 0, the velocity of each layer (for a Newtonian fluid) is given by *u(y)* = *Uy*/*h*, with *y* the coordinate shown in Figure 2(a). Figure 2(b) shows the experimental forces measured for different ratios of *U*/*h*. We observe that the force is proportional to *U*/*h* and *U*/*h* = d*u*(*y*)/d*y*; therefore, we conclude the following (*Newton’s law of viscosity*):

with *u*/d*y* is known as the rate of shear deformation, usually denoted by

A good example of something we may see every day and something that verifies Newton’s law of viscosity is a dashpot. It is used for example as a door closer to prevent it from slamming shut.

### 1.1. Viscoelastic models

The simplest model that considers both viscous and elastic behavior is the linear Maxwell model [1] and can be obtained from a combination in series of a dashpot, *f* and *e* standing for Newtonian fluid and Hookean elastic solid, respectively), as shown in Figure 3.

The total deformation

The three-dimensional version of this model can be easily obtained by considering appropriate tensors instead of the scalar properties of stress and deformation, leading to the following model:

with

where

The Maxwell model is not *observer independent (frame invariant)* and, therefore, the results obtained with this model may not be correct if large deformations are considered (e.g., we may obtain a viscosity that depends directly on the velocity rather than the velocity gradient, which is not correct, and is unphysical). To solve this problem, new models were proposed in the literature that can deal with this non-invariance problem.

Two well-known examples of frame invariant models are the upper-convected Maxwell (UCM) model given by

where

Other well-known example of a frame-invariant but now nonlinear viscoelastic model is the variation of the K-BKZ [2] model proposed by Wagner, Raible and Meissner [3, 4],

where *damping function* [5] (note that it is again assumed that the fluid is at rest for

where

Different differential models were proposed in the literature along the years, with the aim of improving the modeling of complex viscoelastic materials, and with the aim of achieving the same modeling quality of integral models (by only using differential operators). Note that integral models are non-local (in time) operators that take into account all the past deformation of the fluid while differential models ones describe the material response in terms of the rate of change of stress to the *local deformation*, thus influencing the fitting quality of the model and the computational effort to numerically solve them (when performing numerical simulations).

More recently, new models have been proposed in the literature that basically take advantage of the generalization of the exponential function appearing in Eqs. (4), (5), and (7), thus allowing a more broad and accurate description of the relaxation of complex fluids (while the commonly used continuum approach describes the fluid as a whole, with only one relaxation, unless a Prony series is considered, that is, considering a series of the form

## 2. Fractional derivatives

To understand the need and the concept of a fractional derivative and its importance in the context of modeling physical processes, let us start with a simple example (Figure 5).

Imagine a portion of material that is principally formed of two different regions. In these regions, two similar physical processes

Although we have not defined yet what a fractional derivative is, the fact of having the possibility of non-integer derivatives seems quite attractive, allowing the creation of a continuous path between integer-order derivatives that may lead to a better description of the different rates of a certain physical process occurring in the same material. This means that fractional derivatives can *transport* more and more precise local information from the microscopic world to the continuum description.

### 2.1. Riemann-Liouville and Caputo fractional derivatives

Now, to understand a fractional derivative, we start by acknowledging that the *n*-fold integral of a generic function

A generalization to non-integer values of

where we have used

This last fractional derivative is the one chosen to deal with physical processes due to the ease in handling initial and boundary conditions [7].

Next, we present two models that rely on the Mittag-Leffler function (a function closely related to fractional calculus) to improve their modeling and fitting capabilities when describing the behavior of viscoelastic materials. These are the fractional K-BKZ (integral) and the generalized Phan-Thien and Tanner (differential) models.

## 3. Viscoelastic models based on the Mittag-Leffler function

### 3.1. The fractional K-BKZ model

We first note that the Maxwell-Debye relaxation of stress (exponential decay—see Eqs. (4) and (5)) is quite common, but there are many real materials showing different types of fading memory, such as a power law decay

By recognizing that the Caputo fractional derivative of a general function

we obtain a generalized viscoelastic model [10, 11], that can be written in the simple compact form:

This model provides a generalized viscoelastic response, in the sense that when *intermediate* to the spring and dashpot shown in Figure 3 and is thus known as a *spring-pot* [11, 12]. Note that care must be taken when

We can define the fractional Maxwell model (FMM) as a combination of two linear fractional elements (spring-pots) in series. In a series configuration, the stress felt by each spring-pot is the same, that is,

This model allows a much better fit of rheological data, as shown in [12] but it is not frame invariant. However, following the same procedure employed with the Maxwell and K-BKZ model, that is, using the derivative of the relaxation function obtained for the Maxwell model as the memory function of the K-BKZ model, one can also use the derivative of the relaxation function of the FMM and insert it in the K-BKZ model, thus, obtaining a frame-invariant constitutive model, that retains all the good fitting properties of the FMM.

The relaxation function of the FMM can be obtained by solving the fractional differential Eq. (14) considering a constant deformation

where

and a characteristic measure of the relaxation spectrum described by the two spring-pots in series is

This leads to the fractional K-BKZ model proposed by Jaishankar and Mckinley [12, 13], with

Note that here the relaxation modulus

The fractional K-BKZ model is therefore given by:

and we need to ensure that the integral converges (see the Foundations of Linear Viscoelasticity by Coleman and Noll [15]). The main problem seems to be the term

It can be easily shown [1] that a Taylor series expansion of

with

In Refs. [11, 12, 17], the beneficial fitting qualities of this constitutive model framework are discussed in detail. Here, we are interested in determining to what extent the properties of the Mittag-Leffler function can be used to improve the fitting quality of differential models, and this will be discussed in the next subsection.

### 3.2. Generalized Phan-Thien and Tanner model

The previous integral model given by Eq. (18) allows a good fit to experimental rheological data, in flows with defined kinematics where

The original exponential PTT model [19, 20] is given by.

with

The model is then given by

where the factor Γ(β) is used to ensure that

This new model can further improve the accuracy of the description of real data obtained with the original exponential function of the trace of the stress tensor, as shown in [18].

## 4. Parametric study of the GPTT model

We will now present a detailed parametric study on the influence of the new parameters

### 4.1. Steady-state shear flows

As shown in [18], the steady shear viscosity is given by

and

Here

Since we consider a simple plane shear flow aligned with the *x*-axis, we have that

Eqs. (22) and (23) can readily be solved using the Newton-Raphson method (solving first Eq. (23) and then substituting the numerical values obtained for

Figure 6 shows the dimensionless steady shear viscosity obtained for the different parameters of the Mittag-Leffler function,

Figure 7(a) shows the dimensionless steady shear viscosity, now obtained for different values of

Note that (see Figure 7(b)) small variations of the parameter

Figure 7 shows that by setting different combinations of

### 4.2. Steady-state elongational flows

The steady unidirectional extensional viscosity is defined as

with

Figure 8 shows the dimensionless steady elongational viscosity obtained for different parameters of the Mittag-Leffler function. In Figure 8(a), we show the influence of

Note that when we increase

We may conclude that by varying

Figure 9 shows the effect of the parameters used in Figure 7, for the case of elongational viscosity. The results are qualitatively similar to the ones obtained in Figure 7, that is, in terms of changes to the asymptotic slopes at high deformation rates and premature/delayed thinning. It can be observed that the elongational viscosity is more sensitive to changes in the parameters

### 4.3. Steady-state shear and elongational flows

Until now, we have explored generally the influence of the different model parameters on the behavior of the GPTT model for steady flows, but, a more quantitative side-by-side comparison between the shear and elongational flow curves was not performed, and the limited flexibility of the classical exponential PTT model for fitting experimental data (when compared to the GPTT) was not explored. In Figure 10, we try to illustrate the advantages of using the Mittag-Leffler function instead of the classical exponential one. To this end, we present the viscometric predictions obtained for both shear and elongational flows for both models (GPTT and exponential PTT).

Figure 10 illustrates the additional flexibility of using the Mittag-Leffler function, by showing that we can manipulate the magnitude of the increase in the elongational viscosity and at the same time only slightly change the shear viscosity. This allows better fits to rheological data when using the Mittag-Leffler function [18]. Note that in the exponential PTT model, when we increase the

## 5. Conclusions

In this chapter, we have presented a brief introduction to the world of viscoelastic models capable of describing the rheology of complex fluids, and we have summarized some of the well-known classical differential and integral models.

With incorporation of ideas from fractional calculus, most of these models can be further improved, either by changing classical local operators for improved (non-local) fractional versions, or, either using new analytic functions that arise in the realm of fractional differential equations, such as the Mittag-Leffler function.

As an example, we present the fractional K-BKZ model and the recently proposed generalized PTT model. The fractional K-BKZ model allows a better description of fluid flow behavior (when compared to the generalized PTT model), but, increases the need for high computational power. Therefore, the novelty of the present work is our detailed study on the influence of the Mittag-Leffler function in shear and elongational flows of a generalized PTT model.

## Acknowledgments

L.L. Ferrás and A.M. Afonso acknowledge the Project PTDC/EMS-ENE/3362/2014-POCI-010145-FEDER-016665—funded by FEDER funds through COMPETE2020—Programa Operacional Competitividade e Internacionalização (POCI) and by national funds through FCT—Fundação para a Ciência e a Tecnologia; I.P. L.L. Ferrás would also like to thank the funding by FCT through the scholarship SFRH/BPD/100353/2014. M.L. Morgado would like to thank the funding by FCT through Project UID/MULTI/04621/2013 and M. Rebelo would also like to thank the funding by FCT through Project UID/MAT/00297/2013 (Centro de Matemática e Aplicações).