title: Accelerating solutions to the diffusion equation
authors: Felipe A. Asenjo Sergio A. Hojman
text: 
  
  
  
  
  
  felipe.asenjo@uai.cl
  Facultad de Ingeniería y Ciencias,
  Universidad Adolfo Ibáñez, Santiago 7491169, Chile.
  
  sergio.hojman@uai.cl
  Departamento de Ciencias, Facultad de Artes Liberales,
  Universidad Adolfo Ibáñez, Santiago 7491169, Chile.
  Departamento de Física, Facultad de Ciencias, Universidad de Chile,
  Santiago 7800003, Chile.
  Centro de Recursos Educativos Avanzados,
  CREA, Santiago 7500018, Chile.
  
  
  
  We report accelerating diffusive solutions to the diffusion equation with a constant diffusion tensor. The maximum values of the diffusion density evolve in an accelerating fashion described by Airy functions. We show the diffusive accelerating behavior for one–dimensional systems, as well as for a  general three–dimensional case.
   We also construct a modulated modified form of the diffusion solution  that retains the accelerating features.
  
  
  
  
  Accelerating solutions to  diffusion equation
      Sergio A. Hojman
      Aug. 2020
  =============================================
  
  
  
  § INTRODUCTION
  
   
  
  Diffusion is one of the most important problems in physics. Its applications include random walk and heat diffusion <cit.>, economy <cit.>
  and biophysics <cit.>, for instance.  In the most general fashion, the diffusion equation for a three–dimensional diffusing material with density  ϕ(t, x), 
  can be written as
  
      ∂ϕ/∂ t=∂/∂ x_i(D^ij∂ϕ/∂ x_j)  ,
  
  where D^ij=D^ij( x) is the symmetric diffusion tensor, which in principle depends on the position x where diffusion takes place. For Eq. (<ref>) we have used the Einstein summation convention for i, j=1, 2, 3.
  
  Along this work, we show that  accelerating diffusive solutions, described in terms of Airy functions,
  are possible for a general case of Eq. (<ref>) when the  diffusion tensor D^ij has  constant components. 
  The  solutions we present have accelerating properties such that the diffusion evolution of the local density
  has curved trajectories in a space-time plane. In that sense, the position x_M of the maxima for the diffusive density accelerate, i.e, d^2 x_M/dt^2≠ 0, in such a way that the density of the diffusing material has a diffusion which is different from ∼√(t), meaning that it differs from standard solutions. 
  
  The aim of this work is to present these accelerating solutions as robust phenomena, which have not been fully explored in diffusion.
  We start studying the accelerating behavior for one–dimensional systems, to later show the accelerating solution for a general case.
  
  
  
  
  
  § ONE–DIMENSIONAL  DIFFUSION
  
  
  
  We start simply by studying diffusion in a one–dimensional system with constant diffusive coefficient D, such as D^ij=D δ^i x along a preferred direction, say x.
  In this case
  a diffusing material with density  ϕ(t,x), has a dynamics given by 
  
      ∂ϕ/∂ t=D∂^2ϕ/∂ x^2  ,
  
  which is equivalent to the heat equation (where in such case ϕ represents the temperature).
  The standard well–known solution of Eq. (<ref>) has the Gaussian form 
  
      ϕ(t,x)=1/√(D t)exp(-x^2/4D t)  .
  
  The diffusive behavior of solution (<ref>), in a space-time phase space,
  is shown
  in Fig. <ref>(a). 
  
  
  
  
  
   In what follows, we present a different exact solution of Eq. (<ref>) with accelerating properties in the x-t phase space plane.
  Several solutions for the nonlinear version of the diffusion equations have been studied (see for example Refs. <cit.> among others), 
  and also in quantum systems <cit.>. However, the simplest linear model described by Eq. (<ref>) has a solution that has not been fully explored, as yet. 
  This solution corresponds to an initial value condition for the temperature ϕ(0,x) that exhibits an accelerating form of diffusion as it evolves in time.
  
   Several other known solutions of diffusion have accelerating properties. 
  In the context of reaction–difusion equations, different features of initial conditions modify the diffusion behavior <cit.>. On the other hand, for instance,
  Gaussian–like solutions with multiple local maxima as initial conditions satisfying Eq. (<ref>), such as
  
      ϕ(t,x)   =   1/√(D t)exp(-(x-x_0)^2/4D t)
            +λ/√(D t)exp(-(x+x_0/λ)^2/4D t)  ,
   
   have accelerating properties. Here, x_0 and λ are  arbitrary. This solution has the usual diffusive mean square distance ⟨ x^2⟩=∫_-∞^∞ x^2ϕ  dx=4√(π)(1+λ)(t+x_0^2/2λ).
  This density is plotted in Fig. <ref>(b), where its maxima are shown in solid red lines, and the minima is displayed in dashed red line.  The trayectories of the maxima approach to x=0 for longer times, and they display non-constant acceleration that approaches 0 as time grows. 
  This also can be seen in Fig. (<ref>)(a).
  However this kind of behavior is different from the presented solution in this work, as we will show.
  
  
  
  
  
   §.§ Simple accelerating solution
  
  
  Let us start with the  following ansatz, containing the simplest form of acceleration y=x+β t^2 in a one–dimensional system,
  
      ϕ(t,x)=f(k  y)exp(ζ(t,x))  ,
  
  where f and ζ are functions to be determined, β is a constant and 
  k is an arbitrary constant with units of inverse length. A simple solution can be found by assuming that ζ(t,x)=α  t  x+μ(t), in terms of a constant α and a function μ that depends on time only. Using this ansatz in Eq. (<ref>), we find that function f is an Airy function, β=D^2 k^3, α=D k^3, and μ(t)=2 D^3 k^6 t^3/3. In this form, the ansatz (<ref>) acquires the form
  
      ϕ(t,x)=(k  x+ k^4 D^2 t^2 )exp(k^3 D  t  x+2/3k^6 D^3 t^3 )  ,
  
  which can be straightforwardly proved to be a solution of Eq. (<ref>).
  The diffusive density (<ref>) can be also obtained by the Wick rotation in time of the Berry–Balazs's nondiffracting wavepacket  solution <cit.>
  for non–relativistic quantum mechanics.
  
  Solution
  (<ref>) presents local diffusion with acceleration in the x-t plane. This can be explicitly seen in Fig. <ref>(c), where the density (<ref>) is plotted in terms of dimensionless variables k  x and k^2 D  t. 
  This diffusion oscillates (because of the Airy function) representing different local temperature gradients which can be set as initial conditions. In Fig. <ref>(c), 
  we explicitly show the evolution in time of the local maxima  (solid red lines) and local minima (dashed red lines) values of density (<ref>). 
  The maxima and minima have curved trayectories in x-t plane, displaying the  local acceleration features of this solution. 
  
  The dynamics of the local maxima or minima can be found by solving for x_M(t) the following equation
  
      (ξ)/dξ+Dk^2 t (ξ)=0  ,
  
  where ξ=k  x_M+ k^4 D^2 t^2. This equation is solved numerically and shown in Fig. <ref>(b). Cleary d x_M/dt≠ 0 in all cases [different to the Gaussian case of Fig. <ref>(a)]. 
  While the maximum maximorum lobe curves towards the positive values of x, the other maxima and minima curve to its negative values.
  This implies different type of acceleration for different maxima and minima. This can be seen in Fig. <ref>(b),
  where a_M=d^2 x_M/dt^2≠ 0 is shown for the maximum maximorum (solid red line) and the first minimum (dashed red line) of solution (<ref>).
   In particular, the main maxima lobe accelerate always with a_M < 0, decreasing in magnitude as the time grows. The other maxima and minima have a slightly growing negative acceleration. This explains the different path taken by the main lobe and the other ones.
   Compare all these results with diffusion density (<ref>) in  Fig. <ref>(a), where its maxima do not accelerate (they follow a straight line trajectory in the x-t plane), or with diffusion (<ref>) in Fig. <ref>(b), where the accelerating solutions approaches to x=0, changing sign in acceleration in order  
  to form a unique Gaussian–like structure with no acceleration.
  
  In addition to this, we can calculate the trajectories of the zero-levels of density (<ref>). These are the blue lines in Fig. <ref>(c). They display acceleration that can be found analytically, as the zeros of Airy function are known. In this way, for instance, the trajectory of the first zero-level (the first blue line) in Fig. <ref>(c) has a time dependent position x_z given by
  
      x_z(t)=-k^3 D^2 t^2+_z1/k  ,
  
  where _z1≈ -2.33811 is the value for the first zero of Airy function [(_z1)=0]. It presents a change in position with the acceleration a_z=d^2x_z/dt^2=-2 k^3 D^2. This acceleration is shown in Fig. <ref>(b) as the dot-dashed magenta line.
  The procedure to obtain the dynamics of other zero-level curves is straightforward. 
  
  
  
  
  
  
  
  The asymmetry in the accelerating bearing of local maxima and minima in Fig.<ref>(b) can be explained by their asymptotic 
  behavior. For the main lobe, the solution (<ref>) decays in ∞ as ∼|ξ|^-1/4exp(-|ξ|^3/2), which is faster than any other local maxima and minima that decay in -∞ as ∼|ξ|^-1/4  only. These two different kinds of behavior explain also the different local diffusion rates of the solution. We can show that different local solutions of (<ref>) have diffusion rates different to ∼√(t). Gaussian solutions (<ref>) and (<ref>) have diffusion rates (locally and globally) that tend to  √(t). In Fig. (<ref>) we use
  density (<ref>) to display √(t) ϕ
   evaluated in  the maximum maximorum (solid red line) and in the first minimum (dashed redline). This plot indicates that the accelerating solution diffuses in a different way that the Gaussian solution. Its maximum approaches the usual diffusion  ∼√(t). However, the other maxima and minima show a completely different dynamics, diffusing faster as ∼ t^2 as time grows (represented as black dotted lines). In general, all local maxima and minima in (-∞,0) have a ∼ t^2 diffuse behavior.
   This does not occur in Gaussian-like solutions, and it is a direct consequence of the accelerating properties of the Airy function.
  
  
  
  
  
  
  
  
  
   §.§ Generalized accelerating solution
  
  
  Despite of the interesting features of solution (<ref>), it is possible to construct a modulated/attenuated version of the above solution by an arbitrary parameter. Following  Refs. <cit.>, and performing  a Wick rotation in time, the modulated accelerating Airy diffusive solution has  the form
  
      ϕ(t,x)   =(k  x+ k^4 D^2 t^2+2η  k^2 D  t )
         ×exp(k^3 D  t  x+2/3k^6 D^3 t^3+η k x
               
      +2η k^4 D^2t^2 +η^2 k^2 D   t )  ,
  
  which satisfies the diffusion equation (<ref>) for an arbitrary dimensionless constant η. This solution  approaches to a  single profile as η increases, but retaining its accelerating features. This can be seen in Fig. <ref>(d), where the density (<ref>) is shown in terms of k x and k^2 D t, for
   η=1/2.
  
  
   The accelerating trajectory for the diffusion of the maxima of the density (solid red line) shows that the accelerating behavior is present but reduced.  The same occurs for the minima (dashed red lines) and the zero–level solutions (blue lines). In general, as η increases, the acceleration properties  are attenuated.
  
  
  Similar to the simplest case, the time dependent positions x_M of local maxima and minima of  density (<ref>) can be calculated numerically by solving the equation
  
      (ξ̃)/dξ̃+Dk^2 t (ξ̃)=0  ,
  
  for x_M, where now ξ̃=k  x_M+ k^4 D^2 t^2+2η k^2 Dt. Also, the different zero-levels have time dependent positions x_zi (with i=1,2,3,...) given by
  
      x_zi(t)=-k^3 D^2 t^2-2η k D t+_zi/k  ,
  
  where _zi 
  are the different zeros of Airy function. All zero-levels anew move with constant acceleration a_z=-2k^3 D^2.
  
  The accelerating properties of this modulated diffusive solution
  can be seen in Fig. <ref>(b), where the 
  acceleration of the maximum maximorum (solid blue line), of the first minima (dashed blue line), and of the zero-levels (magenta dot-dashed line) are shown. For the maximum maximorum, the acceleration is diminished, whereas for the first minima (and other maxima and minima) the accelerating properties remain the same. Furthermore, the modulated density solution (<ref>) diffuses slower than solution (<ref>), and different to Gaussian solutions. In Fig. <ref>,  we use
  density (<ref>) to calculate √(t) ϕ evaluated in its maximum maximorum (in blue solid line). This approaches to ∼√(t) as time grows. Similarly,  density (<ref>) evaluated in its first minima (blue dashed line) shows a diffusion 
  rate of the order ∼t^2 (in black dotted line) as time grows. 
  
  
  
  
  
  
  
  
  
  
  
  § THREE–DIMENSIONAL GENERAL ACCELERATING  DIFFUSION
  
  
  The above solution can be generalized to a general three-dimensional system, described by Eq. (<ref>) for a general constant diffusion tensor, with components D_ij≠ 0, for i,j=1,2,3.
  Employing
  the same mathematical tools from previous section, it is possible to show that an accelerating solution of Eq. (<ref>) is
  
      ϕ(t,x,y,z)   =   (k_1 x+k_2 y+k_3 z+ γ^2 t^2 )
            ×exp(γ t(k_1 x+k_2 y+k_3 z)+ 2/3γ^3 t^3 )  ,
  
  with γ=D_ijk_ik_j. We emphasize that this solution is for a general constant diffusion with diffusive density
  in every possible direction.
  Notice that  general solution (<ref>) is equivalent to one–dimensional solution (<ref>) by performing the change of variables x'=(k_1 x+k_2 y+k_3 z)/k,
  with k=(k_1^2+k_2^2+k_3^2)^1/2.
  
  
   On the other hand, in the particular case of a constant diffusion tensor with the form D_ij=D_iδ_ij (with a delta function δ_ij) the following accelerating density is a different  solution to Eq. (<ref>)
  
      ϕ(t,x,y,z)   =   (κ_1 x+ ϱ_1^2 t^2 )(κ_2 y+ ϱ_2^2 t^2 )
            ×(κ_3 z+ ϱ_3^2 t^2 )
            ×exp(t(κ_1ϱ_1 x+κ_2ϱ_2 y+κ_3ϱ_3 z)+ς t^3)  ,
  
  with ϱ_i=κ_i^2 D_i,
  and ς=2(ϱ_1^3+ϱ_2^3+ϱ_3^3)/3.
  
  Lastly, to construct a modulated accelerating version of the solutions above, the same scheme 
  followed in previous sections may be used.
  
  
  
  
  
  § DISCUSSION
  
  
  The accelerating diffusion processes presented above have
   interesting  features that depend on the initial conditions of the system. They are the Wick rotated analogues of the accelerating solutions that appear  in quantum mechanics in theory <cit.> and experiments <cit.>.  In this sense, this accelerating diffusion phenomena  can be harnessed  in new processes where now the way it diffuses can be directed in space by proper choices of the free parameters in the solution.
  
  
  Although other accelerating diffusion processes which propagate depending on initial conditions do exist, such as what happens in reaction-diffusion equations with 
  faster diffusion and long distance dispersal <cit.>, the locally accelerating system presented here satisfy an equation which does not present non-linearities. It belongs to the same linear diffusive phenomena as Gaussian densities.
  
  The above accelerating behavior corresponds to local properties of the solutions (<ref>)  and (<ref>). They display different kind of accelerations that involve different dynamics than the ones of Gaussian solutions. However, their global features are analogous to their Gaussian counterparts. For instance, the average position of densities (<ref>) and (<ref>) is ⟨ x⟩=∫_-∞^∞ x ϕ(x,t) dx=, that depends only on η (for η=0, then 
  ⟨ x⟩=0), and thus the global solution behaves in the same way as Gaussian ones do, with d⟨ x⟩/dt=0. 
  On the other hand, both mean square distance of densities 
  (<ref>) and (<ref>) result to be
  
      ⟨ x^2⟩∼ t  ,
   
  as Gaussian solutions. Therefore, the accelerating densities have an analogue global behavior to standard solutions of diffusion.
  
  
  Finally, we expect to have other possible accelerating solutions of this kind for the nonlinear version of the diffusion equation. These endeavors are left for future work.
  
  
  
  
  
  
  
  
  
  17
  
  crank J. Crank, The Mathematics of Diffusion (Oxford University Press, 1975).
  blackscholes F. Black and M. Scholes,  Jour. Political Economy 81, 637 (1973). 
  bioph D. B. Chang et al., Jour. Theo. Bio. 50,  285 (1975).
  burgers J. M. Burgers, The Nonlinear Diffusion Equation:
  Asymptotic Solutions and Statistical Problems (Springer Netherlands, 1974).
  king1 J. R. King, J. Phys. A: Math. Gen. 23, 3681 (1990).
  king2 J. R. King, J. Phys. A: Math. Gen. 24, 3213 (1991).
  chole F. M. Cholewinski and J. A. Reneke, Elec. Jour. Diff. Eqs., 2003, 1 (2003).
  metzler R. Metzler and J. Klafter, Europhys. Lett. 51, 492 (2000).
  chai Chai Hok Eab and S. C. Lim, J. Phys. A: Math. Theor. 45, 145001 (2012).
  said E. A. Saied and M. M. Hussein, J. Phys. A: Math. Gen. 27 4867 (1994).
  ander D. Anderson and M. Lisak, Phys. Rev. A 22, 2761 (1980).
  james J. M. Hill, A. J. Avagliano, M. P. Edwards,  IMA Journal of Applied Mathematics 48, 283 (1992).
  
  
  
  
  
  tsekov R. Tsekov, Phys. Scr. 83, 035004 (2011). 
  
  garnier J. Garnier, SIAM J. Math. Anal. 43, 1955 (2011).
  dstan D. Stan and J. L. Vázquez, SIAM  J. Math. Anal. 46, 3241 (2014).
  hamel F. Hamel and L. Roques, J. Diff. Equations 249, 1726 (2010).
  rking R. King and P. M. McCabe, Proc. R. Soc. A  459, 2529 (2003).
  
  berry M. V. Berry and N. L. Balazs, Am. J. Phys. 47, 264 (1979).
  sivi G. A. Siviloglou and D. N. Christodoulides, Opt. Lett. 32, 979 (2007).
  lekn J. Lekner, Eur. J. Phys. 30 L43 (2009).
  hojmanasenjo S. A. Hojman and F. A. Asenjo, Phys. Lett. A 384, 126913 (2020).
  Broky G. A. Siviloglou, J. Broky, A. Dogariu, and D. N. Christodoulides,
  Phys. Rev. Lett. 99, 213901  (2007).
  
  
  
  