摘要： 基于隐式微积分建模方法，提出分形维空间基本解的概念，从而定义分形维上的微积分算子，用以描述分形材料的各种力学行为.分形微积分算子极大地推广经典的连续介质力学微积分建模方法的使用范围，是分形导数概念的进一步发展.运用奇异边界法成功地数值模拟分形维拉普拉斯算子方程唯象描述的分形材料势问题.
关键词： 分形微积分算子； 隐式微积分方程建模； 唯象模型； 基本解； 分形导数
中图分类号： O39；O241.8文献标志码： A
Definition of calculus operator on fractal and its applications
CHEN Wen， WANG Fajie， YANG Xu
（State Key Laboratory of HydrologyWater Resources and Hydraulic Engineering； College of
Mechanics and Materials， Hohai University， Nanjing 210098， China）
Abstract： Based on the implicit calculus equation modeling approach， the concept of fundamental solution on fractal is defined and the fractal calculus operator is defined to describe the various mechanical behaviors of fractal materials. The fractal calculus operator significantly extends the application scope of the classical calculus modeling approach under the framework of continuum mechanics. It is also a stepforward advance of fractal derivative. By the singular boundary method， the numerical simulation is successfully performed on fractal Laplacian equation for phenomenological modeling potential problems in fractal media.
Key words： fractal calculus operator； implicit calculus equation modeling； phenomenological model； fundamental solution； fractal derivative
收稿日期： 2016[KG*9〗04[KG*9〗04修回日期： 2016[KG*9〗04[KG*9〗25
基金项目： 国家自然科学基金面上项目（11372097）；高等学校学科创新引智计划（B12032）
作者简介： 陈文（1967—），男，教授，博导，博士，研究方向为计算力学和软物质力学，（Email）chenwen@hhu.edu.cn0引言
现代计算科学主要是建立在微积分方程概念和建模方法基础上的，特别是连续介质力学问题的描述离不开微积分方程建模方法，但对于复杂分形结构材料和系统，经典的微积分方程方法面临着巨大的困难.一般的应用策略是直接拓广经典连续介质力学模型，运用非线性项描述分形介质中的复杂力学行为，因此模型中往往含有多个经验参数，且部分人为参数缺乏物理意义.
近年来，分数阶微积分方程建模方法引起广泛关注，成为描述复杂物理力学问题的一个有竞争力的建模方法[1].由于分数阶模型仍然是线性的，能够较好地刻画系统的历史和路径依赖特征，应用在某些问题上比非线性方法有一定的优越性.但是，分数阶微积分和分形几何的数学联系至今还不是很清楚，已有的研究多是定性讨论.[2]
分形几何方法在描述复杂系统的几何特征、统计行为、数据结果的幂律特征等方面取得很多有意义的成果[3]，但其对应的微积分建模方法至今没有完整地建立起来.这极大地限制分形方法在科学和工程问题中的应用.CHEN等[4]首次定义分形维α上分形导数的概念为dg（t）dtα=limt′→tg（t）-g（t′）tα-（t′）α（1）式中：g（t）为所考察的物理量；t为自变量；α为任意实数分形维.
此后，分形导数建模在反常扩散等问题上取得一些有意义的结果.[56]分形导数是局部导数，不同于全域定义的分数阶微积分，因而计算量和内存需求大大减少，但分形导数微分方程的应用目前很不成熟，在多维问题中的应用还很少.
针对多维分形空间问题，本文进一步发展分形导数的概念，定义分形维上的微积分算子.这项研究的关键创新点是拓广经典微分算子的基本解，提出分形维上微分算子基本解的概念.运用陈文等[7]提出的隐式微积分建模方法，根据分形维上的基本解“隐式”地定义分形微积分算子.分形微积分算子可以方便地数值计算和使用，但不一定具有显式表达式或其显示表达式难以得到.
本文以分形维上的拉普拉斯算子为例，详细介绍分形微积分算子的概念和具体应用，主要数值求解技术是奇异边界法[89].该方法以距离为基本变量，不依赖于问题的维数，本质上是无网格无数值积分方法，编程容易，能够计算高维复杂几何形状问题.
首先，引入分形维上微分算子基本解的概念.以拉普拉斯算子为例，比较分形和分数阶导数2种拉普拉斯算子基本解的区别与联系；然后，采用隐式微积分方程建模方法定义分形微积分算子，并给出分形维上拉普拉斯算子、亥姆霍兹算子、修正亥姆霍兹算子、扩散算子的定义；再次，以分形拉普拉斯算子方程为例，采用奇异边界法数值模拟二维和三维分形拉普拉斯算子方程，并对数值结果进行讨论和分析；最后，总结分形微积分算子的特点和建模方法的优势，以及若干有待深入研究解决的问题.
1分形维上微分算子的基本解
为不失一般性，以整数维上的整数阶拉普拉斯方程为例，其数学形式为Δu（x）=0， x∈Rn（2）式中：Δ为Rn上的拉普拉斯算子；n为整数阶空间维数（二维n=2得到的是平凡基本解）；u为待求势函数.相应的基本解[10]为u*n（r）=1（n-2）Sn（1）r2-n（3）式中：Sn（1）=2πn/2/Γ（n/2）；r=||xξ||为点x和ξ的欧氏距离.近年来引起广泛关注的分数阶拉普拉斯算子（-Δ）s/2能够表征物理力学系统的空间非局部性.采用隐式微积分建模方法，从其Riesz分数阶势出发，直接构造出分数阶拉普拉斯算子的基本解[11]为u*s（r）=1（d-s）Sd（1）rs-d（4）式中：s为分数阶数是0～2范围内的任意实数.经典整数阶拉普拉斯算子是一个特例，即s=2；这里s表征材料的非局部性，刻画幂律特征.
推广式（3）和（4）得到整数阶拉普拉斯算子在分形维d上的基本解为u*d（r）=1（d-2）Sd（1）r2-d（5）这里d可以是任意实数.
以三维空间问题为例，比较讨论分形维上的拉普拉斯基本解与分数阶拉普拉斯算子基本解的区别和联系.大部分三维空间问题的分形维在（2， 3]范围内，相应的分形维拉普拉斯算子的距离变量指数（2-d）在[-1， 0）范围内；分数阶拉普拉斯算子基本解的距离变量指数（s-3）在（-3，-1]范围内.由此可见，分形和分数阶拉普拉斯算子有各自不同的适用对象和范围，经典的整数阶拉普拉斯算子基本解1/r是两者的极端特例.
2分形微分算子的定义
根据隐式微积分建模方法，可以用基本解定义微分方程模型，不需要微分方程的显式表达式.基于此，本节运用分形维上的算子基本解，定义分形维上的4类典型微分算子方程.
拉普拉斯方程Δdu（x）=0， x∈Ω（6）亥姆霍兹方程（Δ+k2）du（x）=0， x∈Ω（7）修正亥姆霍兹方程（Δ-k2）du（x）=0， x∈Ω（8）扩散方程αΔdu（x）=u（x）t， x∈Ω，t≥0（9）式（6）～（9）中：下标d为分形维值为d的微分算子，以区别于经典的整数阶和分数阶微分算子.推广相应整数阶基本解[10]，分形维上亥姆霍兹、修正亥姆霍兹以及扩散算子的基本解定义为u*d（r）=12π-ik2πr（d/2）-1K（d/2）-1（-ikr）（10）
u*d（r）=12πk2πr（d/2）-1K（d/2）-1（kr）（11）
u*d（r）=H（t）（4παt）d/2e-r2/4αt（12）式中：K（d/2）-1为第二类修正贝塞尔函数；H（t）为赫维赛德阶跃函数；t=|t2-t1|为时刻到时刻的时间间隔；α为扩散系数；d为分形维数.分形维上的拉普拉斯算子基本解见式（5）.
3分形拉普拉斯势问题的数值模拟
拉普拉斯算子是最重要的椭圆型算子，在物理和力学中有着广泛而重要的应用.本节以拉普拉斯方程为例，数值考察分形维微分算子方程的行为特征.
奇异边界法[89]是一种边界型径向基函数配点法，以基本解作为插值基函数，能够无网格、无数值积分求解高维复杂几何域问题，不需要微分方程的具体表达式.本节基于分形维上拉普拉斯算子的基本解，采用奇异边界法求解分形维拉普拉斯控制方程和相应边界条件的稳态热传导问题.
首先，考虑一个二维正方形域分形介质中的稳态热传导问题，其边界条件见图1：左右边界绝热，热流量q=0，上边界温度u=0 °C，下边界温度u=10 °C.为考察温度变化与分形维数之间的关系，不同分形维数d情况下沿直线x=1.0温度值变化的数值计算结果见图2.由此可见：二维整数维情况下，温度的变化呈线性减小；相比较而言，分形维时温度变化呈指数趋势减小，且维数越小温度变化越剧烈.一般情况下，在不知道分形维上拉普拉斯方程的精确解时，可以通过指定与整数维方程相同的边界条件，考察分形维方程的数值解是否逼近于整数维方程的精确解.在本算例中，考察d趋于2时，方程的解是否逼近d=2整数阶拉普拉斯方程的解.从图2中可以看到，当维数d趋近于2时，分形维拉普拉斯方程的解确实单调趋近于整数维2的解.
Fig.4Variation of temperature u on line {（x， y， z） | x=1， y=1， 0≤z≤2}against fractal dimension d由图4可以看出：在完全相同边界条件下维数d趋近于3时，分形维拉普拉斯方程的解单调趋近于整数维为3的解；另外，三维整数维情形下温度的变化呈线性减小，而当材料具有分形特征时，温度变化在底部附近比整数维的变化缓慢，中间部分比整数维的变化剧烈，接近上顶部时温度的减小趋势又变缓.
4结束语
引入分形微积分算子是分形导数概念的进一步发展，可推广连续介质力学微积分建模方法的使用范围，克服现有分形方法局限于几何描述和数据拟合的瓶颈问题，拓广分形方法的应用范围和深度.
本文提出分形维上基本解的概念，基于隐式微积分建模方法，定义分形维上的微积分算子，微分控制方程表达式本身不再是必要的环节和对象.数学、力学建模和数值建模自然成为一体，极大地简化工程仿真的难度.
从数学上看，分形维上微分算子基本解表达式中的维数d甚至可以是复数或负数，但相关的物理力学意义并不清楚.此外，目前也有多种分形的测量方法和定义，具体到某个应用选择何种定义需要研究.本文提出的分形维微积分算子方法是唯象建模技术，还缺少扎实的数理基础；该方法的适用范围和有效性还有待在科学工程问题中充分验证.
说明：陈文提出本文的基本数学方法和整体研究思路；王发杰负责编程和数值结果整理；杨旭负责收集分形材料的有关数据.参考文献：
[1]陈文， 孙洪广， 李西成. 力学与工程问题的分数阶导数建模[M]. 北京： 科学出版社， 2010.
[2]CARPINTERI A， MAINARDI F. Fractals and fractional calculus in continuum mechanics[M]. SpringerVerlag， 1997.
[3]张济忠. 分形[M]. 北京： 清华大学出版社， 1995.
[4]CHEN W. Timespace fabric underlying anomalous diffusion[J]. Chaos， Solitons & Fractals， 2006， 28（4）： 923929. DOI： 10.1016/j.chaos.2005.08.199.
[5]SUN H G， MEERSCHAERT M M， ZHANG Y， et al. A fractal Richards′equation to capture the nonBoltzmann scaling of water transport in unsaturated media[J]. Advances in Water Resources， 2013， 52： 292295. DOI： 10.1016/j.advwatres.2012.11.005.
[6]CHEN W， ZHANG X D， KOROAK D. Investigation on fractional and fractal derivative relaxationoscillation models[J]. International Journal of Nonlinear Sciences and Numerical Simulation， 2010， 11（1）： 39. DOI： 10.1515/ijnsns.2010.11.1.3.
[7]陈文. 复杂科学与工程问题仿真的隐式微积分建模[J]. 计算机辅助工程， 2014， 23（5）： 16. DOI： 10.13340/j.cae.2014.05.001.
CHEN W. Implicit calculus modeling for simulation of complex scientific and engineering problems[J]. Computer Aided Engineering， 2014， 23（5）： 16. DOI： 10.13340/j.cae.2014.05.001.
[8]陈文. 奇异边界法： 一个新的、简单、无网格、边界配点数值方法[J]. 固体力学学报， 2009， 30（6）： 592599.
CHEN W. Singular boundary method： a novel， simple， meshfree boundary collocation numerical method[J]. Chinese Journal of Solid Mechanics， 2009， 30（6）： 592599.
[9]WEI X， CHEN W， SUN L L， et al. A simple accurate formula evaluating origin intensity factor in singular boundary method for twodimensional potential problems with Dirichlet boundary[J]. Engineering Analysis with Boundary Elements， 2015（58）： 151165. DOI： 10.1016/j.enganabound.2015.04.010.