Python ProblemDefinition.save_state方法代码示例

# 需要导入模块: from sfepy.fem import ProblemDefinition [as 别名]
# 或者: from sfepy.fem.ProblemDefinition import save_state [as 别名]
    def test_solving(self):
        from sfepy.base.base import IndexedStruct
        from sfepy.fem \
             import FieldVariable, Material, ProblemDefinition, \
                    Function, Equation, Equations, Integral
        from sfepy.fem.conditions import Conditions, EssentialBC
        from sfepy.terms import Term
        from sfepy.solvers.ls import ScipyDirect
        from sfepy.solvers.nls import Newton

        u = FieldVariable('u', 'unknown', self.field, self.dim)
        v = FieldVariable('v', 'test', self.field, self.dim,
                          primary_var_name='u')

        m = Material('m', lam=1.0, mu=1.0)
        f = Material('f', val=[[0.02], [0.01]])

        bc_fun = Function('fix_u_fun', fix_u_fun,
                          extra_args={'extra_arg' : 'hello'})

        fix_u = EssentialBC('fix_u', self.gamma1, {'u.all' : bc_fun})
        shift_u = EssentialBC('shift_u', self.gamma2, {'u.0' : 0.1})

        integral = Integral('i', order=3)

        t1 = Term.new('dw_lin_elastic_iso(m.lam, m.mu, v, u)',
                      integral, self.omega, m=m, v=v, u=u)

        t2 = Term.new('dw_volume_lvf(f.val, v)', integral, self.omega, f=f, v=v)

        eq = Equation('balance', t1 + t2)
        eqs = Equations([eq])

        ls = ScipyDirect({})

        nls_status = IndexedStruct()
        nls = Newton({}, lin_solver=ls, status=nls_status)

        pb = ProblemDefinition('elasticity', equations=eqs, nls=nls, ls=ls)
        ## pb.save_regions_as_groups('regions')

        pb.time_update(ebcs=Conditions([fix_u, shift_u]))

        state = pb.solve()

        name = op.join(self.options.out_dir, 'test_high_level_solving.vtk')
        pb.save_state(name, state)

        ok = nls_status.condition == 0
        if not ok:
            self.report('solver did not converge!')

        _ok = state.has_ebc()
        if not _ok:
            self.report('EBCs violated!')

        ok = ok and _ok

        return ok

# 需要导入模块: from sfepy.fem import ProblemDefinition [as 别名]
# 或者: from sfepy.fem.ProblemDefinition import save_state [as 别名]
def main():
    from sfepy import data_dir

    parser = OptionParser(usage=usage, version="%prog")
    parser.add_option("-s", "--show", action="store_true", dest="show", default=False, help=help["show"])
    options, args = parser.parse_args()

    mesh = Mesh.from_file(data_dir + "/meshes/2d/rectangle_tri.mesh")
    domain = Domain("domain", mesh)

    min_x, max_x = domain.get_mesh_bounding_box()[:, 0]
    eps = 1e-8 * (max_x - min_x)
    omega = domain.create_region("Omega", "all")
    gamma1 = domain.create_region("Gamma1", "nodes in x < %.10f" % (min_x + eps))
    gamma2 = domain.create_region("Gamma2", "nodes in x > %.10f" % (max_x - eps))

    field = Field("fu", nm.float64, "vector", omega, space="H1", poly_space_base="lagrange", approx_order=2)

    u = FieldVariable("u", "unknown", field, mesh.dim)
    v = FieldVariable("v", "test", field, mesh.dim, primary_var_name="u")

    m = Material("m", lam=1.0, mu=1.0)
    f = Material("f", val=[[0.02], [0.01]])

    integral = Integral("i", order=3)

    t1 = Term.new("dw_lin_elastic_iso(m.lam, m.mu, v, u)", integral, omega, m=m, v=v, u=u)
    t2 = Term.new("dw_volume_lvf(f.val, v)", integral, omega, f=f, v=v)
    eq = Equation("balance", t1 + t2)
    eqs = Equations([eq])

    fix_u = EssentialBC("fix_u", gamma1, {"u.all": 0.0})

    bc_fun = Function("shift_u_fun", shift_u_fun, extra_args={"shift": 0.01})
    shift_u = EssentialBC("shift_u", gamma2, {"u.0": bc_fun})

    ls = ScipyDirect({})

    nls_status = IndexedStruct()
    nls = Newton({}, lin_solver=ls, status=nls_status)

    pb = ProblemDefinition("elasticity", equations=eqs, nls=nls, ls=ls)
    pb.save_regions_as_groups("regions")

    pb.time_update(ebcs=Conditions([fix_u, shift_u]))

    vec = pb.solve()
    print nls_status

    pb.save_state("linear_elasticity.vtk", vec)

    if options.show:
        view = Viewer("linear_elasticity.vtk")
        view(vector_mode="warp_norm", rel_scaling=2, is_scalar_bar=True, is_wireframe=True)

# 需要导入模块: from sfepy.fem import ProblemDefinition [as 别名]
# 或者: from sfepy.fem.ProblemDefinition import save_state [as 别名]
t2 = FieldVariable('t', 'parameter', field_t, 1,
                   primary_var_name='(set-to-None)')
t2.set_data(t() - T0)

term1 = Term.new('dw_lin_elastic_iso(m.lam, m.mu, v, u)',
                 integral, omega, m=m, v=v, u=u)
term2 = Term.new('dw_biot(m.alpha, v, t)',
                 integral, omega, m=m, v=v, t=t2)
eq = Equation('temperature', term1 - term2)
eqs = Equations([eq])

u_bottom = EssentialBC('u_bottom',
                       bottom, {'u.all' : 0.0})
u_top = EssentialBC('u_top',
                    top, {'u.[0]' : 0.0})

pb.set_equations_instance(eqs, keep_solvers=True)
pb.time_update(ebcs=Conditions([u_bottom, u_top]))

displacement = pb.solve()
out.update(displacement.create_output_dict())

pb.save_state('thermoelasticity.vtk', out=out)

view = Viewer('thermoelasticity.vtk')
view(vector_mode='warp_norm',
     rel_scaling=1, is_scalar_bar=True,
     is_wireframe=True,
     opacity={'wireframe' : 0.1})

# 需要导入模块: from sfepy.fem import ProblemDefinition [as 别名]
# 或者: from sfepy.fem.ProblemDefinition import save_state [as 别名]
def main():
    from sfepy import data_dir

    parser = OptionParser(usage=usage, version='%prog')
    parser.add_option('-s', '--show',
                      action="store_true", dest='show',
                      default=False, help=help['show'])
    options, args = parser.parse_args()

    mesh = Mesh.from_file(data_dir + '/meshes/2d/rectangle_tri.mesh')
    domain = Domain('domain', mesh)

    min_x, max_x = domain.get_mesh_bounding_box()[:,0]
    eps = 1e-8 * (max_x - min_x)
    omega = domain.create_region('Omega', 'all')
    gamma1 = domain.create_region('Gamma1',
                                  'nodes in x < %.10f' % (min_x + eps))
    gamma2 = domain.create_region('Gamma2',
                                  'nodes in x > %.10f' % (max_x - eps))

    field = Field('fu', nm.float64, 'vector', omega,
                  space='H1', poly_space_base='lagrange', approx_order=2)

    u = FieldVariable('u', 'unknown', field, mesh.dim)
    v = FieldVariable('v', 'test', field, mesh.dim, primary_var_name='u')

    m = Material('m', lam=1.0, mu=1.0)
    f = Material('f', val=[[0.02], [0.01]])

    integral = Integral('i', order=3)

    t1 = Term.new('dw_lin_elastic_iso(m.lam, m.mu, v, u)',
         integral, omega, m=m, v=v, u=u)
    t2 = Term.new('dw_volume_lvf(f.val, v)', integral, omega, f=f, v=v)
    eq = Equation('balance', t1 + t2)
    eqs = Equations([eq])

    fix_u = EssentialBC('fix_u', gamma1, {'u.all' : 0.0})

    bc_fun = Function('shift_u_fun', shift_u_fun, extra_args={'shift' : 0.01})
    shift_u = EssentialBC('shift_u', gamma2, {'u.0' : bc_fun})

    ls = ScipyDirect({})

    nls_status = IndexedStruct()
    nls = Newton({}, lin_solver=ls, status=nls_status)

    pb = ProblemDefinition('elasticity', equations=eqs, nls=nls, ls=ls)
    pb.save_regions_as_groups('regions')

    pb.time_update(ebcs=Conditions([fix_u, shift_u]))

    vec = pb.solve()
    print nls_status

    pb.save_state('linear_elasticity.vtk', vec)

    if options.show:
        view = Viewer('linear_elasticity.vtk')
        view(vector_mode='warp_norm', rel_scaling=2,
             is_scalar_bar=True, is_wireframe=True)