world-wide energy transition

Apply for credit in the bank

2010-12-31T20:26:00.001-08:00

Now, The Credit People will talk about how the terms apply for credit in the bank. Like when going to open a savings account at a bank, you will be asked to submit copies of identity such as ID cards, driver's license, or passport. You are also asked to fill out a form that contains the data contributing to opening a personal savings yourself. The goal is that Bank has the right information, so as to identify yourself as a legitimate and entitled to perform transactions from your account. That is if you want to save money in the bank. Now what if you want to borrow money from banks? Here, the bank as the lender of funds called the creditors and those who borrow funds from the bank called the debtor. You can continue read here for further info about Repair Your Bad Credit.

Requirements for applying for loans in banks are not as complicated as people expected. Even the condition is actually quite easy. But of course, there is more data to be furnished rather than when you open a savings account. This is only fair. Let alone the bank. You will of course itself was more cautious and not willing to lend money away to just anyone if you do not believe that your money will be returned. Different if you give it as a donation or gift. Well, to assess whether the prospective borrower deserves credit, the bank must obtain the correct and accurate information, such as the character of the debtor, the funds held today, the influence of current economic conditions on borrowers income, collateral was filed, and much more .

Approximately the same as you, the bank was in accepting credit proposal incoming implement prudential principles in lending money. This is required by laws regulating the banking sector, even in the whole world. Remember that every penny of which is channeled back into the community by the bank is owned by the community as well. For each customer funds in the bank, the bank promised to return it to customers at any time with interest. Therefore, banks always do various kinds of credit analysis to assess the credit worthiness will be given to prospective customers. Need credit repair service?

Anyone can apply for loans to banks as long as eligible. In general, bank debtors divided into two major categories, namely individual debtor and the debtor company (again, the debtor is the party that borrows money from a bank).

Understanding Recovery Factors

2010-10-23T10:17:00.000-07:00

A recent TOD post on reserve growth by Rembrandt Kopelaar motivated this analysis.

The recovery factor indicates how much oil that one can recover from the original estimate. This has important implications for the the ultimately recovery resources, and increases in recovery rate has implications for reserve growth.

First of all, we should acknowledge that we still have uncertainty as to the amount of original oil in place, so that the recovery factor has two factors of uncertainty.

The cumulative distribution of reservoir recovery factor typically looks like the following S-shaped curve. The fastest upslope indicates the region closest to the average recovery factor.

Figure 1: Recovery Factor cumulative distribution function (from)

To understand the spread in the recovery factors, one has to first realize that all reservoirs have different characteristics. Some are more difficult to extract from and others have easier recovery factors. One of the principle first-order effects has to do with the size of the reservoir: bigger reservoirs typically have better recovery factors and as one reservoir engineer mentioned on TOD

"Reserve growth tends to happen in bigger fields because thats where you get the most bang for your buck"

So if we make the simple assumption that cumulative recovery factors (RF) have Maximum Entropy uncertainty or dispersion for a given Size:

P(RF) = 1-exp (-k*RF/Size)

this makes sense as the recovery factor will extend for larger fields.

Add to the mix that reservoir Sizes go approximately as (see here):

Pr(Size)= 1/(1+Median/Size)

Then a simple reduction in these sets of equations (with the key insight that RF ranges between 0 and 1, i.e. between 0 and 100%) gives us

P(RF) = 1 - exp(-k*RF*RF/(1-RF)/Median)

the ratio Median/k indicates the fractional average recovery factor relative to the median field size.

A set of curves for various k/Median values below:

Figure 2: Recovery Factor distribution functions assuming maximum entropy

Rembrandt provided some recovery factor curves originally supplied by Laherrere, and I fit these to the Median/k fractions below.

Figure 3: Recovery factor curves from Rembrandt's TOD post,
alongside the recovery factor model described here.

Laherrere also provided curves for natural gas, where recovery factors turn out much higher.

Figure 4: Recovery Factor distribution functions for natural gas.
Note that the recovery factor is much higher than for oil.
(Note: I had to fix the typo in the graph x-axis naming)

It looks like this derivation has strong universality underlying it. This remains a very simple and parsimonious model as it has only one sliding parameter. The parameter Median/k works in a scale-free fashion because both numerator and denominator have dimensions of size. This means that one can't muck with it that much -- as recovery factors increase, the underlying uncertainty will remain and the curves in Figure 2 will simply slide to the right over time while adjusting their shape. This will essentially describe the future reserve growth we can expect; the uncertainty in the underlying recovery factors will remain and thus we should see the limitations in the smearing of the cumulative distributions.

To reverse the entropic dispersion of nature and thus to overcome the recovery factor inefficiency, we will certainly have to expend extra energy.

Bird Surveys

2010-10-20T18:23:00.000-07:00

This post either points out something pretty obvious or else it reveals something of practical benefit. You can judge for now.

I briefly made a reference to bird survey statistics when I wrote this post on econophysics and income modeling. I took a typical rank histogram of bird species abundance and fit it the best I could to a dispersive growth model, further described here. The generally observed trend follows that many species exist in the middle of abundance and relatively small numbers of species exist at each end of the spectrum -- few species exceedingly common (i.e. starling) and few species exceedingly rare (i.e. whooping crane). Since the bird data comes from a large area in North America, the best fit followed a meta-community growth model. The meta-community adjustment impacts the knee of the histogram curve and broadens the Preston plot, effectively smearing over geological ages that different species have had to adapt.

Figure 1: Preston plot (top) and
rank histogram (bottom) of relative bird species abundance

If we assume that the relative species abundance has a underlying model related to steady-state growth according to p(rate), where rate is the relative advantage for species reproduction and survival, then this should transitively might apply to disturbances to growth as well. Recently, I ran into a paper that actually tried to discern some universality in diverse growth papers, and it coincidentally used the bird survey data along with two economic measures of firm size and mutual fund size.

Schwarzkopf, Yonathan, and J. Doyne Farmer, “The Cause of Universality in Growth Fluctuations.” Santa Fe Institute, (April 2010).
Schwarzkopf, Yonathan, and J. Doyne Farmer. “Supporting Information — The Cause of Universality in Growth Fluctuations.” Santa Fe Institute.

I did the best I could with the figures in the paper but eventually went to the source, ftp://ftpext.usgs.gov/pub/er/md/laurel/BBS/DataFiles/, and used data from 1997 to 2009.

I applied the same abundance distribution as before and came up with the fit below (see blue and red curves below, data and model respectively). That provided a sanity check, but Schwarzkopf and Farmer indicated that the year-to-year relative growth fluctuations should also obey some fundamental behavior through the distribution of this metric:

RelativeGrowth(Year) = n(Year+1) / n(Year)

Sure enough, and for whatever reason, the "growth" in the surveyed data does show as much richness as the steady state averaged abundance distribution. The relative growth in terms of a fractional yearly change sits alongside the relative abundance curve below (in green). Notice right off the bat that the distribution of fractional changes drops off much more rapidly.

Figure 2 : The red meta-model curve smears the median from 200 to 60000

I believe that this has a simple explanation having to do with Poisson counting statistics. When estimating fractional yearly growth, we consider that the rarer bird species having the lowest abundance will contribute most strongly to fluctuation noise on year-to-year survey data. Values flipping from 1 to 2 will lead to 100% growth rates for example. (We have to ignore movements from 1 to 0 and 0 to 1 as these growths become infinite.

I devised a simple algorithm that takes two extreme values (R greater than 1 and R less than 1 ) and the steady state abundance N for each species. The lower bound of:

R1 = R * (1-sqrt(2/N))/(1+sqrt(2/N))

and the upper bound becomes:

R2 = R * (1+sqrt(2/N))/(1-sqrt(2/N))

The term 1.4/sqrt(N) derives from Poisson counting statistics in that the relative changes become inversely related to the size of the sample. We double count in this case because we don't know whether the direction will go up or down, relative to R, a number close to unity.

(This has much similarity to the model I just used in understanding language adoption. Small numbers of adopters experience suppressing fluctuations as 1/sqrt(N))

Expanding on the scale, the results of this algorithm are shown in Figure 3.

Figure 3 : Model of yearly growth fluctuation in terms of a cumulative distribution function

Placing it in terms of a binned probability density function, the results look like the following plot. Note the high counts match closely the data simply because the 1/sqrt(N) is relatively small. Away from these points, you can see the general trend develop even though the data is (understandably) obscured by the same counting noise.

Figure 4 : The probability density function of yearly growth fluctuations.

As an essential argument to take home, consider that a counting statistics argument probably accounts for the yearly growth fluctuations observed. Before you make any other assertions, you likely have to remove this source of noise. Looking at Figure 3 & 4, you can potentially see a slight bias toward positive growth for certain lower abundance species. This comes at the expense of lower decline elsewhere, except for some strong declines in several other low abundance species. This may indicate the natural ebb and flow of attrition and recovery in species populations, with some of these undergoing strong declines. I haven't done this but it makes sense to identify the species or sets of species associated with these fluctuations.

Two puzzling points also stick out. For one, I don't understand why Schwarzkopf and Farmer didn't immediately discern this effect. Their underlying rationale may have some of the same elements but it gets obscured by their complicated explanation. They do use a resampling technique (on 40+ years worth of data) but I didn't see much of a reference to conventional counting statistics, only the usual hand-wavy Levy flight arguments. They did find a power law of around-0.3 instead of the -0.5 we used for Poisson, so they may generate something equivalent to Poisson by drawing from a similar Levy distribution. Overall I find this violates Occam's razor, at least for this set of bird data .

Secondly, it seems that these differential growth curves have real significance in financial applications. More of the automated transactions look for short duration movements and I would think that ignoring counting statistics could lead the computers astray.

Epilogue

As an aside, when I first pulled the data off the USGS server, I didn't look closely at the data sets. It turns out that the years 1994,1995,1996 were included in the data but appeared to have much poorer sampling statistics. From 1994 to 1996, the samples got progressively larger but I didn't realize this when I first collected and processed the data.

Figure 6 : CDF of larger data sample.
Note the strange hitch in the data growth fluctuation curve.

At the time, I figured that the slope had a simple explanation related to uncertainties in the surveying practice. I also saw some similarities to the uncertainties in stock market returns that I blogged about recently in an econophysics posting.

Say the survey delta time has a probability distribution with average time -- the T most likely related to the time between surveys:

pt(time) = (1/T)*exp(-time/T)

then we also assume that a surveyor tries to collect a certain amount of data, x, during the duration of the survey. We could characterize this as a mean, X, or some uniform interval. We don't have any knowledge of higher order moments to we just apply the Maximum Entropy Principle

px(x) = (1/X)*exp(-x/X)

The ratio between these two establishes the relative rate of growth, rate = X/T. We can derive the following cumulative quite easily:

P(rate) = T*rate/(T*rate +X)

The yearly growth rate fluctuations of course turn out as the second derivative of this function. We take one derivative to convert :

dp(rate)/drate = 2*T/X/(1+rate*T/X)^{^3}

On a cumulative plot as in Figure 6, this shows a power-law of order 2 (see the orange curve). Near the knee of the curve, it looks a bit sharper. If we use a uniform distribution of px(x) up to some maximum sample interval, then it matches the knee better (see the dashed curve).

So the simple theory says that much of the observed yearly fluctuation may arise simply due to sampling variations during the surveying interval. Plotting as a binned probability density function, the contrast shows up more clearly in Figure 7. In both cases is fit to X/T = 60. This number is bigger than unity because it looks like every year, the number of samples increases (I also did not divide by 15, the number of years in the survey).

But of course, the reason this maximum entropy model works as well as it does came about from real variation in the sampling techniques. Those years from 1994 to 1996 placed enough uncertainty and thus variance in the growth rates to completely smear the yearly growth fluctuation distribution.

Figure 7 : PDF of larger sample which had sampling variations.
Note that this has a much higher width than Figure 4.

Only in retrospect when I was trying to rationalize why a sampling variation this large would occur in a seemingly standardized yearly survey, did I find the real source of this variation. Clearly, the use of the Maximum Entropy Principle explains a lot, but you still may have to dig out the sources of the uncertainty.

Can we understand the statistics of something as straightforward as a bird survey? Probably, but as you can see, we have to go at it from a different angle than that typically recommended. I will keep an eye out if it has more widespread applicability; for now it obviously requires countable discrete entities.

Tower of Babel, How languages diversify

2010-10-16T16:15:00.000-07:00

One pattern that has evaded linguists and cognitive scientists for some time relates to the quantitative distribution in human language diversity. Much like how plant and animal species diversify in a specific pattern, with very few species dominating within an ecosystem and relatively few species exceedingly rare, the same thing happens with natural languages. You find a few languages spoken by many people, and very few spoken seldomly,with the largest number occupying the middle.

Consider a simple model of language growth whereby adoption of languages occur over time by dispersion. The cumulative probability distribution for the number of languages is

P(n) = 1/(1+1/g(n))

This form derives from the application of the maximum entropy principle to any random variate where one only knows the mean in the growth rate and an assumed mean in the saturation level. I refer to this as entropic dispersion and have used this many applications before so I no longer feel a need to rederive this term every time I bring it up.

The key to applying entropic dispersion is in understanding the growth term g(n). In many cases n will grow linearly with time so the result will assume a hyperbolic shape. In another case, an exponential growth brought up by technology advances will result in a logistic sigmoid distribution. Neither of these likely explains the language adoption growth curve.

Intuitively one imagines that language adoption occurs in fits and starts. Initially a small group of people (at least two for arguments sake) have to convince other people on the utility of the language. But a natural fluctuation arises with small numbers as key proponents of the language will leave the picture and the growth of the language will only sustain itself when enough adopters come along and the law of large numbers starts to take hold. A real driving force to adoption doesn't exist, as ordinary people have no real clue as to what constitutes a "good" language, so that this random walk or Brownian motion has to play an important role in the early stages of adoption.

So with that as a premise, we have to determine how to model this effect mathematically. Incrementally we wish to show that the growth term gets suppressed by the potential for fluctuation in the early number of adopters. A weaker steady growth term will take over once a sufficiently large crowd joins the bandwagon.

dn = dt / (C/sqrt(n) + K)

In this differential formulation, you can see how the fluctuation term which goes as 1/sqrt(n) suppresses the initial growth until it reaches a steady state as the K term becomes more important. Integrating this term once and we get the implicit equation:

2*C*sqrt(n) + K*n = t

Plotting this for C=0.007 and K=0.000004, we get the following growth function.

Figure 1 : Growth function assuming suppression during early fluctuations

This makes a lot of sense as you can see that growth occurs very slowly until an accumulated time at which the linear term takes over. That becomes the saturation level for an expanding population base as the language has taken root.

To put this in stochastic terms assuming that the actual growth terms disperse across boundaries, we get the following cumulative dispersion (plugging the last equation into the first equation to simulate an ergodic steady state):

P(n) = 1/(1+1/g(n)) = 1/(1+1/(2*C*sqrt(n) + K*n))

I took two sets of the distribution of population sizes of languages (DPL) of the Earth’s actually spoken languages from the references below and plotted the entropic dispersion alongside the data. The first reference provides the DPL in terms of a probability density function (i.e. the first derivative of P(n)) and the second as a cumulative distribution function. The values for C and K were as used above. The fit works parsimoniously well and it makes much more sense than the complicated explanations offered up previously for language distribution.

Figure 2 : Language diversity (top) probability density function (below) cumulative. The entropic dispersion model in green.

In summary, the two pieces to the puzzle are assuming dispersion according to the maximum entropy principle, and a suppressed growth rate due to fluctuations during the early adoption. This gives two power law slopes in the cumulative; 1/2 in the lower part of the curve and 1 in the higher part of the curve.

References

Scaling Relations for Diversity of Languages (2008)
Competition and fragmentation: a simple model generating
lognormal-like distributions (2009)
Scaling laws of human interaction activity (2009)
Discussions on the fluctuation term.

NY Math Teacher Howard A. Stern Uses Ingenuity To Overcome Failure Statistics

The public school teacher highlighted in the linked article has this to say:

"So much of math is about noticing patterns," says Stern, who should know. Before becoming a teacher, he was a finance analyst and a quality engineer.

I always try to seek interesting patterns in the data, but more to the point, I try to actually understand the behavior from a fundamental perspective.

One way Stern uses technology is by helping his students visualize his lessons through the use of graphing calculators.

Stern has it exactly right, if we treat knowledge seeking as a game, like a suduko puzzle, we can attract more people to science in general.

I think that the pattern in language distribution has similarities to that of innovation adoption as well, similar to what Rogers describes in his book "Diffusions of Innovations". I will try to look into this further as I think the dispersive arguments holds some promise as an analytical approach.

Stock Market as Econophysics Toy Problem

2010-10-12T23:07:00.000-07:00

Consider a typical stock market. It consists of a number of stocks that show various rates of growth, R. Say that these have an average growth rate, r. Then by the Maximum Entropy Principle, the probability distribution function is:

pr(R) = 1/r*exp(-R/r)

We can solve this for an expected valuation, x, of some arbitrary stock after time, t.

n(x|t) = ∫ pr(R) δ(x-Rt) dR

This reduces to the marginal distribution:

n(x|t) = 1/(rt) * exp(-x/(rt))

In general, the growth of a stock only occurs over some average time, τ, which has its own Maximum Entropy probability distribution:

p(t) = 1/τ *exp(-t/τ)

So when the expected growth is averaged over expected times we get this integral:

n(x) = ∫ n(x|t) p(t) dt

We have almost solved our problem, but this integration reduces to an ugly transcendental function K₀otherwise known as a modified Bessel function of the second kind and order 0.

n(x) = 2/(rτ) * K₀(2*sqrt(x/(rτ) ))

Fortunately, the K₀function is available on any spreadsheet program (Excel, OpenOffice, etc) as the function BESSELK(X;0).

Let us try it out. I took 3500 stocks over the last decade (since 1999), and plotted the histogram of all rates of return below.

The red line is the Maximum Entropy model for the expected rate of return, n(x) where x is the rate of return. This has only a single adjustable parameter, the aggregate value rτ. We line this up with the peak which also happens to coincide with the mean return value. For the 10 year period, rτ = 2, essentially indicating an average doubling in the valuation of the average stock. This doesn't say anything about the stock market as a whole, which turned out pretty flat over the decade, only that certain high-rate-of-return stocks upped the average (much like the story of Bill Gates entering a room of average wage earners).

The following figure shows a Monte Carlo simulation where I draw 3500 samples from a rτ value of 1. This gives an idea of the amount of counting noise we might see.

I should point out that the MaxEnt model shows very little by way of excessively fat tails at high returns. A stock has to both survive a long time and grow at a rapid enough rate to get too far out in the tail. You see that in the data as only a couple of the stocks have returns greater than 100x. I don't rule out the possibility of high-return tails but we would need to put even more disorder in the pr(R) distribution than the MaxEnt provides for a mean return rate. The actual data seems a bit sharper and has more outliers than the Monte Carlo simulation, indicating some subtlety that I probably have missed. Yet, this demonstrates how to use the Maximum Entropy Principle most effectively -- you should only include the parameters that you can defend. From this minimal set of constraints you observe how far this can take you. In this case, I could only defend some concept of mean in rτ and then you get a distribution that reflects the uncertainty you have in the rest of the parameter space.

The stock market with its myriad of players follows an entropic model to first-order. All the agents seem to fill up the state space so that we can get a parsimonious fit to the data with an almost laughably simple econophysics model. For this model, the distribution curve on a log-log plot will always take on exactly that skewed shape (excepting for statistical noise of course) -- it will only shift laterally depending the general direction of the market.

The stock market becomes essentially a toy problem, no different than the explanation of statistical mechanics you may encounter in a physics course.

Has anyone else figured this out?

[EDIT]
Besides the slight fat-tail, which may be due to potential compounding growth similar to that found in individual incomes, the sharper peak may also have a second-order basis. This could result from a behavior called implied correlation which measures the synchronized behavior among stocks in the market. According to recent measurements, the correlation has hit all-time highs (the last around October 5). Qualitatively a high correlation would imply that the average growth rate r would show much less dispersion in that variate, and the dispersion would only apply to the length of time, t, that a stock rides the crest. Correlation essentially removes one of the parameters of variability from the model and the distribution sharpens up. The stock distribution then becomes the following simple damped exponential instead of the Bessel.

n(x) = 1/(rτ) * exp(-x/(rτ))

The figure below shows what happens when about 40% of the stocks would show this correlation (in green). The other 60% show independent variability or dispersion in the rates as per the original model.

I don't think this makes the collective stock behavior and more complex. I think it makes it simpler in fact. Implied correlation actually points to the future in the stock market. Dispersion in stock returns will narrow as all stocks move in unison. It makes it even more of a toy, with computers potentially dictating all movements.

Implied correlation has risen in the last few years (from here)

References
I personally don't deal with the stock market, preferring to watch it from afar. I found a few papers that try to understand this effect, but most just try to brute force fit it to various distributions.

Analysis of same data from Seeking Alpha
This paper is close but no cigar. It looks like they "detrend" the data to get of the skew, which I think misses the point :
"Microscopic origin of non-Gaussian distributions of financial returns" (2007)
This book has info on the Bessel distribution:
"Return distributions in finance", J. Knight and S. Satchell
Interesting from an econophysics perspective.
This book appears worthless:
"Fat-Tailed and Skewed Asset Return Distributions", S.T. Rachev, F.J. Fabozzi, C Menn

Black-Scholes

2010-10-07T13:19:00.000-07:00

Games for suits. This post has no relevance in the greater scheme of things.

As a premise, consider that the financial industry needs instruments of wealth creation that work opposite to that of stocks. For example, when stock prices remain low, then something else else should take up the slack -- otherwise important people won't make money. Wall Street invented derivatives, options, and other hedging methods to serve as an investment vehicle under these conditions.

We can try to show how this works.

If S is the stock price, then V ~ 1/S is an example "derivative" that works as a reciprocal to price. This becomes the normative description and defines the basic objective as to what the investment class wants to achieve -- an alternate form of income that balances swings in stock price, potentially reducing risk.

Further, we make the assumption that the derivative will grow or decline over time.

So we get:
V(S,t) = K/S * exp(a*t)

If a > 0 then the derivative will grow and if a is less than zero than the derivative will damp out over time. The term K is a constant of proportionality.

The infamous Black-Scholes equation supposedly governs the behaviour of derivatives with respect to stock prices (and time) according to this invariant:

The particulars may change but this formulation describes THE equation that Merton, Black, and Scholes devised to aid investors in making hedged investments using options and other derivatives. The way to read this equation is to note that derivatives will drift or diffuse into the space of the stock price, and proportional to the stock price itself. The drift term occurs due to the interest rate r providing a kind of forcing function. The derivative, V, can also grow due to pure interest rate compounding, as seen in the last term. Whether this actually holds or not, I don't really care as I don't participate in these schemes.

So if you look at it from a very neutral perspective you come up with some interesting observations. For one, you can trivially solve this partial differential equation for a generally disordered set of initial conditions. And the solution appears exactly the same as my first expression above:
V(S,t) = K/S * exp(a*t)

To verify this assertion, we test the expression in the B-S equation, substituting the partial derivatives as we go along.

a*K/S* exp(a*t) + 1/2(σS)²*2*K/S²*exp(a*t) - rS*K/S²*exp(a*t) - r*K/S * exp(a*t) = 0

Cancelling out all common factors:

a/S + 1/2(σS)²*2/S²- rS/S² - r/S = 0

Reducing the value of S

a/S + 1/2(σ)²*2/S- r/S - r/S = 0

a + 1/2(σ)²*2- r - r = 0

gets us to:

a = 2*r - σ²

The term r is proportional to interest, and σ is volatility or variance in stock price.

So this simple expression that I just cooked up will obey Black-Scholes as long as we choose the constant a term to correspond to the interest and volatility as shown above, and we get:
V(S,t) = K/S * exp((2*r - σ²)*t)

Note that if the volatility (i.e. diffusion) stays high relative to interest, the exponential will damp out with time. If interest (i.e. drift) goes higher than volatility, the exponential will accelerate, creating a huge amount of paper gains.

At this point someone will argue that this solution does not reflect reality. I beg to differ. When you make your bed of mathematical box-springs, you have to lie in it. This solution to Black-Scholes is perfectly fine as it gives a steady-state picture of the partial differential equation. The diffusional and drift components cancel with the right mix of production vs destruction in derivative wealth. If you don't like it, then come up with something different than that specific B-S equation.

I have a feeling that all the seeming complexity of financial quantitative analysis with its Ito calculus and Wiener processes acts as a shiny facade to a simple reality. The math exists to model the inverse relationship of stocks to derivatives. If this didn't happen -- and the lords of high finance absolutely require this relationship to make money -- the math as formulated would vanish from their toolbox. In other words, the math only exists to justify what the financial operatives want to see happen. Everyone appears to implicitly buy this mathematical artifice hook, line, and sinker.

Quantitative analysis and the "quants" who work it have created a fantasy land, where they do not want you to know how easily their quaint ornate universe reduces to a simple function. If they admitted to the charade, the mystery would all disappear and they would no longer have jobs.

Economics and finance does not constitute a science. In science you may need to use partial differential equations. For example, the Fokker-Planck equation shows up quite often -- which incidentally, the Black-Scholes equation shows some similarity to and the quant proponents of B-S certainly like to play up -- but it typically applies to real, physical systems where you use it to try to understand nature, not trying to model some artificial game-like behavior.

I can edit my solution into the Wikipedia page for Black-Scholes and I will bet that someone will immediately remove it. I harbor no illusions. The financial industry depends on the absence of real knowledge to achieve their objectives.

That explains why economics and finance do not classify as sciences; absolute truth does not matter to economists and financiers, only the art of deconstructing profit and the craft of phantom wealth creation does.

Please address editorial comments to: Postings, Main Incinerator, Department of Sanitation, North River Piers, New York, N.Y. 10019.

Lake Size Distributions

2010-10-02T10:28:00.000-07:00

Our environment shows great diversity in the size and abundance in natural structures. Since we extract oil from our environment, it stands to reason that many of the same mechanisms leading to oil formation could also reveal themselves in more familiar natural phenomena. Take the size distribution of lakes as an example.

Freshwater lakes accumulate their volume in a manner analogous to the way that an underground reservoir accumulates oil. Over geologic time, water drifts into a basin at various rates and over a range in collecting regions. In the context of oil reservoirs, I have talked about this behavior before and the Maximum Entropy prediction of the size distribution leads to the following expression:

P(Size) = 1/(1+Median/Size)

Surveys of lake size show the same reciprocal power law dependence, with the exponent usually appearing arbitrarily close to one. In Figure 1 below, the data plotted on a ranked plot clearly shows this dependence over several orders of magnitude.

Figure 1: Northern Quebec lakes [1]

More revealing, in Figure 2 we can observe the bend in the curve that limits the number of small lakes in exact accordance to the equation. The agreement with such a simple model suggests that a universal behavior links the statistics between environmental phenomena as seemingly distinct as those of lakes and oil reservoirs.

Figure 2: Amazon lakes [2]

This provides other intuitive clues to how to think about reservoir sizing. Consider the fact that very few freshwater lakes reach gigantic portions, the Great Lakes serving as a prime example. Similarly, the rare occurrence of “super-giant” reservoirs follow from the same principles. We clearly won’t find any new huge freshwater lakes, while the future occurrence of super-giant oil reservoirs remains very doubtful just from the statistics of oil reservoirs found so far. Finding substantial numbers of super-giant reservoirs would result in deviations from the size distribution plot, making it very unlikely.

References
[1] K&C Science Report – Phase 1 Global Lake Census

[2] Estimation of the fractal dimension of terrain from Lake Size Distributions

Hydrogeology for Dummies

2010-09-06T22:23:00.000-07:00

A running theme of this blog involves the reduction of seemingly complex behaviors into simple mathematical formulations. It remains a bit of a mystery to me why in many situations that no one has either (a) done this work on their own or (b) uncovered the work of someone else who has done the simplifying analysis years ago.

The majority of scientists practicing mainstream research have furthered the cause by following the lead of others who go down blind alleys and over-complicate the analysis. I suspect that a few complicate matters intentionally, as it demonstrates to other scientists their intellectual prowess. In certain cases, creating a private world of intricate analysis acts as a kind of moat around which they can fortify their specialty discipline.

Of course, this doesn't happen universally. Certainly we run across many scientific and engineering subdisciplines that have gone through years of scrubbing. In these cases, the most salient and simple analyses have emerged and stood the test of time. They often share the same traits of elegance and crystalline transparency so that we can use their patterns to understand the world without a lot of extra effort. To me, that seems a reasonable goal to strive for.

In this post, I will go through the derivation of what I consider a very overlooked and simple argument having to do with the transport of materials in porous media -- much as what you would find in tracing a contaminant though a groundwater basin. Or what may happen if you frac for natural gas and open up new pathways to a drinking water aquifer. Or how oil will migrate to a reservoir over time, feeding the production output of a stripper well for years. Or what happens if you spill oil in a waterway.

Unfortunately, when you pose this kind of problem to a research geologist or hydrologist, you will have to prepare for an onslaught of ornate misdirection. They will either derive some hideous numerical model or possibly run a piece of commercial software. Apparently, they will never resort to plain logic and elementary first-principles considerations.

The Problem

1. Consider a contaminant that enters an aquifer in a single dose
2. Predict how long it will take to pass by a downstream location
3. How do you solve this problem?

A large scale experiment typically looks like this scenario:

from Groundwater Tracing in the Woodville Karst Plain

And you get a result that looks like the following figure. Intuitively, one would expect that the concentrated dose will disperse as it travels downstream and that the original concentration will spread out in time. The red curve that goes through the data gives you a feel for what I will derive via a simple model.

As a main premise, I assume that disorder plays a big role in providing a variety of pathways from source to sink. One can imagine that some paths might occur on the main waterway, providing a maximum speed or path of least resistance. Other paths may follow obstructions or diversions which will either slow down or speed up the flow from the main path.

The main path has a mean velocity v₀ and the other paths have probabilities that range below this, with some mean deviation v_m from v₀. A distribution that maximizes entropy while holding to these two minimal constraints looks like the following graph.

Figure 1: MaxEnt velocity distribution for absolute mean deviation

This illustrates simple dispersion. For this post we won't even consider diffusion, which although important may in fact act as only a second-order effect depending on the speed of the main flow.

The calculation of downstream concentration, n(x,t), drops out of the Fokker-Planck equation if we ignore diffusion. Note the delta function, δ(x-vt), which describes a traveling pulse for each velocity component.

n(x,t) = ∫ p(v) δ(x-vt) dv

Next we apply the Maximum Entropy Principle to generate a velocity distribution as shown in the Figure 1:

p(v) = 1/v_m exp(-|v-v_o|/v_m)

No other distribution has a higher entropy given that mean and an absolute deviation from the mean, so it ranks as the least biased estimator for that set of constraints. (Note that this does not describe the normal or Gaussian distribution as that requires a second-moment, i.e variance, constraint. It turns out that the mean deviation distribution, also known as the Laplace, is actually a smeared Gaussian where we have MaxEnt uncertainty in σ-squared. So Laplace entropy is higher than the Gaussian entropy)

We can trivially solve the integral to generate a concentration at some downstream location x (forget about adding extra dimensions as a one-dimensional result should suffice).

n(x,t) = 1/(v_mt) exp(-|x/(v_mt)-v₀/v_m|)

Let's see how this works in practice.

I pulled data from a pair of papers from 2008, "Non-Fickian dispersion in porous media", T Le Borgne, P Gouze, et al. The scientists created a carefully controlled experiment, which relied on a customized apparatus for making precise measurements of the contaminant, a flourescent dye called uranine. The value of this particular experiment lies in the large dynamic range of the resultant data. The concentration runs over 4-orders of magnitude and the time scale 2-orders. Their own model, although generating a good fit to the data, needed a numerical calculation to solve, violating my assertion that we can model via simpler mechanisms.

The following figure allows for the wide dynamic range by plotting the concentration (also known as a breakthrough curve) on a log-log scale. The red triangles ◊ fit the Maximum Entropy dispersion model, n(x,t), for a fixed value of x and a value of v_m/v₀ = 0.18. By inverting the concentration we can get the probability distribution of velocities in the bottom figure; on a semi-log plot a symmetric two-sided exponential looks like a perfect isosceles triangle. Based on the outstanding fit and symmetric distribution I find it blatantly obvious that entropic mechanisms generate the dispersion observed. You won't get this parsimonious a fit from such a simple model -- with essentially a single parameter v_m/v_{0 --}unless it has some real merit.

Figure 2: Breakthrough curve (top) and
measured velocity distribution (bottom)
for flourescent dye tracer experiment.

I would suggest that any further modeling of these kinds of porous structures makes little sense since we have essentially proved that the multitude of the pathways maximize entropy and thus maximized the disorder of the system. In other words, you could not model a more complex system given those constraints if you tried. Nature will always win out with entropy in its back pocket.

The simplicity of the model also points out how readily fat-tail effects emerge from entropic disorder. The power law drop-off obeys a 1/time behavior that certainly has consequences in terms of how long a contaminant will remain in a groundwater basin. Velocity dispersion with a mean MaxEnt constraint will always lead to a power-law drop-off in time (see more here).

See also these posts:

The hydologists and geologists who ignore entropy in favor of some other fancy model do so based on their own stubborness or ignorance. I have observed the practice of making things too complicated runs rampant among geologists and it really strikes me as kind of sad. We have hydrogeologist hacks like Steven Gorelick writing cornucopian books diminishing the significance of peak oil, when they can't even do the science of their own discipline correctly.

Tasseography

2010-08-20T06:49:00.000-07:00

Oil Watch Monthly

Because of the magnified nature of the production scale I find it interesting to place the data on the real scale, which shows the zeros and the full temporal range. See the short black segment in the following figure, which signifies the range reported on TOD.

I don't really understand this infatuation with what I consider noise riding on top of the more important overall scaled profile. Readers must feel a need to see this magnified view which I don't quite grasp.

Is it because people have become accustomed to using the information for futures trading or anticipating the stock market? I presume that every little glitch provides a chance to make some money.

Or do we suffer from climate change envy where temperature trends get studied to death? That works in a different context because temperatures normally occupy a narrow range and the important signal can get buried in the measurement noise.

Or do people want to anticipate seeing that sudden, precipitous drop that will signal us going over the cliff?

More likely the answer is that we continue to plot the magnified view because we can and it gives us a strawman to argue back and forth over. The term tasseography describes this behavior.

Noise can tell us something but it to first-order it really only tells us what we already know. The fewer the number of independent measurements or actors in the market, the greater the noise and fluctuations.

GOM Maximum Production Rate and Macondo

2010-07-01T20:35:00.000-07:00

I did some analysis based on Berman's post from a few days ago:
(Estimated Oil Flow Rates From the BP Mississippi Canyon Block 252 “Macondo” Well)

I think he messed up the statistics because of his use of a truncated data set from the MMS and the log-normal distribution he used.

I wasn't sure exactly how he got his data but I essentially had to screen scrape the data off of about 18 PDF files giving the Maximum Production Rate (MPR) going back to 1975: http://www.gomr.mms.gov/homepg/pubinfo/repcat/product/MPR.html

I plotted the results histogram against a model of dispersive aggregation for reservoir sizes. The maximum rate is then a simple proportional draw-down from the reservoir size. Bigger reservoirs have a higher rate and smaller reservoirs have a smaller rate -- nothing to argue about here as it is a pretty safe approximation. The way you read this histogram is that the flat regions have the highest frequency.

The integrated underneath the two curves is equal and about 16.5 million barrels per day peak. Don't confuse this with any rate attainable from the GOM; it is high because it sums up the peaks from a span of years. The median value is 200 barrels per day.

The interesting point in the curve is that the model predicts a higher peak rate for the largest reservoirs, the curve goes off the graph to above 400,000 barrels per day. Now, I would think that the operators would never try to have that throughput from a single well. So what do they do? Of course they split it into several wells to extract the maximum amount from that reservoir and essentially throttle that from an individual well.

Since the total amount is conserved between the two curves, the bulge that you see in the data is the extra wells drilled to make up for the excess. My model is totally based on the principle of Maximum Entropy applied to reservoir sizing, and the reordering of the rank histogram is caused by artificial constraints set by human intervention. Notice that all the small reservoirs effectively require no throttling.

The point of this comment is that working wells are likely throttled but the Macondo could conceivably be higher than the maximum of 50,000 barrels per day that Berman suggested. The operators have no way of throttling it until the relief wells are put in place. Of course this kind of throughput is very rare, as at the most a couple of dozen out of 10,000 reservoirs will get this big and generate this potential, but this is the way that nature operates, a big fat-tail effect.

Petroleum Engineering

2010-06-19T18:11:00.000-07:00

With all the discussion on the Gulf Oil disaster going on, lots of petroleum engineers and others from the oil industry have pitched in with their opinions. In which case we can see exactly what they think of their profession.

One commenter, an authority on reservoir engineering apparently had this to say about Peak Oil:

We understand how our business works, certainly. Guys like us, (those IN THE KNOW) have been declaring the end of oil since at least 1886. In Pittsburgh to be specific. Can't say we didn't give the rest of you noobs plenty of warning.

So let me understand this statement. Oil industry types apparently have always known that the end of oil would occur since day one. I wonder why no one thought to just ask them? How did we miss that one?

This same fellow has huge problems with my analysis, because he thinks that what I do amounts to "curve fitting".

I mean seriously, who else would confuse curve fitting with knowledge?

In truth, most of the forecasters who point to continually increasing oil production well into the future base their projections on very little real knowledge. They actually practice curve fitting, i.e. fitting a curve to the production level that we need, because they have no other justification for a realistic outlook.

Bayesian analysis works by using past knowledge to predict future outcomes. We have so much knowledge about previous discoveries, reserve growth mechanisms, and extraction rates that our ability to predict should work very effectively ... if we would just start universally using this kind of approach. The other benefit is that the analysis keeps on getting better and better with time due to the Bayesian updating process. The mathematician Laplace first applied this powerful mode of probabilistic reasoning in the late 1700's to real problems, but we still have holdouts in various disciplines. To top it off, if you have a real model underneath the knowledge, it makes the forecasting that much more powerful.

Let them get through diffy-q, I suppose the only other gang besides engineers forced through that one are the more mathematically inclined....and they are mostly jealous because their theoretical skills don't translate into income very well.

Common knowledge in college that students that went into geology, civil, and petroleum engineering didn't want to get stick in a desk job. Lots of them could not imagine being sedentary for 8 hours a day.

Hubbert peak in Five Easy Pieces

2010-06-16T16:50:00.000-07:00

Based on the increase in spill rate from the leaking Gulf of Mexico oil well, HO at TheOilDrum.com suggested a potential explanation. His post essentially argued that sand particles acting as a strong abrasive driven along by the already high velocity stream of escaping oil leads to increasing in the channeling and thus an even faster leak rate.

HO described a process known as CHOPS (Cold Heavy Oil Production with Sand) which can enlarge a well's streaming throughput by promoting the formation of heavily eroded channels. The TOD post provided the following picture of the possible outcome of the behavior.
Note that the lower curve shows the typical output from a throttled flow. Above that curve, the modulated line shows the results of an accelerated extraction -- note that a peak actually appears which pinpoints the maximum flow rate. In terms of the oil spill, we don't want this behavior because it gives us less time to fix or relieve the problem well. Yet, ordinarily we want this same behavior -- that of fast extraction -- in practical situations because we want and need the oil right now! (so that oil companies can make money, of course)

Which leads me to formulating the following very simple but physically correct model of Hubbert's Peak. You won't find this anywhere else, because this derivation does not jive with how geologists think about oil extraction. They get many of the pieces but they never put them all together.

I will offer up a derivation for this behavior leading to a Hubbert Peak in 5 easy pieces.

Piece 1. The standard assumption of draw-down from a reservoir results in an exponential decline over time. You can consider that the exponential shape results from a law of diminishing returns; in that a constant amount proportional to the remainder draws down per unit time. Or you can say that a maximum entropy range of extraction rates gets applied to the volume. A proportional extraction rate that we call R defines the mean and U₀ is the reservoir size. U(t) gives us the cumulative reserve.

U(t) = U₀*exp(-R*t)

Piece 2. Next, we realize that we have uncertainty over the size of the reservoir; the U₀ we have defined actually only serves as an estimate of the size. This means we have an uncertainty over the rate of proportional extraction as well. This turns into a form of hyperbolic discounting and the cumulative draw-down actually looks like this.

U(t) = U₀ / (1+R*t)

Note the fat-tail.

Piece 3. Next we assert that the constant but uncertain proportional extraction rate undergoes an acceleration starting from the original value, R(t) = R₀ + k*t. This acceleration equates to Newton's law, first-order with time. Then the instantaneous absolute rate of extraction from the remaining reservoir looks like:

RateOfExtraction(t) = -dU(t)/dt = U₀*(R₀ + k*t)/(1+R₀*t+k*t²/2)²

For R₀=0.5 and k=2, it results in this shape

This curve we can scale and overlay on top of the CHOPS curve to validate our thought process.

Piece 4. Over a larger set of reservoirs that experience a technical improvement over time, we can assume that the proportional extraction rate can accelerate even more strongly over time, R(t)=C*exp(k*t). This gives us a Moore's law form of acceleration, doubling every set number of years. Then

RateOfExtraction(t) = -dU(t)/dt = U₀ * R(t) / (1+integral(R(t)dt))²

= U₀*C*exp(k*t)/(1+C/k*(exp(k*t)-1))²

For a small starting rate, the acceleration further accentuates the subtle peak that we observe in piece 3 and it turns into a full-fledged symmetric peak as shown in the next figure:

Piece 5. Congratulations. You haven't broken any rules and you have just derived the famed Hubbert Peak, also known as the Logistic Sigmoid function.

Some Backstory
An alternate derivation exists for the corresponding discovery peak, which I call Dispersive Discovery. There, the uncertainty involves how much volume gets explored and at what rate, otherwise the math turns out exactly the same. Both derivations result from an assumed finite constraint but uncertainty in both rates and subvolumes. The only problem with using the Hubbert peak derivation for extraction is that it premises that each extraction rate started at the same time (globally this would be 1858). We know that this has not happened for global production, as extraction can only start after a discovery, and then some variable hold time. By using dispersive discovery, we get a larger spread in start years, and then The Oil Shock model generates the extraction curve. In general, if the discovery peak precedes the oil production peak by a number of years, I would use Dispersive Discovery, but if the two coincide, then extraction tracks discovery and it doesn't really matter how you interpret the rates. This explains why this particular derivation works well for more localized production areas that have seen significant technology changes. In contrast, the technology of discovery has undergone tremendous technology changes over the years, so that dispersive discovery works very well in terms of global modeling. This is actually not much of a caveat, as the more ways that you can find the same result, the more confidence you have that you have remained on the right track.

The current derivation also points out the huge hole in the technique known as Hubbert Linearization (HL). As defined, HL derives from the observation that

dU(t)/dt = U(t)*(U₀-U(t))

Yet this only works for the one case where we can define R(t) as an exponential function, that of piece 4. The formula does not work for either piece 1, 2, or 3. Therefore, HL only serves as a curious mathematical identity for that one exponential case, which we know does not always occur.

The actual "WebHub" Linearization takes the following form:

dU(t)/dt = -U₀ * R(t) / (1+integral(R(t)dt))²

This may not prove as handy as HL perhaps, but it has the benefit of correctness, and it works well for certain cases.

Like me, Robert Rapier has railed against the inadequacy of HL and this may take up the slack.

GOM Reservoir Size Distributions

2010-06-14T18:58:00.000-07:00

Question:

BigMoose on June 14, 2010 - 6:26pm Permalink | Subthread | Parent | Parent subthread | Comments top
I have heard many unofficial estimates of the magnitude of oil in this formation... 2nd largest in America, 2nd largest in the world...

Does anyone have a credible estimate on the formation reserves?

Some historical data available from the MMS.
http://www.gomr.mms.gov/PDFs/2009/2009-064.pdf

On the basis of proved oil, for 8,014 proved undersaturated oil reservoirs, the median is 0.3 MMbbl, the mean is 1.8 MMbbl.

Peak Oil theory (Entropic Dispersive Aggregation) says the cumulative size distribution of reservoirs (ranked small to large) goes as P(Size)=1/(1+0.3/Size) if we assume a median of 0.3. It doesn't quite follow this exactly because infinite sized reservoirs can not exist.

If you want the raw data it is here:
file:///G:/RE/Shared/EOGR%20Report/2008-034%20Estimated%20Oil%20and%20Ga...

Sorry, that was a joke, the MMS puts the information on a public web server, and the data is retrieved as a local filesystem URL?
http://www.gomr.mms.gov/homepg/pubinfo/freeasci/geologic/estimated2006.html

I placed whatever data I could get into Google Docs, and placed theory next to it.

The MMS is to be split into 3 agencies apparently. Throughout their history, they failed in doing any kind of useful depletion analysis in the GOM. Anyone can collect data; interpreting it is the challenging part.

The Mentaculus

2010-06-12T15:00:00.000-07:00

I saw the Coen brothers movie "A Serious Man" a few months ago. A definite period piece from the 1960's, it contrasted two scientists, one an academic and one a hapless amateur. The main protagonist, Larry Gopnick, a physics professor at what looks like a small liberal arts school in the Twin Cities (Macalester, Hamline maybe?), spends time teaching his students what look like elaborate mathematical derivations on a huge chalkboard. He has trouble dealing with some of his students on occasion:

Clive Park: Yes, but this is not just. I was unaware to be examined on the mathematics.
Larry Gopnik: Well, you can't do physics without mathematics, really, can you?
Clive Park: If I receive failing grade I lose my scholarship, and feel shame. I understand the physics. I understand the dead cat.
Larry Gopnik: You understand the dead cat? But... you... you can't really understand the physics without understanding the math. The math tells how it really works. That's the real thing; the stories I give you in class are just illustrative; they're like, fables, say, to help give you a picture. An imperfect model. I mean - even I don't understand the dead cat. The math is how it really works.

His academic colleagues want Professor Gopnick to publish articles at some point (with the implicit threat of not getting tenure). Gopnick's main problem lies in his rationality:

But his rigid framing of a cause-and-effect universe makes him indignant about lack of apparent cause ...

Gopnick's brother, the minor character of Uncle Arthur, takes the role of an almost savant numerologist, busy at work on a treatise called The Mentaculus. Filled with dense illustrations and symbology, it apparently functions as a "probability map" in what appears to spell out a Theory of Everything. It also apparently works to some extent:

We might guess that it makes no sense, but Arthur's "system" apparently "works" as intended, and he applies it to winning at back room card games.

Based on the events that eventually transpire, the theme of the movie essentially says that if you seek rationality, you will ultimately only land on random chance.

I consider myself a "serious man" as well. But do I have a variation of The Mentaculous buried in the contents of this blog?

I tried to make a probability map of all the applications and blog links that I have worked on relating to what I call entropic dispersion in the following table [full HTML]:

The math is how it really works. Perhaps I should publish. Yet blogging is too much fun. Perhaps I need to take a canoe trip.

Good reads describing The Mentaculus of probability and statistics

"Dawning of the Age of Stochasticity", David Mumford
From its shady beginnings devising gambling strategies and counting corpses in medieval London, probability theory and statistical inference now emerge as better foundations for scientific models, especially those of the process of thinking and as essential ingredients of theoretical mathematics, even the foundations of mathematics itself.
"Probability Theory: The Logic of Science", Edwin T. Jaynes

Our theme is simply: probability theory as extended logic. The ‘new’ perception amounts to the recognition that the mathematical rules of probability theory are not merely rules for calculating frequencies of ‘random variables'; they are also the unique consistent rules for conducting inference(i.e. plausible reasoning) of any kind. and we shall apply them in full generality to that end.
"On Thinking Probabilistically", M.E. McIntyre
"The Black Swan" and "Fooled by Chance", N.N. Taleb

Worst Book on Oil Crisis Written Yet

2010-06-11T06:58:00.000-07:00

Former USGS staffer Steven Gorelick has written a book called "Oil Panic and the Global Crisis: Predictions and Myths". It has to rank as the worst of the neo-cornucopian books out there simply because it actually spreads myths instead of deeming to correct them, as the title implies.

The author acts the role of a somewhat neutral bystander and balanced pseudo-journalist, never giving the appearance of a rabid oil cornucopian, yet slipping in so many groaners that he basically gives away his not-so-hidden agenda. From a scientific context, providing both sides of the story makes no sense when the objective is truth rather than balanced reporting. Excerpts of the book would fit right into a Fox news piece.

To give a taste of how little original research that Gorelick has actually performed and how much he relies on other cornucopians, consider the passage wherein he references geology professor Larry Cathless. On page 128, Gorelick quotes Cathles as saying that we may find as much as "1 trillion barrels of oil and gas in just a portion of the gulf oil sediments".

I found the original statement by Cathles here:

Cathles and his team estimate that in a study area of about 9,600 square miles off the coast of Louisiana, source rocks a dozen kilometers down have generated as much as 184 billion tons of oil and gas — about 1,000 billion barrels of oil and gas equivalent. "That's 30 percent more than we humans have consumed over the entire petroleum era," Cathles says. "And that's just this one little postage stamp area; if this is going on worldwide, then there's a lot of hydrocarbons venting out."

Although not directly implicated as an abiotic oil advocate (unlike his late Cornell University colleague Thomas Gold), former Chevron employee Cathles has close ties to the largely mythical Eugene Island story. Several years ago new discoveries from the previously tapped-out Eugene area had people's hopes up that somehow oil reservoirs could go through a near real-time "replenishment".

"We're dealing with this giant flow-through system where the hydrocarbons are generating now, moving through the overlying strata now, building the reservoirs now and spilling out into the ocean now," Cathles says.

Well, as it turned out, the Eugene Island secondary production turned out just a blip on the radar screen, yet Cathles still gets a mention as a credible source? (Think about it, if this turned out true, then the recent Gulf Oil spill could allow a never-ending release of hydrocarbons from beneath the waters, as this urban legend gets repeated still. How embarrassingly timely for Gorelick.).

Elsewhere, the book becomes safe pablum for a narrowly defined audience. Note the limited depth of Gorelick's analysis and the intentional dumbing down in his writing:

Hubbert used a straightforward formula that yields the curve as illustrated in Figure 1.2. The logistic-curve formula is a simple expression with three adjustable parameters (mathematical knobs) that control the slope, peak, height and time of peak

Now you see what happens when an author keeps it too simple. He ends up never explaining anything about the logistic, apart from providing the functional form in a footnote, and makes it worse by calling the parameters "mathematical knobs". That essentially gives a flavor of the depth of the mathematics.

Gorelick has an entire chapter called "Counter-Arguments to Imminent Oil Depletion". Notwithstanding that oil depletion is imminent by definition (it certainly does not regenerate contrary to the implications), this chapter contains some of the most unscientific assertions that I have come across. Consider this bullet point coming from Gorelick:

- The world has never run out of any significant globally traded, non-renewable Earth resource.

This false equivalency comes somewhere from the list of logical fallacies. I find it bizarre that a reputable scientist would appeal to this kind of argument. Further he bullet points:

- The trends in production of global oil and natural gas have not declined as predicted.

I call a strawman fallacy as no one has really come up with a formal theory for depletion. Instead every oil prediction that I have seen has relied on some sort of ad hoc analysis via heuristics. So to imply that something has not followed as predicted does not prove anything. As I have said before, heuristics do not substitute for theory and Gorelick unfortunately has not contributed any research of his own.

I listed only 2 of the 21 bullet pointed counter-arguments that Gorelick concludes the chapter with. I can understand the need for these bullet points if he wanted to act like an objective journalist wanting to tell both sides of the story. Yet we have all learned from Krugman that real science does not scream headlines that say "Shape of Earth--Views Differ". A scientist should dig deep and try to come up with a model or theory that would confirm or rebut the empirical evidence. You just don't rely on tired worn-out assertions (the world has never run out of a resource, predictions have not come true, etc) from the cornucopian right, put them in a book and consider this an advancement of knowledge.

The book industry likely published Oil Panic because it does not even remotely challenge business as usual and actually condones the cornucopian viewpoint.

End of book review.

Musings

Since Gorelick has propagated half-truths and not resolved any myths at all in the oil depletion realm, I figured I would return the favor in his own research area. From his CV, the "honored and awarded" Gorelick moved on from the USGS and became a professor of hydrogeology and part of the Environmental Earth System Science department at Stanford University. If he can write a book on peak oil and turn back progress on understanding oil depletion, I can opine on hydrogeology.

From his research papers, Gorelick claims to understand how to model principles of hydrogeology and presumably knows about breakthrough curves. It turns out that most of the dispersive transport involved in hydrology applications hinges on some very simple overriding principles. These principles are so obvious to me that I don't understand why the brilliant scientific minds in geology have not figured this out. Consider that Gorelick has expertise in "multiple-rate mass transfer" which I associate this with the simple idea of dispersion applied to material transport. I actually ran across Gorelick's work prior to reviewing his book because of my studies of generalized dispersive transport.

As Gorelick should know, all processes do not proceed at the same rate, and this includes variations in oil discovery rates around the world. This leads directly to the fat-tail effects that I see in oil reserves and to the fat-tails that Gorelick observes in solute transport in his groundwater contamination studies. Not all solute diffuses and drifts at the same rate, so that scientists see these long tails. How Gorelick can publish research on groundwater rates, but see no analogy to the larger issue of oil extraction seems such a waste of intellectual potential.

Should Gorelick ever read this review, I challenge him to read my work on dispersion and the math behind depletion of oil. These models come from solid math and probability underpinnings and simple physical first principles, and lead to the kind of insight that we all need to make sense of our fossil fuel energy situation.

Oil Discovery Simulation Reality

2010-06-09T16:38:00.000-07:00

I should have run this particular simulation long ago. In this exercise, I essentially partitioned the Dispersive Discovery model into a bunch of subvolumes. Each subvolume belongs to a specific prospecting entity, which I have given a short alias. The simulation assigns each one of the entities a random search rate and each one of the subvolumes also has a randomly sized value. The physical analogy equates to the prospector (i.e. the entity is an owner, leaser, company, nation, etc.) given their own subvolume (geographic location) to explore for oil. When they exhaustively search that subvolume, they end up with a cumulative amount of oil. The abstraction for subvolumes allows for the random sizing to directly translate to a proportional amount of oil. In general, bigger subvolumes equates to more oil but this does not have to hold, since the random rates blur this distinction.

Removing the technical mumbo-jumbo, the previous paragraph describes quite simply the context for the dispersive discovery model. Nothing about this description can possibly get misinterpreted as it essentially describes the process of a bunch of people systematically searching through a haystack for needles. Each person has varying ability and owns a varying size to search through, which essentially describes the process of dispersion.

The random number distributions derive from a mean search rate and a mean subvolume based on the principle of maximum entropy (MaxEnt). The number of subvolumes multiplied by the mean subvolume generates an ultimately recoverable resource (URR) total. By building a Monte Carlo simulation of this model, we can see how the discovery process plays out for randomly chosen configurations.

When the simulation executes, the search rates accelerate in unison so that the variance remains the same, maintaining MaxEnt of the aggregate. If I choose an exponential acceleration, the result turns precisely into the Logistic sigmoid, also known as the classic Hubbert Curve..

The entire simulation exists on a Google spreadsheet. Each row corresponds to a prospecting entity/subvolume pairing. The first two cells provide a random starting rate and a randomly assigned subvolume. As you move left to right across the row, you see the fraction of the subvolume searched increase in an accelerating fashion with respect to time. The exponential growth factor resides in cell A2. At some point in time, the accelerating search volume meets the fixed volume constraint and the number stops increasing. At that moment, the prospector has effectively finished his search. That subvolume has essentially ceased to yield newly discovered oil.

I reserve the 4th row for the summed values, the 3rd line generates the time derivative which plots out as a yearly discovery. The simulation "runs" one Monte Carlo frame at a time. We essentially see a full snapshot of one sample for about 150 years of dispersive search.

View Google Spreadsheet

I associated short names for each of the prospecting entities[1]. As I did not to want to make the spreadsheet too large, I limited it to 250 entities (which pushes Google to the limit for data). This of course introduces some noise fluctuations. The non-noisy solid line displays the analytical solution to the dispersive discovery model, which happens to match the derivative of the Logistic sigmoid.

The most important insight that we get from this exercise has to do with generating a BLINDINGLY SIMPLE explanation for deriving the Logistic behavior that most oil depletion analysts assume to exist, yet have no basis for. For crying out loud, I have seen children's board games with more complicated instructions than what I have given in the above paragraphs. Honestly, if you find someone that can't understand what it is going on from what I have written, don't ask them to play Chutes & Ladders either. Common sense Peak Oil theory ultimately reduces to this basic argument.

Contrast the elegance of the dispersive model with the most common alternative derivation for the logistic peak shape. This involves a completely misguided deterministic model that not surprisingly makes ABSOLUTELY NO SENSE. Whoever originally dreamed up the Verhulst derivation for ecological modeling and decided to apply it to Peak Oil must have consumed large quantities of mind-altering drugs prior to putting pencil to paper.

I also want to point out that what I did has nothing to do with multi-cycle Hubbert modeling which adds even less insight to the fundamental process.

I hope that this exercise helps in understanding the mechanism behind dispersive discovery. Seriously, the big intuitive sticking point that people have with the model has to do with the lack of any feedback mechanism in dispersive discovery. I imagine that engineers and most scientists get so used to seeing the feedback-derived Verhulst and LV equations derive the Logistic that they can't believe a simple and correct formulation actually exists!

In real terms, at some point the oil companies will cease to discover much of anything as they exhaust search possibilities. I suggest that they might want to consider making up for lost profit by licensing the oil discovery board game. This would help explain to their customers the reality of the situation.

UPDATE:
Occasionally Google does an underflow or overflow on some calculations so that the aggregate curve won't plot. The following animated GIF shows a succession of curves:

[1] I used shortened versions of TOD commenter names in the spreadsheet to make it a little more entertaining. I probably spent more time on writing the names down and battling the sluggishness of Google spreadsheet than I did on the simulation.

Predictably Unreliable

2010-06-08T21:07:00.000-07:00

I wrote about the unpredictably predictable nature of wind power in a few recent posts.

And of course we have watched the unexpected and unpredicted blow-out of the Deepwater Horizon oil well (the ultra-rare 1 out of 30,000 failure according to conventional wisdom) and hoping for the successful deployment of relief wells.

In the wind situation we know that it will work at least part of the time (given sufficient wind power, that is) without knowing precisely when, while in the second case we can only guess when a catastrophe with such safety-critical implications will occur.

We also have the unnerving situation of knowing that something will eventually blow-out, but with uncertain knowledge of exactly when. Take the unpredictability of popcorn popping as a trivial example. We can never predict the time of any particular kernel but we know the vast majority will pop.

In a recent episode that I went through, the specific failure also did not come as a surprise. I had an inkling that an Internet radio that I frequently use would eventually stop working. From everything I had read on-line, my Soundbridge model had a power-supply flaw that would eventually reveal itself as a dead radio. Previous customers had reported the unit would go bad anywhere from immediately after purchase to a few years later. After about 3 years it finally happened to my radio and the failure mode turned out exactly the same as everyone else's -- a blown electrolytic capacitor and a possible burned out diode.

The part obviously blew out because of some heat stress and power dissipation problem, yet like the popcorn popping, my interest lies in the wide range in failure times. The Soundbridge failure in fact looks like the classic Markov process of a constant failure rate per unit time. In a Markov failure process, the number of expected defects reported per day equate proportionally to how many units remain operational. This turns into a flat line when graphed as failure rate versus time. Customers that have purchased Soundbridges will continue to routinely report the failures for the next few years, with fewer and fewer reports as that model becomes obsolete.

Because of the randomness of the failure time, we know that any failures should follow some stochastic principle and likely that entropic effects play into the behavior as well. When the component goes bad, the unit's particular physical state and the state of the environment governs the actual process; engineers call this the physics of failure. Yet, however specific the failure circumstance, the variability in the component's parameter space ultimately sets the variability in the failure time.

So I see another way to look at failure modes. We can either interpret the randomness from the perspective of the component or from the perspective of the user. If the latter, we might expect that someone would abuse the machine more than another customer, and therefore effectively speed up its failure rate. Except for some occasional power-cycling this likely didn't happen with my radio as the clock stays powered in standby most of the time. Further, many people will treat their machine gingerly. So we have a spread in both dimensions of component and environment.

If we look at the randomness from a component quality-control perspective, certainly manufacturing variations and manual assembly plays a role. Upon internal inspection, I noticed the Soundbridge needed lots of manual labor to construct. Someone posting to the online Roku radio forum noticed a manually extended lead connected to a diode on their unit -- not good from a reliability perspective.

So I have a different way of thinking about failures which doesn't always match the conventional wisdom in reliability circles. In certain cases the result derives as expected, but in other cases the result diverges from the textbook solution.

Fixed wear rate, variable critical point: To model this to first-order, we assume a critical-point (cp) in the component that fails and then assume a distribution of the cp value about a mean. Maximum entropy would say that this distribution would approximate an exponential:

p(x) = 1/cp * exp(-x/cp)

The rate at which we approach the variable cp remains constant at R (everyone uses/abuses it at the same rate). Then the cumulative probability of failure is

P(t) = integral of p(x) from x=0 to x=R*t

This invokes the monotonic nature of failures by capturing all the points on the shortest critical path, and anything "longer" than the R*t threshold won't get counted until it fails later on. The solution to this integral becomes the expected rising damped exponential.

P(t) = 1 - exp(-R*t/cp)

Most people will substitute a value of τ for cp/R to make it look like a lifetime. This is the generally accepted form for the expected lifetime of a component to first-order.

P(t) = 1 - exp(-t / τ)

So even though it looks as if we have a distribution of lifetimes, in this situation we actually have as a foundation a distribution in critical points. In other words, I get the correct result but I approach it from a non-conventional angle.

Fixed critical point, variable rate: Now turn this case on its head and say that we have a fixed critical point and we have a maximum entropy variation in rate assuming some mean value, R.

p(r) = 1/R * exp(-r/R)

Then the cumulative integral looks like:

P(t) = integral of p(r) from r=cp/t to r=∞

Note carefully that the critical path in this case captures only the fastest rates and anything slower than the cp/t threshold won't get counted until later.

The result derives to

P(t) = exp(-cp/(R*t))

This has the characteristics of a fat-tail distribution because time goes into the denominator of the exponent, instead of the numerator. Physically, this means that we have very few instantaneously fast rates and many rates proceed slower than the mean.

Variable wear rate, variable critical point: In a sense, the two preceding behaviors act complementary to each other. So we can also derive P(t) for the situation whereby both the rate and critical point vary.

P(t) = integral of P(t | r)*p(r) over all r

This results in the exponential-free cumulative, which has the form of an entroplet.

P(t) = R*t/cp / (1+ R*t/cp) = t/τ/(1+t/τ)

Plotting the three variations side-by-side and assuming that τ=1, we get the following set of cumulative failure distributions. The full variant nestles in between the two other exponential variants, so it retains the character of more early failures (ala the bathtub curve) yet it also shows a fat-tail so that failure-free operation can extend for longer periods of time.

To understand what happens at a more intuitive level we define the fractional failure rate as

F(t) = dP/dt / (1-P(t))

Analysts use this form since it makes it more amenable to predicting failures on populations of parts. The rate then applies only to how many remain in the population, and the ones that have failed drop out of the count.

Only the first case above gives a failure rate that approaches the Markov ideal of constant rate over time. The other two dip below the constant rate of the Markov simply because the fat-tail cumulative requires a finite integrability over the time scale, and so the rates will necessarily stay lower.

Another post gives a full account of what happens when we generalize the first-order linear growth on the rate term, letting R=g(t). The full variant ultimately gives dg/dt / (1+g(t)), so that if g(t) starts rising we get the complete bathtub curve.

If we don't invoke other time dependencies on the rate function g(t), we see how certain systems never show failures after an initial period. Think about it for a moment -- the fat-tails of the variable rate cases push the effective threshold for failure further and further into the future.

In effect, normalizing the failures in this way explains why some components have predictable unreliability, while other components can settle down and seemingly last forever after the initial transient.

I discovered that this paper by Pandey jives with the way I think about the general problem.

Enjoy your popcorn, it should have popped by now.

Reliability of Relief Wells

2010-06-06T14:04:00.000-07:00

I have seen much discussion on TOD and elsewhere of the effectiveness of adding relief wells to take the pressure off the failed well in the Gulf. Occasionally I have noticed questions on how one would make a kind of reliability prediction given estimated success/failure probability numbers. This turns into the classic redundant configuration reliability prediction problem.

Initially, for pure success probabilities I wouldn't add time to the equation. In the steady-state we just work with basic probability multiplications. If the probabilities of success rates remain independent of each other, then they form a pattern. Say we have three tries for relief wells, each one having a value between 0 and 1. If all three fail then the whole attempt failed:

P(failure) = P1(failure)*P2(failure)*P3(failure)

and

P(success)=1-P(failure)

so if P1=P2=P3=1-0.7=0.3

then P(failure)=0.027

and P(success)=0.973

With time you need to work from the notion of a deadline, i.e. that no failures occur in a certain amount of time. Otherwise you end up using the fixed probabilities above because you have essentially infinite time to work with.

Apart from end-state failure analysis, you can also do a time-averaged effectiveness, where the rates help you do a trade-off analysis between how long it takes before you fix the problem and how much oil gets released in the meantime. Unfortunately, when you look at the optimization criteria, the only valid optimum in most people's minds is to stop the oil leak as quickly as possible. Otherwise it looks like we play dictator rolling dice (at least IMO that is the political response I predict to get).

Given that political issue, you can create a set of criteria with weights on the probabilities of success, the cost, and on the amount of oil leaked (the first and third as Markov models as a function of time). When you combine the three and look for an optimum, you might get a result that gives you a number of relief wells somewhere between 1 and infinity. The hard part remains establishing the weighting criteria. Place a lower weight on cost and you will definitely lower the number of wells. And that's where the politics plays in again, as many people will suggest that cost does not form a limitation. We also have the possibility of a massive blow-out by adding a botched relief well, but that risk may turn out acceptable.

Below I show a state diagram from a Markov-based reliability model. With the Markov process premise you can specify rates of probability flow between various states and then execute it without having to resort to Monte Carlo.

I made this diagram for 3 relief wells drilled in succession, when one doesn't work, then we start the next. The term B1 is the rate for a failure specified as 0.01 (or 1 in 100 days). B2 is a success rate of 0.02 (or 1 in 50 days). The start state is P1, the success state is P3, and the end failure state is P5.

When I execute this for 200 days, the probability of entering state P5 is 3.5% and it will rise to 3.7% after 1000 days. P3 is 95% after 200 days. The sanity check on this gives a success ratio of about 0.02/(0.01+0.02)=0.666 and from the formula this gives a probability of failure at the end state of (1/3)^3 = 0.037 = 3.7%. This sanity checks with the output after 1000 days.

The Markov model allows you to predict the time dependence of success and failure based on the assumptions of the individual non-redundant failure rates. You can thus work the model as a straightforward reliability prediction. Change the success probabilities to 50% individual success rate and we still only need three relief wells if we want to get to 87.5% . Contrast that to 97% average success rate with 3 wells, if we remain on the optimistic side of 50%. So you can see that our confidence grows with the confidence in the success of the individual wells, which makes intuitive sense.

This particular model assumes a serial succession of relief wells. You can also model relief wells constructed in parallel, which I believe remains the current strategy in the Gulf. Or you can model the initial delay a little better. With the model as described, we have success rates that can occur earlier than perhaps expected. An exponential on the success rate per time provides a distribution where the standard deviation equals the mean, which is the most conservative estimator should you have no idea what the standard deviation is. To generate a model with about half the standard deviation, we can turn the exponential into a gamma. Each relief well spends about half its time in a "build" stage where it experiences neither success or failure. Then the next stage of its life-cycle gets spent in testing for success. See the following chart:

The overall result doesn't differ much from the previous model but you do see a much diminished success rate early on -- which makes the model match reality better.

As another possibility, we can repeat an individual relief well several times, backing up and retrying if the last one doesn't work. That models as a state that directs back on itself, with a rate B4. I won't run this one because I don't know the rates of retries, but the general shape of the of the failure/success curve looks similar.

I'm sure some group of analysts somewhere has worked a similar kind of calculation. Whether it pays off or not for a single case, I can't really say. However, this kind of model effectively describes how the probabilities work out and how you can use a state diagram to keep track of failure and success transitions.

By the way, this same math goes into the Oil Shock Model which I use for oil production prediction. In the oil shock model, transitions describe the rates between oil production life-cycle states, such as construction and maturation and extraction. So both the reliability model and the Oil Shock model derive from probability-based data flow models. This kind of model works very well for oil production because we have a huge number of independently producing regions around the world and the law of large numbers makes the probability projections that much more accurate. As a result, I would put more trust in relying on the results of the oil shock model than predicting the success of the recovery of a single failed deep-water production well. Yet, the relief well redundancy model does help to estimate how many extra relief wells to add and adds some quantitative confidence to one's intuition.

Based on the post by Joules Burn (JB) on TOD BP's Deepwater Oil Spill: A Statistical Analysis of How Many Relief Wells Are Needed, I added a few comments:

JB did everything perfectly correctly given the premises. Another way to look at it is that you need to accomplish a sequence of steps, each with a probability rate of entering into the next state. This would simulate the construction of the relief well itself (a sequence of steps). Then you would have a rate into a state where you start testing the well for success. This goes into a state that results in either a success, retry, or failure (the utter failure in JB lingo). The convenient thing is that you can draw the retry as a feedback loop, so the result looks like the following for a single well:

I picked some of the numbers from intuition, but the results have the general shape that JB showed. When you look at a rate like 0.1, inverting it gives a mean transition of 10 days.

This is a state diagram simulation like that used in the Oil Shock model, which I use to project worldwide oil production. I find it interesting to see how well accepted the failure rate approach is for failure analysis, but few seem to accept it for oil depletion analysis. I presume oil depletion is not as mission critical a problem as the Gulf spill is :)

Thermal Entropic Dispersion

2010-06-05T12:31:00.000-07:00

As we learn how to extract energy from disordered, entropic systems such as amorphous photovoltaics and wind power, we can really start thinking creatively in terms of our analysis. Most of the conventional thinking goes out the window as considerations of the impact of disorder requires a different mindset.

In a recent post, I solved the Fokker-Planck diffusion/convection equation for disordered systems and demonstrated how well it applied to transport equations; I gave examples for both amorphous silicon photocurrent response and for the breakthrough curve of a solute. Both these systems feature some measurable particle, either a charged particle for a photovoltaic or a traced particle for a dispersing solute.

Similarly, the conduction of heat also follows the Fokker-Planck equation at its most elemental level. In this case, we can monitor the temperature as the heat flows from regions of high temperature to regions of low temperature. In contrast to the particle systems, we do not see a drift component. In a static medium, not abetted by currents (as an example, mobile ground water) or re-radiation, heat energy will only move around by a diffusion-like mechanism.

We can't argue that the flow of heat shows the characteristics of an entropic system -- after all temperature serves as a measure of entropy. However, the way that heat flows in a homogeneous environment suggests more order than you may realize in a practical siuation. In a perfectly uniform medium, we can propose a single diffusion coefficient, D, to describe the flow or flux. A change of units translates this to a thermal conductivity. This value inversely relates to the R-value that most people have familiraity with when it comes to insulation.

For particles in the steady state, we think of Fick's First Law of Diffusion. For heat conduction, the analogy is Fourier's Law. These both rely on the concept of a concentration gradient, and functionally appear the same, only the physical dimensions of the parameters change. Adding the concept of time, you can generalize to the Fokker-Planck equation (i.e Fick's Second Law or the Heat Equation respectively).

Much as with a particle system, solving the one-dimensional Fokker-Planck equation for a thermal impulse you get a Gaussian packet that widens from the origin as it diffuses outward. See the picture to the right for progressively larger values of time. The cumulative amount collected at some point, x, away from the origin results in a sigmoid-like curve known as an complementery error function or erfc.

Yet in practice we find that a particular medium may show a strong amount of uniformity. For example, earth may contain large rocks or pockets which can radically alter the local diffusivity. Same thing occurs with the insulation in a dwelling; doors and windows will have different thermal conductivity than the walls. The fact that reflecting barriers exist means that the effective thermal conductivity can vary (similarly this arises in variations due to Rayleigh scattering in wind and wireless observations). I see nothing radical about the overall non-uniformity concept, just an acknowledgment that we will quite often see a heterogeneous environment and we should know how to deal with it.

Previously, I solved the FPE for a disordered system assuming both diffusive and drift components. In that solution I assumed a maximum entropy (MaxEnt) distribution for mobilities and then tied diffusivity to mobility via the Einstein relation. The solution simplifies if we remove the mobility drift term and rely only on diffusivity. The cumulative impulse response to a delta-function heat energy flux stimulus then reduces to:

T(x,t) = T1* exp(-x/sqrt(D*t)) + T0

No erfc in this equation (which by the way makes it useful for quick analysis). I show the difference between the two solutions in the graph to the right (for a one-dimensional distance x=1 and a scaled diffusivity of D=1). The uniform diffusivity form (red curve) shows a slightly more pronounced knee as the cumulative increases than the disordered form (blue curve) does. The fixed D also settles to an asymptote more quickly than the MaxEnt disordered D does, which continues to creep upward gradually. In practical terms, this says that things will heat up or slow down more gradually when a variable medium exists between yourself and the external heat source

Because of the variations in diffusivity, some of the heat will also arrive a bit more quickly than if we had a uniform diffusivity. See the figure to the right for small times. Overall the differences appear a bit subtle. This has as much to do with the fact that diffusion already implies disorder, while the MaxEnt formulation simply makes the fat-tails fatter. Again it essentially disperses the heat -- some gets to its destination faster and a sizable fraction later.

Which brings up the question of how we can get some direct evidence of this behavior from empirical data. With drift, the dispersion becomes much more obvious, as systems with uniform mobility with little disorder show very distinct knees (ala photocurrent time-of-flight measurements or solute breakthrough curves for uniform materials) . Adding the MaxEnt variation makes the fat-tail behavior very obvious, as you would observe from the anomalous transport behavior in amorphous semiconductors. With diffusion alone, the knee automatically smears, as you can see from the figure to the right for a typical thermal response measurement.

Evidence
Much of the interesting engineering and scientific work in characterizing thermal systems comes out of Europe. This paper investigating earth-based heat exchangers contains an interesting experiment. As a premise, they wrote the following, where incidentally they acknowledge the wide variation in thermal conductivities of soil:

The thermal properties can be estimated using available literature values, but the range of values found in literature for a specific soil type is very wide. Also, the values specific for a certain soil type need to be translated to a value that is representative of the soil profile at the location. The best method is therefore to measure directly the thermal soil properties as well as the properties of the installed heat exchanger.
This test is used to measure with high accuracy:

The temperature response of the ground to an energy pulse, used to calculate:
the effective thermal conductivity of the ground
the borehole resistance, depending on factors as the backfill quality and heat exchanger construction

The average ground temperature and temperature - depth profile.
Pressure loss of the heat exchanger, at different flows.

The authors of this study show a measurement for the temperature response to a thermal impulse, with the results shown over the course of a couple of days. I placed a solid red and blue line indicating the fit to an entropic model of diffusivity in the figure below. The mean diffusivity comes out to D=1.5/hr (with the red and blue curves +/- 0.1 from this value) assuming an arbitrary measurement point of one unit from the source. This fit works arguably better than a fixed diffusivity as the variable diffusivity shows a quicker rise and a more gradual asymptotic tail to match the data.

The transient thermal response tells us a lot about how fast a natural heat exchanger can react to changing conditions. One of the practical questions concerning their utility arises from how quickly the heat exchange works. Ultimately this has to do with extracting heat from a material showing a natural diffusivity and we have to learn how to deal with that law of nature. Much like we have to acknowledge the entropic variations in wind or cope with variations in CO2 uptake, we have to deal with the variability in the earth if we want to take advantage of our renewable geothermal resources.

Wind Variability in Germany

2010-06-01T18:13:00.000-07:00

By adding more data to the post on wind dispersion, we can observe how dispersion in wind speeds has a universal character. I picked up the previous data set from several years worth of output from Ontario. This new set hails from northwest Germany and this site (thanks to globi for the link). The data consists of wind power collected at 15 minute intervals.

Note that the same entropic dispersion holds as for Ontario (see graph to the right). Both curves display the same damped exponential probability distribution function for frequency of wind power (derived from wind speed). We also see the same qualitative cut-out above a certain power or wind energy level. As I said previously, we don't gain much by drawing from these higher power levels as they occur more sporadically than the nominally rated wind speeds at the upper reaches of the curve.

The following figure gives an explanation for the cutout above the "max" wind speed. Globi also provided this PDF from Vestas, a maker of wind turbines. The end of the document has the complete spec.

Power regulation : pitch regulated with variable speed
Operating data
Rated power : 3,000 kW Cut-in wind speed : 3 m/s
Rated wind speed : 12 m/s
Cut-out wind speed : 25 m/s

Too many people get the idea that the sporadic nature of wind confronts us with some kind of "problem". We will have to get used to a different way of thinking about wind. The entropic dispersion of wind acts much like a variation of the Carnot cycle. In the Carnot cycle of engine efficiency, we have to live with a maximum level of energy conversion based on temperature differences of the input and output reservoirs. With wind, the earth's environment and atmosphere provides the temperature differences which leads directly to the variability over time.

Which leads to the fact that WITH WIND POWER, WE CAN ACHIEVE VERY HIGH USAGE EFFICIENCY GIVEN THE ENTROPIC CHARACTERISTICS OF THE WIND. I put this in upper case because it amounts to a law of nature. We need to talk about efficiencies within the constraints of the physical laws just as with the Carnot cycle. We will observe intermittency as a result of entropic dispersion and we have to get used to it. We should not call it a fundamental "problem", as we cannot change the characteristics of entropy (apart from adding energy, and that just moves us back to square one).

Other people would suggest that the fundamental problem with farming derives from the intermittent nature of the rain. With farming, we adapt -- likewise with wind energy. Instead of a problem, we need to call it an opportunity.

As a blast from the past check out my expose of the forged video editing by the George Bush marketing team against John Kerry. Wind energy advocates will have to watch out for these tactics as the right-wingers will project and frame any way they can to make wind look like a wimpy exercise designed by the elite for the elite.

The Word on Dispersion

2010-05-29T20:15:00.000-07:00

Credit the Gulf oil disaster with allowing the words dispersion and dispersants to enter our common vocabulary. In the context of the spill, the use of dispersants on the oil causes the potentially sticky coagulating oil to split apart into finer granularity drops and somehow make it more amenable to breaking down. Dispersion in terms of a chemical definition simply means spreading out particles in the medium, in this case seawater. So a dispersant breaks it up and dispersion scatters it about.

The BP team apparently wanted to break up the oil up so that it could easily migrate and essentially dilute its strength within a larger volume. So instead of allowing a highly concentrated dose of oil to impact a seashore or the ocean surface, the dispersants would force the oil to remain in the ocean volume, and let the vast expanse of nature take its course. Somebody in the bureaucratic hierarchy made the calculated decision to apply dispersants as a judgment call. I can't comment on the correctness of that decision but I can expound on the topic of dispersion, which no one seems to fully understand, even in a scientific context.

As the media has forced us to listen to made up technical terms such as "top kill", "junk shot", and "top hat" which describe all sorts of wild engineering fixes, I will take a turn toward the more fundamental notions of disorder, randomness, and entropy to explain that which we cannot necessarily control. I always think that if we can understand concepts such as dispersion from first principles, we actually have a good chance of understanding how to apply it to a range of processes besides oil spill dispersal. In other words, well beyond this rather specific interpretation, we can apply the fundamentals to other topics such as green-house gases, financial market fluctuations, and oil discovery and production, amongst a host of other natural or man-made processes. Really, it is this fundamental a concept.

Background

If by the process of dispersion we want the particles to dilute as rapidly as possible, we need to somehow accelerate the rate or kinetics of the interactions. This becomes a challenge of changing the fundamental nature of the process, via a homogeneous change, or by introducing additional heterogeneous pathways that provide alternate pathways to faster kinetics. From this perspective, dispersion describes a mechanism to divergently spread-out the rates and dilute the material from its originally concentrated form. One can analogize in terms of a marathon race; the initial concentration of runners at the starting line rapidly disperses or spreads out as the faster runners move to the front and the slower runners drop to the rear. In a typical race, you see nothing homogeneous about the makeup of the runners (apart from their human qualities); the elites, competitive amateurs, and spur-of-the-moment entrants cause the dispersion. Whether we want to achieve a homogeneous dispersion or not, we have to account for the heterogeneous nature of the material. In other words, we rarely deal with pure environments so have to solve for much more than the limited variability we originally imagined. Generalizing from the rather artificial constraints of a marathon race, dispersion in other contexts (such as crystal growth or reservoir growth) results from an increase of disorder as a direct consequence of entropy and the second law of thermodynamics.

In terms of the spread in dispersion, we might often observe a tight bunching or a wide span in the results. The wider dispersion usually indicates a larger disorder, variability, or uncertainty in the characteristics -- a "fat-tail" to the statistics so to speak. So when we introduce a dispersant into the system, we add another pathway and basically remove order (or introduce disorder) into the system. Dispersion may thus not accelerate a process in a uniform manner, but instead accelerates the differences in the characteristic properties of the material. This again describes an entropic process, and we have to add energy or find exothermic pathways to fight the tide of increasing disorder.

This seems like such a simple concept, yet it rarely gets applied to most scientific discussions of the typical disordered process. Instead, particularly in an academic setting, what one usually reads amounts to pontificating about some abnormal or anomalous kind of random-walk that must occur in the system. The scientists definitely have a noble intention -- that of explaining a fat-tail phenomenon -- yet they don't want to acknowledge the most parsimonious explanation of all. They simply do not want to consider heterogeneous disorder as described by the maximum entropy principle.

Figure 1: Difference between a classical random walk (left) and an anomalous random walk (right). The salient difference is that occasional long jumps (Levy flights) occur in the anomalous random walk. A much simpler approach admits that a heterogeneous nix of random walkers of different rates exists. This will give essentially the same observable outcome without resorting to arcane mathematical modeling.

The complicating factor in discussions about dispersion involves the intuitively related concept of diffusion and convection or drift. Diffusion also derives from the statistics of disorder and describes how particles can spontaneously spread out without a real driving force, apart from the uniform environment, for example from the thermal background. The analysis of a particle undergoing random walk leads directly to the concept of diffusion. Random walk ideas seem to intrigue mathematicians and scientists because it places the concept of diffusion into a real concrete representation. In some sense everyone can relate to the idea of a particles bouncing around, but not necessarily to the idea of a gradient in concentration.

Convection and drift describe the motion of particles under an applied force, say charged particles under the influence of an electric field (Haynes-Shockley), or of solute or suspended particles under the influence of gravity (Darcy's Law). This essentially describes the typical constant velocity, akin to a terminal velocity, that we observe in a pure semiconductor (Haynes-Shockley) or a uniformly porous media (Darcy's).

Dispersion can effect both diffusion and drift, and that establishes the premise for the novel derivation that I came up with.

Breakthrough

The unification of the dispersion and diffusion concepts could have a huge influence on the way we think about practical systems, if we could only factor the mathematics describing the process. I can straightforwardly demonstrate a huge simplification assuming a single somewhat obvious premise. This involves applying the conditions of maximum entropy, by essentially maximizing disorder under known constraints or moments (i.e. mean values, etc).

The obviousness of this unifying solution contrasts with my lack of awareness of of any such similar simplification in the scientific literature. Surprisingly, I can't even confirm that anyone has really looked into the general idea. So far, I can't find any definitive work on this unification and little interest in pursuing this premise. Stating my point-of-view flatly, the result has such a comprehensive and intuitive basis that it should have a far-reaching impact on how we think about dispersion and diffusion. It just needs to gain a foothold of wider acceptance in the marketplace of ideas.

Which brings up a valid point I have heard directed my way. From my postings on TheOilDrum.com, commenters occasionally ask me why I don't publish these results in an academic setting, such as a journal article. To answer that, journals have evidently failed in this case, as I never find any serious discussion of dispersion unification. So consider that even if I submitted these ideas to a journal, it may just sit there and no one would ever apply the analysis in any future topics. This makes it an utterly useless and ultimately futile exercise. I will risk putting the results out on a blog and take my chances. A blog easily has as much archival strength, much more rapid turnaround, the potential for critiquing, and has searchability (believe it or not, googling the term "dispersive transport" yields this blog as the #3 result, out of 16,200,000). The general concepts do not apply to any specific academic discipline apart perhaps applied math, and I certainly won't consider publishing the results in that arena with out risking it disappear without a trace. Eventually, I want to place this information in a Wikipedia entry and see how that plays out. I would call it an experiment in Open Source science.

But that gets a little ahead of the significance of the current result.

The Unification of Diffusion and Drift with Dispersion

As my most recent post described, solving the Fokker-Planck equation (FPE) under maximum entropy conditions provides the fundamental unification between dispersion, diffusion and drift. For fans of Taleb and Mandelbrot, this shows directly how "thin-tail" statistics become "fat-tail" statistics without resorting to fractal arguments.

The Fokker-Planck equation shows up in a number of different disciplines. Really, anything having to do with diffusion or drift has a relation to Fokker-Planck. Thus you will see FPE show up in its various guises: Convection-Diffusion equation, Fick's Second Law of Diffusion, Darcy's Law, Navier-Stokes (kind of), Shockley's Transport Equation, Nernst-Planck; even something as seemingly unrelated as the Black-Scholes equation for finance has applicability for FPE (where the random walk occurs as fractional changes in a metric).

Because of its wide usage, the FPE tends to take the form of a hammer, where everything it applies to acts as the nail. (You don't see this more frequently than in finances, where Black-Scholes played the role of the hammer) Since the solution of FPE results in a probability distribution, it gives the impression that some degree of disorder prevails in the system under study. I find this understandable since the concept of diffusion implies an uncertainty exactly like a random walk shows uncertainty. In other words, no two outcomes will turn out exactly the same. Yet, in mathematical terms, the measurable value associated with diffusion, the diffusion constant D, has a fixed value for random motion in a homogeneous environment. When the parameters actually change, you enter in the world of stochastic differential equations; I won't descend to deeply into this area, only to apply this as a basic concept. The diffusion and mobility parameters have a huge variability that we have yet adequately accounted for in many disordered systems.

For that reason, the FP equation really applies to ordered systems that we can characterize well. Not surprisingly the ordinary solution to FPE gives rise to the conventional ideas of normal statistics and thin-tails.

So for phenomenon that appear to depart from conventional normal diffusion (the so-called anomalous diffusion) we have two distinct camps and corresponding solution paths to choose from. The prevailing wisdom suggests that an entirely different kind of random walk occurs (Camp 1). No longer does the normal diffusion apply, giving rise to normal statistics; instead we get the statistics of fat-tails and random walk trajectories called Levy flights to concretely describe the situation (see Figure 1). The mathematics quickly gets complicated here and most of the results get cast into heuristic power-laws. It takes a leap of faith to follow these arguments.

The question comes down to whether we wish to ascribe anomalous diffusion as a strange kind of random walk (Camp 1) or simply suggest that heterogeneity in diffusional and drift properties adequately describes the situation (Camp 2). I take the stand in the latter category and stand pretty much alone in this regard. Find some academic research article on anything related to anomalous diffusion and very few will accept the most parsimonious explanation -- that a range of diffusion constants and mobilities explain the results. Instead the researcher will punt and declare that some abstract Levy flight describes the motion. Above all I would rather think in practical terms, and simple variability has a very pragmatic appeal to it.

I went through the derivation of the dispersive FPE solution for a disordered semiconductor in the last post, and want to generalize it here. This makes it especially applicable to notions of transport physical transport of material in porous matter. This would include the motion of oil underground, CO2 in the air, and perhaps even spilled oil at sea.

In the one-dimensional model of applying an impulse function of material, the concentration n will disperse according to the following equation:

n(x, z) = (z + sqrt(zL + z^2)/sqrt(zL + z^2)*exp(-2x/(z + sqrt(zL + z^2))

where
z= μFt
L = β/F

The term z takes the place of a time-scaled distance, which can speed up or slow down under the influence of a force F (i.e. gravity, or electric field for a charged particle). The characteristic distance L represents the effect of the stochastic force β (aka Boltzmann's constant) and ties in the diffusional aspects of the system. The specific parameterization of the exponential results in the fat-tail observed.

In the past, I had never gone through the trouble of solving the FPE, simply because intuition would suggest that the dispersive envelope would cancel out most of the details of the diffusion term. In the dispersive transport model that I originally conceived, the dispersion would at most follow the leading wavefront of the drifting diffusional field as "sqrt(Lz+z^2)" as described here or as "sqrt(Lz)+z" here.

I estimated that the diffusion term would follow as the square root of time according to Fick's first law and that drift would follow time linearly, with only an idea of the qualitative superposition of the terms in my mind.

As one might expect, the actual entropic FPE solution borrowed from a little of each of my estimates, essentially averaging between the two:

(z + sqrt(zL + z^2))/2

So the solution to the dispersive FPE form for a disordered system turns out entirely intuitive , and one can almost generate the result from inspection. The difference between the original entropic dispersion derivation and the full FPE treatment amounts to a bit of pre-factor bookkeeping in the first equation above. You can see this by comparing the two approaches for the case of L=1 and unity width for the dispersive transport current model.

Figure 2: Differences between the original entropic dispersive model and the fully quantified FPE solution will converge as L gets smaller.

Dispersive Transport in Porous Media.

The above solved equations can actually apply directly as solutions to Darcy's law when it comes to describing the flow of material in a disordered porous media. I suppose this will irk the petroleum engineers, hydrologists, and geologists out there who have long sought the solution to this particular problem.

Yet we should not act surprised by this result. The actions of multiple processes acting concurrently on a mobile material will generally result in a universal form governed by maximum entropy. It doesn't matter if we model carriers in a semiconductor or particles in a medium, the result will largely look the same. In a hydraulic conductivity experiment, Lange treated the breakthrough curve of a trace element through a natural catchment as a FPE convection-dispersion model, and came up with the same results independent of the fractionation of the media.

By applying the simple dispersion model (blue curve below) to Lange's results, one sees that an excellent fit results with the fat-tail exactly following the hyperbolic decline that reservoir engineers often see in long-term flow behavior. This could includes the time dependent emptying of the currently leaking deep sea Gulf reservoir!

Figure 3: Breakthrough curve of a traced material showing results from an entropic dispersion model in blue.

Moreover, the amount of diffusion that occurs appears quite minimal. Adding a greater proportion of diffusion by increasing L does not improve the fit of the curve (see the chart to the right). Just as in the semiconductor case, the shape has a significant meaning when analyzed from the perspective of maximum entropy.

Nothing complicated about this other than admitting to the fact that heterogeneous disordered systems appear everywhere and we have to use the right models to characterize their behavior.

The details of this experiment are described in the following papers:

D.Haag and M.Kaupenjohann, Biogeochemical Models in the Environmental Sciences: The Dynamical System Paradigm and the Role of Simulation Modeling
H. Lange, Are Ecosystems Dynamical Systems?

The authors of these papers have mixed feelings about the applicability of modeling biogeochemical systems and speculate whether we should use any kinds of models for "ecological risk assessment". They point out that ecological systems obviously can adapt under certain circumstances and no amount of physical modeling can predict which way the system will go. Will spilled oil decompose faster as the environment adapts around it? Will that make dispersion less relevant? Who knows?

Still the work of modeling the physical process alone has enormous value as Haag and Kaupenjohann point out:

Despite not being a ‘real’ thing, "a model may resonate with nature" (Oreskes et al. 1994) and thus has heuristic value, particular to guide further study. Corresponding to the heuristic function, Joergensen (1995) claims that models can be employed to reveal ecosystem properties and to examine different ecological theories. Models can be asked scientific questions about properties. According to Joergensen (1994), examples for ecosystem properties found by the use of models as synthesizing tools are the significance of indirect effects, the existence of a hierarchy, and the ‘soft’ character of ecosystems. However, we agree with Oreskes et al. (1994) who regard models as "most useful when they are used to challenge existing formulations rather than to validate or verify them". Models, as ‘sets of hypotheses’, may reveal deficiencies in hypotheses and the way biogeochemical systems are observed. Moreover, models frequently identify lacunae in observations and places where data are missing (Yaalon 1994).
As an instrument of synthesis (Rastetter 1996), models are invaluable. They are a good way to summarize an individual research project (Yaalon 1994) and they are capable of holding together multidisciplinary knowledge and perspectives on complex systems (Patten 1994).
While models as a product may have heuristic value, we would like to emphasize also the role of the modeling process: "[…] one of the most valuable benefits of modeling is the process itself. These benefits accrue only to participants and seem unrelated to the character of the model produced" (Patten 1994). Model building is a subjective procedure, in which every step requires judgment and decisions, making model development ‘half science, half art’ and a matter of experience (Hoffmann 1997, Hornung 1996). Thus modeling is a learning process in which modelers are forced to make explicit their notions about the modeled system and in which they learn how the analytically isolated components of a system can be ‘glued’ (Paton 1997). As modeling mostly takes place in groups, modeling and the synthesis of knowledge has to be envisaged as a dynamic communication process, in which criteria of relevance, the meaning of terms, the underlying concepts and theories, and so forth are negotiated. Model making may thus become a catalyst of interdisciplinary communication.
In the assessment of environmental risks, however, an exclusively scientific modeling process is not sufficient, as technical-scientific approaches to ‘post-normal’ risks are unsatisfactory (Rosa 1998) and as the predictive capacity and operational validity of models (e.g. for scenario computation) is in doubt. The post-normal science approach (Funtowicz & Ravetz 1991, 1992, 1993) takes account of the stakes and values involved in environmental decision making. Following a ‘post-normal’ agenda, model development and model validation for risk assessment should become a trans-scientific (communication) task, in which "extended peer communities" participate and in which non-equivalent descriptions of complex systems are made explicit, negotiated, and synthesized. In current modeling practice, however, models are highly opaque and can rarely be penetrated even by other scientists (Oreskes, personal communication). As objects of communication, models still are closed systems and black boxes.

We need to really take up the charge on this as our future depends on understanding the role of entropy in nature. For too long, we have not shown the intellectual curiosity to model how much oil we have underground, what size distribution the reservoirs take, and how fast that they can epmty, even though some perfectly acceptable models can describe this statistically, using dispersion no less!

Now that the Macondo oil has discovered an escape hatch and has gone disordered on us and will go who-knows-where, it seems we can really make some headway in our common understanding. Nothing like having your feet in the fire.

Fokker-Planck for Disordered Systems

2010-05-24T18:59:00.000-07:00

To get the cost of photovoltaic (PV) systems down, we will have to learn how to efficiently use crappy materials. By crap I mean that mass-produced PV materials will end up getting rolled or extruded or organically grown. Unless we perfect the process, most everything will turn out non-optimal. We already know the difference between clean-room cultivated single crystal semiconducting material and the defect-ridden and often amorphous materials that nature and entropy drives us to. For performance sensitive applications such as communications and computing we would only rarely consider disordered material as a candidate semiconductor. Certainly, the performance of these materials makes them unlikely candidates for high speed processing -- yet for solar cell applications, they may serve us well. In the end, we just have to learn how to understand and deal with crap.

The following will revisit a couple of previous posts where I outlined a novel way to analyze the behavior of disordered semiconducting material. I know for certain that no one has proposed the particular approach before. If it does exist, I certainly can't find it in the literature. From one perspective, this analysis sets forth a baseline for the characterization of a maximally disordered semiconductor.

Background

The prehistoric 1949 Haynes-Shockley experiment first measured the dynamic behavior of charged carriers in a semiconducting sample. It basically confirmed the solution of the diffusion (Fokker-Planck) equation and it demonstrated diffusion, drift, and recombination in a conceptually simple setup. This animated site gives a very interesting overview of PV electrical behavior.

Figure 1: Apparatus for the Haynes-Shockley experiment

This setup works according to theory for an ordered semiconductor with uniform properties but apparently gets a bit unwieldy for any disordered or non-uniform material sample. I inferred this as conventional wisdom since most scientists either punt or use heuristics partially derived from the inscrutable work of a select group of random-walk theorists (see Scher & Montroll).

I had previously applied a very straightforward interpretation to the problem of carrier transport in disordered material. My dispersion analysis essentially set aside the Fokker-Planck formalism for a mean value approximation where I tactically applied the Maximum Entropy Principle. In particular, I really like the MaxEnt solution because I can recite the solution from memory. It matches intuition in a conceptually simple way once you get into a disordered mind-set.

In the real Haynes-Shockley experiment, a pulse gets injected at one electrode, and a nearly pure time-of-flight (TOF) profile results. The initial pulse ends up spreading out in width a bit, but the detected pulse usually maintains the essential Gaussian sigmoid shape.

Adding Disorder

For the time-of-flight for a disordered system, the Maximum Entropy solution looks like:

q(t) = Q * exp(-w/(sqrt((μEt)² + 2Dt))

This essentially states that the expected amount of charge accumulated at one end of the sample (at a distance w) at time t, follows a maximum entropy probability distribution. The varying rates described by μ and D disperse the speed of the carriers so that a broadened profile results from the initial pulse spike.

The equation above formed the baseline for the interpretation I described initially here.

For completeness, I figured to test my luck and see if I can bull my way through the basic diffusion laws. If I could produce an equivalent solution by applying the Maximum Entropy Principle directly to the Fokker-Planck equation, then this would give a better foundation for the "inspection" result above.

The F-P diffusion equation gets expressed as a partial differential equation with a conservation law constraint:

In this case D1=μ* (carrier mobility) and D2=D* (diffusion coefficient), and f(x,t)=n(x,t) (carrier concentration). With recombination, the solution in one-dimension looks like:

This of course works for well-ordered semiconductors, but D* and μ* will likely vary for disordered material. I made the standard substitution via the Einstein Relation for

D* = V_t μ*

where V_t = β/q stands for the chemical or thermal potential at equilibrium (usually β equals kT where k is Boltzmann's constant and T is absolute temperature). At equilibrium, the stochastic force of diffusion exactly balances the electrostatic force F = qE.

From the basic physics, we can generate a maximum entropy density function for D

p(D*) = 1/D * exp(-D*/D)

then

n(x,t) = Integral p(D*) * n_mean(x,t) over all D*

This looks hairy but the integral comes out straightforwardly as (ignoring the constant factors)

n(x,t) = 1/sqrt(t*(4D+t*(Eμ)²)) * exp(-x*R(t)) / R(t)

where

R(t) = sqrt(1/(Dt) + E/(2V_t)²) - E/(2V_t)

If we evaluate this for carriers that have reached the drain electrode at x=w, the total charge collected q is:

q(t) = Q/sqrt(t*(4D+t*(Eμ)²) * exp(-w*R(t)) / R(t)

The measured current is

I(t) = mean of dq(t)/dt from 0 to w

The simple entropic dispersive expression and the Fokker-Planck result obviously differ in their formulation, yet the two show the same asymptotic trends. For an arbitrary set of parameters, one can't detect a practical difference. Use whichever you feel comfortable with.

I show the dynamics of the carrier profile in the animated GIF to the right. The initial profile starts with a spike at the origin and then the profile broadens as the mean starts drifting and diffusing to the opposing contact. You don't see much from this perspective as it looks completely like mush. Yet, when plotted on a log-log scale, it does take on more character.

The collected current profile looks like the following

Figure 2: Typical photocurrent trace showing the initial diffusional spike, a plateau for relatively constant collection from the active region, and then a power-law tail produced from the entropic drift dispersion.

Organic Semiconductor Applications

The photocurrent profile displayed above came from from Andersson's "Electronic Transport in Polymeric Solar Cells and Transistors" (2007) wherein he analyzed the transport in a specific organic semiconducting material, the polymer APFO.

The blue line drawn through the set of traces follows the entropic dispersion formulation. The upper part of the curve describes the diffusive spike while the lower part generates the fat-tail due to the drift component (this shows an inverse square power law in the tail).

Figure 3: Universal profile generated over a set of applied electric field values. For this set, scaling of transit time with respect to the applied field holds, indicative of a constant mobility. However, carrier diffusion causes the initial transient and this does not scale, as the electric field has no effect on diffusion, as shown in the lower set of blue curves.

As I stated in the previous post, most scientists when discussing this shape have either (1) referred to Scher/Montroll and the vague heuristic α, (2) dismissed these features, or (3) labelled them as uninteresting. Andersson follows suit:

At best this transient, as the high α value indicates, might be possible to evaluate in a meaningful way with a bit of error and at worst it is of no use. Either way the amount of material and effort required is rather large compared to the usefulness of the results. APFO-4 is also the polymer that, among the investigated, gives the ”nicest” transients. The conclusion from this is that if alternative measurement techniques can be used it is not worthwhile to do TOF.

Not to dismiss the hard work that went into Andersson's experiment, but I would beg to differ with his assessment of the worthiness of the approach. When characterizing a novel material, every measurement adds to the body of knowledge, and as the interpretation of the aggregation of data becomes more cohesive, we end up learning much more of the internal structure. As I have learned, if someone does not understand a phenomena, they tend to dismiss it (myself included).

By their very nature, disordered systems contain a huge state space and we really can't afford to throw out any information.

Which brings up another interesting set of TOF experiments that I dug up. These also deal with organic semiconducting materials -- the polymers with the abbreviations ANTH-OXA6t-OC12 and TPA-Cz3d. The following figures show the TOF results for various applied voltages. I superimposed the entropic dispersion equation form as the red line with the derived mobility in the caption below each figure. The original researcher had applied the Scher&Montroll Continuous Time Random Walk (CTRW) heuristic as indicated by the intersecting sloped lines. The CTRW model clearly fails in this situation as the slopes need quite a bit of creative interpretation. Note that we don't observe the diffusive spike; I integrated the charge from 10% to 100% of the width instead of 0% to 100%.

ANTH-OXA6t-OC12 μ = 0.0025	TPA-Cz3d μ = 0.0013
μ = 0.00155	μ = 0.0004
μ = 0.00125	μ = 0.0005
μ = 0.00085	μ = 0.0006

The number of papers I find, especially when dealing with organic semiconductors, that cannot apply the Scher/Montroll theory indicates that it truly lacks any generality. In other words, it works crappily for describing disorderly crap. I will also say the theory has some very serious flaws, including the claim that an α = 1 defines a non-dispersive material. How could a power-law of -2 be anything but dispersive?

The fact that the entropic dispersion formulation works on any disordered material makes it much more general. Several years ago Scher wrote a popular article for Physics Today extolling the wonders of his theory, and how it seemed to fit a variety of disordered systems. He mentioned how well it fit amorphous silicon based on the number of orders of magnitude that his piece-wise line segments matched. Well, the entropic dispersion does just as well:

And nothing mysterious about that slope of 0.5; that results from the diffusion having a square root dependence with time.

Waste Half-Life

2010-05-21T05:46:00.000-07:00

The big Gulf Spill got me thinking about the half-life of the leaking crude oil and the expanding slick. First of all, the oil will biodegrade over time. We don't have the situation as in CO2 where a sizable fraction will wander around the atmosphere trying to find a suitable location to react and form solutes.

Most of the oil will stay on the surface where it will get plenty of attention from aerobic microoganisms. Some of the oil will sink into the ocean and find anaerobic conditions at the bottom and essentially become inert or wash up on shore as sticky globs. Also the composition of crude oil includes many different hydrocarbons, some of which biodegrade at much slower rates, due to their molecular structure.

So I imagine that we can't calculate the half-life of the spilled oil in terms of a single rate constant, k. This kind of first-order kinetics would likely show an exponential decline, which proceeds pretty quickly once you get past the half-lifetime, 1/k . Instead we will get a mix of various rates, with the fast rates occurring initially and the slower rates picking up the slack.

Radioactive waste-dumps also show a mix of decay constants. Nominally, radioactive material will show a single Poisson emission rate, leading to an exponential decline over time. But when the different radioactive materials get combined, the Geiger counter will pick up this mixture of rates, and the decline will turn from an exponential to a fat tail distribution See the red curve below.

A maximum entropy mix of decay rates (where a high decay rate indicates a potentially more energetic state) will generate the following half-life decline profile:

P(t) = 1/(1+k*t)

where k is the average of the individual rates. This looks exactly the same as the hyperbolic decline of reservoirs in my last post.

As you can see, the combined activity shows a much larger equivalent half-life since the tail has so much meat in it. In the limit of a full dispersion of rate constants, the average half-life will actually slowly diverge as the log of infinity. However, it never reaches this because the slowest decay rate will eventually dominate and that will not diverge.

In any case, this gives a good qualitative description of a random waste dump.

If I make the same MaxEnt assumption for crude oil and assume that the most energetic oil (by the bond strength of the hydrocarbon [1]) will likely prove the most difficult to decompose, then the half-life may also show a similar kind of fat-tail as that of a waste dump. It looks like benzene breaks down much slower than diesel oil for example.

As usual, disordered natural phenomena show many of the same dispersive characteristics, driven largely by maximizing entropy.

Notes:

[1] For the derivation, we assume that we have a mean energy E0 and then a probability density function will show many small energies and progressively fewer high energies.

p(E) = exp(-E/Eo)/E0

but the decomposition rate R depends on E, so that

P(t) = integral of P(t|E)p(E) over all E
P(t|E) = exp(-kE*t)

P(t) = 1/(1+tkEo)

(See this for a more detailed derivation.)

Hyperbolic Decline a Fat-Tail Effect

2010-05-18T06:31:00.000-07:00

If the Gulf Oil spill shows results of a hyperbolic decline, the effects can go on for quite some time.

For a typical reservoir, oil depletion goes through either an exponential decline or a hyperbolic decline. Geologists by and large don't realize this, and definitely don't teach this, but hyperbolic decline constitutes a "fat-tail" effect that results from an aggregation of varying exponential declines summed together. As to the behavior of hyperbolic decline, one notices that the effects tend to drag out for a long time. The fast exponential decline finishes more quickly than the slower exponential components. That's where the fat-tail comes from and why the hyperbolic decline can proceed endlessly are at least as long as the longest exponential portion.

Derivation of hyperbolic decline as a one-liner:

The exponential has a rate of x, and x gets integrated over all possible values of r according to an exponential Maximum Entropy probability density function. You can see the fat-tail in the plot below:

This is just entropy at work because nature tends to want to disperse.

EDIT:

JB asked the question on the slope of the two functions. As plotted, these give the cumulatives. If we want to look at the probability density functions, then yes you will see that the hyperbolic gives a mix of these rates more in line with intuition, with a faster initial slope and then the fatter tail later. See the figure below:

Characterizing mobility in disordered semiconductors

2010-05-09T22:26:00.000-07:00

I always look for analogies between physical systems. This often leads to dead ends but sometimes you uncover some interesting parallels that actually add to the knowledge-base of information and ideas for both systems.

As I worked out the problem of CO2 dispersion in the atmosphere, I went back and revisited the work I did on dispersive transport in amorphous semiconductors. Essentially the same math gets used on both analyses, with the same fundamental goal in mind -- that of trying to characterize the annoyingly sluggish response from an input stimulus.

For the climate case, the poor response comes from CO2 molecules wandering around aimlessly trying to find a good resting place. For the disordered semiconductor, the carrier of electricity (the electron or hole) encounters so many trapping states and scattering centers, that it effectively takes much longer for the charge to cross a region. It does have the advantage of the assist of an electric field, but the low effective transport rate makes an amorphous semiconductor such as hydrogenated amorphous silicon (a-Si:H) marginally useful for any time-sensitive applications -- yet eminently usable as a photo-voltaic.

Still, knowing the physical characteristics helps to understand the nature of the material, and could unlock some secrets beneficial to future applications of material such as polycrystalline or amorphous silicon, or any disordered semiconductor. In the future, we will make mass quantities of this material for the PV industry and we won't have the luxury of single crystal material.

The fact that dispersive transport does have the help of an electric field, makes it amenable to experimentation. By applying various electric fields, one can distinguish between a drift component and a diffusive component (of the photoelectric current, for example). With no electric field, any photo-generated carriers will wander around until they recombine. This can take relatively long times, especially in comparison to a piece of single crystal silicon. As the electric field increases, the carriers get swept out faster and the diffusion plays less of a role.

The fact that the atmosphere has no drift role apart from turbulent diffusion, means that CO2 plays the analogous part of a electronic device with generated carriers, but nowhere to remove them (alas, we have no electrodes attached to the atmosphere). So, I wanted to get a bit of insight by looking at the carrier transport problem, and as a goal, perhaps find a way to increase the removal of CO2 by something equivalent to an electric field, and particularly to ask if this could reduce the CO2 mean residence time.

I noticed one detail that I left hanging on the dispersive carrier transport problem. This had to do with the initial diffusion transient often observed. See the figure to the right (from here). You can see the transient near the start time as a quickly declining response from the initial impulse. The particular trace in the tiny inset came from a non-disordered device (perhaps from a commercial-grade photodetector), as the individual regions show sharp delineations. For a disordered material, the regions show more blurring, as shown in the following figure.

Figure 1: Fitting to the dispersive tail from previous posting. Note the missing initial transient in the curve fit in the curves in color.

I did not include the initial spike term in my initial analysis from last year, as I forgot to apply the chain rule to one set of rate equations. I had justified a transport function that had a non-linear component that propagates as the square-root of time, characteristic of diffusion. Yet to generate a current from this, one needs to differentiate this as a simple chain rule. Not too surprisingly, but perhaps non-intuitively, the derivative of a square root generates the reciprocal of the square root, which of course will spike to infinity at times close to zero. However, the accumulated amount of current generated by this spike nowhere approaches infinity, as the transient has very little width to it. Looking at it on a log-log plot, the width appears long but that occurs simply as an optical illusion.

For a pulsed light source, the entire impulse response equation boils down to a simple charge conservation problem. We know that charge builds up as the photons excite the carriers, but we only know the mean rate and we let the Maximum Entropy Principle figure out the rest. The concentrations build up as the following form, with g(t) acting as the transport growth term across a region w:

C(t) = C0 * g(t)* (1 - exp(-w/g(t))

with

g(t) = sqrt(2*D*t) + u*E*t

where D is the diffusivity, w is the active width, u is the charge mobility, and E is the electric field strength (E = Voltage/w). The total number of excited carriers is C0, and this number provides the maximum amount of current that gets collected. Common to all stochastic probability problems, the conservation of probability becomes a strong constraint.

The current derives as:

I(t) = dC(t)/dt = C0 * dg(t)/dt * (1 - exp(-w/g(t) * (1 + w/g(t))

Note the dg(t)/dt term, which I had neglected to derive completely before, keeping only the drift term (note the 1/sqrt(t) term below).

dg(t)/dt = 0.5*sqrt(D/t) + u*E

Re-plotting the original fitted curve trace with the extra chain-rule term, we can actually see the initial transient. Take a close look at the figure below, and observe how well the curve matches all the inflection points, and works over several orders of magnitude. Mystery solved.

Figure 2 : Dispersive transport which includes a term to describe the initial transient. Note the agreement of the dispersive transport model at short durations. Upper curve fits a fixed average mobility sample. For the lower curve, the average mobility depends on applied electric field strength.

You can spend all sorts of time trying to fit the curves; the more time you spend, the better an estimate you can make of the average mobility, u, and diffusivity, D. Suffice to say, no fudge factors play into the equations. If this isn't a textbook ready formula, I don't know what is.

As I said before, no one in the semiconductor industry seems to use this simple dispersive formulation, preferring to hand-wave and heuristically account for the fat-tails of the transient. Importantly, this particular impulse response function both explains the behavior seen, and derives from the most simple particle counting statistics (i.e maximum entropy randomness), so it likely serves as the most canonical model for dispersive transport in disordered materials.

Linking back to CO2

Now a curious fact presents itself. Not many people in science and engineering seem to understand disorder. If they did, somebody would have discovered this dispersion formulation. Yet they haven't (AFAIK). Billions of dollars goes into semiconductor research and I can only find several purely academic papers on anomolous diffusion and Levy flights and fractional random walks. It really is not that complicated to derive the physics behavior, if you simply assume entropic disorder.

So as it turns out, the dispersion math essentially matches that of what happens to CO2 as it enters the atmosphere. The peculiar piece in the transport that provides that initial photo-current spike acts identically to the fast rate of CO2. In other words, a fraction of charged carriers that can diffuse quickly to a recombination site (i.e. an electrode) act precisely the same as CO2 that reacts quickly and removes itself from the atmosphere. Yet the long tails in the dispersion remain, both in the disordered semiconductor, and in the disordered atmosphere. The fat-tails will kill us in atmospheric CO2 build-up, just like the fat-tails in amorphous semiconductors make it useless to use in a fast microprocessor or in a cell-phone receiver.

Now put 2 and 2 together. No wonder no one knows how to simply describe the CO2 buildup problem! Like the scientists and engineers who experiment with dispersive transport can't see the forest for the trees and thus can't come up with a simple derivation that a near layman can understand, the climate scientists also completely miss out on the obvious and have never come up with the equivalent "probability as logic" formulation.

ImpulseResponseCO2(t) = 1/(1+sqrt(t/T))

That is all there is to it.

Wind Energy Dispersion Analysis

2010-05-06T18:13:00.000-07:00

subtitle: Wind is entirely predictable in its unpredictability

A few weeks ago I wrote about how to derive wind speed characteristics from a straightforward maximum entropy analysis: Wind Dispersion and the Renewable Hubble Curve. This assumed only a known mean of wind energy levels (measured as power integrated over a fixed time period).

From this simple formulation, one can get a wind speed probability graph. Knowing the probability of wind speed, you can perform all kinds of interesting extrapolations -- for example, how long it would take to accumulate a certain level of energy.

I received a few comments on the post, with one by BDog pointing out how the wind flow affects the rate of energy transfer, i.e. the load of kinetic energy enclosed by a volume of air gets pushed along at a rate proportional to its speed. I incorporated that modification in a separate calculation and did indeed notice a dispersive effect on the output. I didn't pick up on this at first so I edited the post with BDog's new correction included.

As a fortunate coincidence, Jerome posted a wind-themed article at TheOilDrum and in the comment section LenGould volunteered a spreadsheet of Ontario wind speed data (thanks Len).

In the past 12 months, the max output was 1017 MW, so there's at least that much online, quite widely distributed accross the 500 mile width of the southern part of the province near the great lakes (purportedly excellent wind resource territory).

On April 20th from 8:00 to 10:00 AM, the output averaged 3.5 MW. (0.34%)
On Mar 16th from 11:00AM to 1:00 PM, the output averaged 4.0 MW. (0.39%)
On Mar 9th from 10:00AM to 6:00 PM, the output averaged 6.7 MW. (0.66%)

That's just a few random picks I made in peak demand hours. I've done thorough analysis of this before and found the data to completely contradict your statement. These wind generators aren't anywhere NEAR to baselaod, and look like they never will be, since winds from here to North Dakota all travel in the same weather patterns.

I used LenGould's data set to try to verify the entropic dispersion model.

The data file consisted of about 36,000 sequential hourly measurements in terms of energy (kilowatt-hours). The following chart shows the cumulative probability distribution function of the energy values. This shows the classic damped exponential function, which derives from either the Maximum Entropy Principle (probability) or the Gibbs-Boltzmann distribution (statistics). It also shows a knee in the curve at about 750 KWh, which I assume comes from a regulating governor of some sort designed to prevent the wind turbine from damaging itself at high winds.

I also charted the region around zero energy to see any effect in the air flow transfer regime (which should be strong near zero). In this regime the probability would go as sqrt(E)*exp(-E/E0) instead of exp(-E/E0). Yet only a linearized trend appears courtesy of the Taylor's series expansion of the exponential around E=0.

Remember that this data consists of a large set of independent turbines. You might think that because of the law of large numbers that the distribution might narrow or show a peak. Instead, the mixture of these turbines over a wide variation in the wind speed provides a sufficiently disordered path so that we can apply the maximum entropy principle.

With a gained confidence in the entropic dispersive model, we can test the obvious nagging question behind wind energy -- How long do we have to wait until we get a desired level of energy?

I generated a resampled set of the data (only resampled in the sense that I used a wraparound at the 4 year length of the data to create a set free from any boundary effects). The output of the resampling essentially generated a histogram of years it would take to reach a given energy level. I chose two levels, E(T)=1000 MW-hrs and E(T)=200 MW-hrs. I plotted the results below along with the predetermined model fit next to the data.

The model from the previous post predicts the behavior used in the two fits:

p(t | E>E(T)) = T * exp(-T/ t ) / t²

where T is the average time it will take to reach E(T). From the exponential fit in the first figure, this gives T= 200/178 and T=1000/178, reespectively, for the two charts. As expected we get the fat-tails that fall off as 1/t^2 (not 1/t^1.5 as the velocity flow argument would support).

The models do not work real effectively at the boundary conditions, simply because the wind turbine limiting governors prevent the accumulation of any energy levels above 1000 MWh level; this occurs either in a short amount of time or at long times as a Poisson process of multiple gusts of lower energy. That said, any real deviations likely arise from short-duration correlations between wind energy measurements spaced close together. We do see this as the lower limit of E(200) shows more correlation curvature than E(1000) does. Wind speeds do change gradually so these correlations will occur; yet these seem minor perturbations on the fundamental entropic dispersion model, which seems to work quite well under these conditions.

As a bottom-line, this analysis tells us what we already intuited. Because of intermittency in wind speed, it often takes a long time to accumulate a specific level of energy. Everyone knows this from their day-to-day experience dealing with the elements. However, the principle of maximum entropy allows us to draw on some rather simple probability formula so that we can make some excellent estimates for long-term use.

The derivation essentially becomes the equivalent of a permanent weather forecast. Weathermen perform a useless function in this regard. Only something on the scale of massive global warming will likely effect the stationary results.

How Shock Model Analysis relates to CO2 Rise

2010-05-03T22:21:00.000-07:00

I would rate the graph below as one of the most famous charts in the annals of science, rivalled only by its close kin, the "hockey stick" graph ( the sketch of Hubbert's Peak is an also-ran in this contest):

Figure 0: The classic and frightening atmospheric CO2 build-up.

From just a technical perspective, it has an interesting composition -- a committed research team that has collected data for some 50 years, measurements showing very little noise, the fascinating periodic cycle due to seasonal variations, and Al Gore to present it.

I don't think many people realize how easy one can derive this curve. You only need a historical record of fossil fuel usage, a few parameters and conversion factors, and the knowledge of how to do a convolution. Since I use convolutions heavily in the Oil Shock model, doing this calculation has become second nature to me.

The way I view it, the excess CO2 production becomes just another stage in the set of shock model convolutions, which model how fossil fuel discoveries transition into reserves and then production. The culminating step in oil usage becomes a transfer function convolution from fuel consumption to a transient or persistent CO2 (depending on what you want to look at). Add in the other hydrocarbon sources of coal and natural gas and you have a starting point for generating the Mauna Loa curve.

The Recipe

First of all, we can roughly anticipate what the actual CO2 curve will look like, as it will lie somewhere between the two limits of immediate recapture of CO2 (the fast transient regime hovering just above the baseline) and no recapture (the persistent integrated regime which keeps accumulating). See Figure 1.

Figure 1: The actual CO2 levels fall between the constraints of immediate uptake (red curve) and persistent inertness (orange curve). The latter results from an accumulation or integration of carbon emissions.

Although this transient will show very long persistence and a very fat tail as I described here, we only need an average rate to generate the initial rise curve. (The oscillating part decomposes trivially, and we can safely add that in later)

So the ingredients:

Conversion factor between tons of carbon generated and an equivalent parts-per-million volume of CO2. This is generally accepted as 2.12 Gigatons carbon to 1 ppmv of CO2. Or ~7.8 Gt CO2 to 1 via purely molecular weight considerations.
A baseline estimate of the equilibrium CO2, also known as the pre-industrial level. This ranges anywhere from 270 ppm to 300 ppm, with 280 ppm the most popular (although not necessarily definitive).
A source of historical fossil fuel usage. The further back this goes in time the better. I have two locations: one from the Wikipedia site on atmospheric CO2 (Image) or one from the NOAA site.
A probability density function (PDF) for the CO2 impulse response (see the previous post). If you don't have this PDF, use the first-order reaction rate exponential function, R(t)=exp(-kt).
A convolution function, which you can do on a spreadsheet with the right macro [1].

The convolution of carbon production Pc(t) with the impulse response R(t) generates C(t):

C(t) = k*[Integral of Pc(t-x)*R(x) from x=0 to x=t] + L

Multiplying the result by a conversion factor k; then adding this to the baseline L generates the filtered Mauna Loa curve as a concentration in CO2 parts per million.

I used R(t)=exp(-t/T), where T=42 years and L=1280 ppm baseline for the following curve fit (using data from Figure 3 for Pc(t)) .

Figure 2: Convolution ala the Shock Model of the yearly carbon emission with an impulse response function. An analytical result from a power-law (N=4) carbon emission model is shown as a comparison..

For Figure 2, I also applied a curve fit model of the carbon generated, which followed a Time⁴ acceleration, and which had the same cumulative as of the year 2004 [2]. You can see subtle differences between the two which indicates that the rate function does not completely smooth out all the yearly variations in carbon emission (see Figure 3). So the two convolution approaches show some consistency with each other, but the fit to the Mauna Loa data appears to have a significant level shift. I will address this in a moment.

Figure 3: Carbon emission data used for Figure 2. A power-law starting in the year 1800 generates a smoothed idealized version of the curve useful for generating a closed-form expression.

The precise form of the impulse response function, other than the average rate selected, does not change the result too much. I can make sense out of this since the strongly increasing carbon production wipes out the fat-tails tails of slower order reaction kinetics (see Figure 4). In terms of the math, a Time⁴ power effectively overshadows a weak 1/sqrt(Time) or 1/Time response function. However, you will see start to see this tail if and when we start slowing down the carbon production. This will give a persistence in CO2 above the baseline for centuries.

Figure 4: Widening the impulse response function by dispersing the rates to the maximum entropy amount, does not significantly change the curvature of the CO2 concentration. Dispersion will cause the curve to eventually diverge and more closely follow the integrated carbon curve but we do not see this yet on our time scale.

Once we feel comfortable doing the convolution, we can add in a piecewise extrapolated production curve and we can anticipate future CO2 levels. We need a fat-tail impulse response function to see the long CO2 persistence in this case (unless 42 years is long enough for your tastes).

The Loose End

If you look at Figure 1, you can obviously see an offset of the convolution result from the actual data. This may seem a little puzzling until you realize that the background (pre-industrial) level of CO2 can shift the entire curve up or down. I used the background level of 280 ppm purely out of popularity reasons. More people quote this number than any other number. However, we can always evaluate the possibility that a higher baseline value would fit the convolution model more closely. Let's give that a try.

The following figure (adapted from here) shows a different CO2 data set which includes the Mauna Loa data as well as earlier proxy ice core data. Based on the levels of CO2, I surmised that the NOAA scientist that generated this graph subtracted out the 280ppm value and plotted the resultant offset. I replotted the data convolution as the dotted gray line.

Figure 5: The CO2 data replotted with extra proxy ice core data, assuming a 280ppm baseline (pre-industrial) level. The carbon production curve is also plotted. You can clearly see that the convolution of the impulse response results in a curve that has a consistent shift of between 10 and 20 ppm below the actual data.

Note that my curve consistently shows a shift 14ppm below the actual data (note the log-scale). This indicates to me that the actual background CO2 level sits 14ppm above 280ppm or at approximately 294ppm. When I add this 14ppm to the curve and replot, it looks like:

Figure 6: The convolution model replotted from Figure 5 with a baseline of 294ppm CO2 instead of 280. Note the generally better agreement to the subtle changes in slope

Although the data does not go through a wide dynamic range, I see a rather parsimonious agreement with the two parameter convolution fit.

Just like in the oil shock model, the convolution of the stimulus with an impulse response function will tend to dampen and shift the input perturbations. If you look closely at Figure 6, you can see faint reproductions of the varying impulse, only shifted by about 25 years. I contend that this "delayed ghosting" comes about directly as a result of the 42-year time constant I selected for the reaction kinetics rate. This same effect occurs with the well-known shift between the discovery peak and production peak in peak oil modeling. Even though King Hubbert himself pointed out this effect years ago, no one else has explained the fundamental basis behind this effect, other than through the application of the shock model. That climate scientists most assuredly use this approach as well points out a potential unification between climate science and peak oil theory. I know David Rutledge of CalTech has looked at this connection closely, particularly in relation to future coal usage.

Bottom Line

To believe this model, you have to become convinced that 294 ppm is the real background pre-industrial level (not 280), and that 40 years is a pretty good time constant for CO2 decomposition kinetics. Everything else follows from first-order rate laws and the estimated carbon emission data.

Of course, this simple model does not take into possible positive feedback effects, yet it does give one a nice intuitive framework to think about how hydrocarbon production and combustion leads directly to atmospheric CO2 concentration changes and ultimately climate change. Doing this exercise has turned into an eye-opener for me, as it didn't really occur to me how straightforward one can derive the CO2 results. Gore had it absolutely right.

Update: From the feedback from some astute TOD readers, it has become clear that some other forcing inputs could easily make up the 14 ppm offset. Changing agriculture and forestry patterns, and other human modifications of the biota could alter the forcing function during the 200+ year time-span since the start of the industrial revolution. Although recyclable plant life should eventually become carbon neutral, the fat-tail of the CO2 impulse response function means that sudden changes will persist for long periods of time. A slight rise from time periods from before the 1800's coupled with an extra stimulus on the order of 500 million tons of carbon per year (think large-scale clearcutting and tilling from before and after this period) would easily close the 14 ppm CO2 gap and maintain the overall fit of the curve.

However, we would need to apply the fat-tail response function, g/(g+sqrt(t)), to maintain the offset for the entire period.

Another comment by EoS:

I don't think it is useful to think of an average CO2 lifetime. That implies a lumped linear model with only a single reservoir, hence an exponential decay towards equilibrium. In reality there are lots of different CO2 reservoirs with different capacities and time constants. So any lumped model had better use several reservoirs with widely varying time constants at a minimum, or else it will get the time behavoir seriously wrong.

It turns out that the variation or dispersion in reaction rates makes very little difference in the slope on the climb up. That is fundamental and I addressed that in Figure 4. The reason for this is very simple mathematics -- the climb up in CO2 is generated by power laws on the order of N>3 or by exponential increases. That is the nature of accelerating fossil fuel usage. In contrast the reaction rates of CO2 have exponents that are negative or have inverse power laws of very low order, the so-called fat-tail distributions. When you put these together, the power law increase essentially crushes the long-tails and all you see are the average value of the faster kinetics. I put in the analytical solution so you can see this directly in the convolution results.

Alternately, apply a simple convolution of accelerating growth [exp(at)] with a first-order reaction decline [exp(-kt)] and you will see what I mean. You get this:

C(t) = (exp(at)-exp(-kt)/(a+k)

The accelerating rate a will quickly overtake the decline term k. If you put in a spread in k values as a distributed model, the same result will occur. That essentially demonstrates Figure 4. Climate scientists should realize this as well since they have known about the uses of convolution in the carbon cycle for years (see chapter 16 in "The carbon cycle" by T. M. L. Wigley and David Steven Schimel).

Yet, if we were to stop burning hydrocarbons today, then we would see the results of the fat-tail decline. Again, I think the climate scientists understand this fact as well but that idea gets obscured by layers of computer simulations and the salient point or insight doesn't get through to the layman. This is understandable because these are not necessarily intuitive concepts.

This following figure models CO2 uptake if we suddenly stop growing fossil fuel use after the year 2007. We don't simple stop using oil and coal, we simply keep our usage constant.

Figure 7: Extrapolation of slow kinetics vs fat-tail kinetics

Up to that point in time a dispersive (i.e. variable) set of rate kinetics will be virtually indistinguishable from a single rate (see Figure 4). And you can see that behavior as the curves match for the same average rate. But once the growth increase is cut off, the dispersive/diffusive kinetics takes over and the rise continues. With the first-order kinetics the growth continues but it becomes self-limiting as it reaches an equilibrium. (see http://mobjectivist.blogspot.com/2010/04/fat-tail-in-co2-persistence.html). This works as a plain vanilla rate theory with nothing by the way of feedbacks in the loop. When we include a real positive feedback, that curve can even increase more rapidly.

Recall that this analysis carries over from studying dispersion in oil discovery and depletion. The rates in oil depletion disperse all over the map, yet the strong push of technology acceleration essentially narrows the dispersed elements so that we can get a strong oil production peak or a plateau with a strong decline. In other words, if we did not have the accelerating components, we would have had a long drawn out usage of oil that would reflect the dispersion. That explains why I absolutely hate the classical derivation of the Hubbert Logistics curve, as it reinforces the opinion of peak oil as some "single-rate" model. In fact just like climate science, everything gets dispersed and follows multiple pathways, and we need to use the appropriate math to analyze that kind of situation.

Climate scientists understand convolution, but peak oil people don't, except when you apply the shock model.

That basically outlines why I want to share these ideas with climate scientists and unify the concepts. It will help both camps, simply by dissemination of fresh ideas and unification of the strong ones.

Notes:
[1] Excel VB convolution script

http://www.microsoft.com/communities/newsgroups/list/en-us/default.aspx?dg=microsoft.public.excel.worksheet.functions&tid=933752da-6f86-4af8-9dba-b9edf57f77d9&cat=en_us_b5bae73e-d79d-4720-8866-0da784ce979c&lang=en&cr=us&sloc=&p=1
Copy the function below into a regular codemodule, then use it like

=SumRevProduct(A2:E2,A3:E3)

It will work with columns as well as rows.

HTH,
Bernie
MS Excel MVP

Function SumRevProduct(R1 As Range, R2 As Range) As Variant
Dim i As Integer
If R1.Cells.Count <> R2.Cells.Count Then GoTo ErrHandler
If R1.Rows.Count > 1 And R1.Columns.Count > 1 Then GoTo ErrHandler
If R2.Rows.Count > 1 And R2.Columns.Count > 1 Then GoTo ErrHandler

For i = 1 To R1.Cells.Count
SumRevProduct = SumRevProduct + _
R1.Cells(IIf(R1.Rows.Count = 1, 1, i), _
IIf(R1.Rows.Count = 1, i, 1)) * _
R2.Cells(IIf(R2.Rows.Count = 1, 1, R2.Cells.Count + 1 - i), _
IIf(R2.Rows.Count = 1, R2.Cells.Count + 1 - i, 1))
Next i
Exit Function
ErrHandler:
SumRevProduct = "Input Error"
End Function

[2] Try this with Wolfram Alpha. It gets finicky sometimes but it does symbolic algebra fairly well.

Rayleigh Fading, Wireless Gadgets, and a Global Context

2010-04-30T19:00:00.000-07:00

The intermittent nature of wind power that I recently posted on has a fundamental explanation based on entropy arguments. It turns out that this same entropy-based approach explains some other related noisy and intermittent phenomena that we deal with all the time. The obvious cases involve the use of mobile wireless gadgets such as WiFi devices, cell phones, and GPS navigation aids in an imperfect (i.e. intermittent) situation. The GPS behavior has the most interesting implications which I will get to in a moment.

First of all, consider that we often use these wireless devices in cluttered environments where the supposedly constant transmitted power results in frustrating fade-outs that we have all learned to live with. An example of Rayleigh fading appears to the right. You can find some signal interference-based explanations for why this happens, originating via the same intentional phase cancellations that occur in noise-cancelling headphones. For the headphones, the electronics flip the phase so all interferences turn destructive, but for wireless devices, the interferences turn random, some positive and some negative, so the result gives the random signal shown.

In the limit of a highly interfering environment the amplitude distribution of the signal shows a Rayleigh distribution, the same observed for wind speed.

p(r) = 2kr exp(-kr²)

Because all we know about our signal is an average power, it occurred to me that one can use Maximum Entropy Principles to estimate the amplitude from the energy stored in the signal, just like one can derive it for wind speed. So, as a starting premise, if we know the average power alone, then we can derive the Rayleigh distribution.

The following figure (taken from here) shows the probability density function of the correlated power measured from a GPS signal. Since power in an electromagnetic signal relates to energy as a flow of constant energy per unit time, then we would expect the energy or power distribution to look like a damped exponential, in line with the maximum entropy interpretation. And sure enough, it does exactly match a damped exponential (note that the Std Dev = Mean, a dead giveaway for an exponential distribution).

p(E) = k*exp(-kE)

Yet since power (E) is proportional to Amplitude squared (r²), we can derive the probability density function by invoking the chain rule.

p(A) = p(E)*dE/dr = exp(-kr²) * d(r²)/dr = 2kr * exp(-kr²)

which precisely matches the Rayleigh distribution, implying that Rayleigh fits the bill as a Maximum Entropy (MaxEnt) distribution. So too does the uniformly random phase in the destructive interference process qualify as a MaxEnt distribution, which will range from 0 to 360 degrees (which gives an alternative derivation of Rayleigh). So all three of these distributions, Exponential, Rayleigh, and Uniform all act together to give a rather parsimonious application of the maximum entropy principle.

The most interesting implication of an entropic signal strength environment relates to how we deal with this power variation in our electronic devices. If you own a GPS, you know this when when trying to acquire a GPS signal from a cold-start. The amount of time it takes to acquire GPS satellites can range from seconds to minutes, and sometimes we don't get a signal at all, especially if we have tree cover with branches swaying in the wind.

Explaining the variable delay in GPS comes out quite cleanly as a fat-tail statistic if you understand how the GPS locks into the set of satellite signals. The solution assumes the entropy variations of the signal strength and integrating this against the search space that the receiver needs to lock-in to the GPS satellites.

Since the search space involves time on one axis and frequency in the other, it takes in the limit ~N² steps to decode a solution that identifies a particular satellite signal sequence for your particular unknown starting position [1]. This gets reduced because of the mean number of steps needed on average in the search space. We can use some dynamic programming matrix methods and parallel processing (perhaps using an FFT) to get this to order N, so the speed-up for a given rate is t². So this will take a stochastic amount of time according to MaxEnt of :

P (t | R) = exp(-c*R*t²)

However because of the Rayleigh fading problem we don't know how long it will take to integrate our signal with regard to the rate R. This rate has a density function proportional to the power level distribution :

p (R) = k * exp(-k*R)

then according to Bayes the conditionals line up to give the probability of acquiring a signal within time t:

P (t) = integral of P(t |R) * p(R) over all R

this leads to the entropic dispersion result of:

P (t < T) = 1/(1+(T/a)²)

where a is an empirically determined number derived from k and c. I wouldn't consider this an extremely fat tail because the acceleration of the search by quadrature tends to mitigate very long times.

I grabbed some data from a GPS project that has a goal to speed up wild-fire response times by cleverly using remote transponders : The FireBug project. They published a good chunk of data for cold-start times as shown in the histogram below. Note that the data shows many times that approach 1000 seconds. The single parameter entropic dispersion fit (a=62 seconds) appears as the blue curve, and it fits the data quite well:

Interesting how we can sharpen the tail in a naturally entropic environment by applying an accelerating technology (also see oil discovery). Put this in the context of a diametrically opposite situation where the diffusion limitations of CO2 slow down the impulse response times in the atmosphere, creating huge fat-tails which will inevitably lead to global warming.

If we can think of some way to accelerate the CO2 removal, we can shorten the response time, just like we can speed up GPS acquisition times or speed up oil extraction. Or should we have just slowed down oil extraction to begin with?

How's that for some globally relevant context?

Notes:

[1] If you already know your position and have that stored in your GPS, the search time shrinks enormously. This is the warm or hot-start mode that is currently used by most manufacturers. The cold-start still happens if you transport a "cold" GPS to a completely different location and have to re-ascquire the position based on unknown starting coordinates.

The Fat-Tail in CO2 Persistence

2010-04-27T23:46:00.000-07:00

I constantly read about different estimates for the "CO2 Half-Life" of the atmosphere. I have heard numbers as short as 6 years and others as long as 100 years or more.

ClimateProgress.org -- Strictly speaking, excess atmospheric CO2 does not have a half-life. The distribution has a very long tail, much longer than a decaying exponential.
As an approximation, use 300-400 years with about 25% ‘forever’.
....
ClimateProgress.org -- David is correct. Half-life is an inappropriate way to measure CO2 in the atmosphere. The IPCC uses the Bern Carbon Cycle Model. See Chapter 10 of the WG I report (Physical Basis) or http://www.climate.unibe.ch/ ~joos/ OUTGOING/ publications/ hooss01cd.pdf

This issue has importance because CO2 latency and the possible slow retention has grave implications for rebounding from a growing man-made contribution of CO2 to the atmosphere. A typical climate sceptic response will make the claim for a short CO2 lifetime :

Endangerment Finding Proposal
Lastly; numerous measurements of atmospheric CO2 resident lifetime, using many different methods, show that the atmospheric CO2 lifetime is near 5-6 years, not 100 year life as stated by Administrator (FN 18, P 18895), which would be required for anthropogenic CO2 to be accumulated in the earth's atmosphere under the IPCC and CCSP models. Hence, the Administrator is scientifically incorrect replying upon IPCC and CCSP -- the measured lifetimes of atmospheric CO2 prove that the rise in atmospheric CO2 cannot be the unambiguous result of human emissions.

Not knowing a lot about the specific chemistry involved but understanding that CO2 reaction kinetics has much to do with the availability of reactants, I can imagine the number might swing all over the map, particular as a function of altitude. CO2 at higher altitudes would have fewer reactants to interact with.

So what happens if we have a dispersed rate for the CO2 reaction?

Say the CO2 mean reaction rate is R=0.1/year (or a 10 year half-life). Since we only know this as a mean, the standard deviation is also 0.1. Placing this in practical mathematical terms, and according to the Maximum Entropy Principle, the probability density function for a dispersed rate r is:

p(r) = (1/R) * exp(-r/R)

One can't really argue about this assumption, as it works as a totally unbiased estimator, given that we only know the global mean reaction rate.

So what does the tail of reaction kinetics look like for this dispersed range of half-lifes?

Assuming the individual half-life kinetics act as exponential declines then the dispersed calculation derives as follows

P(t) = integral of p(r)*exp(-rt) over all r

This expression when integrated gives the following simple expression:

P(t) = 1/(1+Rt)

which definitely gives a fat-tail as the following figure shows (note the scale in 100's of years). I can also invoke a more general argument in terms of a mass-action law and drift of materials; this worked well for oil reservoir sizing. Either way, we get the same characteristic entroplet shape.

Figure 1: Drift (constant rate) Entropic Dispersion

For the plot above, at 500 years, for R=0.1, about 2% of the original CO2 remains. In comparison for a non-dispersed rate, the amount remaining would drop to exp(-50) or ~2*10^-20 % !

Now say that R holds at closer to a dispersed mean of 0.01, or a nominal 100 year half-life. Then, the amount left at 500 years sits at 1/(1+0.01*500) = 1/6 ~ 17%.

In comparison, the exponential would drop to exp(-500/100) = 0.0067 ~ 0.7%

Also, 0.7% of the rates will generate a half-life of 20 years or shorter. These particular rates quoted could conceivably result from those volumes of the atmosphere close to the ocean.

Now it gets interesting ...

Climatologists refer to the impulse response of the atmosphere to a sudden injection of carbon as a key indicator of climate stability. Having this kind of response data allows one to infer the steady state distribution. The IPCC used this information in their 2007 report.

Current Greenhouse Gas Concentrations
The atmospheric lifetime is used to characterize the decay of an instanenous pulse input to the atmosphere, and can be likened to the time it takes that pulse input to decay to 0.368 (l/e) of its original value. The analogy would be strictly correct if every gas decayed according to a simple expotential curve, which is seldom the case.
...
For CO2 the specification of an atmospheric lifetime is complicated by the numerous removal processes involved, which necessitate complex modeling of the decay curve. Because the decay curve depends on the model used and the assumptions incorporated therein, it is difficult to specify an exact atmospheric lifetime for CO2. Accepted values range around 100 years. Amounts of an instantaneous injection of CO2 remaining after 20, 100, and 500 years, used in the calculation of the GWPs in IPCC (2007), may be calculated from the formula given in footnote a on page 213 of that document. The above-described processes are all accounted for in the derivation of the atmospheric lifetimes in the above table, taken from IPCC (2007).

Click on the following captured screenshot for the explanation of the footnote.

The following graph shows impulse responses from several sets of parameters using the referenced Bern IPCC model (found in Parameters for tuning a simple carbon cycle model). What I find bizarre about this result is that it shows an asymptotic trend to a constant baseline, and the model parameters reflect this. For a system at equilibrium, the impulse response decay should go to zero. I believe that it physically does, but that this model completely misses the fact that it eventually should decay completely. In any case, the tail shows huge amount of "fatness", easily stretching beyond 100 years, and something else must explain this fact.

Figure 2: IPCC Model for Impulse Response

If you think of what happens in the atmosphere, the migration of CO2 from low-reaction rate regions to high-reaction rate regions can only occur via the process of diffusion. We can write a simple relationship for Fick's Law diffusion as follows:

dG(t)/dt = D (C(0)-C(x))/G(t)

This states that the growth rate dG(t)/dt remains proportional to the gradient in concentration it faces. As a volume gets swept clean of reactants, G(t) gets larger and it takes progressively longer for the material to "diffuse" to the side where it can react. This basically describes oxide growth as well.

The outcome of Fick's Law generates a growth law that goes as the square root of time - t^1/2. According to the dispersion formulation for cumulative growth, we simply have to replace the previous linear drift growth rate shown in Figure 1 with the diffusion-limited growth rate.

P(t) = 1/(1+R*t^1/2)

or in an alternate form where we replace the probability P(t) with a normalized response function R(t):

R(t) = a/(a+t^1/2)

At small time scales, diffusion can show an infinite growth slope, so using a finite width unit pulse instead of a delta impulse will create a reasonable picture of the dispersion/diffusion dynamics.

Remarkably, this simple model reproduces the IPCC-SAR model almost exactly, with the appropriate choice of a and a unit pulse input of 2 years. The IPCC-TAR fit uses a delta impulse function. The analytically calculated points lie right on top of the lines of Figure 2, which actually makes it hard to see the excellent agreement. The window of low to high reaction rates generates a range of a from 1.75 to 3.4, or approximately a 50% variation about the nominal. I find it very useful that the model essentially boils down to a single parameter of entropic rate origin (while both diffusion and dispersion generates the shape) .

Figure 3: Entropic Dispersion with diffusional growth kinetics describes the CO2 impulse response function with a single parameter a. The square of this number describes a characteristic time for the CO2 concentration lifetime.

You don't see it on this scale, but the tail will eventually reach zero, but at a rate asympotically proportional to the square root of time. In 10,000 years, it will reach approximately the 2% level (i.e. 2/sqrt(10000)).

Two other interesting observations grow out of this most parsimonious agreement.

First of all, why did the original IPCC modelers from Bern not use an expression as simple as the entropic dispersion formulation? Instead of using a three-line derivation with a resultant single parameter to model with, they chose an empirical set of 5 exponential functions with a total of 10 parameters and then a baseline offset. That makes no sense unless their model essentially grows out of some heuristic fit to measurements from a real-life carbon impulse (perhaps data from paleoclimatology investigation of an ancient volcanic eruption; I haven't tracked this down yet). I can only infer that they never made the connection to the real statistical physics.

Secondly, the simple model really helps explain the huge discrepancy between the quoted short lifetimes by climate sceptics and the long lifetimes stated by the climate scientists. These differ by more than a magnitude. Yet, just by looking at the impulse response in Figure 3, you can see the fast decline that takes place in less than a decade and distinguish this from the longer decline that occurs over the course of a century. This results as a consequence of the entropy within the atmosphere, leading to a large dispersion in reaction rates, and the rates limited by diffusion kinetics as the CO2 migrates to conducive volumes. The fast slope evolving gradually into a slow slope has all the characteristics of the "law of diminishing returns" characteristic of diffusion, with the precise fit occurring because I included dispersion correctly and according to maximum entropy principles. (Note that I just finished a post on cloud ice crystal formation kinetics which show this same parsimonious agreement).

Think of it this way: if this simple model didn't work, one would have to reason why it failed. I contend that entropy and disorder in physically processes plays such a large role that it ends up controlling a host of observations. Unfortunately, most scientists don't think in these terms; they still routinely rely on deterministic arguments alone. Which gets them in the habit of using heuristics instead of the logically appropriate stochastic solution.

Which leads me to realize that the first two observations have the unfortunate effect of complicating the climate change discussion. I don't really know, but might not climate change deniers twist facts that have just a kernel of truth? Yes, "some" of the CO2 concentrations may have a half-life of 10 years, but that misses the point completely that variations can and do occur. I am almost certain that sceptics that hang around at sites like ClimateAudit.org see that initial steep slope on the impulse response and convince themselves that a 10 year half-life must happen, and then decide to use that to challenge climate change science. Heuristics give the skilled debater ammo to argue their point any way they want.

I can imagine that just having the ability to argue in the context of a simple entropic disorder can only help the discussion along, and relying on a few logically sound first-principles models provides great counter-ammo against the sceptics.

One more thing ...

So we see how a huge fat tail can occur in the CO2 impulse response. What kind of implication does this have for the long term?

Disconcerting, and that brings us to the point that the point that climate scientists have made all along. With a fat-tail, one can demonstrate that a CO2 latency fat-tail will cause the responses to forcing functions to continue to get worse over time.

As this paper notes and I have modeled, applying a stimulus generates a non-linear impulse response which will look close to Figure 3. Not surprisingly but still quite disturbing, applying multiple forcing functions as a function of time will not allow the tails to damp out quickly enough, and the tails will gradually accumulate to a larger and larger fraction of the total. Mathematically you can work this out as a convolution and use some neat techniques in terms of Laplace or Fourier transforms to prove this analytically or numerically.

This essentially explains the 25% forever in the ClimateProgress comment. Dispersion of rates essentially prohibit the concentrations to reach a comfortable equilibrium. The man-made forcing functions keep coming and we have no outlet to let it dissipate quickly enough.

I realize that we also need to consider the CO2 saturation level in the atmosphere. We may asymptotically reach this level and therefore stifle the forcing function build-up, but I imagine that no one really knows how this could play out.

As to one remaining question, do we believe that this dispersion actually exists? Applying Bayes Theorem to the uncertainty in the numbers that people have given, I would think it likely. Uncertainty in people's opinions usually results in uncertainty (i.e. dispersion) in reality.

This paper addresses many of the uncertainties underlying climate change: The shape of things to come: why is climate change so predictable?

The framework of feedback analysis is used to explore the controls on the shape of the probability distribution of global mean surface temperature response to climate forcing. It is shown that ocean heat uptake, which delays and damps the temperature rise, can be represented as a transient negative feedback. This transient negative feedback causes the transient climate change to have a narrower probability distribution than that of the equilibrium climate response (the climate sensitivity). In this sense, climate change is much more predictable than climate sensitivity. The width of the distribution grows gradually over time, a consequence of which is that the larger the climate change being contemplated, the greater the uncertainty is about when that change will be realized. Another consequence of this slow growth is that further eff orts to constrain climate sensitivity will be of very limited value for climate projections on societally-relevant time scales. Finally, it is demonstrated that the e ffect on climate predictability of reducing uncertainty in the atmospheric feedbacks is greater than the eff ect of reducing uncertainty in ocean feedbacks by the same proportion. However, at least at the global scale, the total impact of uncertainty in climate feedbacks is dwarfed by the impact of uncertainty in climate forcing, which in turn is contingent on choices made about future anthropogenic emissions.

In some sense, the fat-tails may work to increase our certainty in the eventual effects -- we only have uncertainty in the when it will occur. People always think that fat-tails only expose the rare events. In this case, they can reveal the inevitable.

Added Info:

Segalstad pulled together all the experimentally estimated residence times for CO2 that he could find, and I reproduced them below. By collecting the statistics for the equivalent rates, it turns out that the standard deviation approximately equals the mean (0.17/year) -- this supports the idea that the uncertainty in rates found by measurement matches the uncertainty found in nature, thus giving the entropic fat tail. These still don't appear to consider diffusion, which fattens the tail even more.

Authors [publication year]	Residence time (years)
Based on natural carbon-14
Craig [1957]	7 +/- 3
Revelle & Suess [1957]	7
Arnold & Anderson [1957]including living and dead biosphere	10
(Siegenthaler,1989)	4-9
Craig [1958]	7 +/- 5
Bolin & Eriksson [1959]	5
Broecker [1963], recalc. by Broecker & Peng [1974]	8
Craig [1963]	5-15
Keeling [1973b]	7
Broecker [1974]	9.2
Oeschger et al. [1975]	6-9
Keeling [1979]	7.53
Peng et al. [1979]	7.6 (5.5-9.4)
Siegenthaler et al. [1980]	7.5
Lal & Suess [1983]	3-25
Siegenthaler [1983]	7.9-10.6
Kratz et al. [1983]	6.7
Based on Suess Effect
Ferguson [1958]	2 (1-8)
Bacastow & Keeling [1973]	6.3-7.0
Based on bomb carbon-14
Bien & Suess [1967]	>10
Münnich & Roether [1967]	5.4
Nydal [1968]	5-10
Young & Fairhall [1968]	4-6
Rafter & O'Brian [1970]	12
Machta (1972)	2
Broecker et al. [1980a]	6.2-8.8
Stuiver [1980]	6.8
Quay & Stuiver [1980]	7.5
Delibrias [1980]	6
Druffel & Suess [1983]	12.5
Siegenthaler [1983]	6.99-7.54
Based on radon-222
Broecker & Peng [1974]	8
Peng et al. [1979]	7.8-13.2
Peng et al. [1983]	8.4
Based on solubility data
Murray< (1992)	5.4
Based on carbon-13/carbon-12 mass balance
Segalstad (1992)	5.4

Dispersive and non-dispersive growth in ice crystals

2010-04-26T17:10:00.000-07:00

I have wanted to do a post on the growth of crystals for a while. Although on a totally different size scale (microns) and time scale (hours) than that of oil reservoir growth (million barrels and eons/ages), the essential behavioral notions behind the two cases of growth remain much the same.

I suggest that in the most disordered environments, the role of entropy overrides other factors enough so that some simple dispersion arguments can explain the size distribution completely.

Take as an example the formation of ice crystals in a cirrus cloud. Depending on the surrounding temperature a crystal nucleates on some foreign particle and then starts growing. The atmospheric conditions have enough variety that the growth rate will disperse to the maximum entropy amount given a mean rate value. The end state for volumetric growth will also show the same amount of variation. Put these two factors together and we end up with the same size distribution as we have for oil reservoir sizes.

p(x) = S/(S + x)²

where x is the size variate and S is the mean size. This becomes the entroplet form of a probability density function. Plotting this on a logarithmic size scale, the envelope looks the following. I roughly plotted how the dispersed velocity components play into the aggregation of the full profile (note that this serves as a mirror image of the other way of thinking about growth, that of varying end-states).

In some sense, the entroplet aggregates the non-dispersed rate functions which show individually much steeper exponential declines. Graphically, these individual curves do not stand-out on their own since the entropic disorder smooths and disperses the curves efficiently. Keep this in mind for a moment.

The following particle size distribution (PSD) graph shows measurements taken from high altitude cloud experiments. The size gets measured along a single length dimension and the density of the particles takes the place of a probability. I assume that they have plotted the density as a mass function along the x-axis. I convert this to a volumetric size growth problem, integrate the growth for a time duration corresponding to the cloud formation period, and plot the entropic dispersion model result below. The value for S is 2 and the integration scale is 16 (see Note [1]) .

Particle Size Distributions
Cloud ice particle number density n(D) vs. the long dimension of particles as observed at the temperature range of -25° C and -30° C [from Platt, 1997]. It shows a bimodal structure in the ice crystal distribution with the second peak at ~500 mm (edit: microns).

The data fits the entropic dispersion model nicely (green line), but notice at low density that an extra mode shows up as the blue line. This clearly has a sharp exponential drop so likely has a non-dispersive origin. In terms of the higher density entropic model, this stands out as an ordered nucleation regime in the midst of a sea of disordered ice crystal growth modes. Compare to the first figure again and one can see how a distinct peak could occur.

Why or how physically can this bump occur? One can understand this in the context of a completely unrelated analysis, that of marathon finishing times. In a race composed of ordinary citizens, mixed in with some elite runners (bribed by prize money), the elites will form a low density bump in a histogram of finishing times.

The incentives and training and good genetics of the elites separate them from the recreational athlete enough so that they generate a statistically measurable deviation from the trend. Anyone who follows sports understands how this can happen.

I assert that an ice crystal growth model could show this same behavior. Some unknown nucleation process has provided an optimal growth environment for these crystals to deviate from the entropic distribution. Hypothesizing, this could take the form of a catalyst or an accommodating growth substrate. With a power tail of -3/2 this might well have a strong diffusive growth component. However, the nuclei occur rarely enough so they do not drown out the much more common random or spontaneously occurring growth centers. It thus shows up as a clear non-dispersive growth mode in a sea of non-uniformity.

On a micro-level, we do have a variety of reproducible structured shapes to bind against. I would predict that nothing like this would ever happen on the scale of oil reservoirs, black swans notwithstanding. On the scale of oil reservoirs, no two substrates will ever have a common origin, so that the entropic trends will dominate.

BTW, I hope this original analysis may prove of some help to those looking at cloud-based climate change forcing functions, or particle size distributions of volcanic ash (see right). It opens up some possibilities to thinking in a different way. From the plot on the Wikipedia page, they use a log-normal fit to the data yet one that uses an entropic dispersion formulation with the appropriate volume/diameter exponent often can work just as well. Below, I use a root 1/2 on volume on dispersive growth which may indicate a diffusion-controlled rate.

Notes

[1] The approach described here, which essentially broadens the peak but retains the power law tail. Convert the linear growth to a time varying parameter (x=kt) and then treat the growth as an evolving pattern which starts from multiple points in time as the cloud develops.

P(x) = Average 1/(1+S/(k*t)) from t=T to t=T+x/k

This is a cumulative so take the derivative with respect to x to get the PDF.

p(x) ~ 1/(S+x) - 1/(S+kT+x)

As described in the reference, this broadening also applies to species diversification and for oil field growth. The only complicating factor in this analysis is that oil reservoir is a volume, yet crystal sizes get reported as a length and we have to convert that to a volume. This means the derivative has to include a chain rule to convert the volume x to a length parameter L, x ~ L³ generates dx/dL ~ L².

Extracting the Learning Curve in Labor Productivity Statistics

2010-04-24T11:37:00.000-07:00

The term ergodicity refers to the uniformity and fairness in occupation of a system of probability-based states. In the way I think about it, an ergodic Tic-Tac-Toe distribution would have gained enough statistics such that the average measured occupancy of squares of an averaged game's end would equal the probability of some theoretically predicted occupancy.

One can get to this state either by capturing statistics over a long period of time or having some process that doesn't have statistical runs that break the stationarity principle. I sometimes confuse the term stationary with ergodic. The first has more to do with suggesting that a particular snapshot in time does not differ from another snapshot at some later time (i.e. independent events). The latter makes certain that the process has enough variety in its trajectories to visit all the possible states.

One of the first mathematical experiments I remember working in school had to do with random number generation. At the time, hand-held calculators still held some excitement and to take advantage of this, the instructor gave the students an assignment to come up with their own uniform random number generator by creating some complicated algorithm on the calculator. The student could then easily test his results.

As I recall, I felt a certain amount of pride in the clever way that I could get results that appeared random. I don't think I understood what pseudo-random meant at the time. In retrospect, I probably had enough knowledge to realize the power of the inverse trig functions and how they could generate numbers between zero and one, and of the idea of truncating to the decimal part of the number (to get a value between 0 and 1). I am sure that my classroom random number generator would not pass any quality tests, but it worked well enough for learning purposes.

I mention this because I also think it first gave me an idea of how easy one can enter into a very disordered state, and one of almost "uniform randomness". It really does not take much in terms of a combination of linear or non-linear calculations before one can obtain a well-travelled set of trajectories and thus an ergodic distribution that also seems to meet the maximum entropy principle.

The post on Japanese labor productivity statistics rests on this same assumption. Enough variety exists in the set of Japanese labor pool for the statistics to reach an ergodic level. And enough complexity exists in the ways that the laborers operate that they will visit all the labor productivity states possible considering the constraints. That partly explains why I can get confidence that some deep fundamental stochastic behavior explains the results.

The labor productivity model that I used to generate the distribution involved the application of a non-linear Lanchester model for competition. The non-linear nature of this model at least contributes to the potential for the states to reach an ergodic limit. Much like a complex expression invoked on a calculator allows one to generate quite easily a pseudo-random distribution.

In the parlance of attractors, these equations naturally translate into a scenario where you can imagine all sorts of possible trajectories to occur:
Unfortunately, this does not help too much with intuition. Yes one can execute a bunch of random calculations and note how easily various states accumulate, but the non-linear math still gets in the way.

For that reason, I will create an alternate model of labor productivity that invokes similar math but relies on the more intuitive concept of a "learning curve". In practice, a learning curve can exist where a worker can pick up much of the rudimentary skills very quickly, yet to get that last level of productivity will often take much more time [1].

Call the labor productivity C(t) and note that it has some minimum (Max) and maximum level (Min) based on the basic minimum requirements of the company and on the maximum technically achievable productivity.

Then over the course of time, a new hire may see productivity rise according to this simple relationship:

dC(t) = k/(C(t) + v) * dt

This mathematical relationship says that the increment of productivity (dC) per unit of time (dt) is inversely proportional to the productivity accumulated so far. In this case productivity equates to skill level. Thus, the worker can learn the early skills very easily, where the accumulated skills have not risen to a high level. (The parameter v prevents the learning curve rate from going to infinity initially). However, as time advances and these skills accumulate, the growth of new skills starts to slow and it reaches something of an asymptotic limit. We can consider this a law of diminishing returns as the weight of the accumulated skills starts to weigh down the worker.

Rearranging the equation, we get:

(C + v) dC = k dt

We can integrate this to get the result:

½C² + vC= t + constant

This results in the same constraint relationship that I used in the previous post in determining P(C) according to the dispersion formulation for C. The equation if given constrained limits (both Max and Min values) has a solution according to the basic quadratic formula which I used for the Monte Carlo simulation. Otherwise we simply use the ergodic view that all the various values of t will get visited over time, and all the possible constant values will show up according to maximum entropy.

The following figure shows a single instance of a labor productivity learning curve. This has a minimum level and a maximum level at which productivity clamps to. Imagine a set of many of these curves, all with different quadratic slope and maximum constraint and that turns into the statistical mechanics of the labor productivity distribution function for P(C) that we see aggregated over all possible learning times.

which leads to this density function and the excellent fit to the data when entropic dispersion gets applied.

This essentially gives an alternate explanation as a learning curve problem for the excellent fit we get for labor productivity. All workers go through a learning curve that shows a minimum proficiency and a maximum productivity that clamps the level, in between we see the quadratic solution growth which shows up as the inverse power law of 2 in the labor productivity distribution function.

To have it make sense with the Lanchester model labor model, the warfare between firms competing for the same resources provide further elements of disorder. Workers switching firms can cause labor productivity to clamp as progress stalls. By the same token, the cross-pollination of worker skills between firms and of compounding growth adds to the dispersion in the growth rates.

[1] Also see Fick's diffusion equation or a random walk, which shows the same quasi-asympotic properties.

Telegraphing Monkeys, Entropy, and 1/f Noise

2010-04-23T22:42:00.000-07:00

Onward with my entropy-based theory of everything.

The observed behavior known as 1/f noise seems to show up everywhere. They call it 1/f noise (also known as flicker or pink noise) because it follows an inverse power law in its frequency spectrum. It shows up both in microelectronic devices as well as emanating from deep space. Its ubiquity gives it an air of mystery and the physicist Bak tried to explain it in terms of self-organized critical phenomena. No need for that level of contrivance, as ordinary entropic disorder will work just as well.

That as an introduction, I have a pretty simple explanation for the frequency spectrum based on a couple of maximum entropy ideas.

The first relates to the origin of the power law in the frequency spectrum.

Something called random telegraph noise (RTS) (or burst or popcorn noise) can occur for a memory-less process. One can describe RTS by simply invoking a square-wave that has a probability of B to switch states at any dt time interval. This turns into a temporal Markov Chain kind of behavior and the typical noise measurement looks like the following monkeys-typing-at-a-telegraph trace and sounds like popcorn popping in its randomness. It pretty much describes an ordinary Poisson process.

The Markov Chain pulse train as described as above has an autocorrelation function that looks like a two-sided damped exponential. The correlation time equals 1/B.
The autocorrelation (or self convolution) of a stochastic statistical has some interesting properties. In this case, it does have maximum entropy content for the two-sided mean of 1/B -- MaxEnt for the positive axis and MaxEnt for the negative axis.

In addition, the Fourier Transform of the autocorrelation gives precisely the frequency power spectrum. This comes out proportionately to:

S(w) = sqrt(2/π) / (B²+w²)

where w is the angular frequency. The figure below shows the B=1 normalized Mathematica Alpha result.

That result only gives one spectrum of the many Markov switching rates that may exist in nature. If we propose that B itself can vary widely, we can solve for the superstatistical RTS spectrum.

Suppose that B ranged from close to zero to some large value R. We don't have a mean but we have these two limits as constraints. Therefore we let maximum entropy generate a uniform distribution for B.

To get the final spectrum, we essentially average the RTS spectrums over all possible intrinsic rates:

Integrate S(w|B) with respect to B from B=0 to B=R.

This generates the following result

S ' (w) = arctan(R/w)/w

If R becomes large enough then the arctan converges to a constant π/2 and reduces to the 1/f spectrum if we convert w=2π*f.

S ' (f) ~ 1/f

If we reduce R in the limit, then we get a regime that has a 1/f component and a 1/f² above the R transition point, where it reverts back to a telegraph noise power-law.

An excellent paper by Edoardo Milotti, titled 1/f noise: a pedagogical review does a very good job of eliminating the mystery behind 1/f noise. He takes a slightly different tact but comes up with the same result that I have above.

In this review we have studied several mechanisms that produce fluctuations with a 1/f spectral density: do we have by now an "explanation" of the apparent universality of flicker noises? Do we understand 1/f noise? My impression is that there is no real mistery behind 1/f noise, that there is no real universality and that in most cases the observed 1/f noises have been explained by beautiful and mostly ad hoc models.

Milotti essentially disproved Bak's theory and said that no universality stands behind the power-law, just some common sense.

From Wikipedia

There are no simple mathematical models to create pink noise. It is usually generated by filtering white noise.

There are many theories of the origin of 1/ƒ noise. Some theories attempt to be universal, while others are applicable to only a certain type of material, such as semiconductors. Universal theories of 1/ƒ noise are still a matter of current research.
A pioneering researcher in this field was Aldert van der Ziel.

I actually took a class from the professor years ago and handed in this derivation as a class assignment (we had to write on some noise topic). I thought I could get his interest up and perhaps get the paper published, but, alas, he muttered a negative with his Dutch accent.

Years later, I finally get my derivation out on a blog.

As Columbo would say, just one more thing.
If the noise represents electromagnetic radiation, then one can perhaps generate an even simpler derivation. The energy of a photon is E(f)=h*f where h=Plank's constant and f is frequency. According to maximum entropy, if energy radiation remains uniform through the frequency spectrum, then we can only get this result if we apply a 1/f probability density function:

E(f) * p(E(f)) = h*f * (1/f) = constant

Wind Dispersion and the Renewable Hubbert Curve

2010-04-21T17:08:00.000-07:00

Most critics of wind energy points to the unpredictability of sustained wind speeds as a potential liability in widespread use of wind turbines. Everyone can intuitively understand the logic behind this statement as they have personally experienced the variability in day-to-day wind statistics.

However comfortably we coexist with the concept of "windiness", people don't seem to understand the mathematical simplicity behind the wind speed variability. Actually the complexity of the earth's climate and environment contributes to this simplicity (see my TOD posts Information and Crude Complexity and Dispersion, Diversity, and Resilience). I suggest that it also could prove useful to understand the dispersion of airborne particulates, such as what occurred in the aftermath of the Icelandic volcano.

Let me go through the derivation of wind dispersion in a few easy steps. I start with the premise that every location on Earth has a mean or average wind speed. This speed has a prevailing direction but assume that it can blow in any direction.

Next we safely assume that the kinetic energy contained in the aggregate speed goes as the square of its velocity.

E ~ v²

This comes about from the Newtonian kinetic energy law ½mv² and it shows up empirically as the aeronautical drag law (i.e. wind resistance) which also goes as the square of the speed. (Note that we can consider E as an energy or modified slightly as a power, since the energy is sustained over time)

As usual, I apply the Principle of Maximum Entropy to the possible states of energy that exist and come up with a probability density function (PDF) that has no constraints other than a mean value (with negative speeds forbidden).

p(E) = k* exp(-kE)

where k is a constant and 1/k defines the mean energy. This describes a declining probability profile, with low energies much more probable than high energies. To convert to a wind dispersion PDF we substitute velocity for energy and simplify:

p(v)dv = p(E)dE

p(v) = p(E) * dE/dv = 2cv * exp(-cv²)

This gives the empirically observed wind speed distribution, showing a peak away from zero wind speeds and a rapid decline of frequency at higher velocity. I would consider this as another variation of an entropic velocity distribution. Many scientists refer to it as a Rayleigh or Weibull distribution. The excellent agreement in the figure below between the model and empirical data occurs frequently.

Wiki Wind Power
The Weibull model closely mirrors the actual distribution of hourly wind speeds at many locations. The Weibull factor is often close to 2 and therefore a Rayleigh distribution can be used as a less accurate, but simpler model.
Distribution of wind speed (red) and energy (blue) for all of 2002 at the Lee Ranch facility in Colorado. The histogram shows measured data, while the curve is the Rayleigh model distribution for the same average wind speed.

The Wiki explanation pulls its punches by using the heuristic family of Weibull curves. The Rayleigh comes out as the simpler model because it derives from first principles and any deviation from the quadratic exponent seems a bit silly. In the following curve, 1.95 is for all practical purposes the same as 2.

Contrary to other distributions, the wind PDF does not qualify as a fat-tail distribution. This becomes obvious if you consider that the power-law only comes about from the reciprocal measure of time, and since we measure speed directly, we invoke the entropic velocity profile directly as well.

So the interesting measure relates to the indirect way that we perceive the variations in wind speed. Only over the span of time do we detect the unpredictability and disorder in speed -- whether by gustiness or long periods of calm.

We can then pose all sorts of questions based on the entropic wind speed law. For example, how long would we have to wait to generate an accumulated amount of energy?

We can answer this analytically by simply equating the steady-state wind speed to a power and then integrating over all possibilities of the distribution that meet the minimum accumulated energy condition over a period of time.

The naive answer is trivial with the time PDF turning into the following fat-tail power-law:

p(t | E>cv²) = (cv²) * exp(-cv²/ t ) / t²

This equation corresponds to the following graphed probability density function (where we set an arbitrary E of cv² = 25). It basically illustrates the likelihood of accumulating a specific energy goal within a given Time.

Because of the scarcity of high wind speeds and the finite time it takes to accumulate energy, one observes a short ramp-up to the peak. Typical wind-speeds round out the peak and the relatively common slow speeds contribute to the long fat-tail.

But since power has an extra velocity factor (Power = Drag*velocity, see Betz' law), it takes longer to integrate the low power values and the exponent changes from a quadratic value of 2 to a value of 5/3, via a chain rule.

p(t ) = (2*d/3) * exp(-d / t^2/3 ) / t^5/3

Thanks to E.Swanson at TOD for pointing this out. As you can see in the graph below, the fat-tail becomes fatter and the ramp-up a bit sooner for roughly the same peak value.
That gives us the power-law and a shape that looks surprisingly close to the time depletion curve for an oil reservoir. In fact, since probabilities have such universal properties, the curvature of this profile has the same fundamental basis as the Hubbert oil depletion profile. The huge distinction lies in the fact that wind energy provides a renewable source of energy, whereas oil depletion results in a dead-end.

So the hopelessness of the Hubbert curve when applied to Peak Oil turns to a sense of optimism when you realize that wind power generates a case of a Renewable Hubbert Curve. In other words, anytime you spin-up the wind-turbine you can always obtain a mini Hubbert cycle, if you have patience, just like you need patience with a train schedule.

Amazing the power of a simple consistent model; both dispersive discovery and the entropic wind dispersion model use the same set of ideas from probability. I find it unfortunate that not enough analysts can see through the complexity and discover the underlying elegance and intuitive power of simple entropy arguments.

Found recently:

http://www.theoildrum.com/node/6386#comment-613465
Once upon a time, long ago and not all that far away, I got an F in math at a famous engineering school because I simply wrote down the answer to a problem on the final instead of going thru all the tedious procedures. I complained to the TA that the answer was obvious and I didn't need the procedure. He said "Ya didn't show that you had learned the procedure, and so the F stays".

I then went down that long soulless hallway to the office of the SUMG (short ugly math genius) who took one glance at the problem and sneered "don't bother me with such trivia, the answer is obviously--". So I went back to the TA with this, and he answered "no procedure, no grade."

So---I look at wind turbines and think "sure this thing works, and it has a energy return of maybe 20 to 40 as is, and there are obvious ways to make it better. That's good enough, so let's go for it".

Right on.

Power Laws and Entrophic Systems

2010-04-20T19:29:00.000-07:00

Every time I look at a distribution of some new set of data I come across, I can't help but notice a power law describing the trend. Also known as fat-tail distributions, NN Taleb popularized these curves with his book The Black Swan. Curiously, both he and Benoit Mandelbrot, his mentor, never spent much time pursuing the rather simple explanation for the classic 1/X cumulative distribution, often known as Zipf's Law or the Mandelbrot/Zipf variant.

That drives me crazy, as I don't consider making the assertion of the fat-tail's origin a huge academic risk. To me, I find that a power law invariably comes about when you invert a measure that follows a Poisson or exponential distribution with respect to its dimension. For example, a velocity distribution that shows a damped exponential profile, when inverted to the time dimension will give a fat-tail power law. I call that entropic velocity dispersion. See the last post on train scheduling delays to appreciate this.

So why doesn't anyone understand this rather simple transform and the general concept of power laws? Let us go through the various disciplines.

The Statistician: Hopes that everything will fall into a normal Gaussian distribution, so they can do their classical frequentist statistical analysis. Many show wariness of Bayesian analysis, of which the Maximum Entropy Principle falls into. This pollutes their pure statistics with the possibility of belief systems or even physics. Likes to think in terms of random walk, which almost always puts them in the time-domain and they can then derive their Poisson or Gaussian profiles with little effort. They need to usually explain fat-tails by Levy walks, which curiously assumes a power law as a first premise. So they never get the obviousness of the natural power law.

The Physicist: Wants to discover new laws of nature and thus directs research to meet this goal. Over the years, physicists have acquired the idea that power-laws somehow connect to critical phenomena, and this has gotten pounded into their skulls. Observing critical phenomena usually holds out hope that some new state of matter or behavior would reveal itself. So associating a power law within a controlled experiment of some phase transition offers up a possible Nobel Prize. They thus become experts at formulating an expression and simplifying it so that a power-law trend shows up in the tail. They want to discover the wondrous when the mundane explanation will just as easily do.

The Mathematician: Doesn't like probability too much because it messes up their idealized universe. If anything, they need problems that can challenge their mind without needing to interface to the real world. Something like Random Matrix Theory, which may explain velocity dispersion but within an N-dimensional state space. Let the Applied Mathematician do the practical stuff.

The Applied Mathematician: Unfortunately, no one ever reads what Applied Mathematicians have to say because they write in journals called Applied Mathematics, where you will find an article on digital image reconstruction calculation next to a calculation of non-linear beam dynamics. How would anyone searching for a particular solution find this stuff in the first place? (outside of a Google search that is)

The Engineer: Treats everything as noise and turbulence which gets in their way. They see 1/f noise and only want to get rid of it. Why explain something if you can't harness its potential? Anyways, if they want to explain something they can just as easily use a heuristic. Its just fuzzy logic at its core. For the computational engineer, anything else that doesn't derive from Fokker-Planck, Poisson's Equation, or Navier-Stokes is suspect. For the ordinary engineer, if it doesn't show up in the Matlab Toolbox it probably doesn't exist.

The Economist and The Quantitative Analyst: Will invoke things like negative probabilities before they use anything remotely practical. Math only serves to promote their economic theory predicting trends or stability, with guidance from the ideas of econometrics. Besides, applying a new approach will force them to sink the costs of years of using normal statistics. And they should know the expense of this because they invented the theory of sunk costs. The economist is willing to allow the econophysicist free reign to explore all the interesting aspects of their field.

I wrote this list at least partly tongue-in-cheek because I get frustrated by the amount of dead-end trails that I see scientists engaged in. They get tantalizingly close with ideas such as superstatistics but don't quite bridge the gap. If I couldn't rationalize why they can't get there, I would end up banging my head on the wall, thinking that I have made some egregious mistake or typo somewhere in my own derivation.

Lastly, many of the researchers that investigate fat-tail distributions resort to invoking ideas like Tsallis Entropy, which kind of looks like entropy, and q-exponentials, which kind of look like exponentials. The ideas of chaos theory, fractals, non-linear dynamics, and complexity also get lumped together with what may amount to plain disorder.

Unfortunately this just gets everyone deeper in the weeds, IMO.

Entrophic

Over at the Oil Drum, somebody accidentally used the non-existent term "entrophic" in a comment. I asked if they meant to say entropic or eutrophic, of which the suffix trophic derives from the Greek for food or feeding.

Yet I think it has a nice ring to it, seeing as it figuratively describes how disorder/ dispersion play in depletion of resources (i.e. eating). So entrophic phenomena describe the smeared out and often fat-tail depletion dynamics that we observe, but that few scientists have tried to explain.

Dispersion and Train Delays

2010-04-19T20:03:00.000-07:00

Since nothing I write about contains any difficult math, I thought to present this post as a word problem.

Say you reside in England. Consider taking a train from point A to B. You have an idea of how long the trip will take based on the current train schedule but have uncertainty on its latency (i.e. possible time delay). The operators of the trains only tell you at best the fraction that arrive within a few minutes of their scheduled time. In other words, you have no idea of the fatness of the tails, and whether you will get delayed by tens of minutes on a seemingly short run.

How would you derive the probability of a specific lateness of time dt based only on the knowledge of the distance X between point A and B, the maximum train speed Vm, and an average train speed (same as distance X divided by the average trip duration t )?

I won't solve this problem in anyway remotely that I would consider a classical approach. Instead I will make an assumption based on the principle of maximum entropy (which I tend to use for everything I run across these days).

Let us see how close we can come to the empirical distribution of train delay times observed based on assuming very limited infomation.

First of all, I have rather limited experiences travelling by train in England. I do know that the speed of a train can vary quite a bit as I have experienced the crawl from Heathrow to Victoria Station. I also know that the train has some maximum speed that it won't exceed. You realize this when you notice that the train rarely arrives early. So the average train speed and maximum train speed provides a pair of constraints that we can use for estimating a train velocity probability density function (PDF).

I assert via maximum entropy (MaxEnt) arguments that the train velocity (or speed for purists) PDF likely looks like this:
In accordance with the constraints, MaxEnt predicts an exponential profile up to the maximum value, which ends at 1.44 miles/minute in the histogram. I would describe it as a dispersive velocity profile following a reverse damped exponential (the damping takes place for smaller velocities, see this post for a similar behavior describing TCP statistics):

p(V) = 1/dv * exp ((V-Vm)/dv)

where V ranges from 0 to Vm. The factor dv relates to the average speed by

dv = Vm - X/t

The smaller that dv becomes, the closer that the average speed approaches the maximum speed. Since we don't know anything more about the distributions of speeds, MaxEnt suggests that this exponential "approximation" will work just fine.

The only remaining step we need to do involves some probability theory to relate this velocity distribution to a time delay distribution. Fortunately, once you get the hang of it, one can just about do this in your sleep.

As the easiest approach, figure out the cumulative probability of velocities that will reach the destination in time T. This becomes the integral over p(v) from a just-in-time speed X/T to the maximum speed Vm:

Integral of p(v) from v=X/T to Vm

This results in the cumulative probability distribution function (upper-case P):

P(T) = 1 - exp ((X/T-Vm)/dv)

To turn this into a probability density function, we simply need to take the derivative with respect to T.

p(T) = Vm*Tm/(dv*T²) * exp(-dt*Vm/(dv*T))

We can also cast it in terms of the time delay by

T = dt + X/Vm

p(dt) = (X/dv) exp (-dt*Vm/(dv*(dt+X/Vm)) / (dt+X/Vm)²

This might seem like a complicated expression but all of the parameters are well known but one, dv. And even in this case, we can estimate dv from some data. The following plot illustrates how the PDF changes with dv. Note that as dv becomes small, the probability of arriving on time becomes closer to 1 (e.g. the Swiss or German train system)

By definition, we have a fat-tail probability distribution because the density follows off as an inverse power law of exponent 2.

So we need some data to check how good this distribution works. My interest in the topic of train scheduling actually came from a paper by an expert in the field of superstatistics [1]:

"Modeling train delays with q-exponential functions"
Keith Briggs and Christian Beck
Physica A: Statistical Mechanics and its Applications
Volume 378, Issue 2, 15 May 2007, Pages 498-504
arxiv.org

Superstatistics may work as an alternate technique to solve the problem but as used by Briggs and Beck, it requires a couple of arbitrary parameters to fit the distribution to. This actually points to the difference between solving a problem as I am trying to do, versus blindly employing arbitrary statistical functions as Briggs and Beck attempt.

In any case, the authors did the hard work of collating statistics of over 2 million train departures covering 23 train stations in the British Rail Network.

Most of the statistics look generally similar but I take the route from Reading to London Paddington station as a trial run.

Both the superstatistics approach and my entropic dispersion solution fit the data remarkably well. One can see clearly that the profile neither follows Poisson statistics (which would give a straight line) nor does it follow normal Gaussian statistics (which would show an upside-down parabolic shape on a semi-log graph).
The shape of the tail points to the real fat-tail probabilities that occur on British rail lines -- as we can see long delays do occur and the power law likely comes from the principle of entropic dispersion.

So to recap, the values I used stem from real observables:

X = 36 miles (exact)
Vm = 1.44 miles/minute (from the schedule, see below)
dv = 0.26 miles/minute (from the curve fit or alternatively from the average latency)

I stated earlier that most people would approach this problem from the perspective of Poisson statistics using time as the varying parameter. Briggs and Beck do this as well, but they use another layer of probability, called superstatistics [1] to "fatten" the tail. Although I appreciate the idea of superstatistics and have used before without realizing it had a name (see my previous post on hyperbolic discounting of an example of a doubly integrated distribution function ), I believe that entropic dispersion of velocities gives a more parsimonious explanation.

Chalk up another success story for entropic dispersion, especially as it shows consistency with the data on the entropic dispersion for all forms of human travel.

Yet, why does no one else use this approach? Does no one understand how to do word problems any longer? I swear that a smart high school student could derive the solution if given the premise.

[1] "Superstatistics", Beck C.; Cohen E.G.D. Physica A, Volume 322, 1 May 2003, pp. 267-275(9)