ForkMaster/netgen
1
0
Fork 0
mirror of https://github.com/NGSolve/netgen.git synced 2024-12-24 04:50:34 +05:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'master' into memory_tracing

Browse Source
This commit is contained in:
Christopher Lackner 2020-11-19 23:05:33 +01:00
commit a394ffedef
6 changed files with 434 additions and 652 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
libsrc/geom2d/genmesh2d.cpp
View File

@ -482,7 +482,8 @@ namespace netgen
for (const Segment & seg : mesh->LineSegments())
{
dom2seg_creator.Add (seg.domin, &seg);
dom2seg_creator.Add (seg.domout, &seg);
if (seg.domin != seg.domout)
dom2seg_creator.Add (seg.domout, &seg);
}
auto dom2seg = dom2seg_creator.MoveTable();

968
libsrc/meshing/boundarylayer.cpp
View File

File diff suppressed because it is too large Load Diff

3
libsrc/meshing/boundarylayer.hpp
View File

@ -14,10 +14,11 @@ public:
// parameters by Philippose ..
Array<int> surfid;
Array<double> heights;
Array<size_t> new_matnrs;
string new_mat;
BitArray domains;
bool outside = false; // set the boundary layer on the outside
bool grow_edges = false;
Array<size_t> project_boundaries;
};
DLL_HEADER void GenerateBoundaryLayer (Mesh & mesh,

56
libsrc/meshing/python_mesh.cpp
View File

@ -839,6 +839,7 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
.def("FaceDescriptor", static_cast<FaceDescriptor&(Mesh::*)(int)> (&Mesh::GetFaceDescriptor),
py::return_value_policy::reference)
.def("GetNFaceDescriptors", &Mesh::GetNFD)
.def("GetNDomains", &Mesh::GetNDomains)
.def("GetVolumeNeighboursOfSurfaceElement", [](Mesh & self, size_t sel)
{
@ -1016,9 +1017,9 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
.def ("BoundaryLayer", [](Mesh & self, variant<string, int> boundary,
variant<double, py::list> thickness,
variant<string, py::list> material,
string material,
variant<string, int> domain, bool outside,
bool grow_edges)
optional<string> project_boundaries)
{
BoundaryLayerParameters blp;
if(int* bc = get_if<int>(&boundary); bc)
@ -1040,7 +1041,7 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
if(dom_pattern)
{
regex pattern(*dom_pattern);
if(regex_match(self.GetMaterial(fd.DomainIn()), pattern) || (fd.DomainOut() > 0 ? regex_match(self.GetMaterial(fd.DomainOut()), pattern) : false))
if((fd.DomainIn() > 0 && regex_match(self.GetMaterial(fd.DomainIn()), pattern)) || (fd.DomainOut() > 0 && regex_match(self.GetMaterial(fd.DomainOut()), pattern)))
blp.surfid.Append(i);
}
else
@ -1048,6 +1049,15 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
}
}
}
blp.new_mat = material;
if(project_boundaries.has_value())
{
regex pattern(*project_boundaries);
for(int i = 1; i<=self.GetNFD(); i++)
if(regex_match(self.GetFaceDescriptor(i).GetBCName(), pattern))
blp.project_boundaries.Append(i);
}
if(double* pthickness = get_if<double>(&thickness); pthickness)
{
@ -1060,34 +1070,13 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
blp.heights.Append(val.cast<double>());
}
auto prismlayers = blp.heights.Size();
auto first_new_mat = self.GetNDomains() + 1;
auto max_dom_nr = first_new_mat;
if(string* pmaterial = get_if<string>(&material); pmaterial)
{
self.SetMaterial(first_new_mat, *pmaterial);
for(auto i : Range(prismlayers))
blp.new_matnrs.Append(first_new_mat);
}
else
{
auto materials = *get_if<py::list>(&material);
if(py::len(materials) != prismlayers)
throw Exception("Length of thicknesses and materials must be same!");
for(auto i : Range(prismlayers))
{
self.SetMaterial(first_new_mat+i, materials[i].cast<string>());
blp.new_matnrs.Append(first_new_mat + i);
}
max_dom_nr += prismlayers-1;
}
blp.domains.SetSize(max_dom_nr + 1); // one based
int nr_domains = self.GetNDomains();
blp.domains.SetSize(nr_domains + 1); // one based
blp.domains.Clear();
if(string* pdomain = get_if<string>(&domain); pdomain)
{
regex pattern(*pdomain);
for(auto i : Range(1, first_new_mat))
for(auto i : Range(1, nr_domains+1))
if(regex_match(self.GetMaterial(i), pattern))
blp.domains.SetBit(i);
}
@ -1096,19 +1085,15 @@ DLL_HEADER void ExportNetgenMeshing(py::module &m)
auto idomain = *get_if<int>(&domain);
blp.domains.SetBit(idomain);
}
// bits for new domains must be set
if(!outside)
for(auto i : Range(first_new_mat, max_dom_nr+1))
blp.domains.SetBit(i);
blp.outside = outside;
blp.grow_edges = grow_edges;
blp.grow_edges = true;
GenerateBoundaryLayer (self, blp);
self.UpdateTopology();
}, py::arg("boundary"), py::arg("thickness"), py::arg("material"),
py::arg("domains") = ".*", py::arg("outside") = false,
py::arg("grow_edges") = false,
py::arg("project_boundaries")=nullopt,
R"delimiter(
Add boundary layer to mesh.
@ -1133,6 +1118,11 @@ outside : bool = False
grow_edges : bool = False
Grow boundary layer over edges.
project_boundaries : Optional[str] = None
Project boundarylayer to these boundaries if they meet them. Set
to boundaries that meet boundarylayer at a non-orthogonal edge and
layer-ending should be projected to that boundary.
)delimiter")
.def ("EnableTable", [] (Mesh & self, string name, bool set)

43
tests/pytest/test_boundarylayer.py
View File

@ -18,7 +18,7 @@ def test_boundarylayer(outside, capfd):
mesh = unit_cube.GenerateMesh(maxh=0.3)
ne_before = mesh.ne
layer_surfacenames = ["right", "top", "left", "back", "bottom"]
mesh.BoundaryLayer("|".join(layer_surfacenames), [0.01, 0.02], "layer", outside=outside, grow_edges=True)
mesh.BoundaryLayer("|".join(layer_surfacenames), [0.01, 0.01], "layer", outside=outside)
should_ne = ne_before + 2 * GetNSurfaceElements(mesh, layer_surfacenames)
assert mesh.ne == should_ne
@ -26,7 +26,7 @@ def test_boundarylayer(outside, capfd):
assert not "elements are not matching" in capture.out
for side in ["front"]:
mesh.BoundaryLayer(side, [0.001, 0.002], "layer", outside=outside, grow_edges=True)
mesh.BoundaryLayer(side, [0.001, 0.001], "layer", outside=outside)
should_ne += 2 * GetNSurfaceElements(mesh, [side])
assert mesh.ne == should_ne
capture = capfd.readouterr()
@ -53,7 +53,42 @@ def test_boundarylayer2(outside, version, capfd):
geo.CloseSurfaces(top, bot, [])
mesh = geo.GenerateMesh()
should_ne = mesh.ne + 2 * GetNSurfaceElements(mesh, ["default"], "part")
layersize = 0.05
mesh.BoundaryLayer("default", [0.5 * layersize, layersize], "layer", domains="part", outside=outside, grow_edges=True)
layersize = 0.025
mesh.BoundaryLayer("default", [layersize, layersize], "part", domains="part", outside=outside)
assert mesh.ne == should_ne
assert not "elements are not matching" in capfd.readouterr().out
import netgen.gui
ngs = pytest.importorskip("ngsolve")
ngs.Draw(ngs.Mesh(mesh))
mesh = ngs.Mesh(mesh)
assert ngs.Integrate(1, mesh.Materials("part")) == pytest.approx(0.5*2.05*2.05 if outside else 0.4*2*2)
assert ngs.Integrate(1, mesh) == pytest.approx(3**3)
@pytest.mark.parametrize("outside", [True, False])
def test_wrong_orientation(outside):
geo = CSGeometry()
brick = OrthoBrick((-1,0,0),(1,1,1)) - Plane((0,0,0), (1,0,0))
geo.Add(brick.mat("air"))
mesh = geo.GenerateMesh()
mesh.BoundaryLayer(".*", 0.1, "air", domains="air", outside=outside)
ngs = pytest.importorskip("ngsolve")
mesh = ngs.Mesh(mesh)
assert ngs.Integrate(1, mesh) == pytest.approx(1.2**3 if outside else 1)
def test_splitted_surface():
geo = CSGeometry()
brick = OrthoBrick((0,0,0), (1,1,1))
slots = OrthoBrick((0.2,0,-1), (0.4, 1, 2)) + OrthoBrick((0.6, 0,-1), (0.8, 1,2))
geo.Add((brick-slots).mat("block"))
geo.Add((brick*slots).mat("slot"))
mesh = geo.GenerateMesh()
mesh.BoundaryLayer(".*", [0.001, 0.001], "block", "block", outside=False)
ngs = pytest.importorskip("ngsolve")
mesh = ngs.Mesh(mesh)
assert ngs.Integrate(1, mesh) == pytest.approx(1)
assert ngs.Integrate(1, mesh.Materials("slot")) == pytest.approx(0.4)

13
tests/pytest/test_geom2d.py Normal file
View File

@ -0,0 +1,13 @@
from netgen.geom2d import SplineGeometry
def test_leftdom_equals_rightdom():
geo = SplineGeometry()
pnts = [(0,0), (1,0), (2,0), (2,1), (1,1), (0,1)]
gp = [geo.AppendPoint(*p) for p in pnts]
lines = [(0,1,0), (1,2,0), (2,3,0), (3,4,0), (4,5,0), (5,0,0), (1,4,1)]
for p1, p2, rd in lines:
geo.Append(["line", p1, p2], leftdomain=1, rightdomain=rd)
mesh = geo.GenerateMesh()