借助仿真 App 评估静态混合器的性能

2016年 7月 8日

静态混合器因其高效、成本低廉、易于安装且维护要求低等优点,成为各类工程领域的常用工具。在评估混合器能否满足某种使用目的时,一个重要的判断指标是得到的混合物是否足够均匀。在本篇博客文章中,我们将介绍如何借助“粒子追踪模块”,开发一款能定量和定性分析静态混合器性能的 App。

层流静态粒子混合器设计器 App 的设计基础

在创建“层流静态粒子混合器设计器”App 前,我们首先前往官网的“案例下载”页面下载层流混合器教程。此模型通过计算贯穿设备的粒子轨迹对静态混合器的混合性能进行了评估。如果您希望深入了解此教程及学习基本的混合器建模知识,我们推荐您阅读下列博客文章:

上述教程使用了与我们的 App 相同的几何模型。如下图所示,模型由一段装有三个螺旋叶片的管道构成,叶片可交替旋转。搅拌叶片的表面呈灰色,周围勾勒出管道的轮廓线。流体运送粒子通过管道时,静态搅拌叶片会使粒子混合在一起。

图像显示了层流静态混合器模型的几何结构。
层流静态混合器模型的几何结构。

用于研究静态混合器性能的 App

借助下图中的层流静态粒子混合器设计器 App,我们可以先计算出通过混合器的粒子轨迹,然后利用一些内置后处理工具定量及定性地对混合器的性能进行评估。

图像展示了层流静态粒子混合器设计器的用户界面。
“层流静态粒子混合器设计器”的用户界面(UI)截图。

此 App 中包含大量的几何参数和材料属性,在创建混合器模型时,可选择使用一个、两个或三个螺旋式搅拌叶片。除此之外,通过 App 的高级设置,您还能更改粒子数量和后处理参数。

为了更准确地绘制静态混合器内不同物质的分布及变化,我们可以释放粒子,并使用流体流动粒子追踪 接口对其轨迹进行计算。粒子位置是通过牛顿运动方程公式计算出来的,而计算位置矢量的分量则需求解一组二阶方程:

(1)

\frac{d}{dt}\left(m_p\frac{d\mathbf{q}}{dt}\right) = \mathbf{F}_t

其中 \mathbf{q}(国际单位:m)表示粒子位置,m_p(国际单位:kg)表示粒子质量,\mathbf{F}_t(国际单位:N)表示粒子上的总力。牛顿公式考虑到了粒子的惯性,并允许粒子横穿速度流线。

在模型中,粒子只受到曳力的作用,并使用 Stokes 曳力定律计算得到:

(2)

\begin{aligned}
\mathbf{F}_D &= \frac{1}{\tau_p} m_p \left(\mathbf{u}-\mathbf{v}\right)\\
\tau_p &= \frac{\rho_p d_p^2}{18 \mu}
\end{aligned}

其中应用了下列物理量:

  • \mathbf{v}(国际单位:m/s)表示粒子速度
  • \mathbf{u}(国际单位:m/s)表示流体速度
  • d_p(国际单位:m)表示粒子直径
  • \rho_p(国际单位:kg/m^3)表示粒子密度
  • \mu(国际单位:kg/(m\cdot s))表示粒子动力粘度

Stokes 曳力定律适用于相对雷诺数远小于 1 的粒子;也就是,

(3)

\textrm{Re}_r = \frac{\rho \parallel {\mathbf{u}-\mathbf{v}}\parallel d_p}{\mu} \ll 1

其中 \rho(国际单位:kg/m^3)表示流体密度。本文的案例符合定律要求。下图描绘了溶液中具有代表性的粒子样本。这些粒子从混合器的右下角释放并流向左上角。其颜色代表粒子在入口处的初始 z 坐标,可直观地对比混合器横截面中粒子的最终位置和初始位置。

绘图显示了层流静态混合器内的粒子轨迹。
静态混合器内粒子轨迹的绘图。

借助 App 开发器量化静态混合器的性能

肉眼观察在一定程度上可以判断混合物的均匀程度,针对本文的案例,混合性能的可视化能够通过创建粒子位置的相图来实现。在相图中,我们可以在任意的二维相空间中绘制粒子轨迹——也就是说,可以通过用户自定义的表达式在二维绘图的坐标轴上表示粒子轨迹。相图常用于绘制特定方向上粒子位置随动量的变化,即相空间 分布。

在下方动画中,相图用于观察每个粒子在通过混合器的过程中横向位置的变化。管道设为 y 方向,因此横向为 x 和 z 方向。我们使用不同的颜色来表示初始状态下每种粒子所占据的象限。从图中可以看出,初始时释放的深蓝色粒子位于 xz 轴的正方向。

相图显示了粒子在通过混合器过程中的横向位置。

相位图定性地表面出口处的粒子没有完全均匀地混合。图中仍存在粒子数量密度偏高或偏低的区域,同时还能看到相同颜色的粒子(来自同一象限)聚集在了一起。

相图存在一个潜在的缺陷,这是因为它绘制的是同一时刻相空间内的粒子,而非基于相同 y 坐标的粒子。生成的混合器可视化绘图也因此带有一定的误导性,这因为一些更靠近叶片的粒子可能会晚于其他粒子到达出口。因此我们提出了另一种方法:创建庞加莱映射,该映射可绘制出粒子轨迹与特定位置上平面的交点。

下图中,每个截面上粒子的颜色取决于释放时该粒子的初始 x 坐标,蓝色表示初始坐标为正方向,红色表示为负方向。接下来我们再一次观察出口处红蓝粒子的聚集情况。

视图显示了庞加莱映射。
庞加莱映射显示了二维绘图中粒子的位置。

我们可以从相图与庞加莱映射中获取大量有关混和器性能的信息,然而对于精密的工业应用而言,多数信息过于主观。肉眼观察可以大概判断不同物质处于完全分离、部分混合或充分混合的状态,但这些说法太模糊,难以量化。举例来说,任何一个观察者都能看出上图中聚集着一小群同色粒子,但要指定数值来描述粒子混合程度就十分困难了。

幸运的是,“App 开发器”和“方法编辑器”提供的工具可以帮助我们创建专业的高端后处理程序,这些程序可以对特定混合器几何模型的性能进行赋值。分布指数是一个常用于评估粒子空间分布均匀性的度量标准,它被定义为方差和平均数之比:

(4)

D=\frac{\sigma^2}{\mu}

我们将出口划分为面积相等的区域或象限,然后计算平均数和方差。因为出口为圆形,所以可以通过绘制多个半径不同的同心圆来将出口划分成面积相等的 N_r 个环状区域

r_i = \sqrt{\frac{i}{N_r}} \hspace{1cm} \textrm{for } i=1,2,3\ldots N_r-1

通过绘制多条夹角相等的直径,每个环状区域又可以分割成 N_{\phi} 个面积相等的域

\phi_j = \frac{2\pi j}{N_{\phi}} \hspace{1cm} \textrm{for } j=0,1,2\ldots N_{\phi}-1

这些细分的区域会产生 N_q=N_r N_{\phi} 个面积相同的象限。x_i 表示第 i 象限内的粒子数量,则每个象限内粒子的平均数量为

\bar{x}=\frac{1}{N_q}\sum_{i=1}^{N_q} x_i

每个象限内粒子数量的方差为

\sigma = \frac{1}{N_q}\sum_{i=1}^{N_q} (x_i-\bar{x})^2

我们在 App 中使用了(p_computeIndexOfDispersion)方法来计算分布指数。

/*
 * p_computeIndexOfDispersion
 * This method computes the index of dispersion at the outlet.
 * The method is called in p_initApplication and in m_compute.
 */

// Get the x- and z-coordinates of the particles at the outlet
// and store them in matrices qx and qz, respectively.
model.result().numerical().create("par1", "Particle");
model.result().numerical("par1").set("solnum", new String[]{"14"});
model.result().numerical("par1").set("expr", "qx");
model.result().numerical("par1").set("unit", "m");
double[][] qx = model.result().numerical("par1").getReal();
model.result().numerical("par1").set("expr", "qz");
model.result().numerical("par1").set("unit", "m");
model.result().numerical("par1").set("solnum", new String[]{"14"});
double[][] qz = model.result().numerical("par1").getReal();

// Use the "at" operator to get the initial x-coordinates of all particles
// and store them in matrix qx0.
model.result().numerical("par1").set("expr", "at(0,qx)");
model.result().numerical("par1").set("unit", "m");
model.result().numerical("par1").set("solnum", new String[]{"14"});
double[][] qx0 = model.result().numerical("par1").getReal();

// The Particle Evaluation is no longer needed.
model.result().numerical().remove("par1");

double Ra = model.param().evaluate("Ra"); // Radius of the outlet
int Np = qx.length;                      // Number of particles
int Nr = nbrR;                            // Number of subdivisions in the radial direction
int Nphi = nbrPhi;                        // Number of subdivisions per quadrant in the azimuthal direction
int nbrQuad = Nr*4*Nphi;                  // Total number of quadrats (regions)
double deltaPhi = Math.PI/(2*Nphi);       // Angular width of each quadrat
int index = 0;
int ir = 0;
int iphi = 0;
int[] x = new int[nbrQuad]; // Array to store number of points per quadrat

// Begin loop over all particles
for (int i = 0; i < Np; ++i) {
  // Determine which quadrat each particle is in.
  ir = (int) Math.floor((Math.pow(qx[i][0], 2)+Math.pow(qz[i][0], 2))*Nr/Math.pow(Ra, 2));
  iphi = (int) Math.floor(Math.atan2(qz[i][0], qx[i][0])*Math.signum(qz[i][0])/deltaPhi);
  if (Math.signum(qz[i][0]) < 0) {
    iphi = (int) Math.floor((2*Math.PI-Math.atan2(qz[i][0], qx[i][0])*Math.signum(qz[i][0]))/deltaPhi);
  }
  index = 4*Nphi*ir+iphi;
  // Consider only half of the particles when evaluating mixer performance.
  if (qx0[i][0] < 0) {
    x[index] = x[index]+1;
  }
}
// compute the mean
double sum = 0;
for (int i = <0; i < nbrQuad; ++i) {
  sum += x[i];
}
// compute the variance
double xmean = sum/nbrQuad;
sum = 0;
for (int i = 0; i < nbrQuad; ++i) {
  sum += Math.pow(x[i]-xmean, 2);
}
indexOfDispersion = sum/xmean;

上述方法的最后一行用于返回分布指数。通常情况下,分布指数减小意味着粒子会分布得更加均匀。根据 App 的默认参数,使用三个叶片时分布指数约为 900,使用两个叶片时为 1200,只用一个叶片时为 1400。因此,分布指数定量地说明了绘图显示的信息:使用的搅拌叶片越多,粒子混合物的均匀度就越高。

借助仿真 App 优化静态混合器的性能分析

今天,我们演示了如何利用“App 开发器”来优化静态混合器的研究。通过创建 App,您可以将仿真能力分享给更广泛的受众,从而优化整个设计工作流,与此同时,还能向不同的混和器几何指定数值,从而更精确地掌握设备的混合性能。

有兴趣深入了解如何设计仿真 App 吗?请浏览下方资源。

指导 App 开发的实用资源


评论 (0)

正在加载...
浏览 COMSOL 博客