【问题标题】:Finding the volume fraction of 3D cells cut by a closed 3D surface查找被封闭 3D 表面切割的 3D 细胞的体积分数
【发布时间】:2015-12-18 18:55:21
【问题描述】:

我有一个 3D 三角封闭曲面。该表面浸入矩形笛卡尔网格中。表面以 STL 格式存储。我的目标是计算被表面切割的细胞的比例。换句话说,当一个单元格与表面有交点,实际上被它切割时,它就被分为内部子单元和外部子单元。我的目标是找到内部子细胞和被表面切割的细胞的体积分数。例如,在下图中为了说明目的显示了一个 2D 案例,我想计算青色区域占总区域的比例(即青色 + 灰色)。

有人可以帮我找到一种可以做到这一点的有效方法/算法吗?我可以在MATLABC 中进行实现。

【问题讨论】:

    标签: geometry computational-geometry geometry-surface


    【解决方案1】:

    精确的方法

    首先考虑将曲面与由平面分隔的半空间相交。

    一些三角形完全在半空间内并保持不变。其他一些完全在外面并被丢弃。其余的被切成两半,形成一个三角形和一个四边形。保留右侧部分,如有必要,对四边形进行三角剖分。

    您还需要将生成的面视为位于曲面内部的平面区域。为此,取所有新边(由三角形与平面的截面产生的边),并将它们链接起来。我的意思是将具有共同顶点的线段联系在一起(注意数值准确性);最后你会有一个或多个循环,形成简单的多边形。对这些多边形进行三角剖分(例如通过剪耳法)。

    通过此过程,您将获得一个新的三角曲面,该曲面描述了所请求的交点。对定义单元的六个平面重复此操作,您将得到单元与曲面中包含的体积之间的交集。

    要计算它的体积,您可以将每个三角形及其在固定平面上的投影所形成的棱柱的体积相加;确保以代数方式(带符号)计算这些体积。

    这是一个 2D 类比。 “新边缘”为绿色。


    由于此过程非常耗时,如果您需要对整个网格的每个单元格重复此过程,成本会很高。您可以通过保留中间结果并逐片工作来简化计算。

    垂直于第一个平面方向的排序将加速外部三角形的拒绝。

    【讨论】:

    • 当为 3D 表面(Sutherland-Hogdman 可重入裁剪)执行此操作时,图中突出显示的段变成了多边形,创建它们很繁琐(需要维护组合),但由于只有体积需要,可以创建任意三角形,请参阅我在alice.loria.fr/index.php/…的文章
    【解决方案2】:

    除了 Yves Daoust 的答案中的精确解之外,如果可以接受近似解,为了获得相当有效的东西,我会将表面“光栅化”为更高分辨率的网格,并计算微单元在每个宏小区。为了“光栅化”表面,如果需要性能,我会使用 GPU [1](但这很痛苦)。如果性能不是真正关心的问题,那么您可以在每个三角形上循环,然后对于落在三角形边界框内的每个微单元,测试它是否与三角形有交集。最后,用洪水填充体积的内部。

    [1]https://developer.nvidia.com/content/basics-gpu-voxelization

    【讨论】:

      【解决方案3】:

      我相信您需要通过 PLIC 方法重建流体界面的体积。以下是 VOF 研究社区通常遵循的步骤 {ref: Rider, W. J., & Kothe, D. B. (1998)。重建音量跟踪。计算物理学杂志,141(2), 112-152.}。

      第 1 步:将方程写入/开发到接口。在 PLIC 中,我编写了 ax+by+c = 0 形式的接口。如果您更喜欢二次/样条方程,请这样做。

      第二步:获取所有网格点的坐标。例如:如果您想开发一个 x = 1:5 和 y = 1:5 的网格,请使用以下

      x = 1:5;
      y = 1:5; 
      [X,Y] = meshgrid(x,y);
      

      ** 第 3 步:设置容差值(我使用 10^-7)来定义界面内/外/上的区域。

      对于流体/液体域内的所有点 (X,Y):

      a*X+b*Y+c >tol
      

      流体/液体域外的所有点 (X,Y)

      a*X+b*Y+c <-tol
      

      界面上的点:

      [-tol<=(a*X+b*Y+c)<=tol]
      

      【讨论】:

        猜你喜欢
        • 2016-11-07
        • 2013-02-04
        • 2010-11-27
        • 2013-12-25
        • 1970-01-01
        • 2017-12-19
        • 2011-01-21
        • 2013-11-07
        • 1970-01-01
        相关资源
        最近更新 更多