如何使用等角三角形实现自适应网格细化

计算科学 有限元 算法 C++ 自适应网格细化
2021-12-08 03:06:19

我正在尝试为有限元代码实现自适应网格细化。代码使用(至少现在)线性三角形,所以当我进行网格细化时,我希望三角形网格保持保形,即没有悬挂节点:

在此处输入图像描述

我计划使用 Rivara [33,34] 的最长边二等分法。不幸的是,我无法访问这些论文,因为我还没有找到可免费下载的 pdf。我查看了描述此方法的其他几篇论文/资源,即fea8.pdf以及 Wolfgang Bangerth 的视频系列讲座 14-18,但是我对如何实际实现该算法的细节有疑问。

到目前为止我所做的是假设一个初始的粗网格,即:

在此处输入图像描述

在我的程序中,我的网格由二叉树数组表示,其中每个原始粗三角形是数组中二叉树之一的根节点。因此,例如在上面的网格中,我将有一个大小为 2 的二叉树数组。从图形上看,我的树看起来像:

在此处输入图像描述

现在说例如三角形e1被标记为细化。我们用最长边平分这个三角形,创建两个新的子三角形:

在此处输入图像描述

和树结构:

在此处输入图像描述

问题是那个三角形f1有一个悬挂节点,最糟糕的是三角形树中的根节点f1不“知道”它沿着该边缘有一个悬挂节点,也不“知道”这个悬挂节点的编号为 5。我想出的一个解决方案是让每个节点代表一个三角形,以使其还有额外的指针指向那个三角形的邻居。然而,这很快使树数组变得非常复杂,到处都是节点之间的指针,试图跟踪哪些三角形是其他三角形的邻居。

例如,在 C++ 中,我的二叉树数据结构类似于:

template <class T>
class Tree
{
  private:
    struct Node
    {
      T element;
      Node *parent;
      Node *left;
      Node *right;

      Node *neigh01;
      Node *neigh12;
      Node *neigh20;

      int hnode01;   //store hanging node numbers for this node
      int hnode12;  
      int hnode20;   

      bool refine;   //refinement flag saying this node is in need of refinement

      Node(T e) : element(e), parent(NULL), left(NULL), right(NULL), neigh01(NULL), neigh12(NULL), neigh20(NULL), hnode01(-1), hnode12(-1), hnode20(-1)
      {
      }
    };

    //rest of class ...constuctors, methods, etc...
};

你可以看到我已经包含了额外的节点指针,它们试图跟踪邻居节点的存储位置。实际上实现这种自适应网格细化的二等分方法似乎相当复杂,我不知道如何去做。我的想法是有额外的指针来跟踪邻居是解决这个问题的正确方法吗?有没有人有免费可用的资源来更详细地解释二分算法,即解释使用的数据结构等?谁能指出我使用二等分进行自适应网格细化的算法?用fea8.pdf写的那个似乎表明我们首先遍历所有需要细化的节点并将它们拆分。然后我们通过(迭代地?)纠正挂起的节点。这是正确的/通常会做什么?我将我的网格表示为二叉树数组的想法是表示将自适应细化的网格的良好/标准方式吗?基本上我需要任何帮助来实现这种自适应网格细化算法。

谢谢

[33] MC里瓦拉。完全自适应多重网格有限元软件的设计和数据结构。ACM 数学软件交易,10:242 264,1984。

[34] MC里瓦拉。基于单纯形的广义二等分的网格细化过程。SIAM 数值分析杂志,21:604 613,1984。

2个回答

通常处理这种情况的方式是在您实际拆分单元格之前运行一个通道,以确定哪些其他单元格也必须进行细化。在伪代码中,这可能如下所示:

for (cell in cells)
  if (cell->is_marked_for_refinement && cell->edge_to_split == unassigned)
  {
     e = longest edge of cell (0...2);
     cell->set_edge_to_split(e);
     n = cell->neighbor(e);
     n->mark_for_refinement();  // if we split this cell along edge e,
                                // then we also need to split cell n
     n_of_n = n->get_neighbor_index(cell); // find out along which edge
                                // of n the current cell is
     n->set_edge_to_split (n_of_n);
  }

此过程可确保您首先标记所有实际需要拆分的单元格,然后再进行拆分。

在沃尔夫冈斯回答的帮助下,我想我已经弄清楚了。有必要在我的 Node 类中包含额外的指针,以便二叉树中的每个 Node 都可以知道谁是它的邻居。例如,我的 Node 和 Tree 类如下所示:

template <class T>
class Node
{
  public:
    T element;
    Node *parent;
    Node *left;
    Node *right;

  private:
    Node *neigh01;
    Node *neigh12;
    Node *neigh20;

    int hnode01;   //store hanging node numbers for this node
    int hnode12;  
    int hnode20;   

    bool refine;   //refinement flag saying this node is in need of refinement

  public:
    Node(T e) : element(e), parent(NULL), left(NULL), right(NULL), neigh01(NULL), neigh12(NULL), neigh20(NULL), hnode01(-1), hnode12(-1), hnode20(-1)
    {
    }

  //rest of methods
};


template <class T>
class Tree
{
  private:
    int numNodes;
    int numLeafNodes
    Node<T> *root;
    Node<T> *iterator;

  public:
    Tree();
    Tree(T e);
    ~Tree();

  void insert(T e);
  void deleteTree(Node<T> *leaf);

  //rest of class methods, etc...
};

三个额外的节点指针——neigh01、neigh12 和 neigh20——允许我们跟踪任何相邻节点。从图形上看,我们的树结构可能如下所示:

在此处输入图像描述 在此处输入图像描述 在此处输入图像描述 在此处输入图像描述

其中红线表示相邻节点连接。当一个三角形被细化(即被最长边二等分)时要记住的重要一点是这些新的子节点只能是邻居:

  1. NULL(在有边界的情况下,即没有邻居)
  2. 彼此(每个左子节点总是有一个相邻的右子节点)
  3. 子节点的父节点的邻居或子节点(作为父节点的邻居或子节点的任何节点都可能是新创建的子节点的邻居)

这有点棘手,但是通过一些浅层递归函数,我能够确保每次精炼节点时都分配了正确的邻居指针。以下是细化前后测试网格的结果:

在此处输入图像描述 在此处输入图像描述

事实上,网格仍然是保形的!