large deflection elgeom() result in empty nodal stress output

[hawkoli1987]

Dear All,

I have a peculiar issue when I enable large deflection elgeom() function. I suspect my problem is generic, and my code is a bit lengthy, and so I didn't paste it here for the time being. Please let me know if otherwise.

What I'm doing here is to model a simple rectangular elastic body in 2D-axisymmetric plane. There is an tensile force stretching the elastic body, and the global mesh size is incremented in steps via a while loop. In each step, I'm measuring the total equivalent strain (and a few other outputs) of the body.

The program is able to produce consistent results when large deflection is disabled.
The program will be able to produce a few steps of results when the element size is large, a quick printout of the key parameters I extracted from the code as shown below:

---------------------------------------------------------------------------
max_stn: 77.1%;           disp: 11.12um;           F:  92.4uN;           Ue: 538.0;           #Node: 475;	 esize = 4.281332398719393
---------------------------------------------------------------------------

However, when element size becomes small enough, e.g. <4um when overall area dimension was 25um*100um, the nodal strain output will become empty, as shown in the error below:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-12-893049cd07ff> in <module>
     18 # run the mesh convergence and output the results on the fly
     19 for esize in esizes:
---> 20     sum_dict, N_nodes = get_Ue(mat, geo, c_r, force, None, AESIZE = esize)
     21     list_N_nodes.append(N_nodes)
     22     print(f'#Nodes = {N_nodes};\t esize = {esize}')

<ipython-input-7-8120267635f2> in get_Ue(mat, geo, c_r, force, disp, AESIZE, refine_lvl)
    131 #     print(mapdl.nlist('ALL'))
    132 #     print(mapdl.post_processing.nodal_total_eqv_strain())
--> 133     max_strain = np.max(mapdl.post_processing.nodal_total_eqv_strain())
    134     max_disp = np.max(mapdl.post_processing.nodal_displacement())
    135 

<__array_function__ internals> in amax(*args, **kwargs)

~\anaconda3\envs\pyAnsys\lib\site-packages\numpy\core\fromnumeric.py in amax(a, axis, out, keepdims, initial, where)
   2704     """
   2705     return _wrapreduction(a, np.maximum, 'max', axis, None, out,
-> 2706                           keepdims=keepdims, initial=initial, where=where)
   2707 
   2708 

~\anaconda3\envs\pyAnsys\lib\site-packages\numpy\core\fromnumeric.py in _wrapreduction(obj, ufunc, method, axis, dtype, out, **kwargs)
     85                 return reduction(axis=axis, out=out, **passkwargs)
     86 
---> 87     return ufunc.reduce(obj, axis, dtype, out, **passkwargs)
     88 
     89 

ValueError: zero-size array to reduction operation maximum which has no identity
---------------------------------------------------------------------------

May I have some suggestions for the possible cause of this error?

Thank you!
hawkoli1987

[mikerife]

Hi @hawkoli1987 please zip up your example and post it. I suspect that the solve errored out for some reason and that there are no results available. Mike

[hawkoli1987]

Thanks Mike! I've attached my code sample here for your kind reference.
sample.zip

[mikerife]

Hi @hawkoli1987 so it turns out your model is not so simple after all. I have some questions on the nature of the model, specifically the material. I find it non-physical...Young's modulus of 0.07 psi seems far too small for a real material. Unless it's really exotic!

That being said a model is not unconditionally stable just because it is linear elastic. At that mesh size where you get no results the solution does not converge. And the unconverged result set does not contain all of the needed results. Hence the zero size array error. Here is what the deformed mesh starts to look like with the smaller mesh sizes:

The red dot indicates where the last UY boundary condition is. There are a few things to do to address this:

1 First, use the NSUBST command to add a reasonable number of allowed sub-steps. This is typical for a non-linear (large deflection in this case) solution. I used mapdl.nsubst(10,1000,10) right before the solve command.

2 If the mesh_refine_tip is not used in get_Ue then all the esizes will run to completion.

3 If you want to keep the mesh refinement then saving the results at every sub-step, then post-processing the next to last will at least keep the zero array issue from happening:

    #adding some reasonable sub-stepping
    mapdl.nsubst(10,1000,10)
    mapdl.outres('all','all')
    mapdl.solve()
    mapdl.finish()

    # Enter post-processor
    mapdl.post1()
    mapdl.set('LAST')  # Select the first load step and (first substep)
    mapdl.set('previous')
    # select all element

4 If you need to use the mesh refinement then we need to look at ways to address the mesh deforming so badly. First I'd look at the material...is it hyperelastic and if so maybe a different material model is in order. May also want to look at non-linear mesh adaptivity. Which is automatic remeshing during the solution if some criteria is met...like some element(s) becoming too distorted.

Mike

[hawkoli1987]

Hi Mike,

Many thanks to your quick and informative explanation! Understanding the root cause (small element ->large deformation->non-convergence) really helps! I will also take a look at the non-linear mesh adaptivity as it seems to be a much more efficient way than my manual meshing refinement.

Just to provide a bit more clarity on my work.

I'm trying to model an microscale elastomeric pillar being pulled, the bottom surface of the pillar is attached to another surface, and will be progressively detached from outside inwards. I will have to create a custom routine to check on each step of detachment, readout the total elastic energy, and determine using an equation if the detachment will continue.

First on the unit system, I didn't find an existing unit system in ANSYS that directly matches my application, so I created one such that lengths are in micrometer (um) and forces are in microNetwon (uN). The matching unit for stress will hence be uN/(um^2). The material elastic modulus is about 500 kPa ~= 75psi, which converts to 0.5 uN/(um^2). I guess that's where it might be confusing. However, if I have made a mistake in the conversion, please yell at me!

The next is on non-linear elasticity. The material I have is indeed hyper-elastic. At this early phase of project, I want to keep things simple first by just applying a linear elastic model. Subsequently I will introduce ultimate tensile stress and hyper-elasticity into the model. By doing so, I believe it will definitely help with the convergence issue. On the other hand, I don't know if it will also help to address the 'stress singularity' as I suspect I'm currently having. Let me explain on my suspicion here:

When I refined the mesh at the tip of detachment while keeping the other parameters (loads and geometry) constant, the strain near the tip increases indefinitely with the refinement level. The total strain energy in the fiber kept constant. I don't know if it fits into the 'stress singularity' definition, as it's not similar to the cases of applying a point load or having sharp corners. If it's indeed a 'stress singularity', does it go away after the large deflection, hyper-elasticity and ultimate tensile stress are applied?

Thank you!
hawkoli1987

[mikerife]

Hi @hawkoli1987 Sorry, my mistake on that units of E! Mike

[mikerife]

@hawkoli1987 I also forgot to mention the solution.converged property: https://mapdl.docs.pyansys.com/api/_autosummary/ansys.mapdl.core.solution.Solution.converged.html#ansys.mapdl.core.solution.Solution.converged

You can add this as a check prior to post processing and bypass the post processing if the solution did not converge.
Mike

[hawkoli1987]

Thanks again Mike!

[mikerife]

Hi @hawkoli1987 can you please contact me at [email protected] if you have an issue with this. But I took the overall model idea and recreated it from scratch (so different dimensions, material properties etc) to use as an example of what can go wrong with respect to initial mesh shape. I'd like to use this in a public-facing material. If that is ok.
Mike

[hawkoli1987]

Hi Mike, that's totally fine with me. Glad that my example can be of some help to others. Thanks again to your help!