Dear Animesh,

Well, with such a small incremental displacement its probable that you're either working in the linear range of the constitutive model or that the update between equilibrium states is at least linear. So after a single Newton iteration you have probably converged system, and given that the relative error norms do not decrease any further it is also probable that its converged down to machine precision. You can check this by printing out the absolute error norms in addition to the relative error norms (see https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L1354-L1355 <https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L1354-L1355>). If you're satisfied that that this the "problem" (i.e. that the system is converged but its not being detected by the nonlinear solver) then you can simply amend the convergence criteria for the NL solver (see https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L1272-L1273 <https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L1272-L1273>) to include checks on the absolute error.

I hope that this helps.
Jean-Paul

On 23.11.20 11:03, Animesh Rastogi IIT Gandhinagar wrote:
Hi Jean-Paul,

Thank you for your detailed response. I checked everything twice and found my code consistent with what you suggested. I think I understood the issue of Newton method convergence. The issue is occurring when I am giving a very small compressive displacement as a boundary condition (near -1e-06).

Is there anyway I can get the solution converged even with a very small compressive displacement?

Thanks!

Animesh

I am getting something like following -

Assembly method: Residual and linearisation are computed manually.
Grid:
     Reference volume: 1200
Triangulation:
     Number of active cells: 4800
     Number ofset degrees of freedom: 9922
    Setting up quadrature point data...

Timestep 1 @ 1s
_______________________________________________________________________________________
    SOLVER STEP     |  LIN_IT   LIN_RES    RES_NORM RES_U     NU_NORM      NU_U
_______________________________________________________________________________________
  0  CST  ASM  SLV  |       1  0.000e+00  1.000e+00 1.000e+00    1.000e+00  1.000e+00   1  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   2  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   3  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   4  CST  ASM  SLV  |       1  0.000e+00  7.223e-09 7.223e-09    2.843e-08  2.843e-08   5  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   6  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   7  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   8  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09   9  CST  ASM  SLV  |       1  0.000e+00  7.223e-09 7.223e-09    2.843e-08  2.843e-08  10  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  11  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  12  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  13  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  14  CST  ASM  SLV  |       1  0.000e+00  7.223e-09 7.223e-09    2.843e-08  2.843e-08  15  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  16  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  17  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  18  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  19  CST  ASM  SLV  |       1  0.000e+00  7.223e-09 7.223e-09    2.843e-08  2.843e-08  20  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  21  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  22  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  23  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  24  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09  25  CST  ASM  SLV  |       1  0.000e+00  7.223e-09 7.223e-09    2.843e-08  2.843e-08  26  CST  ASM  SLV  |       1  0.000e+00  1.696e-09 1.696e-09    6.675e-09  6.675e-09
..................


On Thursday, November 12, 2020 at 1:34:30 PM UTC+5:30 Jean-Paul Pelteret wrote:

    Hi Animesh,


    > In the documentation in code gallery
    
<https://dealii.org/developer/doxygen/deal.II/code_gallery_Quasi_static_Finite_strain_Compressible_Elasticity.html>,
    it is written that -

    This code is based off of step-44, and its documentation is
    largely identical to what's written there. It looks like that
    comment was missed when adjustments were being made to the
    documentation, so you should ignore it. You can see that the
    stress and tangent are computed directly using the kinematic
    variables passed as arguments to the class methods:

    const SymmetricTensor<2,dim,NumberType>
    <https://dealii.org/developer/doxygen/deal.II/classSymmetricTensor.html>
    tau = lqph[q_point]->get_tau(det_F,b_bar);
    const SymmetricTensor<4,dim,NumberType>
    <https://dealii.org/developer/doxygen/deal.II/classSymmetricTensor.html>
    Jc = lqph[q_point]->get_Jc(det_F,b_bar);

    > I would also need to update the quadrature point data to the
    previous step. Could you please let me know, how am I supposed to
    do that?

    You can see from the code snippet , as well as from the way that
    the assembly is structured

    template <int dim,typename NumberType>
    void Solid<dim,NumberType>::assemble_system(const
    BlockVector<double>
    <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html>
    &solution_delta)
    {
    ...
    const BlockVector<double>
    <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html>
    solution_total(get_total_solution(solution_delta));
    ...
    }

    template <int dim,typename NumberType>
    BlockVector<double>
    <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html>
    Solid<dim,NumberType>::get_total_solution(const
    BlockVector<double>
    <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html>
    &solution_delta)const
    {
    BlockVector<double>
    <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html>
    solution_total(solution_n);
    solution_total += solution_delta;
    return solution_total;
    }

    that the only requirement to ensure a consistent state for the
    entire problem is that the addition of "solution_n" with
    "solution_delta" reflects the configuration of the body that
    you're wanting. So if that you're doing in the main "time loop" is
    something like

    while (time.current() <= time.end())
    {
    solution_delta = 0.0;
    solve_nonlinear_timestep(solution_delta);
    solution_n += solution_delta;

    auto tmp_soln_delta = solution_delta;

    solution_delta = 0.0; // get_total_solution() takes in the
    soltuion_delta, so don't solve the eigenvalue problem with the
    delta accounted for twice
    const bool accept_step = solve_eigenvalue_problem();
    if (accept_step == false)
    {
      solution_n -= tmp_soln_delta;
    
adjust_boundary_conditions_or_timestep_size_to_apply_a_smaller_displacement_increment();
    }

    time.increment();
    }

    then that would seem a reasonable approach to me. As I said
    before, the stresses and strains are kept in sync with the
    solution_total within the assembly loop.

    I hope that this helps you to locate the source of your issues.

    Best,
    Jean-Paul


    On 10.11.20 16:06, Animesh Rastogi IIT Gandhinagar wrote:
    Hi Jean-Paul,

    Another question which might be the reason why I am getting the
    inconsistencies in my newton-raphson convergence. I just realised
    this..

    What I am doing is the following -

    I am giving small compressive displacements to the
    film+substrate. I am also calculating the eigenvalues of the
    tangent_stiffness matrix (after the convergence of newton
    method), to understand the stability of the system. The moment I
    encounter a negative smallest eigenvalue, I want to go to a
    previous time step and would like to solve for a small
    compressive displacement. For instance, if the compression is 0.1
    at this time step and the eigenvalue is negative, I go to the
    previous time step and give a compressive displacement of 0.05.

    Here is my concern, I have added a condition that, if I find a
    negative eigenvalue during any time step, I update the
    /solution_n = solution_n - solution_delta/; and then run the time
    step again (with a smaller compressive displacement) after doing
    /solution_delta = 0.0/. Since the next state is also dependent on
    the updated quadrature point data (stresses etc.), I would also
    need to update the quadrature point data to the previous step.
    Could you please let me know, how am I supposed to do that?

    In the documentation in code gallery
    
<https://dealii.org/developer/doxygen/deal.II/code_gallery_Quasi_static_Finite_strain_Compressible_Elasticity.html>,
    it is written that -

    "/The first function is used to create a material object and to
    initialize all tensors correctly: The second one updates the
    stored values and stresses based on the current deformation
    measure //Grad//u//n/" (inside the Point History class). However,
    I only find one function which is setup_lqp, and not the
    update_lqp or something similar.

    Thanks!

    Animesh

    On Tuesday, November 10, 2020 at 5:24:13 PM UTC+5:30 Animesh
    Rastogi IIT Gandhinagar wrote:

        Hi Jean-Paul,

        There seem to be some inconsistencies in the simulations
        results that I am getting, therefore I wanted to check if I
        have applied the material laws properly accounting for the
        heterogeneities. Below, I have addted the edited functions
        that are responsible for that. Could you please let me know
        if there is something wrong here? I have checked the material
        properties (mu1, mu2 etc), and they are being read correctly
        by the code from the parameters file. I have also checked
        that the material ids (1 and 2) are correctly applied to the
        cells where I need them.

        --------------------------------------
        void setup_lqp (const Parameters::AllParameters &parameters,
        const unsigned int &materialid)
            {

            //Material 1 - Thin Film
        if (materialid == 1)
            {
        
std::cout<<"1_id"<<std::endl<<parameters.mu1<<std::endl<<parameters.nu1<<std::endl;
              material.reset(new
        Material_Compressible_Neo_Hook_One_Field<dim,NumberType>(parameters.mu1,
                             parameters.nu1));
            }
            //Material 2 - Substrate
            if (materialid ==2)
            {
              material.reset(new
        Material_Compressible_Neo_Hook_One_Field<dim,NumberType>(parameters.mu2,
                                   parameters.nu2));
            }

            }

        -------------------------------------------

        template <int dim,typename NumberType>
          void Solid<dim,NumberType>::setup_qph()
          {
            std::cout << "    Setting up quadrature point data..." <<
        std::endl;

        quadrature_point_history.initialize(triangulation.begin_active(),
        triangulation.end(),
                                                n_q_points);

            for (typename Triangulation<dim>::active_cell_iterator cell =
                   triangulation.begin_active(); cell !=
        triangulation.end(); ++cell)
              {
                const
        std::vector<std::shared_ptr<PointHistory<dim,NumberType> > >
        lqph =
                  quadrature_point_history.get_data(cell);
                Assert(lqph.size() == n_q_points, ExcInternalError());

                for (unsigned int q_point = 0; q_point < n_q_points;
        ++q_point)
                  lqph[q_point]->setup_lqp(parameters,
        cell->material_id());
              }
          }

        --------------------------------------

        Thanks!

        Animesh
        On Tuesday, November 3, 2020 at 12:09:03 PM UTC+5:30 Animesh
        Rastogi IIT Gandhinagar wrote:

            Hi Jean,

            Sorry for late response. Thank you for your detailed
            answer. Using your suggestions I was able to assign the
            material ids to different area of my domain of interest.

            Thanks again!

            Animesh

            On Tuesday, October 27, 2020 at 3:05:39 AM UTC+5:30
            Jean-Paul Pelteret wrote:

                Dear Animesh,

                Let me preface my suggestion by saying that you can
                introduce some heterogeneity in more than one way.
                That is, you're not bound to using a material ID
                (which is what I'm going to suggest), but you could
                use some geometric arguments to determine which class
                or material properties should dictate the
                constitutive response at a continuum point. You also
                need not use the class-based structure that this code
                galley example has. There are other tutorials (e.g.
                step-8) that use Functions that return the material
                coefficients at any point in space. So you should
                keep that in mind as you read my suggestion.

                In my opinion, the simplest approach would be the
                following:

                1. The PointHistory::setup_lqp() function
                
<https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L720-L724>
                is the point at which the constitutive law at some
                quadrature point is set. So you should extend this
                function in some way to make a decision about which
                material it is in, and what the material law at that
                point is. The first thing that you might do is add an
                argument for the material ID, and then choose
                constiutitve parameters based on the material ID.
                You'll note that the concrete constitutive law is
                stored in a pointer. This makes it easy to extend to
                support other constitutive laws, because you can for
                example abstract the
                Material_Compressible_Neo_Hook_One_Field class into
                some base class and other classes (maybe a
                Material_Compressible_Mooney_Rivlin_One_Field) and
                then, based on the material ID, select which
                constitutive law is applied at that quadrature point.
                2. The PointHistory::setup_lqp() function is called
                from Solid::setup_qph()
                
<https://github.com/dealii/code-gallery/blob/master/Quasi_static_Finite_strain_Compressible_Elasticity/cook_membrane.cc#L1203-L1211>
                , and more specifically inside a loop over all cells.
                Its possible to retrieve the cell->material_id() or
                cell->center() (or, if you create a quadrature rule
                to help you then you can get the actual position of
                the quadrature point if you want such granularity).
                You can then pass this information in to
                PointHistory::setup_lqp() to help decide which
                constitutive law / parameters to apply. If using the
                material ID then you will, of course, want to set the
                material ID when building or reading in your
                triangulation.

                That is, at least, how I typically do things when
                using this code as a basis for extension. I hope that
                this helps you introduce heterogeneities into your
                problem!

                Best,
                Jean-Paul


                On 26.10.20 08:15, Animesh Rastogi IIT Gandhinagar wrote:
                Hello All,

                I am trying to solve the wrinkling problem of a thin
                stiff film attached to a soft substrate in dealii. I
                am using the code - from the code gallery for this
                purpose. For the simulations, I would have to
                consider the film and the substrate with different
                material properties.

                I understand that I would have to do something
                within the class Qudrature point history. However, I
                am getting a little confused on how to approach it.
                It would be great if someone could help me with this.

                This is the code that I am using for this purpose -
                
https://dealii.org/developer/doxygen/deal.II/code_gallery_Quasi_static_Finite_strain_Compressible_Elasticity.html
                
<https://dealii.org/developer/doxygen/deal.II/code_gallery_Quasi_static_Finite_strain_Compressible_Elasticity.html>

                Thanks!

                Animesh

-- The deal.II project is located at
                http://www.dealii.org/ <http://www.dealii.org/>
                For mailing list/forum options, see
                https://groups.google.com/d/forum/dealii?hl=en
                <https://groups.google.com/d/forum/dealii?hl=en>
                ---
                You received this message because you are subscribed
                to the Google Groups "deal.II User Group" group.
                To unsubscribe from this group and stop receiving
                emails from it, send an email to
                dealii+un...@googlegroups.com.
                To view this discussion on the web visit
                
https://groups.google.com/d/msgid/dealii/26db8831-e190-4fcd-8247-b7172ab880fen%40googlegroups.com
                
<https://groups.google.com/d/msgid/dealii/26db8831-e190-4fcd-8247-b7172ab880fen%40googlegroups.com?utm_medium=email&utm_source=footer>.

-- The deal.II project is located at http://www.dealii.org/
    <http://www.dealii.org/>
    For mailing list/forum options, see
    https://groups.google.com/d/forum/dealii?hl=en
    <https://groups.google.com/d/forum/dealii?hl=en>
    ---
    You received this message because you are subscribed to the
    Google Groups "deal.II User Group" group.
    To unsubscribe from this group and stop receiving emails from it,
    send an email to dealii+un...@googlegroups.com.
    To view this discussion on the web visit
    
https://groups.google.com/d/msgid/dealii/f1209ac5-a4b8-46d4-84f3-12ca38f7d8f2n%40googlegroups.com
    
<https://groups.google.com/d/msgid/dealii/f1209ac5-a4b8-46d4-84f3-12ca38f7d8f2n%40googlegroups.com?utm_medium=email&utm_source=footer>.

--
The deal.II project is located at http://www.dealii.org/ <http://www.dealii.org/> For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en <https://groups.google.com/d/forum/dealii?hl=en>
---
You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com <mailto:dealii+unsubscr...@googlegroups.com>. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/4e377c7f-240f-4130-a82e-f2e77fee5fafn%40googlegroups.com <https://groups.google.com/d/msgid/dealii/4e377c7f-240f-4130-a82e-f2e77fee5fafn%40googlegroups.com?utm_medium=email&utm_source=footer>.

--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- You received this message because you are subscribed to the Google Groups "deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/a0af6dbe-d56e-4a88-7448-b1ed7dd687ce%40gmail.com.

Reply via email to