离散元课件二学习教案

会计学1离散离散(lsn)元课件二元课件二第一页，共38页。第1页/共38页第二页，共38页。二二 基本原理基本原理-球形颗粒球形颗粒(kl)元离散元法元离散元法第2页/共38页第三页，共38页。第3页/共38页第四页，共38页。gFFVikjijdijciimtmi1,dd式中， 与 分别为颗粒 的质量和速度。 为时间(shjin)， 为颗粒的重力， 与 分别为颗粒 与 的接触力与粘性接触阻尼力， 为所有与颗粒接触的颗粒总数。 imiVitgimijc,Fijd,Fijik第4页/共38页第五页，共38页。jiijctijcnijc,FFF 同理，粘性接触阻尼力也可分解为法向与切向分量(fn ling)形式，即 ijdtijdnijd,FFF 第5页/共38页第六页，共38页。颗粒间的接触力作用在两个颗粒的接触点上，而不是(b shi)作用在颗粒的中心，所以这些接触力（除法向接触力 外）将会对颗粒产生力矩 ， ijcn,FiT ijdtijctii,FFRT 式中， 为从颗粒 的质心指向接触点的矢量，其幅值为 （颗粒的半径）。iRiiR第6页/共38页第七页，共38页。转动(zhun dng)方程可以表示为 ikjiiitI1ddT式中， 与 分别为颗粒 的转动惯量与角速度，对于球形颗粒 为 iIii252iiiRmI iI第7页/共38页第八页，共38页。关于接触力的计算模型已有大量的研究成果，目前仍旧是一个活跃的研究领域，特别是对于切向力的计算方法。对于理想散体颗粒（无粘连）：采用Hertz理论描述法向作用(zuyng)，而采用Mindlin与Deresiewicz理论描述切向作用(zuyng)；对于存在粘连的散体颗粒：法向接触力根据在Hertz理论基础上考虑粘连力的JKR（Johnson-Kendall-Roberts）理论确定，切向接触力增量则根据把Savkoor和Briggs理论与Mindlin和Deresiewicz理论相结合形成的理论确定。 第8页/共38页第九页，共38页。iiR/R n单位(dnwi)法向向量 单位切向向量 nnnnt ijijijijVVVV单位切向量之所以通过两个颗粒的相对速度来计算，是因为接触力与粘性阻尼力的方向与相对速度的方向相同。 第9页/共38页第十页，共38页。法向相对速度(xin du s d)为 iijjRRVVVijijnnijijn,VV切向相对速度为 nnijijVVVijt,nn ijijt,VV或者写为 第10页/共38页第十一页，共38页。 为颗粒i与j接触时的侵入(qnr)深度n2/3*,34nijcnRE F式中)1 (22*vEE jiRRR11* nijjinRRRR 第11页/共38页第十二页，共38页。nnnijcnk ,F式中， 为法向弹簧(tnhung)刚度。 nk第12页/共38页第十三页，共38页。式中， 为法向粘性(zhn xn)接触阻尼系数。 nnijVFnijdnc,nc第13页/共38页第十四页，共38页。处于接触中的两个颗粒(kl)的切向作用，从本质上讲，是一种摩擦行为，按照摩擦机理，摩擦力包括：滑动摩擦、滚动摩擦与静摩擦，其中滑动摩擦与静摩擦属于切向摩擦力；滚动摩擦是由于法向接触应力的不均匀分布产生的。介绍两个切向接触力模型： Coulomb准则Mindlin与Deresiewicz切向接触力模型 第14页/共38页第十五页，共38页。 在离散元模拟中，一般用Coulomb准则这种简单的形式描述，静摩擦的详细刻画需要涉及切向位移甚至可能要考虑时间(shjin)依赖效应。式中， 为静摩擦系数，切向摩擦力的方向为与相对滑动的趋势相反。 sijcnsijctijcnsijcnsijctijctijct,FFFFFFF第15页/共38页第十六页，共38页。 tmaxt,maxt,ttijcnsijctFF 23,min11式中 为颗粒(kl) 与 间的累积切向位移矢量tij nsmaxt,vv 122 tNtNtnn ijV1第16页/共38页第十七页，共38页。式中， 为切向粘性接触阻尼(zn)系数。 nn ijVFtijdtc,tc第17页/共38页第十八页，共38页。n接触力的计算接触力的计算n法向接触力法向接触力n tmaxt,maxt,ttijcnsijctFF 23,min11gFFVikjijdijciimtmi1,dd ikjiiitI1ddTn2/3*,34nijcnRE Fnnnijcnk ,Fijcnsijctijcnsijcnsijctijctijct,FFFFFFFnnijVFnijdnc, nn ijVFtijdtc,切向接触力切向接触力 第18页/共38页第十九页，共38页。gFFVikjijdijciimtmi1,dd ikjiiitI1ddT一般采用两种方法一般采用两种方法(fngf)求解运动方程：求解运动方程：中心差分法中心差分法Verlet积分法积分法第19页/共38页第二十页，共38页。hXXXnnn 2121hXXXnnn211运动方程运动方程(fngchng)的求解的求解中心差分法中心差分法第20页/共38页第二十一页，共38页。运动方程运动方程(fngchng)的求解的求解Verlet积分法积分法第21页/共38页第二十二页，共38页。hXXXnnn211hXXXnnn211 22242hXXXnnn2112hXXXnnn运动方程运动方程(fngchng)的求解的求解Verlet积分法积分法定义定义(dngy)：从而导出Verlet方程为2112hXXXXnnnn 第22页/共38页第二十三页，共38页。第23页/共38页第二十四页，共38页。第24页/共38页第二十五页，共38页。 当需要判别体系中某个颗粒的邻居数量时，在该粒子周围构建一个球（称之为参考球，称该颗粒为核心颗粒），参考球半径为体系中最大粒子半径的若干倍，那么参考球所包围的所有粒子为该球中心粒子的邻居。参考球半径的选取取决于粒子的运动速度及体系中粒子的密度。对于(duy)每个粒子，都可生成一个邻居粒子的目录。为了得到邻居目录，对每一个粒子而言，所有标号大于该粒子的粒子都必须被检验，判断是否位于该粒子的参考球中，而对于(duy)标号小于该粒子标号的粒子则没有必要被检验，因为邻居是互相的，没有必要对一个邻居对检验两次。 第25页/共38页第二十六页，共38页。Verlet邻居目录法及粒子邻居目录法及粒子(lz)存储目录存储目录接触接触(jich)发现算法发现算法Verlet邻居目录法邻居目录法第26页/共38页第二十七页，共38页。接触发现算法接触发现算法Verlet邻居邻居(ln j)目录法目录法对于n个颗粒组成的体系，用Verlet邻居目录法需要 次计算，也就是说计算次数仍旧为 量级。然而，并不需要在每个时间步长上都对邻居目录进行更新。更新的频率取决于体系中粒子的密度、粒子的运动(yndng)速度以及参考球的尺寸。参考球的半径也可以根据颗粒体系的稠密程度及运动(yndng)速度进行动态调整，并且参考球半径与邻居目录的更新频率呈反比关系，参考球半径越小，邻居目录的更新频率越高；但是参考球半径越大，则有更多的粒子位于球体内，所以判别是否为邻居就需要较长的时间。2/1nn 2no第27页/共38页第二十八页，共38页。将颗粒体系所占据的空间划分成规则的网格，对于三维问题，可以划分为 个立方体单元，对于二维体系，则划分为 个正方形单元，对于颗粒体系所占据空间形状不规则时，也可采用其他形状单元划分。但是所有单元的尺寸必须大于粒子的尺寸。与Verlet邻居目录法的主要区别(qbi)在于：相邻单元法中的单元不依附于粒子，单元不随粒子的运动而运动。如果粒子就当前的位置被分配到某个单元，显然，只有在同一个单元或直接相邻单元内的粒子间才可能发生相互作用，也就是说，只有相邻单元内的粒子才能成为邻居。 mmmmm第28页/共38页第二十九页，共38页。例如，对一个二维体系而言，只有在9个不同的单元内可能包含邻居粒子，对于三维体系，则只有在27个不同单元内可能包含邻居粒子。与前面介绍的Verlet邻居目录法相同，每个粒子对只需要检验一次，这样，没必要对9个单元都进行检验，只需检验中心单元及邻居单元的一半(ybn)即可，即在二维情况下，只需检验5个单元，在三维情况下，只需检验14个单元。 第29页/共38页第三十页，共38页。2D体系体系(tx)中的连接单元法中的连接单元法第30页/共38页第三十一页，共38页。这个方法与前述两种方法不同。首先，在每一个粒子周围构建一个边界(binji)盒，边界(binji)盒的尺寸按这样的方式选取：使每个粒子刚好放进它的边界(binji)盒内。边界(binji)盒的边为直线，并且与体系的坐标轴平行。在判别颗粒的邻居时，把边界(binji)盒投影到体系的坐标轴上。通过边界(binji)盒在坐标轴上投影的起点和终点来判别是否为邻居。 第31页/共38页第三十二页，共38页。 a a 围绕围绕在每在每个粒个粒子周子周围围(zhu(zhuwi)wi)的边的边界盒界盒 b b两个不同两个不同时刻粒子时刻粒子的边界盒的边界盒在轴上的在轴上的投影投影 在判别颗粒的邻居时，把边界盒投影到体系的坐标轴上。例如，图a表明了粒子及边界盒的位置，图b为图a中的边界盒在体系x轴的投影。通过边界盒在坐标轴上投影的起点和终点来判别是否为邻居，出于这个原因，投影的起点和终点序列被存储在目录中。 第32页/共38页第三十三页，共38页。对于3D体系，必须把边界盒子在三个坐标轴上投影，所以生成三个目录。如果一个粒子的边界盒在某个坐标轴上投影的起点和终点间包含(bohn)另外一个粒子边界盒投影的起点、终点或起点和终点，就说明这两个粒子的边界盒在该坐标轴上的投影发生了重叠。如果两个边界盒的投影在每一个坐标轴上都发生重叠，那么就说明这两个边界盒发生了接触。 第33页/共38页第三十四页，共38页。存储目录的更新存储目录的更新检查一个边界盒在每一个坐标轴的投影是否在另外一个边界盒投影的起点和终点之间仍旧需要花费检查一个边界盒在每一个坐标轴的投影是否在另外一个边界盒投影的起点和终点之间仍旧需要花费(hufi)(hufi)很多计算时间。但是，尽管这些目录在每个计算时间步长上都必须被更新，所用的计算时间仍旧可以被减少到与体系中粒子数量成比例的量级，因为在每一个新步长上所需要做的只是对旧目录的更新，并且对旧目录的更新只是对旧目录的重新分类，因而，更新的过程非常简单，只需要对旧目录进行顺序地检查，判断是否在顺序上有新的变化即可。这种变化仅仅是位置的变化。很多计算时间。但是，尽管这些目录在每个计算时间步长上都必须被更新，所用的计算时间仍旧可以被减少到与体系中粒子数量成比例的量级，因为在每一个新步长上所需要做的只是对旧目录的更新，并且对旧目录的更新只是对旧目录的重新分类，因而，更新的过程非常简单，只需要对旧目录进行顺序地检查，判断是否在顺序上有新的变化即可。这种变化仅仅是位置的变化。 第34页/共38页第三十五页，共38页。第35页/共38页第三十六页，共38页。sh)第36页/共38页第三十七页，共38页。只有当选取的时间步长小于临界时间步长时，数值计算才是稳定的，这是模型(mxng)采用显式积分的结果。临界时间步长可以通过计算估计，然后，根据临界时间步长确定计算所采用的时间步长。可以采用质量为m、用刚度k为的弹簧与地面相连接的单自由度体系来估计，这种情况下，临界步长为: )/(kmct 时间时间(shjin)步长的确定步长的确定第37页/共38页第三十八页，共38页。