diff --git a/src/pygmsh/__about__.py b/src/pygmsh/__about__.py
index 54e1cb45..b5f56b0b 100644
--- a/src/pygmsh/__about__.py
+++ b/src/pygmsh/__about__.py
@@ -1 +1 @@
-__version__ = "7.1.16"
+__version__ = "7.1.17"
diff --git a/src/pygmsh/common/geometry.py b/src/pygmsh/common/geometry.py
index 4fd3a7d3..c5bccee9 100644
--- a/src/pygmsh/common/geometry.py
+++ b/src/pygmsh/common/geometry.py
@@ -364,9 +364,6 @@ def generate_mesh(  # noqa: C901
         for entity in self._OUTWARD_NORMALS:
             gmsh.model.mesh.setOutwardOrientation(entity.id)
 
-        if order is not None:
-            gmsh.model.mesh.setOrder(order)
-
         gmsh.option.setNumber("General.Terminal", 1 if verbose else 0)
 
         # set algorithm
@@ -376,6 +373,11 @@ def generate_mesh(  # noqa: C901
 
         gmsh.model.mesh.generate(dim)
 
+        # setOrder() after generate(), see
+        # <https://github.com/nschloe/pygmsh/issues/515#issuecomment-1020106499>
+        if order is not None:
+            gmsh.model.mesh.setOrder(order)
+
         return extract_to_meshio()
 
     def save_geometry(self, filename: str):
diff --git a/tests/built_in/test_hole_in_square.py b/tests/built_in/test_hole_in_square.py
index 56214875..3a8c8db9 100644
--- a/tests/built_in/test_hole_in_square.py
+++ b/tests/built_in/test_hole_in_square.py
@@ -24,6 +24,9 @@ def test():
             xmin, xmax, ymin, ymax, 0.0, lcar, holes=[squareHole.curve_loop]
         )
         mesh = geom.generate_mesh(order=2)
+
+    assert "triangle6" in mesh.cells_dict
+
     # TODO support for volumes of triangle6
     # ref = 16.0
     # from helpers import compute_volume