Skip to content

Commit

Permalink
Merge branch '264-in-amfe-ui-change-all-solve-methods-to-take-a-mecha…
Browse files Browse the repository at this point in the history 
…nical-system-instead-of-a-component' into 'developer'

Resolve "In amfe.ui change all solve methods to take a mechanical system instead of a component"

See merge request AM/AMfe!257
  • Loading branch information
@c-meyer
c-meyer committed May 14, 2020
2 parents b34804e + 8c8c1f1 commit 2bb0b2f
Show file tree
Hide file tree
Showing 3 changed files with 139 additions and 54 deletions.
158 changes: 119 additions & 39 deletions amfe/ui.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,13 +18,9 @@
from amfe.neumann.structural_neumann import FixedDirectionNeumann, NormalFollowingNeumann
from amfe.solver.translators import create_constrained_mechanical_system_from_component
from amfe.solver import SolverFactory, AmfeSolution
from amfe.linalg.linearsolvers import ScipySparseLinearSolver

from amfe.io.mesh import GidJsonMeshReader, AmfeMeshObjMeshReader, GmshAsciiMeshReader
from amfe.io.mesh import AmfeMeshConverter
from amfe.io.postprocessing.reader import AmfeSolutionReader
from amfe.io.postprocessing.writer import Hdf5PostProcessorWriter
from amfe.io.postprocessing.tools import write_xdmf_from_hdf5
from amfe.io.postprocessing import *

import numpy as np
Expand All @@ -38,6 +34,7 @@
__all__ = ['import_mesh_from_file',
'create_structural_component',
'create_material',
'create_mechanical_system',
'assign_material_by_group',
'assign_material_by_elementids',
'set_dirichlet_by_group',
Expand Down Expand Up @@ -70,9 +67,9 @@ def import_mesh_from_file(filename):
"""
_, extension = splitext(filename)
if extension in formats:
reader = formats[extension](filename)
mreader = formats[extension](filename)
converter = AmfeMeshConverter()
reader.parse(converter)
mreader.parse(converter)
return converter.return_mesh()
else:
raise ValueError('No reader available for files with fileextension {}'.format(extension), 'Please use a mesh ',
Expand All @@ -87,8 +84,9 @@ def create_mechanical_system(structural_component, constant_mass=False,
constant_damping=False, all_linear=False, constraint_formulation='boolean',
**formulation_options):
system, formulation = create_constrained_mechanical_system_from_component(structural_component, constant_mass,
constant_damping, all_linear,
constraint_formulation, **formulation_options)
constant_damping, all_linear,
constraint_formulation,
**formulation_options)
return system, formulation


Expand Down Expand Up @@ -177,7 +175,7 @@ def set_dirichlet_by_nodeids(component, nodeids, direction=('ux'), constraint_na


def set_neumann_by_group(component, group_name, direction, following=False, neumann_name='Neumann0',
F=constant_force(1.0)):
f=constant_force(1.0)):
"""
Sets a neumann condition on a component by addressing the group in the mesh.
Expand All @@ -194,12 +192,12 @@ def set_neumann_by_group(component, group_name, direction, following=False, neum
If following is set to false (default), the direction is fixed.
neumann_name : str
name of the condition, defined by user
F : function
f : function
pointer to function with signature float func(float: t)
"""
if direction == 'normal':
if following:
neumann = NormalFollowingNeumann(time_func=F)
neumann = NormalFollowingNeumann(time_func=f)
else:
raise NotImplementedError('There is no implementation for normal direction that is fixed with inertial'
'frame')
Expand All @@ -208,19 +206,35 @@ def set_neumann_by_group(component, group_name, direction, following=False, neum
if following:
raise NotImplementedError('There is no implementation for forces that follow the body frame')
else:
neumann = FixedDirectionNeumann(direction, time_func=F)
neumann = FixedDirectionNeumann(direction, time_func=f)
component.assign_neumann(neumann_name, neumann, [group_name], '_groups')


def set_neumann_by_elementids(component, elementids, direction_vector=np.array([0, 1]), neumann_name='Neumann0',
F=constant_force(1.0)):
neumann = FixedDirectionNeumann(direction_vector, time_func=F)
f=constant_force(1.0)):
neumann = FixedDirectionNeumann(direction_vector, time_func=f)
component.assign_neumann(neumann_name, neumann, elementids, '_eleids')


def solve_linear_static(component):
system, formulation = create_mechanical_system(component, all_linear=True, constraint_formulation='boolean')
def solve_linear_static(system, formulation, component):
"""
Solves a linear, static MechanicalSystem.
Parameters
----------
system : MechanicalSystem
MechanicalSystem object describing the equilibrium forces
formulation : ConstraintFormulation
A ConstraintFormulation object that can recover the nodal unconstrained displacements of the component
from the constrained states of the system
component : StructuralComponent
Structural Component belonging to the system
Returns
-------
solution : AmfeSolution
Simple Container for solutions of AMfe simulations
"""
solfac = SolverFactory()
solfac.set_system(system)
solfac.set_analysis_type('static')
Expand All @@ -246,9 +260,33 @@ def solve_linear_static(component):
return solution_writer


def solve_linear_dynamic(component, t0, t_end, dt, write_timestep=1):
system, formulation = create_mechanical_system(component, all_linear=True, constraint_formulation='boolean')
def solve_linear_dynamic(system, formulation, component, t0, t_end, dt, write_each=1):
"""
Solves a linear, dynamic MechanicalSystem from t0 to t_end with dt as intermediate steps.
Parameters
----------
system : MechanicalSystem
MechanicalSystem object describing the equations of motion
formulation : ConstraintFormulation
A ConstraintFormulation object that can recover the nodal unconstrained displacements of the component
from the constrained states of the system
component : StructuralComponent
StructuralComponent belonging to the system
t0 : float
initial time
t_end : float
final time
dt : float
increment between time steps
write_each : int
Specifies that the solver writes every n'th solution (eg. type 2 to write every second timestep)
Returns
-------
solution : amfe.solution.AmfeSolution
AmfeSolution container that contains the solution
"""
solfac = SolverFactory()
solfac.set_system(system)
solfac.set_analysis_type('transient')
Expand All @@ -267,7 +305,7 @@ def solve_linear_dynamic(component, t0, t_end, dt, write_timestep=1):

def write_callback(t, x, dx, ddx):
cur_ts = int((t-t0) // dt)
if abs(cur_ts % write_timestep) <=1e-7:
if abs(cur_ts % write_each) <= 1e-7:
u, du, ddu = formulation.recover(x, dx, ddx, t)
solution_writer.write_timestep(t, u, du, ddu)

Expand All @@ -281,11 +319,28 @@ def write_callback(t, x, dx, ddx):
return solution_writer


def solve_nonlinear_static(component, load_steps):
def solve_nonlinear_static(system, formulation, component, load_steps):
"""
Solves a MechanicalSystem using a nonlinear static method using newton for the nonlinear solver.
The deformation is divided into an amount of intermediate steps which are defined in the number of load_steps.
system, formulation = create_mechanical_system(component, constant_mass=True, constant_damping=True,
constraint_formulation='boolean')
Parameters
----------
system : MechanicalSystem
Created MechanicalSystem object describing the Structural Component
formulation : ConstraintFormulation
A ConstraintFormulation object that can recover the nodal unconstrained displacements of the component
from the constrained states of the system
component : StructuralComponent
StructuralComponent belonging to the system
load_steps : int
Number of load steps used for the solution
Returns
-------
solution : amfe.solution.AmfeSolution
Simple Container for solutions of AMfe simulations
"""
solfac = SolverFactory()
solfac.set_system(system)
solfac.set_analysis_type('static')
Expand Down Expand Up @@ -319,10 +374,33 @@ def write_callback(t, x, dx, ddx):
return solution_writer


def solve_nonlinear_dynamic(component, t0, t_end, dt, write_timestep=1):
system, formulation = create_mechanical_system(component, constant_mass=True, constant_damping=True,
constraint_formulation='boolean')
def solve_nonlinear_dynamic(system, formulation, component, t0, t_end, dt, write_each=1):
"""
Solves a system
Parameters
----------
system : MechanicalSystem
Created MechanicalSystem object describing the Structural Component
formulation : ConstraintFormulation
A ConstraintFormulation object that can recover the nodal unconstrained displacements of the component
from the constrained states of the system
component : StructuralComponent
StructuralComponent belonging to the system
t0 : float
initial time
t_end : float
final time
dt : float
increment between time steps
write_each : int
Specifies that the solver writes every n'th solution (eg. type 2 to write every second timestep)
Returns
-------
solution : amfe.solution.AmfeSolution
Simple Container for solutions of AMfe simulations
"""
solfac = SolverFactory()
solfac.set_system(system)
solfac.set_analysis_type('transient')
Expand All @@ -341,7 +419,7 @@ def solve_nonlinear_dynamic(component, t0, t_end, dt, write_timestep=1):

def write_callback(t, x, dx, ddx):
cur_ts = int((t - t0) // dt)
if abs(cur_ts % write_timestep) <= 1e-7:
if abs(cur_ts % write_each) <= 1e-7:
u, du, ddu = formulation.recover(x, dx, ddx, t)
strains, stresses = component.strains_and_stresses(u, du, t)
solution_writer.write_timestep(t, u, du, ddu, strains, stresses)
Expand All @@ -356,7 +434,7 @@ def write_callback(t, x, dx, ddx):


def write_results_to_paraview(solution, component, paraviewfilename, displacements_only=True,
plot_strains_and_stresses=True):
write_strains_and_stresses=True):
"""
Writes results to an xdmf paraview file
Expand All @@ -372,6 +450,8 @@ def write_results_to_paraview(solution, component, paraviewfilename, displacemen
In static models, there are no velocities and accelerations, that need plotting. Therefore switch off that
functionality. This is default, because it is rarely needed to plot velocities and accelerations in dynamic
cases as well.
write_strains_and_stresses : bool
Flag if strains and stresses shall be exported to paraview
"""
paraviewfilename = splitext(paraviewfilename)[0]

Expand Down Expand Up @@ -405,23 +485,23 @@ def write_results_to_paraview(solution, component, paraviewfilename, displacemen
}
})

if plot_strains_and_stresses:
if write_strains_and_stresses:
fielddict['strains_normal'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/strains_normal'
}
fielddict['stresses_normal'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/stresses_normal'
}
fielddict['strains_shear'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/strains_shear'
'hdf5path': '/results/strains_normal'
}
fielddict['stresses_shear'] = {'mesh_entity_type': MeshEntityType.NODE,
fielddict['stresses_normal'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/stresses_shear'
'hdf5path': '/results/stresses_normal'
}
fielddict['strains_shear'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/strains_shear'
}
fielddict['stresses_shear'] = {'mesh_entity_type': MeshEntityType.NODE,
'data_type': PostProcessDataType.VECTOR,
'hdf5path': '/results/stresses_shear'
}

with open(xdmfresultsfilename, 'wb') as xdmffp:
with File(hdf5resultsfilename, mode='r') as hdf5fp:
Expand Down
1 change: 1 addition & 0 deletions docs/release/1.1.0-notes.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -89,6 +89,7 @@ amfe.tools improvements
amfe.ui improvements
--------------------

- in all solvers, changed argument from taking a Component to MechanicalSystem
- added option 'normal' to set_neumann_by_group in order to set a neumann condition with a following normal direction
to surface
- UI exports strains and stresses by default
Expand Down
34 changes: 19 additions & 15 deletions examples/ui_example_beam.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
from amfe.forces import constant_force, triangular_force, linearly_increasing_force

###################################################
### Simple linear-elastic cantilever beam ###
# Simple linear-elastic cantilever beam #
###################################################

# Cantilever beam with dimensions 10m x 1m (L x H)
Expand All @@ -22,14 +22,15 @@
ui.set_dirichlet_by_nodeids(model, [1], ('uy'), 'Dirichlet_y')

F = constant_force(5E7)
ui.set_neumann_by_group(model, 'neumann', np.array([0.0, -1.0]), neumann_name='Load', F=F)
ui.set_neumann_by_group(model, 'neumann', np.array([0.0, -1.0]), neumann_name='Load', f=F)

solution_writer = ui.solve_linear_static(model)
system, formulation = ui.create_constrained_mechanical_system_from_component(model, all_linear=True)
solution_container = ui.solve_linear_static(system, formulation, model)

ui.write_results_to_paraview(solution_writer, model, amfe_dir('results/gmsh/ui_example_beam_linear'))
ui.write_results_to_paraview(solution_container, model, amfe_dir('results/gmsh/ui_example_beam_linear'))

###########################################################
### Dynamic linear heterogeneous cantilever beam ###
# Dynamic linear heterogeneous cantilever beam #
###########################################################

model = ui.create_structural_component(mesh)
Expand All @@ -40,15 +41,16 @@
ui.set_dirichlet_by_nodeids(model, [1], ('uy'), 'Dirichlet_y')

F = triangular_force(0, 0.01, 0.02, 1E7)
ui.set_neumann_by_group(model, 'neumann', np.array([0.0, -1.0]), neumann_name='Load', F=F)
ui.set_neumann_by_group(model, 'neumann', np.array([0.0, -1.0]), neumann_name='Load', f=F)

solution_writer = ui.solve_linear_dynamic(model, 0.0, 1.0, 0.0001, 10)
system, formulation = ui.create_constrained_mechanical_system_from_component(model, all_linear=True)
solution_container = ui.solve_linear_dynamic(system, formulation, model, 0.0, 1.0, 0.0001, 10)

ui.write_results_to_paraview(solution_writer, model, amfe_dir('results/gmsh/ui_example_beam_linear_dynamic'),
ui.write_results_to_paraview(solution_container, model, amfe_dir('results/gmsh/ui_example_beam_linear_dynamic'),
displacements_only=False)

###################################################
### Nonlinear heterogeneous cantilever beam ###
# Nonlinear heterogeneous cantilever beam #
###################################################

# Cantilever beam with layers of soft and stiff material and dimensions 5m x 1m (L x H)
Expand All @@ -69,14 +71,15 @@
ui.set_dirichlet_by_group(model, 'xy_dirichlet_point', ('uy'), 'Dirichlet_y')

F = linearly_increasing_force(0, 1.00001, 1.2E5)
ui.set_neumann_by_group(model, 'z_neumann', np.array([0.0, -1.0]), neumann_name='Load', F=F)
ui.set_neumann_by_group(model, 'z_neumann', np.array([0.0, -1.0]), neumann_name='Load', f=F)

solution_writer = ui.solve_nonlinear_static(model, load_steps=10)
system, formulation = ui.create_constrained_mechanical_system_from_component(model)
solution_container = ui.solve_nonlinear_static(system, formulation, model, load_steps=10)

ui.write_results_to_paraview(solution_writer, model, amfe_dir('results/gmsh/ui_example_beam_nonlinear'))
ui.write_results_to_paraview(solution_container, model, amfe_dir('results/gmsh/ui_example_beam_nonlinear'))

###########################################################
### Dynamic nonlinear heterogeneous cantilever beam ###
# Dynamic nonlinear heterogeneous cantilever beam #
###########################################################
model = ui.create_structural_component(mesh)

Expand All @@ -87,9 +90,10 @@
ui.set_dirichlet_by_group(model, 'xy_dirichlet_point', ('uy'), 'Dirichlet_y')

F = triangular_force(0, 0.015, 0.03, 8E5)
ui.set_neumann_by_group(model, 'z_neumann', np.array([0.0, -1.0]), neumann_name='Load', F=F)
ui.set_neumann_by_group(model, 'z_neumann', np.array([0.0, -1.0]), neumann_name='Load', f=F)

solution_writer = ui.solve_nonlinear_dynamic(model, 0.0, 0.1, 0.0001, 10)
system, formulation = ui.create_constrained_mechanical_system_from_component(model)
solution_writer = ui.solve_nonlinear_dynamic(system, formulation, model, 0.0, 0.1, 0.0001, 10)

ui.write_results_to_paraview(solution_writer, model, amfe_dir('results/gmsh/ui_example_beam_nonlinear_dynamic'),
displacements_only=False)

0 comments on commit 2bb0b2f

Please sign in to comment.