On 6/26/23 08:35, Robert Cimrman wrote:
Hi Sebastian
On 6/22/23 11:20, Sebastian Bachmann wrote:
Hi! I have a tetrahedral mesh and know the nodal displacements of the mesh from an external software. I now want to calculate the element' strains and thought I give sfepy a go. However, the learning curve is very steep and I struggle a bit with this - I thought - easy task.
I thought I basically need this function here: https://sfepy.org/doc-devel/primer.html#post-processing but how do I initially setup the problem and set the already solved displacement vector?
This is what I got so far:
import numpy as np from sfepy.discrete.fem import Mesh, FEDomain from sfepy.discrete import Problem # Given: mesh and nodal displacements mesh = Mesh.from_file('the_mesh.vtk') u_vec = np.random.random((mesh.n_nod, 3)) # here I would insert the actual displacements. Use random numbers for testing domain = FEDomain('domain', mesh) omega = domain.create_region('Omega', 'all') # ??? pb = Problem('elasticity', domain=domain, ...) # ??? # ??? strain = pb.evaluate('ev_cauchy_strain.2.Omega(u)', mode='el_avg')Am I on the right track here?
Yes, you can use something like this:
import numpy as np from sfepy.discrete.fem import Mesh, FEDomain, Field from sfepy.discrete import FieldVariable, Integral from sfepy.terms import Term
# Given: mesh and nodal displacements mesh = Mesh.from_file('cylinder.vtk') u_vec = np.random.random((mesh.n_nod, 3)) # here I would insert the actual displacements. Use random numbers for testing
domain = FEDomain('domain', mesh) omega = domain.create_region('Omega', 'all')
field = Field.from_args('fu', np.float64, 'vector', omega, approx_order=1) u = FieldVariable('u', 'parameter', field, primary_var_name='(set-to-None)')
u.set_data(u_vec)
integral = Integral('i', order=2)
strain1 = u.evaluate(mode='cauchy_strain', integral=integral)
term = Term.new('ev_cauchy_strain(u)', integral, omega, u=u) term.setup()
strain2 = term.evaluate(mode='qp')
print(np.linalg.norm(strain1 - strain2))
- It is not necessary to create a Problem when there is nothing to solve.
- Variables can be evaluated in quadrature points in several modes ('val', 'grad', 'div', 'cauchy_strain')
- Terms can be evaluated individually.
Then, to get element strains and save the results to a VTK file:
estrain = term.evaluate(mode='el_avg')
from sfepy.base.base import Struct out = {'strain' : Struct(name='output', mode='cell', data=estrain)}
mesh.write('strain.vtk', out=out)
r.