Regular Functionals & Green’s Theorem in Several Variables: A Mathematical Insight

Exploring Distributions, Derivatives, and Green's Theorem with SageMath: Extended Examples

When dealing with functions that are piecewise smooth or exhibit discontinuities, traditional calculus reaches its limits. Generalized functions (distributions) provide the tools to rigorously define derivatives even when classical ones fail. In this guide, we’ll use SageMath to explore:

• Distributional derivatives
• Green’s Theorem (with singularities)
• Mollifiers and smoothing
• Visualizations of discontinuous and singular behavior

Recall Generalized Functions: Extending Differentiation & Integration and How to Calculate Distributional Derivatives: Step-by-Step Examples of Piecewise Functions and the Dirac Delta Function

Mathematical Background

Let \(f(x_1, x_2))\) be a function defined on a region \(\quad G \subset \mathbb{R}^2\) , bounded by a piecewise smooth curve \(\quad \Gamma\). In distribution theory, the partial derivative of \(f\) is interpreted weakly as: \[ \left( \frac{\partial f}{\partial x_1}, \varphi \right) = - \left( f, \frac{\partial \varphi}{\partial x_1} \right) \] for all test functions \(\varphi\) If \( f \) is discontinuous, its derivative may include singular terms (like Dirac delta functions), which are fundamental in quantum physics and PDEs.

SageMath Setup

      
import sympy as sp

# Define variables
x1, x2 = sp.symbols('x1 x2')

# Define function
f = sp.Function('f')(x1, x2)

# Compute gradient
grad_f = [sp.diff(f, x1), sp.diff(f, x2)]

print("Gradient (abstract):", grad_f)
	
    

Python code correctly computes the symbolic gradient of f(x1, x2). If you'd like to evaluate the gradient at specific points, you can do this using .subs():

      
import sympy as sp

# Define variables
x1, x2 = sp.symbols('x1 x2')

# Define an explicit function (example: f(x1, x2) = x1**2 + x2**3)
f = x1**2 + x2**3

# Compute gradient
grad_f = [sp.diff(f, x1), sp.diff(f, x2)]

# Evaluate at (x1=2, x2=3)
evaluated_grad = [df.subs({x1: 2, x2: 3}) for df in grad_f]

print("Gradient (symbolic):", grad_f)
print("Gradient at (2,3):", evaluated_grad)
	
    

Region Integration Over 𝐺

Define a square region \( G = [-1,1] \times [-1,1] \) is also the distributional derivative.

      
from sage.all import var, integrate

# Define variables
x1, x2 = var('x1 x2')

# Define function (example function)
f = x1**2 + x2**3  # Replace this with your actual function

# Define integration region
G_region = [(x1, -1, 1), (x2, -1, 1)]

# Perform double integration
integral_G = integrate(integrate(f, x1, -1, 1), x2, -1, 1)

print("Integral of f over G:", integral_G)
	
    

Verifying Green’s Theorem

Green’s theorem relates a line integral around a curve to a double integral over the region it bounds: \[ \oint_{\Gamma} \left( P \, dx_1 + Q \, dx_2 \right) = \iint_{G} \left( \frac{\partial Q}{\partial x_1} - \frac{\partial P}{\partial x_2} \right) dx_1 \, dx_2 \]

      
from sage.all import var, function, diff, integrate

# Define variables
x1, x2 = var('x1 x2')

# Define explicit vector field components (Example: P = x1**2, Q = x2**3)
P = x1**2
Q = x2**3

# Compute curl (∂Q/∂x1 - ∂P/∂x2)
curl_G = diff(Q, x1) - diff(P, x2)

# Apply Green's theorem by integrating over the region G = [-1,1] × [-1,1]
integral_Green = integrate(integrate(curl_G, x1, -1, 1), x2, -1, 1)

print("Green's theorem integral:", integral_Green))
	
    

Extensions with Discontinuities and Distributions

1. Jump Discontinuity via Heaviside

      
f_heaviside = heaviside(x1)
df_dx1 = diff(f_heaviside, x1)  # Yields DiracDelta(x1)
df_dx2 = diff(f_heaviside, x2)  # Equals 0
	
    

SageMath correctly returns DiracDelta(x1) to model the singularity.

2. Green’s Theorem with a Discontinuous Vector Field

Let:

• \[ P(x_1, x_2) = x_2 \cdot H(x_1) \]
• \[ Q(x_1, x_2) = x_1^2 \]

      
P_dis = x2 * heaviside(x1)
Q_dis = x1^2
curl_dis = diff(Q_dis, x1) - diff(P_dis, x2)
integral_curl_dis = integrate(integrate(curl_dis, x1, -1, 1), x2, -1, 1)
print("Green’s theorem result with discontinuity:", integral_curl_dis)  # Should return -2
	
    

Even with a jump in 𝑃, integration smooths the singularity, and Green’s theorem still holds.

3. Mollifiers for Smoothing

A mollifier smooths discontinuities. Define a compactly supported mollifier:

      
from sage.all import var, exp, abs, numerical_integral, plot

# Define the compact mollifier (non-normalized)
def compact_mollifier(x, epsilon):
    return exp(-1 / (1 - (x/epsilon)^2)) * (abs(x) < epsilon)

# Define the normalized mollifier using numerical integration
def normalized_mollifier(x, epsilon):
    # Compute the normalization constant C numerically
    C = numerical_integral(lambda t: float(exp(-1 / (1 - (t/epsilon)^2))) if abs(t) < epsilon else 0.0, -epsilon, epsilon)[0]
    # Return the mollifier divided by the constant, zero outside support
    return (1/C) * exp(-1 / (1 - (x/epsilon)^2)) * (abs(x) < epsilon)

# Example usage: plot the normalized mollifier
x = var('x')
epsilon = 0.5
plot(normalized_mollifier(x, epsilon), (x, -1, 1), color='blue', legend_label='Normalized Mollifier').show()
	
    

Visualizing the Rect Function and Mollification

      
plot_rect = plot(heaviside(x1 + 0.5) - heaviside(x1 - 0.5), (x1, -2, 2), color='red', legend_label='Rect Function')
plot_mollifier = plot(compact_mollifier(x1, 0.5), (x1, -2, 2), color='blue', legend_label='Mollifier')
plot_smoothed = plot(tanh((x1+0.5)/0.5) - tanh((x1-0.5)/0.5), (x1, -2, 2), color='green', linestyle='--', legend_label='Smoothed Rect (tanh approx)')
show(plot_rect + plot_mollifier + plot_smoothed)
	
    

Note: The tanh-based smoothed rectangle is an approximation, not an exact convolution with the mollifier.

Visualizations

• Region 𝐺: Unit Disk

        
  region_G = region_plot(x1^2 + x2^2 <= 1, (x1, -1.2, 1.2), (x2, -1.2, 1.2), color='lightblue')
  region_G.show(title="Region G: Unit Disk")
  	
      

• 3D Discontinuous Function Plot

        
  plot3d(heaviside(x1), (x1, -2, 2), (x2, -2, 2), color='purple', opacity=0.8).show(title="f(x1, x2) = H(x1)")
  	
      

• Gradient Field of Smooth Function

        
  f_smooth = x1^2 + x2^2
  vector_field_plot([diff(f_smooth, x1), diff(f_smooth, x2)], (x1, -1, 1), (x2, -1, 1)).show(title="∇f = ∇(x1² + x2²)")
  	
      

• Conceptual Vector Field with a Jump Discontinuity

        
  def jump_field(x, y):
      return [0, heaviside(x)]  # y-component has a jump
  
  vector_field_plot(jump_field(x1, x2), (x1, -1, 1), (x2, -1, 1)).show(title="Jump Discontinuity in Vector Field")
  	
      

Conclusion

Through symbolic differentiation, region integration, and visualization, SageMath proves powerful for exploring generalized functions, singularities, and foundational theorems like Green’s. We’ve:

• Explored distributional derivatives (Heaviside → DiracDelta)
• Verified Green’s Theorem with singular fields
• Simulated smoothing with mollifiers
• Visualized discontinuities and their effects

These tools are essential for advanced applications in mathematical physics, PDEs, and engineering, where idealizations meet real-world irregularities.

👉 Run SageMath Code Here