diff --git a/.gitignore b/.gitignore index feac29bf0dc2a0c0ce9af09a2085e4cd718891ce..e090b1d87d5c4921b7bd12087cccd1b157c9cf07 100644 --- a/.gitignore +++ b/.gitignore @@ -13,3 +13,5 @@ distribute*egg distribute*tar.gz a.out *.vtu + +.cache diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index bc7599d4ce9d16f916e91e1bb109ac4f426023ff..7c2f2c973ec844128cbbe291a6008a1b894344a0 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -34,3 +34,12 @@ Python 3.5 POCL: - pocl except: - tags +Documentation: + script: + - EXTRA_INSTALL="numpy" + - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/build-docs.sh + - ". ./build-docs.sh" + tags: + - python3.5 + only: + - master diff --git a/doc/index.rst b/doc/index.rst index f652d1ba9f72dffe21584f07737ebf1a7b40c4a9..4bde9b27731468d2d65b2fe35c3ee887242debea 100644 --- a/doc/index.rst +++ b/doc/index.rst @@ -8,6 +8,7 @@ Contents: mesh discretization + misc Indices and tables ================== diff --git a/doc/misc.rst b/doc/misc.rst new file mode 100644 index 0000000000000000000000000000000000000000..488262acbc29882a93c25d4ba137d77ac37e8660 --- /dev/null +++ b/doc/misc.rst @@ -0,0 +1,79 @@ +.. _installation: + +Installation +============ + +This command should install :mod:`meshmode`:: + + pip install meshmode + +(Note the extra "."!) + +You may need to run this with :command:`sudo`. +If you don't already have `pip `_, +run this beforehand:: + + curl -O https://raw.github.com/pypa/pip/master/contrib/get-pip.py + python get-pip.py + +For a more manual installation, `download the source +`_, unpack it, and say:: + + python setup.py install + +You may also clone its git repository:: + + git clone --recursive git://github.com/inducer/meshmode + git clone --recursive http://git.tiker.net/trees/meshmode.git + +User-visible Changes +==================== + +Version 2016.1 +-------------- +.. note:: + + This version is currently under development. You can get snapshots from + meshmode's `git repository `_ + +.. _license: + +Licensing +========= + +:mod:`meshmode` is licensed to you under the MIT/X Consortium license: + +Copyright (c) 2014-16 Andreas Klöckner and Contributors. + +Permission is hereby granted, free of charge, to any person +obtaining a copy of this software and associated documentation +files (the "Software"), to deal in the Software without +restriction, including without limitation the rights to use, +copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the +Software is furnished to do so, subject to the following +conditions: + +The above copyright notice and this permission notice shall be +included in all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +OTHER DEALINGS IN THE SOFTWARE. + +Acknowledgments +=============== + +Andreas Klöckner's work on :mod:`meshmode` was supported in part by + +* US Navy ONR grant number N00014-14-1-0117 +* the US National Science Foundation under grant numbers DMS-1418961 and CCF-1524433. + +AK also gratefully acknowledges a hardware gift from Nvidia Corporation. The +views and opinions expressed herein do not necessarily reflect those of the +funding agencies. diff --git a/doc/upload-docs.sh b/doc/upload-docs.sh index 25c97f52894584dfd72f3955c20026f475136e39..2c219a8c84b8f6b1835a9fb6ffb25b2274461616 100755 --- a/doc/upload-docs.sh +++ b/doc/upload-docs.sh @@ -1,3 +1,3 @@ #! /bin/sh -rsync --progress --verbose --archive --delete _build/html/* doc-upload:doc/meshmode +rsync --verbose --archive --delete _build/html/* doc-upload:doc/meshmode diff --git a/meshmode/discretization/__init__.py b/meshmode/discretization/__init__.py index e1ed150f132a516d409757567d4b5219db669235..45fdc22e4b64f93dfd8c04f34d0e5d29cbecf7a6 100644 --- a/meshmode/discretization/__init__.py +++ b/meshmode/discretization/__init__.py @@ -311,7 +311,7 @@ class Discretization(object): knl = lp.split_iname(knl, "i", 16, inner_tag="l.0") knl = lp.tag_inames(knl, dict(k="g.0")) - knl = lp.tag_data_axes(knl, "result", + knl = lp.tag_array_axes(knl, "result", "stride:auto,stride:auto,stride:auto") return knl diff --git a/meshmode/discretization/connection/refinement.py b/meshmode/discretization/connection/refinement.py index 75604ec05dd945cde45e0e6cdd55724f88083c3f..5ef91fe7e41dc464061e73438796afa174756560 100644 --- a/meshmode/discretization/connection/refinement.py +++ b/meshmode/discretization/connection/refinement.py @@ -51,6 +51,8 @@ def _map_unit_nodes_to_children(unit_nodes, tesselation): """ ref_vertices = np.array(tesselation.ref_vertices, dtype=np.float) + assert len(unit_nodes.shape) == 2 + for child_element in tesselation.children: center = np.vstack(ref_vertices[child_element[0]]) # Scale by 1/2 since sides in the tesselation have length 2. @@ -161,7 +163,9 @@ def make_refinement_connection(refiner, coarse_discr, group_factory): coarse_mesh = refiner.get_previous_mesh() fine_mesh = refiner.last_mesh - assert coarse_discr.mesh is coarse_mesh + if coarse_discr.mesh != coarse_mesh: + raise ValueError( + "course_discr must live on the same mesh given to the refiner!") from meshmode.discretization import Discretization fine_discr = Discretization( diff --git a/meshmode/discretization/poly_element.py b/meshmode/discretization/poly_element.py index a5319bc84f47ab6689f502629ff6009893033832..54aaabef97ac5f9bfc97b642ebe1a476608e0603 100644 --- a/meshmode/discretization/poly_element.py +++ b/meshmode/discretization/poly_element.py @@ -144,7 +144,13 @@ class QuadratureSimplexElementGroup(PolynomialSimplexElementGroupBase): @property @memoize_method def unit_nodes(self): - return self._quadrature_rule().nodes + result = self._quadrature_rule().nodes + if len(result.shape) == 1: + result = np.array([result]) + + dim2, nunit_nodes = result.shape + assert dim2 == self.mesh_el_group.dim + return result @property @memoize_method diff --git a/meshmode/discretization/visualization.py b/meshmode/discretization/visualization.py index 956e301533b3222db7d50916ed381222b26e416a..d68ec0c58078e68417721c8205fc3f24021d0bad 100644 --- a/meshmode/discretization/visualization.py +++ b/meshmode/discretization/visualization.py @@ -291,7 +291,7 @@ def write_nodal_adjacency_vtk_file(file_name, mesh, compressor=None,): nconnections = len(adj.neighbors) connections = np.empty((nconnections, 2), dtype=np.int32) - nb_starts = adj.neighbor_starts + nb_starts = adj.neighbors_starts for iel in range(mesh.nelements): connections[nb_starts[iel]:nb_starts[iel+1], 0] = iel diff --git a/meshmode/mesh/__init__.py b/meshmode/mesh/__init__.py index 5c2b5c1cb6ec7be972ed1ca15e0eb2db1c4b852c..da0fdfe242536997dc3dc2f71c465f2eddf07560 100644 --- a/meshmode/mesh/__init__.py +++ b/meshmode/mesh/__init__.py @@ -292,6 +292,42 @@ class SimplexElementGroup(MeshElementGroup): from modepy.tools import unit_vertices return unit_vertices(self.dim) + +class TensorProductElementGroup(MeshElementGroup): + def __init__(self, order, vertex_indices, nodes, + element_nr_base=None, node_nr_base=None, + unit_nodes=None): + """ + :arg order: the mamximum total degree used for interpolation. + :arg nodes: ``[ambient_dim, nelements, nunit_nodes]`` + The nodes are assumed to be mapped versions of *unit_nodes*. + :arg unit_nodes: ``[dim, nunit_nodes]`` + The unit nodes of which *nodes* is a mapped + version. + + Do not supply *element_nr_base* and *node_nr_base*, they will be + automatically assigned. + """ + + if not issubclass(vertex_indices.dtype.type, np.integer): + raise TypeError("vertex_indices must be integral") + + dims = unit_nodes.shape[0] + + if vertex_indices.shape[-1] != 2**dims: + raise ValueError("vertex_indices has wrong number of vertices per " + "element. expected: %d, got: %d" % (dims+1, + vertex_indices.shape[-1])) + + super(TensorProductElementGroup, self).__init__(order, vertex_indices, nodes, + element_nr_base, node_nr_base, unit_nodes) + + def face_vertex_indices(self): + raise NotImplementedError() + + def vertex_unit_coordinates(self): + raise NotImplementedError() + # }}} @@ -454,6 +490,7 @@ class Mesh(Record): face_id_dtype = np.int8 def __init__(self, vertices, groups, skip_tests=False, + node_vertex_consistency_tolerance=None, nodal_adjacency=False, facial_adjacency_groups=False, boundary_tags=None, @@ -464,6 +501,9 @@ class Mesh(Record): :arg skip_tests: Skip mesh tests, in case you want to load a broken 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. :arg nodal_adjacency: One of three options: *None*, in which case this information will be deduced from vertex adjacency. *False*, in which case @@ -540,15 +580,16 @@ class Mesh(Record): ) if not skip_tests: - assert _test_node_vertex_consistency(self) + assert _test_node_vertex_consistency( + self, node_vertex_consistency_tolerance) for g in self.groups: assert g.vertex_indices.dtype == self.vertex_id_dtype if nodal_adjacency: - assert nodal_adjacency.neighbor_starts.shape == (self.nelements+1,) + assert nodal_adjacency.neighbors_starts.shape == (self.nelements+1,) assert len(nodal_adjacency.neighbors.shape) == 1 - assert nodal_adjacency.neighbor_starts.dtype == self.element_id_dtype + assert nodal_adjacency.neighbors_starts.dtype == self.element_id_dtype assert nodal_adjacency.neighbors.dtype == self.element_id_dtype if facial_adjacency_groups: @@ -655,7 +696,7 @@ class Mesh(Record): # {{{ node-vertex consistency test -def _test_node_vertex_consistency_simplex(mesh, mgrp): +def _test_node_vertex_consistency_simplex(mesh, mgrp, tol): if mgrp.nelements == 0: return True @@ -674,7 +715,8 @@ def _test_node_vertex_consistency_simplex(mesh, mgrp): np.sum((map_vertices - grp_vertices)**2, axis=0), axis=-1)) - tol = 1e3 * np.finfo(per_element_vertex_errors.dtype).eps + if tol is None: + tol = 1e3 * np.finfo(per_element_vertex_errors.dtype).eps from meshmode.mesh.processing import find_bounding_box @@ -687,14 +729,17 @@ def _test_node_vertex_consistency_simplex(mesh, mgrp): return True -def _test_node_vertex_consistency(mesh): +def _test_node_vertex_consistency(mesh, tol): """Ensure that order of by-index vertices matches that of mapped unit vertices. """ + if tol is False: + return + for mgrp in mesh.groups: if isinstance(mgrp, SimplexElementGroup): - assert _test_node_vertex_consistency_simplex(mesh, mgrp) + assert _test_node_vertex_consistency_simplex(mesh, mgrp, tol) else: from warnings import warn warn("not implemented: node-vertex consistency check for '%s'" diff --git a/meshmode/mesh/generation.py b/meshmode/mesh/generation.py index 65761586e01c8b0aa60f82d6f8b4bc00941eeed7..ce2a5f236d0ec524a411c78c87b7904abd423db6 100644 --- a/meshmode/mesh/generation.py +++ b/meshmode/mesh/generation.py @@ -175,7 +175,10 @@ def qbx_peanut(t): # {{{ make_curve_mesh -def make_curve_mesh(curve_f, element_boundaries, order): +def make_curve_mesh(curve_f, element_boundaries, order, + unit_nodes=None, + node_vertex_consistency_tolerance=None, + return_parametrization_points=False): """ :arg curve_f: A callable representing a parametrization for a curve, accepting a vector of point locations and returning @@ -183,15 +186,20 @@ def make_curve_mesh(curve_f, element_boundaries, order): :arg element_boundaries: a vector of element boundary locations in :math:`[0,1]`, in order. 0 must be the first entry, 1 the last one. - :returns: a :class:`meshmode.mesh.Mesh` + :arg unit_nodes: if given, the unit nodes to use. Must have shape + ``(dim, nnoodes)``. + :returns: a :class:`meshmode.mesh.Mesh`, or if *return_parametrization_points* + is True, a tuple ``(mesh, par_points)``, where *par_points* is an array of + parametrization points. """ assert element_boundaries[0] == 0 assert element_boundaries[-1] == 1 nelements = len(element_boundaries) - 1 - unodes = mp.warp_and_blend_nodes(1, order) - nodes_01 = 0.5*(unodes+1) + if unit_nodes is None: + unit_nodes = mp.warp_and_blend_nodes(1, order) + nodes_01 = 0.5*(unit_nodes+1) vertices = curve_f(element_boundaries) @@ -200,7 +208,8 @@ def make_curve_mesh(curve_f, element_boundaries, order): # (el_nr, node_nr) t = el_starts[:, np.newaxis] + el_lengths[:, np.newaxis]*nodes_01 - nodes = curve_f(t.ravel()).reshape(vertices.shape[0], nelements, -1) + t = t.ravel() + nodes = curve_f(t).reshape(vertices.shape[0], nelements, -1) from meshmode.mesh import Mesh, SimplexElementGroup egroup = SimplexElementGroup( @@ -210,12 +219,18 @@ def make_curve_mesh(curve_f, element_boundaries, order): np.arange(1, nelements+1, dtype=np.int32) % nelements, ]).T, nodes=nodes, - unit_nodes=unodes) + unit_nodes=unit_nodes) - return Mesh( + mesh = Mesh( vertices=vertices, groups=[egroup], nodal_adjacency=None, - facial_adjacency_groups=None) + facial_adjacency_groups=None, + node_vertex_consistency_tolerance=node_vertex_consistency_tolerance) + + if return_parametrization_points: + return mesh, t + else: + return mesh # }}} diff --git a/meshmode/mesh/io.py b/meshmode/mesh/io.py index e27f8de352016e6f70c86be905fd4e247213d72b..f3cdaea853d0f9441d7d463f2155ba8d95004d1e 100644 --- a/meshmode/mesh/io.py +++ b/meshmode/mesh/io.py @@ -27,18 +27,23 @@ from six.moves import range, zip import numpy as np from meshpy.gmsh_reader import ( # noqa - GmshMeshReceiverBase, FileSource, LiteralSource) + GmshMeshReceiverBase, ScriptSource, FileSource, LiteralSource, + ScriptWithFilesSource, + GmshSimplexElementBase, + GmshTensorProductElementBase) __doc__ = """ +.. autoclass:: ScriptSource .. autoclass:: FileSource -.. autoclass:: LiteralSource +.. autoclass:: ScriptWithFilesSource .. autofunction:: read_gmsh .. autofunction:: generate_gmsh .. autofunction:: from_meshpy .. autofunction:: from_vertices_and_simplices +.. autofunction:: to_json """ @@ -96,6 +101,9 @@ class GmshMeshReceiver(GmshMeshReceiverBase): for el_type in self.element_types: el_type_hist[el_type] = el_type_hist.get(el_type, 0) + 1 + if not el_type_hist: + raise RuntimeError("empty mesh in gmsh input") + groups = self.groups = [] ambient_dim = self.points.shape[-1] @@ -126,7 +134,8 @@ class GmshMeshReceiver(GmshMeshReceiverBase): # }}} - from meshmode.mesh import SimplexElementGroup, Mesh + from meshmode.mesh import (Mesh, + SimplexElementGroup, TensorProductElementGroup) for group_el_type, ngroup_elements in six.iteritems(el_type_hist): if group_el_type.dimensions != mesh_bulk_dim: @@ -154,21 +163,33 @@ class GmshMeshReceiver(GmshMeshReceiverBase): unit_nodes = (np.array(group_el_type.lexicographic_node_tuples(), dtype=np.float64).T/group_el_type.order)*2 - 1 - group = SimplexElementGroup( - group_el_type.order, - vertex_indices, - nodes, - unit_nodes=unit_nodes - ) + if isinstance(group_el_type, GmshSimplexElementBase): + group = SimplexElementGroup( + group_el_type.order, + vertex_indices, + nodes, + unit_nodes=unit_nodes + ) + + if group.dim == 2: + from meshmode.mesh.processing import flip_simplex_element_group + group = flip_simplex_element_group(vertices, group, + np.ones(ngroup_elements, np.bool)) + + elif isinstance(group_el_type, GmshTensorProductElementBase): + group = TensorProductElementGroup( + group_el_type.order, + vertex_indices, + nodes, + unit_nodes=unit_nodes + ) + else: + raise NotImplementedError("gmsh element type: %s" + % type(group_el_type).__name__) # Gmsh seems to produce elements in the opposite orientation # of what we like. Flip them all. - if group.dim == 2: - from meshmode.mesh.processing import flip_simplex_element_group - group = flip_simplex_element_group(vertices, group, - np.ones(ngroup_elements, np.bool)) - groups.append(group) return Mesh( @@ -195,8 +216,9 @@ def read_gmsh(filename, force_ambient_dim=None): return recv.get_mesh() -def generate_gmsh(source, dimensions, order=None, other_options=[], - extension="geo", gmsh_executable="gmsh", 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"): """Run :command:`gmsh` on the input given by *source*, and return a :class:`meshmode.mesh.Mesh` based on the result. @@ -293,4 +315,57 @@ def from_vertices_and_simplices(vertices, simplices, order=1, fix_orientation=Fa # }}} +# {{{ to_json + +def to_json(mesh): + """Return a JSON-able Python data structure for *mesh*. The structure directly + reflects the :class:`Mesh` data structure.""" + + def btag_to_json(btag): + if isinstance(btag, str): + return btag + else: + return btag.__name__ + + def group_to_json(group): + return { + "type": type(group).__name__, + "order": group.order, + "vertex_indices": group.vertex_indices.tolist(), + "nodes": group.nodes.tolist(), + "unit_nodes": group.unit_nodes.tolist(), + "element_nr_base": group.element_nr_base, + "node_nr_base": group.node_nr_base, + "dim": group.dim, + } + + from meshmode import DataUnavailable + + def nodal_adjacency_to_json(mesh): + try: + na = mesh.nodal_adjacency + except DataUnavailable: + return None + + return { + "neighbors_starts": na.neighbors_starts.tolist(), + "neighbors": na.neighbors.tolist(), + } + + return { + "version": 0, + "vertices": mesh.vertices.tolist(), + "groups": [group_to_json(group) for group in mesh.groups], + "nodal_adjacency": nodal_adjacency_to_json(mesh), + # not yet implemented + "facial_adjacency_groups": None, + "boundary_tags": [btag_to_json(btag) for btag in mesh.boundary_tags], + "btag_to_index": dict( + (btag_to_json(btag), value) + for btag, value in six.iteritems(mesh.btag_to_index)), + } + +# }}} + + # vim: foldmethod=marker diff --git a/meshmode/mesh/refinement/__init__.py b/meshmode/mesh/refinement/__init__.py index 2154d45195412c6c012a54ace9b8fbcb07561055..16639442bdbaf10d6dd2def5c87bd8ff3254c25a 100644 --- a/meshmode/mesh/refinement/__init__.py +++ b/meshmode/mesh/refinement/__init__.py @@ -73,17 +73,19 @@ class Refiner(object): # {{{ constructor def __init__(self, mesh): - from meshmode.mesh.tesselate import tesselatetet, tesselatetri + from meshmode.mesh.tesselate import tesselateseg, tesselatetet, tesselatetri self.lazy = False self.seen_tuple = {} self.group_refinement_records = [] + seg_node_tuples, seg_result = tesselateseg() tri_node_tuples, tri_result = tesselatetri() tet_node_tuples, tet_result = tesselatetet() #quadrilateral_node_tuples = [ #print tri_result, tet_result - self.simplex_node_tuples = [None, None, tri_node_tuples, tet_node_tuples] + self.simplex_node_tuples = [ + None, seg_node_tuples, tri_node_tuples, tet_node_tuples] # Dimension-parameterized tesselations for refinement - self.simplex_result = [None, None, tri_result, tet_result] + self.simplex_result = [None, seg_result, tri_result, tet_result] #print tri_node_tuples, tri_result #self.simplex_node_tuples, self.simplex_result = tesselatetet() self.last_mesh = mesh @@ -348,6 +350,10 @@ class Refiner(object): indicating which elements should be split. """ + if len(refine_flags) != self.last_mesh.nelements: + raise ValueError("length of refine_flags does not match " + "element count of last generated mesh") + #vertices and groups for next generation nvertices = len(self.last_mesh.vertices[0]) @@ -590,8 +596,7 @@ class Refiner(object): midpoint_vertices = [] vertex_indices = grp.vertex_indices[iel_grp] #if simplex - if len(grp.vertex_indices[iel_grp]) == len(self.last_mesh.vertices)+1: - + if len(vertex_indices) == grp.dim + 1: # {{{ Get midpoints for all pairs of vertices for i in range(len(vertex_indices)): @@ -845,7 +850,7 @@ class Refiner(object): assert neighbors_starts[-1] == len(neighbors) from meshmode.mesh import NodalAdjacency - return NodalAdjacency(neighbor_starts=neighbors_starts, neighbors=neighbors) + return NodalAdjacency(neighbors_starts=neighbors_starts, neighbors=neighbors) # }}} diff --git a/meshmode/mesh/refinement/utils.py b/meshmode/mesh/refinement/utils.py index 9a78c31d6cd2dbe5bdb2f3729815b06630c96c7f..d34fa14afeeb8edc8d8e1574cdc0dfc8d4f3808d 100644 --- a/meshmode/mesh/refinement/utils.py +++ b/meshmode/mesh/refinement/utils.py @@ -64,7 +64,7 @@ def check_nodal_adj_against_geometry(mesh, tol=1e-12): for igrp, grp in enumerate(mesh.groups): for iel_grp in range(grp.nelements): iel_g = group_and_iel_to_global_iel(igrp, iel_grp) - nb_starts = nadj.neighbor_starts + nb_starts = nadj.neighbors_starts for nb_iel_g in nadj.neighbors[nb_starts[iel_g]:nb_starts[iel_g+1]]: connected_to_element_connectivity[iel_g].add(nb_iel_g) @@ -80,7 +80,7 @@ def check_nodal_adj_against_geometry(mesh, tol=1e-12): nearby_origin_vertex = mesh.vertices[ :, nearby_grp.vertex_indices[nearby_iel][0]] # noqa transformation = np.empty( - (len(mesh.vertices), nvertices_per_element-1)) + (mesh.ambient_dim, mesh.ambient_dim)) vertex_transformed = vertex - nearby_origin_vertex for inearby_vertex_index, nearby_vertex_index in enumerate( diff --git a/meshmode/mesh/tesselate.py b/meshmode/mesh/tesselate.py index 5d0e03aef4b9b2a8ac9ead1571dbdea1290ea2a7..c9527802526b4a5e6f34660a2d3200229b889ec5 100644 --- a/meshmode/mesh/tesselate.py +++ b/meshmode/mesh/tesselate.py @@ -6,6 +6,12 @@ def add_tuples(a, b): return tuple(ac+bc for ac, bc in zip(a, b)) +def tesselateseg(): + node_tuples = [(0,), (1,), (2,)] + result = [(0, 1), (1, 2)] + return [node_tuples, result] + + def tesselatetri(): result = [] diff --git a/meshmode/mesh/visualization.py b/meshmode/mesh/visualization.py index caa9c1beb3cdb7123a3c9798ce04a3c882b14f93..601086c1570e1b2a87f57b68361010774a78aeb7 100644 --- a/meshmode/mesh/visualization.py +++ b/meshmode/mesh/visualization.py @@ -1,6 +1,4 @@ -from __future__ import division -from __future__ import absolute_import -from six.moves import range +from __future__ import division, absolute_import __copyright__ = "Copyright (C) 2014 Andreas Kloeckner" @@ -25,6 +23,7 @@ THE SOFTWARE. """ import numpy as np +from six.moves import range # {{{ draw_2d_mesh diff --git a/meshmode/version.py b/meshmode/version.py index 1da70a1e09bc16f4b7e677505af6cf4dd4502b86..d26fbc2f9341a880b1119e7e6079bd51e59e11b9 100644 --- a/meshmode/version.py +++ b/meshmode/version.py @@ -1,2 +1,2 @@ -VERSION = (2015, 1) +VERSION = (2016, 1) VERSION_TEXT = ".".join(str(i) for i in VERSION) diff --git a/requirements.txt b/requirements.txt index 8e99cd7fdd25e578a5dc11c56f4c43506d18974f..6ed102e6d1309b0e65dcfec56e15aac08448a7d4 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,10 +1,11 @@ numpy -git+git://github.com/inducer/modepy -git+git://github.com/pyopencl/pyopencl -git+git://github.com/inducer/islpy -git+git://github.com/inducer/loopy +git+https://github.com/inducer/meshpy.git +git+https://github.com/inducer/modepy.git +git+https://github.com/pyopencl/pyopencl.git +git+https://github.com/inducer/islpy.git +git+https://github.com/inducer/loopy.git # required by pytential, which is in turn needed for some tests -git+git://github.com/inducer/boxtree -git+git://github.com/inducer/sumpy -git+git://github.com/inducer/pytential +git+https://github.com/inducer/boxtree.git +git+https://github.com/inducer/sumpy.git +git+https://github.com/inducer/pytential.git diff --git a/test/test_meshmode.py b/test/test_meshmode.py index 7624fbe93dee63e467f04c6f296726b12ace1e37..01e6c35cd95cd5d87270d43962d57bd25b605514 100644 --- a/test/test_meshmode.py +++ b/test/test_meshmode.py @@ -802,6 +802,29 @@ def test_nd_quad_submesh(dims): # }}} +# {{{ test_quad_mesh + +def test_quad_mesh(): + from meshmode.mesh.io import generate_gmsh, ScriptWithFilesSource + print("BEGIN GEN") + mesh = generate_gmsh( + ScriptWithFilesSource( + """ + Merge "blob-2d.step"; + Mesh.CharacteristicLengthMax = 0.05; + Recombine Surface "*" = 0.0001; + Mesh 2; + Save "output.msh"; + """, + ["blob-2d.step"]), + force_ambient_dim=2, + ) + print("END GEN") + print(mesh.nelements) + +# }}} + + if __name__ == "__main__": import sys if len(sys.argv) > 1: diff --git a/test/test_refinement.py b/test/test_refinement.py index 9c36b9750ebe87a899a8dd5c5f17f4191ed1b0b4..fca3aa9381eac1ec2e392c0e6d0e5e4f307499d3 100644 --- a/test/test_refinement.py +++ b/test/test_refinement.py @@ -31,7 +31,7 @@ from pyopencl.tools import ( # noqa pytest_generate_tests_for_pyopencl as pytest_generate_tests) from meshmode.mesh.generation import ( # noqa - generate_icosahedron, generate_box_mesh) + generate_icosahedron, generate_box_mesh, make_curve_mesh, ellipse) from meshmode.mesh.refinement.utils import check_nodal_adj_against_geometry from meshmode.mesh.refinement import Refiner @@ -80,6 +80,15 @@ def uniform_refine_flags(mesh): # partial(random_refine_flags, 0.4), # 3), + ("3_to_1_ellipse_unif", + partial( + make_curve_mesh, + partial(ellipse, 3), + np.linspace(0, 1, 20), + order=2), + uniform_refine_flags, + 4), + ("rect2d_rand", partial(generate_box_mesh, ( np.linspace(0, 1, 3), @@ -138,6 +147,7 @@ def test_refinement(case_name, mesh_gen, flag_gen, num_generations): PolynomialEquidistantGroupFactory ]) @pytest.mark.parametrize(("mesh_name", "dim", "mesh_pars"), [ + ("circle", 1, [20, 30, 40]), ("blob", 2, [1e-1, 8e-2, 5e-2]), ("warp", 2, [4, 5, 6]), ("warp", 3, [4, 5, 6]), @@ -171,7 +181,12 @@ def test_refinement_connection( for mesh_par in mesh_pars: # {{{ get mesh - if mesh_name == "blob": + if mesh_name == "circle": + assert dim == 1 + h = 1 / mesh_par + mesh = make_curve_mesh( + partial(ellipse, 1), np.linspace(0, 1, mesh_par + 1), order=1) + elif mesh_name == "blob": assert dim == 2 h = mesh_par mesh = gen_blob_mesh(h)