Understanding Delta Function Approximations: Lorentzian Delta Sequence (Cauchy Kernel)

Delta-Convergent Sequences — Refined Blog with SageMath Symbolics, Physics Insights, and Cleaner Code

The Dirac delta function isn’t a “normal” function — it’s an idealization used to represent a point source. It's infinitely narrow, infinitely tall, and yet integrates to 1. We approximate it using delta-convergent sequences: real functions depending on a parameter that becomes increasingly peaked at zero as the parameter vanishes.
This post explores the three most common delta-approximating sequences using SageMath, including plots, integration checks, and real-world meaning.

Why Study These Approximations?

Delta functions are central in many fields:

• Signal Processing: Ideal impulse, filter response
• Physics: Point charges/masses, Green's functions
• Spectral Theory: Lorentzian profiles in resonance
• Diffusion Models: Gaussians arise from the heat equation
• Numerics: Regularizing singular integrals

Each kernel has a story to tell.

Lorentzian Delta Sequence

Formula \[ f_{\epsilon}(x) = \frac{1}{\pi}. \frac{\epsilon}{x^2+ {\epsilon}^2} \]

• Sharp peak at 𝑥=0
• Always integrates to 1
• Smooth and rational

# Check the Function Definition
var('x epsilon')
f_lorentz(x, epsilon) = (1/pi) * (epsilon / (x^2 + epsilon^2))
f_lorentz(x, epsilon) 

# Check Symbolic Integration
var('xi')
assume(epsilon > 0)  # Ensure SageMath understands epsilon is positive
integral(f_lorentz(xi, epsilon), xi, -oo, oo).simplify_full()

#Check Limit at x → 0
limit(f_lorentz(x, epsilon), epsilon=0)

#  Integral Test (Distributional Behavior)
var('a b')
assume(a < 0, b > 0)  # Ensure a < 0 < b to match delta behavior
integral(f_lorentz(xi, epsilon), xi, a, b).simplify_full()

#  Numerical Evaluation
# To see how the integral converges to 1 as ( \epsilon \to 0 ), run:

import numpy as np
import matplotlib.pyplot as plt
import sage.all as sage

def lorentzian_integral(epsilon, a=-1, b=1):
    return (np.arctan(b/epsilon) - np.arctan(a/epsilon)) / np.pi

# Test for different epsilon values
epsilons = np.logspace(-3, 0, 50)  # Log-spaced values from 0.001 to 1
integral_values = [lorentzian_integral(eps) for eps in epsilons]

# Plotting
plt.figure(figsize=(8, 5))
plt.plot(epsilons, integral_values, marker='o', linestyle='-', color='blue')
plt.axhline(y=1, color='r', linestyle='--', label="Expected Limit (1)")
plt.xscale("log")
plt.xlabel(r"$\epsilon$")
plt.ylabel(r"Integral Value")
plt.title("Numerical Verification: Lorentzian Integral Convergence")
plt.legend()
plt.grid(True)
plt.show()

# Lorentzian delta approximation
var('x epsilon xi')
assume(epsilon > 0)  # Ensure epsilon is positive for proper symbolic handling
f_lorentz(x, epsilon) = (1/pi) * (epsilon / (x^2 + epsilon^2))

# Symbolic normalization check
integral(f_lorentz(xi, epsilon), xi, -oo, oo).simplify_full()

# Limit at x = 0 to verify δ-behavior
limit(f_lorentz(x, epsilon), epsilon=0)

# Plot Lorentzian for different ε
p1 = plot(f_lorentz(x, 0.5), (x, -5, 5), color='red', legend_label='ε = 0.5') + \
     plot(f_lorentz(x, 0.2), (x, -5, 5), color='blue', legend_label='ε = 0.2') + \
     plot(f_lorentz(x, 0.05), (x, -5, 5), color='green', legend_label='ε = 0.05')
p1.show(title='Lorentzian Approximation to δ(x)', ymin=0, ymax=3)

# Numerical Limits Instead of Symbolic
epsilon_vals = [0.1, 0.01, 0.001, 0.0001]
[f_lorentz(0, eps).n() for eps in epsilon_vals]
# Check Symbolic Integration
# Confirm normalization:

var('xi')
assume(epsilon > 0)  # Ensure SageMath understands epsilon is positive
integral(f_lorentz(xi, epsilon), xi, -oo, oo).simplify_full()

# First & Second Derivative Computation
f_lorentz_prime(x, epsilon) = diff(f_lorentz(x, epsilon), x)
f_lorentz_double_prime(x, epsilon) = diff(f_lorentz_prime(x, epsilon), x)

f_lorentz_prime(x, epsilon), f_lorentz_double_prime(x, epsilon)

#Plot First & Second Derivative
p1 = plot(f_lorentz_prime(x, 0.5), (x, -5, 5), color='red', legend_label="ε=0.5") + \
     plot(f_lorentz_prime(x, 0.2), (x, -5, 5), color='blue', legend_label="ε=0.2") + \
     plot(f_lorentz_prime(x, 0.05), (x, -5, 5), color='green', legend_label="ε=0.05")

p1.show(title="First Derivative of Lorentzian Approximation")

p2 = plot(f_lorentz_double_prime(x, 0.5), (x, -5, 5), color='red', legend_label="ε=0.5") + \
     plot(f_lorentz_double_prime(x, 0.2), (x, -5, 5), color='blue', legend_label="ε=0.2") + \
     plot(f_lorentz_double_prime(x, 0.05), (x, -5, 5), color='green', legend_label="ε=0.05")

p2.show(title="Second Derivative of Lorentzian Approximation")

#Integration of the Lorentzian Sequence
# Compute symbolic integral
var('a b')
assume(a < 0, b > 0)
integral(f_lorentz(x, epsilon), x, a, b).simplify_full()

p1=plot(integral(f_lorentz(x, 0.1), x, -5, 5), (x, -5, 5), color='blue', legend_label="Lorentzian")
p1.show(title="Integrated Delta Approximations")

# Plotting the Integrated Sequences
# To visualize how the cumulative integrals approach step-like functions:

p1=plot(integral(f_lorentz(x, 0.1), x, -5, 5), (x, -5, 5), color='blue', legend_label="Lorentzian")
p1.show(title="Integrated Delta Approximations")

💡 Try It Yourself! Now You can copy and paste directly into here Run SageMath Code Here

Physics Note
Used in resonance and spectral line broadening (e.g. Lorentzian profile in spectroscopy).