Python Mesh.order方法代码示例

# 需要导入模块: from dolfin import Mesh [as 别名]
# 或者: from dolfin.Mesh import order [as 别名]
def convert_and_create_facet_mesh_function ( ifilename, ofilename ):
    # First convert the gmsh mesh
    meshconvert.convert2xml ( ifilename, ofilename )
    
    # Now load the created mesh and initialise the required connectivity information
    mesh = Mesh ( ofilename )
    mesh.order()
    
    File ( ofilename ) << mesh
    
    D = mesh.topology().dim()
    mesh.init(D-1, 0)
    
    # read the data from the gmsh file once again
    dim_count, vertices_used, tags = process_gmsh_elements( ifilename, D-1 )
    # Get the facet-node connectivity information (reshape as a row of node indices per facet)
    facets_as_nodes = mesh.topology()(D-1,0)().reshape ( mesh.num_facets(), D )
    
    # Create and initialise the mesh function
    facet_mark_function = MeshFunction ( 'uint', mesh, D-1 )
    facet_mark_function.set_all( 0 )
    
    # set the relevant values of the mesh function
    facets_to_check = range( mesh.num_facets() )
    for i in range(len(tags)):
        nodes = np.sort(np.array(vertices_used[2*i:(2*i+D)]))
        value  = tags[i][0]
        
        if value != 0:
            found = False
            for j in range(len(facets_to_check)):
                index = facets_to_check[j]
                if np.array_equal(facets_as_nodes[index,:], nodes):
                    found = True;
                    facets_to_check.pop(j)
                    # set the value of the mesh function
                    facet_mark_function[index] = value
                    break;
                
            if not found:
                raise Exception ( "The facet (%d) was not found to mark: %s" % (i, nodes) )
        
    # save the mesh function to file
    fname = os.path.splitext(ofilename)[0]
    mesh_function_file = File("%s_%s.xml" % (fname, "facet_regions"))
    
    mesh_function_file << facet_mark_function

# 需要导入模块: from dolfin import Mesh [as 别名]
# 或者: from dolfin.Mesh import order [as 别名]
def polyhedron_mesh(data):
    """
    Build polyhedral mesh. Must be strlike with respect to the origin.

    Input:
        data[0] - list of vertices
        data[1] - list of faces
        data[2] - optional other starlike point, instead of the origin
    """
    # TODO: Center of mass of the vertices as origin
    vertex_data = np.array(data[0], dtype='double')
    lowest = np.min(flatten(data[1]))
    face_data = [list(np.array(d) - lowest) for d in data[1]]
    numv = len(vertex_data)  # will be the index of the origin
    if len(data) > 2:
        origin = np.array(data[2], dtype='double')
    else:
        origin = [0.0, 0.0, 0.0]
    mesh = Mesh()
    editor = DynamicMeshEditor()
    editor.open(mesh, "tetrahedron", 3, 3, numv + 1, len(face_data))
    for i, vert in enumerate(vertex_data):
        editor.add_vertex(i, *vert)
    editor.add_vertex(numv, *origin)
    newv = numv + 1  # next vertex index
    newf = 0  # next face index
    for face in face_data:
        if len(face) == 3:
            # triangular face, no splitting
            editor.add_cell(newf, numv, *face)  # face + origin
            newf += 1
        else:
            # split face into triangles using center of mass
            # average face vertices to get the center
            vert = list(np.mean(vertex_data[np.array(face)], axis=0))
            editor.add_vertex(newv, *vert)  # new vertex: face center
            face.append(face[0])
            for i in zip(face[:-1], face[1:]):
                # pairs of vertices
                editor.add_cell(newf, numv, newv, *i)  # + face center + origin
                newf += 1
            newv += 1
    editor.close()
    mesh.order()
    return mesh

# 需要导入模块: from dolfin import Mesh [as 别名]
# 或者: from dolfin.Mesh import order [as 别名]
def write_fenics_file(dim, ofilename):
    ofile  = File(ofilename + '.xml')
    mesh = Mesh()
    editor = MeshEditor()
    editor.open(mesh, dim, dim)
    editor.init_vertices(nodes.shape[1])
    editor.init_cells(len(cell_map))    
    for i in range(nodes.shape[1]):
        if dim == 2:
            editor.add_vertex(i, nodes[0, i], nodes[1, i])
        else:
            editor.add_vertex(i, nodes[0, i], nodes[1, i], nodes[2, i])
            
    for i in range(1, len(cell_map)+1):
        if dim == 2:
            editor.add_cell(i-1, cell_map[i][0]-1, cell_map[i][1]-1, cell_map[i][2]-1)
        else:
            editor.add_cell(i-1, cell_map[i][0]-1, cell_map[i][1]-1, cell_map[i][2]-1, cell_map[i][3]-1)
    
    mesh.order()
    mvc = mesh.domains().markers(dim-1)
    for zone, cells in boundary_cells.iteritems():
        for cell, nds in cells.iteritems():
            dolfin_cell = Cell(mesh, cell-1)
            nodes_of_cell = dolfin_cell.entities(0)
            #print cell
            #print nodes_of_cell
            nodes_of_face = nds - 1
            #print nodes_of_face
            for jj, ff in enumerate(facets(dolfin_cell)):
                facet_nodes = ff.entities(0)
                #print facet_nodes
                if all(map(lambda x: x in nodes_of_face, facet_nodes)):
                    local_index = jj
                    break
            mvc.set_value(cell-1, local_index, zone)
        
    ofile << mesh        
    from dolfin import plot
    plot(mesh, interactive=True)
    print 'Finished writing FEniCS mesh\n'

# 需要导入模块: from dolfin import Mesh [as 别名]
# 或者: from dolfin.Mesh import order [as 别名]
# Rossby radius.
LR=c/params["f"]

class InitialConditions(Expression):
    def __init__(self):
        pass
    def eval(self, values, X):
        r=(X[0]**2+X[1]**2)**0.5
        if r>0.0001:
            values[0]=-0.05*c*exp((r-r0)/LR)*X[0]/r*X[1]/r
            values[1]= 0.05*c*exp((r-r0)/LR)*X[0]/r*X[0]/r
            values[2]= 0.05*exp((r-r0)/LR)*X[0]/r
        else:
            values[0]=0.
            values[1]=0.
            values[2]=0.
    def value_shape(self):
        return (3,)

try:
    mesh=Mesh("basin.xml")
except RuntimeError:
    import sys
    import os.path

    mesh=Mesh(os.path.dirname(sys.argv[0]) + os.path.sep + "basin.xml")

mesh.order()
mesh.init()