Python Mesh.refine_all_elements方法代码示例

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_07():
    from hermes2d.examples.c07 import set_bc, set_forms

    set_verbose(False)

    # The following parameters can be changed:
    P_INIT = 2             # Initial polynomial degree of all mesh elements.
    INIT_REF_NUM = 4       # Number of initial uniform refinements

    # Load the mesh
    mesh = Mesh()
    mesh.load(get_07_mesh())

    # Perform initial mesh refinements.
    for i in range(INIT_REF_NUM):
        mesh.refine_all_elements()

    # Create an H1 space with default shapeset
    space = H1Space(mesh, P_INIT)
    set_bc(space)

    # Initialize the weak formulation
    wf = WeakForm()
    set_forms(wf)

    # Initialize the linear system.
    ls = LinSystem(wf)
    ls.set_spaces(space)

    # Assemble and solve the matrix problem
    sln = Solution()
    ls.assemble()
    ls.solve_system(sln)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_plot_mesh3e():
    mesh = Mesh()
    mesh.load(domain_mesh)
    mesh.refine_all_elements()
    mesh.refine_all_elements()

    plot_mesh_mpl(mesh.nodes_dict, mesh.elements)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_04():
    from hermes2d.examples.c04 import set_bc

    set_verbose(False)

    # Below you can play with the parameters CONST_F, P_INIT, and UNIFORM_REF_LEVEL.
    INIT_REF_NUM = 2         # number of initial uniform mesh refinements
    P_INIT = 2               # initial polynomial degree in all elements

    # Load the mesh file
    mesh = Mesh()
    mesh.load(get_example_mesh())

    # Perform initial mesh refinements
    for i in range(INIT_REF_NUM):
        mesh.refine_all_elements()

    # Create an H1 space with default shapeset
    space = H1Space(mesh, P_INIT)
    set_bc(space)

    # Initialize the weak formulation
    wf = WeakForm()
    set_forms(wf)

    # Initialize the linear system
    ls = LinSystem(wf)
    ls.set_spaces(space)

    # Assemble and solve the matrix problem
    sln = Solution()
    ls.assemble()
    ls.solve_system(sln)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_03():
    from hermes2d.examples.c03 import set_bc

    set_verbose(False)

    P_INIT = 5                # Uniform polynomial degree of mesh elements.

    # Problem parameters.
    CONST_F = 2.0

    # Load the mesh file
    mesh = Mesh()
    mesh.load(get_example_mesh())

    # Sample "manual" mesh refinement
    mesh.refine_all_elements()

    # Create an H1 space with default shapeset
    space = H1Space(mesh, P_INIT)
    set_bc(space)

    # Initialize the weak formulation
    wf = WeakForm(1)
    set_forms(wf)

    # Initialize the linear system
    ls = LinSystem(wf)
    ls.set_spaces(space)

    # Assemble and solve the matrix problem.
    sln = Solution()
    ls.assemble()
    ls.solve_system(sln)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_08():
    from hermes2d.examples.c08 import set_bc, set_forms

    set_verbose(False)

    # The following parameter can be changed:
    P_INIT = 4

    # Load the mesh file
    mesh = Mesh()
    mesh.load(get_sample_mesh())

    # Perform uniform mesh refinement
    mesh.refine_all_elements()

    # Create the x- and y- displacement space using the default H1 shapeset
    xdisp = H1Space(mesh, P_INIT)
    ydisp = H1Space(mesh, P_INIT)
    set_bc(xdisp, ydisp)

    # Initialize the weak formulation
    wf = WeakForm(2)
    set_forms(wf)

    # Initialize the linear system.
    ls = LinSystem(wf)
    ls.set_spaces(xdisp, ydisp)

    # Assemble and solve the matrix problem
    xsln = Solution()
    ysln = Solution()
    ls.assemble()
    ls.solve_system(xsln, ysln, lib="scipy")

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_plot_mesh3c():
    mesh = Mesh()
    mesh.load(domain_mesh)
    mesh.refine_all_elements()
    mesh.refine_all_elements()

    plot_mesh_mpl_simple(mesh.nodes_dict, mesh.elements, plot_nodes=False)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_plot_mesh3d():
    mesh = Mesh()
    mesh.load(domain_mesh)
    mesh.refine_all_elements()
    mesh.refine_all_elements()

    view = MeshView("Solution")
    view.show(mesh, lib="mpl", method="orders", show=False)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_09():
    from hermes2d.examples.c09 import set_bc, temp_ext, set_forms

    # The following parameters can be changed:
    INIT_REF_NUM = 4      # number of initial uniform mesh refinements
    INIT_REF_NUM_BDY = 1  # number of initial uniform mesh refinements towards the boundary
    P_INIT = 4            # polynomial degree of all mesh elements
    TAU = 300.0           # time step in seconds

    # Problem constants
    T_INIT = 10           # temperature of the ground (also initial temperature)
    FINAL_TIME = 86400    # length of time interval (24 hours) in seconds

    # Global variable
    TIME = 0;

    # Boundary markers.
    bdy_ground = 1
    bdy_air = 2

    # Load the mesh
    mesh = Mesh()
    mesh.load(get_cathedral_mesh())

    # Perform initial mesh refinements
    for i in range(INIT_REF_NUM):
        mesh.refine_all_elements()
    mesh.refine_towards_boundary(bdy_air, INIT_REF_NUM_BDY)

    # Create an H1 space with default shapeset
    space = H1Space(mesh, P_INIT)
    set_bc(space)

    # Set initial condition
    tsln = Solution()
    tsln.set_const(mesh, T_INIT)

    # Initialize the weak formulation
    wf = WeakForm()
    set_forms(wf)

    # Initialize the linear system.
    ls = LinSystem(wf)
    ls.set_spaces(space)

    # Time stepping
    nsteps = int(FINAL_TIME/TAU + 0.5)
    rhsonly = False;

    # Assemble and solve
    ls.assemble()
    rhsonly = True
    ls.solve_system(tsln, lib="scipy")

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_02():
    set_verbose(False)    
    P_INIT = 3

    # Load the mesh file
    domain_mesh = get_example_mesh()    # Original L-shape domain
    mesh = Mesh()
    mesh.load(domain_mesh)

    # Refine all elements (optional)
    mesh.refine_all_elements()

    # Create a shapeset and an H1 space
    space = H1Space(mesh)

    # Assign element orders and initialize the space
    space.set_uniform_order(P_INIT)    # Set uniform polynomial order

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_mesh_create1():
    mesh = Mesh()
    mesh.create([
           [0, 0],
           [1, 0],
           [1, 1],
           [0, 1],
       ], [
           [2, 3, 0, 1, 0],
       ], [
           [0, 1, 1],
           [1, 2, 1],
           [2, 3, 1],
           [3, 0, 1],
       ], [])
    assert compare(mesh.nodes, [[0, 0], [1, 0], [1, 1], [0, 1]])
    assert mesh.elements_markers == [[2, 3, 0, 1, 0]]
    assert mesh.elements == [[2, 3, 0, 1]]
    # not yet implemented:
    #assert mesh.boundaries == [[0, 1, 1], [1, 2, 1], [2, 3, 1], [3, 0, 1]]
    #assert mesh.nurbs is None
    mesh.refine_all_elements()
    assert compare(mesh.nodes, [(0, 0), (1, 0), (1, 1), (0, 1),
        (1, 0.5), (0.5, 0), (0, 0.5), (0.5, 1), (0.5, 0.5)])
    assert mesh.elements == [[2, 7, 8, 4], [7, 3, 6, 8], [8, 6, 0, 5],
            [4, 8, 5, 1]]
    mesh.refine_all_elements()
    assert mesh.nodes_dict == {0: (0.0, 0.0), 1: (1.0, 0.0), 2: (1.0, 1.0), 3:
            (0.0, 1.0), 4: (1.0, 0.5), 5: (0.5, 0.0), 6: (0.0, 0.5), 7: (0.5,
                1.0), 8: (0.5, 0.5), 9: (0.5, 0.75), 10: (0.25, 1.0), 12:
            (0.75, 1.0), 13: (0.25, 0.5), 14: (0.0, 0.75), 15: (0.5, 0.25), 16:
            (0.25, 0.0), 17: (0.0, 0.25), 18: (0.75, 0.25), 19: (1.0, 0.25),
            20: (0.75, 0.0), 21: (0.75, 0.5), 22: (1.0, 0.75), 23: (0.75,
                0.75), 36: (0.25, 0.75), 47: (0.25, 0.25)}
    assert mesh.elements == [[2, 12, 23, 22], [12, 7, 9, 23], [23, 9, 8, 21],
            [22, 23, 21, 4], [7, 10, 36, 9], [10, 3, 14, 36], [36, 14, 6, 13],
            [9, 36, 13, 8], [8, 13, 47, 15], [13, 6, 17, 47], [47, 17, 0, 16],
            [15, 47, 16, 5], [4, 21, 18, 19], [21, 8, 15, 18], [18, 15, 5, 20],
            [19, 18, 20, 1]]

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_example_06():
    from hermes2d.examples.c06 import set_bc, set_forms

    set_verbose(False)

    # The following parameters can be changed:

    UNIFORM_REF_LEVEL = 2;   # Number of initial uniform mesh refinements.
    CORNER_REF_LEVEL = 12;   # Number of mesh refinements towards the re-entrant corner.
    P_INIT = 6;              # Uniform polynomial degree of all mesh elements.

    # Boundary markers
    NEWTON_BDY = 1

    # Load the mesh file
    mesh = Mesh()
    mesh.load(get_example_mesh())

    # Perform initial mesh refinements.
    for i in range(UNIFORM_REF_LEVEL):
        mesh.refine_all_elements()
    mesh.refine_towards_vertex(3, CORNER_REF_LEVEL)

    # Create an H1 space with default shapeset
    space = H1Space(mesh, P_INIT)
    set_bc(space)

    # Initialize the weak formulation
    wf = WeakForm()
    set_forms(wf)

    # Initialize the linear system.
    ls = LinSystem(wf)
    ls.set_spaces(space)

    # Assemble and solve the matrix problem
    sln = Solution()
    ls.assemble()
    ls.solve_system(sln)

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
#! /usr/bin/env python

# This example shows how to load a mesh, perform various types
# of "manual"  element refinements.

# Import modules
from hermes2d import Mesh, MeshView
from hermes2d.examples import get_example_mesh

# Load the mesh file
domain_mesh = get_example_mesh()
mesh = Mesh()
mesh.load(domain_mesh)

# Perform some sample initial refinements
mesh.refine_all_elements();           # Refines all elements.
#mesh.refine_towards_vertex(3, 4);    # Refines mesh towards vertex #3 (4x).
#mesh.refine_towards_boundary(2, 4);  # Refines all elements along boundary 2 (4x).

# Display the mesh
mesh.plot()

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def test_mesh_create2():
    m = Mesh()
    m.create([
            [0, -1],
            [1, -1],
            [-1, 0],
            [0, 0],
            [1, 0],
            [-1, 1],
            [0, 1],
            [0.707106781, 0.707106781],
        ], [
            [0, 1, 4, 3, 0],
            [3, 4, 7, 0],
            [3, 7, 6, 0],
            [2, 3, 6, 5, 0],
        ], [
            [0, 1, 1],
            [1, 4, 2],
            [3, 0, 4],
            [4, 7, 2],
            [7, 6, 2],
            [2, 3, 4],
            [6, 5, 2],
            [5, 2, 3],
        ], [
            [4, 7, 45],
            [7, 6, 45],
        ])
    assert compare(m.nodes, [
            [0, -1],
            [1, -1],
            [-1, 0],
            [0, 0],
            [1, 0],
            [-1, 1],
            [0, 1],
            [0.707106781, 0.707106781],
        ], eps=1e-4)
    assert m.elements_markers == [
            [0, 1, 4, 3, 0],
            [3, 4, 7, 0],
            [3, 7, 6, 0],
            [2, 3, 6, 5, 0],
        ]
    assert m.elements == [
            [0, 1, 4, 3],
            [3, 4, 7],
            [3, 7, 6],
            [2, 3, 6, 5],
        ]
    # This is not yet implemented:
    #assert m.boundaries == [
    #        [0, 1, 1],
    #        [1, 4, 2],
    #        [3, 0, 4],
    #        [4, 7, 2],
    #        [7, 6, 2],
    #        [2, 3, 4],
    #        [6, 5, 2],
    #        [5, 2, 3],
    #    ]
    #assert mesh.nurbs == ...
    m.refine_all_elements()
    assert m.elements == [
            [0, 11, 20, 19], [11, 1, 9, 20], [20, 9, 4, 8], [19, 20, 8, 3],
            [3, 8, 34], [8, 4, 33], [34, 33, 7], [33, 34, 8], [3, 34, 10],
            [34, 7, 12], [10, 12, 6], [12, 10, 34], [2, 18, 16, 15],
            [18, 3, 10, 16], [16, 10, 6, 17], [15, 16, 17, 5]
            ]

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
    def calc(
        threshold=0.3,
        strategy=0,
        h_only=False,
        error_tol=1,
        interactive_plotting=False,
        show_mesh=False,
        show_graph=True,
    ):
        mesh = Mesh()
        mesh.create(
            [[0, 0], [1, 0], [1, 1], [0, 1]], [[2, 3, 0, 1, 0]], [[0, 1, 1], [1, 2, 1], [2, 3, 1], [3, 0, 1]], []
        )

        mesh.refine_all_elements()

        shapeset = H1Shapeset()
        pss = PrecalcShapeset(shapeset)

        space = H1Space(mesh, shapeset)
        set_bc(space)
        space.set_uniform_order(1)

        wf = WeakForm(1)
        set_forms(wf)

        sln = Solution()
        rsln = Solution()
        solver = DummySolver()

        selector = H1ProjBasedSelector(CandList.HP_ANISO, 1.0, -1, shapeset)

        view = ScalarView("Solution")
        iter = 0
        graph = []
        while 1:
            space.assign_dofs()

            sys = LinSystem(wf, solver)
            sys.set_spaces(space)
            sys.set_pss(pss)
            sys.assemble()
            sys.solve_system(sln)
            dofs = sys.get_matrix().shape[0]
            if interactive_plotting:
                view.show(sln, lib=lib, notebook=True, filename="a%02d.png" % iter)

            rsys = RefSystem(sys)
            rsys.assemble()

            rsys.solve_system(rsln)

            hp = H1Adapt([space])
            hp.set_solutions([sln], [rsln])
            err_est = hp.calc_error() * 100

            err_est = hp.calc_error(sln, rsln) * 100
            print "iter=%02d, err_est=%5.2f%%, DOFS=%d" % (iter, err_est, dofs)
            graph.append([dofs, err_est])
            if err_est < error_tol:
                break
            hp.adapt(selector, threshold, strategy)
            iter += 1

        if not interactive_plotting:
            view.show(sln, lib=lib, notebook=True)

        if show_mesh:
            mview = MeshView("Mesh")
            mview.show(mesh, lib="mpl", notebook=True, filename="b.png")

        if show_graph:
            from numpy import array

            graph = array(graph)
            import pylab

            pylab.clf()
            pylab.plot(graph[:, 0], graph[:, 1], "ko", label="error estimate")
            pylab.plot(graph[:, 0], graph[:, 1], "k-")
            pylab.title("Error Convergence for the Inner Layer Problem")
            pylab.legend()
            pylab.xlabel("Degrees of Freedom")
            pylab.ylabel("Error [%]")
            pylab.yscale("log")
            pylab.grid()
            pylab.savefig("graph.png")

# 需要导入模块: from hermes2d import Mesh [as 别名]
# 或者: from hermes2d.Mesh import refine_all_elements [as 别名]
def schroedinger_solver(n_eigs=4, iter=2, verbose_level=1, plot=False,
        potential="hydrogen", report=False, report_filename="report.h5",
        force=False, sim_name="sim", potential2=None):
    """
    One particle Schroedinger equation solver.

    n_eigs ... the number of the lowest eigenvectors to calculate
    iter ... the number of adaptive iterations to do
    verbose_level ...
            0 ... quiet
            1 ... only moderate output (default)
            2 ... lot's of output
    plot ........ plot the progress (solutions, refined solutions, errors)
    potential ... the V(x) for which to solve, one of:
                    well, oscillator, hydrogen
    potential2 .. other terms that should be added to potential
    report ...... it will save raw data to a file, useful for creating graphs
                    etc.

    Returns the eigenvalues and eigenvectors.
    """
    set_verbose(verbose_level == 2)
    set_warn_integration(False)
    pot = {"well": 0, "oscillator": 1, "hydrogen": 2, "three-points": 3}
    pot_type = pot[potential]
    if report:
        from timeit import default_timer as clock
        from tables import IsDescription, UInt32Col, Float32Col, openFile, \
                Float64Col, Float64Atom, Col, ObjectAtom
        class Iteration(IsDescription):
            n = UInt32Col()
            DOF = UInt32Col()
            DOF_reference = UInt32Col()
            cpu_solve = Float32Col()
            cpu_solve_reference = Float32Col()
            eig_errors = Float64Col(shape=(n_eigs,))
            eigenvalues = Float64Col(shape=(n_eigs,))
            eigenvalues_reference = Float64Col(shape=(n_eigs,))
        h5file = openFile(report_filename, mode = "a",
                title = "Simulation data")
        if hasattr(h5file.root, sim_name):
            if force:
                h5file.removeNode(getattr(h5file.root, sim_name),
                        recursive=True)
            else:
                print "The group '%s' already exists. Use -f to overwrite it." \
                        % sim_name
                return
        group = h5file.createGroup("/", sim_name, 'Simulation run')
        table = h5file.createTable(group, "iterations", Iteration,
                "Iterations info")
        h5eigs = h5file.createVLArray(group, 'eigenvectors', ObjectAtom())
        h5eigs_ref = h5file.createVLArray(group, 'eigenvectors_reference',
                ObjectAtom())
        iteration = table.row

    mesh = Mesh()
    mesh.load("square.mesh")
    if potential == "well":
        # Read the width of the mesh automatically. This assumes there is just
        # one square element:
        a = sqrt(mesh.get_element(0).get_area())
        # set N high enough, so that we get enough analytical eigenvalues:
        N = 10
        levels = []
        for n1 in range(1, N):
            for n2 in range(1, N):
                levels.append(n1**2 + n2**2)
        levels.sort()

        E_exact = [pi**2/(2.*a**2) * m for m in levels]
    elif potential == "oscillator":
        E_exact = [1] + [2]*2 + [3]*3 + [4]*4 + [5]*5 + [6]*6
    elif potential == "hydrogen":
        Z = 1 # atom number
        E_exact = [-float(Z)**2/2/(n-0.5)**2/4 for n in [1]+[2]*3+[3]*5 +\
                                    [4]*8 + [5]*15]
    else:
        E_exact = [1.]*50
    if len(E_exact) < n_eigs:
        print n_eigs
        print E_exact
        raise Exception("We don't have enough analytical eigenvalues.")
    #mesh.refine_element(0)
    mesh.refine_all_elements()
    #mesh.refine_all_elements()
    #mesh.refine_all_elements()
    #mesh.refine_all_elements()

    #mview = MeshView()
    #mview.show(mesh)

    shapeset = H1Shapeset()
    space = H1Space(mesh, shapeset)
    space.set_uniform_order(2)
    space.assign_dofs()

    pss = PrecalcShapeset(shapeset)
    #bview = BaseView()
    #bview.show(space)
#.........这里部分代码省略.........