一种关于混凝土中非凸形颗粒的三维数值构造方法与流程

本发明属于混凝土力学性能分析领域，具体涉及一种关于混凝土中非凸形颗粒的三维数值构造方法。

背景技术

混凝土作为一种颗粒增强复合材料，它的承载能力很大程度上受制于骨料，特别是它的形状和大小。长期以来，人们主要通过实验来分析骨料特征。但是，由于实验过程中出现的一系列不可控因素，经常会导致实验结果离散、可重复性差。近几年，随着计算机的迅猛发展，计算机仿真模拟骨料的几何形状已经在混凝土研究中得到广泛地应用。

关于混凝土骨料形状的相关数值构造工作，前期研究大部分是集中在球体、椭球体、凸面体等规则的几何凸体构造。实际混凝土中的骨料的形状并非完全是完美的凸体颗粒，而是复杂的非凸形颗粒组成。因此，构建非凸形颗粒的几何形貌是实施精确分析混凝土和其他颗粒材料的物理力学性能的基础和保证。相关研究人员提出了球拼装、像素扫描等方法来构造复杂的非凸形颗粒，但是它们存在着一些缺陷。对于球拼装法，很大程度上受到球体数量的限制，不能够精确的模拟出所要构造的非凸形颗粒，而像素扫描法则很大程度上受制于成像设备的分辨率。

技术实现要素：

针对现有技术的不足，本发明的目的在于提供一种关于混凝土中非凸形颗粒的三维数值构造方法，用此方法模拟混凝土中的骨料有效克服了骨料与现实的差距，为混凝土的力学性能提供更有利的保障。

为解决现有技术问题，本发明采取的技术方案为：

一种关于混凝土中非凸形颗粒的三维数值构造方法，包括以下步骤：

步骤1，在o-xyz全局坐标系中随机生成点p0(x0,y0,z0)，利用公式(1)平移矩阵

将原点(0，0，0)平移至p0(x0,y0,z0)形成局部坐标系o1-x1y1z1，p0(x0,y0,z0)即为局部坐标系o1-x1y1z1的原点，将局部坐标系x1o1y1面内(0,0)点作为极坐标的极心，若生成多边形的边数为n，则在区间[0,2π)内随机的取n个不同的极角值，可以通过公式(2)生成第i个顶点所对应的极角，将所取的极角按照从小到大的顺序排列,

θi＝ηi×2π

(2)，其中，θi是第i个顶点所对应的极角，ηi是[0,1)内的一个随机数；

在生成多边形各个顶点对应的极角的同时，随机地生成对应的极径，极径随机的分布在(a0-a1)到(a0+a1)之间，通过公式(3)可以随机地生成第i个顶点所对应的极径ri

ri＝a0+αi×a1(3)，

其中，a0是颗粒的平均粒径；αi是[-1,1]内的一个随机数，a1是颗粒的最小粒径，ri是第i个顶点所对应的极径；

根据每个顶点所对应的极角和极径，通过公式(4)可以计算出在局部直角坐标系o1-x1y1z1下的第i个顶点的顶点坐标p′i(x′i，y′i，z′i)，i＝1，2，...，n；

xi＝x0+ricos(θi)

yi＝y0+risin(θi)

zi＝z0

(4)，

其中，ri和θi分别是第i个顶点对应的极径和极角；pi′(x′i,y′i,z′i)是在局部直角坐标系下的第i个顶点的顶点坐标,i＝1，2…n；

步骤2，将步骤1生成的多边形通过平移变换和欧拉变换扩展至三维空间中得多边形，并判断所得多边形中所有顶点的凹凸性

根据步骤1所得顶点坐标，通过平移变换和欧拉变换，将其在各自范围内随机地生成三个欧拉角α、β和γ。其中，α和γ的取值区间在[-π,π]，β的取值区间在[0,π]，计算其余弦和正弦值，得到如公式(5)所示的旋转矩阵，顺时针旋转的角度为正值，逆时针时旋转的角度为负值；

其中，r是旋转矩阵，通过公式(6)，将在局部坐标o1-x1y1z1求得的第i个顶点的坐标pi′(x′i,y′i,z′i)，变换成在全局坐标系o-xyz下的顶点坐标pi(xi,yi,zi)

x＝trx′

(6)，其中，x′是局部坐标系下的坐标，x是全局坐标系下的坐标，

通过多边形各个顶点转向来判断顶点的凹凸性，可以通过公式(7)计算出过顶点a(x1,y1,z1)、点b(x2,y2,z2)和点c(x3,y3,z3)的多边形的法向量n0，通过公式(8)来判断多边形的顶点a(x1,y1,z1)的凹凸性，其中，公式(7)中的点b(x2,y2,z2)和点c(x3,y3,z3)是与顶点a(x1,y1,z1)相邻的两个顶点，依次判断完多边形的每个顶点的凹凸性。若多边形的所有顶点为凸顶点，此时生成的多边形为无效多边形，重新生成；

其中，n0是过顶点(x1,y1,z1)的多边形法向量，v0是以(x1,y1,z1)，(x0,y0,z0)为端点的多边形的棱边的方向向量，v1是以(x1,y1,z1)，(x2,y2,z2)为端点的非凸多边形的棱边的方向向量；

其中，n是非多边形的法向量，n0是过顶点(x1,y1,z1)的多边形法向量。λ是判定顶点凹凸性的参数：当λ>0时，该顶点为凸顶点；当λ<0时，该顶点为凹顶点；

步骤3，确定扩展球半径和球心的位置

根据粒子的厚度确定扩展球的半径，根据步骤2得到的非凸多边形的顶点坐标，任意取一个顶点作为扩展球体的球心；

步骤4，确定扩展球的滚动方向，令扩展球沿非凸多边形的边滚动

根据步骤2得到的顶点坐标，确定非凸多边形各个边的方向，以任一顶点为扩展球体的球心,让扩展球体沿着经过该点的边滚动，滚动至球心位于邻近顶点位置，滚动的轨迹构成一个球柱体，按照相同的滚动方式，依次滚完非凸多边形的其余的边；

步骤5，确定新的表面顶点坐标

根据步骤4滚动所得的球柱体，首先通过公式(9)计算得出初始状态的法向量n(x0,y0,z0)，基于该法向量通过公式(10)计算得出新的表面顶点坐标(x,y,z)

x0＝y1z2+y2z3+y3z1-y1z3-y2z1-y3z2

y0＝x1z3+x2z1+x3z2-x1z2-x2z3-x3z1

z0＝x1y2+x2y3+x3y1-x1y3-x2y1-x3y2(9)，

其中，(x1,y1,z1)、(x2,y2,z2)、(x3,y3,z3)为非凸多边形平面内任意不共线的三点坐标，通常随机地取其三个顶点，计算出初始状态的法向量；(x0,y0,z0)是初始状态的法向量n，即非凸多边形的法向量；

其中，(x,y,z)为扩展球的球心坐标；r是扩展球体的半径；(x0,y0,z0)是初始状态的法向量n。

有益效果

与现有技术相比，本发明在计算机图形学基础上，构造出三维状态下的非凸形颗粒，克服了模拟出的颗粒形状和现实相差很大的问题。

本发明构造方法简单，操作简便，且对硬件要求低，成本低，为混凝土性能的计算提供了更有利的证据。

附图说明

图1为本发明的三维非凸形颗粒构造的流程图；

图2为模拟生成的在三维空间中的一个非凸多边形单元示意图；

图3为扩展球体所在非凸多边形单元某个顶点的示意图；

图4为扩展球体沿着非凸多边形单元的棱边滚动示意图；

图5为生成的新的表面顶点坐标示意图；

图6为实现的边数为8的三维非凸形颗粒的示例图。

具体实施方式

下面结合附图和具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围。在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

本发明的思路是在现有的计算机图形学理论的基础上，实现了三维非凸形颗粒的构造。具体而言，本发明采用以下技术方案解决上述技术问题。

对于需要构造的三维非凸形颗粒，将在三维空间中的非凸多边形作为基本单元体。扩展单元体设为半径为r的球体，且球心在基本单元体表面上。具体的构造方法图如流程图1所示。

一种关于混凝土中非凸形颗粒的三维数值构造方法，包括以下步骤：

步骤1，在o-xyz全局坐标系中随机生成点p0(x0,y0,z0)，利用公式(1)平移矩阵(1)，将原点(0，0，0)平移至p0(x0,y0,z0)形成局部坐标系o1-x1y1z1，p0(x0,y0,z0)即为局部坐标系o1-x1y1z1的原点，将局部坐标系x1o1y1面内(0,0)点作为极坐标的极心，若生成多边形的边数为n，则在区间[0,2π)内随机的取n个不同的极角值，可以通过公式(2)生成第i个顶点所对应的极角，将所取的极角按照从小到大的顺序排列,

θi＝ηi×2π

(2)，其中，θi是第i个顶点所对应的极角，ηi是[0,1)内的一个随机数；

在生成多边形各个顶点对应的极角的同时，随机地生成对应的极径，极径随机的分布在(a0-a1)到(a0+a1)之间，通过公式(3)可以随机地生成第i个顶点所对应的极径ri

ri＝a0+αi×a1(3)，

其中，a0是颗粒的平均粒径；αi是[-1,1]内的一个随机数，a1是颗粒的最小粒径，ri是第i个顶点所对应的极径；

根据每个顶点所对应的极角和极径，通过公式(4)可以计算出在局部直角坐标系o1-x1y1z1下的第i个顶点的顶点坐标p＇i(x＇i，y＇i，z＇i)，i＝1，2，...，n；

其中，ri和θi分别是第i个顶点对应的极径和极角；pi′(x′i,y′i,z′i)是在局部直角坐标系下的第i个顶点的顶点坐标,i＝1，2…n；

步骤2，，将步骤1生成的多边形通过平移变换和欧拉变换扩展至三维空间中得多边形，并判断所得多边形中所有顶点的凹凸性

根据步骤1所得顶点坐标，通过平移变换和欧拉变换，将其在各自范围内随机地生成三个欧拉角α、β和γ。其中，α和γ的取值区间在[-π,π]，β的取值区间在[0,π]，计算其余弦和正弦值，得到如公式(5)所示的旋转矩阵。顺时针旋转的角度为正值，逆时针时旋转的角度为负值；

其中，r是旋转矩阵，通过公式(6)，将在局部坐标o1-x1y1z1求得的第i个顶点的坐标p′i(x′i,y′i,z′i)，变换成在全局坐标系o-xyz下的顶点坐标pi(xi,yi,zi)

x＝trx′(6)，

其中，x′是局部坐标系下的坐标，x是全局坐标系下的坐标，

通过多边形各个顶点转向来判断顶点的凹凸性，可以通过公式(7)计算出过顶点a(x1,y1,z1)、点b(x2,y2,z2)和点c(x3,y3,z3)的多边形的法向量n0，通过公式(8)来判断多边形的顶点a(x1,y1,z1)的凹凸性，其中，公式(7)中的点b(x2,y2,z2)和点c(x3,y3,z3)是与顶点a(x1,y1,z1)相邻的两个顶点，依次判断完多边形的每个顶点的凹凸性。若多边形的所有顶点为凸顶点，此时生成的多边形为无效多边形，重新生成；

其中，n0是过顶点(x1,y1,z1)的多边形法向量，v0是以(x1,y1,z1)，(x0,y0,z0)为端点的多边形的棱边的方向向量，v1是以(x1,y1,z1)，(x2,y2,z2)为端点的非凸多边形的棱边的方向向量；

其中，n是非多边形的法向量，n0是过顶点(x1,y1,z1)的多边形法向量。λ是判定顶点凹凸性的参数：当λ>0时，该顶点为凸顶点；当λ<0时，该顶点为凹顶点，如图2所示；

步骤3，设扩展球体的半径为r，以某个顶点为扩展球体的球心,如图3所示；

步骤4，扩展球体以该顶点位置为出发点，沿着棱边滚动至球心位于另一顶点位置，滚动的轨迹构成一个球柱体，球柱体的上、下顶点分别对应于该棱边的两个顶点，则该圆柱体表露在颗粒外围的区域即三维非凸颗粒的边，如图4所示；

步骤5，根据所得的球柱体，计算出新的表面顶点位置。具体细节如下:假设，非凸多边形的三个顶点坐标为(x1,y1,z1)、(x2,y2,z2)、(x3,y3,z3)，通过公式(1)可以求出，基本单元非凸多边形的平面法向量n(x0,y0,z0)

x0＝y1z2+y2z3+y3z1-y1z3-y2z1-y3z2

y0＝x1z3+x2z1+x3z2-x1z2-x2z3-x3z1

z0＝x1y2+x2y3+x3y1-x1y3-x2y1-x3y2(9)

其中，(x1,y1,z1)、(x2,y2,z2)、(x3,y3,z3)是非凸多边形的三个顶点坐标；(x0,y0,z0)是基本单元的法向量，即非凸多边形的法向量。

然后通过公式(10)可以计算出现的表面顶点位置，从而构造出三维的非凸形颗粒。

其中，(x,y,z)为扩展球的球心坐标；r是扩展球体的半径；(x0,y0,z0)是初始状态的法向量n。

实施例2

图6展示了，基于边数n＝8的非凸边形构造出三维的非凸形颗粒。将扩展球体的半径r设为1，其球心为非凸多边形的某个顶点；让扩展球沿着非凸多边形的棱边滚动；而后计算其新的表面的顶点坐标。最终得到图6最右侧的三维非凸形颗粒。