**Journal of Modern Physics**

Vol.08 No.12(2017), Article ID:80594,25 pages

10.4236/jmp.2017.812120

*Dedicating this work to my mother Purabi Das on her 83^{rd} birthday.

Revisiting the Curie-Von Schweidler Law for Dielectric Relaxation and Derivation of Distribution Function for Relaxation Rates as Zipf’s Power Law and Manifestation of Fractional Differential Equation for Capacitor^{*}

Shantanu Das^{1,2}^{ }

^{1}Scientist-Reactor Control System, Design Section (E & I Group), Bhabha Atomic Research Centre (BARC), Mumbai, India

^{2}Honorary Senior Research Professor, Condensed Matter Physics Research Centre (CMPRC), Department of Physics, Jadavpur University, Kolkata, India

Copyright © 2017 by author and Scientific Research Publishing Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

Received: October 11, 2017; Accepted: November 24, 2017; Published: November 27, 2017

ABSTRACT

The classical power law relaxation, i.e. relaxation of current with inverse of power of time for a step-voltage excitation to dielectric―as popularly known as Curie-von Schweidler law is empirically derived and is observed in several relaxation experiments on various dielectrics studies since late 19^{th} Century. This relaxation law is also regarded as “universal-law” for dielectric relaxations; and is also termed as power law. This empirical Curie-von Schewidler relaxation law is then used to derive fractional differential equations describing constituent expression for capacitor. In this paper, we give simple mathematical treatment to derive the distribution of relaxation rates of this Curie-von Schweidler law, and show that the relaxation rate follows Zipf’s power law distribution. We also show the method developed here give Zipfian power law distribution for relaxing time constants. Then we will show however mathematically correct this may be, but physical interpretation from the obtained time constants distribution are contradictory to the Zipfian rate relaxation distribution. In this paper, we develop possible explanation that as to why Zipfian distribution of relaxation rates appears for Curie-von Schweidler Law, and relate this law to time variant rate of relaxation. In this paper, we derive appearance of fractional derivative while using Zipfian power law distribution that gives notion of scale dependent relaxation rate function for Curie-von Schweidler relaxation phenomena. This paper gives analytical approach to get insight of a non-Debye relaxation and gives a new treatment to especially much used empirical Curie-von Schweidler (universal) relaxation law.

**Keywords:**

Power Law, Relaxation Rate Distribution, Fractional Derivative, Fractional Integration, Curie-Von Schweidler Law, Time-Constants, Laplace Integral, Zipf’s Law, Integral Representation, Time Dependent Relaxation Rate, Scale Dependent Relaxation Rate, Non-Debye Relaxation

1. Introduction

The Curie-von Schweidler law relates to relaxation current in dielectric when a step DC voltage is applied and is given by
$i\left(t\right)\sim {t}^{-n}$ , where
$t>0$ and the power (exponent) i.e. n is called relaxation constant or decay constant, where
$0<n<1$ [1] [2] [3] [4] . We note that n is non-integer. This relaxation law is taken as universal law, at least for dielectric relaxations. Whereas we are used to Debye type of relaxation i.e. exponential decay law given by
$i\left(t\right)\sim {\text{e}}^{-t/{\tau}_{0}}$ or
$i\left(t\right)\sim {\text{e}}^{-{\lambda}_{0}t}$ where
${\tau}_{0}$ is the relaxation time constant while
${\lambda}_{0}$ denotes the relaxation rate of the process with
${\lambda}_{0}={\tau}_{0}^{-1}$ . The radioactive decay is example of ideal Debye law where the exponential decay is governed by “one-lumped” decay constant i.e.
${\lambda}_{0}$ . The Curie-von Schweidler behavior has been observed in many instances, since late 19^{th} Century, such as those shown in dielectric studies and experiments [3] - [10] .

This power law relaxation of the non-Debye type i.e. $i\left(t\right)\sim {t}^{-n}$ has been interpreted as a many-body problem but can also be formulated as an infinite number of independent relaxing bodies meaning infinite number of time constants $\tau $ or relaxation rates $\lambda $ varying from near zero to infinity [11] [12] . The observations of power law relaxation are also made in the experiments and studies with super-capacitors [13] [14] [15] [16] [17] . These studies also indicate the fractional calculus is used as constituent expression to describe super-capacitors. The use of empirical power law i.e. Curie-von Schweidler Law of relaxation of current to a step input of voltage to get constituent relation with fractional derivative was proposed in [5] . Apart from relaxation of current decay in dielectrics and super-capacitors, the power law type or non-Debye relaxation is observed in visco-elastic experiments strain relaxation in [18] [19] [20] [21] .

In this paper, we are giving the derivation of the distribution of relaxation rates ( $\lambda $ ) particularly for Curie-von Schweidler law and we observe the distribution nature as Zipf’s distribution [22] - [27] . We try to reason out as to why this distribution of relaxation rates takes Zipfian nature. We also show that Curie-von Schweidler law has time varying rate of relaxation. This paper will not deal with the mathematics of Zipfian distribution (or power law distribution) like probability density function, cumulative probability density function, and the conditions of finding finite mean, variance or standard deviation for power law distribution. This paper describes finding the distribution function of relaxation rates (or histogram) by formulating Laplace integral, and show that the distribution thus obtained is a Zipf’s power law.

We extend this mathematical approach to get the distribution function for time constants ( $\tau $ ). We observe that time constants $\tau $ are also distributed as Zipf’s power law; but this observation points to a contrary physical interpretation derived from this obtained power law distribution for the relaxation rate ( $\lambda $ ) distributions. Thus we can conclude that this method developed by Laplace integral approach is restricted to get only distribution of relaxation rates i.e. $\lambda $ and not to get the distribution of time constants i.e. $\tau $ . Though we are discussing especially Curie-von Schweidler law, yet we will tabulate relaxation rate distributions obtained for some other relaxation functions which are obtained via this Laplace integral method.

We shall demonstrate the formation of fractional derivative in the expression relating current and voltage considering the relaxation rates as Zifian distribution; and thus forming a scale dependent power law for relaxation rates as the scale varies from zero to infinity. Though by experiments one cannot make histogram directly for the rates of relaxation for any non-Debye processes, yet this mathematical procedure that we develop helps in extracting this information from the observations relaxation function. This is new treatment, and much more research is required, across various dynamic processes.

2. Obtaining Fractional Derivative Directly from Curie-Von Schweidler Law for Capacitor

Practically on applying a step input voltage $v\left(t\right)={V}_{BB}$ Volts at $t=0$ to a capacitor which is initially uncharged; we get a power-law decay of current given by empirical Curie-von Schweidler as $i\left(t\right)\sim {t}^{-n};0<n<1$ [5] . That we write in following way as indicated by experimental studies [5] - [10] :

$i\left(t\right)={K}_{n}\frac{{V}_{BB}}{{t}^{n}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}t>0$ (1)

The parameter
${K}_{n}$ is constant. This is from observation and the evaluation of order of power-law function is
$0.5<n<1$ [5] [13] [14] [15] [16] [17] . Let the capacitor be excited by a step input of V_{BB} Volts, i.e. written as
$v\left(t\right)={V}_{BB}\left(u\left(t\right)\right)$ , where
$u\left(t\right)$ is unit step function at time
$t=0$ . The Laplace transform of step input is
$V\left(s\right)=\mathcal{L}\left\{v\left(t\right)\right\}=\mathcal{L}\left\{{V}_{BB}\left(u\left(t\right)\right)\right\}={V}_{BB}/s$ . Then taking Laplace transform of above power-law decay current (1), we obtain

$I\left(s\right)=\mathcal{L}\left\{i\left(t\right)\right\}=\mathcal{L}\left\{{K}_{n}{V}_{BB}{t}^{-n}\right\}={K}_{n}{V}_{BB}\left(\frac{\left(-n\right)!}{{s}^{-n+1}}\right)$ [28] . Then using the formula

for generalization of factorial i.e. $\left(\alpha -1\right)!=\Gamma \left(\alpha \right)$ [28] [29] [30] , we get the following expressions

$I\left(s\right)={K}_{n}\frac{\Gamma \left(1-n\right){V}_{BB}}{{s}^{1-n}}={K}_{n}\frac{\Gamma \left(1-n\right)}{{s}^{-n}}\left(\frac{{V}_{BB}}{s}\right)$ (2)

We get Transfer function [28] of capacitor as following expression

$\begin{array}{c}G\left(s\right)=\frac{I\left(s\right)}{V\left(s\right)}=\frac{{K}_{n}\frac{\Gamma \left(1-n\right)}{{s}^{-n}}\left(\frac{{V}_{BB}}{s}\right)}{\left(\frac{{V}_{BB}}{s}\right)}\\ ={K}_{n}\left(\Gamma \left(1-n\right)\right){s}^{n}={C}_{n}{s}^{n},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{C}_{n}={K}_{n}\left(\Gamma \left(1-n\right)\right)\end{array}$ (3)

This Expression (3) i.e. $G\left(s\right)$ is admittance expression in complex frequency s-domain of a capacitor. From here we write impedance expression for capacitor as following

$Z\left(s\right)=\frac{1}{{C}_{n}{s}^{n}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1$ (4)

From the obtained Expressions (3) (4) i.e. $I\left(s\right)={C}_{n}{s}^{n}\left(V\left(s\right)\right)$ and by taking inverse Laplace transform by using the identity ${\mathcal{L}}^{-1}\left\{{s}^{n}F\left(s\right)\right\}={}_{0}D{}_{t}^{n}\left[f\left(t\right)\right]$ i.e. fractional derivative operation [12] [31] , we get the constituent relation for capacity as following

$i\left(t\right)={C}_{n}\left({}_{0}D{}_{t}^{n}\left[v\left(t\right)\right]\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1$ (5)

This fractional derivative Expression (5) gives a new capacitor theory [5] and we utilize this above Formula (5) to find characteristics of super-capacitors, like the variation of n with the current excitation, and the efficiency of the energy discharged to the energy stored [13] [14] [15] [16] [17] . Classically the expression of capacitor is $i\left(t\right)=C\left({D}_{t}^{\left(1\right)}\left[v\left(t\right)\right]\right)$ i.e. with integer “one-whole” order classical derivative. Therefore Curie-von Schweidler law gives a different approach based on fractional calculus [12] [31] . In experimental observations we find that capacitor has fractional order impedance [5] - [10] [13] [14] [15] [16] , [17] . The impedance $Z\left(\omega \right)\sim {\omega}^{-n},0<n<1$ (is obtained by writing Laplace variable s in (4) as $i\omega ;i=\sqrt{-1}$ , i.e. considering steady state analysis [5] [28] . This fractional impedance observed in [5] - [10] [13] [14] [15] [16] [17] , has implication in dissipation [5] theory of di-electrics that we will not cover here. We state that while classical capacitor unit is in Farad, the ${C}_{n}$ the “fractional capacity” is in units of $Farad/{\mathrm{sec}}^{1-n}$ [5] .

This section gives us the understanding that this Curie-von Schweidle law i.e. the empirical law gives a relation of voltage and current of capacitor by using fractional derivative. In this paper we will show how we get the same relation (i.e. via use of fractional derivative) by considering Zipfian distribution of relaxation rates ( $\lambda $ ) that we get for Curie-von Schweidler relaxation law.

3. About the Zipf’s Power Law Distribution and Probable Hypothesis for Its Mechanism

The Zipf’s law is widely referred in linguistic studies, economics studies, population studies [22] - [27] . We use this for a dielectric relaxation law (i.e. Curie-von Schweidler law), which is observed as
$i\left(t\right)\sim {t}^{-n},0<n<1$ , since late 19^{th} century. We derived histogram of relaxation rates for relaxation function
$i\left(t\right)\sim {t}^{-n}$ and show that it follows Zipf’s power law. We try to give possible reasons as to why Zipfian distribution is observed for the distribution of relaxation rates. The histogram function of Zipf’s law is
$H\left(x\right)\sim {x}^{-m}$ a power law type. In this section we assume that relaxation rates l’s follow a Zipfian histogram say
${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{-m}$ . This we will derive in subsequent section. The l’s are relaxation rates of infinite number of relaxing bodies, simultaneously relaxing as per Debye law i.e.
$\sim {\text{e}}^{-\lambda t}$ .

Zipf’s Law is an empirical law formulated using mathematical statistics that refers to the fact that many types of data studied in the physical and social sciences can be approximated with a Zipfian distribution. This distribution is one of a family of related discrete power law probability distributions [22] - [27] . This power law distribution help to describe phenomena where large events are rare, but small ones are quite common. For example, there are few large earthquakes but many small ones. There are a few mega-cities, but many small towns. There are few words, such as “and” and “the” that occur very frequently, but many which occur rarely.

The emergence of a complex language is one of the fundamental events of human evolution, and several remarkable features suggest the presence of fundamental principles of organization. These principles seem to be common to all languages. The best known is the so-called Zipf’s law, which states that the frequency of a word decays as a (universal) power law of its rank. The possible origins of this law have been controversial, and its meaningfulness is still an open question. One of the early hypotheses of Zipf of a principle of least effort for explaining the law is shown to be sound [26] [27] . But still the exact mechanism how the Zipf’s distribution manifests is debated.

Many of the things that we measure have a typical size or “scale”. We ask ourselves why the relaxation rates $\lambda $ cannot be arranged as simple “normal distribution”. Like while we plot the height of person in X-axis and the percentage of occurrence of that particular height in Y-axis, we get a “normal distribution” peaked around mean height with a spread both ways, that is a histogram. We find that ratio of maximum height and minimum height of a person is finite (or relatively low value). For example as per Guinness book of records tallest person was having height 272 cm and shortest person was having the height of 57 cm, making this ratio 4.8. This ratio is relatively low value. We see the most adults are about 170 cm tall-there is some variations around this figure notably depending on sex, but we never measure persons having height of 1 cm or 1000 cm.

But not all things we measure are peaked around a typical value. Some may over a very large dynamic range, sometimes many orders of magnitude. For example the ratio of population of largest town to population of smallest town is about 250,000. The histogram if plotted for X-axis with population of cities and Y-axis with percentage of cities having that population; the distribution will not show the “normal-distribution”. The histogram of cities & population is highly “right-skewed”, meaning that while the bulk of distribution occurs for fairly small sizes―i.e. most cities have small population-there is small number of cities with population much higher than a said typical value, producing the long tail to the right of histogram. This “right skewed” form is qualitatively quite different from histogram of person’s height. That is because we know that there is large dynamic range from smallest to largest city sizes, we can immediately infer that there can only be a small number of very large cities. The histogram of this sort is like a function i.e. ${H}_{x}\left(x\right)\sim {x}^{-m}$ . The distribution of this nature is called Zip’s power law distribution.

The same we observe when relaxation rates call them $\lambda $ having a large (ideally infinite) spreads follow Zipfian distribution, we call that ${H}_{\lambda}\left(\lambda \right)$ and will show that the histogram follows the function i.e. ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{-m}$ . Thus one reason that this non-Debye relaxation (explained in subsequent sections) of Curie-von Schweidler Law ( $i\left(t\right)\sim {t}^{-n}$ ) in dielectric is having infinite spread of relaxation rates of l’s-thus forming a Zipfian power law.

4. Zipfian Power Law Distribution as a Result of Connected Exponential Processes―A Postulate

Having discussed the formation of a histogram as power law type, when there is very large dynamic spreads amongst the relaxation rates of a complex relaxing process we move to a probable postulate of explaining this process via exponentially distributed processes. A much more common distribution than power law is the exponential distribution. In this complex relaxation mechanism i.e. $i\left(t\right)\propto {t}^{-n}$ that we are discussing we consider infinite number of bodies relaxing simultaneously, in different time scales ( $T$ ). We consider that a complex relaxation mechanism and a quantity $T$ say survival time of a relaxing body, has exponential distribution of probability $p\left(T\right)\sim {\text{e}}^{-aT}$ . This means that a probability for a body having very large survival time (age) is very low; and vice-versa. Then $p\left(T\right)dT$ indicates the fraction of survival numbers of bodies between survival time $T$ and $T+dT$ . Now suppose that the real quantity that we are interested is not $T$ but other quantity $\lambda $ , say the relaxation rate of discharge which is exponentially related to $T$ ; thus $\lambda \sim {\text{e}}^{-bT}$ . That implies the surviving bodies with very large time of survival (age) have a very low rate of relaxation. This also states that $d\lambda =-dT$ . Then if probability distribution of $\lambda $ is $p\left(\lambda \right)$ ; then we have $p\left(\lambda \right)d\lambda =-p\left(T\right)dT$ (effect of conservation of probability [25] ). The negative sign indicates opposite movement, as $T$ is increased from $T$ to $T+dT$ , then $\lambda $ is decreased from $\lambda $ to $\lambda +\left(-d\lambda \right)$ . This means that number of discharging units having relaxation rates between $\lambda $ and $\lambda +d\lambda $ is equal to number of surviving bodies having survival time between $T$ and $T+dT$ . Thus we write following steps:

$\begin{array}{c}p\left(\lambda \right)=-p\left(T\right)\frac{dT}{d\lambda}=-\frac{p\left(T\right)}{\left(\frac{d\lambda}{dT}\right)}\sim \frac{{\text{e}}^{-aT}}{b{\text{e}}^{-bT}}=\frac{{\text{e}}^{\frac{-a}{-b}\left(-bT\right)}{\text{e}}^{bT}}{b}\\ =\frac{{\lambda}^{\frac{a}{b}}{\lambda}^{-1}}{b}\sim {\lambda}^{-\left(1-\frac{a}{b}\right)}\sim {\lambda}^{-m};\text{\hspace{0.17em}}\text{\hspace{0.17em}}m=1-\frac{a}{b}\end{array}$ (6)

The above discussion in steps (6) gives a power law distribution for relaxation rates l’s where there is combination of exponential processes. Thus we expect that in our complex relaxation process governed by Curie-von Scweidler Law $i\left(t\right)\propto {t}^{-n}$ which is having infinite number of simultaneously discharging bodies will have a power law distribution for relaxation rates as a histogram ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{-m}$ . This we will derive subsequently. We proceed with this explanation and hypothesis. This could be one explanation in physical sense, in line with exponential distribution in the Boltzmann distribution of energies in statistical mechanics.

5. Complex Relaxation of Non-Debye Type Composing with Several Exponentially Relaxing Functions of Debye Type

We call the Curie-von Shweidler relaxation law $i\left(t\right)\sim {t}^{-n};0<n<1$ as complex process, of non-Debye type. Where a Debye type relaxation is a decay function given by function of exponential type as $\sim {\text{e}}^{-{\lambda}_{0}t};{\lambda}_{0}>0$ . We mention that Curie-von Shweidler relaxation law is not exponential relaxation process of Debye type.

In this section we formulate the method to extract the histogram of the relaxation rates call it
${H}_{\lambda}\left(\lambda \right)$ , for a complex non-Debye relaxation process
$i\left(t\right)$ , which we assume to be composed of several Debye type exponential relaxation functions
${\text{e}}^{-\lambda t}$ , with
$\lambda $ varying from zero to infinity. The complex decay may be expressed as following with several rate constants
${\lambda}_{1},{\lambda}_{2},{\lambda}_{3},\cdots $ with weights
${a}_{1},{a}_{2},{a}_{3},\cdots $ , where
${\lambda}_{k}$ is having units in sec^{−}^{1} i.e. “per second”, and is equal to inverse of time constant i.e.
${\lambda}_{k}={\left({\tau}_{k}\right)}^{-1};k=1,2,3,\cdots $ . We write following composite relaxation expression as sum of several “discrete” relaxations of Debye type i.e.

$\begin{array}{l}i\left(t\right)={a}_{1}{\text{e}}^{-{\lambda}_{1}t}+{a}_{2}{\text{e}}^{-{\lambda}_{2}t}+{a}_{3}{\text{e}}^{-{\lambda}_{3}t}+\cdots ={\displaystyle \sum {a}_{k}{\text{e}}^{-{\lambda}_{k}t}}\\ i\left(0\right)={a}_{1}+{a}_{2}+{a}_{3}+\cdots \end{array}$ (7)

The coefficients ${a}_{k}$ ’s in (7) can be positive or negative that we will elucidate in later section. In the continuum limit, we may write the above discrete expression as following integral equation

$i\left(t\right)={\displaystyle \underset{0}{\overset{\infty}{\int}}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-\lambda t}d\lambda}$ (8)

Where the function i.e. ${H}_{\lambda}\left(\lambda \right)$ is the distribution-function of the rate of the relaxation ( $\lambda $ ) of the process, or we may call it as histogram of relaxation rates. The function ${H}_{\lambda}\left(\lambda \right)$ can be positive or negative that we will elucidate in later section.

While for the case with discrete set of relaxation rates i.e. ${\lambda}_{1},{\lambda}_{2},{\lambda}_{3},\cdots $ the rate distribution function would be having discrete delta functions ( $\delta \left(t-{\lambda}_{k}\right),k=1,2,3,\cdots $ ) at points ${\lambda}_{1},{\lambda}_{2},{\lambda}_{3},\cdots $ ; which we write like following expression

$\begin{array}{c}{H}_{\lambda}\left(\lambda \right)={a}_{1}\left(\delta \left(\lambda -{\lambda}_{1}\right)\right)+{a}_{2}\left(\delta \left(\lambda -{\lambda}_{2}\right)\right)+{a}_{3}\left(\delta \left(\lambda -{\lambda}_{3}\right)\right)+\cdots \\ ={\displaystyle \sum {a}_{k}\left(\delta \left(\lambda -{\lambda}_{k}\right)\right)}\end{array}$ (9)

From above formulation (9) we infer that if we have only one single Debye relaxation i.e. having only one rate constant say ${\lambda}_{0}$ i.e. $i\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ then ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda -{\lambda}_{0}\right)$ . This is verified in the following expression

$i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-\lambda t}d\lambda}={\displaystyle {\int}_{0}^{\infty}\left(\delta \left(\lambda -{\lambda}_{0}\right)\right){\text{e}}^{-\lambda t}d\lambda}={\text{e}}^{-{\lambda}_{0}t}$ (10)

In above derivation (10) we used the property of delta function [29] [30] [32] i.e. $\int \left(\delta \left(x-{x}_{0}\right)\right)\left(f\left(x\right)\right)dx}=f\left({x}_{0}\right)$ .

6. Extraction of Rate Distribution Function by Formulating Laplace Integral

In this section we formulate Laplace integral of the complex decay given by Curie-von Shweidler relaxation law $i\left(t\right)\sim {t}^{-n};0<n<1$ and then getting by inverse Laplace transform of time domain response i.e. ${\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}$ we get relaxation rate distribution function i.e. ${H}_{\lambda}\left(\lambda \right)$ . Conventionally we are used to get inverse Laplace transform of a frequency domain function to time domain function; we note here we will be inverting a time domain function i.e. ${\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}$ .

The Laplace transform $F\left(s\right)$ of a function in time domain $f\left(t\right)$ is defined as following integral transform relation [28] [29] [30] , i.e. called Laplace integral

$\begin{array}{l}F\left(s\right)\stackrel{def}{=}{\displaystyle \underset{0}{\overset{\infty}{\int}}\left(f\left(t\right)\right){\text{e}}^{-st}}dt,\text{\hspace{0.17em}}\text{\hspace{0.17em}}F\left(s\right)=0\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}for\text{\hspace{0.17em}}s<0\\ F\left(s\right)=\mathcal{L}\left\{f\left(t\right)\right\}\end{array}$ (11)

This Expression (11) is standard integral transform of a function $f\left(t\right)$ from a time domain (t) to a complex frequency domain i.e. $s=\mathrm{Re}\left\{s\right\}+i\omega ;i=\sqrt{-1}$ ; where real part is significant in the transient response and the imaginary part of the frequency corresponds to “steady-state” response; in classical “Control Science” [28] . Here $f\left(t\right)$ is “inverse Laplace transform” of $F\left(s\right)$ , and we write ${\mathcal{L}}^{-1}\left\{F\left(s\right)\right\}=f\left(t\right)$ and $\mathcal{L}\left\{f\left(t\right)\right\}=F\left(s\right)$ .

We have in earlier Section 8 derived $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-\lambda t}d\lambda}$ . Compare this with defined Laplace integral expression as follows

$i\left(t\right)={\displaystyle \underset{0}{\overset{\infty}{\int}}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-t\lambda}}d\lambda ,\text{\hspace{0.17em}}\text{\hspace{0.17em}}F\left(s\right)={\displaystyle \underset{0}{\overset{\infty}{\int}}\left(f\left(t\right)\right){\text{e}}^{-st}}dt$ (12)

Both expressions in (12) above are Laplace transform expressions, (or Laplace integrals). The first expression is transforming the function ${H}_{\lambda}\left(\lambda \right)$ from $\lambda $ domain to “complex” t time domain; while the second one is transforming $f\left(t\right)$ from t domain to “complex” s frequency domain. Thus both expressions are Laplace integral expressions with change of variable and symbol. Therefore we can say ${H}_{\lambda}\left(\lambda \right)$ is inverse Laplace Transform of $i\left(t\right)$ in the first expression, i.e. ${H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}$ ; as $f\left(t\right)$ is inverse Laplace of $F\left(s\right)$ in the second expression, i.e. $f\left(t\right)={\mathcal{L}}^{-1}\left\{F\left(s\right)\right\}$ .

Therefore in order to get the rate distribution-function ${H}_{\lambda}\left(\lambda \right)$ from the decay curve (or relaxation-function $i\left(t\right)$ ), we need to perform inverse Laplace Transform of the time function $i\left(t\right)$ . The definition of inverse Laplace Transform is described as following integral expressions

$f\left(t\right)=\frac{1}{2\text{\pi}i}{\displaystyle \underset{x-i\infty}{\overset{x+i\infty}{\int}}\left(F\left(s\right)\right){\text{e}}^{st}}ds,\text{\hspace{0.17em}}\text{\hspace{0.17em}}{H}_{\lambda}\left(\lambda \right)=\frac{1}{2\text{\pi}i}{\displaystyle \underset{x-i\infty}{\overset{x+i\infty}{\int}}\left(i\left(t\right)\right){\text{e}}^{t\lambda}}dt$ (13)

In the above Expression (13) $x$ is real number larger than ${x}_{0}$ , where ${x}_{0}$ being such that $i\left(t\right)$ has some form of singularity on the real line $\mathrm{Re}\left\{t\right\}={x}_{0}$ but is analytic in the complex plane to the right of that line, i.e. for $\mathrm{Re}\left\{t\right\}>{x}_{0}$ , [28] [29] [30] . Thus in this formulation we treat time variable as complex quantity say $t=x+iy$ in the Expression (13) of inverse Laplace Transform i.e.

${H}_{\lambda}\left(\lambda \right)=\frac{1}{2\text{\pi}i}{\displaystyle {\int}_{x-i\infty}^{x+i\infty}\left(i\left(t\right)\right){\text{e}}^{t\lambda}d\lambda}$ . Though we cannot explain presently physical

meaning of concept an “imaginary time” in the expression of complex time quantity $t=x+iy$ , yet mathematically there is no restriction in assuming time to be complex number. We thus proceed in mathematical sense to invert a function in complex time variable, by techniques of Laplace inversion.

7. Extracting the Rate Distribution Function by Performing Inverse Laplace Transform on Relaxation Function of Time Variable

Table 1 lists several types’ relaxation functions $i\left(t\right)$ and its inverse Laplace ${H}_{\lambda}\left(\lambda \right)$ describing the rate distribution function; mostly got from standard Laplace transform tables [28] . The integral representations of ${H}_{\lambda}\left(\lambda \right)$ , shown in Table 1 i.e. for entries 12 to 16 is got via Berberan-Santos method [33] [34] . The entry 12 is for stretched exponential decay function and entry 13 is Becquerel’s compressed hyperbolic radioactive decay function; the entry 15 and 16 is for Mittag-Leffler function and the entry 14 is general power law relaxation. These integral representations of ${H}_{\lambda}\left(\lambda \right)$ are difficult to solve but are easy to plot via use of numerical integration techniques.

We have observed in the previous section that a Debye relaxation of $i\left(t\right)\sim {\text{e}}^{-{\lambda}_{0}t}$ has rate distribution as ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda -{\lambda}_{0}\right)$ i.e. it is given by a delta function at point $\lambda ={\lambda}_{0}$ (10). This we verify with known Laplace relation i.e. $\mathcal{L}\left\{f\left(t-{t}_{0}\right)\right\}={\text{e}}^{-s{t}_{0}}F\left(s\right)$ [28] , where $\mathcal{L}\left\{f\left(t\right)\right\}=F\left(s\right)$ . Also we have $\mathcal{L}\left\{\delta \left(t\right)\right\}=1$ ; thus we can write $\mathcal{L}\left\{\delta \left(t-{t}_{0}\right)\right\}={\text{e}}^{-s{t}_{0}}$ . From here we can write with change of variable for $i\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ the inverse Laplace of this time domain function in $\lambda $ domain we get as ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda -{\lambda}_{0}\right)$ , i.e. the rate distribution function.

If there is no decay then say $i\left(t\right)=1={\text{e}}^{-{\lambda}_{0}t};{\lambda}_{0}=0$ ; the rate distribution function is delta function at origin i.e. ${H}_{\lambda}\left(\lambda \right)=\delta \left(t\right)$ .

If the relaxation function is of say $i\left(t\right)=1-{\text{e}}^{-{\lambda}_{0}t}$ ; then we have

${H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{1-{\text{e}}^{-{\lambda}_{0}t}\right\}$ ; giving ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda \right)-\delta \left(\lambda -{\lambda}_{0}\right)$ . From this observation we say that for our earlier Expressions (7) and (9) i.e. ${H}_{\lambda}\left(\lambda \right)={\displaystyle \sum {a}_{k}\left(\delta \left(\lambda -{\lambda}_{k}\right)\right)}$ ,

Table 1. Several relaxation functions and corresponding rate distribution functions.

the coefficients ${a}_{k}$ ’s can have negative values as well for some type of relaxation function.

For example, if $i\left(t\right)=t{\left(t+1\right)}^{-1}$ is a relaxation function that initially grows to a maximum value and then starts falling as time increases, it has rate distribution function as ${H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{t{\left(t+1\right)}^{-1}\right\}=\mathrm{cos}\lambda $ , a oscillatory one. Thus in this case the distribution function i.e. ${H}_{\lambda}\left(\lambda \right)$ can take positive as well as negative values (8).

One interesting observation is for a relaxation function $i\left(t\right)={t}^{-1}$ the relaxation function is ${H}_{\lambda}\left(\lambda \right)=1$ -a “uniform distribution”, for $\lambda \ge 0$ . All these are listed in Table 1.

The inverse Laplace transformation is usually carried out by contour integration. But the very modern technique of Berberan-Santos [33] [34] method is the analytical Laplace inversion without the usual contour integration. We describe this now briefly.

Our aim is to evaluate Laplace inverse ${H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}$ which is given as Laplace inversion (13) integral expression i.e.

${H}_{\lambda}\left(\lambda \right)=\frac{1}{2\text{\pi}i}{\displaystyle \underset{x-i\infty}{\overset{x+i\infty}{\int}}\left(i\left(t\right)\right){\text{e}}^{t\lambda}dt}$ (14)

Here we describe Berberan-Santos method formulas for evaluation of the Laplace inversion without going for contour integration. First is change of variable i.e. from “real time variable” to “complex time variable” as $t=x+iy$ ; with $i=\sqrt{-1}$ . Here the real part i.e. $x$ is constant as a vertical line calls it $x={x}_{0}$ a constant. The formulas are following [33] [34]

${H}_{\lambda}(\lambda )=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle \underset{0}{\overset{\infty}{\int}}\left(\mathrm{Re}\left\{i({x}_{0}+iy)\right\}\mathrm{cos}(\lambda y)-\mathrm{Im}\left\{i({x}_{0}+iy)\right\}\mathrm{sin}(\lambda y)\right)}dy$ (15)

Consider a very simple case of decay function $i\left(t\right)={\left(t-a\right)}^{-1}$ and convert to complex time by putting $t={x}_{0}+iy$ as $i\left({x}_{0}+iy\right)={\left(\left({x}_{0}-a\right)+iy\right)}^{-1}$ [33] [34] . We know from standard Laplace pair that is ${\mathcal{L}}^{-1}{\left(s\pm a\right)}^{-1}={\text{e}}^{\mp at}$ . Thus, for $i\left(t\right)={\left(t-a\right)}^{-1}$ we should get via inverse Laplace the rate distribution functions as ${H}_{\lambda}\left(\lambda \right)={\text{e}}^{a\lambda}$ . The application of the Berberan-Santros formula [17] [35] with ${x}_{0}>a$ yields the following steps

$\begin{array}{l}{H}_{\lambda}(\lambda )=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle \underset{0}{\overset{\infty}{\int}}\left(\mathrm{Re}\left\{i({x}_{0}+iy)\right\}\mathrm{cos}(\lambda y)-\mathrm{Im}\left\{i({x}_{0}+iy)\right\}\mathrm{sin}(\lambda y)\right)}dy\\ \begin{array}{ccc}& & \end{array}=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{({x}_{0}-a)cos(\lambda y)dy}{{({x}_{0}-a)}^{2}+{y}^{2}}}+\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{y\mathrm{sin}(\lambda y)dy}{{({x}_{0}-a)}^{2}+{y}^{2}}}\\ \begin{array}{ccc}& & \end{array}=\frac{{e}^{{x}_{0}\lambda}({x}_{0}-a)}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{cos(\lambda y)dy}{{({x}_{0}-a)}^{2}+{y}^{2}}}+\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{y\mathrm{sin}(\lambda y)dy}{{({x}_{0}-a)}^{2}+{\omega}^{2}}}\\ \begin{array}{ccc}& & \end{array}=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{({x}_{0}-a)cos(\lambda y)+y\mathrm{sin}(\lambda y)}{{({x}_{0}-a)}^{2}+{y}^{2}}}dy\end{array}$ (16)

Here we say that ${\text{e}}^{a\lambda}$ has integral representation as

${e}^{a\lambda}=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{({x}_{0}-a)cos(\lambda y)+y\mathrm{sin}(\lambda y)}{{({x}_{0}-a)}^{2}+{y}^{2}}}dy$ .

Particularly for $a=-1$ , we have $i\left(t\right)={\left(t+1\right)}^{-1}$ . The condition ${x}_{0}>-1$ enables us to choose ${x}_{0}=0$ we get following integral representation for ${\text{e}}^{-\lambda}$ [33] [34] which is also rate distribution function ${H}_{\lambda}\left(\lambda \right)$ is following

${H}_{\lambda}(\lambda )={e}^{-\lambda}=\frac{1}{\pi}{\displaystyle {\int}_{0}^{\infty}\frac{cos(\lambda y)+y\mathrm{sin}(\lambda y)}{{1}^{2}+{y}^{2}}}dy$ (17)

8. Derivation of Rate Distribution Function for Curie-Von Schweidler Law

For the Curie-von Schweidler relaxation of type function i.e. $i\left(t\right)\sim {t}^{-n}$ then rate distribution function is ${H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{{t}^{-n}\right\}$ . With the known Laplace pair

i.e. ${\mathcal{L}}^{-1}\left\{{s}^{-\left(\alpha +1\right)}\right\}=\frac{1}{\alpha !}{t}^{\alpha}$ , we can write the following steps

$\begin{array}{c}{H}_{\lambda}\left(\lambda \right)={\mathcal{L}}^{-1}\left\{{t}^{-n}\right\}\\ =\frac{1}{\left(n-1\right)!}{\lambda}^{n-1}=\frac{1}{m!}{\lambda}^{m};\text{\hspace{0.17em}}\text{\hspace{0.17em}}m=n-1;\text{\hspace{0.17em}}\text{\hspace{0.17em}}\Gamma \left(\alpha \right)=\left(\alpha -1\right)!,\text{\hspace{0.17em}}\alpha \in \mathbb{R}\\ =\frac{1}{\Gamma \left(n\right)}{\lambda}^{n-1}=\frac{1}{\Gamma \left(m+1\right)}{\lambda}^{m};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda >0\end{array}$ (18)

Therefore above discussion suggests that for a power law type relaxation, i.e. Curie-von Schweidler law i.e. $i\left(t\right)\propto {t}^{-n};\text{\hspace{0.17em}}0<n<1$ , the relaxation rates l’s are having a power law distribution of type i.e. ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{m}$ , $m=n-1$ , $-1<m<0$ , $\lambda >0$ . This is Zipf’s power law with $m<0$ .

For dielectric relaxation as observed that $0<n<1$ in Curie-von Schweidler relaxation $i\left(t\right)\sim {t}^{-n}$ , the rate relaxation distribution function ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{m}$ has exponent in power in the range i.e. $-1<m<0$ . Considering graph of ${H}_{\lambda}\left(\lambda \right)={\lambda}^{m};m<0$ as histogram, we infer that for Curie-von Schweidler relaxation function i.e. $i\left(t\right)\sim {t}^{-n};0<n<1$ there are very large number of relaxations with small $\lambda $ i.e. large number of slower decay takes place, compared to fewer faster decay rates-and the histogram ${H}_{\lambda}\left(\lambda \right)={\lambda}^{m};m<0$ is highly right skewed with long tail.

From the above discussion (18) and using our Laplace integral (8) i.e. $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-t\lambda}d\lambda}$ we write for Curie-von Schweidler relaxation function the following

${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}\left({\lambda}^{n-1}\right){\text{e}}^{-t\lambda}d\lambda}$ (19)

The above Expression (19) is integral representation of the ${t}^{-n}$ shows weighted averaging of infinite Debye relaxations i.e. ${\text{e}}^{-\lambda t}$ with weight ${\lambda}^{n-1}$ applied for all $\lambda $ from zero to infinity.

Using Berberan-Santos method [33] [34] we get the Laplace inversion of $i\left(t\right)={t}^{-n}$ . In reality of decay functions, we can take ${x}_{0}=0$ ; in complex time variable i.e. $t={x}_{0}+iy$ as the decay function is not expected to have singularity at time $t>0$ . Choosing ${x}_{0}=0$ in $t={x}_{0}+iy$ we have $i\left(iy\right)={\left(iy\right)}^{-n}$ , that is

$i\left(iy\right)={\left(iy\right)}^{-n}={y}^{-n}{\left(\mathrm{cos}\left(\frac{n\text{\pi}}{2}\right)+i\mathrm{sin}\left(\frac{n\text{\pi}}{2}\right)\right)}^{-1}$ where we used

${i}^{n}={\text{e}}^{\left(in\text{\pi}/2\right)}=\mathrm{cos}\left(\frac{n\text{\pi}}{2}\right)+i\mathrm{sin}\left(\frac{n\text{\pi}}{2}\right)$ . The real part of the complex function is $\mathrm{Re}\left\{i\left(iy\right)\right\}={y}^{-n}\mathrm{cos}\left(\frac{n\text{\pi}}{2}\right)$ . Now using the Berberan-Snatos formula [33] [34] , we

get following steps

$\begin{array}{l}{H}_{\lambda}(\lambda )=\frac{{e}^{{x}_{0}\lambda}}{\pi}{\displaystyle \underset{0}{\overset{\infty}{\int}}\left(\mathrm{Re}\left\{i({x}_{0}+iy)\right\}\mathrm{cos}(\lambda y)-\mathrm{Im}\left\{i({x}_{0}+iy)\right\}\mathrm{sin}(\lambda y)\right)}dy\\ \begin{array}{ccc}& & \end{array}=\frac{1}{\pi}{\displaystyle \underset{0}{\overset{\infty}{\int}}\left({y}^{-n}\mathrm{cos}\left({\scriptscriptstyle \frac{n\pi}{2}}\right)\mathrm{cos}(\lambda y)+{y}^{-n}\mathrm{sin}\left({\scriptscriptstyle \frac{n\pi}{2}}\right)\mathrm{sin}(\lambda y)\right)}dy\\ \begin{array}{ccc}& & \end{array}=\frac{1}{\pi}{\displaystyle {\int}_{0}^{\infty}{y}^{-n}}\mathrm{cos}\left(\lambda y-{\scriptscriptstyle \frac{n\pi}{2}}\right)dy\end{array}$ (20)

From Laplace transform tables [28] (and Table 1) we have ${H}_{\lambda}\left(\lambda \right)=\frac{1}{\Gamma \left(n\right)}{\lambda}^{n-1}$

i.e. from inverse Laplace transformed of $i\left(t\right)={t}^{-n}$ , therefore we write following representation

${H}_{\lambda}(\lambda )=\frac{1}{\Gamma (n)}{\lambda}^{n-1}=\frac{1}{\pi}{\displaystyle {\int}_{0}^{\infty}{y}^{-n}}\mathrm{cos}\left(\lambda y-{\scriptscriptstyle \frac{n\pi}{2}}\right)dy$ (21)

From (21) we write the following

${\lambda}^{n-1}=\frac{\Gamma (n)}{\pi}{\displaystyle {\int}_{0}^{\infty}{y}^{-n}}\mathrm{cos}\left(\lambda y-{\scriptscriptstyle \frac{n\pi}{2}}\right)dy$ (22)

Now by changing variable $\lambda $ to $t$ , $n-1$ to $-n$ , and rearrange above

expression (22) to get ${t}^{-n}=\frac{2\Gamma \left(1-n\right)}{\text{\pi}}\mathrm{cos}\left(\frac{\left(1-n\right)\text{\pi}}{2}\right){\displaystyle {\int}_{0}^{\infty}{u}^{\left(n-1\right)}}\left(\mathrm{cos}\left(tu\right)\right)\left(du\right)$ .

Considering now $u$ as $\lambda $ we write another integral representation of ${t}^{-n}$ as follows

${t}^{-n}=\frac{\Gamma (1-n)}{\pi}{\displaystyle {\int}_{0}^{\infty}{\lambda}^{n-1}}\mathrm{cos}\left(\lambda t-{\scriptscriptstyle \frac{(1-n)\pi}{2}}\right)d\lambda $ (23)

This means that if we chose basic relaxation function as $\mathrm{cos}\left(\lambda t\right)$ , then Curie-von Schweidler relaxation $i\left(t\right)={t}^{-n}$ is weighted sum of all $\mathrm{cos}\left(\lambda t\right)$ ’s with weights ${\lambda}^{n-1}$ , as $\lambda $ is varied from zero to infinity.

9. Zipf’s Distribution for Relaxation Time Constants for Curie-Von Schweidler Law―A Contradiction

Now converting to $\tau ={\lambda}^{-1}$ , we assume the Distribution of time-constants call it ${H}_{\tau}\left(\tau \right)\sim {\tau}^{-m}$ , is the Zipf’s power law distribution. However direct taking of reciprocal of obtained inversion of ${H}_{\lambda}\left(\lambda \right)$ i.e. the rate distribution function got via Laplace inversion of $i\left(t\right)$ is not possible. This we demonstrate in this section.

As we have formulated Laplace integral (8) i.e. $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-t\lambda}d\lambda}$ , just

by replacing $\lambda ={\tau}^{-1}$ , $d\lambda =-{\left(\tau \right)}^{-2}d\tau $ we will get

$i\left(t\right)={\displaystyle {\int}_{\infty}^{0}\left({H}_{\lambda}\left(\lambda \right)\right){\left(-\tau \right)}^{-2}{\text{e}}^{-t/\tau}d\tau}$ that is $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\left(\tau \right)}^{-2}{\text{e}}^{-t/\tau}d\tau}$ . This is not Laplace integral. Now we do the following steps, for $i\left(t\right)={t}^{-n}$ and obtained

${H}_{\lambda}\left(\lambda \right)=\frac{1}{\Gamma \left(n\right)}{\lambda}^{n-1}$

$\begin{array}{c}{t}^{-n}={\displaystyle {\int}_{0}^{\infty}\left(\frac{1}{\Gamma \left(n\right)}{\lambda}^{n-1}\right){\left(\tau \right)}^{-2}{\text{e}}^{-t/\tau}d\tau};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda ={\tau}^{-1}\\ =\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{\tau}^{1-n}{\tau}^{-2}{\text{e}}^{-t/\tau}d\tau}\\ =\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{\tau}^{-\left(n+1\right)}{\text{e}}^{-t/\tau}d\tau}\end{array}$ (24)

We write the two representations of ${t}^{-n}$ as following integrals

${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}\left({\lambda}^{\left(n-1\right)}\right){\text{e}}^{-t\lambda}d\lambda},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}\left({\tau}^{-\left(n+1\right)}\right){\text{e}}^{-t/\tau}}d\tau $ (25)

Thus we have ${H}_{\tau}\left(\tau \right)\sim {\tau}^{-\left(n+1\right)}$ , as we have ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1}$ . Now we verify the above obtained result in the subsequent discussion.

By the logic that we had constructed $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-t\lambda}d\lambda}$ which is Laplace integral (8); we will similarly get the integral $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\tau}\left(\tau \right)\right){\text{e}}^{-t/\tau}d\tau}$ which is not a direct Laplace Transform formula. Following steps will convert this expression into the Laplace Transform formula, and from there we will extract ${H}_{\tau}\left(\tau \right)$ :

$\begin{array}{c}i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\tau}\left(\tau \right)\right){\text{e}}^{-t/\tau}d\tau};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\tau ={\lambda}^{-1},\text{\hspace{0.17em}}d\tau =-{\lambda}^{-2}d\lambda \\ ={\displaystyle {\int}_{\infty}^{0}\left(-{\lambda}^{-2}\left({H}_{\tau}\left(\tau \right)\right)\right)}\left({\text{e}}^{-t\lambda}\right)d\lambda ;\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}F={\lambda}^{-2}\left({H}_{\tau}\left(\tau \right)\right)\\ ={\displaystyle {\int}_{0}^{\infty}\left(F\right)}\left({\text{e}}^{-t\lambda}\right)d\lambda \end{array}$ (26)

Proceeding further we obtain following result

$\begin{array}{l}F={\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}=\frac{1}{\Gamma \left(n\right)}{\lambda}^{\left(n-1\right)}\\ {\lambda}^{-2}\left({H}_{\tau}\left(\tau \right)\right)=\frac{1}{\Gamma \left(n\right)}{\lambda}^{\left(n-1\right)}\\ {H}_{\tau}(\tau )=\frac{{\lambda}^{n-1}{\lambda}^{2}}{\Gamma (n)};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda ={\tau}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.05em}}=\frac{{\tau}^{-\left(n+1\right)}}{\Gamma \left(n\right)}\end{array}$ (27)

Now we take different approach to verify the above (26) (27) obtained expression for ${H}_{\tau}\left(\tau \right)$ . Let us have set of relaxation functions with various time constants $\tau $ ranging from 0 to infinity that is $\left\{{\text{e}}^{-t/{\tau}_{1}},{\text{e}}^{-t/{\tau}_{2}},{\text{e}}^{-t/{\tau}_{3}},\cdots \right\}$ , comprising of infinite number of functions, in continuum in $\tau $ . The relaxation function varies from very-very quick decay (when $\tau \approx 0$ ) to very-very slow decay curve (when $\tau \approx \infty $ ). We construct a weighted decay function as ${\tau}^{-m}{\text{e}}^{-t/\tau}$ . This shows that we are multiplying by weight ${\tau}^{-m};m>0$ the decay function ${\text{e}}^{-t/\tau}$ . We are assuming Zipf’s type distribution of $\tau $ , in form of ${H}_{\tau}\left(\tau \right)\sim {\tau}^{-m}$ , meaning the lowest time constant i.e. fastest decay occurs more frequent than slow decay i.e. large time constant. The time constant parameter $\tau $ let vary from 0 to infinity and construct the following integral $I$ , i.e.

$I={\displaystyle {\int}_{0}^{\infty}{\tau}^{-m}}{\text{e}}^{-t/\tau}d\tau $ (28)

The integral of (28) i.e. $I={\displaystyle {\int}_{0}^{\infty}{\tau}^{-m}}{\text{e}}^{-t/\tau}d\tau $ gives notion of weighted average of

infinite relaxation functions. We do the substitution i.e. $\frac{t}{\tau}=y$ , i.e. $\tau =\frac{t}{y}$ and

$d\tau ={\left(-y\right)}^{-2}t\left(dy\right)$ in the above integral to get following steps

$\begin{array}{c}I={\displaystyle {\int}_{0}^{\infty}{\tau}^{-m}}{\text{e}}^{-t/\tau}d\tau \\ ={\displaystyle {\int}_{\infty}^{0}{\left(\frac{t}{y}\right)}^{-m}}{\text{e}}^{-y}{\left(-y\right)}^{-2}t\left(dy\right)\\ ={t}^{-m+1}{\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-y}}{y}^{m-2}dy\end{array}$ (29)

By using the definition of the Gamma function [29] [32] in integral form i.e. $\Gamma \left(\alpha \right)={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-y}}{y}^{\alpha -1}dy$ , we write the above integral as following

$\begin{array}{l}I={t}^{-\left(m-1\right)}{\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-y}}{y}^{m-2}dy={t}^{-\left(m-1\right)}\left(\Gamma \left(m-1\right)\right)\\ {\displaystyle {\int}_{0}^{\infty}{\tau}^{-m}}{\text{e}}^{-t/\tau}d\tau =\frac{\Gamma \left(m-1\right)}{{t}^{m-1}}\end{array}$ (30)

Putting $m-1=n$ in above we get integral representation of the power law ${t}^{-n}$ and we represent this by time constant distribution function ${H}_{\tau}\left(\tau \right)$ in following expressions

$\begin{array}{l}{t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{\tau}^{-\left(n+1\right)}{\text{e}}^{-t/\tau}d\tau}\\ {t}^{-n}={\displaystyle {\int}_{0}^{\infty}\left({H}_{\tau}\left(\tau \right)\right){\text{e}}^{-t/\tau}d\tau}\end{array}$ (31)

Earlier in (19) we have obtained ${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}\left({\lambda}^{n-1}\right){\text{e}}^{-t\lambda}d\lambda}$ ; where we called

rate distribution function as ${H}_{\lambda}\left(\lambda \right)=\frac{1}{\Gamma \left(n\right)}{\lambda}^{n-1}$ . Now from above (31) weighted average logic we get ${H}_{\tau}\left(\tau \right)=\frac{1}{\Gamma \left(n\right)}{\tau}^{-\left(n+1\right)}$ . We note that these two are not

reciprocal of each other.

Therefore we can conclude that Curie-von Schweidler law ( $\sim {t}^{-n}$ ) relates to weighted averaging of several classical Debye relaxations (of type ${\text{e}}^{-t/\tau}$ ) over several time constants from zero to infinity, that is having Zipf’s power-law with time constant distribution as ${H}_{\tau}\left(\tau \right)\sim {\tau}^{-\left(n+1\right)}$ . What does it say for $0<n<1$ , that ${H}_{\tau}\left(\tau \right)\sim {\tau}^{-\left(n+1\right)}$ is also a right-skewed distribution, where the lower time constants (faster decay) appear more than larger time constant (slower decay). This is contradiction to what we inferred for ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1}$ .

10. Demonstration of Contradiction of Obtained Zipf’s Power Law Distribution for Relaxation Rates and Time Constant Distributions

This contradiction we demonstrate now. For $i\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ we got ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda -{\lambda}_{0}\right)$ ; we expect ${H}_{\tau}\left(\tau \right)=\delta \left(\tau -{\tau}_{0}\right)$ for $i\left(t\right)={\text{e}}^{-t/{\tau}_{0}}$ ; let us see what happens in following steps

$\begin{array}{l}i\left(t\right)={\text{e}}^{-{\lambda}_{0}t};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}F={\mathcal{L}}^{-1}\left\{i\left(t\right)\right\}\\ F={\mathcal{L}}^{-1}\left\{{\text{e}}^{-{\lambda}_{0}t}\right\}=\delta \left(\lambda -{\lambda}_{0}\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}F={\lambda}^{-2}\left({H}_{\tau}\left(\tau \right)\right)\\ {\lambda}^{-2}\left({H}_{\tau}\left(\tau \right)\right)=\delta \left(\lambda -{\lambda}_{0}\right)\\ {H}_{\tau}\left(\tau \right)={\lambda}^{2}\left(\delta \left(\lambda -{\lambda}_{0}\right)\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda ={\tau}^{-1}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}=\frac{\delta \left(\tau -{\tau}_{0}\right)}{{\tau}^{2}}\ne \delta \left(\tau -{\tau}_{0}\right)\end{array}$ (32)

Though mathematically we can get integral representation for any relaxation function as $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\tau}\left(\tau \right)\right){\text{e}}^{-t/\tau}}d\tau $ but physically it will be contradictory to Laplace integral i.e. $i\left(t\right)={\displaystyle {\int}_{0}^{\infty}\left({H}_{\lambda}\left(\lambda \right)\right){\text{e}}^{-\lambda t}}d\lambda $ . Hence we will deal with the relaxation rate distribution function that we extracted as ${H}_{\lambda}\left(\lambda \right)$ from $i\left(t\right)$ via our devised method of Laplace inversion.

11. Experimental Validation of Range of Relaxation Exponent in Curie-Von Schweidler Law

The Curie-von Schweidler empirical law of power law relaxation, i.e. $i\left(t\right)\propto {t}^{-n}$ states that $0<n<1$ . This is validated via experiments on dielectric relaxations. A 100V step input applied to a completely discharged capacitor of 0.47 μF having metalized paper dielectric, and the current decay is recorded with time. The graphs of log-log plot i.e. $\mathrm{log}\left(i\left(t\right)\right)$ vs. $\mathrm{log}\left(t\right)$ show a straight line of average slope −0.86 [5] - [10] . This experiment indicates a Curie-von Schweidler law, with $i\left(t\right)\propto {t}^{-n}$ , having $n=0.86$ . This gives ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{-0.14}$ . The exponent n is in the range of $0.85<n<1$ in several di-electric relaxation experiments [5] - [10] . The experiments with super-capacitors [13] [14] , show range as $0.5<n<1$ . A very low value of exponent n is found in relaxation of Laponite studies averagely $n=0.09$ [35] . Thus in case of Laponite studies we have relaxation rate distribution function as ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{-0.91}$ . In this Laponite study [35] though the exponent n was obtained on “self-discharge” curves with various charging time history-showing memory effect, the expression obtained for self-discharge decay of voltage assumes fractional capacity-that in turn assumes Curie-von Schweidler law as current relaxation function.

12. Time Variant Relaxation Rate for Non-Debye Relaxation & Curie-Von Schweidler Law

Any decay function $i\left(t\right)$ is written as a general formulation in following way

$i\left(t\right)=\mathrm{exp}\left(-{\displaystyle {\int}_{0}^{t}\left(\lambda \left(\xi \right)\right)d\xi}\right)$ (33)

where $\lambda \left(\xi \right)$ is the time ( $\xi $ ) dependent rate coefficient. When the relaxation is pure exponential, one has $\lambda \left(\xi \right)$ as constant say ${\lambda}_{0}$ described as $\lambda \left(t\right)={\lambda}_{0}$ , expressed in following steps

$\begin{array}{c}i\left(t\right)=\mathrm{exp}\left(-{\displaystyle {\int}_{0}^{t}\left(\lambda \left(\xi \right)\right)d\xi}\right)=\mathrm{exp}\left(-{\displaystyle {\int}_{0}^{t}\left({\lambda}_{0}\right)d\xi}\right)\\ =\mathrm{exp}\left(-{\lambda}_{0}{\xi |}_{0}^{t}\right)={\text{e}}^{-{\lambda}_{0}t}\end{array}$ (34)

Thus we get a Debye relaxation for a system having constant rate of relaxation. To extract $\lambda \left(t\right)$ that is time dependent rate coefficient we have to follow the following steps by taking logarithm of (33) and then differentiating both sides

$\mathrm{ln}\left(i\left(t\right)\right)=-{\displaystyle {\int}_{0}^{t}\left(\lambda \left(\xi \right)\right)d\xi},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda \left(t\right)=-\frac{d}{dt}\mathrm{ln}\left(i\left(t\right)\right)$ (35)

We use the above rule for Mittag-Leffler relaxation function i.e. $i\left(t\right)={E}_{\alpha}\left(-t\right)$ ,

this is defined as [29] [30] [31] ${E}_{\alpha}\left(-t\right)={\displaystyle {\sum}_{k=0}^{\infty}\frac{{\left(-1\right)}^{k}}{\Gamma \left(\alpha k+\beta \right)}{t}^{k}}$ . For $i\left(t\right)={E}_{\alpha}(-t)$

the $\lambda \left(t\right)$ is extracted as in following steps

$\begin{array}{c}\lambda \left(t\right)=-\frac{d}{dt}ln\left({E}_{\alpha}\left(-t\right)\right)=-\frac{1}{{E}_{\alpha}\left(-t\right)}\frac{d{E}_{\alpha}\left(-t\right)}{dt}\\ =-\frac{1}{{E}_{\alpha}\left(-t\right)}\frac{d}{dt}\left({\displaystyle \underset{k=0}{\overset{\infty}{\sum}}\frac{{\left(-1\right)}^{k}{t}^{k}}{\Gamma \left(\alpha k+1\right)}}\right)\\ =-\frac{1}{{E}_{\alpha}\left(-t\right)}\left({\displaystyle \underset{k=0}{\overset{\infty}{\sum}}\frac{\left(k+1\right){\left(-t\right)}^{k}}{\Gamma \left(1+\alpha +\alpha k\right)}}\right)\end{array}$ (36)

For Curie-von Schweilder relaxation law $i\left(t\right)\sim {t}^{-n}$ , we have time dependent rate relaxation rate as

$\begin{array}{l}\lambda \left(t\right)=-\frac{d}{dt}ln\left({t}^{-n}\right)=-\frac{d}{dt}\left(-n\mathrm{ln}t\right)=n\left({t}^{-1}\right)\\ \tau \left(t\right)=\frac{1}{\lambda \left(t\right)}=\frac{t}{n}\end{array}$ (37)

Therefore we have two observations that Curie-von Schweidler relaxation $i\left(t\right)\sim {t}^{-n}$ has time constant distributed as Zipfian power law ${H}_{\lambda}\sim {\lambda}^{n-1};0<n<1$ , while the relaxation rate constant is variable in time as a function $\lambda \left(t\right)=n/t$ . Thus implying the relaxation starts with very-very fast relaxation at a very-very high rate ( $\lambda $ ) and as the time goes the rate constant decreases indicating slow rate of current decay. This is a case of “equivalent” single body relaxation where the rate is varying with time; whereas the multi-body relaxation gives simultaneous relaxations with rates distributed as Zipf’s power law.

13. Scale Dependence Relaxation Rates Give Capacitors Charging Current as per Curie-Von Schweidler Law

Let a uncharged capacitor $C$ be connected to a voltage source $v\left(t\right)$ Volts, at time $t=0$ ; obviously this capacitor will get charged to the battery voltage. Let this capacitor is uncharged at $t<0$ , thus there is no charge held by it, therefore the voltage across the capacitor is zero at $t<0$ , and the circuit current is $i\left(t\right)=0;t<0$ . The voltage balance equation assuming $R$ be the total resistance of the circuit (including internal resistance of Capacitor) at $t>0$ is the following

$\left(\frac{1}{C}{\displaystyle \underset{0}{\overset{t}{\int}}i\left(x\right)d}x\right)+Ri\left(t\right)=v\left(t\right)$ (38)

where $i\left(t\right)$ is the charging current flowing into the capacitor. The above integral Equation (38) is differentiated and is put as following, for $t>0$

$\begin{array}{l}\frac{di\left(t\right)}{dt}+{\lambda}_{0}i\left(t\right)=\left(\frac{1}{R}\right)\frac{dv\left(t\right)}{dt},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\lambda}_{0}={\left(RC\right)}^{-1},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\tau}_{0}=RC\\ {i}^{\left(1\right)}\left(t\right)+{\lambda}_{0}i\left(t\right)=f\left(t\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}f\left(t\right)=\frac{1}{R}\left({v}^{\left(1\right)}\left(t\right)\right)\end{array}$ (39)

The RHS of above Expression (39) first order system indicates “forcing

function” which is $f\left(t\right)$ . The forcing function is $f\left(t\right)=\frac{1}{R}{v}^{\left(1\right)}\left(t\right)$ in this case.

If we take $v\left(t\right)={V}_{BB}$ i.e. a constant, then considering $u\left(t\right)=1;t\ge 0$ and $u\left(t\right)=0;t<0$ , i.e. “unit-step function”, we have for RHS of the above Equation (39) as following

$\begin{array}{l}\frac{dv\left(t\right)}{dt}={v}^{\left(1\right)}\left(t\right)\\ {v}^{\left(1\right)}\left(t\right)={V}_{BB}\frac{d\left(u\left(t\right)\right)}{dt}={V}_{BB}\left(\delta \left(t\right)\right)\end{array}$ (40)

Thus substituting the above (40) into (39) we get following equation for $v\left(t\right)={V}_{BB}$ for $t\ge 0$ the following

$\begin{array}{l}\frac{di\left(t\right)}{dt}+{\lambda}_{0}i\left(t\right)=\frac{{V}_{BB}}{R}\delta \left(t\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\lambda}_{0}={\left(RC\right)}^{-1},\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\tau}_{0}=RC\\ \frac{di\left(t\right)}{dt}+{\lambda}_{0}i\left(t\right)={I}_{0}\delta \left(t\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{I}_{0}=\frac{{V}_{BB}}{R}\end{array}$ (41)

We see that forcing function of above first order Equation (41) is by delta function $f\left(t\right)={I}_{0}\delta \left(t\right)$ . The solution to the above equation gives Debye relaxation function i.e.

$i\left(t\right)={I}_{0}{\text{e}}^{-{\lambda}_{0}t}$ (42)

This solution $i\left(t\right)\sim {\text{e}}^{-{\lambda}_{0}t}$ is the “impulse response” of the circuit equation. The relaxation current of the above system (41) follows Debye’s relaxation, with one relaxation rate ${\lambda}_{0}$ (also termed as Debye law). The rate distribution function is ${H}_{\lambda}\left(\lambda \right)=\delta \left(\lambda -{\lambda}_{0}\right)$ ; that we discussed in previous sections. With this $i\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ as Green’s function call it $g\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ i.e. solution of differential equation with “unit” impulse excitation ( ${I}_{0}=1$ ) or say Homogeneous solution, i.e.

$\frac{di\left(t\right)}{dt}+{\lambda}_{0}i\left(t\right)=\delta \left(t\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}i\left(t\right)=g\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ (43)

Now we find if the input is step function at time $t=0$ , call it ${I}_{0}u\left(t\right)$ , where $u\left(t\right)=1$ for $t\ge 0$ and $u\left(t\right)=0$ for $t<0$ ; then we get relaxation function for current as convolution integral, i.e. depicted as in following steps

$\begin{array}{l}\frac{di\left(t\right)}{dt}+{\lambda}_{0}i\left(t\right)={I}_{0}u\left(t\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}g\left(t\right)={\text{e}}^{-{\lambda}_{0}t}\\ i\left(t\right)=\left({I}_{0}u\left(t\right)\right)\ast \left(g\left(t\right)\right)={\displaystyle {\int}_{-\infty}^{t}\left({I}_{0}u\left(x\right)\right)}\left({\text{e}}^{-{\lambda}_{0}\left(t-x\right)}\right)dx={I}_{0}{\displaystyle {\int}_{0}^{t}{\text{e}}^{-{\lambda}_{0}\left(t-x\right)}dx}\\ i\left(t\right)=\frac{{I}_{0}}{{\lambda}_{0}}\left(1-{\text{e}}^{-{\lambda}_{0}t}\right)\end{array}$ (44)

We saw in earlier sections the relaxation rates ( $\lambda $ ) distribution, for a Curie-von Schweidler relaxation law, i.e. $i\left(t\right)\sim {t}^{-n}$ is ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1};0<n<1$ ; for relaxations in dielectrics. This is histogram of rates. It says that the relaxation of current is with several relaxation rates, which are distributed as discussed in Zipf’s law fashion with right-skewed-histogram. Thus if we represent the equivalent relaxation rate say ${\lambda}_{eq}\sim {\lambda}^{1/q};0<q<1$ with $\lambda $ as scale of relaxation where the scale $\lambda $ varies from zero to infinity; we will not be incorrect in assuming this. That is as we slide from a low scale $\lambda $ to high scale $\lambda $ the equivalent relaxation rate ${\lambda}_{eq}$ will be different at different scales of relaxation. If the index parameter i.e. $q=1$ then we have single rate constant system given by ${\lambda}_{eq}=\lambda $ always at all scales of relaxation i.e. $\lambda ={\lambda}_{0}={\tau}_{0}^{-1}={\left(RC\right)}^{-1}$ , and with solution as $i\left(t\right)={\text{e}}^{-{\lambda}_{0}t}$ , i.e. Debye relaxation function.

We thus modify the capacitor discharge current equation, with ${\lambda}_{eq}=\lambda ={\lambda}_{0}$ i.e. with one relaxation rate at any scale of relaxation (λ), i.e. ${i}^{\left(1\right)}\left(t\right)+{\lambda}_{0}i\left(t\right)=\delta \left(t\right)$ to following i.e. variable ${\lambda}_{eq}={\lambda}^{1/q}$ at any scale of relaxation rate(λ)

$\frac{di\left(t\right)}{dt}+\left({\lambda}_{eq}\right)i\left(t\right)=\delta \left(t\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\frac{di\left(t\right)}{dt}+{\left(\lambda \right)}^{1/q}i\left(t\right)=\delta \left(t\right)$ (45)

The initial condition is given as $i\left(t\right)=0$ for $t<0$ . The above equation is having a free “scale” parameter $\lambda $ varying from zero to infinity. The solution of the above is $i\left(t\right)={\text{e}}^{-{\lambda}_{eq}t}$ . We call this $i\left(t\right)=\mathrm{exp}\left(-{\lambda}_{eq}t\right)=\mathrm{exp}\left(-{\lambda}^{1/q}t\right)$ as “impulse response function” at a particular scale $\lambda $ , i.e. we call it $h\left(\lambda ,t\right);\lambda \in \left(0,\infty \right)$

$i\left(t\right)=h\left(\lambda ,t\right)={\text{e}}^{-\left({\lambda}^{1/q}t\right)};\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<q<1$ (46)

The above Expression (46) actually is valid for all scale $\lambda $ varying from zero to infinity. Thus on integrating this “impulse response function” $h\left(\lambda ,t\right)$ on the free variable ( $\lambda $ ) from 0 to $\infty $ , we get the function of time and that is called “impulse response” or the Green’s function $g\left(t\right)$ as depicted in following derivation

$g\left(t\right)={\displaystyle \underset{0}{\overset{\infty}{\int}}h\left(\lambda ,t\right)d\lambda}={\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-\left({\lambda}^{1/q}t\right)}d\lambda}=\frac{\Gamma \left(1+q\right)}{{t}^{q}}$ (47)

To get above Expression (47), we substitute in $g\left(t\right)={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-\left({\lambda}^{1/q}t\right)}d\lambda}$ , ${\lambda}^{1/q}t=x$

that makes following changes:

$\begin{array}{l}\lambda ={\left(\frac{x}{t}\right)}^{q};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\left(\frac{x}{t}\right)={\lambda}^{1/q}\\ d\lambda =q{x}^{q-1}{\left(\frac{1}{t}\right)}^{q}dx=\left(\frac{q}{t}\right){\left(\frac{x}{t}\right)}^{q-1}dx=\left(\frac{q}{t}\right){\left({\lambda}^{1/q}\right)}^{q-1}dx\\ d\lambda ={\lambda}^{1-\left(1/q\right)}\left(\frac{q}{t}\right)dx\end{array}$ (48)

Then by using definition of Gamma function i.e. $\Gamma \left(\alpha \right)={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-y}{y}^{\alpha -1}dy}$ , and its property $\alpha \left(\Gamma \left(\alpha \right)\right)=\Gamma \left(1+\alpha \right)$ the following steps are followed to get the desired

expression i.e. $g\left(t\right)=\frac{\Gamma \left(1+q\right)}{{t}^{q}}$

$\begin{array}{c}g\left(t\right)={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-\left({\lambda}^{1/q}t\right)}d\lambda}={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-x}{\lambda}^{1-\left(1/q\right)}\left(\frac{q}{t}\right)dx}\\ ={\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-x}\left(\frac{q}{t}\right)\lambda {\lambda}^{-\left(1/q\right)}dx},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\lambda ={\left(\frac{x}{t}\right)}^{q}\\ ={\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-x}}\left(\frac{q}{t}\right){\left(\frac{x}{t}\right)}^{q}{\left(\frac{x}{t}\right)}^{-1}dx=\left(\frac{q}{t}\right){\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-x}}{\left(\frac{x}{t}\right)}^{q}{\left(\frac{x}{t}\right)}^{-1}dx\\ =\left(\frac{q}{t}\right){\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-x}}\frac{{x}^{q-1}}{{t}^{q-1}}dx=\left(\frac{q}{{t}^{q}}\right){\displaystyle \underset{0}{\overset{\infty}{\int}}{\text{e}}^{-x}{x}^{q-1}dx}=\frac{q\left(\Gamma \left(q\right)\right)}{{t}^{q}}=\frac{\Gamma \left(1+q\right)}{{t}^{q}}\end{array}$ (49)

By changing q to n we get integral representation of ${t}^{-n}$ as following

${t}^{-n}=\frac{1}{n\left(\Gamma \left(n\right)\right)}{\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-\left({\lambda}^{1/n}t\right)}d\lambda}$ (50)

For q = 1 case we have scale invariance λ thus $i\left(t\right)=h\left({\lambda}_{0},t\right)=\mathrm{exp}\left(-{\lambda}_{0}t\right);\text{\hspace{0.17em}}q=1$ , where ${\lambda}_{eq}={\lambda}_{0}$ at all scales. For this case $q=1$ “impulse response” or Green’s function is $g\left(t\right)=h\left({\lambda}_{0},t\right)={\text{e}}^{-{\lambda}_{0}t}$ same as “impulse response function” i.e. $h\left({\lambda}_{0},t\right)$ .

We find that for a system where the equivalent relaxation rate is ${\lambda}_{eq}={\lambda}^{1/n}$ ; similar to a distribution function that we obtained as ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1}$ gives relaxation current as $i\left(t\right)\sim {t}^{-n}$ . We write the two currents expressions $i\left(t\right)\sim {t}^{-n}$ obtained as following for $0<n<1$

${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{\lambda}^{n-1}{\text{e}}^{-t\lambda}d\lambda},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}{t}^{-n}=\frac{1}{n\left(\Gamma \left(n\right)\right)}{\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-\left({\lambda}^{1/n}t\right)}d\lambda}$ (51)

Therefore we infer that the Curie-von Schweidler relaxation current for dielectric excited by a step voltage that follows the relation $i\left(t\right)\sim {t}^{-n};0<n<1$ has distribution function ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1}$ a power law or Zipfian distribution, with scale dependent relaxation rate described as ${\lambda}_{eq}={\lambda}^{1/n}$ .

Here we mention that in this study, we have derived mathematically several integral representations of Curie-von Schweidler relaxation function $i\left(t\right)\sim {t}^{-n};0<n<1$ , (19), (23), (24), (50). Those we list as

${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{y}^{n-1}{\text{e}}^{-ty}dy}$ , ${t}^{-n}=\frac{\Gamma (1-n)}{\pi}{\displaystyle {\int}_{0}^{\infty}{y}^{n-1}}\mathrm{cos}\left(yt-{\scriptscriptstyle \frac{(1-n)\pi}{2}}\right)dy$ ,

${t}^{-n}=\frac{1}{\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{y}^{-\left(n+1\right)}}{\text{e}}^{-t/y}dy$ and ${t}^{-n}=\frac{1}{n\Gamma \left(n\right)}{\displaystyle {\int}_{0}^{\infty}{\text{e}}^{-\left({y}^{1/n}t\right)}dy}$ . These are mathemati-

cally equivalent representations of the function $i\left(t\right)\sim {t}^{-n};0<n<1$ , yet using for physical explanations one needs care. However, these are formulations for few definite integrals giving the same result as ${t}^{-n}$ .

14. Appearance of Fractional Derivative―In the System Having Zipfian Power Law Distribution in Relaxation Rates, where the Equivalent Relaxation Rate Is Scale Dependent

The delta-function for excitation as shown in above section gives homogeneous system with solution as $g\left(t\right)={t}^{-n}\left(\Gamma \left(1+n\right)\right)$ i.e. described as following

$\frac{di\left(t\right)}{dt}+{\left(\lambda \right)}^{1/n}i\left(t\right)=\delta \left(t\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1;\text{\hspace{0.17em}}\text{\hspace{0.17em}}i\left(t\right)=g\left(t\right)=\frac{\Gamma \left(1+n\right)}{{t}^{n}}$ (52)

Now let the system described above be excited by a signal proportional to $f\left(t\right)\sim {v}^{\left(1\right)}\left(t\right)$ , a derivative of voltage excitation function $v\left(t\right)$ ; so we write this as following

$\frac{di\left(t\right)}{dt}+{\left(\lambda \right)}^{1/n}i\left(t\right)={v}^{\left(1\right)}\left(t\right)$ (53)

Note that if $v\left(t\right)=u\left(t\right)$ , that is unit-step-function at time $t=0$ then ${v}^{\left(1\right)}\left(t\right)=\delta \left(t\right)$ , we recover the above homogeneous differential Equation (52). Then the response to this new excitation function ${v}^{\left(1\right)}\left(t\right)$ is convolution of

Green’s function obtained i.e. $g\left(t\right)=\frac{\Gamma \left(1+n\right)}{{t}^{n}}$ above (52), with the forcing

function $f\left(t\right)$ i.e. now $\sim {v}^{\left(1\right)}\left(t\right)$ . We write the following steps to get $i\left(t\right)$ for a forcing function $f(t)$

$\begin{array}{c}i\left(t\right)=\left(g\left(t\right)\right)\ast \left(f\left(t\right)\right)=\left(g\left(t\right)\right)\ast \left({v}^{\left(1\right)}\left(t\right)\right)\\ ={\displaystyle \underset{0}{\overset{t}{\int}}\left(g\left(t-x\right)\right)\left({v}^{\left(1\right)}\left(x\right)\right)dx};\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}g\left(t\right)=\frac{\Gamma \left(1+n\right)}{{t}^{n}}\\ =\Gamma \left(1+n\right){\displaystyle \underset{0}{\overset{t}{\int}}\frac{{v}^{\left(1\right)}\left(x\right)}{{\left(t-x\right)}^{n}}dx},\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1\end{array}$ (54)

Multiplying and dividing the above Expression (54) with $\Gamma \left(1-n\right)$ and using the definition of fractional integral [6] [34] that is

${}_{0}\mathcal{I}{}_{t}^{\alpha}\left(f\left(t\right)\right)={}_{0}D{}_{t}^{-\alpha}\left(f\left(t\right)\right)={\displaystyle \underset{0}{\overset{t}{\int}}\frac{1}{\Gamma \left(\alpha \right)}{\left(t-x\right)}^{\alpha -1}\left(f\left(x\right)\right)dx},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\alpha >0$ (55)

we get the following derivation

$\begin{array}{c}i\left(t\right)=\Gamma \left(1+n\right)\Gamma \left(1-n\right){\displaystyle \underset{0}{\overset{t}{\int}}\frac{{\left(t-x\right)}^{\left(-n\right)}}{\Gamma \left(1-n\right)}\left({v}^{\left(1\right)}\left(x\right)\right)dx}\\ =\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}\mathcal{I}{}_{t}^{\left(1-n\right)}\left[{v}^{\left(1\right)}\left(t\right)\right]\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\left(1-n\right)>0\\ =\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}D{}_{t}^{-\left(1-n\right)}\left[{v}^{\left(1\right)}\left(t\right)\right]\right)\\ =\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}D{}_{t}^{n}{}_{0}D{}_{t}^{-1}\left[{v}^{\left(1\right)}\left(t\right)\right]\right)\\ =\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}D{}_{t}^{n}\left[v\left(t\right)\right]\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}n<1\end{array}$ (56)

In above derivation (56) we have used the relation ${}_{0}\mathcal{I}{}_{t}^{\alpha}={}_{0}D{}_{t}^{-\alpha};\alpha >0$ and the composition rule i.e. ${}_{0}D{}_{t}^{\alpha +\beta}={}_{0}D{}_{t}^{\alpha}{}_{0}D{}_{t}^{\beta};\alpha >0;\beta <0$ [12] [31] . This derivation implies the appearance of fractional derivative for cases where several relaxation rates (ideally infinite of them) define a relaxation process; which are having a scale dependence behavior, i.e. ${\lambda}_{eq}={\lambda}^{1/n}$ with histogram distributed as Zipf’s power law i.e. ${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1}$ , and the relaxation is by Curie-von Schweidler law i.e. $i\left(t\right)\sim {t}^{-n},0<n<1$ . Thus we have current through a system (having a complex relaxation process with several rate distributed as power law excited by a voltage $v\left(t\right)$ as fractional derivative of it, i.e. $i\left(t\right)\propto {}_{0}D{}_{t}^{n}\left[v\left(t\right)\right]$ .

Let this system with scale dependent relaxation rates with $i\left(t\right)=0$ for $t<0$ i.e.

$\frac{di\left(t\right)}{dt}+{\lambda}^{1/n}i\left(t\right)={v}^{\left(1\right)}\left(t\right);\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1$ (57)

be excited by a source which is a delta function say ${v}^{\left(1\right)}\left(t\right)=\left(\left({I}_{0}\right)\delta \left(t\right)\right)$ ; at $t=0$ . This means $v\left(t\right)=\left({I}_{0}\right)u\left(t\right)$ ; where $u\left(t\right)$ is unit step function at $t=0$ . With this excitation the relaxation current would be fractional integral of the input excitation that is from as depicted in above derivation (56) i.e. $i\left(t\right)=\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}\mathcal{I}{}_{t}^{\left(1-n\right)}\left[{I}_{0}\delta \left(t\right)\right]\right)$ . We have fractional integration of delta

function [12] [31] as ${}_{0}\mathcal{I}{}_{x}^{\alpha}\left[\delta \left(x\right)\right]=\frac{1}{\Gamma \left(\alpha \right)}{x}^{\alpha -1}$ ; and using this formula we get

following

$\begin{array}{l}i\left(t\right)=\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}\mathcal{I}{}_{t}^{\left(1-n\right)}\left[{I}_{0}\delta \left(t\right)\right]\right)\\ i\left(t\right)={I}_{0}\left(\frac{\Gamma \left(1+n\right)}{{t}^{n}}\right)\end{array}$ (58)

This $i\left(t\right)={I}_{0}\left(\Gamma \left(1+n\right)\right){t}^{-n}$ was what was derived in (49) above as impulse response (where ${I}_{0}=1$ ) i.e. $g\left(t\right)={t}^{-n}\left(\Gamma \left(1+n\right)\right)$ .

If the excitation source of (57) is a step function as ${v}^{\left(1\right)}\left(t\right)={I}_{0}\left(u\left(t\right)\right)$ at $t=0$ ; meaning $v\left(t\right)=\left({I}_{0}\right)t;t\ge 0$ where the unit step function is $u\left(t\right)=1,t\ge 0$ ; $u\left(t\right)=0,t<0$ then the relaxation current is fractional integration of order $\left(1-n\right)$ ; that is $i\left(t\right)=\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}\mathcal{I}{}_{t}^{\left(1-n\right)}\left[{I}_{0}u\left(t\right)\right]\right)$ . Using the formula for

fractional integration of a constant i.e. ${}_{0}\mathcal{I}{}_{x}^{\alpha}C=\frac{C}{\Gamma \left(1+\alpha \right)}{x}^{\alpha}$ [12] [31] we have;

the relaxation current as following

$\begin{array}{l}i\left(t\right)=\Gamma \left(1+n\right)\Gamma \left(1-n\right)\left({}_{0}\mathcal{I}{}_{t}^{\left(1-n\right)}\left[{I}_{0}u\left(t\right)\right]\right)\\ i\left(t\right)=\left(\Gamma \left(1+n\right)\Gamma \left(1-n\right)\right)\left({I}_{0}\right)\frac{{t}^{1-n}}{\Gamma \left(2-n\right)}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}=\frac{\Gamma \left(1+n\right)}{1-n}\left({I}_{0}\right){t}^{1-n},\text{\hspace{0.17em}}\text{\hspace{0.17em}}0<n<1\end{array}$ (59)

15. Conclusion

The empirical law that is Curie-von Schweidler law, which is a type of non-Debye relaxation, (that is also considered to be universal law of dielectric relaxation of current, since late 19^{th} century), states when dielectric is stressed with a constant voltage, gives relaxation current as
$i\left(t\right)\sim {t}^{-n},0<n<1$ . We derived the rate distribution-function (the histogram function) for current relaxation as Zipf’s power law distribution. The histogram function we found out to be of a function of type
${H}_{\lambda}\left(\lambda \right)\sim {\lambda}^{n-1},0<n<1$ . We infer the Curie-von Schweidler relaxation
$i\left(t\right)\sim {t}^{-n},0<n<1$ is simultaneous multi-body relaxations which has a distribution for relaxation rates i.e. right-skewed. That is with large number of relaxations with lower value of rate (slow rates) followed with long tail of small number of relaxations with faster relaxation rates, relaxing simultaneously. We noted that the possibility of having Zipfian distribution arises due to very-very large ratio of maximum to minimum spread in the relaxation rates l’s, and possibility of connected exponential distribution of many body simultaneous relaxations. The method we obtained for getting rate distributions of relaxation rates via formation of Laplace integral. When this method is extended for finding distribution of time constants though mathematically correct yet gave contrary physical interpretation. Thus we carried out the entire discussion with rate distribution functions i.e.
${H}_{\lambda}\left(\lambda \right)$ and not the time constant distribution function i.e.
${H}_{\tau}\left(\tau \right)$ . We also showed that Curie-von Schweidler law gives constituent of current and voltage of capacitor via use of fractional derivative, i.e.
$i\left(t\right)\sim {}_{0}D{}_{t}^{n}\left[v\left(t\right)\right];0<n<1$ , unlike classical capacitor relation i.e.
$i\left(t\right)\sim {D}_{t}^{\left(1\right)}\left[v\left(t\right)\right]$ . We verified by using obtained by Zipf’s distribution as power law for Curie-von Schwidler current relaxation law, assuming the scale dependence equivalent relaxation rate in the classical charging equation of capacitor with scale of relaxation varying from zero to infinity, i.e.
${\lambda}_{eq}\sim {\lambda}^{1/n}$ . We also related the Curie-von Schweidler relaxation law which gives a time varying rate i.e.
$\lambda \left(t\right)=n/t$ , indicating that the relaxation starts with very-very high rate, and becomes slower and slower with elapse of time. The paper gives a possible foundation for further studies in obtaining the rate relaxation distribution functions for other non-Debye type relaxation functions, and new type of explanation regarding reasons of Zipfian distributions.

Acknowledgements

I sincerely acknowledge support received from Prof Sujata Tarafdar, Prof Tapati Dutta, Research Scholars Smoasri Hazra, Moutushi Dutta Choudhury, Tania Basu, Simantini Mazumdar (Dept. of Physics, CMPRC Jadavpur University) students Adreja Mondol Rivu Gupta of St Xaviers Univ. Kolkata to have put to use concept of fractional calculus in relaxation experiments in visco-elastic studies, cooling law dynamics, impedance spectroscopy, studies in crack formation under electric stress. I also acknowledge Dr N C Pramanik Scientist CMET Thrissur, Prof Subhojit Ghosh and Research Scholar Manoranjan Kumar of Dept. of EE NIT Raipur for using fractional calculus in super-capacitor characterization experiments.

Cite this paper

Das, S. (2017) Revisiting the Curie-Von Schweidler Law for Dielectric Relaxation and Derivation of Distribution Function for Relaxation Rates as Zipf’s Power Law and Manifestation of Fractional Differential Equation for Capacitor. Journal of Modern Physics, 8, 1988- 2012. https://doi.org/10.4236/jmp.2017.812120

References

- 1. Jaques, C. (1889) Annales de Chimie et de Physique, 17, 384-434.
- 2. Schweidler, E.R. (1907) Annalen der Physik, 329, 711-770.
- 3. Jonscher, A.K. (1983) Dielectric Relaxation in Solids. Chelsea Dielectrics Press Limited.
- 4. Jameson, N.J., Azarian, M.H. and Pecht, M. (2017) Thermal Degradation of Polyimide Insulation and Its Effect on Electromagnetic Coil Impedance. Proceedings of the Society for Machinery Failure Prevention Technology 2017 Annual Conference.
- 5. Westerlund, S. and Ekstam, L. (1994) IEEE Trans on Dielectrics and Insulation, 1, 826. https://doi.org/10.1109/94.326654
- 6. Chiniba, S. and Tobazeon, R. (1985) Long Term Break-Down of Polypropylene Films in Different Media. International Conference on Properties and Application of Dielectric Materials, Xian, 24-29 June 1985, 566-584.
- 7. Borough, J.W., Brammer, W.G. and Burnham, J. (1980) Degradation of PVDF Capacitors during Accelerated Tests. Report Huges Aircraft Company, Culver City, 219-223.
- 8. Thorborg, K. (1985) Power Electronics. Chalmer Tech. Univ, Goteborg.
- 9. Jonscher, A.K. (1983) Dielectric Relaxation in Solids. Chelsea Dielectric Press, London.
- 10. Blythe, A.R. (1979) Electrical Properties of Polymer. Cambridge University Press.
- 11. Das, S. and Pramanik, N.C. (2013) International Journal of Mathematics and Computation, 20, 94-113.
- 12. Das, S. (2011) Functional Fractional Calculus. 2nd Edition, Springer-Verlag, Berlin. https://doi.org/10.1007/978-3-642-20545-3
- 13. Kumar, M., Ghosh, S. and Das, S. (2015) Journal of Power Sources, 335, 98-104. https://doi.org/10.1016/j.jpowsour.2016.10.024
- 14. Kumar, M., Ghosh, S. and Das, S. (2017) International Journal of Electronics & Communications, 78, 2714-280. https://doi.org/10.1016/j.aeue.2017.05.011
- 15. Elwakil, A.S., Allagui, A., Maundy, B.J. and Psycchalinos, C.A. (2016) International Journal of Electronics and Communications, 70, 970-973. https://doi.org/10.1016/j.aeue.2016.03.020
- 16. Freebon, T.J., Maundy, B. and Elwakil, A.S. (2015) Mater Renew Sustain Energy, 4, 1-7.
- 17. Freebon, T.J., Maundy, B. and Elwakil, A.S. (2013) IEEE Journal on Emerging and Selected Topics in Circuits and Systems, 3, 367-376. https://doi.org/10.1109/JETCAS.2013.2271433
- 18. Choudhury, M.D., Chandra, S., Nag, S., Das, S. and Tarafdar, S. (2011) Spreading of Non-Linear and Newtonian Fluids on a Solid Substrate under Pressure. Journal of Physics: Conference Series 319, IOP Publishing.
- 19. Choudhury, M.D., Chandra, S., Nag, S., Das, S. and Tarafdar, S. (2016) Colloids and Surfaces A: Physicochemical and Engineering Aspects, 492, 47-53. https://doi.org/10.1016/j.colsurfa.2015.12.007
- 20. Choudhury, M.D., Chandra, S., Nag, S., Das, S. and Tarafdar, S. (2012) Colloids and Surfaces A: Physicochemical and Engineering Aspects, 407, 64-70. https://doi.org/10.1016/j.colsurfa.2012.05.008
- 21. Majumdar, S., Hazra, S., Choudhury, M.D., Sinha, S.D., Das, S., Middya, T.R., Tarafdar, S. and Dutta, T. (2017) Colloids and Surfaces A: Physicochemical and Engineering Aspects, 516, 181-189. https://doi.org/10.1016/j.colsurfa.2016.12.019
- 22. Powers, D.M.W. (1998) Application and Explanation of Zipf’s Law. Association for Computational Linguistics, 151-160.
- 23. Bak, P. (1996) How Nature Works: The Science of Self-Organized Criticality. Springer-Verlag, New York. https://doi.org/10.1007/978-1-4757-5426-1
- 24. Zipf, G.K. (1949) Human Behavior and the Principle of Least Effort. Addison-Wesley.
- 25. Sornette, D. (2001) International Journal of Modern Physics C, 13, 133-136. https://doi.org/10.1142/S0129183102003036
- 26. Cancho, R.F. and Sole, R.V. (2003) Proceedings of the National Academy of Sciences of USA, 100, 788-791. https://doi.org/10.1073/pnas.0335980100
- 27. Easley, D. and Kleinberg, J. (2010) Power Laws and Rich Get Richer Phenomena. (Chapter-18): Networks Crowds and Markets: Reasoning about a Highly Connected World. Cambridge University Press. https://doi.org/10.1017/CBO9780511761942
- 28. Nagrath, I.J. and Gopal, M. (1982) Control System Engineering. 2nd Edition, Wiley Eastern Limited New Age International Limited.
- 29. Polyanin, A.D. and Manzhirov, A.V. (2008) Handbook of Integral Equations. 2nd Edition, Chapman & Hall/CRC, Boca Raton. https://doi.org/10.1201/9781420010558
- 30. Abramowitz, M. and Stegun, I.A. (1964) Handbook of Mathematical Functions. Nat. Bur. Stand. Appl. Math. Ser 55, U.S. Govt. Printing Office, Washington DC.
- 31. Oldham, K.B. and Spanier, J. (1974) The Fractional Calculus. Academic Press.
- 32. Fourier, J.B. (1822) The Analytical Theory of Heat-English. The University Press, 408.
- 33. Berberan-Santos, M.N. (2005) Journal of Mathematical Chemistry, 38, No. 2.
- 34. Berberan-Santos, M.N. (2005) Journal of Mathematical Chemistry, 38, No. 4.
- 35. Hazra, S., Dutta, T., Das, S. and Tarafdar, S. (2017) Memory of Electric Field in Laponite and How It Affects Crack Formation: Modeling through Generalized Calculus. Langmuir. https://doi.org/10.1021/acs.langmuir.7b02034