/usr/lib/python2.7/dist-packages/dolfin/mesh/ale.py is in python-dolfin 1.3.0+dfsg-2.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 | # Copyright (C) 2009-2011 Anders Logg
#
# This file is part of DOLFIN.
#
# DOLFIN is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# DOLFIN is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
#
# First added: 2009-02-12
# Last changed: 2012-09-28
__all__ = ["compute_vertex_map", "compute_edge_map", "init_parent_edge_indices"]
from dolfin.cpp import Mesh, Vertex, error, edges, info
from dolfin import cpp
def compute_vertex_map(mesh0, mesh1):
"""
Compute map from vertices of mesh0 to vertices of mesh1.
*Arguments*
mesh0
a :py:class:`Mesh <dolfin.cpp.Mesh>`.
mesh1
a :py:class:`Mesh <dolfin.cpp.Mesh>`.
It is assumed that both meshes have a :py:class:`MeshFunction
<dolfin.cpp.MeshFunction>` over the vertices named
"parent_vertex_indices" which contain a mapping from the
local vertices to a common parent vertex numbering.
"""
# Check arguments
if not isinstance(mesh0, Mesh):
raise TypeError, "expected 'Mesh' as argument"
if not isinstance(mesh1, Mesh):
raise TypeError, "expected 'Mesh' as argument"
# Get parent vertex numbers
vertices0 = mesh0.data().array("parent_vertex_indices", 0)
vertices1 = mesh1.data().array("parent_vertex_indices", 0)
# Check mappings
if vertices0 is None or vertices1 is None:
cpp.dolfin_error("ale.py",
"compute vertex map",
"Parent vertex indices are missing")
# Compute parent-to-local mapping for mesh1
parent_to_local_mesh1 = {}
for i in range(len(vertices1)):
parent_to_local_mesh1[vertices1[i]] = i
# Compute local-to-local mapping
vertex_map = {}
for i in range(vertices0.size):
parent_vertex = vertices0[i]
if parent_vertex in parent_to_local_mesh1:
vertex_map[i] = parent_to_local_mesh1[parent_vertex]
return vertex_map
def compute_edge_map(mesh0, mesh1):
"""
Compute map from edges of mesh0 to vertices of mesh1.
*Arguments*
mesh0
a :py:class:`Mesh <dolfin.cpp.Mesh>`.
mesh1
a :py:class:`Mesh <dolfin.cpp.Mesh>`.
It is assumed that both meshes have a :py:class:`MeshFunction
<dolfin.cpp.MeshFunction>` over the vertices named
"parent_vertex_indices" which contain a mapping from the
local vertices to a common parent vertex numbering.
"""
# Check arguments
if not isinstance(mesh0, Mesh):
raise TypeError, "expected 'Mesh' as argument"
if not isinstance(mesh1, Mesh):
raise TypeError, "expected 'Mesh' as argument"
# Get parent vertex numbers
vertices0 = mesh0.data().array("parent_vertex_indices", 0)
vertices1 = mesh1.data().array("parent_vertex_indices", 0)
# Check mappings
if len(vertices0) == 0 or len(vertices1) == 0:
cpp.dolfin_error("ale.py",
"compute edge map",
"Parent vertex indices are missing")
# Initialize edges
mesh0.init(1)
mesh1.init(1)
# Build parent to local map from vertex pair to local edge for mesh0
parent_to_local_mesh0 = {}
for edge in edges(mesh0):
v = [vertices0[int(i)] for i in edge.entities(0)]
v.sort()
parent_to_local_mesh0[tuple(v)] = edge.index()
# Build parent to local map from vertex pair to local edge for mesh1
parent_to_local_mesh1 = {}
for edge in edges(mesh1):
v = [vertices1[int(i)] for i in edge.entities(0)]
v.sort()
parent_to_local_mesh1[tuple(v)] = edge.index()
# Get common edges
common_edges = set(parent_to_local_mesh0.iterkeys()).intersection(set(parent_to_local_mesh1.iterkeys()))
# Compute map
edge_map = {}
for edge in common_edges:
edge_map[parent_to_local_mesh0[edge]] = parent_to_local_mesh1[edge]
return edge_map
def init_parent_edge_indices(submesh, mesh):
"Initialize data 'parent_edge_indices' for submesh."
# Check arguments
if not isinstance(submesh, Mesh):
raise TypeError, "expected 'Mesh' as argument"
if not isinstance(mesh, Mesh):
raise TypeError, "expected 'Mesh' as argument"
# Check if edge map has already been computed
if not submesh.data().exists("parent_edge_indices", 1):
info("Edge map 'parent_edge_indices' already computed, not computing again.")
# Get parent vertex numbers
parent_vertex_indices = submesh.data().array("parent_vertex_indices", 0)
if parent_vertex_indices is None:
cpp.dolfin_error("ale.py",
"initialize parent edge indices",
"Parent vertex indice are missing")
# Make sure that we have edges for both meshes
submesh.init(1)
mesh.init(1)
# Make sure we have vertex-edge connectivity for parent mesh
mesh.init(0, 1)
# Create the edge map
parent_edge_indices = submesh.data().create_array("parent_edge_indices", 1)
# Iterate over the edges and figure out their parent number
for local_edge in edges(submesh):
# Get parent indices for edge vertices
v0, v1 = local_edge.entities(0)
V0 = Vertex(mesh, parent_vertex_indices[int(v0)])
V1 = Vertex(mesh, parent_vertex_indices[int(v1)])
# Get outgoing edges from the two parent vertices
edges0 = set(V0.entities(1))
edges1 = set(V1.entities(1))
# Check intersection
common_edges = edges0.intersection(edges1)
if not len(common_edges) == 1:
cpp.dolfin_error("ale.py",
"initialize parent edge indices",
"Parent vertices do not share exactly one common edge")
parent_edge_index = list(common_edges)[0]
# Set value
parent_edge_indices[local_edge.index()] = parent_edge_index
|