[sfepy] 04/06: Add autopkgtest.
Anton Gladky
gladk at moszumanska.debian.org
Tue May 3 20:00:49 UTC 2016
This is an automated email from the git hooks/post-receive script.
gladk pushed a commit to branch master
in repository sfepy.
commit 3667f7f466a223d776281b9324fdd4a815336015
Author: Anton Gladky <gladk at debian.org>
Date: Tue May 3 21:57:58 2016 +0200
Add autopkgtest.
---
debian/tests/build1 | 159 ++++++++++++++++++++++++
debian/tests/build2 | 232 +++++++++++++++++++++++++++++++++++
debian/tests/build3 | 334 +++++++++++++++++++++++++++++++++++++++++++++++++++
debian/tests/control | 2 +
4 files changed, 727 insertions(+)
diff --git a/debian/tests/build1 b/debian/tests/build1
new file mode 100755
index 0000000..1a8c644
--- /dev/null
+++ b/debian/tests/build1
@@ -0,0 +1,159 @@
+#!/bin/sh
+# autopkgtest check
+# (C) 2016 Anton Gladky
+
+set -e
+
+WORKDIR=$(mktemp -d)
+trap "rm -rf $WORKDIR" 0 INT QUIT ABRT PIPE TERM
+cd $WORKDIR
+
+cat <<EOF > box.py
+#!/usr/bin/env python
+"""
+Laplace equation with shifted periodic BCs.
+
+Display using::
+
+ ./postproc.py laplace_shifted_periodic.vtk --wireframe -b -d'u,plot_warp_scalar,rel_scaling=1'
+
+or use the --show option.
+"""
+import sys
+sys.path.append('.')
+from optparse import OptionParser
+import numpy as nm
+
+from sfepy.base.base import output
+from sfepy.discrete import (FieldVariable, Integral, Equation, Equations,
+ Function, Problem)
+from sfepy.discrete.fem import FEDomain, Field
+from sfepy.terms import Term
+from sfepy.discrete.conditions import (Conditions, EssentialBC,
+ LinearCombinationBC)
+from sfepy.solvers.ls import ScipyDirect
+from sfepy.solvers.nls import Newton
+from sfepy.mesh.mesh_generators import gen_block_mesh
+import sfepy.discrete.fem.periodic as per
+
+def run(domain, order):
+ omega = domain.create_region('Omega', 'all')
+ bbox = domain.get_mesh_bounding_box()
+ min_x, max_x = bbox[:, 0]
+ min_y, max_y = bbox[:, 1]
+ eps = 1e-8 * (max_x - min_x)
+ gamma1 = domain.create_region('Gamma1',
+ 'vertices in (x < %.10f)' % (min_x + eps),
+ 'facet')
+ gamma2 = domain.create_region('Gamma2',
+ 'vertices in (x > %.10f)' % (max_x - eps),
+ 'facet')
+ gamma3 = domain.create_region('Gamma3',
+ 'vertices in y < %.10f' % (min_y + eps),
+ 'facet')
+ gamma4 = domain.create_region('Gamma4',
+ 'vertices in y > %.10f' % (max_y - eps),
+ 'facet')
+
+ field = Field.from_args('fu', nm.float64, 1, omega, approx_order=order)
+
+ u = FieldVariable('u', 'unknown', field)
+ v = FieldVariable('v', 'test', field, primary_var_name='u')
+
+ integral = Integral('i', order=2*order)
+
+ t1 = Term.new('dw_laplace(v, u)',
+ integral, omega, v=v, u=u)
+ eq = Equation('eq', t1)
+ eqs = Equations([eq])
+
+ fix1 = EssentialBC('fix1', gamma1, {'u.0' : 0.4})
+ fix2 = EssentialBC('fix2', gamma2, {'u.0' : 0.0})
+
+ def get_shift(ts, coors, region):
+ return nm.ones_like(coors[:, 0])
+
+ dof_map_fun = Function('dof_map_fun', per.match_x_line)
+ shift_fun = Function('shift_fun', get_shift)
+
+ sper = LinearCombinationBC('sper', [gamma3, gamma4], {'u.0' : 'u.0'},
+ dof_map_fun, 'shifted_periodic',
+ arguments=(shift_fun,))
+
+ ls = ScipyDirect({})
+
+ pb = Problem('laplace', equations=eqs, auto_solvers=None)
+
+ pb.time_update(ebcs=Conditions([fix1, fix2]), lcbcs=Conditions([sper]))
+
+ ev = pb.get_evaluator()
+ nls = Newton({}, lin_solver=ls,
+ fun=ev.eval_residual, fun_grad=ev.eval_tangent_matrix)
+
+ pb.set_solver(nls)
+
+ state = pb.solve()
+
+ return pb, state
+
+usage = '%prog [options]\n' + __doc__.rstrip()
+
+helps = {
+ 'dims' :
+ 'dimensions of the block [default: %default]',
+ 'centre' :
+ 'centre of the block [default: %default]',
+ 'shape' :
+ 'numbers of vertices along each axis [default: %default]',
+ 'show' : 'show the results figure',
+}
+
+def main():
+ parser = OptionParser(usage=usage, version='%prog')
+ parser.add_option('-d', '--dims', metavar='dims',
+ action='store', dest='dims',
+ default='[1.0, 1.0]', help=helps['dims'])
+ parser.add_option('-c', '--centre', metavar='centre',
+ action='store', dest='centre',
+ default='[0.0, 0.0]', help=helps['centre'])
+ parser.add_option('-s', '--shape', metavar='shape',
+ action='store', dest='shape',
+ default='[11, 11]', help=helps['shape'])
+ parser.add_option('', '--show',
+ action="store_true", dest='show',
+ default=False, help=helps['show'])
+ (options, args) = parser.parse_args()
+
+ dims = nm.array(eval(options.dims), dtype=nm.float64)
+ centre = nm.array(eval(options.centre), dtype=nm.float64)
+ shape = nm.array(eval(options.shape), dtype=nm.int32)
+
+ output('dimensions:', dims)
+ output('centre: ', centre)
+ output('shape: ', shape)
+
+ mesh = gen_block_mesh(dims, shape, centre, name='block-fem')
+ fe_domain = FEDomain('domain', mesh)
+
+ pb, state = run(fe_domain, 1)
+ pb.save_state('laplace_shifted_periodic.vtk', state)
+
+ if options.show:
+ from sfepy.postprocess.viewer import Viewer
+ from sfepy.postprocess.domain_specific import DomainSpecificPlot
+
+ view = Viewer('laplace_shifted_periodic.vtk')
+ view(rel_scaling=1,
+ domain_specific={'u' : DomainSpecificPlot('plot_warp_scalar',
+ ['rel_scaling=1'])},
+ is_scalar_bar=True, is_wireframe=True,
+ opacity=0.3)
+
+if __name__ == '__main__':
+ main()
+
+EOF
+
+python box.py
+echo "run: OK"
+ls
diff --git a/debian/tests/build2 b/debian/tests/build2
new file mode 100755
index 0000000..e13bd29
--- /dev/null
+++ b/debian/tests/build2
@@ -0,0 +1,232 @@
+#!/bin/sh
+# autopkgtest check
+# (C) 2016 Anton Gladky
+
+set -e
+
+WORKDIR=$(mktemp -d)
+trap "rm -rf $WORKDIR" 0 INT QUIT ABRT PIPE TERM
+cd $WORKDIR
+
+cat <<EOF > box.py
+#!/usr/bin/env python
+from optparse import OptionParser
+import sys
+sys.path.append('.')
+
+import numpy as nm
+
+def define():
+ """Define the problem to solve."""
+ from sfepy.discrete.fem.meshio import UserMeshIO
+ from sfepy.mesh.mesh_generators import gen_block_mesh
+
+ def mesh_hook(mesh, mode):
+ """
+ Generate the block mesh.
+ """
+ if mode == 'read':
+ mesh = gen_block_mesh([2, 2, 3], [2, 2, 4], [0, 0, 1.5], name='el3',
+ verbose=False)
+ return mesh
+
+ elif mode == 'write':
+ pass
+
+ filename_mesh = UserMeshIO(mesh_hook)
+
+ options = {
+ 'nls' : 'newton',
+ 'ls' : 'ls',
+ 'ts' : 'ts',
+ 'save_steps' : -1,
+ }
+
+ functions = {
+ 'linear_tension' : (linear_tension,),
+ 'linear_compression' : (linear_compression,),
+ 'empty' : (lambda ts, coor, mode, region, ig: None,),
+ }
+
+ fields = {
+ 'displacement' : ('real', 3, 'Omega', 1),
+ }
+
+ # Coefficients are chosen so that the tangent stiffness is the same for all
+ # material for zero strains.
+ # Young modulus = 10 kPa, Poisson's ratio = 0.3
+ materials = {
+ 'solid' : ({
+ 'K' : 8.333, # bulk modulus
+ 'mu_nh' : 3.846, # shear modulus of neoHookean term
+ 'mu_mr' : 1.923, # shear modulus of Mooney-Rivlin term
+ 'kappa' : 1.923, # second modulus of Mooney-Rivlin term
+ 'lam' : 5.769, # Lame coefficients for LE term
+ 'mu_le' : 3.846,
+ },),
+ 'load' : 'empty',
+ }
+
+ variables = {
+ 'u' : ('unknown field', 'displacement', 0),
+ 'v' : ('test field', 'displacement', 'u'),
+ }
+
+ regions = {
+ 'Omega' : 'all',
+ 'Bottom' : ('vertices in (z < 0.1)', 'facet'),
+ 'Top' : ('vertices in (z > 2.9)', 'facet'),
+ }
+
+ ebcs = {
+ 'fixb' : ('Bottom', {'u.all' : 0.0}),
+ 'fixt' : ('Top', {'u.[0,1]' : 0.0}),
+ }
+
+ integrals = {
+ 'i' : 1,
+ 'isurf' : 2,
+ }
+ equations = {
+ 'linear' : """dw_lin_elastic_iso.i.Omega(solid.lam, solid.mu_le, v, u)
+ = dw_surface_ltr.isurf.Top(load.val, v)""",
+ 'neo-Hookean' : """dw_tl_he_neohook.i.Omega(solid.mu_nh, v, u)
+ + dw_tl_bulk_penalty.i.Omega(solid.K, v, u)
+ = dw_surface_ltr.isurf.Top(load.val, v)""",
+ 'Mooney-Rivlin' : """dw_tl_he_neohook.i.Omega(solid.mu_mr, v, u)
+ + dw_tl_he_mooney_rivlin.i.Omega(solid.kappa, v, u)
+ + dw_tl_bulk_penalty.i.Omega(solid.K, v, u)
+ = dw_surface_ltr.isurf.Top(load.val, v)""",
+ }
+
+ solvers = {
+ 'ls' : ('ls.scipy_direct', {}),
+ 'newton' : ('nls.newton', {
+ 'i_max' : 5,
+ 'eps_a' : 1e-10,
+ 'eps_r' : 1.0,
+ }),
+ 'ts' : ('ts.simple', {
+ 't0' : 0,
+ 't1' : 1,
+ 'dt' : None,
+ 'n_step' : 101, # has precedence over dt!
+ }),
+ }
+
+ return locals()
+
+##
+# Pressure tractions.
+def linear_tension(ts, coor, mode=None, **kwargs):
+ if mode == 'qp':
+ val = nm.tile(0.1 * ts.step, (coor.shape[0], 1, 1))
+ return {'val' : val}
+
+def linear_compression(ts, coor, mode=None, **kwargs):
+ if mode == 'qp':
+ val = nm.tile(-0.1 * ts.step, (coor.shape[0], 1, 1))
+ return {'val' : val}
+
+
+def store_top_u(displacements):
+ """Function _store() will be called at the end of each loading step. Top
+ displacements will be stored into \`displacements\`."""
+ def _store(problem, ts, state):
+
+ top = problem.domain.regions['Top']
+ top_u = problem.get_variables()['u'].get_state_in_region(top)
+ displacements.append(nm.mean(top_u[:,-1]))
+
+ return _store
+
+def solve_branch(problem, branch_function):
+ displacements = {}
+ for key, eq in problem.conf.equations.iteritems():
+ problem.set_equations({key : eq})
+
+ load = problem.get_materials()['load']
+ load.set_function(branch_function)
+
+ time_solver = problem.get_time_solver()
+ time_solver.init_time()
+
+ out = []
+ for _ in time_solver(save_results=False, step_hook=store_top_u(out)):
+ pass
+ displacements[key] = nm.array(out, dtype=nm.float64)
+
+ return displacements
+
+usage = ''
+
+helps = {
+ 'no_plot' : 'do not show plot window',
+}
+
+def main():
+ from sfepy.base.base import output
+ from sfepy.base.conf import ProblemConf, get_standard_keywords
+ from sfepy.discrete import Problem
+ from sfepy.base.plotutils import plt
+
+ parser = OptionParser(usage=usage, version='%prog')
+ parser.add_option('-n', '--no-plot',
+ action="store_true", dest='no_plot',
+ default=False, help=helps['no_plot'])
+ options, args = parser.parse_args()
+
+ required, other = get_standard_keywords()
+ # Use this file as the input file.
+ conf = ProblemConf.from_file(__file__, required, other)
+
+ # Create problem instance, but do not set equations.
+ problem = Problem.from_conf(conf, init_equations=False)
+
+ # Solve the problem. Output is ignored, results stored by using the
+ # step_hook.
+ u_t = solve_branch(problem, linear_tension)
+ u_c = solve_branch(problem, linear_compression)
+
+ # Get pressure load by calling linear_*() for each time step.
+ ts = problem.get_timestepper()
+ load_t = nm.array([linear_tension(ts, nm.array([[0.0]]), 'qp')['val']
+ for aux in ts.iter_from(0)],
+ dtype=nm.float64).squeeze()
+ load_c = nm.array([linear_compression(ts, nm.array([[0.0]]), 'qp')['val']
+ for aux in ts.iter_from(0)],
+ dtype=nm.float64).squeeze()
+
+ # Join the branches.
+ displacements = {}
+ for key in u_t.keys():
+ displacements[key] = nm.r_[u_c[key][::-1], u_t[key]]
+ load = nm.r_[load_c[::-1], load_t]
+
+
+ if plt is None:
+ output('matplotlib cannot be imported, printing raw data!')
+ output(displacements)
+ output(load)
+ else:
+ legend = []
+ for key, val in displacements.iteritems():
+ plt.plot(load, val)
+ legend.append(key)
+
+ plt.legend(legend, loc = 2)
+ plt.xlabel('tension [kPa]')
+ plt.ylabel('displacement [mm]')
+ plt.grid(True)
+
+ plt.gcf().savefig('pressure_displacement.png')
+
+if __name__ == '__main__':
+ main()
+
+
+EOF
+
+python box.py
+echo "run: OK"
+ls
diff --git a/debian/tests/build3 b/debian/tests/build3
new file mode 100755
index 0000000..118417d
--- /dev/null
+++ b/debian/tests/build3
@@ -0,0 +1,334 @@
+#!/bin/sh
+# autopkgtest check
+# (C) 2016 Anton Gladky
+
+set -e
+
+WORKDIR=$(mktemp -d)
+trap "rm -rf $WORKDIR" 0 INT QUIT ABRT PIPE TERM
+cd $WORKDIR
+
+cat <<EOF > box.py
+#!/usr/bin/env python
+"""
+Modal analysis of a linear elastic block in 2D or 3D.
+
+The dimension of the problem is determined by the length of the vector
+in \`\`--dims\`\` option.
+
+Optionally, a mesh file name can be given as a positional argument. In that
+case, the mesh generation options are ignored.
+
+The default material properties correspond to aluminium in the following units:
+
+- length: m
+- mass: kg
+- stiffness / stress: Pa
+- density: kg / m^3
+
+Examples
+--------
+
+- Run with the default arguments, show results (color = strain)::
+
+ python examples/linear_elasticity/modal_analysis.py --show
+
+- Fix bottom surface of the domain, show 9 eigen-shapes::
+
+ python examples/linear_elasticity/modal_analysis.py -b cantilever -n 9 --show
+
+- Increase mesh resolution::
+
+ python examples/linear_elasticity/modal_analysis.py -s 31,31 -n 9 --show
+
+- Use 3D domain::
+
+ python examples/linear_elasticity/modal_analysis.py -d 1,1,1 -c 0,0,0 -s 8,8,8 --show
+
+- Change the eigenvalue problem solver to LOBPCG::
+
+ python examples/linear_elasticity/modal_analysis.py --solver="eig.scipy_lobpcg,i_max:100,largest:False" --show
+
+ See :mod:\`sfepy.solvers.eigen\` for available solvers.
+"""
+import sys
+sys.path.append('.')
+from optparse import OptionParser
+
+import numpy as nm
+import scipy.sparse.linalg as sla
+
+from sfepy.base.base import assert_, output, Struct
+from sfepy.discrete import (FieldVariable, Material, Integral, Integrals,
+ Equation, Equations, Problem)
+from sfepy.discrete.fem import Mesh, FEDomain, Field
+from sfepy.terms import Term
+from sfepy.discrete.conditions import Conditions, EssentialBC
+from sfepy.mechanics.matcoefs import stiffness_from_youngpoisson
+from sfepy.mesh.mesh_generators import gen_block_mesh
+from sfepy.solvers import Solver
+
+usage = '%prog [options] [filename]\n' + __doc__.rstrip()
+
+helps = {
+ 'dims' :
+ 'dimensions of the block [default: %default]',
+ 'centre' :
+ 'centre of the block [default: %default]',
+ 'shape' :
+ 'numbers of vertices along each axis [default: %default]',
+ 'bc_kind' :
+ 'kind of Dirichlet boundary conditions on the bottom and top surfaces,'
+ ' one of: free, cantilever, fixed [default: %default]',
+ 'axis' :
+ 'the axis index of the block that the bottom and top surfaces are related'
+ ' to [default: %default]',
+ 'young' : "the Young's modulus [default: %default]",
+ 'poisson' : "the Poisson's ratio [default: %default]",
+ 'density' : "the material density [default: %default]",
+ 'order' : 'displacement field approximation order [default: %default]',
+ 'n_eigs' : 'the number of eigenvalues to compute [default: %default]',
+ 'ignore' : 'if given, the number of eigenvalues to ignore (e.g. rigid'
+ ' body modes); has precedence over the default setting determined by'
+ ' --bc-kind [default: %default]',
+ 'solver' : 'the eigenvalue problem solver to use. It should be given'
+ ' as a comma-separated list: solver_kind,option0:value0,option1:value1,...'
+ ' [default: %default]',
+ 'show' : 'show the results figure',
+}
+
+def main():
+ parser = OptionParser(usage=usage, version='%prog')
+ parser.add_option('-d', '--dims', metavar='dims',
+ action='store', dest='dims',
+ default='[1.0, 1.0]', help=helps['dims'])
+ parser.add_option('-c', '--centre', metavar='centre',
+ action='store', dest='centre',
+ default='[0.0, 0.0]', help=helps['centre'])
+ parser.add_option('-s', '--shape', metavar='shape',
+ action='store', dest='shape',
+ default='[11, 11]', help=helps['shape'])
+ parser.add_option('-b', '--bc-kind', metavar='kind',
+ action='store', dest='bc_kind',
+ choices=['free', 'cantilever', 'fixed'],
+ default='free', help=helps['bc_kind'])
+ parser.add_option('-a', '--axis', metavar='0, ..., dim, or -1', type=int,
+ action='store', dest='axis',
+ default=-1, help=helps['axis'])
+ parser.add_option('--young', metavar='float', type=float,
+ action='store', dest='young',
+ default=6.80e+10, help=helps['young'])
+ parser.add_option('--poisson', metavar='float', type=float,
+ action='store', dest='poisson',
+ default=0.36, help=helps['poisson'])
+ parser.add_option('--density', metavar='float', type=float,
+ action='store', dest='density',
+ default=2700.0, help=helps['density'])
+ parser.add_option('--order', metavar='int', type=int,
+ action='store', dest='order',
+ default=1, help=helps['order'])
+ parser.add_option('-n', '--n-eigs', metavar='int', type=int,
+ action='store', dest='n_eigs',
+ default=6, help=helps['n_eigs'])
+ parser.add_option('-i', '--ignore', metavar='int', type=int,
+ action='store', dest='ignore',
+ default=None, help=helps['ignore'])
+ parser.add_option('', '--solver', metavar='solver',
+ action='store', dest='solver',
+ default="eig.scipy,method:'eigh',tol:1e-5,maxiter:1000",
+ help=helps['solver'])
+ parser.add_option('', '--show',
+ action="store_true", dest='show',
+ default=False, help=helps['show'])
+ options, args = parser.parse_args()
+
+ aux = options.solver.split(',')
+ kwargs = {}
+ for option in aux[1:]:
+ key, val = option.split(':')
+ kwargs[key.strip()] = eval(val)
+ eig_conf = Struct(name='evp', kind=aux[0], **kwargs)
+
+ output('using values:')
+ output(" Young's modulus:", options.young)
+ output(" Poisson's ratio:", options.poisson)
+ output(' density:', options.density)
+ output('displacement field approximation order:', options.order)
+ output('requested %d eigenvalues' % options.n_eigs)
+ output('using eigenvalue problem solver:', eig_conf.kind)
+ output.level += 1
+ for key, val in kwargs.iteritems():
+ output('%s: %r' % (key, val))
+ output.level -= 1
+
+ assert_((0.0 < options.poisson < 0.5),
+ "Poisson's ratio must be in ]0, 0.5[!")
+ assert_((0 < options.order),
+ 'displacement approximation order must be at least 1!')
+
+ if len(args) == 1:
+ filename = args[0]
+
+ mesh = Mesh.from_file(filename)
+ dim = mesh.dim
+ dims = nm.diff(mesh.get_bounding_box(), axis=0)
+
+ else:
+ dims = nm.array(eval(options.dims), dtype=nm.float64)
+ dim = len(dims)
+
+ centre = nm.array(eval(options.centre), dtype=nm.float64)[:dim]
+ shape = nm.array(eval(options.shape), dtype=nm.int32)[:dim]
+
+ output('dimensions:', dims)
+ output('centre: ', centre)
+ output('shape: ', shape)
+
+ mesh = gen_block_mesh(dims, shape, centre, name='mesh')
+
+ output('axis: ', options.axis)
+ assert_((-dim <= options.axis < dim), 'invalid axis value!')
+
+ eig_solver = Solver.any_from_conf(eig_conf)
+
+ # Build the problem definition.
+ domain = FEDomain('domain', mesh)
+
+ bbox = domain.get_mesh_bounding_box()
+ min_coor, max_coor = bbox[:, options.axis]
+ eps = 1e-8 * (max_coor - min_coor)
+ ax = 'xyz'[:dim][options.axis]
+
+ omega = domain.create_region('Omega', 'all')
+ bottom = domain.create_region('Bottom',
+ 'vertices in (%s < %.10f)'
+ % (ax, min_coor + eps),
+ 'facet')
+ bottom_top = domain.create_region('BottomTop',
+ 'r.Bottom +v vertices in (%s > %.10f)'
+ % (ax, max_coor - eps),
+ 'facet')
+
+ field = Field.from_args('fu', nm.float64, 'vector', omega,
+ approx_order=options.order)
+
+ u = FieldVariable('u', 'unknown', field)
+ v = FieldVariable('v', 'test', field, primary_var_name='u')
+
+ mtx_d = stiffness_from_youngpoisson(dim, options.young, options.poisson)
+
+ m = Material('m', D=mtx_d, rho=options.density)
+
+ integral = Integral('i', order=2*options.order)
+
+ t1 = Term.new('dw_lin_elastic(m.D, v, u)', integral, omega, m=m, v=v, u=u)
+ t2 = Term.new('dw_volume_dot(m.rho, v, u)', integral, omega, m=m, v=v, u=u)
+ eq1 = Equation('stiffness', t1)
+ eq2 = Equation('mass', t2)
+ lhs_eqs = Equations([eq1, eq2])
+
+ pb = Problem('modal', equations=lhs_eqs)
+
+ if options.bc_kind == 'free':
+ pb.time_update()
+ n_rbm = dim * (dim + 1) / 2
+
+ elif options.bc_kind == 'cantilever':
+ fixed = EssentialBC('Fixed', bottom, {'u.all' : 0.0})
+ pb.time_update(ebcs=Conditions([fixed]))
+ n_rbm = 0
+
+ elif options.bc_kind == 'fixed':
+ fixed = EssentialBC('Fixed', bottom_top, {'u.all' : 0.0})
+ pb.time_update(ebcs=Conditions([fixed]))
+ n_rbm = 0
+
+ else:
+ raise ValueError('unsupported BC kind! (%s)' % options.bc_kind)
+
+ if options.ignore is not None:
+ n_rbm = options.ignore
+
+ pb.update_materials()
+
+ # Assemble stiffness and mass matrices.
+ mtx_k = eq1.evaluate(mode='weak', dw_mode='matrix', asm_obj=pb.mtx_a)
+ mtx_m = mtx_k.copy()
+ mtx_m.data[:] = 0.0
+ mtx_m = eq2.evaluate(mode='weak', dw_mode='matrix', asm_obj=mtx_m)
+
+ try:
+ eigs, svecs = eig_solver(mtx_k, mtx_m, options.n_eigs + n_rbm,
+ eigenvectors=True)
+
+ except sla.ArpackNoConvergence as ee:
+ eigs = ee.eigenvalues
+ svecs = ee.eigenvectors
+ output('only %d eigenvalues converged!' % len(eigs))
+
+ output('%d eigenvalues converged (%d ignored as rigid body modes)' %
+ (len(eigs), n_rbm))
+
+ eigs = eigs[n_rbm:]
+ svecs = svecs[:, n_rbm:]
+
+ omegas = nm.sqrt(eigs)
+ freqs = omegas / (2 * nm.pi)
+
+ output('number | eigenvalue | angular frequency '
+ '| frequency')
+ for ii, eig in enumerate(eigs):
+ output('%6d | %17.12e | %17.12e | %17.12e'
+ % (ii + 1, eig, omegas[ii], freqs[ii]))
+
+ # Make full eigenvectors (add DOFs fixed by boundary conditions).
+ variables = pb.get_variables()
+
+ vecs = nm.empty((variables.di.ptr[-1], svecs.shape[1]),
+ dtype=nm.float64)
+ for ii in xrange(svecs.shape[1]):
+ vecs[:, ii] = variables.make_full_vec(svecs[:, ii])
+
+ # Save the eigenvectors.
+ out = {}
+ state = pb.create_state()
+ for ii in xrange(eigs.shape[0]):
+ state.set_full(vecs[:, ii])
+ aux = state.create_output_dict()
+ strain = pb.evaluate('ev_cauchy_strain.i.Omega(u)',
+ integrals=Integrals([integral]),
+ mode='el_avg', verbose=False)
+ out['u%03d' % ii] = aux.popitem()[1]
+ out['strain%03d' % ii] = Struct(mode='cell', data=strain)
+
+ pb.save_state('eigenshapes.vtk', out=out)
+ pb.save_regions_as_groups('regions')
+
+ if len(eigs) and options.show:
+ # Show the solution. If the approximation order is greater than 1, the
+ # extra DOFs are simply thrown away.
+ from sfepy.postprocess.viewer import Viewer
+ from sfepy.postprocess.domain_specific import DomainSpecificPlot
+
+ scaling = 0.05 * dims.max() / nm.abs(vecs).max()
+
+ ds = {}
+ for ii in xrange(eigs.shape[0]):
+ pd = DomainSpecificPlot('plot_displacements',
+ ['rel_scaling=%s' % scaling,
+ 'color_kind="tensors"',
+ 'color_name="strain%03d"' % ii])
+ ds['u%03d' % ii] = pd
+
+ view = Viewer('eigenshapes.vtk')
+ view(domain_specific=ds, only_names=sorted(ds.keys()),
+ is_scalar_bar=False, is_wireframe=True)
+
+if __name__ == '__main__':
+ main()
+
+EOF
+
+python box.py
+echo "run: OK"
+ls
diff --git a/debian/tests/control b/debian/tests/control
new file mode 100644
index 0000000..238d0d1
--- /dev/null
+++ b/debian/tests/control
@@ -0,0 +1,2 @@
+Tests: build1 build2 build3
+Depends: python-sfepy
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/sfepy.git
More information about the debian-science-commits
mailing list