Fractional Generalizations of Maxwell and Kelvin-Voigt Models for Biopolymer Characterization

The paper proposes a fractional generalization of the Maxwell and Kelvin-Voigt rheological models for a description of dynamic behavior of biopolymer materials. It was found that the rheological models of Maxwell-type do not work in the case of modeling of viscoelastic solids, and the model which significantly better describes the nature of changes in rheological properties of such media is the modified fractional Kelvin-Voigt model with two built-in springpots (MFKVM2). The proposed model was used to describe the experimental data from the oscillatory and creep tests of 3% (w/v) kuzu starch pastes, and to determine the values of their rheological parameters as a function of pasting time. These parameters provide a lot of additional information about structure and viscoelastic properties of the medium in comparison to the classical analysis of dynamic curves G’ and G” and shear creep compliance J(t). It allowed for a comprehensive description of a wide range of properties of kuzu starch pastes, depending on the conditions of pasting process.


Introduction
Biopolymers produced by living organisms can be divided into three main groups [1]: polysaccharides (cellulose, starch, pectin, chitin, glycogen, inulin), polypeptides (proteins) and nucleic acids (DNA and RNA). Biopolymers are characterized by significant sensitivity to physical and chemical factors, therefore, careful examination of their structure and rheological properties requires the use of non-invasive measurement methods. One of them is the oscillation technique. It consists in subjecting the material sample to sinusoidal strain or stress and recording its reaction. This measurement does not affect the structure of the medium-as long as the amplitude is not too large-and can be used to control the processes taking place in time [2].
The obtained experimental data from dynamic tests in the form of storage modulus G' (responsible for elastic properties of the material) and loss modulus G" (representing viscous characteristics of medium) can be described by rheological models. The classical phenomenological rheological models are composed of two types of elements-spring and dashpot. The spring element behaves in accordance with Hooke's law and represents ideal elastic response of a material to applied stress; while the viscous element obeys the law of Newton and describes the energy losses due to viscous dissipation [3][4][5]. Depending on the number and manner of combining of these basic elements, in the literature there have been proposed a number rheological models [3,[6][7][8]. However, their practical use is limited to media with not very complex rheological properties [9][10][11], which are certainly not biopolymers. Much more possibilities in the description of the behavior of such materials provide fractional rheological models, created on the basis of differential calculus of fractional order 0 α 1 [12]. Fractional model in terms of construction differs from a classical model in that the certain standard components as a spring or dashpot are replaced by Scott-Blair elements (springpots). Each springpot can be understood as a component having intermediate properties between a purely elastic element (for which α = 0) and a perfectly viscous element (for which α = 1) [13,14].
Fractional rheological models allow to describe the dynamic behavior of a medium with a single constitutive equation which contains a certain number of parameters that are the constants determining viscoelastic properties of a given material [15]. The identification of these values is so-called reverse problem-in the first place, approximation of the experimental data with trigonometric functions is made, and then rheological parameters of the applied model are determined. Obtained quantities allow for a comprehensive assessment of the medium structure [16][17][18][19].
In the literature, there are two basic phenomenological models that describe the rheological behaviors of viscoelastic materials. These are: Maxwell model and Kelvin-Voigt model [3,20,21]. This division results from the way of connecting the elastic Hooke and viscous Newton elements (serial or parallel). Two-element models take into account only a single relaxation time and a single elastic modulus, which practically excludes the possibility of their use for the characterization of viscoelastic media in a wide range of oscillation frequencies [6,22,23].
The aim of the study was to propose and compare the new modified forms of Maxwell and Kelvin-Voig fractional rheological models to improve the description of the dynamic behavior of biopolymeric materials in the area of viscoelastic plateau and within the scope of the smallest oscillation frequencies ω. Finding the optimal rheological model whose parameters are simultaneously the material constants, allows for a comprehensive assessment of the structure and viscoelastic properties of the medium.

Kuzu starch pastes
The rheological studies involved 3% (w/v) Japanese kuzu starch pastes (Terrasana, Netherlands). Pasting process of aqueous starch suspensions was carried out at 90°C for 15, 30, 45, 60 or 75 min, while stirring with a magnetic stirrer at a constant rotation speed of 300 rpm. After 24 hours, the obtained pastes (Fig 1) were subjected to oscillatory and creep tests by means of rotary rheometer Physica MCR 301 (Anton Paar, Austria) with a cone-plate configuration (cone diameter-60 mm, cone slope-1°, gap width-117 μm). Rheological measurements conducted at a constant temperature of 25°C included the determination of storage modulus G' and loss modulus G" for oscillation frequencies ω in the range of 6.3Á10 −4 to 450 s -1 , and at a given sinusoidal strain with 3% amplitude, within the linear viscoelastic region. In addition, the 2500-second measurements of shear creep compliance J(t) at a specified constant shear stress value of 1 Pa, also within the range of linear viscoelasticity, were performed.

Maxwell-type models
The classical Maxwell model (CMM) is composed of serially connected Hooke and Newton elements (Fig 2a). The total shear stress is equal to shear stresses acting on both elements Eq (1) and the total shear strain is the sum of the deformations of spring and dashpot Eq (2)   Fractional Generalizations of Maxwell and Kelvin-Voigt Models [3,10].
where σ tot is the total shear stress; σ s is the shear stress acting on spring, expressed by Eq (3); σ d is the shear stress acting on dashpot, expressed by Eq (4); γ tot is the total shear strain; γ s is the shear strain of spring; and γ d is the shear strain of dashpot.
sðtÞ ¼ G e gðtÞ ð 3Þ where σ(t) is the shear stress; γ(t) is the shear strain;η is the viscosity; τ 0 is the characteristic relaxation time; G e is the equilibrium modulus; G 0 N is the plateau modulus; and t is time. By differentiating Eq (2) and substituting Eqs (3) and (4), the constitutive equation of the model (5) is obtained. The stress-strain relation contains a single relaxation time τ 0 and two different elastic moduli G e and G 0 N : Carrying out Fourier transform of Eq (5) according to the Eqs (6) and (7) z Á f ðtÞ! and knowing that the ratio of obtained stress and strain transforms defines complex modulus G Ã Eq (8) [24], the equation describing the value of complex modulus G Ã as a function of oscillation frequency ω for the classical Maxwell model (CMM) is obtained Eq (9).
GÃ ¼ŝ ðoÞ gðoÞ ð8Þ where f(t) is the original function;f ðoÞ is the Fourier transform of the function f(t); i is the imaginary unit; m is the order of derivative; and z is the constant. Separating real and imaginary parts of Eq (9), the equations describing storage modulus G' Eq (10) and loss modulus G" Eq (11), respectively, are obtained.
G@ðoÞ  loss angle δ and model curves resulting from the classical Maxwell model (CMM) for kuzu starch pastes which were pasted at 90°C for 30 min. Presented model contains only three rheological parameters (τ 0 , G e , G 0 N ) and is completely inadequate to describe the obtained experimental data. The simplest classical form of Maxwell model provides a completely different trend for modeling curves and experimental points.
In order to improve the description of the experimental data, it was proposed to introduce the fractional Maxwell model with one built-in springpot (FMM1). In this model, the dashpot has been replaced with a Scott-Blair element (Fig 2b). The behavior of the springpot-type element is expressed by Eq (12) [13,14]: where α is the fractional exponent. By differentiating Eq (2) and substituting Eqs (3) and (12), the constitutive equation of the model (13) is obtained. The stress-strain relation contains a single relaxation time τ 0 , two different elastic moduli G e ,G 0 N , and fractional exponent α: Using Eqs (6), (7) and (8), the equation describing complex modulus G Ã as a function of oscillation frequency ω for the fractional Maxwell model with one springpot (FMM1) is obtained: Separating real and imaginary parts of Eq (14)-with the use of Eq (15) [21]-the equations describing storage modulus G' Eq (16) and loss modulus G" Eq (17), respectively, are obtained.
where n is the exponent of imaginary unit. Fig 3b shows the experimental values of storage modulus G', loss modulus G", tangent of loss angle δ and model curves resulting from the fractional Maxwell model with one springpot (FMM1). Presented four-parameter model (τ 0 , G e , G 0 N , α) describes the experimental data with an error smaller by several magnitude orders in comparison to the classical Maxwell model (CMM). However, accurate description of the experimental data is still unsatisfactory.
Further expansion of the model by replacing the elastic Hook component with a Scott-Blair element led to the creation of the fractional Maxwell model with two built-in springpots (FMM2) (Fig 2c). Additional springpot is associated with equilibrium modulus G e , characteristic (the longest) relaxation time τ m and fractional exponent β by the following relation: By differentiating Eq (2) according to the rule Eq (19) [10] and substituting Eqs (12) and (18), the constitutive equation of the model (20) is obtained. The stress-strain relation contains two relaxation times τ 0 , τ m , two different elastic moduli G e ,G 0 N and two fractional exponents α, β.
Using Eqs (6), (7) and (8), the equation describing complex modulus G Ã as a function of oscillation frequency ω for the fractional Maxwell model with two springpots (FMM2) is obtained: Separating real and imaginary parts of Eq (21)-with the use of Eq (15)-the equations describing storage modulus G' Eq (22) and loss modulus G" Eq (23), respectively, are obtained.
G@ðoÞ ¼ where: The fractional Maxwell model with two springpot-type elements (FMM2) contains six parameters (τ 0 , τ m , G e , G 0 N , α, β) which leads to a further increase in the quality of model fitting (Fig 3c). In the case of this model, relatively good description of loss modulus G" in the range of intermediate oscillation frequencies ω was achieved. The improvement in the description of storage modulus G' for large values of the oscillation frequencies ω was also obtained. On the other hand, the trend of experimental data and model predictions concerning the tangent of loss angle δ are fundamentally different.

Kelvin-Voigt-type models
The classical Kelvin-Voigt model (CKVM) is composed of parallel-connected Hooke and Newton elements (Fig 4a). In this case, the total shear stress is the sum of shear stresses acting on both elements Eq (26), and the shear strain is the same for spring and dashpot Eq (27) [3,10].
Substituting Eqs (3) and (4) into Eq (26), the constitutive equation of the model (28) is obtained. The stress-strain relation contains a single relaxation time τ 0 and two different elastic moduli G e ,G 0 N : Using Eqs (6), (7) and (8), the equation describing complex modulus G Ã as a function of oscillation frequency ω for the classical Kelvin-Voigt model (CKVM) is obtained: Separating real and imaginary parts of Eq (29), the equations describing storage modulus G' Eq (30) and loss modulus G" Eq (31), respectively, are obtained.   Therefore, there was an attempt to modify the classical Kelvin-Voigt model (CKVM) by replacing the dashpot with a Scott-Blair element (Fig 4b). Substituting Eqs (3) and (12) into Eq (26), the constitutive equation of the fractional Kelvin-Voigt model with one built-in springpot (FKVM1) Eq (32) is obtained. The stress-strain relation contains a single relaxation time τ 0 , two different elastic moduli G e ,G 0 N and fractional exponent α: Using Eqs (6), (7) and (8), the equation describing complex modulus G Ã as a function of oscillation frequency ω for the fractional Kelvin-Voigt model with one springpot (FKVM1) is obtained: Separating real and imaginary parts of Eq (33)-with the use of Eq (15)-the equations describing storage modulus G' Eq (34) and loss modulus G" Eq (35), respectively, are obtained. G@ðoÞ Replacement of the dashpot with a Scott-Blair element significantly improves a description of the experimental data by means of the fractional Kelvin-Voigt model with one springpot (FKVM1) (Fig 5b). Presented four-parameter fractional model (τ 0 , G e , G 0 N , α) very well captures the course of storage modulus G' in the whole investigated range of oscillation frequencies ω. While in case of loss modulus G" and tangent of loss angle δ a substantial improvement in description of the experimental data in the range of intermediate values of oscillation frequencies ω, was obtained.
Further modification of the Kelvin-Voigt model by replacing the elastic Hooke element with another Scott-Blair component allowed to obtain the fractional Kelvin-Voigt model with two built-in springpots (FKVM2) (Fig 4c). Substituting Eqs (12) and (18) into Eq (26), the constitutive equation of the model Eq (36) is obtained. The stress-strain relation contains two relaxation times τ 0 , τ m , two different elastic moduli G e ,G 0 N and two fractional exponents α, β: Using Eqs (6), (7) and (8), the equation describing complex modulus G Ã as a function of oscillation frequency ω for the fractional Kelvin-Voigt model with two springpots (FKVM2) is obtained: Separating real and imaginary parts of Eq (37)-with the use of Eq (15)-the equations describing storage modulus G' Eq (38) and loss modulus G" Eq (39), respectively, are obtained. Presented six-parameter fractional Kelvin-Voigt model (τ 0 , τ m , G e , G 0 N , α, β) with two springpot-type elements (FKVM2) very well describes the course of storage modulus G' over the entire range of oscillation frequencies ω (Fig 5c). Trends in loss modulus G" and tangent of loss angle δ are also very well described by the model curves, except for a narrow range of the smallest oscillation frequencies ω.
The above applicability analysis of various forms of Maxwell and Kelvin-Voigt models indicates that it is necessary to propose modifications to these models that allow for a description of the experimental data on viscoelastic solids in broad range of oscillation frequencies ω.Modification of fractional rheological models The two-element rheological models are not able to describe the dynamic behavior of real materials in a wide range of oscillation frequencies ω with an acceptable accuracy [9]. In order to improve the description quality of the loss peak on loss modulus curve G" (Figs 3 and 5), the paper proposes a modification of the fractional Maxwell and Kelvin-Voigt models with two built-in springpots (FMM2 and FKVM2) by adding to the Eqs (23) and (39) a new component called the network durability. Rayleigh dimensional analysis indicates that a function describing the dependence of loss modulus G" from the rest of parameters associated with the phenomenon of energy dissipation must fulfill the following condition: where η 0 is the Newtonian steady state shear viscosity. According to the Fourier principle, all the laws of physics are expressed by dimensionally homogeneous equations [25]. This makes it possible to write Eq (41), which after substituting the appropriate units for all variables, takes the form Eq (42).
The dimensional compatibility on both sides of the Eq (42) occurs when: Presented problem has three possible solutions under the assumption that the exponent of oscillation frequency ω takes the value of a = −1 and the remaining exponents are natural numbers: The analysis of experimental data has allowed to determine that the best fit to the loss peak on G" curve is achieved when b = 2 and c = 0.
The network durability D thus obtained characterizes flowing abilities of the imaginary elementary cells of biopolymer network-units locked by a minimal number of nodes and having the individual movement ability: Substituting Eq (45) into Eqs (23) and (39), the equations describing the value of loss modulus G" as a function of oscillation frequency ω for the modified fractional Maxwell model with two springpots (MFMM2) Eq (46) and the modified fractional Kelvin-Voigt model with two springpots (MFKVM2) Eq (47), respectively, are obtained. The proposed modified fractional models increase the ability to describe both fast and slow processes of energy dissipation. Fig 6 shows the experimental values of storage modulus G', loss modulus G", tangent of loss angle δ and model curves resulting from the seven-parameter (τ 0 , τ m , G e , G 0 N , α, β, η 0 ) modified fractional models with two built-in springpots: Maxwell-type (MFMM2) (Fig 6a) and Kelvin-Voigt-type (MFKVM2) (Fig 6b). In the case of MFKVM2, a modification in the form of additional component called the network durability D allows for a significant improvement in a description of the experimental data over the entire range of oscillation frequencies ω.

Parameters of fractional rheological models and their limitations
Presented fractional rheological models of Maxwell and Kelvin-Voigt-Eqs (10) and (11); (16) and (17); (22) and (23); (30) and (31) (Table 1), which describe a number of viscoelastic properties of a given material: • plateau modulus G 0 N -represents the power of biopolymer network and its resistance to aging during time; value of the parameter was determined using the minimum method [26]-it is based on the assumption that the plateau modulus G 0 N corresponds to the value of storage modulus G' at the oscillation frequency ω in which the tangent of loss angle δ reaches a minimum: • equilibrium modulus G e -illustrates the total elasticity of biopolymer network; value of the parameter was taken as equal to the inverse of the intercept in the equation of line tangent to the shear creep compliance curve J(t) at endpoint [16]; • characteristic relaxation times τ 0 and τ m -represent the shortest and the longest time required to complete stress relaxation in the biopolymer network; values of the parameters correspond to the inverses of oscillation frequencies ω at which the dynamic curves G' and G" intersect themselves; • fractional exponents α and β-indicate which properties dominate in the material: these parameters take values from 0 (for perfectly elastic solid) to 1 (for ideal Newtonian fluid); • Newtonian steady state shear viscosity η 0 -characterizes flowing abilities of the elementary units of biopolymer network; value of the parameter was taken as equal to the inverse of the slope in the equation of line tangent to the shear creep compliance curve J(t) at endpoint [16].
Thanks to the correlations available in the literature [16,18,22,[26][27][28][29][30][31], it is possible to determine further indicators characterizing the rheological properties of biopolymer, such as: • plateau compliance J 0 Nrepresenting the force with which the entanglements of the biopolymer network suppress any kind of long-range configurational rearrangements: • steady state compliance J e -being the measure of stored energy in the steady-state flow under the influence of low stresses: • dispersion modulus f-characterizing the molecular weight distributions of the biopolymer: • coefficient of the network vibration damping k-representing the resistance of biopolymer network to oscillatory deformations: • width of the viscoelastic plateau L-specifying the polydispersity index of the biopolymer: • cross-linking density ω 0 -characterizing the structure of biopolymer network: • gel stiffness S-representing the degree of fragility, brittleness of gel:  • average entanglement molecular weight M e -being the average molecular weight between topological constraints: • average molecular weight between cross-links M c -being the average molecular weight of biopolymer chains between two consecutive junctions: where ρ is density of biopolymer; R is universal gas constant; and T is temperature. According to the second law of thermodynamics, the dynamically deformable real solids must be characterized by positive values of energy dissipation and internal work [32]. It means that the proposed rheological models have physical meaning only when model values of storage modulus G' and loss modulus G" are positive in the entire analyzed range of oscillation frequencies ω [33]. This condition is fulfilled when: Statistical evaluation of rheological models In order to determine the quality of the experimental data description by means of the proposed Maxwell and Kelvin-Voigt rheological models, the statistical evaluation referred to tangent of loss angle δ Eq (59) was performed.
The analysis was carried out using the statistical indicators, such as [34,35]: • mean percentage error (MPE): • mean bias error (MBE): ðtan d mod; j À tan d exp; j Þ ð 61Þ • root mean square error (RMSE): • modelling efficiency (EF): • chi-square test (χ 2 ): where tan δ exp,j is the experimental value of tangent of loss angle δ; tan δ mod,j is the model value of tangent of loss angle δ; tan δ exp,ave is the average experimental value of tangent of loss angle δ; j is the index of experimental point; N is the number of experimental points; n is the number of parameters in the model.

Results and Discussion
In order to analyze the accuracy of the experimental data description by the proposed modified Maxwell and Kelvin-Voigt models, the new rheological studies for kuzu starch pastes were carried out. The obtained results (Fig 6b) revealed that in the oscillation frequency range from 6.3Á10 −4 to 450 s -1 the biopolymer was located in a hyperelastic physical state called the viscoelastic plateau region. This area is associated with changes in the position of chain segments of biopolymer (rotary and sliding movement) in an absence of movement of all macroparticles. It may result in significant deformations of the material even in the case of small external stresses [16,22]. In the entire analyzed range of oscillation frequencies ω, the storage modulus G' (representing elastic properties of biopolymer) was greater than the loss modulus G" (characterizing viscous features of medium). Table 2 lists the goodness-of-fit indicators for presented in the work rheological models. The analysis showed that the optimal model for describing the dynamic behavior of biopolymers, such as kuzu starch pastes, is modified fractional Kelvin-Voigt model with two built-in springpots (MFKVM2). In this case, for pasting time t = 15 min, the mean percentage error was equal to MPE = 2.12%, and the modeling efficiency has reached the value of EF = 0.995. It was also confirmed that the rheological models of Maxwell-type do not work in the case of modeling of viscoelastic solids, that is the media in which elastic properties dominate over viscous properties [36].
The rheological parameters of Japanese kuzu starch pastes determined on the basis of the modified fractional Kelvin-Voigt model with two springpots-Eqs (38) and (47)-provide a lot of additional information about the structure and viscoelastic properties of the biopolymer in comparison to the classical analysis of dynamic curves G' and G" and shear creep compliance J (t). This is particularly important in the field of materials science to design the utility and functional characteristics of products.

Conclusions
The analysis of various types of Maxwell and Kelvin-Voigt rheological models for the description of viscoelastic properties of biopolymers indicated that: • Maxwell-type rheological models do not work in the case of modeling of viscoelastic solids.
The optimal model for describing the dynamic behavior of kuzu starch pastes is modified fractional Kelvin-Voigt model with two built-in springpots (MFKVM2). This model can be  Fractional Generalizations of Maxwell and Kelvin-Voigt Models an important tool for specialists in the field of materials engineering, to design the structure and rheological properties of the media, which have a direct impact on utility and functional characteristics of products.
• The proposed modification of the fractional Kelvin-Voigt model with two springpots (MFKVM2), consists in introducing an additional component called the network durability D into the equation for loss modulus G". It allowed for a significant improvement in the quality of experimental data description in the range of the lowest oscillation frequencies ω, which correspond to slow dissipative processes ( Table 2). The proposed modulus D Eq (45) characterizes the impact of cross-linking on the flowing abilities of imaginary cells of the biopolymer network.