diff --git a/meshmode/mesh/__init__.py b/meshmode/mesh/__init__.py
index 89e97a137fd84b1ed291d17f80192afb6854b09e..b8935ca2e54574732f3737baee4e4b158aeb8ce0 100644
--- a/meshmode/mesh/__init__.py
+++ b/meshmode/mesh/__init__.py
@@ -537,6 +537,7 @@ class Mesh(Record):
 
     def __init__(self, vertices, groups, skip_tests=False,
             node_vertex_consistency_tolerance=None,
+            skip_element_orientation_test=False,
             nodal_adjacency=None,
             facial_adjacency_groups=None,
             boundary_tags=None,
@@ -550,7 +551,10 @@ class Mesh(Record):
             mesh anyhow and then fix it inside of this data structure.
         :arg node_vertex_consistency_tolerance: If *False*, do not check
             for consistency between vertex and nodal data. If *None*, use
-            the (small, near FP-epsilon) tolerance.
+            the (small, near FP-epsilon) default tolerance.
+        :arg skip_element_orientation_test: If *False*, check that
+            element orientation is positive in volume meshes
+            (i.e. ones where ambient and topological dimension match).
         :arg nodal_adjacency: One of three options:
             *None*, in which case this information
             will be deduced from vertex adjacency. *False*, in which case
@@ -673,7 +677,7 @@ class Mesh(Record):
             from meshmode.mesh.processing import \
                     test_volume_mesh_element_orientations
 
-            if self.dim == self.ambient_dim:
+            if self.dim == self.ambient_dim and not skip_element_orientation_test:
                 # only for volume meshes, for now
                 assert test_volume_mesh_element_orientations(self), \
                         "negatively oriented elements found"
diff --git a/meshmode/mesh/io.py b/meshmode/mesh/io.py
index c663a3eb6ac863aceb89b4f5d4ed1ea5f572a246..d19f041df65de8f9bf6715dbc93437ae7dbbe919 100644
--- a/meshmode/mesh/io.py
+++ b/meshmode/mesh/io.py
@@ -51,7 +51,7 @@ __doc__ = """
 # {{{ gmsh receiver
 
 class GmshMeshReceiver(GmshMeshReceiverBase):
-    def __init__(self):
+    def __init__(self, mesh_construction_kwargs=None):
         # Use data fields similar to meshpy.triangle.MeshInfo and
         # meshpy.tet.MeshInfo
         self.points = None
@@ -60,6 +60,11 @@ class GmshMeshReceiver(GmshMeshReceiverBase):
         self.element_markers = None
         self.tags = None
 
+        if mesh_construction_kwargs is None:
+            mesh_construction_kwargs = {}
+
+        self.mesh_construction_kwargs = mesh_construction_kwargs
+
     def set_up_nodes(self, count):
         # Preallocate array of nodes within list; treat None as sentinel value.
         # Preallocation not done for performance, but to assign values at indices
@@ -212,22 +217,25 @@ class GmshMeshReceiver(GmshMeshReceiverBase):
 
         return Mesh(
                 vertices, groups,
-                is_conforming=is_conforming)
+                is_conforming=is_conforming,
+                **self.mesh_construction_kwargs)
 
 # }}}
 
 
 # {{{ gmsh
 
-def read_gmsh(filename, force_ambient_dim=None):
+def read_gmsh(filename, force_ambient_dim=None, mesh_construction_kwargs=None):
     """Read a gmsh mesh file from *filename* and return a
     :class:`meshmode.mesh.Mesh`.
 
     :arg force_ambient_dim: if not None, truncate point coordinates to
         this many dimensions.
+    :arg mesh_construction_kwargs: *None* or a dictionary of keyword
+        arguments passed to the :class:`meshmode.mesh.Mesh` constructor.
     """
     from gmsh_interop.reader import read_gmsh
-    recv = GmshMeshReceiver()
+    recv = GmshMeshReceiver(mesh_construction_kwargs=mesh_construction_kwargs)
     read_gmsh(recv, filename, force_dimension=force_ambient_dim)
 
     return recv.get_mesh()
@@ -235,7 +243,7 @@ def read_gmsh(filename, force_ambient_dim=None):
 
 def generate_gmsh(source, dimensions=None, order=None, other_options=[],
         extension="geo", gmsh_executable="gmsh", force_ambient_dim=None,
-        output_file_name="output.msh"):
+        output_file_name="output.msh", mesh_construction_kwargs=None):
     """Run :command:`gmsh` on the input given by *source*, and return a
     :class:`meshmode.mesh.Mesh` based on the result.
 
@@ -243,8 +251,10 @@ def generate_gmsh(source, dimensions=None, order=None, other_options=[],
         :class:`LiteralSource`
     :arg force_ambient_dim: if not None, truncate point coordinates to
         this many dimensions.
+    :arg mesh_construction_kwargs: *None* or a dictionary of keyword
+        arguments passed to the :class:`meshmode.mesh.Mesh` constructor.
     """
-    recv = GmshMeshReceiver()
+    recv = GmshMeshReceiver(mesh_construction_kwargs=mesh_construction_kwargs)
 
     from gmsh_interop.runner import GmshRunner
     from gmsh_interop.reader import parse_gmsh
diff --git a/test/test_meshmode.py b/test/test_meshmode.py
index 3ee4953f2d2876c625ebc4819d08d22eb4a82ff2..7173fc29d359b9a6e7ef65f889f819e3290a455a 100644
--- a/test/test_meshmode.py
+++ b/test/test_meshmode.py
@@ -793,7 +793,10 @@ def test_sanity_balls(ctx_getter, src_file, dim, mesh_order,
         if visualize:
             vol_vis.write_vtk_file("volume-h=%g.vtu" % h, [
                 ("f", vol_one),
-                ("area_el", bind(vol_discr, sym.area_element())(queue)),
+                ("area_el", bind(
+                    vol_discr,
+                    sym.area_element(mesh.ambient_dim, mesh.ambient_dim))
+                    (queue)),
                 ])
             bdry_vis.write_vtk_file("boundary-h=%g.vtu" % h, [("f", bdry_one)])