学习中心

Course:

Modeling Speaker Drivers in COMSOL Multiphysics^®

返回学习中心

有限元瞬态分析

在扬声器驱动器仿真系列课程的第 1 部分，我们讨论了适用于频域和瞬态分析的情况，让我们来回顾一下什么时候需要运行瞬态模型：

激励为一个短脉冲，或包括一个以上的频率分量。
问题变为非线性，可能原因包括：大信号激励和高次谐波生成、非线性磁化模型，或大变形导致的几何或拓扑变化。
关注输出是对激励的瞬态响应。

在这部分课程中，我们将使用与第 2 部分（侧重于频域研究）相同的动态扬声器驱动器模型来讨论运行瞬态分析时的具体注意事项。本模型在二维轴对称空间维度上求解，主要研究描述驱动器因运动部件的大变形而产生的非线性行为，以及确定系统产生的声学信号的总谐波失真（THD）。

几何结构

用于运行瞬态分析的几何结构如下图所示。用蓝色突出显示的层是为了应用 完美匹配层 特征而创建，以避免为形成有限的计算域而添加的边界产生非物理反射。第 2 部分中使用的 完美匹配边界 特征仅适用于频域压力声学模型，对于时域仿真，需要在连接到开放边界的虚拟域上应用类似特征。

The transient analysis geometry in 2D.

用于瞬态分析的几何结构。为了应用完美匹配层特征吸收外波，必须在开放边界添加一个虚拟域。

几何结构中完美匹配层（PML）的几何厚度需要进行合理的设置，可以通过以下两种网格划分方法来确定：

使用与相邻物理域相同的网格单元大小
使用至少有 6 层网格的结构化网格进行有理数缩放，或使用至少有 8 层网格的结构化网格进行多项式缩放

我们推荐您在物理场和完美匹配层域使用相同的网格大小，是因为完美匹配层是在时域中实现的：在完美匹配层中，相对于实际坐标没有拉伸，因此波在完美匹配层中的传播波长与在物理场中的波长相同。

需要使用的物理场接口

运行瞬态分析需要使用 AC/DC 模块中的磁场接口和声学模块中的 声-固相互作用，瞬态 多物理场接口。请注意，COMSOL^® 中内置了不同的物理场接口，可针对不同类型的流体分析求解波动方程的特定公式。因此，必须确保选择正确的物理场接口以匹配研究类型。下图显示了模型向导中的 选择物理场 窗口中 声-固相互作用，瞬态 接口的位置。

The Select Physics window with the Acoustic-Solid Interaction, Transient interface selected and a description of the feature.

声-固相互作用，瞬态 多物理场场接口位于 模型向导选择物理场 窗口的 声学>声-固相互作用 分支。

该特征将 压力声学，瞬态 接口（瞬态声学模型）和 固体力学 接口进行耦合。与我们在第二部分中进行频域分析的方法类似，声-结构边界 耦合特征会自动添加到模型中，以连接 压力声学，瞬态 和 固体力学 接口。随后，您应该手动将 洛伦兹耦合 特征添加到模型树中的 多物理场 节点，以耦合磁场和 固体力学 接口（与本课程系列的第 2 部分所述方法相同）。

在瞬态模型中使用动网格

洛伦兹耦合 特征的瞬态分析与频域分析的公式有所不同。在频域分析中，计算公式明确包括反向电磁场的贡献：将洛伦兹力从磁场接口传递到 固体力学 接口，并将引起感应电场的结构性速度从 固体力学 接口传递到磁场接口。然而，在时域中，只有洛伦兹力从磁场接口传递到 固体力学 接口，感应电流并不包括在计算中。您需要使用动网格来表征反向电磁场的贡献，这种对瞬态分析的特殊处理是为了避免重复计算反向电磁场。当您想使用动网格来表征由于大变形和几何或拓扑变化引起的非线性时，例如线圈部分或全部在磁隙中移动，就需要使用动网格。由于动网格会自动表征感应电流，因此在耦合公式中加入感应电流会得出错误的结果。

在瞬态研究中，即使位移很小，也必须使用动网格，因为只有这样才能表征感应电流。在本文的示例中，我们还想模拟结构变形对磁场和声场的影响，因此求解的问题又增加了一个非线性源。

如下图所示，单击 动网格 并选择定义工具栏中的 域 >变形域，即可添加 动网格 特征。

The top half of the COMSOL UI, with the Moving Mesh feature selected in the toolbar and the resulting options with Deforming Domain selected.

可从定义 工具栏中添加 动网格 特征。

完成上述操作后，模型树中将添加一个 动网格 接口，请按照以下步骤操作进行正确的设置：

在 变形域 1 子节点中，为 域选择 选择运动部件上方和下方的空气域，并从 网格平滑类型 列表中为平滑选择 拉普拉斯，如下图所示。指定变形域后，使用固体位移分量的指定网格位移将自动添加到固体附近的所有边界，而所有其他外部边界将使用固定边界。

The COMSOL Multiphysics UI showing the Model Builder with the Deforming Domain feature selected and the corresponding Settings window with Laplace selected for the mesh smoothing type.

在变形域中使用 拉普拉斯 网格平滑类型。

与其他非线性平滑方法相比，拉普拉斯 平滑法是计算量最少的方法，因为它是线性的，每个坐标方向只使用一个方程，不需要耦合在一起。这种方法适用于网格变形不太大的情况，尤其是与 自动重新划分网格 特征一起使用时，可确保在网格质量低于给定限值时创建新网格。

添加 对称/辊支承 功能，以在 r = 0 处指定一个对称的边界条件。

The COMSOL Multiphysics UI with the Symmetry/Roller feature selected and the corresponding settings and model.

在 r = 0 处应用 对称/辊支承 边界条件。

时域完美匹配层的设置

多项式 拉伸类型的完美匹配层比例因子和比例曲率参数推荐值分别为 1 和 3，有理数 拉伸类型的推荐值为 1 和 1。在示例中，我们使用 多项式 选项，如下图所示。

The COMSOL Multiphysics UI showing the Model Builder with the Perfectly Matched Layer feature selected, the corresponding settings, and the geometry with the outer layer highlighted in blue to show where the PML feature is applied.

用于瞬态分析的 完美匹配层 设置。

瞬态分析的信号输入

对扬声器进行瞬态分析，需要向音圈施加瞬态电压信号 V(t)，具体输入信号类型的选择取决于所研究的特性。在本例中，输入的是谐波电压信号：

输入信号的频率和振幅分别选择 f₀=70Hz 和 V₀=10V。

在求解瞬态模型时，较好的做法是逐渐增加激励载荷，以尽量减少数值误差。由于系统的惯性效应，瞬时施加的载荷会产生大量高频信号，当这些信号未充分解析时，就会出现数值误差，甚至导致求解器发散。逐渐增加载荷也有其物理意义，因为大多数激励不可能随时间立刻产生对应变化，总是会有一定的延迟。为此，您可以使用有适当平滑功能的阶跃或斜坡函数。

在本示例的分析中，我们使用驱动信号周期 T₀=1/f₀, 来定义斜率为 1/T₀ 的斜坡函数 rm1(t)。如下图所示，输入斜坡函数信号 V₀*sin(2πf₀t)*rm1(t)，定义线圈电压，这样就可以使载荷需要一个周期达到其全幅值。

The Ramp function selected in the Model Builder, the corresponding Parameters and Smoothing settings, and the function plot.

A close-up of the Coil feature selected in the Model Builder and the Coil section of the corresponding Settings window, with a close-up of the model in the Graphics window.

利用驱动信号周期定义 斜坡 函数（上），然后对输入信号进行平滑处理（下）。

物理场设置中控制时间步长的设置

求解瞬态模型时，首先应该确定要求解的最大频率（即，f_max）。请在 压力声学，瞬态 物理场接口设置窗口的 瞬态求解器和网格设置 栏的 要求解的最大频率 框中输入最大频率值 3f₀，如下图所示。此栏还可以进行时间步长的设置，即，可以将 时间步进 方法设置为 固定（首选） 或自由。在本例中，建议使用最适合处理波传播的 固定（首选） 方法。根据此处显示的设置，自动生成的求解器在大多数情况下都能满足要求，前提是网格能够解析模型的频率组成。在这个示例中，我们将 f_max 设置为 3f₀ ，以模拟高达 3 次的谐波。

A close-up of the Pressure Acoustics, Transient feature selected in the Model Builder and the corresponding Settings window showing the Maximum frequency to resolve and Time stepping settings.

物理场的瞬态求解器设置。

研究步骤

与频域分析类似，瞬态研究也执行两个研究步骤。首先，稳态步骤仅求解电磁部分的问题，用于计算驱动器静止时永磁体的磁场。然后，完整的瞬态研究步骤将处理运动扬声器的所有相关多物理场的相互作用。

The Stationary Study selected in the Model Builder and the corresponding Settings window showing the tables in the Physics and Variables Selection setting.

The Time Dependent study selected in the Model Builder and the corresponding Settings window showing settings such as the output times and physics and variables selection.

瞬态分析的研究步骤设置。

瞬态研究步骤在 0 到 T_end = 4T₀ 时间段内求解。完整输出时间范围有两个子区间：从 0 到 3T₀，从 3T₀ 到 T_end。第二个子区间最值得关注，因为假定从时间 3T₀ 开始趋于稳定。因此，为获得时间分辨率更高的解，第二个子区间的输出时步更加精细。

在瞬态研究步骤中启用 自动重新划分网格 选项，可避免使用高度扭曲的网格，因为这可能会导致数值问题。在本模型中，一旦当前网格单元的变形超过预先指定的特定水平，就会创建一个新网格。

A close-up of the Time Dependent study selected in the Model Builder and the Study Extensions settings in the corresponding Settings window.

在瞬态设置 窗口的 研究扩展 栏启用 自动重新划分网格 选项。

通过以下步骤指定重新划分网格的条件：

单击研究工具栏中的 显示默认求解器。
展开求解器序列，单击 瞬态求解器 1 节点下的 自动重新划分网格。
如下图所示，在 自动重新划分网格 设置窗口的 重新划分网格的条件 栏，选择失真作为 条件类型，并在 超出失真范围时停止 文本框中输入 2.5。

The Automatic Remeshing feature selected in the Model Builder and the corresponding Settings window showing the Condition for Remeshing settings, such as Distortion set as the condition type.

重新划分网格的条件 设置 窗口。

这个条件可以确保求解器在最大单元变形的平方根（存储在名为 comp1.spatial.I1isoMax 的变量中）大于 2.5 时重新划分网格。

如下图所示，在 时间步进 设置中，默认瞬态求解器使用 最大频率 的输入值来进行设定以解析最大频率分量。

The Time-Dependent Solver selected in the Model Builder and a close-up of the Time Stepping settings in the corresponding Settings window.

自动生成求解器使用的时间步进。

初始网格

开始时需要先创建一个初始网格，每当达到网格变形容差时，就会自动生成一个新网格。

对于瞬态模型，网格需要能够求解您想要表征的最大频率，例如 f_max。将其转换为最小波长，λ_min = c/f_max，进而获得最大单元大小，h_max < λ_min/5。使用此标准确定每个域的最大单元大小。

初始网格如下图所示。使用映射网格来处理薄实体域，包括线圈、前体、喇叭、音盆和悬架，并使用 自由三角形 网格来处理其余物理域。含 8 层单元分布的映射网格用于离散完美匹配层域。详细的网格划分顺序请参阅随附的 MPH 文件。

A close-up of the initial mesh and the mesh of the half-circle.

用于瞬态分析的初始网格。

此模型在相对较低的频率（包括谐波）下运行，因此不需要对极片中的涡流（集肤深度）进行解析。因此，与频域分析相比，去除了边界层网格。

瞬态解和分析

下方动画显示了音圈间隙内和周围的磁场、音圈（橙色）、前体、喇叭和音盆（粉色）的运动，以及动网格如何在 3T₀ 到 4T₀ 的时间段内使网格适应系统拓扑结构的改变。

在 3T₀ 到 4T₀ 时间段内，音圈间隙内和周围的磁场，以及运动部件的运动和啮合变形。

下图描述了通过计算从原始位置到音圈域内的平均垂直位移得到的音圈的相对位置，如下图所示。在 0.044 s 时，音圈位于最底部（绿色标记）。在 0.048 s（红色标记）和 0.052 s（青色标记）时，音圈分别处于向上和向下运动状态。

A graph plotting the position of the coil as it goes up and then down.

音圈的相对位置。

在 0.044 s（左）、0.048 s（中）和 0.052 s（右）时步，扬声器驱动器周围声压的分布如下图所示。您还可以看到完美匹配层如何衰减由扬声器音盆运动产生的声波。

Three side-by-side images showing the acoustic pressure at different time steps.

3 个不同时步的声压分布。

计算总谐波失真

瞬态分析的一个重要结果是系统产生的声音信号的总谐波失真。失真分析，包括总谐波失真（THD）和互调失真（IMD），通常用于判断扬声器驱动器的质量。这些非线性效应会产生泛音和额外的非谐波频率分量，失真的程度可通过博客：如何对扬声器驱动器进行非线性失真分析中讨论的方法进行量化和计算。要计算总谐波失真，首先需要选择一个感兴趣的点作为聆听点，并获取该聆听点稳态声压的频谱。总谐波失真值的计算公式为

式中，H{sub}N 是 N^th 阶谐波的谐波响应，H₁ 是基波响应。

在此，我们选择物理域中位于声轴的最高点 r = 0 和 z = 0.115 m 处作为聆听点。假设扬声器在 t=3T₀ 时达到稳态，因此在这个特定的聆听点，从 t=3T₀ 到 t=4T₀ 的压力信号与时间的函数关系就是达到稳定状态的解，可用于计算总谐波失真，如下图所示。您可以看到，该曲线与完美的正弦波略有不同，表示输出信号还包含除 f₀ 之外的频率分量。

A graph plotting the pressure at listening point.

聆听点的声压。

请注意，时间间隔越长，快速傅里叶变换（FFT）计算的精度就越高。因此，为了获得更精确的单周期信号频谱，我们对随时间变化的数据进行周期性扩展。具体步骤如下：

使用 点计算 特征，以 T₀/50 为时步，从 3T₀ 到 4T₀ 计算聆听点的压力。这会创建一个表格，存储声压随时间变化的函数。详情请参见随附模型文件中 结果 > 导出值 下的 点压力 节点。
将 结果表 作为数据源，并使用步骤 1 中创建的表格，定义一个插值函数。通过三次样条曲线对信号进行插值。详情请参见 全局定义 下的 插值 1 (p_point) 节点。
使用步骤 2 中指定的插值函数定义一个解析函数，并通过 周期性扩展 使其具有周期性。详情请参见 全局定义 下的 Analytic 1 (p_periodic) 节点。扩展后的压力信号如下图所示。

A graph plotting the period signal as it goes up and down several times.

聆听点压力信号的周期性扩展。

如下图所示，扩展信号的频谱分析由 组件 2 和第二次研究中定义的 全局常微分方程 和 微分代数方程 接口进行。本文示例的研究包括一个瞬态和 时域到频域的 FFT 步骤。第一个研究步骤求解一个全局方程，将周期性扩展信号的值传递给一个名为 P 的因变量；第二个研究步骤计算 P 的 FFT。在示例中，我们选择将信号扩展到 10 个声波周期，以便进行更精确的 FFT 计算。

A close-up of Global Equations selected in the Model Builder and the corresponding Settings window with the Global Equations table for variable P.

压力信号的周期性扩展被传递到一个因变量 P，使用全局常微分方程和代数微分方程 接口求解。

The Model Builder with the Time Dependent study selected and the Settings window showing the Study Settings and Physics and Variables Selection table.

The Model Builder with the Time to Frequency FFT study selected and the Settings window showing the Study Settings and the Physics and Variables Selection table.

用于扩展信号频谱分析的研究步骤。

下图显示了通过 FFT 研究获得的声学信号频谱，即 SPL 随频率变化的值。最高峰值出现在 f₀ 的倍数频率上（红点），值得注意的是奇次谐波的声压级高于偶次谐波的声压级，这表示对称系统的非线性特性优先于非对称系统的非线性特性。

A graph showing the 1/3 octave bands (PSD) as green bar graphs and the continuous PSD as a blue line, with both lines peaking after 50 Hz.

聆听点的声学信号频谱（用 SPL 表示频率）。

在 派生值 > THD 计算 中，使用以下表达式进行 全局计算 计算，得出总谐波失真：

sqrt(sum(with(11 + 10*k, abs(comp2.P)^2), k, 1, 9))/with(11, abs(comp2.P))

The Model Builder with Global Evaluation selected and the corresponding Settings window showing the dataset data and the calculated expression.

使用全局计算 计算总谐波失真。

此表达式使用 with 算子获取特定频率分量的解，使用 sum 算子计算所有索引的表达式之和。具体地说，就是 with(11, abs(comp2.P)) 提供了基频的声压振幅（第 11 次输出频率为 70 Hz），而 sum(with(11 + 10*k, abs(comp2.P)^2), k, 1, 9)（其中 k 为求和指数）计算了前 9 次谐波的压力振幅平方和。如果查看 FFT 研究的输出频率，就会看到基频及其谐波对应的索引。这个包括 9 次谐波计算的总谐波失真等于 2.5%。