一维多孔介质域污染物迁移解析模拟：傅里叶-快速傅里叶变换(Fourier-FFT)方法

《Geosciences》：Analytical Modelling of Contaminant Transport in One-Dimensional Porous Medium Domains: The Fourier-FFT Approach

　　

多孔介质中污染物迁移的解析解对于理解地下过程及验证数值模型具有重要意义。然而，传统的基于拉普拉斯变换(Laplace-transform)的方法往往难以处理实际的瞬态边界条件(transient boundary conditions)，且通常导致复杂的逆变换形式，需要计算量大的卷积积分。为解决这些局限，本文提出一种傅里叶-快速傅里叶变换(Fourier-FFT)解析框架，用于求解均质与非均质一维多孔介质中经典的一维对流-弥散-反应(Advection–Dispersion–Reaction, ADR)方程。该傅里叶-FFT方法具有数学推导简洁、计算速度快的优点，并能方便地引入超越理想阶跃或脉冲输入的真实瞬态边界条件。研究人员通过对均质域中基于拉普拉斯变换的解析解以及双渗透率(dual-permeability)域中有限元解的对比验证，表明二者吻合良好，证实了方法的准确性。参数分析进一步证明该框架能在不同水文地质及反应条件下捕捉污染物迁移的预期物理行为。

论文解读：一维多孔介质域污染物迁移解析模拟——Fourier-FFT方法

研究背景与意义

地下水中污染物迁移是环境与健康领域的重大关切问题。准确预测污染物在多孔介质中的行为是风险评估与修复设计的基础，通常通过求解描述对流(advection)、弥散(dispersion)及反应(reaction)过程的偏微分方程即一维对流-弥散-反应方程(One-Dimensional Advection–Dispersion–Reaction Equation, ADR方程)来实现。现有解析方法多采用拉普拉斯变换(Laplace transform)或格林函数法，但其逆变换常涉及复杂半无限积分与卷积运算，且通常仅适用于理想的狄拉克脉冲或阶跃边界输入，难以匹配实际工程中任意时变(transient/time-dependent)的边界条件，数值求逆亦存在稳定性与收敛性问题。相比之下，傅里叶变换结合快速傅里叶变换(Fast Fourier Transform, FFT / IFFT)可在频域将时间导数转化为代数项，避免繁琐的逆变换过程，天然适配任意时序边界信号。本文在《Geosciences》发表，旨在建立一套系统化的Fourier-FFT解析求解框架，克服传统Laplace方法的局限，并验证其在均质及双重渗透率(dual-permeability)非均质多孔介质中的适用性与精度。

主要关键技术方法

研究人员针对一维ADR方程分别在时间域施加连续Fourier变换（将?/?t转换为iω·?(ω,x)），空间域保留常系数ODE求解，得到频域解析解?(ω,x)后采用Inverse Fast Fourier Transform(IFFT)反变换回时间域获得c(t,x)。研究涵盖：(1)半无限及有限长均质域，考虑无源汇、零级产污(zero-order production/γ)、一级衰变(first-order decay/λ)及吸附(retardation factor/Rd = 1+ρ_bK_d/θ)情形，搭配Dirichlet、Neumann及Robin三类边界条件；(2)双渗透率(dual-permeability/fracture-matrix)耦合域的一维半无限问题，引入域间传质系数(transfer coefficient/α)。模型通过van Genuchten & Alves基于Laplace的近似解析解（均质域）及COMSOL Multiphysics有限元解（双渗透率域）进行验证，参数分析考察产污、衰变及分配系数(K_d)影响。时间离散采用N=2ⁿ点FFT采样（Δt=0.1 d）。

研究结果

2. 一维均质多孔介质中污染物迁移(One-Dimensional Contaminant Transport in Homogeneous Porous Domains)

研究人员推导了含吸附修正的ADR方程：R_d·?c/?t = D·?²c/?x²- v·?c/?x + γ - λc，其中R_d=1+ρ_bK_d/θ为阻滞因子(retardation factor)。

2.1 无产污与衰变的污染物迁移(Contaminant Transport with No Production or Decay)

• •
  2.1.1 半无限一维域Dirichlet边界(Semi-Infinite 1D Domain with Dirichlet Boundary Condition)：对半无限域x≥0施加c(0,t)=c₀(t)，Fourier变换后得通解?(ω,x)=A(ω)exp[((v/D)-√((v/D)²+4iωR_d/D)/2)·x]，由有界性排除正向指数项，A(ω)=?₀(ω)，经IFFT还原c(t,x)。结论：该方法可处理任意时变入口浓度。
• •
  2.1.2 有限一维域Dirichlet–Neumann边界(Finite 1D Domain with a Dirichlet–Neumann Boundary Condition)：域长L，c(0,t)=c₀(t)，?c/?x|_x=L=0（隔水端）。频域解为?(ω,x)=?₀(ω)·cosh[λ(ω)·(L-x)]/cosh[λ(ω)L]，λ(ω)=√((v/D)²+4(iωR_d/D))/2 - v/(2D)（取实部保证衰减）。IFFT得时域解。结论：能正确模拟封闭端无通量条件。
• •
  2.1.3 有限一维域Robin–Neumann边界(Finite 1D Domain with Robin–Neumann Boundary Conditions)：入口处-D?c/?x+v c = h·(c_s(t)-c)（第三类Robin边界代表与环境交换），出口?c/?x|_x=L=0。代入频域通解并由边界条件联立解得系数。结论：可模拟具边界交换反应的受限系统。

2.2 含零级产污的污染物迁移(Mass Transport with Zero-Order Production)

非齐次ODE分解为齐次解+特解（常数特解?_p=γ/(iωR_d+λ_decay→0+...)，本节λ暂为0）。

• •
  2.2.1 有限域零级产污Dirichlet–Neumann(Finite 1D Domain with Zero-Order Production: Dirichlet–Neumann Boundary Conditions)：频域解含齐次部分与特解γ/(iωR_d)，由两端边界定常数。半无限极限L→∞退化为?(ω,x)=?₀(ω)exp[μ(ω)x]+γ/(iωR_d)，其中μ(ω)为负实部特征值。结论：零级源项可通过频域特解法自然纳入。
• •
  2.2.2 有限域零级产污Robin–Neumann(Finite 1D Domain with Zero-Order Production: Robin–Neumann Boundary Conditions)：同理叠加Robin条件，给出含β=(vh+Dα_Robin)/D组合的闭式频域解。结论：框架兼容混合边界与内源项。

2.3 含零级产污与一级衰变的污染物迁移(Contaminant Transport with Zero-Order Production and First-Order Decay)

• •
  2.3.1 有限域稳态Dirichlet–Neumann(Steady-State Contaminant Transport in a Finite 1D Domain with Dirichlet-Neumann Boundary Conditions)：令?c/?t=0得D·c''-v·c'+γ-λc=0，解析解c(x)=γ/λ + [c₀-γ/λ]·(exp(αx)-exp(αL))/(1-exp(αL))，α=v/D·(1-√(1+4λD/v²))/2（取合适根）。L→∞时退化为c(x)=(γ/λ)+[c₀-(γ/λ)]exp(Pe·x/L·(1-√(1+4Da))/... )简化形式。结论：稳态解与经典结果一致。
• •
  2.3.2 有限域Robin–Neumann瞬态及算例(Finite 1D Domain with Robin–Neumann Boundary Conditions)：频域非齐次方程特解?_p=γ/(iωR_d+λ)，齐次解特征根含λ_ADR=iωR_d+λ，由Robin及Neumann边界确定常数。算例（L=150 m, ρ_b=1600 kg/m³, θ=0.3, K_d=0.001或0.01 m³/kg, D=1.5 m²/d, v=0.5 m/d, γ=0.1 mg/(m³·d), λ=0.001 d^-1, 钟形入口c₀(t), N=2¹⁴FFT, Δt=0.1 d）显示：产污提高域内浓度并延长持留；衰变降低浓度；大K_d（强吸附）增强相分配导致溶质滞后与展宽。结论：Fourier-FFT框架能合理复现产污、衰变及吸附对羽状污染物时空分布的影响。

2.4 模型验证(Model Verification)

均质有限域对比van Genuchten & Alves (1982) Laplace近似解：无产污/衰变情况完全吻合；含γ=0.1 mg/(m³·d)及λ=0.001 d^-1时基本吻合，微小偏差源于van Genuchten模型中γ在全域瞬时激活之假设差异。结论：Fourier-FFT结果与经典Laplace基准解一致性良好，验证了方法准确性。

3. 双渗透率域一维污染物迁移(One-Dimensional Contaminant Transport in Dual Permeability Domains)

裂隙(fracture, 下标f)与基质(matrix, 下标m)分别满足R_d,f?c_f/?t = D_f?²c_f/?x²- v_f?c_f/?x - λc_f- α(c_f-c_m) 与 R_d,m?c_m/?t = D_m?²c_m/?x²- v_m?c_m/?x - λc_m+ α(c_f-c_m)，α为域间传质系数(mass transfer coefficient)。Fourier变换后得四阶特征多项式，由有界性剔除正实部指数，通过特征向量关联两域振幅，由x=0处c_f(0,t)=c₀(t)及c_m(0,t)=c₀(t)（或指定）定解。验证对比COMSOL（TDS-PM接口, N=2¹⁵, Δt=0.1 d）：(i)无吸附无衰变、(ii)有衰变无吸附、(iii)有吸附(K_d=0.005 m³/kg)及衰变三种情形Fourier-FFT与有限元结果均吻合良好。参数分析表明：增大K_d强化阻滞使裂隙与基质突破延迟且浓度降低；增大λ加速衰减且下游效应更显著。结论：所提框架成功扩展至裂隙-基质双渗透率耦合系统，能刻画优势通道流、基质慢迁及域间质量交换引起的羽流衰减与延迟。

讨论与结论(Conclusions)

研究人员指出，与传统Laplace变换方法相比，Fourier-FFT方法避免了复杂的逆Laplace积分与卷积运算，IFFT算法成熟高效，数学处理更为简捷；且该方法可直接接受任意实测时序边界信号，适合真实污染场景。经均质域Laplace基准解与双渗透率域有限元解交叉验证，精度可靠。框架适用于均质及非均质（dual-permeability）一维ADR问题，能纳入吸附、零级产污、一级衰变及各类(Dirichlet/Neumann/Robin)边界条件。参数分析证实方法可正确反映阻滞(retardation)、产污(production)与衰变(decay)对污染物迁移与突破曲线(breakthrough curve)的物理影响。该Fourier-FFT解析框架为地下水污染物迁移的快速解析预测、数值模型验证及教学演示提供了实用工具，降低了高阶数学门槛，便于工程师与应用科学家针对具体边界条件定制求解。