matlab ceil函数_径向基函数生成的有限差分法(RBF-FD)原理及程序实现(2)——Halton节点分布

参考文献

Fornberg, B. and Flyer, N., 2015.A primer on radial basis functions with applications to the geosciences. Society for Industrial and Applied Mathematics.

Fornberg, B., Lehto, E. and Powell, C., 2013. Stable calculation of Gaussian-based RBF-FD stencils. Computers & Mathematics with Applications, 65(4), pp.627-637.

Larsson, E., Lehto, E., Heryudono, A. and Fornberg, B., 2013. Stable computation of differentiation matrices and scattered node stencils based on Gaussian radial basis functions. SIAM Journal on Scientific Computing, 35(4), pp.A2096-A2119.

Fornberg, B. and Flyer, N., 2015. Solving PDEs with radial basis functions. Acta Numerica, 24, p.215.

Flyer, N., Fornberg, B., Bayona, V. and Barnett, G.A., 2016. On the role of polynomials in RBF-FD approximations: I. Interpolation and accuracy.Journal of Computational Physics,321, pp.21-38.

Flyer, N., Wright, G.B. and Fornberg, B., 2014. Radial basis function-generated finite differences: A mesh-free method for computational geosciences. Handbook of geomathematics, pp.1-30.

地球物理局 地震波动力学实验室 无网格组

声明:
# 系列文章优先满足个人研究需求
# 欢迎批评指正,禁止转载

目 录

石中居士:径向基函数生成的有限差分法(RBF-FD)原理及程序实现——目录​zhuanlan.zhihu.com
8f5c65b7215a7f20e543d80084c6c0e3.png

基本概念

在基于网格的常规有限差分近似的情况中,可以在所有节点上(重新)使用相同的模板(stencil)形状和权重,可能在边界处进行较小的修改。在当前的分散节点的情况中,每一个模版将变得不同。当整个域中有

个节点时,一种常见的策略是将RBF-FD模版置于每个节点中央,并使其扩展到其
个最近的邻点,使每个模版中总共有
个节点。直到最近,数值稳定性问题已将
通常限制为相对较小的值,例如在二维中,
之间。“高粘度”方法允许在某些情况下,将点数增加至少一个数量级,从而实现高精度近似。从概念上讲,RBF-FD可以看作是重叠域分解的一种极端情况,每个节点周围都有一个独立的域。

一些节点分布策略

  • Halton节点

在创建用于求解PDE的节点分布时,完全随机分布(例如由MATLAB中的“rand”生成)和高度规则的分布(例如笛卡尔网格)都具有缺点或局限性。在需要准均匀分散的节点的情况下,一个易于实现的选择是使用Halton型节点分布。为了获得最佳的精度,它们通常很不规则,但是它们比完全随机分布更方便并且更“健壮(robust)”,因为可以避免出现节点偶然的巧合。例如,图1对比了单位正方形中

的Halton节点与相同数量的完全随机节点。

7fa193821e0fc87a2eb47b3885078e42.png
图1 N=2000的Halton节点(左)和随机节点(右)的对比(Fornberg and Flyer, 2015)

Halton序列最早是在1960年引入的[1],从那时起,它已成为各种应用(如蒙特卡洛模拟和高维数值积分)中的准随机数的常规来源。它们的数学论证已被反复描述[2][3]。Fornberg and Flyer(2015)给出了一个特别短而且快速的MATLAB代码。在高达6维的情况下,此代码产生的Halton序列与MATLAB Central中的haltonseq算法和MATLAB统计工具箱中的haltonset算法相同(但是缺少后面这些代码的某些额外选项,例如生成指定的子序列)。

halton.m

%   Fornberg, B. and Flyer, N., 2015.
%   A primer on radial basis functions with applications to the geosciences. 
%   Society for Industrial and Applied Mathematics.
%   P.191

%   地球物理局
%   地震波动力学实验室
%   无网格组

%   仅供学习与参考,请勿用于商业目的

function H = halton(numpts,ndims)
    %   该子例程创建与haltonseq和haltonset相同的序列
    %   numpts (标量) :在Halton序列中生成的点数
    %   ndims (标量) :维数,应小于等于6
    if ndims > 6 
        error('ndims > 6');     %   当维数大于6时显示错误
    end
    p = [2 3 5 7 11 13]; 
    H = zeros(numpts,ndims);    %   初始化序列
    for k = 1:ndims
        N = p(k); 
        v1 = 0; 
        v2 = 0:N-1; 
        lv1 = 1;    %   v1的长度
        while lv1 <= numpts
           v2 = v2(1:max(2,min(N,ceil((numpts+1)/lv1))))/N;
           [x1,x2] = meshgrid(v2,v1); 
           v1 = x1+x2; 
           v1 = v1(:);
           lv1 = length(v1);
        end
        H(:,k) = v1(2:numpts+1);
    end
end

test.m

%   地球物理局
%   地震波动力学实验室
%   无网格组

%   仅供学习与参考,请勿用于商业目的

clear;
clc;

H = halton(2000,2);

figure(1)
scatter(H(:,1),H(:,2),1,[0,0,0]);
axis equal;
axis([0 1 0 1]);
set(gca,'Fontname','Times New Roman','FontSize',10);
xlabel('x','Fontname','Times New Roman','FontSize',10);
ylabel('y','Fontname','Times New Roman','FontSize',10);

图像

634a49bee3ac88870e52e1160c3d85f9.png
N=2000的Halton节点

商业推广


参考

  1. ^J. H. Halton, On the efficiency of certain quasi-random sequences of points in evaluating multi- dimensional integrals, N u m e r . M a t h . 2 (1960), 84–90. (Cited on p. 191)
  2. ^G. E. Fasshauer, Meshfree Approximation Methods with MATLAB, Interdisciplinary Mathematical Sci- ences 6, World ScientificPublishers, Singapore, 2007. (Cited on pp. 39, 42, 44, 45, 47, 56, 58, 69, 88, 191)
  3. ^L. Kocis and W. J. Whiten, Computational investigations of low-discrepancy sequences, A C M Trans. Math. Softw. 23 (1997), 266–294. (Cited on p. 191)
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值