Moving least square(MLS),移动最小二乘法,被广泛地用在了曲面拟合和函数逼近问题上。它通过在采样点周围拟合出一个具有local support性质的函数,来构筑连续光滑的曲面。最近在做的布料仿真就用到了MLS,中间遇到了个坑,花费了几天时间才解决,特在此记录一下。本文侧重讨论于MLS的实现,关于原理部分可参见以下文献:
1. Garimella, R.V., 2017. A simple introduction to moving least squares and local regression estimation (No. LA-UR-17-24975). Los Alamos National Lab.(LANL), Los Alamos, NM (United States). 详细介绍了如何一步步地从global least square逐步推导出moving least square.,但本文不涉及MLS导数的计算。
2.Yaw, L., 2009. Introduction to moving least squares (MLS) shape functions. Walla Walla University, College Place. 给出了MLS approximation及其导数的具体公式。
3.Kushawaha, J.S. and Misra, B.K., 2013. MESHFREE SHAPE FUNCTION AND DERIVATIVES IN ONE DIMENSION WITH ARBITRARY NODAL DISTRIBUTION: MATLAB PROGRAMMING. Suranaree Journal of Science & Technology, 20(4). 讨论了MLS的一些性质并用画出了图表。
4.Li, W., Nguyen-Thanh, N., Huang, J. and Zhou, K., 2020. Adaptive analysis of crack propagation in thin-shell structures via an isogeometric-meshfree moving least-squares approach. Computer Methods in Applied Mechanics and Engineering, 358, p.112613. Chapter 2.1 给出了二维MLS的公式,方便实现。
1. c++实现
![](https://i-blog.csdnimg.cn/blog_migrate/506f9419e577ef4605ac7c47ee6a7404.png)
Fig.1 是二维MLS系数计算时用到的变量,包括A,B和P。其中, P是采样点的多项式向量, 而A和B则是所有控制点多项式向量组成的矩阵。下面给出代码求解这些变量。
2. 代码
#include <Eigen/Core>
#include <Eigen/Eigenvalues>
// cubic weight function
double cubicWeight(double sA, double s, double radius)
{
double result = 0.0;
double r = std::abs(sA - s) / radius;
if (r <= 1.0)
{
result = 1.0 - 6.0 * r * r + 8.0 * r * r * r - 3.0 * r * r * r * r;
}
else
{
result = 0.0;
}
return result;
}
// first order derivative
double cubicWeight_d1(double sA, double s, double radius)
{
double result = 0.0;
double r = std::abs(sA - s) / radius;
if (r <= 1.0)
{
result =