学习笔记：透明电子纸的粒子运动仿真模型

学习笔记：透明电子纸的粒子运动仿真模型

---

前言

透明电子纸，又称透明墨水屏，其显示机理是通过TFT的电场调控实现墨水颜色颗粒的移动，进而产生显示画面或透明区域。在国内外研究进展中，透明电子纸技术相关的仿真模型的研究较少，这主要是因为目前还没有较完整以及较精确的粒子运动仿真模型。而构建一个符合实际情况的仿真模型，对于深入研究电泳电子纸技术，或者是电泳、电渗流等物理现象具有很强的实际意义。因此本文将根据所学知识以及相关经验，提出了一种基于Python编写的透明电子纸的粒子运动模型。本人水平有限，请看官们多多包涵。

---

一、粒子运动模型的基本物理背景

电子纸墨水中的颜色颗粒的粒径一般在百纳米级别，根据界面与胶体学相关知识，其颗粒表面会形成“电双层（electric double layer）”结构。该结构的作用是：稳定了颜色粒子在墨水溶剂中的稳定性及分散性，使其不会轻易发生颗粒之间的团聚及吸附现象。“电双层”这一物理结构对于颜色颗粒的运动行为有着至关重要的影响，因此透明电子纸的粒子运动仿真模型也主要是围绕“电双层”这一结构而展开。

根据Peter Debye 和 Erich Hückel等人的相关研究，电子墨水中颜色颗粒的“电双层”结构可由如下表达式表达：
$${\kappa }^{-1}=\sqrt{\frac{{\epsilon }_0{\epsilon }_r{k}_{B}T}{2e^{2}I}}$$
其中${\epsilon }_0{\epsilon }_r$为溶剂的介电常数，$k_B$是玻尔兹曼常数，$T$是绝对温度，$e$是基本电荷，$I$是以量纲$m^{-3}$为单位的离子浓度。颜色颗粒的电泳迁移率和其移动速度由如下表达式表述：

$${\mu }_e=\frac{{2\epsilon }_0{\epsilon }_r\zeta }{3\eta }{f}_1\left(\kappa a\right)$$
$$\vec{v}={\mu }_e\vec{E}$$
其中$a$为粒子半径，$\zeta$为粒子的Zeta电势，$\eta$是溶液的黏度系数，$f_1\left(ka\right)$是一个在1至1.5之间变化的系数，相关的知识及资料均可以在谷歌学术或者专业书籍中找到，这里不再赘述。

二、粒子运动仿真模型

我是用 Python 构建透明电子纸的粒子运动模型，原因是 Python 代码编写自由度大，并且社区里有很多的开源代码及拓展包。同时，用 Python 进行仿真计算所占用的资源少，计算速度快（相对于某些专业多物理场仿真软件而言）。当然，用 Python 建模仿真对使用者的代码能力要求较高，并且还要对电子纸的物理模型有深入了解。那么，接下来本文将详细介绍使用 Python 构建电子纸粒子运动模型的思路及步骤。

1. 导入Python库

导入的Python库大致如下：

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from matplotlib import animation
from mpl_toolkits.mplot3d.art3d import Poly3DCollection
from parameters import *
import matplotlib.pylab as py
from mpl_toolkits.mplot3d import Axes3D
import time

当然，后续如果有其他的需求，如：图形化界面设计、神经网络搭建等，也可以导入相关的包。

2. 修改相关参数及输入仿真文件

在上一篇博客：“学习心得：使用 Comsol 5.6 with Matlab 构建自动仿真模型”中，我介绍了如何自动化批量生成大量电场仿真数据。在这一步骤，我们将这些电场仿真数据（csv格式）导入模型之中，相关代码如下：

df_1 = pd.read_csv('15um方块空间电场0321.csv')
df_2 = pd.read_csv('15um方块空间电场电场翻转.csv')
df = df_1

对 csv 文件中的 x,y,z 电场模分量进行检索：

        self.idx = np.argmin((df[['x', 'y', 'z']] - [x0, y0, z0]).pow(2).sum(1))
        self.E_x = df.loc[self.idx, 'es.Ex (V/m)']
        self.E_y = df.loc[self.idx, 'es.Ey (V/m)']
        self.E_z = df.loc[self.idx, 'es.Ez (V/m)']

同时，下列代码中列举了典型的电子纸墨水的器件材料参数，后续所有的粒子运动计算都将基于这里给出的参数而开展。使用者可以根据自己的需求修改其中的参数。相关参数如下：

import pandas as pd
import numpy as np

# Parameters (All SI units)
# =============================================================================
k_B = 1.38e-23  # Boltzmann constant
T = 300  # Temperature (K)
epsilon_0 = 8.854187817e-12  # 真空介电常数F/m
epsilon_r = 2.01  # 溶剂介电常数
e = 1.6e-19
k_electricity = 9.0e9
I = 0.41153e18  # 分散液电荷浓度
g = 9.8  # m/s2

eta = 1.4e-3  # Dynamic viscosity of water (PaS)
rou_m = 750  # 溶剂密度 kg/m3

# 电泳粒子参数：
a = 400e-9  # 粒子直径
zeta_p = -60e-3  # 白粒子zeta 电位
rou_p = 1120  # 粒子密度 kg/m3
m_p = 4 / 3 * np.pi * (a / 2) ** 3 * rou_p

# 带电反胶束CIMs参数：
a_c = 10e-9  # 粒子直径
q_c = e
zeta_c = q_c / (6 * np.pi * a_c / 2 * epsilon_0 * epsilon_r)

# 非极性溶剂粒子带电量公式
q = 6 * np.pi * epsilon_0 * epsilon_r * a/2 * zeta_p

k = np.sqrt((2 * e ** 2 * I) / (epsilon_0 * epsilon_r * k_B * T))
debye_length = 1/k
print(f'德拜长度为：{debye_length}')

fH_ka = 1/(3*(1+2.5/k*a*(1+2*np.exp(-k*a))))
print(f'亨利修正系数为：{fH_ka}')

mu_e = ((2 * epsilon_0 * epsilon_r * zeta_p)/(3 * eta)) * fH_ka
print(f'电泳迁移率为：{mu_e}')

rho = 1055  # density of polystyrene beads in kg/m3

3. 粒子三维运动速度计算

在透明电子纸的像素单元中，颜色颗粒的运动是沿各个方向的（三维方向），因此需要对粒子的三维方向的运动速度进行分别计算，计算基于了 Debye–Hückel 电泳迁移公式。如下是相关的计算代码：

 # Debye–Hückel迁移率公式，v_x,v_y,v_z的 随着当前坐标的空间电场分量进行变化，
        self.v_x = -self.E_x * (2 * epsilon_0 * epsilon_r * zeta_p * (
                1 + np.sqrt((2 * e ** 2 * I) / (epsilon_0 * epsilon_r * k_B * T)) * a)) / (3 * eta)
        self.v_y = -self.E_y * (2 * epsilon_0 * epsilon_r * zeta_p * (
                1 + np.sqrt((2 * e ** 2 * I) / (epsilon_0 * epsilon_r * k_B * T)) * a)) / (3 * eta)
        self.v_z = -self.E_z * (2 * epsilon_0 * epsilon_r * zeta_p * (
                1 + np.sqrt((2 * e ** 2 * I) / (epsilon_0 * epsilon_r * k_B * T)) * a)) / (3 * eta)

然后，再根据以上计算结果对每颗粒子的空间位置进行更新：

        self.x += self.v_x * dt
        self.y += self.v_y * dt
        self.z += self.v_z * dt

4. 其它物理模块的添加

以上三个步骤仅仅是满足了基本的粒子运动的仿真，但是这样离符合真实情况的粒子运动仿真模型还相距甚远。差距主要在于：

1. 颜色粒子运动至靠近电极界面的时候，会受到来自其它带电物质的“电场屏蔽”作用；
2. 当粒子数量众多（紧密堆积）时，粒子与粒子间会存在库伦排斥等作用，影响粒子的运动轨迹；
3. 颜色粒子的粒径大小不一，其运动速度及运动轨迹也有差别；
4. 电子墨水中的其它带电物质的影响；
   …

本文根据上述第2点，添加了一个粒子间距离判定的模块，其基本思路就是强行让各个粒子间距离保持在两者半径之和之外，倘若发生距离小于两者半径之和，则发生粘弹性碰撞反弹。相关代码演示如下：

    # # 增加电泳粒子间的碰撞检测，0.9为碰撞衰减系数
    def collision_detect(self, other):
        for other in particles:
            if other is not self:
                distance = np.sqrt((self.x - other.x) ** 2 +
                                   (self.y - other.y) ** 2 +
                                   (self.z - other.z) ** 2)

                if distance < (self.r + other.r) * 1:
                    overlap_distance = (other.r + self.r) * 2
                    overlap_direction_x = (other.x - self.x) / distance
                    overlap_direction_y = (other.y - self.y) / distance
                    overlap_direction_z = (other.z - self.z) / distance
                    self.v_x, other.v_x = -self.v_x * 0.9, -other.v_x * 0.9
                    self.v_y, other.v_y = -self.v_y * 0.9, -other.v_y * 0.9
                    self.v_z, other.v_z = -self.v_z * 0.9, -other.v_z * 0.9
                    other.x += overlap_distance * overlap_direction_x
                    other.y += overlap_distance * overlap_direction_y
                    other.z += overlap_distance * overlap_direction_z

以上仅做了一个复杂物理模块的展示，针对上述其余几点的物理现象，使用者可以根据自己的需求进行修改。

5. 粒子运动轨迹动画展示

下列代码展示了使用 python 自带的动画模块 animation 绘制粒子的运动轨迹。使用者亦可根据自己的需求，调整绘图框的格式以及动画的相应参数。

def animate(t):
    ax.clear()
    ax.set_xlim(0, 100e-6)
    ax.set_ylim(0, 100e-6)
    ax.set_zlim(0, 40e-6)
    ax.ticklabel_format(style='sci', axis='both', scilimits=(0, 0))

    for p in particles:
        p.move_p(0.005)  # 更新时间
        for i in range(len(particles)):
            particles[i].move_p(dt)

            for j in range(i + 1, len(particles)):
                particles[i].collision_detect(particles[j])
        ax.scatter(p.x, p.y, p.z, s=100, c='r', marker='o')
         # 时间步长显示
    time_string = ax.text(0.2e-3, 0.5e-3, 0.4e-3, f'Time = {t * 0.005}s', size=20, transform=ax.transAxes)

    # 电极图形：顺时针的四个顶点
    verts = [
        [0.0425e-3, 0.0425e-3, 0.000],
        [0.0425e-3, 0.0575e-3, 0.000],
        [0.0575e-3, 0.0575e-3, 0.000],
        [0.0575e-3, 0.0425e-3, 0.000]]

    collection = Poly3DCollection([verts], alpha=.25)
    collection.set_facecolor('#9400D3')
    ax.add_collection(collection)

    return ax, time_string

ni = animation.FuncAnimation(fig, animate, frames=100)  # 动画总共有100帧

# ani.save('15um方块电极.gif', writer='Pillow', dpi=100)  # 默认和代码存在一个路径下
plt.show()

6. 粒子运动仿真结果展示

最后，经过以上步骤，我们得到了透明电子纸三维的粒子运动仿真结果，粒子运动仿真示意图如下图所示：

通过对上述仿真模型的绘图和动画部分进行稍许修改，在输入合适的参数以及电场仿真文件后，我们便可以得到一组粒子运动轨迹仿真的可视化动图，以及其相对应的结果（如：粒子运动响应时间、粒子三维坐标值等）。根据这些仿真结果，我们可以判定不同的器件结构、墨水参数、驱动波形对粒子运动行为的影响关系。下图展示了一组粒子运动仿真的三视图动画：

---

总结

以上就是关于透明电子纸的粒子运动仿真模型的全部介绍。本模型由 Python 代码编写，粒子的运动计算基于 Debye–Hückel 电泳迁移理论，空间电场需要由 Comsol with matlab 导出。使用者可以根据自己的需求，灵活修改其中电子墨水的材料参数，及仿真结果的判定标准，同时也可以根据自己审美，修改动画绘图的格式，以适应不同的需求。

该粒子运动仿真模型不仅可以模拟电泳粒子的运动轨迹，还可以计算电子纸器件的响应时间、粒子聚集情况、粒子三维坐标等，为电子纸显示屏的设计制造提供一定的参考依据。同时也可为研究电泳、电渗流等物理现象提供一定的仿真基础。

最后，如果对电泳理论以及胶体界面理论的研究足够深入后，粒子运动仿真模型可以较好地模拟真实的粒子运动情况，或许可以考虑搭配 Comsol with Matlab 进行大批量自动化的建模仿真，构建庞大且较为准确AI训练数据集！