# Nonlinearities in Drug Release Process from Polymeric Microparticles: Long-Time-Scale Behaviour

Authors: Elena Simona Bacaita; Costica Bejinariu; Borsos Zoltan; Catalina Peptu; Gabriela Andrei; Marcel Popa; Daniela Magop; Maricel Agop
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

A theoretical model of the drug release process from polymeric microparticles (a particular type of polymer matrix), through dispersive fractal approximation of motion, is built. As a result, the drug release process takes place through cnoidal oscillations modes of a normalized concentration field. This indicates that, in the case of long-time-scale evolutions, the drug particles assemble in a lattice of nonlinear oscillators occur macroscopically, through variations of drug concentration. The model is validated by experimental results.


---

## Body

## 1. Introduction

Polymer matrices can be produced in one of the following forms: micro/nanoparticles, micro/nanocapsules, hydrogels, films, and patches. Due to the multitude of biocompatible polymers in the experimental protocol, drug proper delivery via many administration routes occurs. No matter what their form might be, drug carrier polymer matrices should have the following characteristics: biocompatibility, biodegradability, and controlled release capacity. The last one refers to the relationship between the efficient, nontoxic drug administration and therapeutic window type concentration, that is, minimum concentration is required to produce the wanted effect, but in the case of high levels, a toxic barrier occurs.Given the importance of the released drug concentration, numerous studies have been performed with the purpose of identifying the mathematical function that describes time dependence. Many papers show how various factors, such as polymer molecular weight [1, 2], polymer chemistry, monomer ratios [3, 4], pH of release media, additives to the release media [5, 6], and particle size [7], affect the release kinetics. At the same time, certain phenomena appearing in the release process have been studied. Of these, we mention (in approximate order of their occurrence) polymer swelling and degradation [8–11], drug dissolution and diffusion [12, 13], and above all, permanent chemical and physical interaction among components (drug, polymer, and release medium). Since all these phenomena are not independent, their analysis becomes complicated; consequently, it will not be possible to treat them separately and cumulate the effects. For example, microparticle morphology changes due to polymer degradation, their surfaces becoming highly porous. This will lead to increased diffusion coefficients and hence certain connected phenomena, such as polymer degradation and drug diffusion [7].The multitude of phenomena and dependencies occurring in drug release process as well as numerous structural entities (polymer, drug, and release medium) will turn the system into a complex one. Consequently, the complete theoretical analysis becomes more difficult in terms of performing.Nevertheless, significant amount of work has been accomplished in mathematical modeling, with the purpose of predicting the concentration of the released drug and providing the analysis of fundamental processes that govern release. Higuchi [14] was among the first who produced a drug release model from nonswelling and nondissolving polymer matrices, assuming that such phenomenon is purely controlled through diffusion. A number of other models have also been proposed in order to predict drug release in the case of erosion [9, 11], swelling [8], and dissolution [13] influenced processes. These mathematical patterns have chosen only two phenomena, with the purpose of simplifying mathematical modeling, which, otherwise, proves to be quite difficult.That is the main reason why it is necessary to use alternative approaches with reduced number of the approximations. One of such possible approaches is the fractal one [15, 16]. Its use is justified by natural and synthetic polymers that have been included in the category of fractional-dimensioned objects whose structures and behaviour can be described by means of fractal geometry [17, 18]. Moreover, it has been observed that the dynamics of drug release systems is a fractal one, because, in spite of complex phenomena and factors, mathematical expressions describing drug release kinetics from a variety of polymer matrices are power type laws (Higuchi [14] for nonswelling and nondissolving polymer, Ritger and Peppas [19] for nonswellable polymer in the form of slabs, spheres, cylinders, or discs, Peppas Sahlin [20] for solute release, Alfrey et al. [21] for diffusion in glassy polymers, etc.) specific for the fractal system evolution [22]. At the same time, it is quite important to emphasize that correlation of experimental data with the above-mentioned laws revealed good correlation in the first part (approximate 60%) of the release kinetics, the correlation coefficient decreasing according to time evolution.Studies on the release from different types of systems (HPMC matrices [23], inert porous matrices [24], and sponges [25, 26]) have been performed. Such approach analyzes drug release kinetics through Monte Carlo simulation. In this perspective, release systems are considered as three-dimensional lattices with leak sites located at the boundaries of the lattice pattern. Particles are free to move inside the porous network according to the random walk model of the Fickian diffusion (the moving particles act as hard spheres colliding with each other and having no possibility to mutually penetrate).The first studies started with simplifying approximations. Kosmidis et al. [23] consider that porosity has a constant value. Later on, Villalobos et al. [24] improved the model, assumed that network porosity behaves dynamically, and considered the effects of drug spatial distribution and initial drug concentration. All these approaches proved the validity of Weibull function (a continuous probability distribution function) for the entire release kinetics and consequently eliminate Peppas’ temporal limitation of the equation and criticism lacking kinetic basis and physical nature of parameters [27].Our new approach considers the entire system (drug-loaded polymer matrix in the release environment) as a type of “fluid” totally lacking interaction or neglecting physical interactions among particles. At the same time, the induced complexity is replaced by fractality. This will lead to particles moving on certain trajectories called geodesics within fractal space. This assumption represents the basis of the fractal approximation of motion in scale relativity theory (SRT) [28, 29], leading to a generalized fractal “diffusion” equation that can be analyzed in terms of two approximations (dissipative and dispersive).The comparison between dissipative approximation (with dominant convective and dissipative processes) and the dispersive one allows theoretical demonstration of Weibull function that best describes the behaviour of drug release systems at short time scales. These phenomena will be the object of subsequent analysis since they are responsible for certain types of behaviour and characterized by  high degree of nonlinearity in drug release systems.This paper is structured as follows: theoretical model (Section2), experimental results that validate the theoretical model (Section 3), and conclusions (Section 4).


## 2. Theoretical Model

### 2.1. Consequences of Nondifferentiability

We suppose that the drug release process takes place on continuous, but nondifferentiable curves (fractal curves). Then, nondifferentiability implies [28–30] the following.(i)
A continuous and a nondifferentiable curve (or almost nowhere differentiable) is explicitly scale dependent, and its length tends to infinity, when the scale interval tends to zero. In other words, a continuous and nondifferentiable space is fractal, and in the general meaning Mandelbrot used this concept [15];(ii)
Physical quantities will be expressed through fractal functions, namely, through functions that are dependent both on coordinate field and resolution scale. The invariance of the physical quantities in relation with the resolution scale generates special types of transformations, called resolution-scale transformations. In what follows, we will explain the above statement.LetF(x) be a fractal function in the interval x∈[a,b], and let the sequence of values for x be
(2.1)xa=x0,x1=x0+ε,…,xk=x0+kε,xn=x0+nε=xb.Let us note byF(x,ε) the broken line that connects the points
(2.2)F(x0),…,F(xk),…,F(xn).We can now say thatF(x,ε) is a ε-scale approximation.Let us now considerF(x,ε-) as a ε¯-scale approximation of the same function. Since F(x) is everywhere almost self-similar, if ε and ε¯ are sufficiently small, both approximations F(x,ε) and F(x,ε-) must lead to the same results; in the particular case, a fractal phenomenon is studied through approximation. By comparing the two cases, one notices that scale expansion is related to the increase dε of ε, according to an increase dε¯ of ε¯ (see Figure 1). But, in this case, we have
(2.3)dεε=dε¯ε¯=dρ,
a situation in which we can consider the infinitesimal-scale transformation as being
(2.4)ε′=ε+dε=ε+εdρ.Figure 1
Dilatation-scale invariance.Such transformation in the case of functionF(x,ε) leads to
(2.5)F(x,ε′)=F(x,ε+εdρ),
respectively, if we limit ourselves to a first-order approximation:
(2.6)F(x,ε′)=F(x,ε)+∂F(x,ε)∂ε(ε'-ε)=F(x,ε)+∂F(x,ε)∂εεdρ.Moreover, let us notice that for an arbitrary but fixedε0, we obtain
(2.7)∂ln(ε/ε0)∂ε=∂(lnε-lnε0)∂ε=1ε,
a situation in which (2.6) can be written as follows:
(2.8)F(x,ε′)=F(x,ε)+∂F(x,ε)∂ln(ε/ε0)dρ=[1+∂∂ln(ε/ε0)dρ]F(x,ε).Therefore, we can introduce the dilatation operator(2.9)D^=∂∂ln(ε/ε0).At the same time, relation (2.9) shows that the intrinsic variable of resolution is not ε, but ln(ε/ε0).The fractal function is explicitly dependant on the resolution(ε/ε0); therefore, we have to solve the differential equation
(2.10)dFdln(ε/ε0)=P(F),
where P(F) is now an unknown function. The simplest explicit suggested form for P(F) is linear dependence [29]
(2.11)P(F)=A+BF,A,B=const.,
in which case the differential equation (2.10) takes the form
(2.12)dFdln(ε/ε0)=A+BF.Hence, by integration and substituting(2.13a)B=-τ,(2.13b)-AB=F0,we obtain
(2.14)F(εε0)=F0[1+(ε0ε)τ].This solution is independent as compared to parameterization on fractal curve.We can now generalize the previous result by considering thatF is dependent on parameterization of the fractal curve. If p characterizes the position on the fractal curve, then, following the same algorithm as above, the solution will be as a sum of two terms, that is, both classical and differentiable (depending only on position) and fractal, nondifferentiable (depending on position and, divergently, on ε/ε0)
(2.15)F(p,εε0)=F0(p)[1+ξ(p)(ε0ε)τ(p)],
where ξ(p) is a function depending on parameterization of the fractal curve.The following particular cases are to be considered.(
ii
1
) In asymptotic small-scale regimeε≪ε0, τ is constant (with no scale dependence) and power-law dependence on resolution is obtained:(2.16a)F(p,εε0)=T(p)(ε0ε)τ,(2.16b)T(p)=F0(p)Q(p).At this stage, some power laws should also be considered, namely, those equations describing drug release kinetics from a different type of polymer matrix [14, 19–21]. Consequently, through the appropriate correspondence among quantities from (2.16a) and (2.16b) and those from drug release processes, we will obtain the following: (a)
Higuchi law:(2.17)MtM∞=kH·t1/2,whereMt and M∞ are cumulative amounts of drug release at time t and infinity, respectively, and kH is a constant characteristic of the system [14];(b)
Peppas law:(2.18)MtM∞=k·tn,where k is an experimentally obtained parameter, and n is a real number geometrically  related to the system and to drug release mechanism. The n value is used to characterize different release mechanisms, that is, n=0.5 indicates a Fickian diffusion. In their turn, different from 0.5  n values refer to mass transport according to non-Fickian model [31]. This equation is a generalization of a square-root time law and an approximation for short times of Weibull function.In these expressions, we recognize the standard form of self-similar fractal behaviour with fractal dimensionDF=DT+τ, which has already been used for accurately describing many physical and biological systems [15]. The topological dimensions are here DT=1, since we deal with length, but this can be easily generalized to surfaces (DT=2) and volumes (DT=3). Therefore, such result is not a consequence of postulation or deduction, but an aftermath of first principle theoretical analysis.Considering that the resolutionε is a length, ε=δX, the scale-dependent length is given, by definition, by the law
(2.19)X(p,δX)=X0(p)·(λδX)DF-1,
where λ is a scale characteristic length, and the exponent is identified with τ=DF-1.Now, in the above solution, one may use timet as parameter, and if one constantly moves along the curve, one obtains X0(t)=at. Then, a differential version of the above relation will be
(2.20)δX=aδt·(λδX)DF-1,
so that the following fundamental relation among space and time elements on a fractal curve or function is obtained:
(2.21)δXDF∝δt.
In other words, they are differential elements of different orders.(
ii
2
) In the asymptotic big-scale regime ε≫ε0,τ is constant (with no scale dependence), and, in terms of resolution, one obtains an independent law
(2.22)F(p,εε0)→F0(p).Particularly, ifF(p,ε/ε0) are the coordinates in given space, we can write
(2.23)X(p,εε0)=x(p)[1+ξ(p)(ε0ε)τ].In this situation,ξ(p) becomes a highly fluctuating function which can be described by stochastic process, while τ represents (according to previous description) the difference between fractal and topological dimensions. The result is a sum of two terms, a classical, differentiable one (dependent only on the position) and a fractal, nondifferentiable one (dependent both on the position and, divergently, on ε/ε0).  This represents the importance of the above analysis.By differentiating these two parts, we obtain(2.24)dX=dx+dξ,
where dx is the classical differential element, and dξ is a differential fractal one.(iii) There is infinity of fractal curves (geodesics) relating to any couple of points (or starting from any point) and applied for any scale. The phenomenon can be easily understood at the level of fractal surfaces, which, in their turn, can be described in terms of fractal distribution of conic points of positive and negative infinite curvature. As a consequence, we have replaced velocity on a particular geodesic by fractal velocity field of the whole infinite ensemble of geodesics. This representation is similar to that of fluid mechanics [32] where the motion of the fluid is described in terms of its velocity field v=(x(t),t), density ρ=(x(t),t), and, possibly, its pressure. We will, indeed, recover the fundamental equations of fluid mechanics (Euler and continuity equations), but we will write them in terms of a density of probability (as defined by the set of geodesics) instead of a density of matter and adding an additional term of quantum pressure (the expression of fractal geometry).(iv) The local differential time invariance is broken, so the time derivative of the fractal fieldQ can be written as twofold:(2.25a)d+Qdt=limΔt→0+Q(t+Δt)-Q(t)Δt,(2.25b)d-Qdt=limΔt→0-Q(t)-Q(t-Δt)Δt.Both definitions are equivalent in the differentiable casedt→-dt. In the nondifferentiable situation, these definitions are no longer valid, since limits are not defined anymore. Fractal theory defines physics in relationship with the function behavior during the “zoom” operation on the time resolution δt, here identified with the differential element dt (substitution principle), which is considered an independent variable. The standard field Q(t) is therefore replaced by fractal field Q(t,dt), explicitly dependent on time resolution interval, whose derivative is not defined at the unnoticeable limit dt→0. As a consequence, this leads to the two derivatives of the fractal field Q as explicit functions of the two variables t and dt,(2.26a)d+Qdt=limΔt→0+Q(t+Δt,Δt)-Q(t,Δt)Δt,(2.26b)d-Qdt=limΔt→0-Q(t,Δt)-Q(t-Δt,Δt)Δt.Notation “+” corresponds to the forward process, while “−” to the backward one.(v) LetP(x1,x2) be a point of the fractal curve, and let us consider a line which starts from this point. Let Mbe the first intersection of this line with the fractal curve. By dX+i, we denote the components of the vector PM, to the right of the line (d), and by dX-i the components of the vector PM′, to the left of the line (d)—see Figure 2.Figure 2
The continuous curves which are not fractals but have certain points where they are not differentiable.If we consider all the lines (segments) starting fromP, we denote the average of these vectors by dx±i, that is, (2.27a)〈dX+i〉=dx+i,i=1,2,(2.27b)〈dX-i〉=dx-i,i=1,2.Since, according to (2.24), we can write(2.28a)dX+i=dx+i+dξ+i,(2.28b)dX-i=dx-i+dξ-i,it results that(2.29a)〈dξ+i〉=0,(2.29b)〈dξ-i〉=0.(vi) The differential fractal part satisfies, according to (2.21), the fractal equation(2.30a)d+ξi=λ+i(dt)1/DF,(2.30b)d-ξi=λ-i(dt)1/DF,whereλ+i and λ-i are some constant coefficients, and DF is a constant fractal dimension. We note that the use of any Kolmogorov or Hausdorff [15, 28, 33–35] definitions can be accepted for fractal dimension, but once a certain definition is admitted, it should be used until the end of analyzed dynamics.(vii) The local differential time reflection invariance is recovered by combining the two derivatives,d+/dt and d-/dt, in the complex operator
(2.31)d^dt=12(d++d-dt)-i2(d+-d-dt).Applying this operator to the “position vector,” a complex speed yields(2.32)V^=d^Xdt=12(d+X+d-Xdt)-i2(d+X-d-Xdt)=V++V-2-iV+-V-2=V-iU,
with(2.33a)V=V++V-2,(2.33b)U=V+-V-2.The real part,V, of the complex speed V^ represents the standard classical speed, which does not depend on resolution, while the imaginary part, U, is a new quantity coming from  resolution-dependant fractal.


### 2.2. Covariant Total Derivative in Drug Release Mechanism

Let us now assume that curves describing drug release (continuous but nondifferentiable) are immersed in a 3-dimensional space, and thatX of components Xi(i=1,3¯) is the position vector of a point on the curve. Let us also consider the fractal field Q(X,t) and expand its total differential up to the third order(2.34a)d+Q=∂Q∂tdt+∇Q·d+X+12∂2Q∂Xi∂Xjd+Xid+Xj+16∂3Q∂Xi∂Xj∂Xkd+Xid+Xjd+Xk,(2.34b)d-Q=∂Q∂tdt+∇Q·d-X+12∂2Q∂Xi∂Xjd-Xid-Xj+16∂3Q∂Xi∂Xj∂Xkd-Xid-Xjd-Xk,where only the first three terms were used in Nottale’s theory (i.e., second-order terms in the motion equation). Relations (2.34a) and (2.34b) are valid in any point both for the spatial manifold and for the points X on the fractal curve (selected in relations (2.34a) and (2.34b)). Hence, the forward and backward average values of these relations take the form
(2.35)〈d±Q〉=〈∂Q∂tdt〉+〈∇Q·d±X〉+12〈∂2Q∂Xi∂Xjd±Xid±Xj〉+16〈∂3Q∂Xi∂Xj∂Xkd±Xid±Xjd±Xk〉,(2.36)〈d-Q〉=〈∂Q∂tdt〉+〈∇Q·d-X〉+12〈∂2Q∂Xi∂Xjd-Xid-Xj〉+16〈∂3Q∂Xi∂Xj∂Xkd-Xid-Xjd-Xk〉.The following aspects should be mentioned: the mean value of functionf and its derivatives coincide with themselves, and the differentials d±Xi and dt are independent; therefore, the average of their products coincides with the product of averages. Consequently, (2.35) and (2.36) become(2.37a)d+Q=∂Q∂tdt+∇Q〈d+X〉+12∂2Q∂Xi∂Xj〈d+Xid+Xj〉+16∂3Q∂Xi∂Xj∂Xk〈d+Xid+Xjd+Xk〉,(2.37b)d-Q=∂Q∂tdt+∇Q〈d-X〉+12∂2Q∂Xi∂Xj〈d-Xid-Xj〉+16∂3Q∂Xi∂Xj∂Xk〈d-Xid-Xjd-Xk〉,or more, using (2.28a) and (2.28b) with characteristics (2.29a) and (2.29b),(2.38a)d+Q=∂Q∂tdt+∇Q·d+X+12∂2Q∂Xi∂Xj(d+xid+xj+〈d+ξid+ξj〉)+16∂3Q∂Xi∂Xj∂Xk(d+xid+xjd+xk+〈d+ξid+ξjd+ξk〉),(2.38b)d-Q=∂Q∂tdt+∇Q·d-X+12∂2Q∂Xi∂Xj(d-xid-xj+〈d-ξid-ξj〉)+16∂3Q∂Xi∂Xj∂Xk(d-xid-xjd-xk+〈d-ξid-ξjd-ξk〉).Even if the average value of the fractal coordinated±ξi is null (see (2.29a) and (2.29b)), for higher order of fractal coordinate average, the situation can still be different. Firstly, let us focus on the averages 〈d+ξid+ξj〉 and 〈d-ξid-ξj〉. If i≠j, these averages are zero due to the independence of d±ξi and d±ξj. So, using (2.30a) and (2.30b), we can write(2.39a)〈d+ξid+ξj〉=λ+iλ+j(dt)(2/DF)-1dt,(2.39b)〈d-ξid-ξj〉=λ-iλ-j(dt)(2/DF)-1dt.Then, let us consider the averages〈d+ξid+ξjd+ξk〉 and 〈d-ξid-ξjd-ξk〉. If i≠j≠k, these averages are zero due to independence of d±ξi on d±ξj and d±ξk. Now, using (2.30a) and (2.30b), we can write(2.40a)〈d+ξid+ξjd+ξk〉=λ+iλ+jλ+k(dt)(3/DF)-1dt,(2.40b)〈d-ξid-ξjd-ξk〉=λ-iλ-jλ-k(dt)(3/DF)-1dt.Then, (2.38a) and (2.38b) may be written as follows:(2.41a)d+Q=∂Q∂tdt+d+x·∇Q+12∂2Q∂Xi∂Xjd+xid+xj+12∂2Q∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1dt+16∂3Q∂Xi∂Xj∂Xkd+xid+xjd+xk+16∂3Q∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1dt,(2.41b)d-Q=∂Q∂tdt+d-x·∇Q+12∂2Q∂Xi∂Xjd-xid-xj+12∂2Q∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1dt+16∂3Q∂Xi∂Xj∂Xkd-xid-xjd-xk+16∂3Q∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1dt.If we divide bydt and neglect the terms containing differential factors (for details on the method, see [36, 37]), (2.41a) and (2.41b) are reduced to(2.42a)d+Qdt=∂Q∂t+V+·∇Q+12∂2Q∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1+16∂3Q∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1,(2.42b)d-Qdt=∂Q∂t+V-·∇Q+12∂2Q∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1+16∂3Q∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1.These relations also allow us to define the operator(2.43a)d+dt=∂∂t+V+·∇+12∂2∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1+16∂3∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1,(2.43b)d-dt=∂∂t+V-·∇+12∂2∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1+16∂3∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1.Under these circumstances, let us calculate(∂^Q/∂t). Taking into account (2.43a), (2.43b), (2.31), and (2.32), we will obtain
(2.44)∂Q∧∂t=12[d+Qdt+d-Qdt-i(d+Qdt-d-Qdt)]=12∂Q∂t+12V+·∇Q+λ+iλ+j14(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ+iλ+jλ+k112(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk+12∂Q∂t+12V-·∇Q+λ-iλ-j14(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ-iλ-jλ-k112(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk-i2∂Q∂t-i2V+·∇Q-λ+iλ+ji2(dt)(2/DF)-1∂2Q∂Xi∂Xj-λ+iλ+jλ+ki12(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk+i2∂Q∂t+i2V-·∇Q+λ-iλ-ji2(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ-iλ-jλ-ki12(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk=∂Q∂t+(V++V-2-iV+-V-2)·∇Q+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2Q∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3Q∂Xi∂Xj∂Xk=∂Q∂t+V∧·∇Q+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2Q∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3Q∂Xi∂Xj∂Xk.This relation also allows us to define the fractal operator(2.45)∂^∂t=∂∂t+V^·∇+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3∂Xi∂Xj∂Xk.Particularly, by choosing(2.46a)λ+iλ+j=-λ-iλ-j=2𝒟δij,(2.46b)λ+iλ+jλ+k=-λ-iλ-jλ+k=22𝒟3/2δijk,the fractal operator (2.45) takes the usual form
(2.47)∂^∂t=∂∂t+V^·∇-i𝒟(dt)(2/DF)-1Δ+23𝒟3/2(dt)(3/DF)-1∇3.We now apply the principle of scale covariance and postulate that the passage from classical (differentiable) to “fractal” mechanics can be implemented by replacing the standard time derivative operator,d/dt, with the complex operator ∂^/∂t (this results in a generalization of Nottale’s [28, 29] principle of scale covariance). Consequently, we are now able to write the diffusion equation in its covariant form
(2.48)∂^Q∂t=∂Q∂t+(V^·∇)Q-i𝒟(dt)(2/DF)-1ΔQ+23𝒟3/2(dt)(3/DF)-1∇3Q=0.This means that at any point on a fractal path, the local temporal∂tQ, the nonlinear (convective), (V^·∇)Q, the dissipative, ΔQ, and the dispersive, ∇3Q, terms keep their balance.The dissipative approximation was applied for the drug release processes, and the result was a Weibull type function that was analyzed in [38, 39].  In what follows, we will focus on dispersive approximation.


### 2.3. The Dispersive Approximation

Let us now consider that, in comparison with dissipative processes, convective and dispersive processes are dominant ones. Consequently, we are now able to write the diffusion equation in its covariant form, as a Korteweg de Vries type equation(2.49)∂^Qdt=∂Q∂t+(V^·∇)Q+23D3/2(dt)(3/DFD)-1∇3Q=0.If we separate the real and imaginary parts from (2.49), we will obtain(2.50a)∂Q∂t+V·∇Q+23D3/2(dt)(3/DF)-1∇3Q=0,(2.50b)-U·∇Q=0.By adding them, the fractal diffusion equation is(2.51)∂Q∂t+(V-U)·∇Q+23D3/2(dt)(3/DF)-1∇3Q=0.From (2.50b), we see that, at fractal scale, there will be no Q field gradient.Assuming that|V-U|=σ·Q with σ=constant (in systems with self-structuring processes, the speed fluctuations induced by fractal/nonfractal are proportional with the concentration field [22]), in the particular one-dimensional case, (2.51) will lack parameters(2.52a)τ¯=ωt,(2.52b)ξ¯=kx,(2.52c)Φ=QQ0,and normalizing conditions(2.53)σQ0k6ω=23D3/2(dt)(3/DF)-1k3ω=1
take the form
(2.54)∂τ¯ϕ+6ϕ∂ξ¯ϕ+∂ξ-ξ-ξ-ϕ=0.In relations (2.52a), (2.52b), (2.52c), and (2.53), ω corresponds to a characteristic pulsation, k to the inverse of a characteristic length, and Q0 to balanced concentration.Through substitutions(2.55a)w(θ)=ϕ(τ¯,ξ¯),(2.55b)θ=ξ¯-uτ¯,(2.54), by double integration, becomes
(2.56)12w'2=F(w)=-(w3-u2w2-gw-h),
with g, h are two integration constants, and u is the normalized phase velocity. If F(w) has real roots, (2.54) has the stationary solution
(2.57)ϕ(ξ¯,τ¯,s)=2a(E(s)K(s)-1)+2a·cn2[as(ξ¯-u2τ¯+ξ0¯);s],
where cn is Jacobi’s elliptic function of s modulus [39], a is  the amplitude, ξ0¯ is a constant of integration, and(2.58a)K(s)=∫0π/2(1-s2sin2φ)-1/2dφ,(2.58b)E(s)=∫0π/2(1-s2sin2φ)1/2dφ,are the complete elliptic integrals [40].Parameters represents measure characterizing the degree of nonlinearity in the system. Therefore, the solution (2.57) contains (as subsequences for s=0) one-dimensional harmonic waves, while for, s→0 one-dimensional wave packet. These two subsequences define the nonquasiautonomous regime of the drug release process [22, 36, 37], that is, the system should receive external energy in order to develop. For s=1, the solution (2.57) becomes one-dimensional soliton, while for s→1, one-dimensional soliton packet will be generated. The last two imply a quasiautonomous regime (self-evolving and independent [22]) for drug particle release process [22, 36, 37].The three-dimensional plot of solution (2.57) shows one-dimensional cnoidal oscillation modes of the concentration field, generated by similar trajectories of the drug particles (see Figure 3). We mention that cnoidal oscillations are nonlinear ones, being described by the elliptic function cn, hence the name (cnoidal).Figure 3
One-dimensional cnoidal oscillation modes of the fieldΦ.It is known that in nonlinear dynamics, cnoidal oscillation modes are associated with nonlinear lattice of oscillators (the Toda lattice [41]). Consequently, large-time-scale drug particle ensembles can be compared to a lattice of nonlinear oscillators which facilitates drug release process.


## 2.1. Consequences of Nondifferentiability

We suppose that the drug release process takes place on continuous, but nondifferentiable curves (fractal curves). Then, nondifferentiability implies [28–30] the following.(i)
A continuous and a nondifferentiable curve (or almost nowhere differentiable) is explicitly scale dependent, and its length tends to infinity, when the scale interval tends to zero. In other words, a continuous and nondifferentiable space is fractal, and in the general meaning Mandelbrot used this concept [15];(ii)
Physical quantities will be expressed through fractal functions, namely, through functions that are dependent both on coordinate field and resolution scale. The invariance of the physical quantities in relation with the resolution scale generates special types of transformations, called resolution-scale transformations. In what follows, we will explain the above statement.LetF(x) be a fractal function in the interval x∈[a,b], and let the sequence of values for x be
(2.1)xa=x0,x1=x0+ε,…,xk=x0+kε,xn=x0+nε=xb.Let us note byF(x,ε) the broken line that connects the points
(2.2)F(x0),…,F(xk),…,F(xn).We can now say thatF(x,ε) is a ε-scale approximation.Let us now considerF(x,ε-) as a ε¯-scale approximation of the same function. Since F(x) is everywhere almost self-similar, if ε and ε¯ are sufficiently small, both approximations F(x,ε) and F(x,ε-) must lead to the same results; in the particular case, a fractal phenomenon is studied through approximation. By comparing the two cases, one notices that scale expansion is related to the increase dε of ε, according to an increase dε¯ of ε¯ (see Figure 1). But, in this case, we have
(2.3)dεε=dε¯ε¯=dρ,
a situation in which we can consider the infinitesimal-scale transformation as being
(2.4)ε′=ε+dε=ε+εdρ.Figure 1
Dilatation-scale invariance.Such transformation in the case of functionF(x,ε) leads to
(2.5)F(x,ε′)=F(x,ε+εdρ),
respectively, if we limit ourselves to a first-order approximation:
(2.6)F(x,ε′)=F(x,ε)+∂F(x,ε)∂ε(ε'-ε)=F(x,ε)+∂F(x,ε)∂εεdρ.Moreover, let us notice that for an arbitrary but fixedε0, we obtain
(2.7)∂ln(ε/ε0)∂ε=∂(lnε-lnε0)∂ε=1ε,
a situation in which (2.6) can be written as follows:
(2.8)F(x,ε′)=F(x,ε)+∂F(x,ε)∂ln(ε/ε0)dρ=[1+∂∂ln(ε/ε0)dρ]F(x,ε).Therefore, we can introduce the dilatation operator(2.9)D^=∂∂ln(ε/ε0).At the same time, relation (2.9) shows that the intrinsic variable of resolution is not ε, but ln(ε/ε0).The fractal function is explicitly dependant on the resolution(ε/ε0); therefore, we have to solve the differential equation
(2.10)dFdln(ε/ε0)=P(F),
where P(F) is now an unknown function. The simplest explicit suggested form for P(F) is linear dependence [29]
(2.11)P(F)=A+BF,A,B=const.,
in which case the differential equation (2.10) takes the form
(2.12)dFdln(ε/ε0)=A+BF.Hence, by integration and substituting(2.13a)B=-τ,(2.13b)-AB=F0,we obtain
(2.14)F(εε0)=F0[1+(ε0ε)τ].This solution is independent as compared to parameterization on fractal curve.We can now generalize the previous result by considering thatF is dependent on parameterization of the fractal curve. If p characterizes the position on the fractal curve, then, following the same algorithm as above, the solution will be as a sum of two terms, that is, both classical and differentiable (depending only on position) and fractal, nondifferentiable (depending on position and, divergently, on ε/ε0)
(2.15)F(p,εε0)=F0(p)[1+ξ(p)(ε0ε)τ(p)],
where ξ(p) is a function depending on parameterization of the fractal curve.The following particular cases are to be considered.(
ii
1
) In asymptotic small-scale regimeε≪ε0, τ is constant (with no scale dependence) and power-law dependence on resolution is obtained:(2.16a)F(p,εε0)=T(p)(ε0ε)τ,(2.16b)T(p)=F0(p)Q(p).At this stage, some power laws should also be considered, namely, those equations describing drug release kinetics from a different type of polymer matrix [14, 19–21]. Consequently, through the appropriate correspondence among quantities from (2.16a) and (2.16b) and those from drug release processes, we will obtain the following: (a)
Higuchi law:(2.17)MtM∞=kH·t1/2,whereMt and M∞ are cumulative amounts of drug release at time t and infinity, respectively, and kH is a constant characteristic of the system [14];(b)
Peppas law:(2.18)MtM∞=k·tn,where k is an experimentally obtained parameter, and n is a real number geometrically  related to the system and to drug release mechanism. The n value is used to characterize different release mechanisms, that is, n=0.5 indicates a Fickian diffusion. In their turn, different from 0.5  n values refer to mass transport according to non-Fickian model [31]. This equation is a generalization of a square-root time law and an approximation for short times of Weibull function.In these expressions, we recognize the standard form of self-similar fractal behaviour with fractal dimensionDF=DT+τ, which has already been used for accurately describing many physical and biological systems [15]. The topological dimensions are here DT=1, since we deal with length, but this can be easily generalized to surfaces (DT=2) and volumes (DT=3). Therefore, such result is not a consequence of postulation or deduction, but an aftermath of first principle theoretical analysis.Considering that the resolutionε is a length, ε=δX, the scale-dependent length is given, by definition, by the law
(2.19)X(p,δX)=X0(p)·(λδX)DF-1,
where λ is a scale characteristic length, and the exponent is identified with τ=DF-1.Now, in the above solution, one may use timet as parameter, and if one constantly moves along the curve, one obtains X0(t)=at. Then, a differential version of the above relation will be
(2.20)δX=aδt·(λδX)DF-1,
so that the following fundamental relation among space and time elements on a fractal curve or function is obtained:
(2.21)δXDF∝δt.
In other words, they are differential elements of different orders.(
ii
2
) In the asymptotic big-scale regime ε≫ε0,τ is constant (with no scale dependence), and, in terms of resolution, one obtains an independent law
(2.22)F(p,εε0)→F0(p).Particularly, ifF(p,ε/ε0) are the coordinates in given space, we can write
(2.23)X(p,εε0)=x(p)[1+ξ(p)(ε0ε)τ].In this situation,ξ(p) becomes a highly fluctuating function which can be described by stochastic process, while τ represents (according to previous description) the difference between fractal and topological dimensions. The result is a sum of two terms, a classical, differentiable one (dependent only on the position) and a fractal, nondifferentiable one (dependent both on the position and, divergently, on ε/ε0).  This represents the importance of the above analysis.By differentiating these two parts, we obtain(2.24)dX=dx+dξ,
where dx is the classical differential element, and dξ is a differential fractal one.(iii) There is infinity of fractal curves (geodesics) relating to any couple of points (or starting from any point) and applied for any scale. The phenomenon can be easily understood at the level of fractal surfaces, which, in their turn, can be described in terms of fractal distribution of conic points of positive and negative infinite curvature. As a consequence, we have replaced velocity on a particular geodesic by fractal velocity field of the whole infinite ensemble of geodesics. This representation is similar to that of fluid mechanics [32] where the motion of the fluid is described in terms of its velocity field v=(x(t),t), density ρ=(x(t),t), and, possibly, its pressure. We will, indeed, recover the fundamental equations of fluid mechanics (Euler and continuity equations), but we will write them in terms of a density of probability (as defined by the set of geodesics) instead of a density of matter and adding an additional term of quantum pressure (the expression of fractal geometry).(iv) The local differential time invariance is broken, so the time derivative of the fractal fieldQ can be written as twofold:(2.25a)d+Qdt=limΔt→0+Q(t+Δt)-Q(t)Δt,(2.25b)d-Qdt=limΔt→0-Q(t)-Q(t-Δt)Δt.Both definitions are equivalent in the differentiable casedt→-dt. In the nondifferentiable situation, these definitions are no longer valid, since limits are not defined anymore. Fractal theory defines physics in relationship with the function behavior during the “zoom” operation on the time resolution δt, here identified with the differential element dt (substitution principle), which is considered an independent variable. The standard field Q(t) is therefore replaced by fractal field Q(t,dt), explicitly dependent on time resolution interval, whose derivative is not defined at the unnoticeable limit dt→0. As a consequence, this leads to the two derivatives of the fractal field Q as explicit functions of the two variables t and dt,(2.26a)d+Qdt=limΔt→0+Q(t+Δt,Δt)-Q(t,Δt)Δt,(2.26b)d-Qdt=limΔt→0-Q(t,Δt)-Q(t-Δt,Δt)Δt.Notation “+” corresponds to the forward process, while “−” to the backward one.(v) LetP(x1,x2) be a point of the fractal curve, and let us consider a line which starts from this point. Let Mbe the first intersection of this line with the fractal curve. By dX+i, we denote the components of the vector PM, to the right of the line (d), and by dX-i the components of the vector PM′, to the left of the line (d)—see Figure 2.Figure 2
The continuous curves which are not fractals but have certain points where they are not differentiable.If we consider all the lines (segments) starting fromP, we denote the average of these vectors by dx±i, that is, (2.27a)〈dX+i〉=dx+i,i=1,2,(2.27b)〈dX-i〉=dx-i,i=1,2.Since, according to (2.24), we can write(2.28a)dX+i=dx+i+dξ+i,(2.28b)dX-i=dx-i+dξ-i,it results that(2.29a)〈dξ+i〉=0,(2.29b)〈dξ-i〉=0.(vi) The differential fractal part satisfies, according to (2.21), the fractal equation(2.30a)d+ξi=λ+i(dt)1/DF,(2.30b)d-ξi=λ-i(dt)1/DF,whereλ+i and λ-i are some constant coefficients, and DF is a constant fractal dimension. We note that the use of any Kolmogorov or Hausdorff [15, 28, 33–35] definitions can be accepted for fractal dimension, but once a certain definition is admitted, it should be used until the end of analyzed dynamics.(vii) The local differential time reflection invariance is recovered by combining the two derivatives,d+/dt and d-/dt, in the complex operator
(2.31)d^dt=12(d++d-dt)-i2(d+-d-dt).Applying this operator to the “position vector,” a complex speed yields(2.32)V^=d^Xdt=12(d+X+d-Xdt)-i2(d+X-d-Xdt)=V++V-2-iV+-V-2=V-iU,
with(2.33a)V=V++V-2,(2.33b)U=V+-V-2.The real part,V, of the complex speed V^ represents the standard classical speed, which does not depend on resolution, while the imaginary part, U, is a new quantity coming from  resolution-dependant fractal.


## 2.2. Covariant Total Derivative in Drug Release Mechanism

Let us now assume that curves describing drug release (continuous but nondifferentiable) are immersed in a 3-dimensional space, and thatX of components Xi(i=1,3¯) is the position vector of a point on the curve. Let us also consider the fractal field Q(X,t) and expand its total differential up to the third order(2.34a)d+Q=∂Q∂tdt+∇Q·d+X+12∂2Q∂Xi∂Xjd+Xid+Xj+16∂3Q∂Xi∂Xj∂Xkd+Xid+Xjd+Xk,(2.34b)d-Q=∂Q∂tdt+∇Q·d-X+12∂2Q∂Xi∂Xjd-Xid-Xj+16∂3Q∂Xi∂Xj∂Xkd-Xid-Xjd-Xk,where only the first three terms were used in Nottale’s theory (i.e., second-order terms in the motion equation). Relations (2.34a) and (2.34b) are valid in any point both for the spatial manifold and for the points X on the fractal curve (selected in relations (2.34a) and (2.34b)). Hence, the forward and backward average values of these relations take the form
(2.35)〈d±Q〉=〈∂Q∂tdt〉+〈∇Q·d±X〉+12〈∂2Q∂Xi∂Xjd±Xid±Xj〉+16〈∂3Q∂Xi∂Xj∂Xkd±Xid±Xjd±Xk〉,(2.36)〈d-Q〉=〈∂Q∂tdt〉+〈∇Q·d-X〉+12〈∂2Q∂Xi∂Xjd-Xid-Xj〉+16〈∂3Q∂Xi∂Xj∂Xkd-Xid-Xjd-Xk〉.The following aspects should be mentioned: the mean value of functionf and its derivatives coincide with themselves, and the differentials d±Xi and dt are independent; therefore, the average of their products coincides with the product of averages. Consequently, (2.35) and (2.36) become(2.37a)d+Q=∂Q∂tdt+∇Q〈d+X〉+12∂2Q∂Xi∂Xj〈d+Xid+Xj〉+16∂3Q∂Xi∂Xj∂Xk〈d+Xid+Xjd+Xk〉,(2.37b)d-Q=∂Q∂tdt+∇Q〈d-X〉+12∂2Q∂Xi∂Xj〈d-Xid-Xj〉+16∂3Q∂Xi∂Xj∂Xk〈d-Xid-Xjd-Xk〉,or more, using (2.28a) and (2.28b) with characteristics (2.29a) and (2.29b),(2.38a)d+Q=∂Q∂tdt+∇Q·d+X+12∂2Q∂Xi∂Xj(d+xid+xj+〈d+ξid+ξj〉)+16∂3Q∂Xi∂Xj∂Xk(d+xid+xjd+xk+〈d+ξid+ξjd+ξk〉),(2.38b)d-Q=∂Q∂tdt+∇Q·d-X+12∂2Q∂Xi∂Xj(d-xid-xj+〈d-ξid-ξj〉)+16∂3Q∂Xi∂Xj∂Xk(d-xid-xjd-xk+〈d-ξid-ξjd-ξk〉).Even if the average value of the fractal coordinated±ξi is null (see (2.29a) and (2.29b)), for higher order of fractal coordinate average, the situation can still be different. Firstly, let us focus on the averages 〈d+ξid+ξj〉 and 〈d-ξid-ξj〉. If i≠j, these averages are zero due to the independence of d±ξi and d±ξj. So, using (2.30a) and (2.30b), we can write(2.39a)〈d+ξid+ξj〉=λ+iλ+j(dt)(2/DF)-1dt,(2.39b)〈d-ξid-ξj〉=λ-iλ-j(dt)(2/DF)-1dt.Then, let us consider the averages〈d+ξid+ξjd+ξk〉 and 〈d-ξid-ξjd-ξk〉. If i≠j≠k, these averages are zero due to independence of d±ξi on d±ξj and d±ξk. Now, using (2.30a) and (2.30b), we can write(2.40a)〈d+ξid+ξjd+ξk〉=λ+iλ+jλ+k(dt)(3/DF)-1dt,(2.40b)〈d-ξid-ξjd-ξk〉=λ-iλ-jλ-k(dt)(3/DF)-1dt.Then, (2.38a) and (2.38b) may be written as follows:(2.41a)d+Q=∂Q∂tdt+d+x·∇Q+12∂2Q∂Xi∂Xjd+xid+xj+12∂2Q∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1dt+16∂3Q∂Xi∂Xj∂Xkd+xid+xjd+xk+16∂3Q∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1dt,(2.41b)d-Q=∂Q∂tdt+d-x·∇Q+12∂2Q∂Xi∂Xjd-xid-xj+12∂2Q∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1dt+16∂3Q∂Xi∂Xj∂Xkd-xid-xjd-xk+16∂3Q∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1dt.If we divide bydt and neglect the terms containing differential factors (for details on the method, see [36, 37]), (2.41a) and (2.41b) are reduced to(2.42a)d+Qdt=∂Q∂t+V+·∇Q+12∂2Q∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1+16∂3Q∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1,(2.42b)d-Qdt=∂Q∂t+V-·∇Q+12∂2Q∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1+16∂3Q∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1.These relations also allow us to define the operator(2.43a)d+dt=∂∂t+V+·∇+12∂2∂Xi∂Xjλ+iλ+j(dt)(2/DF)-1+16∂3∂Xi∂Xj∂Xkλ+iλ+jλ+k(dt)(3/DF)-1,(2.43b)d-dt=∂∂t+V-·∇+12∂2∂Xi∂Xjλ-iλ-j(dt)(2/DF)-1+16∂3∂Xi∂Xj∂Xkλ-iλ-jλ-k(dt)(3/DF)-1.Under these circumstances, let us calculate(∂^Q/∂t). Taking into account (2.43a), (2.43b), (2.31), and (2.32), we will obtain
(2.44)∂Q∧∂t=12[d+Qdt+d-Qdt-i(d+Qdt-d-Qdt)]=12∂Q∂t+12V+·∇Q+λ+iλ+j14(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ+iλ+jλ+k112(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk+12∂Q∂t+12V-·∇Q+λ-iλ-j14(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ-iλ-jλ-k112(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk-i2∂Q∂t-i2V+·∇Q-λ+iλ+ji2(dt)(2/DF)-1∂2Q∂Xi∂Xj-λ+iλ+jλ+ki12(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk+i2∂Q∂t+i2V-·∇Q+λ-iλ-ji2(dt)(2/DF)-1∂2Q∂Xi∂Xj+λ-iλ-jλ-ki12(dt)(3/DF)-1∂3Q∂Xi∂Xj∂Xk=∂Q∂t+(V++V-2-iV+-V-2)·∇Q+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2Q∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3Q∂Xi∂Xj∂Xk=∂Q∂t+V∧·∇Q+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2Q∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3Q∂Xi∂Xj∂Xk.This relation also allows us to define the fractal operator(2.45)∂^∂t=∂∂t+V^·∇+(dt)(2/DF)-14[(λ+iλ+j+λ-iλ-j)-i(λ+iλ+j-λ-iλ-j)]∂2∂Xi∂Xj+(dt)(3/DF)-112[(λ+iλ+jλ+k+λ-iλ-jλ-k)-i(λ+iλ+jλ+k-λ-iλ-jλ-k)]∂3∂Xi∂Xj∂Xk.Particularly, by choosing(2.46a)λ+iλ+j=-λ-iλ-j=2𝒟δij,(2.46b)λ+iλ+jλ+k=-λ-iλ-jλ+k=22𝒟3/2δijk,the fractal operator (2.45) takes the usual form
(2.47)∂^∂t=∂∂t+V^·∇-i𝒟(dt)(2/DF)-1Δ+23𝒟3/2(dt)(3/DF)-1∇3.We now apply the principle of scale covariance and postulate that the passage from classical (differentiable) to “fractal” mechanics can be implemented by replacing the standard time derivative operator,d/dt, with the complex operator ∂^/∂t (this results in a generalization of Nottale’s [28, 29] principle of scale covariance). Consequently, we are now able to write the diffusion equation in its covariant form
(2.48)∂^Q∂t=∂Q∂t+(V^·∇)Q-i𝒟(dt)(2/DF)-1ΔQ+23𝒟3/2(dt)(3/DF)-1∇3Q=0.This means that at any point on a fractal path, the local temporal∂tQ, the nonlinear (convective), (V^·∇)Q, the dissipative, ΔQ, and the dispersive, ∇3Q, terms keep their balance.The dissipative approximation was applied for the drug release processes, and the result was a Weibull type function that was analyzed in [38, 39].  In what follows, we will focus on dispersive approximation.


## 2.3. The Dispersive Approximation

Let us now consider that, in comparison with dissipative processes, convective and dispersive processes are dominant ones. Consequently, we are now able to write the diffusion equation in its covariant form, as a Korteweg de Vries type equation(2.49)∂^Qdt=∂Q∂t+(V^·∇)Q+23D3/2(dt)(3/DFD)-1∇3Q=0.If we separate the real and imaginary parts from (2.49), we will obtain(2.50a)∂Q∂t+V·∇Q+23D3/2(dt)(3/DF)-1∇3Q=0,(2.50b)-U·∇Q=0.By adding them, the fractal diffusion equation is(2.51)∂Q∂t+(V-U)·∇Q+23D3/2(dt)(3/DF)-1∇3Q=0.From (2.50b), we see that, at fractal scale, there will be no Q field gradient.Assuming that|V-U|=σ·Q with σ=constant (in systems with self-structuring processes, the speed fluctuations induced by fractal/nonfractal are proportional with the concentration field [22]), in the particular one-dimensional case, (2.51) will lack parameters(2.52a)τ¯=ωt,(2.52b)ξ¯=kx,(2.52c)Φ=QQ0,and normalizing conditions(2.53)σQ0k6ω=23D3/2(dt)(3/DF)-1k3ω=1
take the form
(2.54)∂τ¯ϕ+6ϕ∂ξ¯ϕ+∂ξ-ξ-ξ-ϕ=0.In relations (2.52a), (2.52b), (2.52c), and (2.53), ω corresponds to a characteristic pulsation, k to the inverse of a characteristic length, and Q0 to balanced concentration.Through substitutions(2.55a)w(θ)=ϕ(τ¯,ξ¯),(2.55b)θ=ξ¯-uτ¯,(2.54), by double integration, becomes
(2.56)12w'2=F(w)=-(w3-u2w2-gw-h),
with g, h are two integration constants, and u is the normalized phase velocity. If F(w) has real roots, (2.54) has the stationary solution
(2.57)ϕ(ξ¯,τ¯,s)=2a(E(s)K(s)-1)+2a·cn2[as(ξ¯-u2τ¯+ξ0¯);s],
where cn is Jacobi’s elliptic function of s modulus [39], a is  the amplitude, ξ0¯ is a constant of integration, and(2.58a)K(s)=∫0π/2(1-s2sin2φ)-1/2dφ,(2.58b)E(s)=∫0π/2(1-s2sin2φ)1/2dφ,are the complete elliptic integrals [40].Parameters represents measure characterizing the degree of nonlinearity in the system. Therefore, the solution (2.57) contains (as subsequences for s=0) one-dimensional harmonic waves, while for, s→0 one-dimensional wave packet. These two subsequences define the nonquasiautonomous regime of the drug release process [22, 36, 37], that is, the system should receive external energy in order to develop. For s=1, the solution (2.57) becomes one-dimensional soliton, while for s→1, one-dimensional soliton packet will be generated. The last two imply a quasiautonomous regime (self-evolving and independent [22]) for drug particle release process [22, 36, 37].The three-dimensional plot of solution (2.57) shows one-dimensional cnoidal oscillation modes of the concentration field, generated by similar trajectories of the drug particles (see Figure 3). We mention that cnoidal oscillations are nonlinear ones, being described by the elliptic function cn, hence the name (cnoidal).Figure 3
One-dimensional cnoidal oscillation modes of the fieldΦ.It is known that in nonlinear dynamics, cnoidal oscillation modes are associated with nonlinear lattice of oscillators (the Toda lattice [41]). Consequently, large-time-scale drug particle ensembles can be compared to a lattice of nonlinear oscillators which facilitates drug release process.


## 3. Experimental Results

Most of the experimental data in the literature reveal that, on average, drug release from polymeric matrices takes place according to a power law in the first 60% of the release curve and/or to exponential Weibull function on the entire drug release curve, reaching an average constant balanced value.  The majority of these experimental results are carried out on relatively short time intervals, dissolution and diffusion being the dominant systems. The system exhibits a “burst effect” due to highly concentrated gradient. The phenomenon is followed by linear evolution on a constant value that corresponds to the balanced state with diluted gradient.Nevertheless, some experimental results with long enough time intervals allow complete evolution of the process (polymer degradation stage is included here) and show unusually strong fluctuating behaviour.Experimental data of drug release, at short and long time scales, for polymeric microparticles (as polymeric matrices) are presented below.


### 3.1. Materials

The following materials were used: low-molecular-weight chitosan-CS, deacetylation degree 75–80% (Aldrich), type B gelatin-GEL (Aldrich), glutaraldehyde-GA (Aldrich)-25% aqueous solution, sodium tripolyphosphate-TPP (Sigma), Levofloxacin-LEV (Sigma), Tween 80 (Aldrich), and Span 80 (Aldrich).


### 3.2. Preparation of Microparticles by Ionic Gelation and Covalent Cross-Linking in O/W/O Emulsion

Microparticles were prepared using an original double cross-linking method of a CS-GEL mixture. Different weight ratios CS/GEL (in terms of amino groups of both polymers) were dissolved in acetic acid solution 2%, and then Tween 80 was added to make a 2% (w/w) surfactant in the solution. The mixture was magnetically stirred until the surfactant was completely dissolved. Two different solutions of 2% Span 80 in toluene were prepared according to O1/W (v/v) = 1/4 and O2/(O1 + W) = 4/1. The organic phase O1 was dripped within the aqueous polymer phase, W under homogenization with an Ultraturax device at 9000 rpm. The primary emulsion was transformed into a double one through dripping in the second organic phase O2, according to the same hydrodynamic regime. The emulsion was then gelled by slowly adding a TPP solution at a rate of 2 mL/min with continuous stir for extra 10 min.The suspension was then transferred to a round-bottom flask and mechanically stirred at 500 rpm. A certain amount of a saturated solution of GA in toluene was consecutively added and stirred for 60 min. The particles were separated by centrifugation (6000 rpm) and repeatedly washed with acetone and water in order to eliminate residual compounds. After hexane wash, the particles were dried at room temperature.


### 3.3. Preparation Parameters

A two-step solidification method was used. The first step, which has critical influence over the subsequent particle shape and size, included ionic cross-linking with TPP effect through phosphate bridges among amino functionalities in both types of polymeric chains. The GA covalent cross-linking (also taking place in NH2 groups) was performed with the purpose of stabilizing gel capsules. Our study analyzes the influence of the following cross-linking reaction parameters on the levofloxacin release kinetics:(i)
concentration of the ionic cross-linker,(ii)
ratio among amino functionalities of the two polymers and the ionic crosslinker,(iii)
polymer composition of the polymer mixture.Table1 shows the variable parameters in the preparation protocol that have been grouped according to the variable parameters.Table 1
The variable parameters in the preparation protocol.
Sample code
CS/GEL (w/w)
Conc. of TPP sol. (%)
NH2/TPP (mols/mols)
C3
1/1
1
2.4/1
C1
5
C4
10
C2
15
C7
1/1
5%
1.17/1
C1
2.4/1
C5
4.8/1
C6
11.7/1
C5
1/1
5%
4.8/1
C8
1/0
C9
3/1


### 3.4. Levofloxacin Release Kinetics

#### 3.4.1. Levofloxacin Release Kinetics at Small-Time Scales

If the experimental time scale is of minutes order, the evolution of the released drug concentration will be described by Peppas law. In this case, the correlated factor ranges between 0.8413 and 0.9983. Experimental and Peppas curves can be observed in Figure4 (the Peppas parameters and the correlation coefficient R2 for each sample are given in Table 2). The plots group according to variable preparation parameters (for a better observation of the first points, time scale is 500 min, although the fitting was made on the points up to 1440 min (one day)). Relative errors range between 1% and 5%, with no important influence on release kinetic evolution.Table 2
The variable parameters of the preparation protocol.
Sample code
Conc. of TPP sol. (%)
k
n
Fractal dimension
R
2
C3
1
0.0142
0.3225
6.20
0.9884
C1
5
0.0176
0.2879
6.95
0.9645
C4
10
0.0069
0.4326
4.62
0.9932
C2
15
0.0261
0.2357
8.49
0.8413
Sample code
NH2/TPP (mols/mols)
k
n
Fractal dimension
R
2
C7
1.17/1
0.0032
0.5131
3.90
0.9983
C1
2.4/1
0.0176
0.2879
6.95
0.9645
C5
4.8/1
0.0431
0.128
15.63
0.9787
C6
11.7/1
0.0116
0.3529
5.67
0.9879
Sample code
CS/GEL (w/w)
k
n
Fractal dimension
R
2
C5
1/1
0.0431
0.128
15.63
0.9787
C8
1/0
0.0256
0.1854
10.79
0.9539
C9
3/1
0.0201
0.2994
6.68
0.9895Levofloxacin release kinetics (experimental and Peppas fitting), at time scale of hours order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)Previous works have shown the form dependence [38, 39] between the value of parameter n in Peppas equation (considered as short-time approximation of Weibull function) and the fractal dimension of the drug particle during the release process (Df)
(3.1)n=2Df.Thus, according to experimental data, the following values were obtained in Table2.One first observation refers to the proportional dependence among experimental variable, on one hand, and Peppas parameters, on the other, in the particular case of the third sample group, that is,n increases with the chitosan/gelatin ratio. This proves to be experimentally useful if we want to obtain a Fickian diffusion. At the same time, the concentration of the released drug proves to be very low. This could be explained by drug crystallization inside the microparticle and the dependence of its release (dissolution followed by diffusion) on polymer degradation.In our opinion, the value of the fractal dimension is important as long as its values are unusually high and indicate that either fractal dimension must be considered as function of structure “classes,” or drug release processes (implicitly drug particle trajectories) have high degrees of complexity and nonlinearities, implying many freedom degrees in the phase space [42].This analysis (small concentration of the released drug and high fractal dimensions) made us continue the experiment until the system reached a stationary state.


#### 3.4.2. Levofloxacin Release Kinetics at Large-Time Scales

The experiments at large-time scales (of days order) revealed unusual behavior characterized by large variations. The release kinetics of levofloxacin is plotted in Figure5. The relative errors range between 1% and 5%, and, for better visualization, the error bars are plotted in Figure 6.Experimental release kinetics of levofloxacin, at time scale of days order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)The best correlations among experimental and theoretical curves (blue line—experimental curve, red line—theoretical curve).
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Experiments have been performed for 28 days, the concentration of the released drug being measured daily, at the same hour. The general characteristic of the above kinetics refers to strong variations of concentration in time, approximately at the same moment.In the following section, we will explain the evolution of these systems through the theoretical model (developed in Section2) based on fractal approximation of motion.


### 3.5. The Correspondence between Theoretical Model and Experimental Results

In what follows, we identify the fieldΦ from relation (2.57) with normalized concentration field of the released drug from microparticles.For best correlation between experimental data and the theoretical model (for each sample), we used a planar intersection of the graph in Figure3, where the two variables are y=(ξ¯-τ¯u)/2 and x=s. With these variables, (2.57) becomes(3.2)ϕ1(x,y)=2a(E(x)K(x)-1)+2a·cn2[ax(y+ξ0);x].Thus, in order to find the one-dimensional equation for a planar intersection, perpendicular to planexOy, we used y=mx+n (linear function equation), where m and n are two parameters. This equation is transformed into a parametric equation by means of the following substitutions:(3.3a)x=lm2+1,(3.3b)y=n+mlm2+1,in (3.2).Afterwards, we obtain one-dimensional function(3.4)ϕ2(t,m,n)=ϕ1(tm2+1,n+mtm2+1).The highest value of the correlation coefficient (for two vectors: one obtained from this very function, the other from experimental data) for different values ofm and n (in the particular experimental case) will represent the best approximation of experimental data with the theoretical model.Our goal was to find the right correlation coefficient which should be higher than0.6-0.7, in order to demonstrate the relevance of the model we had in view. Figure 6 shows experimental and theoretical curves that were obtained through our method, where R2 represents the correlation coefficient and η a normalized variable which is simultaneously dependent on normalized time and on nonlinear degree of the system (s parameter). Geometrically, η represents the congruent angle formed by the time axis and the vertical intersection plane.Parametersm and n of the planar intersections for the above theoretical curves are shown in Table 3.Table 3
Parameters of the planar intersections.
Sample code
Conc. of TPP sol. (%)
n
m
C3
1
3.342
9.17
C1
5
6.902
8.024
C4
10
8.125
13.486
C2
15
9.479
12.25
Sample code
NH2/TPP (mols/mols)
n
m
C7
1.17/1
7.322
4.297
C1
2.4/1
6.902
8.024
C5
4.8/1
2.414
8.665
C6
11.7/1
4.24
12.747
Sample code
CS/GEL (w/w)
n
m
C5
1/1
2.414
8.665
C8
1/0
8.303
5.941
C9
3/1
8.678
3.738We must mention that for each sample the fitting process was an independent one. The corresponding intersection plane that offers best correlation factors had to be identified by each sample in turn.A first observation refers to the correlation among plane and variable parameters (within experimental protocol) differ from Peppas small-time-scale fitting.We consider that this could be a starting point in establishing dependence among experimental parameters involved in the protocol. The purpose of this analysis is to obtain polymer matrices together with characteristics of release kinetics, taking into account that until now, this type of dependence had to pass through intermediary stages of the physical and chemical characterization of polymer matrices.The few experimental data could not sustain a general conclusion on the existing dependence among plane and experimental parameters, but this will be the purpose of a next paper.


## 3.1. Materials

The following materials were used: low-molecular-weight chitosan-CS, deacetylation degree 75–80% (Aldrich), type B gelatin-GEL (Aldrich), glutaraldehyde-GA (Aldrich)-25% aqueous solution, sodium tripolyphosphate-TPP (Sigma), Levofloxacin-LEV (Sigma), Tween 80 (Aldrich), and Span 80 (Aldrich).


## 3.2. Preparation of Microparticles by Ionic Gelation and Covalent Cross-Linking in O/W/O Emulsion

Microparticles were prepared using an original double cross-linking method of a CS-GEL mixture. Different weight ratios CS/GEL (in terms of amino groups of both polymers) were dissolved in acetic acid solution 2%, and then Tween 80 was added to make a 2% (w/w) surfactant in the solution. The mixture was magnetically stirred until the surfactant was completely dissolved. Two different solutions of 2% Span 80 in toluene were prepared according to O1/W (v/v) = 1/4 and O2/(O1 + W) = 4/1. The organic phase O1 was dripped within the aqueous polymer phase, W under homogenization with an Ultraturax device at 9000 rpm. The primary emulsion was transformed into a double one through dripping in the second organic phase O2, according to the same hydrodynamic regime. The emulsion was then gelled by slowly adding a TPP solution at a rate of 2 mL/min with continuous stir for extra 10 min.The suspension was then transferred to a round-bottom flask and mechanically stirred at 500 rpm. A certain amount of a saturated solution of GA in toluene was consecutively added and stirred for 60 min. The particles were separated by centrifugation (6000 rpm) and repeatedly washed with acetone and water in order to eliminate residual compounds. After hexane wash, the particles were dried at room temperature.


## 3.3. Preparation Parameters

A two-step solidification method was used. The first step, which has critical influence over the subsequent particle shape and size, included ionic cross-linking with TPP effect through phosphate bridges among amino functionalities in both types of polymeric chains. The GA covalent cross-linking (also taking place in NH2 groups) was performed with the purpose of stabilizing gel capsules. Our study analyzes the influence of the following cross-linking reaction parameters on the levofloxacin release kinetics:(i)
concentration of the ionic cross-linker,(ii)
ratio among amino functionalities of the two polymers and the ionic crosslinker,(iii)
polymer composition of the polymer mixture.Table1 shows the variable parameters in the preparation protocol that have been grouped according to the variable parameters.Table 1
The variable parameters in the preparation protocol.
Sample code
CS/GEL (w/w)
Conc. of TPP sol. (%)
NH2/TPP (mols/mols)
C3
1/1
1
2.4/1
C1
5
C4
10
C2
15
C7
1/1
5%
1.17/1
C1
2.4/1
C5
4.8/1
C6
11.7/1
C5
1/1
5%
4.8/1
C8
1/0
C9
3/1


## 3.4. Levofloxacin Release Kinetics

### 3.4.1. Levofloxacin Release Kinetics at Small-Time Scales

If the experimental time scale is of minutes order, the evolution of the released drug concentration will be described by Peppas law. In this case, the correlated factor ranges between 0.8413 and 0.9983. Experimental and Peppas curves can be observed in Figure4 (the Peppas parameters and the correlation coefficient R2 for each sample are given in Table 2). The plots group according to variable preparation parameters (for a better observation of the first points, time scale is 500 min, although the fitting was made on the points up to 1440 min (one day)). Relative errors range between 1% and 5%, with no important influence on release kinetic evolution.Table 2
The variable parameters of the preparation protocol.
Sample code
Conc. of TPP sol. (%)
k
n
Fractal dimension
R
2
C3
1
0.0142
0.3225
6.20
0.9884
C1
5
0.0176
0.2879
6.95
0.9645
C4
10
0.0069
0.4326
4.62
0.9932
C2
15
0.0261
0.2357
8.49
0.8413
Sample code
NH2/TPP (mols/mols)
k
n
Fractal dimension
R
2
C7
1.17/1
0.0032
0.5131
3.90
0.9983
C1
2.4/1
0.0176
0.2879
6.95
0.9645
C5
4.8/1
0.0431
0.128
15.63
0.9787
C6
11.7/1
0.0116
0.3529
5.67
0.9879
Sample code
CS/GEL (w/w)
k
n
Fractal dimension
R
2
C5
1/1
0.0431
0.128
15.63
0.9787
C8
1/0
0.0256
0.1854
10.79
0.9539
C9
3/1
0.0201
0.2994
6.68
0.9895Levofloxacin release kinetics (experimental and Peppas fitting), at time scale of hours order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)Previous works have shown the form dependence [38, 39] between the value of parameter n in Peppas equation (considered as short-time approximation of Weibull function) and the fractal dimension of the drug particle during the release process (Df)
(3.1)n=2Df.Thus, according to experimental data, the following values were obtained in Table2.One first observation refers to the proportional dependence among experimental variable, on one hand, and Peppas parameters, on the other, in the particular case of the third sample group, that is,n increases with the chitosan/gelatin ratio. This proves to be experimentally useful if we want to obtain a Fickian diffusion. At the same time, the concentration of the released drug proves to be very low. This could be explained by drug crystallization inside the microparticle and the dependence of its release (dissolution followed by diffusion) on polymer degradation.In our opinion, the value of the fractal dimension is important as long as its values are unusually high and indicate that either fractal dimension must be considered as function of structure “classes,” or drug release processes (implicitly drug particle trajectories) have high degrees of complexity and nonlinearities, implying many freedom degrees in the phase space [42].This analysis (small concentration of the released drug and high fractal dimensions) made us continue the experiment until the system reached a stationary state.


### 3.4.2. Levofloxacin Release Kinetics at Large-Time Scales

The experiments at large-time scales (of days order) revealed unusual behavior characterized by large variations. The release kinetics of levofloxacin is plotted in Figure5. The relative errors range between 1% and 5%, and, for better visualization, the error bars are plotted in Figure 6.Experimental release kinetics of levofloxacin, at time scale of days order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)The best correlations among experimental and theoretical curves (blue line—experimental curve, red line—theoretical curve).
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Experiments have been performed for 28 days, the concentration of the released drug being measured daily, at the same hour. The general characteristic of the above kinetics refers to strong variations of concentration in time, approximately at the same moment.In the following section, we will explain the evolution of these systems through the theoretical model (developed in Section2) based on fractal approximation of motion.


## 3.4.1. Levofloxacin Release Kinetics at Small-Time Scales

If the experimental time scale is of minutes order, the evolution of the released drug concentration will be described by Peppas law. In this case, the correlated factor ranges between 0.8413 and 0.9983. Experimental and Peppas curves can be observed in Figure4 (the Peppas parameters and the correlation coefficient R2 for each sample are given in Table 2). The plots group according to variable preparation parameters (for a better observation of the first points, time scale is 500 min, although the fitting was made on the points up to 1440 min (one day)). Relative errors range between 1% and 5%, with no important influence on release kinetic evolution.Table 2
The variable parameters of the preparation protocol.
Sample code
Conc. of TPP sol. (%)
k
n
Fractal dimension
R
2
C3
1
0.0142
0.3225
6.20
0.9884
C1
5
0.0176
0.2879
6.95
0.9645
C4
10
0.0069
0.4326
4.62
0.9932
C2
15
0.0261
0.2357
8.49
0.8413
Sample code
NH2/TPP (mols/mols)
k
n
Fractal dimension
R
2
C7
1.17/1
0.0032
0.5131
3.90
0.9983
C1
2.4/1
0.0176
0.2879
6.95
0.9645
C5
4.8/1
0.0431
0.128
15.63
0.9787
C6
11.7/1
0.0116
0.3529
5.67
0.9879
Sample code
CS/GEL (w/w)
k
n
Fractal dimension
R
2
C5
1/1
0.0431
0.128
15.63
0.9787
C8
1/0
0.0256
0.1854
10.79
0.9539
C9
3/1
0.0201
0.2994
6.68
0.9895Levofloxacin release kinetics (experimental and Peppas fitting), at time scale of hours order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)Previous works have shown the form dependence [38, 39] between the value of parameter n in Peppas equation (considered as short-time approximation of Weibull function) and the fractal dimension of the drug particle during the release process (Df)
(3.1)n=2Df.Thus, according to experimental data, the following values were obtained in Table2.One first observation refers to the proportional dependence among experimental variable, on one hand, and Peppas parameters, on the other, in the particular case of the third sample group, that is,n increases with the chitosan/gelatin ratio. This proves to be experimentally useful if we want to obtain a Fickian diffusion. At the same time, the concentration of the released drug proves to be very low. This could be explained by drug crystallization inside the microparticle and the dependence of its release (dissolution followed by diffusion) on polymer degradation.In our opinion, the value of the fractal dimension is important as long as its values are unusually high and indicate that either fractal dimension must be considered as function of structure “classes,” or drug release processes (implicitly drug particle trajectories) have high degrees of complexity and nonlinearities, implying many freedom degrees in the phase space [42].This analysis (small concentration of the released drug and high fractal dimensions) made us continue the experiment until the system reached a stationary state.


## 3.4.2. Levofloxacin Release Kinetics at Large-Time Scales

The experiments at large-time scales (of days order) revealed unusual behavior characterized by large variations. The release kinetics of levofloxacin is plotted in Figure5. The relative errors range between 1% and 5%, and, for better visualization, the error bars are plotted in Figure 6.Experimental release kinetics of levofloxacin, at time scale of days order, at different concentration of TPP solution (a), NH2: TPP mol ratios (b), and NH2 mols (c).
(a)
(b)
(c)The best correlations among experimental and theoretical curves (blue line—experimental curve, red line—theoretical curve).
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Experiments have been performed for 28 days, the concentration of the released drug being measured daily, at the same hour. The general characteristic of the above kinetics refers to strong variations of concentration in time, approximately at the same moment.In the following section, we will explain the evolution of these systems through the theoretical model (developed in Section2) based on fractal approximation of motion.


## 3.5. The Correspondence between Theoretical Model and Experimental Results

In what follows, we identify the fieldΦ from relation (2.57) with normalized concentration field of the released drug from microparticles.For best correlation between experimental data and the theoretical model (for each sample), we used a planar intersection of the graph in Figure3, where the two variables are y=(ξ¯-τ¯u)/2 and x=s. With these variables, (2.57) becomes(3.2)ϕ1(x,y)=2a(E(x)K(x)-1)+2a·cn2[ax(y+ξ0);x].Thus, in order to find the one-dimensional equation for a planar intersection, perpendicular to planexOy, we used y=mx+n (linear function equation), where m and n are two parameters. This equation is transformed into a parametric equation by means of the following substitutions:(3.3a)x=lm2+1,(3.3b)y=n+mlm2+1,in (3.2).Afterwards, we obtain one-dimensional function(3.4)ϕ2(t,m,n)=ϕ1(tm2+1,n+mtm2+1).The highest value of the correlation coefficient (for two vectors: one obtained from this very function, the other from experimental data) for different values ofm and n (in the particular experimental case) will represent the best approximation of experimental data with the theoretical model.Our goal was to find the right correlation coefficient which should be higher than0.6-0.7, in order to demonstrate the relevance of the model we had in view. Figure 6 shows experimental and theoretical curves that were obtained through our method, where R2 represents the correlation coefficient and η a normalized variable which is simultaneously dependent on normalized time and on nonlinear degree of the system (s parameter). Geometrically, η represents the congruent angle formed by the time axis and the vertical intersection plane.Parametersm and n of the planar intersections for the above theoretical curves are shown in Table 3.Table 3
Parameters of the planar intersections.
Sample code
Conc. of TPP sol. (%)
n
m
C3
1
3.342
9.17
C1
5
6.902
8.024
C4
10
8.125
13.486
C2
15
9.479
12.25
Sample code
NH2/TPP (mols/mols)
n
m
C7
1.17/1
7.322
4.297
C1
2.4/1
6.902
8.024
C5
4.8/1
2.414
8.665
C6
11.7/1
4.24
12.747
Sample code
CS/GEL (w/w)
n
m
C5
1/1
2.414
8.665
C8
1/0
8.303
5.941
C9
3/1
8.678
3.738We must mention that for each sample the fitting process was an independent one. The corresponding intersection plane that offers best correlation factors had to be identified by each sample in turn.A first observation refers to the correlation among plane and variable parameters (within experimental protocol) differ from Peppas small-time-scale fitting.We consider that this could be a starting point in establishing dependence among experimental parameters involved in the protocol. The purpose of this analysis is to obtain polymer matrices together with characteristics of release kinetics, taking into account that until now, this type of dependence had to pass through intermediary stages of the physical and chemical characterization of polymer matrices.The few experimental data could not sustain a general conclusion on the existing dependence among plane and experimental parameters, but this will be the purpose of a next paper.


## 4. Conclusions

If the particle moves on fractal curves, a new model for drug release mechanism from polymer matrix (namely, polymeric particles) is obtained. This model offers new alternatives for the theoretical study of drug release process (on large time scale) in the presence of all phenomena and considering a highly complex and implicitly nonlinear system. Consequently, the concentration field has cnoidal oscillation modes, generated by similar trajectories of drug particles. This means that the drug particle ensemble (at time large scale) works in a network of nonlinear oscillators, with oscillations around release boundary. Moreover, the normalized concentration field simultaneously depends on normalized time nonlinear system (throughs parameter).The fitting procedure among experimental and theoretical curves revealed the existing correlation of some characteristics of the release kinetics (the parameters of the intersection plane) with variable experimental parameters.At the same time, we consider that this could be a starting point in establishing dependence among experimental parameters, taking into account that until now, this type of dependence had to pass through intermediary stages of physical and chemical characteristics of polymer matrices.



---

*Source: 653720-2012-09-16.xml*