电泳迁移率拓展猜想-优快云博客

⚠️ 声明：猜想类模型，尚未严格证明，需要数值模拟标定参数。

简介

传统电泳迁移率模型有一组常用简化前提：颗粒为标准球形、德拜层极薄、惯性效应忽略、表面电势均匀，形状扰动仅保留二阶勒让德项。
这套拓展模型针对传统框架五个薄弱点做补充：

放弃表面电荷均匀假设，支持表面电荷分异分布；
不直接抹除高阶形状勒让德项，给出高阶项衰减猜想；
不只使用无限薄德拜层极限，给出有限德拜层显式修正；
纳入低雷诺数下被丢弃的惯性小量；
引入收敛判据，判断级数展开是否发散。

核心公式

E1 表面电荷分域积分

(\mu_\mathrm{charge}= \frac{\varepsilon\zeta}{\eta}\int_{\partial\Omega}\omega_q(\boldsymbol{x}),\hat{\boldsymbol{n}}(\boldsymbol{x}),\mathrm{d}S)

当表面电荷密度(\omega_q)取常数时，退化为经典 Henry 迁移率形式。

E2 高阶形状项

(\mu_\ell = c_\ell(\kappa a)\cdot\varepsilon\cdot\langle P_\ell(\cos\theta)\rangle_\mathrm{surface},\quad \ell\ge0)

猜想：(\ell=2)项保留；(\ell\ge3)项在(\kappa a\to\infty)条件下趋于 0；(c_\ell(\kappa a))是依赖德拜数的系数，需要数值扫描确定衰减行为。

E3 有限德拜层渐近展开

(\mu(\kappa a)=\mu_\infty\left[1-\frac{\alpha_1}{\kappa a}-\frac{\alpha_2}{(\kappa a)^2}+O\left(\frac{1}{(\kappa a)^3}\right)\right])

(\kappa a\to\infty)回归 Smoluchowski 极限(\mu_\infty)；(\alpha_1,\alpha_2)为待标定修正系数。

E4 颗粒惯性修正

(\mu_\mathrm{inertia}= \beta\cdot Re_p\cdot \mu_0,\quad Re_p=\frac{\rho_f a^2 E}{\eta})

(Re_p\to0)惯性修正归零，回到 Stokes 低雷诺数假设。

E5 收敛判据

(N_\mathrm{core}= \omega_1X_1+\omega_2X_2+\omega_3(\kappa-X_3)+\omega_4X_4-\omega_5X_5-\omega_6|X_6|,\quad \kappa=\frac{4}{\pi})

(N_\mathrm{core}>0)：级数展开自洽收敛；(N_\mathrm{core}\le0)：展开发散，高阶效应不可忽略。

通俗比喻

可以把带电颗粒电泳想象成小船在河里漂流：

E1：船体表面不是整块均匀带电，有的区域电荷密、有的区域疏，不同区域对漂流速度贡献不一样；
E2：船体不是完美圆球，有凹凸形变；二阶形变会实实在在改变漂流速度，但更高阶细碎凹凸的影响会逐步衰减；
E3：船体周围带电的 “影响水层” 不是无限薄，水层有厚度，厚度会拖慢航行速度；
E4：水流很慢的时候惯性几乎看不见，但颗粒尺寸、电场拉高之后，惯性会带来微小偏移；
E5：综合所有误差、耦合程度，判断整套计算方式还能不能够用，什么时候这套近似会彻底失效。

#include <stdio.h>
#include <math.h>
.

double legendre(int l, double x){
    if(l==0) return 1.0;
    if(l==1) return x;
    double p0=1.0, p1=x, p2;
    for(int k=2;k<=l;k++){
        p2=((2*k-1)*x*p1-(k-1)*p0)/k;
        p0=p1; p1=p2;
    }
    return p1;
}

double surface_avg_Pl(int l, double eps, int n){
    double sum=0.0, dtheta=PI/n;
    for(int i=0;i<n;i++){
        double th=(i+0.5)*dtheta;
        double r=1.0+eps*legendre(2,cos(th));
        double w=r*r*sin(th)*dtheta;
        sum += legendre(l,cos(th))*w;
    }
    double area=0.0;
    for(int i=0;i<n;i++){
        double th=(i+0.5)*dtheta;
        double r=1.0+eps*legendre(2,cos(th));
        area += r*r*sin(th)*dtheta;
    }
    return sum/area;
}

double mu_finite_Debye(double mu_inf, double kappa_a, double a1, double a2){
    if(kappa_a<=0) return 0.0;
    return mu_inf*(1.0 - a1/(kappa_a) - a2/(kappa_a*kappa_a));
}

double mu_inertia(double mu0, double Rep, double beta){
    return beta*Rep*mu0;
}

double N_core(double X[6], double w[6]){
    return w[0]*X[0] + w[1]*X[1] + w[2]*(KAPPA - X[2]) + w[3]*X[3]
         - w[4]*X[4] - w[5]*fabs(X[5]);
}

int main(){
    printf("=== EPH‑E1~E5 验证 ===\n\n");
    int n=20000; double eps=0.1;
    printf("[E2] 表面平均勒让德 ⟨P_l⟩\n");
    for(int l=0;l<=5;l++){
        double avg=surface_avg_Pl(l,eps,n);
        printf("  l=%d  ⟨P_l⟩ = %+.10f\n", l, avg);
    }

    printf("\n[E3] 有限Debye层展开\n");
    double mu_inf=1.0, a1=1.0, a2=0.1;
    for(int k=0;k<=4;k++){
        double ka=10.0*pow(10.0,k);
        printf("  κa=%.0f  μ/μ∞ = %.8f\n", ka, mu_finite_Debye(mu_inf,ka,a1,a2));
    }

    printf("\n[E4] 惯性修正\n");
    for(int k=0;k<=4;k++){
        double Rep=1e-6*pow(10.0,k);
        printf("  Re=%.2e  δμ/μ0 = %.3e\n", Rep, mu_inertia(mu_inf,Rep,0.1));
    }

    printf("\n[E5] N_core收敛判据\n");
    double X[6]={0.88,0.50,0.55,0.85,0.65,0.45};
    double w[6]={0.15,0.25,0.20,0.15,0.15,0.10};
    double N=N_core(X,w);
    printf("  N_core = %.6f\n", N);
    printf("  判定: %s\n", N>0?"收敛可行":"展开发散");
    return 0;
}