Simulation of Cyclic Loading in Tensile Test

[MMTUser1]

TL;DR Summary

Tips and suggestions regarding FE simulation of cyclic loading in a tensile test which should display hysteresis trend in the force-displacement curve

I am currently working on a FE simulation project where a tensile test of DP800 steel is subjected to cyclic loading.

It utilizes Yoshida Uemori model (YUM) to formulate the modulus of elasticity (E) to simulate the hysteresis caused by cyclic loading.

We are using a USDFLD subroutine in ABAQUS CAE to formulate the E according to the cyclic loading i.e. using a FORTRAN code, we save the ratio of 'current E' (using YUM) to 'initial E' as a state variable to simulate the process using this ratio. The code recognizes the loading and unloading stages in the curve correctly after checking the PRINT statements in it.

To check the FORTRAN code, we used a RVE (a simple cube element). Original YUM formula mentioned in the research paper works correctly. The modified YUM formula gives incorrect result.

[Mentor Note: link deleted]

I would highly appreciate any tips or suggestions as how to simulate this process correctly. Or how I can test/check the various parameters to ensure the correct simulation of the process.

Thank you.

 

[berkeman]

Welcome to PF.

Paging @FEAnalyst

 

[MMTUser1]

Welcome to PF.

Paging @FEAnalyst

Thank you!

 

[MMTUser1]

Further information:

To check the FORTRAN code, we used a RVE (a simple cube element). Original YUM formula mentioned in the research paper works correctly. The modified YUM formula gives negative E values at certain time points which can be seen as the oscillations on the Force-Displacement curve of Modified YUM formula. (Please refer Image1 below)

I modified the FORTRAN code to convert the negative E values to absolute values so it is always positive. This gave us almost the correct trend in the Force-Displacement curve. (Please refer Image2 below)

I wanted to test this code on the tensile test sample which is a notched tensile test specimen but it does not display the hysteresis as expected. Also, the simulation crashes so the image does not show the complete process (Please refer Image3 below)

I would appreciate any tips or suggestions as to how to simulate the hysteresis caused by cyclic loading on a tensile test specimen.

 

[MMTUser1]

Following is the FORTRAN code snippet:

CALL GETVRM('SINV',ARRAY,JARRAY,FLGRAY,JRCD,JMAC,JMATYP,MATLAYO,
     1     LACCFLA)
               SvM = ARRAY(1)
        
        current_von_mises_stress = SvM
        
        IF (KINC .EQ. 1) THEN
           ! For the first time step, current and previous Von Mises stresses are equal
           previous_von_mises_stress = current_von_mises_stress
           peak_stress = 0
           stress_ratio = 1
        ELSE IF (KINC .EQ. 2) THEN
           previous_von_mises_stress = STATEV(11)
           peak_stress = current_von_mises_stress
           stress_ratio = 1
        ELSE
           ! For subsequent steps, retrieve the previous value using STATEV
           previous_von_mises_stress = STATEV(11)
           peak_stress = STATEV(12)
        ENDIF
        
        IF (current_von_mises_stress .GT. previous_von_mises_stress) THEN       
           IF ((peak_stress .LT. current_von_mises_stress) .AND. (peak_stress .LT. previous_von_mises_stress)) THEN
           ELSE
              IF (KINC .EQ. 2) THEN
                peak_stress = current_von_mises_stress
              ELSE
                peak_stress = previous_von_mises_stress
              ENDIF
           ENDIF
          
           PRINT *, 'Loading stage'
           ! Print the peak equivalent stress of the current loading stage
           PRINT *, 'Peak equivalent stress for current loading stage:', peak_stress
          
           ! Print the current increment number where loading begins
           PRINT *, 'Loading increment number:', KINC
          
           IF (KINC .EQ. 1) THEN
              stress_ratio = 1
           ELSEIF (KINC .EQ. 2) THEN
              stress_ratio = 1
           ELSE
              stress_ratio = current_von_mises_stress/peak_stress
           ENDIF
          
           PRINT *, 'Stress Ratio: ', stress_ratio
          
           ! Modified YUM formula for loading
           E = E_init*(1-STATEV(1)*c1) - (stress_ratio)*(E_init-E_saturation)*(1-exp(-psi*PEEQ))
           PRINT *, 'Modulus of elasticity: ', E
          
           PRINT *, '' ! This line prints a blank line
        ELSE IF (current_von_mises_stress .LT. previous_von_mises_stress) THEN
           IF ((peak_stress .GT. current_von_mises_stress) .AND. (peak_stress .GT. previous_von_mises_stress)) THEN
           ELSE
              peak_stress = previous_von_mises_stress
           ENDIF
          
           PRINT *, 'Unloading stage'
           ! Print the peak equivalent stress of the current unloading stage
           PRINT *, 'Peak equivalent stress for current unloading stage:', peak_stress
          
           ! Print the current increment number where unloading begins
           PRINT *, 'Unloading increment number:', KINC
          
           stress_ratio = current_von_mises_stress/peak_stress
           PRINT *, 'Stress Ratio: ', stress_ratio
          
           ! Modified YUM formula for loading
           E = E_init*(1-STATEV(1)*c1) - (1-stress_ratio)*(E_init-E_saturation)*(1-exp(-psi*PEEQ))
           PRINT *, 'Modulus of elasticity: ', E
          
           PRINT *, '' ! This line prints a blank line
        ELSE
           E = E_init
           PRINT *, 'Modulus of elasticity: ', E
        ENDIF
        
        ! Original YUM formula
        !E = E_init*(1-STATEV(1)*c1) - (E_init-E_saturation)*(1-exp(-psi*PEEQ))
        
        PRINT *, 'Peak stress: ', peak_stress
        PRINT *, 'Previous equivalent stress: ', previous_von_mises_stress
        PRINT *, 'Current equivalent stress: ', current_von_mises_stress
        PRINT *, 'Percentage: ', E/E_init
        PRINT *, '' ! This line prints a blank line
      
       STATEV(9) = E/E_init
       STATEV(10) = E
       STATEV(11) = current_von_mises_stress
       STATEV(12) = peak_stress
       FIELD(1) = STATEV(9)

 

[Tom.G]

Try using Double Precision variables & calculations. Those jagged falling edges suggest there could be a precision problem. (At least it's easy to try!)

 

[FEAnalyst]

I'm a bit late but it can still be of interest. First of all, user-defined material behavior (UMAT, USDFLD or other) should always be tested using a single-element model like in Abaqus Verification Guide. I guess that's what you mean by RVE but actual RVEs are much more complex with periodic BCs and typical inhomogeneity.

Once it works as expected on the single-element model, you can proceed to your actual sample. Now it's best to gradually increase the complexity - try without the subroutine (using only built-in material models) and maybe with constant loading first, for instance. This way it's much easier to diagnose any issues (also with convergence since you've mentioned that the analysis fails - this should be investigated). You have to take many factors into account when analyzing the stress-strain or force-displacement response of the model - output frequency, measurement approach and location, incrementation, presence of artificial damping, mesh density and quality, whether interactions/constraints (if any) work correctly or not and so on.

I know that this reply is rather general but it's really hard to say anything without knowing ALL the details about the model and being able to run the tests ourselves. Only simple issues can be solved without that and this case doesn't look simple.

 

[MMTUser1]

I'm a bit late but it can still be of interest. First of all, user-defined material behavior (UMAT, USDFLD or other) should always be tested using a single-element model like in Abaqus Verification Guide. I guess that's what you mean by RVE but actual RVEs are much more complex with periodic BCs and typical inhomogeneity.

Once it works as expected on the single-element model, you can proceed to your actual sample. Now it's best to gradually increase the complexity - try without the subroutine (using only built-in material models) and maybe with constant loading first, for instance. This way it's much easier to diagnose any issues (also with convergence since you've mentioned that the analysis fails - this should be investigated). You have to take many factors into account when analyzing the stress-strain or force-displacement response of the model - output frequency, measurement approach and location, incrementation, presence of artificial damping, mesh density and quality, whether interactions/constraints (if any) work correctly or not and so on.

I know that this reply is rather general but it's really hard to say anything without knowing ALL the details about the model and being able to run the tests ourselves. Only simple issues can be solved without that and this case doesn't look simple.

Thank you for the reply and all the information! You are right about the RVE, we used a single element model to test the FORTRAN code in the simulation and called it an RVE.

We are now focusing on the single element model to rectify the issue. The original YUM model is giving correct results as far as FORTRAN code is concerned but the issue seems to be in amplitude applied to the load. In case of single element model, the results are shown in Image1.

Here we can see the forces go further negative after unloading which seem to contribute to the plastification under unloading conditions. We found this to be incorrect as there is no compression action in a tensile test, rather just unloading. To correct it, we are tweaking the amplitude applied to the load and testing it. Following is the graph (original - incorrect) of amplitude applied on load (Image4):

These points (original data) are the displacement values taken from the tensile test experiment. This data is reduced and fit to the original data to reduce the number of points using a code. Does it seem accurate to replicate the tensile test according to you? If you have any suggestions, please let me know. I am grateful for the help provided on this forum.

Thank you.

 

[FEAnalyst]

So that cyclic load is modeled as displacement boundary condition with amplitude ? Is incrementation and output frequency adjusted to the resolution of the amplitude ? Abaqus won't do it automatically, it could even solve the analysis in one increment if convergence allows, ignoring the amplitude resolution.

 

[MMTUser1]

So that cyclic load is modeled as displacement boundary condition with amplitude ? Is incrementation and output frequency adjusted to the resolution of the amplitude ? Abaqus won't do it automatically, it could even solve the analysis in one increment if convergence allows, ignoring the amplitude resolution.

Yes, amplitude is applied to the displacement BC to recreate the cyclic loading. The time step is given as 0.001 and there are 1000 points in the amplitude data. That should be the resolution needed for the simulation, if I am correct.

The displacement values are directly taken from the tensile test experiment and made to fit the original data while maintaining the sample size = 1000.

Should we have used a Load BC and apply the amplitude derived from its data instead?

If you have any tips/suggestions, please let me know.

 

[FEAnalyst]

The time step is given as 0.001

I guess that you mean increment size. But is it set for both initial and maximum increment size or are you using direct (fixed) incrementation (not recommended) ? Because you should make sure that the incrementation actually follows the amplitude resolution and that Abaqus is not using larger increments than 0.001.

Also, are the results written every increment ? For history output that's usually the case but not for field output and people sometimes use it to create XY data (plots) in postprocessing.

Displacement control is fine. It should be consistent with the real-life setup.

 

[MMTUser1]

I guess that you mean increment size. But is it set for both initial and maximum increment size or are you using direct (fixed) incrementation (not recommended) ? Because you should make sure that the incrementation actually follows the amplitude resolution and that Abaqus is not using larger increments than 0.001.

Also, are the results written every increment ? For history output that's usually the case but not for field output and people sometimes use it to create XY data (plots) in postprocessing.

Displacement control is fine. It should be consistent with the real-life setup.

Yes, the increment size is set to 0.001 for both initial and maximum increment size. We are using automatic increment and not fixed increments.

The field output request is set for whole model (domain) for frequency = every n increments and n = 0.001.

Thank you for confirming these inputs. Now I am sure that the problem lies somewhere else and have to investigate further.

It is tough to check the issue manually because the log file generated is also too huge to be opened by any text editor. Do you have any recommendations on how to "debug" the simulation step-by-step?

 

[MMTUser1]

I tried to smoothen the amplitude of the loading curve obtained from the experiment. The simulation curve now resembles the experimental curve more.

 

[FEAnalyst]

Do you have any recommendations on how to "debug" the simulation step-by-step?

As I've mentioned in the first reply, the best way is to build up the complexity gradually but if the model is already complete then you can do it the opposite way - simplify or remove features one by one and check the results. Thanks to this approach you have much bigger chances of finding the issue. This is especially important when subroutines are involved.

I tried to smoothen the amplitude of the loading curve obtained from the experiment. The simulation curve now resembles the experimental curve more.

What is the step type here ? Just static or dynamic ? Smooth load applications can be very important in dynamic analyses, especially in the case of quasi-static explicit simulations but also with implicit dynamics.

 

[MMTUser1]

simplify or remove features one by one and check the results

I will try this approach. So far, we built a simple model and have successfully tested till we could apply Original YUM formula through FORTRAN code. I would verify it again.

What is the step type here ? Just static or dynamic ?

Static step.

Smooth load applications can be very important in dynamic analyses, especially in the case of quasi-static explicit simulations but also with implicit dynamics.

A tensile test with cyclic loading would belong to implicit analysis. I wil try with smooth step amplitude and check.

Thank you.