Turek

.github/workflows/test.yaml

@@ -183,12 +183,18 @@ jobs:
         with:
           name: python-package
           path: dist/
+      - name: Install Gmsh
+        # install gmsh via pip on windows and macos and via apt on linux, as
+        # the latter version is dynamically linked and requires libgl etc.
+        run: |
+          ${{ matrix.os == 'ubuntu' && 'sudo apt install gmsh' || 'python -um pip install gmsh' }}
+          echo "NUTILS_TESTING_REQUIRES=$NUTILS_TESTING_REQUIRES app:gmsh" >> $GITHUB_ENV
       - name: Install Nutils and dependencies
         id: install
         run: |
           python -um pip install --upgrade --upgrade-strategy eager wheel
           # Install Nutils from `dist` dir created in job `build-python-package`.
-          python -um pip install "$_wheel[matrix_scipy,export_mpl]"
+          python -um pip install "$_wheel[import-gmsh,matrix-scipy,export-mpl]"
       - name: Test
         run: python -um unittest discover -b -q -t . -s examples
   test-sphinx:

devtools/container/build_base.py

@@ -19,7 +19,7 @@
     container.run('apt', 'update')
     # Package `libtbb12` is required when using Intel MKL with environment
     # variable `MKL_THREADING_LAYER` set to `TBB`, which is nowadays the default.
-    container.run('apt', 'install', '-y', '--no-install-recommends', 'python3', 'python3-pip', 'python3-wheel', 'python3-ipython', 'python3-numpy', 'python3-scipy', 'python3-matplotlib', 'python3-pil', 'libmkl-rt', 'libomp-dev', 'libtbb12', 'python3-gmsh', env=dict(DEBIAN_FRONTEND='noninteractive'))
+    container.run('apt', 'install', '-y', '--no-install-recommends', 'python3', 'python3-pip', 'python3-wheel', 'python3-ipython', 'python3-numpy', 'python3-scipy', 'python3-matplotlib', 'python3-pil', 'libmkl-rt', 'libomp-dev', 'libtbb12', 'python3-gmsh', 'gmsh', env=dict(DEBIAN_FRONTEND='noninteractive'))
     container.run('apt', 'clean')
     container.add_label('org.opencontainers.image.url', 'https://github.com/evalf/nutils')
     container.add_label('org.opencontainers.image.source', 'https://github.com/evalf/nutils')

examples/turek.geo

@@ -0,0 +1,80 @@
+// Geometry file for the turek.py example.
+//
+// This is a generalized description of the setup defined by Turek and Hron,
+// which requires the following lenghts to be supplied externally, using the
+// numbers argument of mesh.gmsh or the -setnumber switch when invoking the
+// gmsh application directly:
+//
+// - channel_length: length of the fluid domain
+// - channel_height: height of the fluid domain
+// - x_center: horizontal position of the cylinder measured from the left edge
+// - y_center: vertical position of the cylinder measured from the bottom edge
+// - cylinder_radius: radius of the cylinder
+// - structure_length: length of the elastic structure measured from the cylinder wall
+// - structure_thickness: thickness of the elastic structure
+// - min_elemsize: mesh element size at the solid/fluid interface
+// - max_elemsize: mesh element size at the channel wall and far field
+//
+// The parameterization matches largely that of Table 1 of Turek and Hron 2006,
+// with the main difference that reference points A and B cannot be
+// independently placed but are always located at the tip of the elastic
+// structure and the leading edge of the cylinder, respectively.
+
+SetFactory("OpenCASCADE");
+
+Rectangle(1) = {0, 0, 0, channel_length, channel_height};
+Rectangle(2) = {x_center, y_center - structure_thickness/2, 0, cylinder_radius + structure_length, structure_thickness, 0};
+Disk(3) = {x_center, y_center, 0, cylinder_radius};
+BooleanDifference(4) = { Surface{2}; }{ Surface{3}; };
+BooleanDifference(5) = { Surface{1}; }{ Surface{2,3}; };
+A = newp; Point(A) = {x_center + cylinder_radius + structure_length, y_center, 0};
+B = newp; Point(B) = {x_center - cylinder_radius, y_center, 0};
+
+// At this point surface 3 (cylinder), 4 (solid domain) and 5 (fluid domain) are
+// non-overlapping. Gmsh promises that the boolean fragments operation with
+// deletion will reuse the surface IDs for the new objects.
+
+_() = BooleanFragments{ Surface{3,4,5}; Point{A,B}; Delete; }{};
+
+// Fragments deduplicates boundary segments, which means that we can now
+// perform boolean operations on the index sets.
+
+bnd_cylinder() = Abs(Boundary{ Surface{3}; });
+bnd_structure() = Abs(Boundary{ Surface{4}; });
+tmp = bnd_structure;
+bnd_structure -= bnd_cylinder();
+bnd_cylinder -= tmp();
+bnd_fluid() = Abs(Boundary{ Surface{5}; });
+bnd_fluid -= bnd_cylinder();
+bnd_fluid -= bnd_structure();
+
+// After subtracting the inner boundaries, only the four boundary segments of
+// rectangle 1 remain in bnd_fluid, and we are going to assume that they are
+// ordered bottom, right, top, left.
+
+Physical Surface("fluid") = {5};
+Physical Surface("solid") = {4};
+Physical Line("inlet") = {bnd_fluid(3)};
+Physical Line("outlet") = {bnd_fluid(1)};
+Physical Line("wall") = {bnd_fluid(0), bnd_fluid(2)};
+Physical Line("structure") = {bnd_structure()};
+Physical Line("cylinder") = {bnd_cylinder()};
+Physical Point("A") = {A};
+Physical Point("B") = {B};
+
+// The element size is set to be uniformly min_elemsize inside the elastic
+// structure, and grow linearly to max_elemsize in the fluid domain over a
+// distance of half the channel height.
+
+Mesh.MeshSizeFromPoints = 0;
+Mesh.MeshSizeFromCurvature = 0;
+Mesh.MeshSizeExtendFromBoundary = 0;
+Field[1] = Distance;
+Field[1].SurfacesList = {3,4};
+Field[2] = Threshold;
+Field[2].InField = 1;
+Field[2].DistMin = 0;
+Field[2].DistMax = channel_height/2;
+Field[2].SizeMin = min_elemsize;
+Field[2].SizeMax = max_elemsize;
+Background Field = 2;