Rigidity of GAP elements in contact
GAP element rigidity will depend on the material of the parts in contact... and also on the mesh size! Learn how to calculate it!12 December 2022
I must confess, that I haven’t always pay attention to mesh size. At the beginning of my career, I simply build my models and was happy when they converged. I was vaguely (and I mean vaguely!) aware, that the element size does “something”. But I was absolutely happy with whatever default options were set! Now I’m a bit wiser, and I will show you a few important things about mesh convergence!
Deciding on the correct mesh size is a rather difficult thing. Sure, you can base that on your experience in any given analysis type. But this experience has to come from somewhere, right? The most accurate approach would be doing a mesh convergence check. I will show you an example of how to do it.
If you are interested in the general principles of how this works, you may want to read my older post. Here, I want to show you how to apply this method, in a simple example (so you can follow along).
I admit, that what I will do here, is a trick to some degree. Mostly, because I will do the example where we can actually know the correct answer (I will calculate it later). I know that you most likely won’t have the luxury like that in your “actual” work. And this is pretty obvious too! If you could easily calculate the outcome by hand… why do FEA in the first place?!
My advice would be, to try doing the same thing, with a different example anyway. Preferably, one that you can easily solve by hand! Maybe even something that is more similar to things you usually solve in FEA? This way, you can check for yourself if this works. I won’t hide I’m a big fan of “check if works first” approach to FEA stuff. So doing simple examples at the start is actually how I learn stuff.
Of course, you can always say: “yea, I saw Łukasz did this, and I’m sure it will work, so let’s just do it”. I can’t (and won’t!) stop you. But I think that learning is best when you are doing simple examples. Mostly because you can easily and quickly “tweak” them just to see “what happens if…”. The “what happens if…” part and simply “experimenting” with various things, will teach you the most!
This is why I’ve selected two simple cantilevers for our problem. We will do them one, by one as I want to show you something neat!
The first one is the “classical” cantilever with the load on the tip. I will use plate elements to model this problem. It’s basically a cantilever plate (0.2 x 0.5m) with only 10mm thickness.
To make this more interesting, I will be doing QUAD4 and QUAD8 elements (so the 4-node and 8-node quad elements). This way, we can see how fast the answer converges in both cases. This will also give us an opportunity to learn some other interesting things!
The model is of course super simple, but there is a small thing I did to make my life easier. You see if I would model the support “only” as fixed at the supported end, the stress distribution I would get would look like this:
This is not what we would expect right? I mean, with the load on the “tip”, one would expect a nice “smooth” stress increase toward the support. And here, we get some weird “stress concentration” near the support, seemingly for no actual reason.
Well, there is a reason of course, and if you already took my awesome FEA Quiz, you know the answer already! This is caused by the fact, that our cantilever wants to change its width under the load as well (Hooke’s law), and near support, it “cannot”, since it is supported. Without going into details, we can fix this quickly by assuming that our “plate” is actually “infinite” in width. This way, we will get a nice “smooth” stress distribution like this:
All of the stresses I will display here are NOT von Misses stress. Instead, I’m displaying the normal stress in the “length direction” of the cantilever. This normal stress will “resist” the bending! Please note, that such stress is usually called “normal X” or “normal Y” or something similar… BUT this is in your Finite Element orientation, so it doesn’t have to follow “global coordinate system direction”.
In my case, all finite elements have their Coordinate System with Y axis “to the right” on the picture above. So I’m displaying “Normal Stress Y”. Even thoug the global UCS would suggest that it should be “X” not “Y’. As I wrote before It’s the local element axis, not the global UCS, deciding which normal stress is in which direction!
It may seem counterintuitive at first, I know. But imagine the same cantilever, that is oriented 45 deg between axis X and Y. Global UCS would not be aligned with elements at all! but still you could display the “normal stress” along and across the cantilever!
With this small trick, we are ready to rumble!
With the model set up already, we can finally start. Let’s just say, that we do not know the maximal stress in our element yet. So we want to answer a question: What is the maximal normal stress in our element?
So, what we could do, is to model out a piece, use the auto mesher and press calculate. This is precisely what I did! Below you can see the mesh size Femap would assume automatically based on the model size. This is what we would get for QUAD4 elements:
As you can see on the right, the max normal stress is 283MPa. There are 2 questions you need to answer in such a case:
Those two questions are the cornerstone of mesh convergence. In essence, if you are certain about your answers in FEA, and you know that you are using “optimal” mesh size… then you are set! However, if you have doubts, you are in the right place! Let’s go deeper into this subject together to find those answers!
We don’t know if the above answer (from the “default mesh”) is correct or not. This means that we need to assume 2 things (at once!) and consider a third:
Knowing that I’ve decided to test QUAD4 against QUAD8 elements for this problem. This is a super small task, so it will be hard to measure time. Most runs are “well below 1sec”, and at such times, it’s hard to measure time accurately. But we can see how the model converges anyway!
Since it won’t cost me much, I can start from the single element (!) and move toward smaller and smaller mesh. Of course, you don’t have to be so “precise” in mesh convergence studies. You can comfortably start from the mesh size you feel “makes sense”.
Let’s deal with the QUAD4 elements first. I won’t show you all of the models I did (this is quite boring, to be honest). Below you can see several of those:
If we would put all of those outcomes on a simple chart, it would look like this:
You can clearly see, that the answer we are getting nicely converges. What I mean by this is, that with the increase in the number of elements answer is “stabilizing”. Seeing the chart, we could also say, that it “looks” like the answer is around 300MPa.
For me, that is enough. But if need be, you could do some math “hurt” to the chart to see what is the precise value this converges to. But since I know that the difference between 297MPa (50 elements) and 299.7MPa (500 elements) is only 0.9%. I’m just “fine” with the outcome from the model with 50 elements.
Personally, I wouldn’t dig deeper unless perhaps I would be doing fatigue analysis. 1-2% of accuracy always satisfies me. Heck! One could argue that the outcome from a model with 10 elements along the length (285MPa) that produces around 5% error is still reasonable! It only depends on how big a model you really have. Waiting an additional 3 weeks for computing to increase accuracy from 5% to 1% may not be the best idea. But of course, it can be a great idea as well… this really depends on what you are doing!
Now it will go quicker since we already know what’s going on! Below, you can see the outcomes for QUAD8 elements. Those are, of course, the same conditions and all.
I’m not sure if you noticed… but the answer is always 300MPa! I mean, starting from the 1 element, up to the 500 elements along the length… it’s constantly 300MPa! I guess we don’t need to plot the chart in this case, as it would be an awesome horizontal like! This clearly shows that the answer “converged”, and that the correct result is 300MPa.
I must admit, that even knowing this, and understanding why… I still would not feel comfortable with 1-3 elements along the length!
This is a place, where we can stop for a second, and wonder what we already got (and why). Also, I would like to clear, why QUAD8 elements “always knew” the correct answer!
Firstly, let’s calculate the stress by hand. I guess you already know it should be 300MPa, but it never hurts to check right?
Awesome! So, now we know, how much stress we should get. And, both QUAD4 and QUAD8 converged to a proper solution (yea!). But it was absolutely clear, that QUAD8 elements converged instantly, while QUAD4 elements struggled quite a bit. Let’s take a closer look.
You see, QUAD8 elements can “easily” do linear stress distribution along their length. This is thanks to the additional nodes they have on their edges. You may think about those as “data points”. Those points allow us to better calculate what is going on inside the element (thanks to shape functions). If this topic fascinates you, you can read more here, and also take a look at one of my YouTube Live sessions about this topic. What this means, is that the QUAD8 element has no issues with the fact that the out-of-plane bending causes linear changes in the stress distribution on the element.
On another hand, QUAD4 elements are “blind” to such actions. If there is an out-of-plane bending causing stress, the QUAD4 element wants this stress to be “constant” on its surface. This is not always the case for QUAD4 elements. For example, in their plane, they would “allow” linear stress distribution. But this is a completely different topic (and not connected to our problem).
The cantilever is set up, to have a nice, linear change in the bending moment. This means, that the stress caused by this bending will be “perfectly” linear as well. You can see this below:
This means, that even a single QUAD8 element can “handle this”. Simply put, it is “well suited” to deal with the linear change of stress across the element. QUAD 4 elements on another side struggled. This is because they want to have a constant stress value across their width. The problem can be easily shown this way:
As you can see, each QUAD4 element simply “gathers up” all the stress from its area, and averages it. Or at least, this is what happens in a super simplistic look. This is why reducing the size of the element changed outcomes. The smaller area of the element you have, the less impact this effect has. To the point that if the element was “infinitely small” the averaged outcome across its minimal area, would actually be precise.
From this example, you could read, that using QUAD8 elements means that you don’t have to worry about mesh convergence. This is of course not the case… and also a reason why there is a second example in this post!
All right, you already know how to solve a mesh convergence problem! So perhaps you may want to try doing the second example on your own before reading along? This is a good practice for sure!
This time, our cantilever is loaded on its surface. Seemingly a small change, but it will provide us with all the emotions we will need for today’s post! Take a look:
I will use the same boundary conditions (along with symmetry assumption) as in the previous case. And I will use the same finite elements, with the same setup. So I will measure how many elements “along the length” there are in the model. This will be an indicator of the mesh size. Not much more to say here, maybe apart from this:
In this example, I’m using the amount of elements along a selected edge as an “indicator” of how small mesh I have. You can always use other mesaures. I admit, that I’m not a fan of element “size” used for this. You know the “width od the element edge” in milimeters. But, you can always go with other “amounts of elements” for mesh convergence. You don’t even have to measure amount of elements in entire model. This can easily be limited to the amount of elements on any particular part of your model. Just use whatever in your mind nicely represent how small elements you have. Preferably in the place where you are searching for an outcome, not necesairly in the entire model.
This also means, that you don’t need to have the mesh of the same size everywhere. It’s completely fine to have bigger elements in a lot of places, and small elements “where it counts”!.
Just as in the last example, we will start with QUAD4 elements. Again, below you can see selected cases (I removed element division for 500 FE, so you can actually see something):
If you have an awesome memory (and you pay attention like crazy!) you may notice that the outcomes actually converge “slower”. I mean, in the previous example with 3 QUAD4 the answer was 250MPa, while wow it’s “only” 216.7MPa. For 10 elements, the difference is still there (285 vs 271.5MPa). Heck! Even for 500FE along the length, the difference is still there (299.7 vs 299.4MPa). Clearly, something in this task is more “difficult” to solve. And yes… you guessed it. The final answer is the same in both cases.
The difference can be nicely shown on the convergence chart. The blue chart is the “old one” from the previous case with the load at the tip. The current case is marked in orange.
I think you can easily see, that the orange curve converges “slower” to the proper solution. It still does, but it’s just not as simple. Interestingly in the previous case, 5 QUAD4 elements were more or less as accurate as 10 QUAD4 elements now (around 270MPa). Previously, I’ve mentioned that 10 elements along the edge could be considered decent (5% accuracy). Now it’s unthinkable, at least for me. In this case, the accuracy of 10 elements along the long edge is only 10%. That is a bit too much for me. This shows us, that mesh size doesn’t only depend on the geometry – it also depends on the load. Or, to be more accurate, on the stress state in the region we are calculating.
I really hope you are already curious why that is the case! But first…
Before we move on, let’s take a look at QUAD8 elements in our small challenge! This time, they won’t perform as awesome as previously:
It’s clear, that they converge nicely, but it’s not a “from the start” thing that happened in the previous example. Here, we would need 3FE along the length to be “decent”, which obviously is way fewer than in the case of QUAD4 elements. Convergence plot can be seen below:
Just as in the previous case, firstly let’s take a look at what the actual answer is:
But knowing the max value in itself isn’t a solution. Just as before, we should also take a look at stress distribution. This is where all the magic happens. You see, since the load is uniformly distributed, the bending moment distribution will be parabolic, just as you can see below:
Knowing the stress distribution, let’s do the trick we did last time. We already know, those QUAD4 elements want to have constant stress on their surface. Well… at least in our case, since it is caused by the out-of-element-plane bending. I guess you can easily imagine what happens next!
I think that after seeing the above schematic, you already know what is going on. QUAD4 elements want to have a constant value of stress in our case. So they work in a way to “make this happen”. Please note, that it’s not as if QUAD4 element “knows” the correct answer, and simply calculate its own value backward (as I showed you above). Instead, the “internal math” of the element is limited, and this limitation produced the outcomes in the elements. And it so happens, that those outcomes are in relation to the “actual” outcome in a way I’ve shown you above.
This also explains, why now it’s more difficult for the answer to converge. The quadratic curve is more difficult to approximate with blocks of “constant values”, than the line from the previous example. A certain amount of elements may be “good enough” to approximate stress that is distributed linearly. But in the same model, if the stress is “quadratic”, the same elements may not suffice! You know, like this:
Funny enough, you won’t know what stress distribution will be in your model. I mean, sure… you can suspect stuff. In simple cases (like this one) you may even be sure. But normally, especially in big and complex models, you just won’t know. So it’s usually better to assume “worse” than “better” and take smaller elements “just in case”.
The only thing that remains, is what “happens” with QUAD8 elements! Of course, this is a similar problem. I’ve already told you that QUAD8 is “well suited” to have the stress change linearly on their surface. Here, the stress change is “quadratic”, so even the QUAD8 doesn’t get the “perfect” answer from the start. But they have a huge advantage: they can have a linear distribution of stress (rather than constant stress as in QUAD4). This makes it way easier for them to “fit” into the actual stress distribution compared to QUAD4.
While the stress changed in a linear way, QUAD8 offered a “perfect fit”. With quadratic distribution, there are some “errors” but, as you can surely imagine, way smaller than in the case of QUAD4. In fact, small enough to make the drawing below a bit awkward…
This is why, QUAD8 had an easier time converging the proper answer, even with the use of a very small amount of elements.
This may suggest that QUAD8 elements are simply better, but there is something else, we need to consider…
When all is said and done, it’s not the “smallest amount of elements needed” that wins! After all, you can get a similar accuracy between the 2 models. One would have fewer QUAD8 elements, the other more of the QUAD4. Since accuracy is similar… the question is, which one is the “better one”?
For me, the victory goes to the model that computes faster! As I wrote at the start, this example is too small to measure time in a reliable way. Most models are computed in under 1s. However, 500 FE along the edge of QUAD4 took 12s… while for QUAD8 it was 52s. Significantly more. Also, those things will “compound” and with bigger models, the difference is even bigger.
What I try to say here is, that there is no easy answer. Sadly, you won’t learn here if it is better for you to use QUAD4 or QUAD8… simply put, I just don’t know!
Checking waht is best for you:
There is a way you could check, which elements are best in most of your cases!
Just select a model, that will represent a “typical problem” you will have to deal with in your work. Then, perform a mesh convergence check for both QUAD4 elements and QUAD8 elements. Heck! You can even throw TRI elements there if you feel fancy.
Measure the time it takes each run to compute.
After all, decide what accuracy level suits you… and simply check which type of elements give you the required accuracy with the shortest computing time!
I admit, you will have to do mesh convergence several times for the same problem (for different element types). This will take some time, and surely is impractical in your daily work. But doing this to 1-2 typical examples, can show you which elements are best for you. The time this will save you in the future on computing, is well worth the effort!
But I feel, I should mention something else. Raw time is a great measure, but don’t be too fixated on it as well.
There are other considerations, you may wish to take into account. The one that comes to my mind is, how many small details you have in your model. Usually, the mesh is small near small features you wish to have in your model. But, on the other hand, “better” elements (like QUAD8) can “allow you” to use super big elements in “general”. In such a case, the transition between small elements near details, and super-big elements elsewhere will cause the mesh to be “ugly” with poor quality. I would think about this as well when deciding on the element type.
Well, this is it! I hope that you enjoyed this mesh convergence example and that it helped you learn something about how FEA works! Let’s try to recap what we talked about here:
If you like my teaching style, definitely check my FREE FEA course below.
10 Lessons I’ve Learned in 10 Years!