Sunday, August 21, 2022

使用WELSIM生成MFEM初始网格文件

 最近几年比较活跃的开源有限元求解器要数MFEM (mfem.org)了。此开源包虽然定位于轻量级可扩展的C++程序求解库,但提供了高阶有限元空间,支持混合单元、非连续伽辽金单元、等几何分析方法等等。尤其在高性能计算上有巨大优势,不仅支持消息传递机制并行(MPI)和共享内存并行(OpenMP),同时对GPU并行计算也有很好的支持。自带的后处理显示GLVis程序可以方便读取显示结果文件。BSD协议对于开发者也是极为友好。最近MFEM支持了Python编程,对各种类型的科研人员显得更加友好。

笔者10多年前就开始关注MFEM这个开源包,当时只能进行一些简单的电磁场计算。随着时间的发展,功能越来越稳定,支持的线性求解器也越来越丰富。如hypre和PETSc,特征值求解器SLEPc,各种并行的迭代求解器Ginkgo都已很好支持。支持众多不同类型的单元和网格自适应功能。程序也提供了大量的实例,可以帮助使用者快速理解并运用此程序。

和很多偏微分方程求解器一样,MFEM需要用户提供网格文件,初始网格文件是不可或缺的输入文件。MFEM的方便之处在于,用户可以选择使用简单粗糙的初始网格,通过求解中的网格再划分功能,可以在计算中生成更加优化的网格。MFEM的网格格式清晰易懂,包含了dimension, elements, nodes, boundary 等信息块。其官网提供了详细的网格格式注解。

尽管MFEM有着众多优异的特性,并且网格文件格式清晰易懂。但对于复杂几何模型,很难通过手动编写完成。目前社区也没有一个款好用的初始网格生成工具,这也一定程度上限制了其在复杂几何模型上的应用。WELSIM不仅具有简单好用的前端界面,同时支持几何到网格的生成,并导出各种格式的网格文件。用户可以快速获得复杂模型的网格,并直接用于MFEM的计算求解。只需要在WELSIM中选择输出网格文件为MFEM格式即可。

在很多科学计算或求解偏微分过程中,网格文件只含有节点和单元信息是不够的。我们还需要提供边界条件所选择的边界网格信息。WELSIM也提供了很好的解决方案,用户只需要在GUI中将分析类型设置为电磁场(Electromagnetics),添加一些边界条件如电压等,同时选择好对应的几何模型边界。选择导出输入文件(Export Input Scripts)。

此时就会生成多个求解器相关的文件,其中网格文件中不仅含有节点和单元信息,还有边界条件中所需要的边界单元信息。

通过这种方式,用户可以快速生成用于MFEM求解器的复杂网格文件,从而进行更为复杂的仿真研究,或者更贴近实际工程应用的计算。

目前WelSim支持自动化生成Tet4/10和Tri3/6网格。

注:WELSIM不隶属于MFEM,与MFEM开发团队没有直接关联。这里引用“MFEM”仅用作技术博客文章与软件使用的参考。

Thursday, August 11, 2022

曲线拟合的核心计算方法与CurveFitter的更新

 随着大数据技术与计算科学的发展,数据挖掘与分析工具发挥着越来越重要的工作。对于常见的二维曲线和三维曲面数据,通过给定的曲线/曲面方程和参数,可以将复杂的数据用精简抽象的函数表达出来,迅速且正确得出相关重要信息。目前曲线拟合已经大量应用于工程与科学研究,涉及社会学,医学,工学,生物学等各个领域。

WELSIM早在2年前就发布了免费的CurveFitter工具软件,获得了很多好评与热度。最早的CurveFitter版本是为了计算结构力学中超弹材料的测试数据-曲线拟合,以及磁性材料中损耗的数据-曲线拟合而开发。之后按用户的要求添加了许多常用曲线,如多项式,幂函数,指数函数等等。这些曲线的大多都是非线性的,只有通过计算机才能快速的准确的得到相应参数。

最小二乘法(Least Squares)

在计算机求解曲线拟合时,最常用的方法是最小二乘法。这种计算方法广泛应用于曲线拟合以及优化计算中。算法的核心是在给定的数据区间内找到泛函的全局最小值。通过构建雅克比矩阵,以及每个步长下的线性化,获得以参数为自变量的线性方程组,求解收敛后得到我们需要的参数值。

对于求解非线性问题的通用方法是用线性化方法将其转化为一系列假设并求解。在每次迭代中,假设被求解来决定正确的步长。为了更快更准确地收敛,步长的控制尤为重要。非线性优化算法(又称最小化方法)主要可以分为两类。1. 信任区域法(Trust Region)。2. 线性搜索法(Line Search)。这两种方法有很多地方相似。信任区域法先选择一个步长,再决定步长方向。而线性搜索法相反,先选择一个方向,再决定步长。通常对于计算量不大的曲线拟合问题,信任区域法是一个更好的选择。

Levenberg-Marquardt算法

Levenberg-Marquardt算法(以下简称LM法)是求解非线性最小二乘法的最常用的方法。LM算法是在 1960 年代初开发的,属于信任区域法,用于解决非线性最小二乘问题。LM算法结合了两种数值最小化算法:梯度下降法(Gradient Descent)和高斯-牛顿法(Gauss-Newton)。在梯度下降法中,通过更新最速下降方向的参数来减小平方误差之和。在高斯-牛顿法中,通过假设最小误差函数在参数中是局部二次的来减小平方误差的总和,并且找到这个二次方的最小值。LM方法当参数远离其最佳值时,像是梯度下降法,当参数接近其最优值时,更像高斯-牛顿法。

除了LM法,狗腿法(Dogleg)也是一种常用计算信任区域优化问题的方法。和LM不同的是,狗腿法计算两个关于步长的向量。狗腿法的优点是不需要从头计算,且信任区域半径小。狗腿法一般只能用直接法(分解因式法)求解线性方程组Ax=b。

CurveFitter最近更新

最近,WELSIM开发的免费软件CurveFitter进行了很多升级。其中包含:

1. 支持多核并行计算。对于大量测试数据,可以迅速得到曲线参数。

2. 支持用户输入求解器参数,如最大迭代次数,函数收敛公差值等。目前,前端支持的求解控制参数有:

最大迭代次数:求解器最大的迭代次数,超过此数值,判定不收敛。默认为100。

函数阈值:用于判断目标函数值变化的临界值。小于此值时判定收敛。默认为1e-9。

梯度阈值:用于判断最大范数的临界值。默认为1e-10。

参数阈值:用于判断步长参数的临界值。默认为1e-8。

3. 新增 1–6阶的Schulz-Flory函数。

Schulz-Flory函数已经广泛应用于描述各种复杂材料的聚合与结晶过程。的控制方程如下:

可以看出此方程具有高度非线性的特征,对于曲线拟合算法有着很高的要求。

4. 大量的软件优化与增强。

目前CurveFitter已经支持多种非线性曲线拟合,随着使用者的需求,将来会添加更多曲线仿真。目前可以从WELSIM的官网下载并免费使用。