2015, Vol. 47
随着空间技术的进步，人类对太空的探索程度不断深入，小行星探测已经成为本世纪深空探测的主要发展方向之一. 小行星上有丰富的资源，保存着太阳系形成初期的原始成分，是研究太阳系起源和演化历史的``活化石''；同时，部分小行星对地球具有潜在撞击威胁，通过深空探测任务了解其组织结构对研究开展规避任务具有重要意义.
目前，各航天大国已先后提出并实施了基于不同对象的深空探测计划，旨在占领空间高新技术战略制高点，从而在未来的空间 资源争夺和综合国力竞争中保持优势. 中国也提出了以月球探测为突破口、稳步发展深空探测技术的策略，相关单位也积极开展小行星探测理论研究并启动了总体方案论证工作. 2012年12月13日，中国嫦娥二号卫星飞抵距离地球约700万千米的深空，并成功对4179号小行星Toutatis进行近距飞越探测，在国 际上首次实现了对小行星Toutatis的近距探测.
随着深空探测技术的发展，国际小行星探测任务从初期的近距飞越探测(1991年美国Galileo号探测器)发展到近期的小行星环绕 着陆探测(2000年美国NEAR号探测器)，再到小行星采样返回探测(2005年日本隼鸟号(Hayabusa)探测器). 其中，对小行星开展采样返回探测将是未来几十年内航天界的研究热点，其涵盖了深空探测的几大关键技术[1,2,3,4]：自主导航技术、轨道设计与控制、自主采样技术. 而对小行星不规则形状、结构及引力场的建模则是实现这3个技术的根基.
现有文献中，对小行星不规则引力场的建模主要采用球谐系数展开和多面体的方法，其在小行星附近轨道动力学以及表面运动的研究中已 取得大量成功应用[5,6,7,8,9,10,11,12,13,14,15,16,17]. 但这两种方法不能直接应用于小行星表面物质在探测活动中的动力特性的研究. 根据已有的小行星探测数据，大多数近地小行星是由几毫米至几十米的碎石颗粒在引力的作用下聚合构成的[18,19](嫦娥 二号飞越探测的4179号小行星Toutatis就是这种结构[20]). 这种结构的小行星被形象地称为``碎石堆''小行星，其内部疏松多孔、易脆，同时由于小行星的引力场较弱，碎石堆结构较容易在强 力的作用下破碎瓦解[21]. 例如日本隼鸟号(Hayabusa)探测器曾经造访的25143号小行星Itokawa也是一颗典型的碎石堆小行星，其质量为3.5×10kg，表面逃逸速度约为0.17 m/s，这意味着只要有颗粒运动速度高于0.17 m/s就可以飞离该小行星. 因此，在对碎石堆小行星进行探测时，需要非常小心以确保航天器的安全和任务的成功. 假若锚定或开采时产生的力或振动导致小行星表面物质飞溅脱离，就很有可能对探测器造成撞击威胁；而探测器在没有锚定的情 况下也非常容易飞离小行星. 这就需要对碎石堆小行星力学特性进行建模和仿真分析，为安全可靠地完成探测任务(如锚定、着陆、挖掘、采样收集、起飞等)提供 理论依据.
理解颗粒物质行为特征对于开展小行星探测活动至关重要，将散体动力学中的数值模拟方法应用到其中也可在一定程度上揭示小行星的形成机制(碰撞、聚合)和演化机理，从而促进太阳系起源和演化的研究. 本文对散体动力学的研究和发展情况进行了调研，并分析了几种具有代表性的模型和方法；再根据国内外的相关研究，对目前碎石堆小行星动力学模型与数值仿真方法的研究进行了分析；最后结合小行星探测的应用背景，总结了小行星散体动力学模型中需解决的关键问题.
1 散体动力学的研究和发展现状散体在我们的生活中普遍存在，如自然界中的砂砾、谷粒以及工业上的固体化工原料、粉末等都属于散体模型. 从20世纪80年代以来，散体的力学性能研究日益受到国内外相关领域研究机构的重视，其研究成果也在Science和Nature等具有重要影响力的期刊上不断涌现，研究背景涉及机械设计、物料输送、化工合成、登月工程等多个方向.
散体的性质与现在研究得较为广泛的固体和流体有很大不同. 原静止的散体在外力作用下没有发生流动时表现出的宏观性质与固体类似；开始流动后其表现出的宏观性质又和流体类似. 但散体在静止时的力学性质又与一般固体不同，流动时的性质也与一般流体不同，现今还没有找到一种很好的理论来描述散体的这种性质.
1.1 散体动力学的实验研究由于散体动力学基础理论不够完善，目前主要采用实验和数值模拟的方法对散体系统进行研究. 国内外相关领域的研究人员开展了大量实验来研究散体的动力学过程，对其在外力作用下的自组织、分层结构和聚集现象有了一定的认识. Cates等[22]通过实验研究发现当散体材料负载变动时，即使变动很微小，粒子群的构型都会重新排列直到能够维持力和力矩的平衡；而在冲击载荷足够大的情况下，散体材料将呈现短暂的流体化状态，并瞬间完成大量的碰撞过程[23,24]. 利用X射线断层摄影术，Mètayer 等[25]研究了不同密度的粒子在剪切载荷下的响应，研究发现疏散排列的散体材料在剪切力的作用下粒子间相互填补空隙使得整体体积压缩，而紧密排列的散体材料在剪切力的作用下则会膨胀. Murdoch 等[26]研究了不同程度的重力驱动下的散体对流行为，发现在微重力情况下散体粒子几乎不发生二次流现象，随着重力的增大二次流现象会增强，而主流不受重力影响，由此认为重力的大小会改变粒子间和粒子与接触面间的摩擦作用力从而影响二次流.
1.2 散体动力学的数值模拟由于实验方法在参数选取上具有不可避免的局限性，同时也不能完全解释散体动力学机理，研究一般的散体动力学过程还须用数值模拟的方法. 散体数值模拟方法中运用的最为广泛的就是离散单元法(distinct element method，DEM). 该方法以牛顿第二定律为基础，最早由Cundall[27]提出和发展. 所研究的对象是大量具有特定形状的离散颗粒(通常是球体)，每个颗粒具有一定的自由度，颗粒间通过接触力和摩擦力相互作用，同时也可施加外力、外力场或位移的边界条件进行限制.
接触模型是DEM的核心，对于干燥球形颗粒目前主要有三种接触模型：硬球模型，软球模型和基于接触力学的模型[28,29]. 其中，硬球模型[30]是将研究的颗粒视为刚体，忽略刚体碰撞的过程，直接根据恢复系数和碰撞刚体的运动状态来求解碰撞后的结果，适合于处理快速稀疏颗粒流中的二体碰撞. 由于不用处理碰撞过程中接触力对时间的积分，硬球模型能够快速处理碰撞. 软球模型[27]采用弹簧阻尼、Hertz接触等力的模型根据粒子交叠的程度刻画相互间的接触力和摩擦力，因其能够处理多体碰撞且可模拟静态颗粒流，接近客观的颗粒碰撞过程，因此得到了大量实际应用和研究[31,32,33,34]，众多学者也针对不同的问题背景提出了各种法向力和切向力模型[35,36,37,38,39]，文献[40]中对这些模型进行了详细地比较分析. 而基于接触力学的模型[41]则是在软球模型的基础上考虑应力加载历史，利用颗粒接触理论和颗粒材料的实际物理参数推导接触力增量和位移增量. 由于基于接触力学的模型不对接触过程做任何简化，计算效率不高，在较大规模的数值研究中运用不多. 实际应用中应根据不同的问题特点选择模型.
虽然球形颗粒模型能够高效地进行大规模模拟(便于检测和处理碰撞，数量可达到105量级)，但球形假设离真实颗粒形状差别较大，对结果有显著影响. 因此近几年来，一些研究人员还发展了非球集合的散体接触模型， 可对椭球[42,43]、圆盘和球的混合模型[44,45]进行模拟. 姜世平等[46]用射线交叉法判断颗粒间的接触实现了对少量凸多面体组成的散体系统接触碰撞动力学的数值仿真. 但在对任意形状的散体粒子建立接触模型时，无论是对接触的判断还是对大量粒子的模拟都有困难. 目前大多数关于颗粒物质的数值模拟研究还是采用圆形或球形模型进行模拟.
2 碎石堆小行星散体动力学建模方法小行星碎石堆模型作为散体动力学模型之一，其受力特征与宏观粒子类似，主要经历以下经典力学过程：受力平衡，碰撞，摩擦，非弹性接触等. 同时，对于碎石堆小行星而言，其粒子间还存在万有引力作用，因此我们在研究其动力学行为时需要综合考虑散体和引力N体两方面的问题.
2.1 引力N体问题引力N体问题是把N(>2)个天体视为质点，研究它们在相互引力作用下各自运动情况. 由于对N>2的N体系统已无解析解，只能 用数值方法在计算机上模拟天体和天体系统的形成和演化. 目前这方面的研究较为成熟，形成了一些常用的算法，例如：低阶时间中心的蛙跳法(leapfrog integrator)[47]，高阶均差单时步多项式法(force polynomials integrator)[48]，Ahmad-Cohen邻域法[49,50]，高阶RK法(Runge-Kutta algorithm)及其衍生方法(如RKF，RKN等)[51]，辛算法(Symplectic algorithm)[52,53]和厄米算法(Hermite algorithm)[54,55]等.
根据现有研究比较几种数值模拟方法，得知：蛙跳法简单易行且能维持能量守恒，但其在处理短距交会问题时会在短时间内产生较大积分误差，需要很小的时步才能够减小误差；多项式法在时步选取正当时，能在较长时间保持能量守恒，但多项式法所需存储空间较大，耗费机时较长，通常仅用于N<100的情况；在多项式法基础上发展的AC邻域法，采用了将近距引力作用与远距引力作用分开计算的方案，效率相比于多项式法有很大提高，也不损失精度，可处理N~1 000的情况；RK高阶算法精度很高，但对于多粒子N体问题，长期积分不能保持能量守恒，使轨道误差迅速增长；辛算法能够严格保持Hamilton系统的几何结构，没有人工耗散，保持能量守恒，但不能高效地模拟非开普勒轨道运动，而且难于采用变步长进行积分，不适用于存在碰撞的系统；Hermite算法效率较高，构造简单，但精度比其他算法低. 这些算法各有所长，但不能直接运用到碎石堆小行星的模拟中，需要结合散体动力学的数值方法进行改进.
2.2 散体-引力N体耦合模型近几年来，国外的一些研究人员将引力N体问题的数值积分方法与DEM结合，发展出了适用于研究碎石堆小行星演化的数值模拟方法.
Richardson 等[56]运用硬球模型的离散元数值方法，并采用四阶预估校正方法进行引力积分，研究了碎石堆小行星在木星潮汐力作用下的 演化情况，所模拟的碎石堆模型由247个球体组成. 在此基础上，Leinhardt等[57]采用并行蛙跳法N-body程序代码来搜寻碰撞和计算球体的轨迹，进一步提高了算法的效率， 模拟了两个碎石堆小行星低速撞击的过程，其模拟个数达到了103量级. 此后，采用此方法的相关研究大量涌现. Richardson[58]等采用该方法模拟了碎石堆小行星旋转角速度与形状的关系. Walsh等[59]采用相同的模型，仿真了一颗具有一定形状的碎石堆小行星在Yarkovsky效应的作用下高速旋转时在赤道处飞射 出物质的现象，并得出结论：若此飞射物质因碰撞而产生能量耗散，则会被吸进它的卫星中，从而演化成双星. Leinhardt 和Stewart[60]将光滑粒子流体动力学方法(smoothed particle hydrodynamics， SPH)与该模型结合，研究比较了碎石堆小行星在不同撞击速度、撞击角度及撞击者尺寸的情况下的撞击结果.
但由于硬球方法本质上不适合于模拟稠密颗粒流的行为，将其应用在模拟聚集的颗粒物质(即碎石堆小行星)时需要进行一些脱离实际的假设以避免发生多球碰撞的情况[61]. 这种简化假设忽略了颗粒物质的表面变形，不适合于处理具有复杂接触现象的颗粒流. 因此，虽然硬球模型具有较高的积分效率(时步为颗粒自由程相应时间)，但在碎石堆的模拟中无论是仿真程度上还是算法的稳定性上都有一定缺陷. 软球模型则可弥补这一缺陷. 软球模型虽然在效率上远小于硬球模型(时步为变形时间的1/50)，但可以有效地解决多球碰撞和变形的情况，近几年来开始在小行星领域得到了应用[62,63]. Schwartz等[64]在软球模型接触力的基础上加入了黏合力，并用进行了实验验证，实验结果与数值模拟结果较为吻合. Yu等[65]用软球模型模拟了小行星Apophis接近地球时表面风化层在潮汐作用下的运动情况，发现近地交会将引发其表面``风化层''的崩塌，揭示了小行星表面地貌重塑的一种重要机制.
软球模型给小行星结构的模拟和风化层的研究注入了新的力量，但球体的模型仍是不足以描述太空中颗粒物质的真实行为. 同时，计算效率始终是颗粒物质数值模拟中的瓶颈(在引入引力场的积分后更为复杂)，如何提高模拟的效率和颗粒的数量也是未来亟待解决的问题.
3 小行星探测中散体动力学建模的关键问题现有关于小行星散体模型的文献都主要关注小行星自身尺度的动力学特点，与小行星探测活动中的着陆、采样、起飞等行为的关系不大. 而现有的关于探测活动的研究中[66- 68]一般只考虑小行星表面颗粒物质在微重力情况下的运动，很少把研究置于整个小行星的引力场环境中，也没有考虑探测活动对于小行星整体的影响，这对于探测碎石堆小行星却是必不可少的[69]. 本文将对小行星探测活动中的散体建模的应用和关键问题进行论述分析.
3.1 小行星散体模型应用 3.1.1 小行星着陆与返回目前，人类关于小行星地质结构、表面强度、质量分布等性质的知识仍然非常匮乏，在实施小行星探测活动之前很难完全掌握小行星及着陆点的相关信息，需要实时施加控制以安全着陆. 并且由于引力场非常弱，在小行星表面着陆器很容易反弹飞离目标，还需要设计锚定系统. 着陆和锚定过程在一定程度上会造成小行星表面风化层的演变，甚至引发内部的振动，给探测任务带来风险. 小行星散体模型可对着陆装置和着陆性能进行仿真模拟，为探测器的安全着陆返回的特征参数提供理论基础.
3.1.2 小行星采样小行星采样过程是衡量采样返回任务成败的重要标志，有极大的不确定性. 地面模拟试验不能反应出空间的力学特性，也无法 进行整星级的试验. 应用小行星散体模型，可对小行星采样过程进行全系统仿真，从而对可能出现的异常情况进行预估，指导采样区域和采样方案的确定.
3.1.3 小行星自身演化影响探测活动往往会对小行星自身的性质带来一定影响，甚至可以改变小行星的演化情况. 例如，小行星防御理论中所提出的利用Yarkovsky 效应改变小行星运行轨道的方案，就是通过改变小行星表面性质来实施的. 探测活动中可以进行这样的试验，在一定程度上改变小行星的运行情况，并与小行星散体模型进行比较，从而促进对小行星物理性质的认识，进而能够改善模型的参数，为今后的小行星防御任务打下基础.
3.2 需解决的关键问题 3.2.1 散体颗粒接触模型碎石堆小行星由大量的散体颗粒构成，颗粒间通过接触力、摩擦力以及引力相互作用，研究其力学行为首先需要对接触模型进行分析. 一般情况下构成碎石堆小行星的碎片(尤其是大块的碎片)不可能是标准的球形，与球体的受力差别较大，对最终结果有一定影响. 为了更真实地反映小行星碎石堆的运动状态，模拟中需要考虑非球形结构的颗粒(如多面体)的接触模型和数值模拟方法. 但实际操作中对多面体散体元的计算和接触点搜索较为复杂，对数量巨大的颗粒群进行实用模拟较为困难，需要对算法实现进行深入研究.
3.2.2 不规则碎片引力场的精确模型小行星及其碎片通常具有不规则的形状，产生的弱引力场变化复杂，必须精确考察这些特点导致的复杂动力学环境对碎石堆构型及其运动的影响，及对探测器的影响. 因而，需要建立高精度的描述不规则碎片的引力场模型，再研究多个碎片之间引力、接触力与碰撞力相互作用的耦合影响，从而分析小行星姿态和构型在外力作用下的变化. 目前，一些学者已经建立了双多面体引力场的数值积分方法[70]，但其在处理效率和接触作用上仍有一定困难，需要修改才能应用于碎石堆模型.
3.2.3 颗粒材料性质的获取小行星散体接触模型中需要的颗粒材料的相关参数(如杨氏模量、泊松比、剪切模量等)往往都是难以获取的，为保证模拟的仿真程度，需要与实验进行对比确定. 实验中可利用与小行星物理观测数据相似的物体(如花岗岩、冰块等)，通过静力学实验、撞击试验、滑坡实验、光弹性实验等对相关量进行测量，并可将其颗粒流动力学行为与模拟结果进行分析比对，优化模型及参数.
4 总结与展望小行星散体动力学模型融合多学科研究方法，模型假设贴近实际，对开展小行星探测任务具有重要指导作用. 本文总结了散体动力学的发展情况，并对其常用的各种模型进行了分析；继而，综述了散体模型在碎石堆小行星建模方面的国内外研究现状以及发展趋势；最后分析了其在小行星探测中的应用优势，指出目前仍需解决的关键问题.
对小行星开展深空探测将是未来几十年内航天界的研究热点. 目前，各航天大国已分别提出相应的小行星采样返回任务的计划，如：JAXA的Hayabusa2号任务(计划2014年发射)；NASA的OSIRIS-Rex探测器(计划2016年发射)； ESA的MarcoPolo-R任务(计划2022年发射)等. 中国作为航天大国，将来一定会开展小行星探测、太阳系科学探索等项目. 建立小行星散体动力学模型，是检验小行星探测任务可行性的基础，具有广泛的应用前景.
| [1] | 李俊峰,崔文,宝音贺西. 深空探测自主导航技术综述. 力学与实践,2012,34(2): 1-9 (Li Junfeng, Cui Wei, Baoyin Hexi. A survey of autonomous navigation for deep space exploration.Mechanics in Engineering, 2012, 34(2): 1-9 (in Chinese)) |
| [2] | 李俊峰,宝音贺西. 深空探测中的动力学与控制. 力学与实践,2007,29(4): 1-8 (Li Junfeng, Baoyin Hexi. The dynamics and control for deep space exploration.Mechanics in Engineering, 2007, 29(4): 1-8 (in Chinese)) |
| [3] | Baoyin HX, Li JF. A survey on orbital dynamics and navigation of asteroid missions.Acta Mechanica Sinica, 2014, 30(3): 282-293 |
| [4] | Barucci MA, Cheng AF, Michel P, et al. MarcoPolo-R near earth asteroid sample return mission. Experimental Astronomy, 2012, 33(2-3): 645-684 |
| [5] | Scheeres DJ, Ostro SJ, Hudson RS, et al. Orbits close to asteroid 4769 Castalia.Icarus, 1996, 121(1): 67-87 |
| [6] | Scheeres DJ, Ostro SJ, Hudson RS, et al. Dynamics of orbits close to asteroid 4179 Toutatis.Icarus, 1998, 132(1): 53-79 |
| [7] | Hu WD, Scheeres DJ. Numerical determination of stability regions for orbital motion in uniformly rotating second degree and order gravity fields.Planetary and Space Science, 2004, 52(8): 685-692 |
| [8] | Hu WD, Scheeres DJ. Averaging analyses for spacecraft orbital motions around asteroids.Acta Mechanica Sinica, 2014, 30(3): 294-300 |
| [9] | Hu WD, Scheeres DJ. Spacecraft motion about slowly rotating asteroids.Journal of Guidance, Control, and Dynamics, 2002, 25(4): 765-775 |
| [10] | Tricarico P, Sykes MV. The dynamical environment of Dawn at Vesta.Planetary and Space Science, 2010, 58(12): 1516-1525 |
| [11] | Liu XD, Baoyin HX, Ma XR. Equilibria, periodic orbits around equilibria, and heteroclinic connections in the gravity field of a rotating homogeneous cube.Astrophysics and Space Science, 2011, 333(2): 409-418 |
| [12] | Yu Y, Baoyin HX. Orbital dynamics in the vicinity of asteroid 216 Kleopatra.The Astronomical Journal, 2012, 143(3): 62 |
| [13] | Yu Y, Baoyin HX. Generating families of 3D periodic orbits about asteroids.Monthly Notices of the Royal Astronomical Society, 2012, 427(1): 872-881 |
| [14] | Yu Y, Baoyin HX. Routing the asteroid surface vehicle with detailed mechanics.Acta Mechanica Sinica, 2014, 30(3): 301-309 |
| [15] | Li XY, Qiao D. Modeling method and preliminary model of asteroid Toutatis from Chang'e-2.Acta Mechanica Sinica, 2014, 30(3): 310-315 |
| [16] | Wang XY, Gong SP, Li JF. A method for classifying orbits near asteroids.Acta Mechanica Sinica, 2014, 30(3): 316-325 |
| [17] | Jiang Y, Baoyin HX, Li JF, et al. Orbits and manifolds near the equilibrium points around a rotating asteroid.Astrophysics and Space Science, 2014, 349(1): 83-106 |
| [18] | Fujiwara A, Kawaguchi J, Yeomans DK, et al. The rubble-pile asteroid Itokawa as observed by Hayabusa.Science, 2006, 312(5778): 1330-1334 |
| [19] | Scheeres DJ, Hartzell CM, Sánchez P, et al. Scaling forces to asteroid surfaces: The role of cohesion.Icarus, 2010, 210(2): 968-984 |
| [20] | Huang JC, Ji JH, Ye PJ, et al. The ginger-shaped asteroid 4179 Toutatis: new observations from a successful flyby of Chang'e-2. Scientific Reports, 2013, 3 |
| [21] | Britt DT, Consolmagno SJ. Modeling the structure of high porosity asteroids.Icarus, 2001, 152(1): 134-139 |
| [22] | Cates ME, Wittmer JP, Bouchaud JP, et al. Jamming, force chains, and fragile matter.Physical Review Letters, 1998, 81(9): 1841 |
| [23] | Daniels KE, Coppock JE, Behringer RP. Dynamics of meteor impacts.Chaos, 2004, 14(4): S4 |
| [24] | Clark AH, Kondic L, Behringer RP. Particle scale dynamics in granular impact.Physical Review Letters, 2012, 109(23): 238302 |
| [25] | Métayer JF, Suntrup DJ III, Radin C, et al. Shearing of frictional sphere packings.Europhysics Letters, 2011, 93(6): 64003 |
| [26] | Murdoch N, Rozitis B, Nordstrom K, et al. Granular convection in microgravity.Physical Review Letters, 2013, 110(1): 018307 |
| [27] | Cundall PA, Strack ODL. A discrete numerical model for granular assemblies.Geotechnique, 1979, 29(1): 47-65 |
| [28] | Haile JM. Molecular Dynamics Simulation: Elementary Methods 1st. New York: John Wiley & Sons, Inc., 1992 |
| [29] | 孙其诚,王光谦. 颗粒流动力学及其离散模型评述. 力学进展,2008,38(1): 87-100 (Sun Qicheng, Wang Guangqian. A review on the granular dynamics and discrete model.Advances in Mechanics, 2008, 38(1): 87-100 (in Chinese)) |
| [30] | Hong DC, McLennan JA. Molecular dynamics simulations of hard sphere granular particles.Physica A: Statistical Mechanics and Its Applications, 1992, 187(1-2): 159-171 |
| [31] | Biswas P, Sanchez P, Swift MR, et al. Numerical simulations of air-driven granular separation. Physical Review E, 2003, 68(5): 050301 |
| [32] | Campbell CS, Cleary PW, Hopkins M. Large-scale landslide simulations: Global deformation, velocities and basal friction. Journal of Geophysical Research: Solid Earth, 1995, 100(B5): 8267-8283 |
| [33] | Moakher M, Shinbrot T, Muzzio FJ. Experimentally validated computations of flow, mixing and segregation of non-cohesive grains in 3D tumbling blenders. Powder Technology, 2000, 109(1): 58-71 |
| [34] | Richard P, Valance A, Métayer JF, et al. Rheology of confined granular flows: Scale invariance, glass transition, and friction weakening.Physical Review Letters, 2008, 101(24): 248002 |
| [35] | Brilliantov NV, Spahn F, Hertzsch JM, et al. Model for collisions in granular gases.Physical Review E, 1996, 53(5): 5382 |
| [36] | Walton OR, Braun RL. Viscosity, granular-temperature, and stress calculations for shearing assemblies of inelastic, frictional disks.Journal of Rheology, 1986, 30(5): 949-980 |
| [37] | Haff PK, Werner BT. Computer simulation of the mechanical sorting of grains.Powder Technology, 1986, 48(3): 239-245 |
| [38] | Thornton C. Coefficient of restitution for collinear collisions of elastic-perfectly plastic spheres.Journal of Applied Mechanics, 1997, 64(2): 383-386 |
| [39] | Shäfer J, Dippel S, Wolf DE. Force schemes in simulations of granular materials.Journal de Physique I, 1996, 6(1): 5-20 |
| [40] | Stevens AB, Hrenya CM. Comparison of soft-sphere models to measurements of collision properties during normal impacts. Powder Technology, 2005, 154(2): 99-109 |
| [41] | 孙其诚,王光谦. 颗粒物质力学导论. 北京:科学出版社, 2009 (Sun Qicheng, Wang Guangqian. An Introduction to Granular Material Mechanics. Beijing: Science press, 2009 (in Chinese)) |
| [42] | Mailman M, Schreck CF, O'Hern CS, et al. Jamming in systems composed of frictionless ellipse-shaped particles.Physical Review Letters, 2009, 102(25): 255501 |
| [43] | Zeravcic Z, Xu N, Liu AJ, et al. Excitations of ellipsoid packings near jamming.Europhysics Letters, 2009, 87(2): 26001 |
| [44] | Phillips CL, Anderson JA, Huber G, et al. Optimal filling of shapes.Physical Review Letters, 2012, 108(19): 198304 |
| [45] | Gravish N, Franklin SV, Hu DL, et al. Entangled granular media.Physical Review Letters, 2012, 108(20): 208001 |
| [46] | 姜世平,芮筱亭,洪俊,等. 散粒体系统动力学仿真. 岩土力学,2011,32(8): 2529-2532 (Jiang Shiping, Rui Xianting, Hong Jun, et al. Dynamic simulation of granular system.Rock and Soil Mechanics, 2011, 32(8): 2529-2532 (in Chinese)) |
| [47] | Hockney RW, Eastwood JW. Computer Simulation Using Particles. CRC Press, 1988 |
| [48] | Aarseth SJ. Direct methods for N-body simulations. Multiple Time Scales, 1985, 1: 377-418 |
| [49] | Ahmad A, Cohen L. A numerical integration scheme for the N-body gravitational problem.Journal of Computational Physics, 1973, 12(3): 389-402 |
| [50] | 刘步林. 改进的 Ahmad-Cohen 方法——引力 N 体问题的一个积分方法. 天文学进展,1985,3(4): 363-370 (Liu Bulin. The improved Ahmad-Cohen scheme: a numerical intergration scheme for the gravitation N-body problem.Progress in Astronomy, 1985, 3(4): 363-370 (in Chinese)) |
| [51] | Dormand JR, Prince PJ. New Runge-Kutta algorithms for numerical simulation in dynamical astronomy.Celestial Mechanics, 1978, 18(3): 223-232 |
| [52] | Chin SA, Chen CR. Forward symplectic integrators for solving gravitational few-body problems. Celestial Mechanics and Dynamical Astronomy, 2005, 91(3-4): 301-322 |
| [53] | Xu J, Wu X. Several fourth-order force gradient symplectic algorithms.Research in Astronomy and Astrophysics, 2010, 10(2): 173-188 |
| [54] | Makino J, Aarseth SJ. On a Hermite integrator with Ahmad-Cohen scheme for gravitational many-body problems. Publications of the Astronomical Society of Japan, 1992, 44(2): 141-151 |
| [55] | Nitadori K, Makino J. Sixth- and eighth-order Hermite integrator for N-body simulations.New Astronomy, 2008, 13(7): 498-507 |
| [56] | Richardson DC, Bottke WF Jr, Love SG. Tidal distortion and disruption of Earth-crossing asteroids.Icarus, 1998, 134(1): 47-76 |
| [57] | Leinhardt ZM, Richardson DC, Quinn T. Direct N-body Simulations of Rubble Pile Collisions.Icarus, 2000, 146(1): 133-151 |
| [58] | Richardson DC, Elankumaran P, Sanderson RE. Numerical experiments with rubble piles: equilibrium shapes and spins.Icarus, 2005, 173(2): 349-361 |
| [59] | Walsh KJ, Richardson DC, Michel P. Rotational breakup as the origin of small binary asteroids.Nature, 2008, 454(7201): 188-191 |
| [60] | Leinhardt ZM, Stewart ST. Collisions between gravity-dominated bodies. I. Outcome regimes and scaling laws. The Astrophysical Journal, 2012, 745(1): 79 |
| [61] | Richardson DC, Walsh KJ, Murdoch N, et al. Numerical simulations of granular dynamics: I. Hard-sphere discrete element method and tests. Icarus, 2011, 212(1): 427-437 |
| [62] | Schwartz SR, Richardson DC, Michel P. An implementation of the soft-sphere discrete element method in a high-performance parallel gravity tree-code.Granular Matter, 2012, 14(3): 363-380 |
| [63] | Sánchez P, Scheeres DJ. Simulating asteroid rubble piles with a self-gravitating soft-sphere distinct element method model.The Astrophysical Journal, 2011, 727(2): 120 |
| [64] | Schwartz SR, Michel P, Richardson DC. Numerically simulating impact disruptions of cohesive glass bead agglomerates using the soft-sphere discrete element method.Icarus, 2013, 226(1): 67-76 |
| [65] | Yu Y, Richardson DC, Michel P, et al. Numerical predictions of surface effects during the 2029 close approach of asteroid 99942 Apophis.Icarus, 2014, 242: 82-96 |
| [66] | Zhao ZJ, Zhao JD, Liu H. Development of a landing mechanism for asteroids with soft surface. International Journal of Aerospace Engineering, 2013, 2013. Article ID 873135. DOI:10.1155/2013/873135 |
| [67] | Zhao ZJ, Zhao JD, Liu H. Landing dynamic and key parameter estimations of a landing mechanism to asteroid with soft surface.International Journal of Aeronautical and Space Sciences, 2013, 14(3): 237-246 |
| [68] | 凌道盛,蒋祝金,钟世英 等. 着陆器足垫冲击模拟月壤的数值分析. 浙江大学学报 (工学版),2013,47(7): 1171-1177 (Ling Daosheng, Jiang Zhujin, Zhong Shiying, et al. Numerical study on impact of lunar lander footpad against simulant lunar soil.Journal of Zhejiang University (Engineering Science), 2013, 47(7): 1171-1177 (in Chinese)) |
| [69] | Gómez LR, Turner AM, van Hecke M, et al. Shocks near jamming.Physical Review Letters, 2012, 108(5): 058001 |
| [70] | Fahnestock EG. The full two-body-problem: Simulation, analysis, and application to the dynamics, characteristics, and evolution of binary asteroid systems. [PhD Thesis]. Michigan: University of Michigan. 2009" |