For comparison I was trying to run simple.py on the attached. It fails with the following error:
[grassy@myhost CFRAC]$ ~/sfepy/simple.py one.py
sfepy: left over: ['lam', 'Dijkl', '__builtins__', '_filename', '__file__', 'stiffness_tensor_lame', '__name__', 'verbose', 'nm',
'mesh_file', 'youngpoisson_to_lame', 'young', 'mu', '__package__', 'traction', 'stress_strain', 'output_dir', 'poisson', '__doc__']
sfepy: guessing abaqus
sfepy: reading mesh (/home/grassy/Dropbox/sfepy/msh/cyl_02.inp)...
sfepy: ...done in 0.01 s
sfepy: setting up domain edges...
sfepy: ...done in 0.01 s
sfepy: setting up domain faces...
sfepy: ...done in 0.01 s
sfepy: creating regions...
sfepy: Omega
sfepy: YZ-Plane
sfepy: Top
sfepy: XZ-Plane
sfepy: XY-Plane
sfepy: ...done in 0.01 s
sfepy: equation "equilibrium":
sfepy: dw_lin_elastic_iso.2.Omega(Asphalt.lam, Asphalt.mu, v, u) = 0
sfepy: setting up dof connectivities...
sfepy: ...done in 0.00 s
sfepy: using solvers:
ts: ts
nls: newton
ls: ls
sfepy: ====== time 0.000000e+00 (step 1 of 101) =====
sfepy: updating variables...
sfepy: ...done
Traceback (most recent call last):
File "/home/grassy/sfepy/simple.py", line 113, in <module>
main()
File "/home/grassy/sfepy/simple.py", line 110, in main
app()
File "/home/grassy/sfepy/sfepy/applications/application.py", line 29, in call_basic
return self.call( **kwargs )
File "/home/grassy/sfepy/sfepy/applications/simple_app.py", line 121, in call
pre_process_hook=self.pre_process_hook)
File "/home/grassy/sfepy/sfepy/solvers/generic.py", line 217, in solve_direct
nls_status=nls_status)
File "/home/grassy/sfepy/sfepy/solvers/generic.py", line 137, in solve_evolutionary_op
for ts, state in time_solver( state0 ):
File "/home/grassy/sfepy/sfepy/solvers/ts.py", line 128, in __call__
state = step_fun( self.ts, state0, *step_args )
File "/home/grassy/sfepy/sfepy/solvers/generic.py", line 71, in time_step_function
problem.time_update( ts )
File "/home/grassy/sfepy/sfepy/fem/problemDef.py", line 460, in time_update
functions, create_matrix)
File "/home/grassy/sfepy/sfepy/fem/problemDef.py", line 422, in update_equations
self.equations.time_update(self.ts, ebcs, epbcs, lcbcs, functions)
File "/home/grassy/sfepy/sfepy/fem/equations.py", line 215, in time_update
self.variables.equation_mapping(ebcs, epbcs, ts, functions)
File "/home/grassy/sfepy/sfepy/fem/variables.py", line 239, in equation_mapping
var.equation_mapping(bcs, var_di, ts, functions)
File "/home/grassy/sfepy/sfepy/fem/variables.py", line 1372, in equation_mapping
self.eq_map.map_equations(bcs, self.field, ts, functions, warn=warn)
File "/home/grassy/sfepy/sfepy/fem/dof_info.py", line 316, in map_equations
if vv is not None: val_ebc[eq] = vv
TypeError: array cannot be safely cast to required type
I believe the input is correct. Any ideas?
a
On Tue, Jan 25, 2011 at 2:52 AM, Robert Cimrman <
cimr...@ntc.zcu.cz> wrote:
It works for me *with a tweak, see below) even for fvar = 0.1, the other settings left as they were, for both compression and
extension.
One thing I noted is, that you compare floating point numbers using the '==' operator - this is not safe due to floating point
approximation. Use
is_min(stiffness0, smin)
instead of
stiffness0 == smin
otherwise you may get wrong indices...
So the current stiffness update logic is:
1. mark elements where stress > strengh using
stiffness = nm.where(abs(szz) > abs(fc), smin, smix)
2. using
stiffness = nm.where(is_min(stiffness0, smin), smin, stiffness)
ensures, that already failed elements remain failed. The "already failed" are taken from the previous load step, not the
sub-iterations.
I think that makes sense.
So good luck with the actual data fitting :)
r.
On Mon, 24 Jan 2011, Andre Smit wrote:
Think I've found the error - we need to increase fvar to 100 say and
increase the number of n_steps = 100.
I believe the problem is now solved. Next step is to compare the output
to actual data and code the tensile failure checks.
a
On Mon, Jan 24, 2011 at 12:38 PM, Andre Smit <
freev...@gmail.com>
wrote:
Robert - sorry for all these posts but I'm making corrections
as I go along.
Please disregard the previous code link, I've corrected as
shown at [1]. The force-displacement curve is now non-linear
as one would expect. Still bothering me is the change in
moduli for the same rate. If the displacement interval (dt)
is -0.1 and the rate is abs(dt), one should get the same
strain rate if dt is -0.01, but it increases. Still trying to
track down this error.
[1]
http://paste.pocoo.org/show/326129/
On Mon, Jan 24, 2011 at 11:26 AM, Andre Smit
<
freev...@gmail.com> wrote:
I've revised the code as shown in [1] to illustrate
what I mean by assigning smin at convergence. It is
interesting to see that the problem still converges
when the elements start to fail!
[1]
http://paste.pocoo.org/show/326090/
On Mon, Jan 24, 2011 at 9:08 AM, Robert Cimrman
<
cimr...@ntc.zcu.cz> wrote:
On Mon, 24 Jan 2011, Andre Smit wrote:
Thanks
works like a charm! Two comments:
1) Since the model is 1/8 of actual
the deformation rate is half of
actual test. So the strain rate in
the Fc, Ft and Ec functions should be
doubled:
e *= 2 # With 1/8 model the applied
deformation rate is half of actual
test
Yes, modify it as needed :)
2) Don't you think failed elements (smin)
should only be assigned at the
end of an iteration run and not during
sub-iterations.
I do not know - if you do not change the stiffness in
sub-iterations, then nothing changes after the first
one, and there is no need for them, right?
r.
On Mon, Jan 24, 2011 at 2:05 AM, Robert
Cimrman <
cimr...@ntc.zcu.cz>
wrote:
ok, so the patch is in... the script
should work
out-of-the-box.
r.
On Sat, 22 Jan 2011, Robert Cimrman wrote:
I have overhauled the script, so that
there are no
globals, subiterations
can be ended as needed etc. The
traction() function is
removed, the
loading displacement is given simply
by value. I hope
it does what you
intended.
It needs a patch to TimeStepper, that
I will be able to
upload first on
Monday, so either just look at it, or
fix it
yourself... :)
r.
----- Reply message -----
From: "Andre Smit"
<
freev...@gmail.com>
To: <
sfepy...@googlegroups.com>
Subject: Conference
Date: Fri, Jan 21, 2011 02:47
Robert
I believe I've sorted out the
oscillation problem,
which was a bug in the
code. I've pasted the latest version
at [1] that now
converges tightly
after 3 or 4 iterations. I have a few
comments.
As you will see:
a) I'm runnng the solver repeatedly
with increasing
displacements. This
was mainly because I'm unable to
change the
displacement variable (disp)
in traction() outside of the function.
Globals and
lists didn't work for
me here.
b) To overcome this obstacle, the
displacement for each
solver/displacement cycle is saved to
file and read
from disk in
traction().
c) I iterate the problem until the
force converges and
then increment the
displacement.
d) The element stiffness matrix (Ele)
is saved to file
and loaded at the
start of the subsequent displacement
interval.
Comments:
1) I considered adding a variable to
ts, say ts.disp
and adding the
displacement to that, making it
available in traction.
I'm sure there is
a better way though.
2) I also cannot terminate the ts loop
and am forced to
run through all
ts.n_steps. This is a problem since
the forces converge
at different
number of steps at which point I would
like the loop to
end.
[1]
http://paste.pocoo.org/show/324169/
Andre
On Wed, Jan 19, 2011 at 10:22 AM,
Robert Cimrman
<
cimr...@ntc.zcu.cz>
wrote:
On Wed, 19 Jan 2011, Andre Smit
wrote:
Robert - I'm trying to change the
displacement
(disp)
in the traction
function when the force
stabilizes but for some
reason
the scope of disp
within traction() appears
different to that within
calc_force(), even
though disp is defined as global.
Any ideas?
old_force = 0
def calc_force(pb, ts, state):
global disp,old_force
d = state()
pb.remove_bcs()
f =
pb.evaluator.eval_residual(d)
pb.time_update()
f.shape =
(pb.domain.mesh.n_nod, 3)
p = []
for n in
pb.domain.regions['Top'].vertices[0]:
p.append(f[n][2])
p = nm.array(p)
if abs(old_force-p.sum()) <
fvar:
my_output(disp,",",p.sum())
disp -= ddisp
old_force=p.sum()
def traction(ts, coors, bc=None):
global disp
#disp = -(ts.dt*(ts.step+1))
print "traction disp: ", disp
val =
nm.empty_like(coors[:,0])
val.fill(disp)
return val
I would use brute force:
globals()['disp'] = ...
You cannot really change value of a
global that is
immutable
(IMHO)...
A common idiom is to use list instead,
like:
disp = [0]
def traction(ts, coors, bc=None):
disp[0] = ...
What is the meaning of
'smix' stiffness, for
strain rate
zero?
I'm not sure.
That might be the problem causing the
oscillations...
r.
--
You received this message because you are subscribed to the Google Groups "sfepy-devel" group.
To post to this group, send email to
sfepy...@googlegroups.com.
To unsubscribe from this group, send email to
sfepy-devel...@googlegroups.com.
For more options, visit this group at
http://groups.google.com/group/sfepy-devel?hl=en.
--
Andre