# Implementation of Z^2 error estimator in Abaqus for adaptive mesh refinement

Currently, I am working on a remeshing routine for my simulations (Abaqus 6.14-1) using python scripts. The simulation deals with the Brinell indentation test and as the remeshing software I use Gmsh 4.5.2 (although the newest probably also works fine).

The main script is given by

main.py:
createIndentationTest.py
runJob.py                    #(only one small step)
odbExport.py
writeGeo.py                  #(.geo files are Gmsh readable geometry files)
remeshGeo.py                 #(here occurs the change from quadrangular elements to triangular ones)
importInAbaqus.py
runJob.py                    #(only one small step)
while totalTime-stepTime>0:
odbExport.py
writeGeo.py
remeshGeo.py
importInAbaqus.py
runJob.py                #(only one small step)
totalTime-=stepTime


As an example, my last routine created 6 separate jobs, where after every step the mesh was exported to Gmsh, remeshed and imported back in Abaqus to continue the simulation. The first job starts with quadrangular elements (element size of 0.25) for one step. Then the remaining 5 jobs use triangular elements with a custom density function, as can be seen in the following figure (The reason for switching from quadrangular elements to triangular ones during the remeshing process is that it was easier for me to work with in Gmsh).

Now, part of my thesis is to determine the error during the whole simulation and remeshing process and I have been given the task to use the $$Z^2$$ error estimator. When compared to my reference simulation, where I ran a whole simulation from start to finish with quadrangular elements with the element size of 0.5, I can see there is a slight discrepancy when analyzing the stresses, as shown in Figure 2.

For now, I assume that the results on the left-hand side are more accurate, therefore I can assume that there is an error (or deviation) in my simulation results on the right-hand side where the remeshing was implemented. For the past 2 weeks, I have studied some works on the $$Z^2$$ error estimator, but I fail to see how to implement this in my script. Therefore, I wanted to ask, if anyone knows how I can apply the error estimator in my python code.

I know that I have to export some data from the COORD fieldOutput and use these to solve this equation

$$F(\underline{a})=\sum_{i=0}^n\,\left(\sigma_h(x_i,y_i)-\sigma_p^*(x_i,y_i)\right)^2=\sum_{i=0}^n\,\left(\sigma_h(x_i,y_i)-\underline{P}(x_i,y_i)\,\underline{a}\right)^2$$

where $$\underline{a}$$ is a set of unknowns, $$\sigma_h$$ is the approximate solution to $$\sigma$$, $$\underline{P}$$ is a set of appropriate polynomial terms and $$\sigma_p^*=\underline{P}\,\underline{a}$$.

Another question I have stems from the fact that I change element types during the simulation process. To estimate the error, the $$Z^2$$ method uses so-called superconvergent points. But I have no clue if there will arise a problem when using these since I change my element type once during the whole simulation process.

I hope, someone can help me with this.

• I think that it's a good idea to use the same type of elements and subdivide them. That is, having a sequence of nested meshes. That way you can think that each coarser mesh represents a subspace of a finer one. – nicoguaro Apr 29 '20 at 4:17
• Thanks for the reply. So, now with some minor changes to my script, I managed to use only triangles in my meshes (especially the first one). I have a question though. I not that experienced in meshes and FEM, therefore I don't really understand what you mean by subdividing them or nested meshes. – fruitiest Punch May 1 '20 at 12:09
• That means that you take each triangle and subdivide it into 4 triangles, so the nodes from the previous are part of the new mesh. In Gmsh this is five with the option "Refine by splitting". – nicoguaro May 1 '20 at 16:22