Введение
В предыдущем посте мы рассмотрели особенности использования сторонних библиотек с открытым кодом Freefem++ и NetGen в программе моделирования аэродинамических процессов. Речь шла о возможности включения этих библиотек в коммерческий проект с позиции лицензирования, об особенностях выполнения функций и включения в программную архитектуру. Данная статья является дополнением к предыдущей, в ней более детально рассмотрим библиотеку NetGen. Интерес представляют функции генерации 3D-сетки конечных элементов с заданными регионами на поверхности модели.
Подключение NetGen в MS Visual Studio
Покажем, как можно подключить NetGen к программе на C++. С официального сайта проекта NetGen скачиваем архив. В нашем случае был доступен NetGen версии 4.9.13. Для подключения библиотеки нужны три папки из архива: libsrc, nglib, windows. В Visual Studio создаем проект и решение для демонстрационной программы и подключаем существующий файл проекта nglib.vcxproj из папки windows. Так как в нашем эксперименте проект nglib был создан в более ранней версии Visual Studio, то попутно выполняем конвертацию проекта под новую версию. В настройках проекта nglib добавляем include-папку libsrcinclude, добавляем символ NO_PARALLEL_THREADS, удаляем pthreadVC2.lib из дополнительных зависимостей компоновщика, удаляем пост-событие компоновщика и настраиваем целевые папки так, чтобы проект собирался без ошибок.
Включение заголовочного файла nglib.h надо записывать следующим не очень привычным способом:
namespace nglib { #include "nglib.h" } using namespace nglib;
Тестовая модель
Пространство, которое должно быть заполнено конечными элементами, может задаваться для NetGen двумя способами:
- В виде операторов конструктивной блочной геометрии (CSG).
- В виде описания границы пространства в файле формата STL.
Здесь рассмотрим только второй вариант описания поверхности – описание трехмерной модели в текстовом формате STL. Графическое изображение тестовой модели показано на рис. 1. Модель представляет собой параллелепипед размером 10 Х 10 Х 5. На одной из граней располагается прямоугольный регион размером 5 Х 2 со смещением (4;2), границы которого должны оставаться неизменными
Рис. 1. Триангуляция поверхности для формирования STL файла. Красным цветом выделен регион, границы которого должны оставаться без изменений
Ниже приводится содержимое STL-файла, описывающего тестовую модель:
solid facet normal 1 0 0 outer loop vertex 10 4 2 vertex 10 6.5 3 vertex 10 4 4 endloop endfacet facet normal 1 0 0 outer loop vertex 10 6.5 3 vertex 10 9 4 vertex 10 4 4 endloop endfacet facet normal 1 0 0 outer loop vertex 10 9 2 vertex 10 6.5 3 vertex 10 4 2 endloop endfacet facet normal 1 0 0 outer loop vertex 10 9 2 vertex 10 9 4 vertex 10 6.5 3 endloop endfacet facet normal 1 0 0 outer loop vertex 10 0 5 vertex 10 0 0 vertex 10 4 2 endloop endfacet facet normal 1 0 0 outer loop vertex 10 4 2 vertex 10 4 4 vertex 10 0 5 endloop endfacet facet normal 1 0 0 outer loop vertex 10 10 5 vertex 10 0 5 vertex 10 4 4 endloop endfacet facet normal 1 0 0 outer loop vertex 10 9 4 vertex 10 10 5 vertex 10 4 4 endloop endfacet facet normal 1 0 0 outer loop vertex 10 9 2 vertex 10 10 5 vertex 10 9 4 endloop endfacet facet normal 1 0 0 outer loop vertex 10 10 0 vertex 10 10 5 vertex 10 9 2 endloop endfacet facet normal 1 0 0 outer loop vertex 10 10 0 vertex 10 9 2 vertex 10 4 2 endloop endfacet facet normal 1 0 0 outer loop vertex 10 10 0 vertex 10 4 2 vertex 10 0 0 endloop endfacet facet normal -1 0 0 outer loop vertex 0 10 5 vertex 0 10 0 vertex 0 0 5 endloop endfacet facet normal -1 0 0 outer loop vertex 0 0 0 vertex 0 0 5 vertex 0 10 0 endloop endfacet facet normal 0 -1 0 outer loop vertex 0 0 5 vertex 0 0 0 vertex 10 0 5 endloop endfacet facet normal 0 -1 0 outer loop vertex 10 0 0 vertex 10 0 5 vertex 0 0 0 endloop endfacet facet normal 0 1 0 outer loop vertex 10 10 5 vertex 10 10 0 vertex 0 10 5 endloop endfacet facet normal 0 1 0 outer loop vertex 0 10 0 vertex 0 10 5 vertex 10 10 0 endloop endfacet facet normal 0 0 1 outer loop vertex 0 0 5 vertex 10 0 5 vertex 0 10 5 endloop endfacet facet normal 0 0 1 outer loop vertex 10 10 5 vertex 0 10 5 vertex 10 0 5 endloop endfacet facet normal 0 0 -1 outer loop vertex 10 0 0 vertex 0 0 0 vertex 10 10 0 endloop endfacet facet normal 0 0 -1 outer loop vertex 0 10 0 vertex 10 10 0 vertex 0 0 0 endloop endfacet endsolid
Информация в STL-файле имеет следующий формат. В начале и конце файла указываются ключевые слова solid и endsolid, между которыми перечисляются треугольники, задающие требуемую поверхность в формате:
facet normal nx ny nz outer loop vertex v1x v1y v1z vertex v2x v2y v2z vertex v3x v3y v3z endloop endfacet
где nx, ny, nz – составляющие вектора нормали по осям X, Y, Z, направленного наружу из моделируемого объекта; (v1x, v1y, v1z), (v2x, v2y, v2z), (v3x, v3y, v3z) – координаты трех вершин треугольника.
Пример программы, генерирующей сетку конечных элементов
Ниже приведен текст программы, выполняющей все действия, необходимые для получения сетки конечных элементов. Данная программа читает STL-файл OneRegionBar.stl и записывает сетку конечных элементов в файл OneRegionBar.vol.
#include "stdafx.h" #include <iostream> #include <fstream> namespace nglib { #include "..NetGennglibnglib.h" } using namespace std; using namespace nglib; void main() { const char *stlFileName = "..\Data\OneRegionBar.stl"; const char *volFileName = "..\Data\OneRegionBar.vol"; // Ребра на поверхности const int EDGES_NUM = 4; const int POINTS_NUM = 4; typedef double Point[3]; Point points[POINTS_NUM] = {{10, 4, 2}, {10, 9, 2}, {10, 9, 4}, {10, 4, 4}}; typedef pair<int, int> Edge; Edge edges[EDGES_NUM] = {Edge(0,1), Edge(1, 2), Edge(2,3), Edge(3, 0)}; // Параметры генерации сетки Ng_Meshing_Parameters ngMeshParameters; ngMeshParameters.maxh = 10; ngMeshParameters.fineness = 0.4; ngMeshParameters.secondorder = 0; ngMeshParameters.quad_dominated = 0; ngMeshParameters.grading = 0.0; ngMeshParameters.optsurfmeshenable = false; ngMeshParameters.optsteps_2d = 0; ngMeshParameters.closeedgeenable = false; // Результат выполнения операций Ng_Result ng_res; // // Начало построения сетки // // Чтение STL-файла, создание объекта "STL-геометрии" Ng_STL_Geometry *stlGeometry = Ng_STL_LoadGeometry(stlFileName); // Добавление в "геометрию" неизменяемых ребер for(int i = 0; i < EDGES_NUM; i++) { double* p1 = points[edges[i].first]; double* p2 = points[edges[i].second]; Ng_STL_AddEdge(stlGeometry, p1, p2); } // Инициализация объекта с описанием STL ng_res = Ng_STL_InitSTLGeometry(stlGeometry); if(ng_res != NG_OK) { cout << "Error in Ng_STL_InitSTLGeometry" << endl; } // Пока что пустой объект сетки элементов Ng_Mesh *mesh = Ng_NewMesh(); // Формировать ребра на поверхности модели ng_res = Ng_STL_MakeEdges(stlGeometry, mesh, &ngMeshParameters); if(ng_res != NG_OK) { cout << "Error in Ng_STL_MakeEdges" << endl; } // Генерировать поверхность пространства ng_res = Ng_STL_GenerateSurfaceMesh(stlGeometry, mesh, &ngMeshParameters); if(ng_res != NG_OK) { cout << "Error in Ng_STL_GenerateSurfaceMesh" << endl; } // Генерировать сеть конечных элементов ng_res = Ng_GenerateVolumeMesh(mesh, &ngMeshParameters); if(ng_res != NG_OK) { cout << "Error in Ng_GenerateVolumeMesh" << endl; } // Сохранить сетку конечных элементов в файл Ng_SaveMesh(mesh,volFileName); }
Далее дадим подробные комментарии к тексту программы:
- Координаты точек, являющихся границами ребер, находятся в массиве points[POINTS_NUM]. В массиве edges[EDGES_NUM] заданы ребра региона поверхности в виде соответствующих пар номеров точек. Эти ребра надо будет указать NetGen, чтобы они сохранились на поверхности модели.
- Управление процессом генерации сетки осуществляется при помощи параметров, собранных в объекте ngMeshParameters класса Ng_Meshing_Parameters. Один из параметров – maxh ограничивает максимальный размер элементов сетки. С помощью других параметров можно делать тонкую настройку процесса генерации: ограничить минимальный размер элементов, инициировать оптимизацию сетки и т.д.
- STL-файл считывается функцией Ng_STL_LoadGeometry(), а возвращается объект STL-геометрии:
Ng_STL_Geometry *stlGeometry = Ng_STL_LoadGeometry(stlFileName);
где Ng_STL_Geometry следует понимать как тип указателя на объект с STL-геометрией. На самом деле Ng_STL_Geometry, как и некоторые другие типы, является синонимом void*. Авторы NetGen таким способом полностью скрыли особенности реализации сложных объектов.
- Далее в объект STL-геометрии надо добавить ребра. Они добавляются в цикле при помощи функции Ng_STL_AddEdge():
Ng_STL_AddEdge(stlGeometry, p1, p2);
где stlGeometry – объект STL-геометрии; p1 и p2 – вершины ребра. Каждая вершина задается массивом из компонент по осям X, Y, Z. Эти вершины выбираются из массива ребер edges, который был сформирован ранее.
- Выполнить вызов функции Ng_STL_InitSTLGeometry() для инициализации объекта STL-геометрии:
ng_res = Ng_STL_InitSTLGeometry(stlGeometry);
- Создать пустой объект сетки элементов, на который указывает mesh. Для этого выполняется функция Ng_NewMesh():
Ng_Mesh *mesh = Ng_NewMesh();
- Формировать ребра на поверхности модели в соответствии с STL-геометрией и параметрами генерации:
ng_res = Ng_STL_MakeEdges(stlGeometry, mesh, &ngMeshParameters);
- Выполнить триангуляцию поверхности модели функцией Ng_STL_GenerateSurfaceMesh(). Здесь уже учитываются ранее заданные ребра, которые останутся неискаженными, и будут учтены параметры генерации:
ng_res = Ng_STL_GenerateSurfaceMesh(stlGeometry, mesh, &ngMeshParameters);
- И, наконец, выполняется непосредственно генерация сетки элементов функцией Ng_GenerateVolumeMesh() с учетом параметров генерации:
ng_res = Ng_GenerateVolumeMesh(mesh, &ngMeshParameters);
- В данной программе результат генерации сетки выводится в VOL-файл посредством функции Ng_SaveMesh(), которой передается указатель на сетку и имя файла:
Ng_SaveMesh(mesh,volFileName);
Результат генерации сетки
Результатом генерации сетки является программный объект по указателю типа Ng_Mesh, в демонстрационной программе это указатель mesh. Имеется набор функций, обеспечивающих программный доступ к свойствам сетки элементов: получить количество точек в сетке Ng_GetNP(), получить количество треугольников на поверхности Ng_GetNSE(), получить количество конечных элементов Ng_GetNE(), получить координаты точки Ng_GetPoint(), получить элемент поверхности Ng_GetSurfaceElement(), получить конечный элемент Ng_GetVolumeElement(), сохранить сетку в файл Ng_SaveMesh.
VOL-файл, сохраненный в программе, может быть открыт в визуализаторе netgen.exe. Результат просмотра файла показан на рис. 2. Как видно на рисунке, границы региона поверхности остались неизменными.
Рис. 2. Поверхность модели после генерации сетки элементов
VOL-файл, формируемый NetGen, имеет следующий формат:
- В начале файла указывается код формата (mesh3d), размерность модели (dimension 3), код геометрии (geomtype 0):
mesh3d dimension 3 geomtype 0
- Далее после ключевого слова surfaceelements указывается количество и перечисляются элементы поверхности. Для каждого элемента записаны номер поверхности, номер «материала» поверхности, зарезервированное целочисленное поле, количество точек, которыми описывается элемент поверхности (для треугольников — 3), номера трех точек вершин элемента поверхности. Здесь и далее точки нумеруются от 1.
# surfnr bcnr domin domout np p1 p2 p3 surfaceelements 546 2 1 1 0 3 3 4 13 2 1 1 0 3 4 5 107 #...
- В файле перечисляются конечные элементы, для каждого из которых записаны номер «материала», количество вершин одного элемента (для тетраэдров их четыре) и номера точек, являющихся вершинами:
# matnr np p1 p2 p3 p4 volumeelements 1009 1 4 213 85 153 214 1 4 298 307 301 305 #...
- Указывается информация о ребрах на поверхности:
# surfid 0 p1 p2 trignum1 trignum2 domin/surfnr1 domout/surfnr2 ednr1 dist1 ednr2 dist2 edgesegmentsgi2 220 1 0 1 2 1 1 0 0 1 0 1 2 2 0 2 1 6 6 0 0 1 2 1 0 #...
- Перечисляются точки, являющиеся узлами сетки:
# X Y Z points 333 10.0000000000000000 4.0000000000000000 4.0000000000000000 10.0000000000000000 4.0000000000000000 2.0000000000000000 #...
- В конце файла записаны компоненты цвета для поверхностей, идентифицируемых номерами. Эти цвета используются визуализатором для раскрашивания поверхности. Здесь для поверхности с номером 2 задан цвет так, чтобы выделить прямоугольный регион на рис. 2.
# Surfnr Red Green Blue face_colours 7 2 1.00000000 1.00000000 0.00000000 3 0.00000000 1.00000000 0.00000000 #...
Формирование 2D-сетки при помощи NetGen
В NetGen имеется функция построения 2D-сетки элементов. В частном случае эта функция может быть полезна, например, для формирования STL файла для 3D-модели, имеющей плоские грани. Геометрия модели описывается во входном файле в виде границ областей пространства. Границы подобластей разбиения задаются либо с помощью отрезков прямых, либо с помощью квадратических сплайнов. При описании границы надо указывать номер области слева и справа от границы. Область вне сетки кодируется номером 0.
Рис. 3. Кодирование точек и областей для генерации 2D-сетки элементов
Приведем пример входного файла, соответствующего рис. 3.
splinecurves2dv2 2 points 1 0 0 2 3 0 3 3 2 4 0 2 5 1 0 6 2 0 7 2 1 8 1 1 segments 1 0 2 1 5 -bc=1 2 0 2 5 6 -bc=1 1 0 2 6 2 -bc=1 1 0 2 2 3 -bc=1 1 0 2 3 4 -bc=1 1 0 2 4 1 -bc=1 2 1 2 6 7 -bc=1 2 1 2 7 8 -bc=1 2 1 2 8 5 -bc=1 materials 1 domain1 -maxh=1 2 domain2 -maxh=1
В начале этого файла указывается ключевое слово splinecurves2dv2. Следом за ним коэффициент перестроения сетки (здесь – 2). После ключевого слова points перечисляются точки описания модели: номер точки и координаты по осям X, Y.
Далее после ключевого слова segments следует перечень сегментов границ в следующем формате:
- Номер области слева от границы.
- Номер области справа от границы.
- Количество точек для описания сегмента границы (в примере 2).
- Номер начальной точки границы.
- Номер конечной точки границы.
- Флаги управления генерацией. Флаг bc (номер граничного условия) приходится указывать даже тогда, когда он не используется.
После ключевого слова materials перечисляются области («материалы») в следующем формате:
- Номер области.
- Название материала.
- Флаги управления генерацией. Здесь может быть указан флаг maxh, ограничивающий максимальный размер сетки элементов.
Ниже приведена программа, генерирующая 2D-сетку:
#include "stdafx.h" namespace nglib { #include "..NetGennglibnglib.h" } using namespace nglib; void main(){ const char *in2DFileName = "..\Data\triangulation.in2d"; const char *volFileName = "..\Data\triangulation.vol"; Ng_Geometry_2D * geom; Ng_Init(); geom = Ng_LoadGeometry_2D(in2DFileName); Ng_Meshing_Parameters mp; mp.maxh = 10000; mp.fineness = 1; mp.secondorder = 0; Ng_Mesh * mesh; Ng_GenerateMesh_2D (geom, &mesh, &mp); Ng_SaveMesh(mesh, volFileName); }
Для генерации 2D-сетки используется объект описания 2D-геометрии по указателю geom. Результат генерации сетки находится в объекте по указателю mesh. В примере для вывода результата генерации используется та же самая функция Ng_SaveMesh(), что и в примере генерации 3D-сетки. Поэтому в данном случае выходной файл имеет структуру для 3D-сетки, но значение координаты Z в нем везде равно нулю. Программный интерфейс для работы с 2D-сеткой содержит функции Ng_GetNP_2D() – получить количество узлов сетки, Ng_GetNE_2D() – получить количество элементов в сетке, Ng_GetPoint_2D() – получить координаты узла сетки по номеру узла, Ng_GetElement_2D() – получить координаты вершин элемента по номеру элемента. Визуализация результата генерации 2D-сетки показана на рис. 4. Как видно из рисунка, границы региона не искажены.
Рис. 4. Результат генерации 2D-сетки
Заключение
В данной статье дано краткое представление NetGen, из которого следует, что эта библиотека удобным образом может быть использована в приложениях, требующих генерацию 3D-сетки элементов. Функции библиотеки могут быть изучены и использованы и без написания кода. В поставку NetGen входит приложение netgen.exe, которое реализует визуальный интерфейс ко всем функциям библиотеки. Рисунки в данной статье, изображающие сетки элементов, получены с помощью этого приложения.
Автор: InrecoLAN