Relaxation Functions and Relaxation Spectra 

 

                                                Univ.-Prof. Dr.-Ing. habil. Josef  BETTEN  

                                                           RWTH Aachen University  

                          Mathematical Models in Materials Science and Continuum Mechanics  

                                                                  Augustinerbach 4-20  

                                                           D-52080  A a c h e n  ,  Germany  

 

                                                           <betten@mmw.rwth-aachen.de>  

 

 

Abstract  

 

This worksheet is concerned with the relaxation behaviour of  viscoelastic materials. Many                                       materials exibit both features of  elastic solids and characteristics of  viscous fluids. Such                                         materials are called viscoelastic. One of  the main features of elastic behaviour is the capacity                                            for materials to store mechanical energy when deformed by loading, and to set free this energy                                completely after removing the load. In viscous flow, however, mechanical energy is continuously                                       and totally dissipated.  

The phenomenological behaviour of  viscoelastic materials can be illustrated by spring-dashpot  

models consisting of several elastic springs and viscous dashpots in parallel or in series.  

  In the following some examples should be discussed in more detail based upon the fundamental  

rheological model of  MAXWELL. Furthermore, nonlinearities are taken into account. 

 

Keywords:    Relaxation Functions and Relaxation Spectra; POYNTING - THOMSON Model;                                                  Distributions by POISSON, MAXWELL, Chi-square, WEIBULL; Nonlinearities  

 

 

Relaxation Function of the POYNTING - THOMSON Model  

 

This model consists of  one linear spring (HOOKE) and one MAXWELL  model in parallel, where                                      the  MAXWELL element is characterized by one linear spring (HOOKE) and one linear dashpot                                              (NEWTON) connected in series. Thus, the POYNTING - THOMSON model has three parameters,                                      two elastic constants and one viscous parameter and is characterized by the following relaxation                                      function:  

>

restart:

 

>	E(t):=sigma(t)/epsilon[0]=E[infinity]+E[1]*exp(-t/lambda[1]);

 

(1)

 

containing the three parameters  

>	E[infinity]:=sigma(infinity)/epsilon[0];  E[1]:=(sigma(0)-sigma(infinity))/epsilon[0]; lambda[1]:=eta[1]/E[1];

 

(2)

 

(2)

 

(2)

 

Substituting the constant strain epsilon[0] in the above relaxation function by epsilon[0]*H(t), we                                   arrive at the representation  

>

restart:

 

>	sigma(t):=epsilon[0]*H(t)*E(t);

 

(3)

 

For a constant strain epsilon[0] within an interval  t = [a, b]  the response can be expressed as  

>	sigma(t):=epsilon[0]*(H(t-a)-H(t-b))*E(t-a);

 

(4)

 

where  H(t)  or  H(t - a)  is the HEAVISIDE unit step function. This response function is illustrated  

in the following Figure.  

>	alias(H=Heaviside,th=thickness,co=color):

 

>	p[1]:=plot(H(t-1)-H(t-6),t=0..7,0..1,co=black,                      title="Relaxation of the POYNTING - THOMSON Model"):

 

>	p[2]:=plot((H(t-1)-H(t-6))*(0.2+0.8*exp(-(t-1)/2)),        t=0..7,th=3,co=black):

 

>	p[3]:=plots[textplot]({[3.6,0.95,`constant  strain`],            [4,0.54,`response`]},co=black):

 

>	plots[display](seq(p[k],k=1..3));

 

 

The POYNTING - THOMSON model is subjected to a constant strain of  height one in  t = [1, 6].                                      Parameters are: epsilon[0] = 1,  E[infinity] = 0.2,  E[1] = 0.8, so that  sigma[0] = 1.  

  Supplementing the POYNTING - THOMSON  model by further MAXWELL elements we arrive                                            at the following  generalized relaxation  function:  

>

restart:

 

>	E(t):=E[infinity]+Sum(E[k]*exp(-t/lambda[k]),k=1..n);

 

(5)

 

 

where E[k] are the spring constants of several MAXWELL elements, while                                                                lambda[k] = eta[k] / E[k] represents the corrosponding relaxation  times.                                                                  Analogous to the discrete creep spectrum (BETTEN, 2008), the values                                                                    (lambda[k],  E[k])  form a discrete relaxation spectrum. The parameter                                                                            E(infinity)  can be interpreted as the long-time equilibrium value of  the                                                                    relaxation modulus  E(t).                                                                                                                                                       The number of  MAXWELL  elements may increase indefinitely (n --> infinity).                                                                  Then the discrete relaxation spectrum proceeds to a continuous one, h = h(lambda).                                                          The relaxation  function results for  n --> infinity in the integral form  

>	restart:  macro(la=lambda,inf=infinity):

 

>	E(t):=E[inf]+beta*Int(h(la)*exp(-t/la),la=0..inf);

 

(6)

 

where the parameter beta can be determined from  E(0)  by assuming normalized relaxation spectra,  

>	Int(h(la),la=0..inf)=1;

 

(7)

 

according to 

>	beta:=E[0]-E[inf];

 

(8)

 

Hence, the relaxation  function is represented in the following dimensionless form with  R(0) = 1: 

>	restart:  macro(la=lambda,inf=infinity):

 

>	R(t):=(E(t)-E[inf])/(E[0]-E[inf])= Int(h(la)*exp(-t/la),la=0..inf);

 

(9)

 

>

 

For the sake of convenience, the integral on the right hand side can be solved by utilizing  

a  LAPLACE  transformation. Thus, we introduce the substitution  1 / lambda  =  xi  and  

arrive at the result: 

>	restart:  macro(la=lambda,inf=infinity):

 

>	R(t):=Int((h(la=1/xi)/xi^2)*exp(-t*xi),xi=0..inf)=L(H(xi));

 

(10)

 

where  L{H(xi)}  is the LAPLACE  transform of  the function 

>

restart:

 

>	H(xi):=(1/xi^2)*h(lambda=1/xi);

 

(11)

 

>

 

In the following some examples should be discussed. At first we select the POISSON                                          distributions as relaxation spectra:  

>	restart:  macro(la=lambda,inf=infinity):

 

>	h[1](la,n):=la^n*exp(-la)/n!;   # n = 0,1,2,...

 

(12)

 

>	for i from 0 to 3 do h[1](la,i):=subs(n=i,h[1](la,n)) od:

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({seq(h[1](la,n),n=0..3)},la=0..5,0..1,              th=3,co=black,axes=boxed,                                          title="POISSON Distributions"):

 

>	p[2]:=plots[textplot]({[0.75,0.8,`n = 0`],[1.5,0.4,`n = 1`],             [2.5,0.3,`n = 2`],[4,0.25,`n = 3`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

Because of   

>	Diff(h[1],la)=simplify(diff(h[1](la,n),la));

 

(13)

 

the curves in the above Figure possess a maximum at  lambda = n > 0.  

  From the relaxation spectra, illustrated in the above Figure, we find the relaxation  functions                                            E(t, n) in dimensionless form  

>	R(t,n):=(E(t,n)-E[inf])/(E[0]-E[inf]);

 

(14)

 

>

 

by utilizing the following MAPLE program including the LAPLACE  transform.  

>	with(inttrans):

 

>	H[1](xi,n):=subs(la=1/xi,h[1](la,n))/xi^2;

 

(15)

 

>	R(t,n):=laplace(H[1](xi,n),xi,t);

 

(16)

 

>	for i in [0,1,2,3] do R(t,n=i):=subs(n=i,R(t,n)) od:

 

>	p[1]:=plot({seq(R(t,n=i),i=0..3)},t=0..5,0..1,                 th=3,axes=boxed,co=black,                                             title="Relaxation  Functions  R(t, n)"):

 

>	p[2]:=plots[textplot]({[1.5,0.78,`n = 3`],                             [1.5,0.1,`n = 0`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

The POISSON  distribution is normalized as we can immediately see by considering  

the gamma function  

>	Gamma(n):=Int(la^(n-1)*exp(-la),la=0..inf); Gamma(n+1):=n!;

 

(17)

 

(17)

 

If  we substitute  (n - 1)  for  n , the POISSON  distribution proceeds to the gamma function.  

  The second example is concerned with the MAXWELL distribution function  

>	h[2](la,a):=(a/2/sqrt(Pi)/(la^(3/2))*exp(-(a^2/4/la)));

 

(18)

 

which is normalized and therefore admissible. This spectrum has been calculated by the  

following MAPLE program.  

>	for i in [3/4,1,3/2] do h[2](la,i):=subs(a=i,h[2](la,a)) od:

 

>	p[1]:=plot({h[2](la,3/4),h[2](la,1),h[2](la,3/2)}, la=0..1,0..1.75,th=3,co=black,axes=boxed,                          title="MAXWELL Distributions as Relaxation Spectra"):

 

>	p[2]:=plots[textplot]({[0.3,1.2,`a = 3/4`],[0.15,0.6,`a = 1`],         [0.3,0.3,`a = 3/2`]},ytickmarks=4,co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

Because of  

>	Diff(h[2],la)=simplify(diff(h[2](la,a),la));

 

(19)

 

the MAXWELL distributions possess a maximum at the position  lambda = (1/6)*a^2.  

The associated  relaxation  functions E(t, a)  in dimensionless form have been calculated  

by utilizing the following MAPLE program including the LAPLACE transform.  

>	with(inttrans):

 

>	H[2](xi):=subs(la=1/xi,h[2](la,a))/xi^2;

 

(20)

 

>	R(t,a):=laplace(H[2](xi),xi,t);

 

(21)

 

>	for i in [3/4,1,3/2] do R(t,i):=subs(a=i,R(t,a)) od:

 

>	p[1]:=plot({R(t,3/4),R(t,1),R(t,3/2)},t=0..5,0..1,                 th=3,co=black, title="Relaxation Functions  R(t, a)"):

 

>	p[2]:=plots[textplot]({[1,0.76,`a = 3/2`],                             [1.1,0.5,`1`],[1,0.25,`a = 3/4`]},co=black,axes=boxed):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

  The third example takes into consideration the Chi-square distripution used as relaxation                                      spectra. This distribution is also normalized and therefore admissible.  

>	restart:  macro(la=lambda,inf=infinity):

 

>	with(Statistics):

 

>	X:=RandomVariable(ChiSquare(n)):

 

>	h[3](la,n):=PDF(X,la) assuming la >= 0;                                         # Probability Density Function

 

(22)

 

>	Int(h[3],la=0..inf)=int(h[3](la,n),la=0..inf);

 

(23)

 

>

 

The Chi-square distributions are used as relaxation spectra and illustrated for several   

parameters  n  in the following Figure.  

>	for i in [2,3,4] do h[3](la,i):=subs(n=i,h[3](la,n)) od:

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({seq(h[3](la,k),k=2..4)},la=0..5,0..0.5,     th=3,co=black,axes=boxed,                                       title="Chi-square Distributions used as relaxation Spectra"):

 

>	p[2]:=plots[textplot]({[1,0.4,`n = 2`],                                [2,0.25,`n = 3`],[4,0.18,`n = 4`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

The curves in this Figure with  n > 2  possess a maximum at the position  lambda = n - 2.  

From the above relaxation spectra we arrive at the relaxation functions  E(t, n) in                                              dimensioless form   

>	R(t,n):=(E(t,n)-E[inf])/(E[0]-E[inf]);

 

 

(24)

 

by utilizing the following MAPLE program including the LAPLACE  transform .  

>	with(inttrans):

 

>	h[3](xi,n):=subs(la=1/xi,h[3](la,n))/xi^2;

 

(25)

 

>	R(t,n):=laplace(h[3](xi,n),xi,t);

 

(26)

 

>	for i in [2,3,4] do R(t,i):=subs(n=i,R(t,n)) od:

 

>	p[1]:=plot({seq(R(t,i),i=2..4)},t=0..5,0..1,            th=3,axes=boxed,co=black,                                            title="Relaxation Functions based upon the Chi-square Distributions"):

 

>	p[2]:=plots[textplot]({[1.5,0.7,`n = 4`],                                [1.5,0.25,`n = 2`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

>	R(inf,n):=Limit(R(t,n),t=inf)=limit(R(t,n),t=inf);

 

(27)

 

>

 

The normalized relaxation spectra lead to dimensionless relaxation functions                                                                 with  R(infinity, n) = 0. 

>

 

  The next example is concerned with relaxation functions based upon the WEIBULL                                        distributions assumed as relaxation spectra. This distribution is characterized by two                                           parameters (a, b), where a is the scale parameter, while b forms the shape of  the curves.  

>	restart: macro(la=lambda,inf=infinity): with(Statistics):

 

>	X:=RandomVariable(Weibull(a,b)):

 

>	h[Weibull](la,a,b):=PDF(X,la);                                                  # PDF = Probability Density Function

 

(28)

 

>	h[4](la,a,3):=subs(b=3,%)  assuming la>=0;

 

(29)

 

>	Int(h[4],la=0..inf)=1; # normalized

 

(30)

 

>	for i in [1,2,3,4] do h[4](la,i,3):=                       subs(a=i,h[4](la,a,3)) od:

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({seq(h[4](la,i,3),i=1..4)},la=0..5,0..1.2, ytickmarks=3,th=3,axes=boxed,co=black,                             title="WEIBULL  Distributions"):

 

>	p[2]:=plots[textplot]({[1.5,1,`a = 1`],[2,0.65,`a = 2`],                  [3,0.45,`a = 3`],[4,0.35,`a = 4`],[4,1,`b = 3`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

Because of   

>	Diff(h[4],la)=simplify(diff(h[4](la,a,3),la));

 

(31)

 

the curves in the above Figure possess a maximum at  lambda^3 = (2/3)*a^3 , i.e.:  

>	for i in [1,2,3,4] do la[opt](a=i):=             evalf(subs(a=i,(2/3)^(1/3)*a),3) od;

 

(32)

 

(32)

 

(32)

 

(32)

 

>

 

  From the relaxation spectra illustrated in the above Figure we find the relaxation  functions  

E(t, a, 3) in dimensionless form  

>	R(t,a,3):=(E(t,a,3)-E[inf])/(E[0]-E[inf]);

 

(33)

 

>	R(t,a,3):=int(h[4](la,a,3)*exp(-t/la),la=0..inf);

 

(34)

 

>	for i from 1 to 4 do R(t,i,3):=subs(a=i,R(t,a,3)) od:

 

>	p[1]:=plot({seq(R(t,i,3),i=1..4)},t=0..5,0..1,   thickness=3,axes=boxed,co=black,                                     title="Relaxation  Functions  R(t, a, 3)"):

 

>	p[2]:=plots[textplot]({[1,0.85,`a = 4`],                                 [1,0.2,`a = 1`],[3.5,0.85,`b = 3`]},co=black):

 

>	plots[display](p[1],p[2]);

 

 

>

 

>

 

 

 

>

 

 

Nonlinearities  

 

In the previous part several relaxation functions have been calculated by using linear rheological  

elements due to MAXWELL. The relaxation of  the MAXWELL  body subjected to  a constant strain  

can be expressed by the relation  

>

restart:

 

>	sigma[MAXWELL](t):=epsilon[0]*E*exp(-t/lambda);

 

(35)

 

In contrast, the following part of  this worksheet is concerned with the modified form 

>	sigma(t):=sigma[modified]=epsilon[0]*E*exp(-c*sqrt(t));

 

(36)

 

>

 

Because of  the good agreement between the modified form, called sqrt(t)-law, and experimental  

results, we assume that creep can be interpreted as a  diffusion controlled process as has been                                 discussed by BETTEN (2008) in more detail. 

  The sqrt(t)-law describes very well the relaxation of  many materials, too. It can be successfully                                     applied to several polymers, e.g., EVE-copolymers at room temperature. The use of  PRONY-series                                  was less successful as has been pointed out in more detail by BETTEN (2008). Furthermore, the                              relaxation of  glass can be expressed by the sqrt(t)-law. AIMEDIEU (2004) has investigated the                              nonlinear relaxation of  brain tissue and found very good agreement between the sqrt(t)-law and                                        his own experiments. The use of  PRONY-series was again less successful.  

  The number of  modified MAXWELL elements may increase indefinitely (n --> infinity). Then the                                     relaxation  function in dimensionless form is given by  

>	restart:  macro(la=lambda,inf=infinity):

 

>	N(t):=(E(t)-E[inf])/(E[0]-E[inf])=           Int(h(la)*exp(-sqrt(t)/la),la=0..inf);

 

(37)

 

where  h(lambda)  are normalized relaxation spectra dicussed in more detail in the previous part  

of  this worksheet. In the following some examples should be discussed. At first we select the  

POISSON distributions as relaxation spectra and find the following modified relaxation functions: 

>	restart:  macro(la=lambda,inf=infinity):

 

>	h[1](la,n):=la^n*exp(-la)/n!;

 

(38)

 

>	N(t,n):=int(h[1](la,n)*exp(-sqrt(t)/la),la=0..inf);

 

(39)

 

>	for i in [0,1,2,3] do N(t,n=i):=subs(n=i,N(t,n)) od:

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({seq(N(t,n=i),i=0..3)},t=0..5,0..1,                th=3,axes=boxed,co=black,                                         title="Modified Relaxation Functions N(t, n)"):

 

>	p[2]:=plots[textplot]({[1.55,0.78,`n = 3`],                           [1.55,0.17,`n = 0`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

In contrast to the preveous relaxation functions  R(t, n)  the modified functions  N(t, n)  have  

a vertical tangent at the beginning t = 0, which can often be observed in experiments on 

several materials (BETTEN , 2008).  

  In the following Figure the relaxation functions  R(t, n)  and  N(t, n) have been compared.  

>

restart:

 

>	R(t,n):=2*t^(1/2+n/2)*BesselK(n+1,2*sqrt(t))/n! ;                               # linear theory

 

(40)

 

>	for i in [0,3] do R(t,n=i):=subs(n=i,R(t,n)) od;

 

(41)

 

(41)

 

>	N(t,n):=2*t^(1/4+n/4)*BesselK(n+1,2*t^(1/4))/n! ;                               # non-linear theory

 

(42)

 

>	for i in [0,3] do N(t,n=i):=subs(n=i,N(t,n)) od;

 

(43)

 

(43)

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({R(t,n=0),R(t,n=3)},t=0..5,0..1,  linestyle=4,th=3,axes=boxed,co=black):

 

>	p[2]:=plot({N(t,n=0),N(t,n=3)},t=0..5,0..1,th=3,co=black,                 title=" Linear and Non-linear Theory"):

 

>	p[3]:=plots[textplot]({[1,0.9,`n = 3`],                                 [1,0.45,`n = 0`]},co=black):

 

>	plots[display](seq(p[k],k=1..3));

 

 

>

 

In this Figure the linear theory is characterized by dashed lines, while the solid lines belong  

to the non-linear theory.  

>

 

The next example is concerned with finding the modified relaxation functions based upon the  

WEIBULL distributions assumed as relaxation spectra.  

>	restart:  macro(la=lambda,inf=infinity):

 

>	h[4](la,a,3):=3*la^2*exp(-la^3/a^3)/a^3;

 

(44)

 

>	N(t,a,3):=int(h[4](la,a,3)*exp(-sqrt(t)/la),la=0..inf);                         # non-linear theory

 

(45)

 

>	for i from 1 to 4 do N(t,i,3):=subs(a=i,N(t,a,3)) od:

 

>	alias(th=thickness,co=color):

 

>	p[1]:=plot({seq(N(t,i,3),i=1..4)},t=0..5,0..1,          th=3,axes=boxed,co=black,                                         title="Modified Relaxation Functions N(t, a, 3)"):

 

>	p[2]:=plots[textplot]({[1,0.82,`a = 4`],[1,0.45,`a = 2`],                [1,0.2,`a = 1`],[4,0.9,`b = 3`]},co=black):

 

>	plots[display](seq(p[k],k=1..2));

 

 

>

 

In contrast to the previous relaxation functions  R(t, a, 3)  the modified functions  N(t, a, 3)  

have a vertical tangent at the beginning  t = 0, which can often be observed in experiments  

on several materials, e.g.  EVE copolymers (BLOK, 2006) and many other materials                                                               (BETTEN, 2008). 

  In the following Figure the relaxation functions R(t, a, 3) and N(t, a, 3) have been compared. 

>

 

>	R(t,a,3):=t^2*sqrt(3)*MeijerG([[],[]],[[1/3,0,-1/3,-2/3],[]], t^3/27/a^3)/18/a^3/Pi/(1/a^3)^(1/3); # linear theory

 

(46)

 

>	for i in [1,4] do R(t,i,3):=subs(a=i,R(t,a,3)) od;

 

(47)

 

(47)

 

 

non-linear theory: 

>	for i in [1,4] do N(t,i,3):=subs(a=i,N(t,a,3)) od;                               # non-linear theory

 

(48)

 

(48)

 

>	p[1]:=plot({R(t,1,3),R(t,4,3)},t=0..5,0..1,   linestyle=4,th=3,axes=boxed,co=black):

 

>	p[2]:=plot({N(t,1,3),N(t,4,3)},t=0..5,0..1,th=3,co=black,           title="Linear and Non-linear Theory"):

 

>	p[3]:=plots[textplot]({[1,0.85,`a = 4`],                                [1,0.5,`a = 1`],[4,0.9,`b = 3`]},co=black):

 

>	plots[display](seq(p[k],k=1..3));

 

 

>

 

 

Experiments 

 

Based upon a lot of experiments on glas, SCHERER (1986) has shown that both the stress  

and structural relaxation in glass can be predicted by the relation 

>

restart:

 

>	r(t):=exp(-(t/lambda)^b);

 

(49)

 

>

 

often called KOHLRAUSCH function or b-function, which is a modified form of  the  

MAXWELL relaxation function. The exponent  b  in the above relation was found to be  

near the value of  b = 0.5, so that the assumption of  the sqrt(t)-law is justified. However,  

the above relation is valid only for stabilized glass, i.e., the glass is held at a given  

temperature until its properties do no longer change with time. Then the load can be  

applied.  

  In unstabilized glasses the viscosity and other typical properties, e.g. the density,  

vary with time. Then, the above relaxation function should be replaced by the formular 

>

restart:

 

>	r(t):=exp(-(Int(G[0]/eta(tau),tau=0..t)^b));

 

(50)

 

where, in agreement with experimental results, the exponent  b  can again be assumed  

near to  b = 0.5 , as has been discussed in more detail by SCHERER (1986).  

  AIMEDIEU (2004) has investigated the nonlinear relaxation of  brain tissue and found  

very good agreement between the sqrt(t)-law and his own experiments. The use of  PRONY  

series was less succesful. 

 

Conclusion  

 

The examples in this worksheet are concerned with the representation of  relaxation curves                                                for both the linear and non-linear theory. The linear theory is based upon linear rheological                                       elements due to MAXWELL , while the non-linear theory is characterized by the sqrt(t)-law.  

Besides very good agreement with experimental results, the sqrt(t)-law has a physical meaning,  

as has been discussed by BETTEN (2008) in more detail. 

 

References  

 

AIMEDIEU, P. (2004). Contributation ? la biom?canique de tissues mous intr?craniens,                                                          PhD Thesis of  the Universit? de Picardie Jules Vernes, Faculte de medicine,                                                                       Amiens Cedex, France.                                                                                                                                                                                                                                                                                                                                                         BETTEN, J.(2008). Creep Mechanics, 3rd Edition, Springer-Verlag, Berlin / Heidelberg . 

BLOK, A. (2006). Vergleich rheologischer Modelle zum Kriechverhalten eines                                            EVA-Copolymers, Presentation on the occation of  doctorial examination,                                                                  RWTH Aachen University.  

SCHERER, G. (1986). Rlaxation in Glass and Composites, John Wiley & Sohn, New York... 

>