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


What is the meaning of 'smix' stiffness, for strain rate zero?

I'm not sure.




On Wed, Jan 19, 2011 at 8:59 AM, Robert Cimrman <cimr...@ntc.zcu.cz> wrote:
On Wed, 19 Jan 2011, Andre Smit wrote:

Robert

See [1] to illustrate the initial instability. I've fixed the
displacement and the loading time (the time steps are really just
providing the loop). The output shows displacement, force. If you plot
the force you see that the load appears to converge after a few steps but
then becomes very variable gain.

Yes, I can see it.

I am not sure how it works, though: how is it possible, that the force can actually grow? I assumed that once an element fails, it's stiffness should be reduced for ever. Then the force should be lower and lower, right?

What is the meaning of 'smix' stiffness, for strain rate zero?


r.

I've also attached the new abaqus generated mesh.

[1] http://paste.pocoo.org/show/323259/

On Wed, Jan 19, 2011 at 8:04 AM, Robert Cimrman <cimr...@ntc.zcu.cz>
wrote:
     Oh, this output! Then you can try also putting

     'output_format' : 'h5',

     into the options. All results will then be saved into a
     single HDF5 file with many time steps.

     hth,
     r.


On Wed, 19 Jan 2011, Andre Smit wrote:

     Thanks Robert - the code works great. Regarding the
     output - I was really
     concerned about the large number of vtk files being
     generated when
     running thru simple.py with the extra iterations
     applied to stabilize the
     force. Now as a standalone I can suppress this and the
     script runs faster
     too! I like the idea of gradually reducing modulus for
     stabilization.



     On Wed, Jan 19, 2011 at 3:09 AM, Robert Cimrman
     <cimr...@ntc.zcu.cz>
     wrote:
          One more thing: there was a bug in the Output()
     class
          constructor arguments, so you need the latest
     sources for the
          code below.

          r.


     On Wed, 19 Jan 2011, Robert Cimrman wrote:

          The zero forces are because you have changed the
          traction() function to the form used by
          material-defining functions. It should look (like
     in
          your previous versions) as:

          def traction(ts, coors, bc=None):
            disp = -(ts.dt*(ts.step+1))
            val = nm.empty_like(coors[:,0])
            val.fill(disp)
            return val

          i.e. as a bc-defining function.

          As for too much output printed, yes, we should
     address
          the issue 16 one day in a systematic way. I do not
     know
          yet, how to control the output flexibly enough yet
          unobtrusively.

          Right now you can do several things:

          - you can use the fact, that all messages that
     SfePy
          outputs are printed using the output() function -
     it's
          not just a function but an instance of the Output
     class
          (see sfepy/base/base.py), and as such can be
          controlled:

          So after adding:

            from sfepy.base.base import output

            output.set_output(filename='output_log.txt')

          to the beginning of your main() function, all the
          output is redirected into the specified file. At
     the
          same time, messages you print using 'print' are
     still
          printed.

          - Note that you could use an Output instance to
     log
          your messages too:

          from sfepy.base.base import Output

          my_output = Output('strainer:',
     filename='my_log.txt',
          combined=True)

          my_output('hello!')

          I have modified your script in this way, see [2].

          - the above method can be used together with using
          simple.py - you get the standard, and your custom
          outputs well-split.

          Just put:

          from sfepy.base.base import output, Output

          my_output = Output('strainer:',
     filename='my_log.txt',
          combined=True)
          output.set_output(filename='output_log.txt')

          at the top of the file, and run it with simple.py

          The example at [2] already has this - it can be
     run as
          both standalone and via simple.py...

          - the deprecation warnings coming from numpy/scipy
     can
          be suppressed by:

          import warnings
          warnings.filterwarnings("ignore")

          (also done at [2])

          - the Log class in sfepy.base.log can be used to
     plot
          and log data, see
          examples/standalone/live_plot/live_plot.py

          As I think of it now, it seems that the issue 16
     can be
          tackled by documenting the above features :)

          Does it help?

          r.
          PS: I used the old cylgeo.vtk mesh, as I do not
     have
          the new one, so you have to fix the paths and
     region
          definitions in the paste [2].

          [2] http://paste.pocoo.org/show/323128/

          On Wed, 19 Jan 2011, Andre Smit wrote:

                Robert

                Using simple.py is generating too much
                output. To have better control
                over this I decided to convert to a
                standalone app as shown at [1] using
                one of the sfepy examples as reference.
                Something is not working as the
                forces at each time step are zero. Can you
                please have a look.

                Thanks  


                [1] http://paste.pocoo.org/show/323091/


                On Tue, Jan 18, 2011 at 3:42 PM, Andre Smit
                <freev...@gmail.com>
                wrote:
                     Robert - I had a good meeting with
                Yannis this afternoon. He
                     is concerned that the strain condition
                in the specimen after
                     each time step as currently coded is
                not stable and suggested
                     that a few iterations be run at each
                displacement interval to
                     ensure that the element moduli
                converge before moving onto
                     the next time step or displacement
                interval. To test his idea
                     I ran a few iterations at a constant
                displacement and found
                     that the forces do fluctuate quite a
                bit but appear to
                     converge (somewhat) after about 10
                iterations. I presume the
                     moduli of the elements would also
                stabilize.



                On Tue, Jan 18, 2011 at 10:48 AM, Robert
                Cimrman
                <cimr...@ntc.zcu.cz> wrote:
                     On Tue, 18 Jan 2011, Andre Smit wrote:
                     Robert

                     I've revised the code loops as shown
                in [1].
                     Seems to work. Yes, I still
                     find myself thinking in FORTRAN :(


                You will get used to it :)

                The code works for me (and is much shorter
                and readable now).

                As a matter of fact, the functions in
                sfepy.mechanics.matcoefs could be easily
                vectorized as well,
                and subsequently your glame() - it's
                another EasyToFix issue
                topic (new issue 150).

                     I'm looking at the hyperelastic option
                now. I've
                     also asked Yannis
                     Tassoulas to consider joining as a
                co-author. He
                     may have some ideas
                     regarding plasticity options.


                OK. Plasticity would come handy - I still
                have not got to
                implementing it (first, I would have to
                study it :]).


                r.

                     [1]
                http://paste.pocoo.org/show/322755/


                     On Tue, Jan 18, 2011 at 8:14 AM, Andre
                Smit
                     <freev...@gmail.com>
                     wrote:
                          Thanks Robert!

                          I'm working on the revisions and
                updates and
                     will get back.
                          The abstract is due at the end of
                January
                     but I hope to get
                          the bulk of the paper done by
                then as well,
                     which is due in
                          Feb sometime i.e. if abstract is
                accepted.


                     On Mon, Jan 17, 2011 at 7:13 AM,
                Robert Cimrman
                     <cimr...@ntc.zcu.cz> wrote:
                          At [1] I have posted the code
                with "correct"
                     avgqp()
                          function. Yours was ok, as long
                as all the
                     quadrature
                          point weights were the same.

                          The new implementation defines an
                auxiliary
                     material
                          with the quadrature point data,
                and calls
                     the
                          de_volume_average_mat term. Note
                that the
                     order of the
                          quadrature has to be the same as
                the one
                     used for
                          evaluating the strain.

                          Then it seems to me, that you
                could avoid
                     all those
                          loops in
                          functions Ft(), Fc(), Ec(), by
                using numpy
                     vectorized
                          operations, i.e. numpy.exp()
                instead of
                     math.exp().
                          Also the failure loops (around
                line 193)
                     could be IMHO
                          vectorized using numpy.where()
                and fancy
                     indexing.
                          Try it as a NumPy excercise :) -
                you get big
                     speed gain
                          by this, and also the code would
                be shorter
                     and more
                          readable. I can give you a hand
                if you get
                     stuck, of
                          course.

                          cheers,
                          r.

                          [1]
                http://paste.pocoo.org/show/322183/


                     On Mon, 17 Jan 2011, Robert Cimrman
                wrote:

                          Hi Andre,

                          (sorry for the delayed response -
                I'm
                          accommodating in Fribourg again,
                and have
                     not yet
                          wifi at home...)

                          On Fri, 14 Jan 2011, Andre Smit
                wrote:

                                Robert

                                I think I have something to
                share.
                                How would you feel about
                co-authoring
                                a paper to this conference
                - mainly
                                to further expose SfePy but
                also to
                                check that I'm not applying
                it
                                incorrectly?


                          That would be great, thanks for
                proposing
                     it!

                                Using the strain rate
                dependency
                                example you provided I
                recoded an
                                example looking at the
                failure of a
                                cylindrical specimen under
                                compression. The cylinder
                is made up
                                of elements that I treat
                                individually or discretely,
                checking
                                the strain rates in each
                and using a
                                relationship established in
                the lab,
                                I calculate the compressive
                                strengths in the elements
                based on
                                these strain rates. If the
                stress in
                                an element exceeds its
                strength then
                                the element fails and I
                assign a low
                                stiffness to that element.
                The code
                                I'm using is available
                here:
                              
                 http://paste.pocoo.org/show/320646/


                          I will look at the code and try
                it ASAP.
                     What is
                          the deadline, anyway?
                          A few notes straight from my head
                follow.

                          What happens with the linear
                elastic
                     assumption
                          when you reduce stiffness in some
                elements?
                     The
                          code has some hyperelastic terms,
                in case
                     they
                          would be needed.

                                Would you mind looking and
                commenting
                                on this code. I'm a bit
                skeptical
                                about the avgqp function
                which
                                averages the strains at the
                quad
                                points.
                                Note that I use a normal
                distribution
                                to vary the moduli of the
                elements
                                initially.  


                          To get averaged element values of
                something
                          defined in quadrature points, one
                can use
                     either
                          directly de_volume_average_mat
                term (a fake
                          material definition would then be
                needed, I
                     can
                          show you), or use directly
                similar code to
                     the
                          one in the body of that term.
                           
                                I've attached the mesh:
                cylgeo.vtk.

                                What I aim to do with the
                paper is
                                demonstrate the influence
                of specimen
                                geometry on strength. For
                simple
                                compression and tension
                tests (on
                                cylinders) this is trivial
                but with
                                indirect tensile tests, for
                example,
                                there is both compression
                and tension
                                failure, and I hope to show
                using
                                FEM that the strengths of
                materials
                                determined using these
                tests are
                                significantly influenced by
                the
                                geometry of the specimen. 
                   


                          Interesting!

                                The force-displacement
                curves
                                currently generated by the
                code look
                                promising. I'm not sure the
                analysis
                                can be refined by using
                smaller
                                elements?


                          The mesh is not very fine, yes.
                Another way
                     to
                          refine the analysis might be to
                use another
                          material relation
                (hyperelastic?).

                          r.

                                On Fri, Jan 14, 2011 at
                3:19 AM,
                                Robert Cimrman
                <cimr...@ntc.zcu.cz>
                                wrote:
                                     I see :) Good luck
                then!

                                     r.


                                On Thu, 13 Jan 2011, Andre
                Smit
                                wrote:

                                     No plans currently -
                still
                                working on an analysis.
                                     Maybe if I can get it
                                     done before the end of
                the month
                                :)

                                     On Thu, Jan 13, 2011
                at 2:43 AM,
                                Robert Cimrman
                                     <cimr...@ntc.zcu.cz>
                                     wrote:
                                          Hi Andre,

                                          are you planning
                to go
                                there?

                                          r.


                                     On Mon, 10 Jan 2011,
                Andre Smit
                                wrote:

                                          FYI:

                                        
                                   
                              
                   
              
          http://congress.cimne.com/CFRAC2011/frontal/default.asp


                                          --
                                          Andre


                                --
                                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

                                --
                                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.



                          --
                          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.



                     --
                     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




                     --
                     Andre

                     --
                     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.



                --
                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




                --
                Andre

                --
                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.



          --
          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.



     --
     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

     --
     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.



--
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

--
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.



--
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