Multiple fluctuations - improved version

This is a revision of stochasflow6.nb, but here a simple approximation for the single fluctuation enhancement factor A1 is found so that the need for separate consideration of various ranges of |τ0| is avoided. Also a factor 4/3 mistake in eq 19.74 is corrected here.

Previous formulas

g[x_, σ_] := 
      (1/(σ (Sqrt[ 2  π]))) Exp[-((x^2)/(2  (σ^2)))]

a2 = (4 (α^2) Δ + 
              (T^4) (((α + Δ))^2) + 
              (T^2) α (-2 Δ (1 + Δ) + α (-1 - 
                          4 Δ + (Δ^2))))/((T^2) (-1 + 
                  (T^2)) (((α + Δ))^2));

or the following alternative expression

Simplify[
      a2 == 1 - ((4 (α^2) Δ)/((T^2)  (((α + Δ))^2))) + (((((1 - α))^2) (Δ^2))/( ((T^2) - 1) (((α + Δ))^2)))]

b0[
          T_] = (((T^2) (1 - Δ) (α + Δ) (Sqrt[(1/((T^2) (((-1 + Δ))^2) + 4 Δ))]))/(Sqrt[((4 (α^2) Δ + (T^4) (((α + Δ))^2) + (T^2) α (-2 Δ (1 + Δ) + α (-1 - 4 Δ + (Δ^2))))/(-1 + (T^2)))]));

b1[
          T_] = -((T^2) β (((-1 + Δ))^2) Δ ((((1/((T^2) (((-1 + Δ))^2) + 4 Δ))))^(3/2)) ((T^2) Δ (-8 Δ + (T^2) (1 + 
                                10 Δ + (Δ^2))) + 
                    2 α Δ (-4 (-3 - 
                                2 Δ + (Δ^2)) - 
                          (T^2) (23 + 9 Δ - 
                                5 (Δ^2) + (Δ^3)) + (T^4) (10 - Δ - 
                                2 (Δ^2) + (Δ^3))) + (α^2) (-8 (-3 - 
                                2 Δ + (Δ^2)) - 
                          2 (T^2) (17 + 5 Δ - 
                                7 (Δ^2) + (Δ^3)) + (T^4) (10 - 5 Δ - 4 (Δ^2) + 
                                3 (Δ^3)))))/(2 (Sqrt[2 π]) (1 + Δ) (4 (α^2) Δ + (T^4) (((α + Δ))^2) + 
                  (T^2) α (-2 Δ (1 + Δ) + α (-1 - 
                              4 Δ + (Δ^2)))));

b2[
          T_] = -((T^2) (β^2) (((-1 + Δ))^3) Δ ((((1/((T^2) (((-1 + Δ))^2) + 4 Δ))))^(5/2)) √((1/(-1 + 
                            (T^2))) (4 (α^2) Δ + 
                          (T^4) (((α + Δ))^2) + 
                          (T^2) α (-2 Δ (1 + Δ) + α (-1 - 
                                      4 Δ + (Δ^2))))) (64 (α^2) Δ (((α + Δ))^2) (3 - Δ - 
                          2 (Δ^2) + 2 (Δ^3)) + 
                    32 (T^2) α (-4 (Δ^4) (3 - Δ + (Δ^2)) + α (Δ^2) (-12 - 55 Δ + 13 (Δ^2) - 
                                13 (Δ^3) - 
                                7 (Δ^4) + 
                                2 (Δ^5)) + 
                          2 (α^2) Δ (-12 - 
                                34 Δ + 
                                13 (Δ^2) - 
                                4 (Δ^3) - 
                                7 (Δ^4) + 
                                2 (Δ^5)) + (α^3) (-12 - 25 Δ + 
                                17 (Δ^2) + (Δ^3) - 
                                7 (Δ^4) + 
                                2 (Δ^5))) + 
                    8 (T^4) α (2 (Δ^3) (21 + 
                                87 Δ - 5 (Δ^2) + 
                                17 (Δ^3) - 
                                4 (Δ^4)) + α (Δ^2) (187 + 456 Δ - 147 (Δ^2) + 
                                140 (Δ^3) - 
                                7 (Δ^4) - 
                                10 (Δ^5) + (Δ^6)) + 2 (α^2) Δ (127 + 
                                197 Δ - 
                                134 (Δ^2) + 
                                91 (Δ^3) + 
                                8 (Δ^4) - 
                                10 (Δ^5) + (Δ^6)) + (α^3) (109 + 106 Δ - 
                                149 (Δ^2) + 
                                58 (Δ^3) + 
                                9 (Δ^4) - 
                                10 (Δ^5) + (Δ^6))) - 
                    4 (T^6) (24 (Δ^4) (1 + Δ) + α (Δ^3) (167 + 
                                314 Δ - 
                                28 (Δ^2) + 
                                56 (Δ^3) - 
                                23 (Δ^4) + 
                                2 (Δ^5)) + (α^2) (Δ^2) (418 + 581 Δ - 294 (Δ^2) + 
                                244 (Δ^3) - 
                                66 (Δ^4) + 
                                3 (Δ^5) + 
                                2 (Δ^6)) + (α^4) (148 + 39 Δ - 206 (Δ^2) + 
                                128 (Δ^3) - 
                                28 (Δ^4) + 
                                5 (Δ^5) + 
                                2 (Δ^6)) + (α^3) Δ (417 + 324 Δ - 460 (Δ^2) + 
                                328 (Δ^3) - 
                                77 (Δ^4) + 
                                4 (Δ^6))) + 
                    (T^8) (3 (Δ^3) (1 + 
                                35 Δ + 
                                35 (Δ^2) + (Δ^3)) + 4 (α^3) Δ (98 - 
                                22 Δ - 
                                97 (Δ^2) + 
                                107 (Δ^3) - 
                                49 (Δ^4) + 
                                11 (Δ^5)) + (α^4) (104 - 78 Δ - 94 (Δ^2) + 
                                159 (Δ^3) - 
                                91 (Δ^4) + 
                                27 (Δ^5) - 
                                3 (Δ^6)) + 
                          4 α (Δ^2) (3 + 
                                92 Δ + 
                                71 (Δ^2) - 
                                37 (Δ^3) + 
                                20 (Δ^4) - 
                                7 (Δ^5) + 
                                2 (Δ^6)) + 
                          2 (α^2) Δ (3 + 
                                277 Δ + 
                                89 (Δ^2) - 
                                227 (Δ^3) + 
                                170 (Δ^4) - 
                                62 (Δ^5) + 
                                14 (Δ^6)))))/(32 (1 + Δ) (α + Δ) (4 (α^2) Δ + (T^4) (((α + Δ))^2) + 
                    (T^2) α (-2 Δ (1 + Δ) + α (-1 - 
                                4 Δ + (Δ^2))))^2);

a2f[tt_] := a2 //. T -> tt

Z2f[T_, 
        Step_] := (((-(1/T) + T) (1 + Δ))/((Sqrt[2]) (β - β Δ))) /. 
        Step

u[x_, 
        n_, σ_] := (( (2^(-n/2)) (Sqrt[π]) ((((x/σ)))^n) g[x, σ])/Gamma)[((1 + n)/2)]

ups = {α -> 41. , β -> 
            0.780869, Δ -> 0.2};

downs = {α -> 41. , β -> 
            0.780869, Δ -> -0.333333};

Cu[z_, T_, 
        Step_] := (b0[T] u[z, 0, (1/(Sqrt[2 a2f[T]]))] + 
            b1[T] u[z, 1, (1/(Sqrt[2 a2f[T]]))] + 
            b2[T] u[z, 2, (1/(Sqrt[2 a2f[T]]))]) /. Step

<< Graphics`Graphics`

Names["Graphics`Graphics`*"]

Stir[x_] := (((x - 1))^(x - 1)) 
        Exp[1 - x] (Sqrt[2  π (x - 1)])

ΔL = {Δ, ((-2 Δ)/(1 + (Δ^2))), Δ};

am1[m_, i_] := 
      Simplify[(((T^2) (((α[m, i] + ΔL[[i]]))^2))/((T^2) (((α[m, i] + ΔL[[i]]))^2) - 4  α[m, i]^2 ΔL[[i]])) //. {T -> (Sqrt[α])[m + If[i == 3, 1, 0], If[i == 3, 1, i + 1]]}]

A1[m_] := (am1[m, 1]*am1[m, 2]*am1[m, 3])

θL := 
        Together[(t0/τ)0 {ϵ1, ϵ1 + ((1 - Δ)/2), ϵ1 + 1}];

αL[m_] = 
      Simplify[{(((m - 1)*θL[[3]] + θL[[1]] - t0)/((m - 1)*θL[[3]] - t0)), (((m - 1)*θL[[3]] + θL[[2]] - t0)/((m - 1)*θL[[3]] + θL[[1]] - t0)), ((m*θL[[3]] - t0)/((m - 1)*θL[[3]] + θL[[2]] - t0))}]

α[m_, i_] := Simplify[αL[m][[i]]]

A1[m_] = Together[A1[m] /. ϵ1 -> (1/2)]

Revision of the fluctuation enhancement factors

The behaviour of the single fluctuation factor with τ0 and the order m:

Plot3D[
      Evaluate[A1[
            m] /. {Δ -> 
              0.1}], {τ0, -0.1, -10}, {m, 1, 20}, 
      PlotRange -> {1, 1.05}, BoxRatios -> {1, 1, 1}]

To see it better, just look at the excess above 1:

Plot3D[
      Evaluate[A1[m] - 
            1 /. {Δ -> 
              0.1}], {τ0, -0.1, -10}, {m, 1, 20}, 
      PlotRange -> {0, Automatic}, 
      BoxRatios -> {1, 1, 1}]

Behaviour still looks simple. To really see what happens, it is necessary to take a large enough τ0 value that there is a substantial range of the integer m values below and above it. Take τ0 = 1000:

LogLinearPlot[
      Evaluate[{A1[m] - 
              1, (10^10) (m + τ0)} //. {Δ -> 
              0.1, τ0 -> -1000}], {m, 1, 40000}, 
      PlotRange -> {0, 0.00006}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

Notice now that in the lower range, A1 stays relatively constant, but beyond m = -τ0 it goes down like 1/m to 1. The behaviour is demonstrated by the following series expansions:

Asermi = Simplify[Series[A1[m] - 1, {m, ∞, 2}]]

Simplify[Series[A1[m] - 1, {τ0, -∞, 2}]]

The striking thing here is that in each of the two limits, the leading term has a coefficient that depends only on Δ, and in fact the two coefficients have a simple numerical ratio. Define the coefficient in the first as Q (Q1 in stochasflow6.nb, but here we only use the first order coefficient throughout), then the second is -3 Q/2 :

QQ[Δ_] = SeriesCoefficient[Asermi, 1]

As shown in a plot further below, using just the leading term in each limit as was done in stochasflow6.nb, gives a reasonable account of the qualitative behaviour, but overestimates A1 in the intermediate region near m = -τ0. To improve on that, look also at the expansion around that point:

Simplify[
      Series[A1[m] - 1, {m, -τ0, 2}, {τ0, -∞, 
          2}]]

Again, the leading order has a coefficient just a simple multiple of Q. 
By trial and error, the following rational expression has the same sort of expansions in terms of m, τ0 and (m+τ0):

Simplify[
      Series[((A m - B τ0)/(C m τ0 + (((m + τ0))^2))), {m, ∞, 2}]]

Simplify[
      Series[((A m - B τ0)/(C m τ0 + (((m + τ0))^2))), {τ0, -∞, 2}]]

Simplify[
      Series[((A m - B τ0)/(C m τ0 + (((m + τ0))^2))), {m, -τ0, 2}, {τ0, -∞, 2}]]

Obtain A,B and C in that order by putting leading orders in the expansions the same; then we have A=Q, B=(3/2) Q and (A+B)/C=-(3/5)Qwhich gives C=-(25/6). Substituting these values into the trial expression simplifies to the following  approximation to A1:

Simplify[((Q m - (3/2) Q τ0)/(-(25/6) m τ0 + (((m + τ0))^2)))]

((3 Q)/(3 m - 2 τ0))

Apx[m_, τ0_, Δ_] := 
      QQ[Δ] (3/(3 m - 2 τ0))

LogLinearPlot[
      Evaluate[{A1[m] - 1, Apx[m, τ0, Δ], 
            QQ[Δ]/m, -(3/2) 
              QQ[Δ]/τ0, (10^10) (m + τ0)} //. {Δ -> 0.1, τ0 -> -1000}], {m, 1, 40000}, 
      PlotRange -> {0, 0.00006}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

The approximation is far simpler than A1(m) itself  (especially noticeable that the dependence on Δ is factored out completely)  and yet captures the behaviour almost perfectly, far better than the two limitng expressions that were used in stochasflow6.nb and plotted as the turquoise and purple curves above.

Look how the approximation performs for low values of τ0:

Plot[
      Evaluate[{A1[m] - 1, 
            Apx[m, τ0, Δ]} //. {Δ -> 
              0.15, τ0 -> -10}], {m, 1, 50}, 
      PlotRange -> {0, Automatic}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

Plot[
      Evaluate[{A1[m] - 1, 
            Apx[m, τ0, Δ]} //. {Δ -> 
              0.15, τ0 -> -1}], {m, 1, 50}, 
      PlotRange -> {0, Automatic}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

Plot[
      Evaluate[{A1[m] - 1, 
            Apx[m, τ0, Δ]} //. {Δ -> 
              0.15, τ0 -> -0.01}], {m, 1, 50}, 
      PlotRange -> {0, Automatic}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

To summarise, the approximation works excellently for large τ0. i.e. small scale fluctuations. They are the most important ones, as large scale fluctuations only come into play when the scale of the experiment approaches the relatively large length Λ determinied by the Peclet number. For such large distances, fluctuations with small τ0 occurr, and these are still well represented by the approximation over most of the m range, only giving underestimates for small m values. Such underestimates imply a constant correction factor and do not affect the trend as m -> ∞.
. So continue by finding the enhancement factor for multiple fluctuations.

QQ[0.17]

Am[M_] = (∏_(m = 1)^M)(1 + ((3 Q)/(3 m - 2 τ0)))

Notice that while in the product expression M has to be an integer, and hence Am rises in discrete steps, when expressed in terms of Pchhammer and hence Gamma-functions, we can consider M a continuous variable (later to be replaced by continuous time or position variables). In so doing Am becomes a continuous function, can be seen as an interpolation of the discrete steps that physically result from the fact that there are discrete fluctuations encountered as the plume moves across the porous medium.

Plot[Am[M] /. {Q -> 0.1, τ0 -> -100}, {M, 1, 10000}, 
      PlotRange -> {1, Automatic}]

Am[1] /. {Q -> 0.1, τ0 -> -100}

Once more, Am rises fast in the beginning, then more slowly but does not appear to saturate. Investigate by taking large M limits.

For τ0 near 1, but M >> 1, this can be simplified in the following way:

Simplify[
      Am[M] /. {Pochhammer[a_, n_] -> 
            Stir[a + n]/Gamma[a], n_! -> Stir[n + 1]}]

But consider first the  more important case where |τ0| >> 1 and we can use the Stirling approximation in both numerator and denominator:

As = 
      Simplify[Am[
            M] /. {Pochhammer[a_, n_] -> 
              Stir[a + n]/Stir[a]}]

Combine and rearrange the factors to simplify the expression:

As1 = (((1 + ((3  Q)/(3  M - 2  τ0))))^((1/2) + 
                M - ((2 τ0)/3))) (((M + 
                  Q - ((2 τ0)/3)))^Q) (((1 - ((3  Q)/(2  τ0))))^(-(1/2) + ((2 τ0)/3))) (((Q - ((2 τ0)/3)))^(-Q));

As2 = (((1 + ((3  Q)/(3  M - 2  τ0))))^((1/2) + 
                M - ((2 τ0)/3))) (((1 + ((3 M)/(3 Q - 2 τ0))))^Q) (((1 - ((3  Q)/(2  τ0))))^(-(1/2) + ((2 τ0)/3)));

Apart[(((1 - ((3 Q)/(2 τ0))))^(-1)) (1 + ((3 Q)/(3 M - 2 τ0)))]

Together[((3 M Q)/((-M + Q) (-3 M + 2 τ0))) - ((3 M Q)/((-M + Q) (-3 Q + 2 τ0)))]

As3 = (((1 + ((3  Q)/(3  M - 2  τ0))))^M) (((1 + ((3 M)/(3 Q - 2 τ0))))^Q) (((1 - ((9 M Q)/((3 M - 2 τ0) (3 Q - 2 τ0)))))^((1/2) - ((2 τ0)/3)));

Check the algebra. Mathematica won't do it with fractional exponents, so help it along a little:

Simplify[(As == As3) /. τ0 -> -3/
            4, {M ∈ Integers, Q ∈ Integers}]

Plot AS3 as a function of M for typical Q=0.1, τ0=-1000:

LogLinearPlot[
      As3 /. {Q -> 0.1, τ0 -> -1000}, {M, 1, 10000}, 
      PlotRange -> {1, Automatic}]

It is easy to confirm that the factors with powers M and Q resepctively reduce to the limiting forms found in stochasflow6.nb for M<<|τ0| and M>>|τ0| respectively, while the remaining factors are close to 1 in both cases.

Accept the AS3 expression for the big |τ0| range:

Abigm[
        M_] := (((1 + ((3  Q)/(3  M - 2  τ0))))^M) (((1 + ((3 M)/(3 Q - 2 τ0))))^Q) (((1 - ((9 M Q)/((3 M - 2 τ0) (3 Q - 2 τ0)))))^((1/2) - ((2 τ0)/3)))

Abigm[0]

This limit, M=0, means that before the first fluctuation is encountered, the enhancement factor is 1 as it should be. On the logarithmic plot above it looks like this value is approached as an asymptote (gradient 0) while the linear plot shows a non-zero gradient at low M. Check that analytically:

Simplify[(∂_M) Abigm[M]] /. 
      M -> 0
Simplify[% /. {Q -> 0.1, τ0 -> -1000}]

So indeed the gradient is positive (although small for Q << |τ0|).

Dispersivity as function of traversal length

Now find the dispersivity, replacing M with a time variable as before and also replacing τ0 by the ratio of length scales (19.67a):

Amt[t_] = 
      Abigm[((-2 Λ t)/(3  L t0))] /. τ0 -> ((-Λ)/L)

varf[t_] := (t - t0) Amt[t]

Simplify[(∂_t) varf[t]]

Now that the time differentiation is done, it is simpler to switch directly to the traversal length Γ as the independent variable:

Simplify[(∂_t) 
          varf[t] //. {t -> -((Γ t0)/Λ)}]

Compare that with the Am expression in the same variables:

Amt[-((Γ t0)/Λ)]

Simplify[
      Amt[-((Γ t0)/Λ)]]

The 3 n-th power brackets in Am reappeared in the derivative. Split them off leaving the rest as a prefactor .

pref[Γ_] := (1/(3 L (Γ + Λ))) ((2^(-1)) (((((3 L Q + 2 (Γ + Λ))/(Γ + Λ))))^(-1)) (3 L (3 L Q + 4 (Γ + Λ)) + 
                4 (Γ + Λ) (3 L Q + 
                      2 (Γ + Λ)) Log[
                    1 + ((3 L Q)/(2 (Γ + Λ)))]))

Rewrite the argument of the log to allow further simplification:

pref[Γ_] := (Λ  (((((2 Λ (6 L Q + 3 Γ + 4 Λ))/((3 L Q + 2 Λ) (3 Γ + 4 Λ)))))^(-1)) (2 L ((((3 Γ + 4 Λ))^2) + 
                      3 L Q (3 Γ + 
                            5 Λ)) + (6 L Q + 
                      3 Γ + 
                      4 Λ) (3 (Γ^2) + 
                      7 Γ Λ + 
                      4 (Λ^2)) Log[
                    1 + ((6 L Q)/(3 Γ + 4 Λ))]))/(L (3 L Q + 
                2 Λ) (((3 Γ + 4 Λ))^2))

In most cases L << Δ and/or L<<Λ and also Q<<1. So unless L is really big, we can approximate the logarithm to get

apref[Γ_] = 
      Simplify[pref[Γ] //. {Log[1 + x_] -> x}]

This prefactor is generally close to 1, and see how well the approximation performs, even for L already as large as Λ:

LogLinearPlot[
      Evaluate[{pref[Γ], 
            apref[Γ]} /. {Λ -> 1000, 
            Q -> 0.1, L -> 1000}], {Γ, 2, 50000}, 
      PlotRange -> {0, Automatic}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

Here is a case where L = 10 Λ:

LogLinearPlot[
      Evaluate[{pref[Γ], 
            apref[Γ]} /. {Λ -> 1000, 
            Q -> 0.1, L -> 10000}], {Γ, 2, 50000}, 
      PlotRange -> {0, Automatic}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

It is only when L >> 10 Λ that the approximation breaks down  A further simplification is found as

Simplify[apref[Γ] - 1 - Q]

So we finally settle on the prefactor as the expression

pref[Γ_] := 
        1 + Q - ((3 L Q)/(6 L Q + 4 (Γ + Λ)));

To recap, the dispersivity factor is the product of pref with the 3 factors that make up Amt[-Γ t0/Λ]. Inspection of the factors reveal that they all have the structure that can be simplified by using either L or Λ to rescale all lengths. Usin L as the scaling length would simplify the exponents best but the next step will be to investigate the effect of combining different fluctuation lengths so it is better to keep the explicit appearance of L in the result. For this reason choose Λ as the scaling length and reduce the variables to just two scaled lengths: the fluctuation length ℓ = (L/Λ) and the traversal length \(TraditionalForm\`ℒ = Γ\/Λ\). In terms of these the final expression for the dispersivity enhancement factor is the product of the 3 factors:

Den0[ℒ_, ℓ_] := (1 + 
            Q - ((3 Q ℓ)/(3 Q ℓ + 4 ℒ + 4))) (((1 - ((6 Q ℓ ℒ)/((3 Q ℓ + 2 ) (2 ℒ + 2 )))))^((1/2) + ((2 )/(3 ℓ))))

Den1[ℒ_, ℓ_] := (((1 + ((3 Q ℓ)/(2  ℒ + 2))))^(((2  ℒ)/(3 ℓ))))

Den2[ℒ_, ℓ_]:=(((1 + ((2 ℒ )/(3 Q ℓ + 2))))^Q)

Den[ℒ_, ℓ_] := 
      Den0[ℒ, ℓ]*Den1[ℒ, ℓ]*
        Den2[ℒ, ℓ]

Here, the pref factor has been grouped together with the 3rd factor in Amt, because they both have a similar algebraic structure in that the additive term in the bracket is of order 1 as ℒ->∞, and at smaller ℒ values they behave oppositely so that their product shows the following simple behaviour:

LogLinearPlot[
      Evaluate[{Den0[ℒ, 0.0001], 
            Den0[ℒ, 0.1], Den0[ℒ, 1]} /. 
          Q -> 0.1], {ℒ, 0.001, 100}, 
      PlotRange -> {0.8, 1.6}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

The behaviour of the other individual terms and the final enhancement factor is as follows. 

In the graph above as well as those below, curves for a range of ℓ values are shown, but bear in mind that only values of ℒ >= ℓ are possible, so some curves need to be truncated at the lower end.

LogLinearPlot[
      Evaluate[{Den1[ℒ, 0.0001], 
            Den1[ℒ, 0.1], Den1[ℒ, 1], 
            Exp[Q]} /. Q -> 0.1], {ℒ, 0.001, 100}, 
      PlotRange -> {1, 1.6}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

LogLinearPlot[
      Evaluate[{Den2[ℒ, 0.0001], 
            Den2[ℒ, 0.1], Den2[ℒ, 1]} /. 
          Q -> 0.1], {ℒ, 0.001, 100}, 
      PlotRange -> {1, 1.6}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

LogLinearPlot[
      Evaluate[{Den[ℒ, 0.0001], Den[ℒ, 0.1], 
            Den[ℒ, 1]} /. Q -> 0.1], {ℒ, 
        0.001, 100}, PlotRange -> {1, 1.6}, 
      PlotStyle -> {Hue[0.0], Hue[0.3], Hue[0.5], Hue[0.8]}]

{Den0[0.01, 0.001], Den1[0.01, 0.001], Den2[0.01, 0.001]} /. 
      Q -> 0.1

Generally, Den1 and Den2 are quite insensitive to the fluctuation length, and Den0 shows similar slowly decreasing behaviour but shifted downwards as the fluctuation length increases. The overall effect is that long fluctuations still increase the dispersivity, but not as much as shorter ones. Values of ℓ >1 are in fact not realistic, in view of the large value that is expected for Λ in a real experiment. Generally, it is clear that the behaviour of Den is dominated for ℒ << 1 by the behaviour of Den1 but for ℒ > 1 by Den2 and essentially rises as (ℒ^Q) .