Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

GEC Monopore Characteristic Function

This notebook derives the GEC (Giddings–Eyring–Carmichael) monopore characteristic function from first principles using SymPy.

Goal: Show how the CF ϕtS(ω)=exp[nˉ(11iωτˉ1)]\phi_{t_S}(\omega) = \exp\left[\bar{n}\left(\frac{1}{1-i\omega\bar{\tau}} - 1\right)\right] emerges from:

  1. Exponential egress time distribution

  2. Poisson number of adsorptions

  3. Independence assumption rM{τS,i}r_M \perp \{\tau_{S,i}\}

Step 1: CF of Exponential Distribution

For a single egress time τExp(1/τˉ)\tau \sim \text{Exp}(1/\bar{\tau}):

ϕτ(ω)=E[eiωτ]=0eiωt1τˉet/τˉdt=11iωτˉ\phi_\tau(\omega) = \mathbb{E}[e^{i\omega\tau}] = \int_0^\infty e^{i\omega t} \frac{1}{\bar{\tau}} e^{-t/\bar{\tau}} dt = \frac{1}{1 - i\omega\bar{\tau}}
from sympy import symbols, exp, factorial, simplify, I, oo, Sum
import sympy as sp

# Define symbols
w, n, t_bar, n_bar = symbols('w n t_bar n_bar', real=True, positive=True)

Step 2: CF of Sum of n i.i.d. Exponentials

For Sn=i=1nτiS_n = \sum_{i=1}^n \tau_i where τi\tau_i are i.i.d. exponential:

ϕSn(ω)=[ϕτ(ω)]n=(11iωτˉ)n\phi_{S_n}(\omega) = [\phi_\tau(\omega)]^n = \left(\frac{1}{1-i\omega\bar{\tau}}\right)^n

This uses independence: E[eiω(X+Y)]=E[eiωX]E[eiωY]\mathbb{E}[e^{i\omega(X+Y)}] = \mathbb{E}[e^{i\omega X}]\mathbb{E}[e^{i\omega Y}]

# CF of single exponential random variable
phi_tau = 1 / (1 - I*w*t_bar)
print("CF of single egress time:")
phi_tau
CF of single egress time:
Loading...

Step 3: Poisson Distribution

For rMPoisson(nˉ)r_M \sim \text{Poisson}(\bar{n}):

P(rM=n)=enˉnˉnn!P(r_M = n) = \frac{e^{-\bar{n}}\bar{n}^n}{n!}
# CF of sum of n i.i.d. exponentials
phi_sum_n = phi_tau**n
print("CF of sum of n egress times:")
phi_sum_n
CF of sum of n egress times:
Loading...

Step 4: Law of Total Expectation (Requires Independence!)

For random sum tS=i=1rMτit_S = \sum_{i=1}^{r_M} \tau_i where rM{τi}r_M \perp \{\tau_i\}:

ϕtS(ω)=ErM[E[eiωtSrM]]=n=0P(rM=n)[ϕτ(ω)]n\phi_{t_S}(\omega) = \mathbb{E}_{r_M}\left[\mathbb{E}[e^{i\omega t_S} | r_M]\right] = \sum_{n=0}^\infty P(r_M=n) \cdot [\phi_\tau(\omega)]^n

Critical: This step requires independence between rMr_M and {τi}\{\tau_i\}!

# Poisson probability mass function
P_n = exp(-n_bar) * n_bar**n / factorial(n)
print("Poisson PMF P(r_M = n):")
P_n
Poisson PMF P(r_M = n):
Loading...

Step 5: Evaluate the Sum (Poisson Generating Function)

We need to evaluate:

n=0enˉnˉnn!znwhere z=ϕτ(ω)\sum_{n=0}^\infty \frac{e^{-\bar{n}}\bar{n}^n}{n!} \cdot z^n \quad \text{where } z = \phi_\tau(\omega)

This is the probability generating function of Poisson:

E[zrM]=enˉ(z1)\mathbb{E}[z^{r_M}] = e^{\bar{n}(z-1)}
# Compound Poisson CF: sum over all possible n values
# φ_tS(ω) = Σ P(n) * [φ_τ(ω)]^n
summand = P_n * phi_sum_n
print("Summand P(r_M=n) * [φ_τ(ω)]^n:")
summand
Summand P(r_M=n) * [φ_τ(ω)]^n:
Loading...

Step 6: Final Result (Dondi Eq. 43)

Substituting z=11iωτˉz = \frac{1}{1-i\omega\bar{\tau}} into enˉ(z1)e^{\bar{n}(z-1)}:

ϕtS(ω)=exp[nˉ(11iωτˉ1)]\phi_{t_S}(\omega) = \exp\left[\bar{n}\left(\frac{1}{1-i\omega\bar{\tau}} - 1\right)\right]

This can also be written as:

ϕtS(ω)=exp[nˉiωτˉ1iωτˉ]\phi_{t_S}(\omega) = \exp\left[\frac{\bar{n} \cdot i\omega\bar{\tau}}{1-i\omega\bar{\tau}}\right]
# The sum Σ (e^(-n_bar) * n_bar^n / n!) * z^n = exp(n_bar*(z-1))
# We substitute z = phi_tau
z = symbols('z')

# Poisson generating function
pgf = exp(n_bar * (z - 1))
print("Poisson generating function E[z^r_M]:")
display(pgf)

# Substitute z = phi_tau to get GEC CF
gec_cf_derived = pgf.subs(z, phi_tau)
print("\nGEC CF (substitute z = φ_τ(ω)):")
gec_cf_derived
Poisson generating function E[z^r_M]:
Loading...

GEC CF (substitute z = φ_τ(ω)):
Loading...

Verification: Check Against Original Form

Now we verify this matches the form used in later cells.

# Simplify the exponent
exponent = simplify(phi_tau - 1)
print("Exponent φ_τ(ω) - 1:")
display(exponent)

# Full GEC CF
gec_cf_final = simplify(gec_cf_derived)
print("\nFinal GEC monopore CF:")
gec_cf_final
Exponent φ_τ(ω) - 1:
Loading...

Final GEC monopore CF:
Loading...
from sympy import symbols, exp, simplify, I, cancel
w = symbols('w')
n1, t1 = symbols('n1, t1', positive=True, real=True)

# Original form from Dondi Eq. 43
gec_monopore_cf_original = exp(n1*(1/(1 - I*w*t1) - 1))

print("Original form (Dondi Eq. 43):")
display(gec_monopore_cf_original)

# Manually simplify the exponent to show equivalence
# 1/(1 - iωτ) - 1 = [1 - (1 - iωτ)]/(1 - iωτ) = iωτ/(1 - iωτ)
exponent_original = n1*(1/(1 - I*w*t1) - 1)
exponent_simplified = simplify(exponent_original)

print("\nOriginal exponent:")
display(exponent_original)
print("\nSimplified exponent:")
display(exponent_simplified)

# This should be: n1*I*w*t1/(1 - I*w*t1)
# Verify by manual calculation
numerator = 1 - (1 - I*w*t1)
denominator = 1 - I*w*t1
manual_simplified = n1 * numerator / denominator
manual_simplified = simplify(manual_simplified)

print("\nManual simplification [1 - (1 - iωτ)]/(1 - iωτ):")
display(manual_simplified)

# Alternative form
gec_alternative = exp(n1*I*w*t1/(1 - I*w*t1))
print("\nAlternative form exp[n̄·iω·τ̄/(1 - iω·τ̄)]:")
display(gec_alternative)

# Verify they're the same by checking the difference
print("\nDifference between original and alternative:")
display(simplify(gec_monopore_cf_original - gec_alternative))

print("\n✓ Both forms are mathematically identical!")
print("  Form 1: exp[n̄(1/(1-iωτ̄) - 1)]")
print("  Form 2: exp[n̄·iωτ̄/(1-iωτ̄)]")
Original form (Dondi Eq. 43):
Loading...

Original exponent:
Loading...

Simplified exponent:
Loading...

Manual simplification [1 - (1 - iωτ)]/(1 - iωτ):
Loading...

Alternative form exp[n̄·iω·τ̄/(1 - iω·τ̄)]:
Loading...

Difference between original and alternative:
Loading...

✓ Both forms are mathematically identical!
  Form 1: exp[n̄(1/(1-iωτ̄) - 1)]
  Form 2: exp[n̄·iωτ̄/(1-iωτ̄)]

Moments from Characteristic Functions

The characteristic function ϕ(ω)\phi(\omega) encodes all statistical properties of a distribution. We can extract moments using derivatives:

Raw Moments (about the origin)

The kk-th raw moment is μk=E[Xk]\mu'_k = \mathbb{E}[X^k], obtained via:

μk=(i)kdkϕ(ω)dωkω=0\mu'_k = (-i)^k \left.\frac{d^k \phi(\omega)}{d\omega^k}\right|_{\omega=0}
  • k=1k=1: Mean (first raw moment)

  • k=2k=2: Second raw moment E[X2]\mathbb{E}[X^2]

Central Moments (about the mean)

The kk-th central moment is μk=E[(Xμ)k]\mu_k = \mathbb{E}[(X-\mu)^k] where μ=E[X]\mu = \mathbb{E}[X]:

  • k=2k=2: Variance σ2\sigma^2

  • k=3k=3: Skewness component (normalized: γ1=μ3/σ3\gamma_1 = \mu_3/\sigma^3)

  • k=4k=4: Kurtosis component (normalized: κ=μ4/σ4\kappa = \mu_4/\sigma^4)

Method: Central moments can be computed from the shifted CF:

ϕ~(ω)=ϕ(ω)eiωμμk=(i)kdkϕ~(ω)dωkω=0\tilde{\phi}(\omega) = \phi(\omega) \cdot e^{-i\omega\mu} \quad \Rightarrow \quad \mu_k = (-i)^k \left.\frac{d^k \tilde{\phi}(\omega)}{d\omega^k}\right|_{\omega=0}

This is implemented below by dividing the CF by eiωμe^{i\omega\mu} before differentiation.

How to use SymPy to symbolically compute moments

We can extract statistical moments from the characteristic function using SymPy as follows.

from sympy import symbols, exp, diff, simplify, I
w = symbols('w')

def raw_moment(cf, k):
    return simplify((-I)**k * diff(cf, w, k).subs(dict(w=0)))

def central_moment(cf, k):
    m = raw_moment(cf, 1)
    return simplify(raw_moment(cf/exp(I*w*m),k))
n1, t1 = symbols('n1, t1')
gec_monopore_cf = exp(n1*(1/(1 - I*w*t1) - 1))
gec_monopore_cf
Loading...
raw_moment(gec_monopore_cf, 1)
Loading...
central_moment(gec_monopore_cf, 2)
Loading...
central_moment(gec_monopore_cf, 3)
Loading...