Week - 2 - Explicit and Implicit Analysis
INTRODUCTION
The finite element method (FEM) is a numerical problem-solving methodology commonly used across multiple engineering disciplines for numerous applications such as structural analysis, fluid flow, heat transfer, mass transport, and anything existing as a real-world force. This practice systematically yields equations and attempts to approximate the values of the unknowns. This method subdivides the overall problem into simpler sub-issues that are easier to solve. In turn, these sub-issues called finite elements require either an implicit or an explicit analysis.
Image Source: Solid Mechanics Classroom
Solvers utilise the finite element method (also called finite element analysis) to simulate naturally or artificially occurring phenomenons. This numerical technique is the foundation of simulation software in order for engineers (like civil and mechanical engineers) to assess their designs for tension, weak spots, etc., before the prototyping or implementation stages.
When running a finite element analysis on a given problem, it is important to identify whether it is a static or a dynamic problem, in order to achieve a stable and accurate solution in a time-efficient manner. But that is not enough to decide on the type of analysis or method that would be best suited to solve the problem.
The key deciding factor of selecting the best analysis method (explicit or implicit) is time and time plays an immense role in situations where kinetic energy is a key player in the analysis. With kinetic energy, there are a lot of quick movements and as a result, accelerations and decelerations occur. In these cases, we can say the analysis is time-dependent. Examples of this would be a car crash or a ballistic simulation.
On the other hand, if the load is gradual, the effect of inertia in the analysis is basically non-existent. In such cases, time has no effect and the model in the analysis is always in a state of static equilibrium. Such an analysis would be considered time-independent.
A general rule of thumb is to typically consider implicit analysis for time-independent problems since it can solve bigger chunks of the problem in one go and is ‘implicitly’ stable whereas explicit would be more appropriate for time-dependent ones.
Image Source: Implicit vs Explicit Approach in FEM by Yasin Çapar
IMPLICIT ANALYSIS
Essentially, implicit analysis tries to achieve static equilibrium through multiple iterations (within a certain tolerance) with respect to the applied loading or displacement. In non-linear analyses, the solution is obtained over several tiny increments and the solution for one ‘step’ would be based on the result acquired in the previous step.
This method ends up being highly time-consuming in complex problems due to the fact that the solution needs to attain equilibrium (also known as ‘converging’) for each step before moving on to the next. So, in large models or models with significant nonlinearity, inverting the matrix (the basis of implicit analysis) is highly computationally expensive and sometimes impossible.
Although, if the load application is gradual and the boundary conditions are well-defined (such as in cases of static loading), this can be a very efficient method. After global equilibrium is converged, the solver calculates all the local finite element variables (stresses, strain, etc.) for each increment.
An example where an implicit solver would be more appropriate:
Stresses generated in a bevel gear pair
EXPLICIT ANALYSIS
Explicit analysis on the other hand does not require convergence in each step. There are no iterations involved like in the implicit method. It works on the assumption that the problem is in equilibrium, instead of trying to ensure ‘global equilibrium’. This is a large chunk of the computational effort behind implicit solving. Instead of that, the explicit solver focuses on solving the local finite element variables for a particular increment and simply moves on to the next one. The solver does it very quickly compared to the implicit method.
The speed of solving does have a drawback though and that is ensuring the time step is extremely small. Hence, even an analysis that is capturing fractions of a second will be broken into hundreds of thousands of increments, which is most of them since explicit solvers focus on cases where a large amount of energy is expended in a short period of time or in other words, there are large strain rates/velocities involved within the simulation, like in a car crash.
An example where an explicit solver is used:
Bullet penetrating a bucket
IMPLICIT OR EXPLICIT?
As discussed previously, the general rule of thumb is to go for the explicit solver in time-dependent problems and the implicit solver where problems are time-independent. But the term 'time-dependency' really depends on other factors specific to each method (implicit or explicit).
So, implicit solvers are used when:
- Sum of forces = 0 (all the way through the analysis)
- Strain rates are very low (< 0.001 per second)
- Velocity is irrelevant
- An increase in stress can be defined easily as a function of strain only.
Examples: static forces, slowly applied displacements, self-weight, constant pressure loading.
And explicit solvers are used in cases where:
- Sum of forces = mass x acceleration
- Strain rates are above 0.001 per second
- Velocity is above 10m per second
- An increase in stress is a function of both strain and strain rate
Examples: automotive crash, ballistic events, metal stamping.
Image source: What is explicit dynamics in ANSYS (mechead.com)
OBJECTIVE
To solve the given non-linear static problem (defined by the equation `F(u) = u^3+9u^2+4u`) using implicit and explicit methods, with a residual tolerance of 0.001.
In the equation `f(u) = u^3+9u^2+4u`, F is the force in the model and u is the displacement. Also, an incremental load scheme is to be used in both cases.
The stiffness can be defined using the equation as `k(u) = (df)/(du) = 2u^2 + 18u + 4`
The two equations for force and stiffness will play big roles in both methods.
SOLVING USING THE EXPLICIT METHOD
If the nonlinear force versus displacement relationship is not known, but the stiffness is known, it is possible to construct the graph of force versus displacement numerically. One method is to use an incremental explicit load control scheme. For this example, three equal load increments (or steps) of ∆F = 1 unit each are used to load the tension bar. The variables used are 'u' for displacement, 'f' for internal force in the subject, 'F' for the external force applied to the subject and k for stiffness. Incremental displacements or incremental externally applied forces are represented as ∆u or ∆F, respectively. The relationship ∆F = k∆u is utilised. The analysis is as follows:
`f(u) = u^3+9u^2+4u`
Stiffness, k(u) =`(df)/(du) = 3u^2 +18u + 4`
Step 1
Starting from rest, we can assign a displacement value to obtain the stiffness of the
At `u = 0 = u_0`
`k(u_0) = 3(0) + 18(0) + 4 = 4`
With increments of `DeltaF = 1`,
`Deltau_1 = (DeltaF)/(k(u_0)) = 1/4`
Therefore, `u_1 = u_0 + Deltau_1 = 0 + (1/4) = 0.25`
and
`f_(i nt) = (u_1)^3 + 9(u_1)^2 + 4(u_1) = (0.25)^3 + 9(0.25)^2 + 4(0.25)` = 1.58
Step 2
Similarly,
`k(u_1) = 3(0.25)^2 + 18(0.25) + 4 = 8.6875`
Hence,
`Deltau_2 = (DeltaF)/(k(u_1)) = 1/8.6875 = 0.115`
`u_2 = u_1 + Deltau_2 = 0.25 + 0.115 = 0.365`
Also,
`k(u_2) = 3(0.365)^2 + 18(0.365) + 4 = 10.97`
and
`f_(i nt) = (u_2)^3 + 9(u_2)^2 + 4(u_2) = (0.365)^3 + 9(0.365)^2 + 4(0.365)` = 2.71
Step 3
Hence,
`Deltau_3 = (DeltaF)/(k(u_2)) = 1/10.97 = 0.091`
`u_3 = u_2 + Deltau_3 = 0.365 + 0.091 = 0.456`
`F_(ext) = 3(DeltaF) = 3(1) = 3`
`f_(i nt) = (u_3)^3 + 9(u_3)^2 + 4(u_3) = (0.456)^3 + 9(0.456)^2 + 4(0.456)` = 3.79
We can see that the internal and external forces are unequal and hence are not in equilibrium
Summarizing the results obtained:
SOLVING USING THE IMPLICIT METHOD
Newton-Raphson iterations are used to ensure equilibrium is achieved before moving on to the next step. Basically, an incremental force is applied to advance the solution forward at the beginning of each step. Although, internal and external forces won't be equal to each other unless the stiffness varies linearly for the particular step. Hence, in order to achieve equilibrium, corrections must be made to the displacement. This is accomplished by using Newton-Raphson iterations to minimize the residual, `R(u) = f_(i nt) − F_(ext)`.
Expanding the residual as a Taylor series about the current displacement `u^j` gives
`R^(j+1) = R^j(u^j) + (dR(u^j))/(du) δu^(j+1) + ... = R^j + k(u^j)δu^(j+1) + ...`
By neglecting higher order terms in the series, setting `R^(j+1) = 0` and solving for `δu^(j+1)` the
following correction is obtained
`δu^(j+1) = −[k(u^j)]^(−1).R^j`
The iteration variable is j and it may take several iterations for the norm of the residual, `||R||`, to be reduced below the chosen tolerance criteria (for the current problem, a tolerance of 0.001 is chosen). The implicit analysis proceeds as follows:
Step 1
We start from u = 0
At `u = 0 = u_0`
`k(u_0) = 3(0)^2 + 18(0) + 4 = 4`
`DeltaF = 1`
`Deltau_1 = (DeltaF)/(k(u_0)) = 1/4`
Therefore, `u_1 = u_0 + Deltau_1 = 0 + (1/4) = 0.25`
Checking the residual,
`F_(ext) = DeltaF = 1`
`f_(i nt) = u_1^3+9u_1^2+4u_1 = 0.25^3 + 9(0.25)^2 +4(0.25) = 1.578`
Therefore, `R^(0) = f_(i nt) - F_(ext) = 1.578 - 1 = 0.578` which is greater than 0.001
As a result, Newton-Raphson iterations are necessary.
Calculating the new total correction to `u_1 = u_1^((0))`,
`δu^((1))=−[K(u_1^((0)))]^(−1)⋅R^((0))`
We previously established that `K(u_1) = 8.6875`, which is to be taken as `K(u_1^((0)))`
`δu^((1)) = -(8.6875)^(-1) . 0.578 = -0.0665`
Therefore, the updated value of `u_1^((1)) = u_1^((0)) + δu^((1)) = 0.25 - 0.0665 = 0.1835`
We can check the residual again for `F_(ext) = 1`:
`f_(i nt) = (u_1^((1)))^3+9(u_1^((1)))^2+4(u_1^((1))) = 0.1835^3 + 9(0.1835)^2 +4(0.1835) = 1.0432`
Therefore, the residual `R^((1)) = f_(i nt) - F_(ext) = 1.0432 -1 = 0.0432`
Again, `||R^((1))|| > 0.001` => This means further iterations are required
Calculating the new total correction to `u_1 = u_1^((1))`
`δu^((2))=δu_1−[K(u_1^((1)))]^(−1)⋅R^((1))`
We know that `u_1^((1)) = 0.1835`
Therefore, `δu^((2))= -0.0665 −[K(u_1^((1)))]^(−1)⋅R^((1))`
= `-0.0665 - [(3(0.1835)^2 + 18(0.1835) + 4)^(−1) . 0.0432]`= -0.0723
Therefore, the updated value of `u_1^((2)) = u_1^((0)) + δu^((2)) = 0.25 - 0.0723 = 0.1777`
Checking the residual again for `F_(ext) = 1`:
`f_(i nt) = (u_1^((2)))^3+9(u_1^((2)))^2+4(u_1^((2))) = 0.1777^3 + 9(0.1777)^2 +4(0.1777) = 1.00061`
Therefore, the residual `R^((2)) = f_(i nt) - F_(ext) = 1.00061 -1 = 0.00061`
`||R^((2))|| = 0.00061` which is lesser than the specified 0.001
This means another Newton-Raphson iteration is not required.
Therefore, final value of `u_1 = u_1^((2)) = 0.1777`
Step 2
`k(u_1) = 3(0.1777)^2 + 18(0.1777) + 4 = 7.2933`
`DeltaF = 1`
`Deltau_2 = (DeltaF)/(k(u_1)) = 1/7.2933 = 0.1371`
Therefore, `u_2 = u_1 + Deltau_2 = 0.1777 + 0.1371 = 0.3148`
Checking the residual,
`F_(ext) = DeltaF = 2`
`f_(i nt) = u_2^3+9u_2^2+4u_2 = 0.3148^3 + 9(0.3148)^2 +4(0.3148) = 2.1824`
Therefore, `R^(0) = f_(i nt) - F_(ext) = 2.1824 - 2 = 0.1824` which is unfortunately greater than 0.001
As a result, Newton-Raphson iterations are necessary.
Calculating the new total correction to `u_2 = u_2^((0))`,
`δu^((1))=−[K(u_2^((0)))]^(−1)⋅R^((0))`
`K(u_2) = 3(u_2)^2 + 18(u_2) + 4 = 3(0.3148)^2 + 18(0.3148) + 4 = 9.9639`
`K(u_2)` is to be taken as `K(u_2^((0)))`
`δu^((1)) = -(9.9639)^(-1) . 0.1824 = -0.0183`
Therefore, the updated value of `u_2^((1)) = u_2^((0)) + δu^((1)) = 0.3148 - 0.0183 = 0.2965`
We can check the residual again for `F_(ext) = 2`:
`f_(i nt) = (u_2^((1)))^3+9(u_2^((1)))^2+4(u_2^((1)))`
`= 0.2965^3 + 9(0.2965)^2 +4(0.2965) = 2.0032`
Therefore, `R^(1) = f_(i nt) - F_(ext) = 2.0032 - 2 = 0.0032` which is still, unfortunately, greater than 0.001. Further iterations are required
Calculating the new total correction to `u_2 = u_2^((0))`
`δu^((2))=δu^((1))−[K(u_2^((1)))]^(−1)⋅R^((1))`
We know that `u_2^((1)) = 0.2965
Therefore, `δu^((2)) = -0.0665 −[K(u_2^((1)))]^(−1)⋅R^((1))`
= `-0.0183 - [(3(0.2965)^2 + 18(0.2965) + 4)^(−1) . 0.0032]= -0.01863`
Therefore, the updated value of `u_2^((2)) = u_2^((0)) + δu^((2)) = 0.3148 - 0.01863 = 0.2962`
Checking the residual again for `F_(ext) = 1`:
`f_(i nt) = (u_2^((2)))^3+9(u_2^((2)))^2+4(u_2^((2)))`
`= 0.2962^3 + 9(0.2962)^2 +4(0.2962) = 2.0003`
Therefore, the residual `R^((2)) = f_(i nt) - F_(ext) = 2.0003 -2 = 0.0003` which is less than 0.001.
No further iterations are required.
Therefore, final value of `u_2 = u_2^((2)) = 0.2962`
STEP 3
`k(u_2) = k(0.2962) = 3(0.2962)^2 + 18(0.2962) + 4 = 9.5928`
`DeltaF = 1`
`Deltau_3 = (DeltaF)/(k(u_2)) = 1/9.5928 = 0.1042`
Therefore, `u_3 = u_2 + Deltau_3 = 0.2961 + 0.1042 = 0.4003`
Checking the residual,
`F_(ext) = DeltaF = 3`
`f_(i nt) = u_3^3+9u_3^2+4u_3 = 0.4003^3 + 9(0.4003)^2 +4(0.4003) = 3.1075`
Therefore, `R^(0) = f_(i nt) - F_(ext) = 3.1075 - 3 = 0.1075` which is unfortunately greater than 0.001
As a result, Newton-Raphson iterations are necessary.
Calculating the new total correction to `u_3 = u_3^((0))`,
`δu^((1))=−[K(u_3^((0)))]^(−1)⋅R^((0))`
`K(u_3) = 3(u_3)^2 + 18(u_3) + 4 = 3(0.3148)^2 + 18(0.3148) + 4 = 11.6861`
`K(u_3)` is to be taken as `K(u_2^((0)))`
`δu^((1)) = -(11.6861)^(-1) . 0.1075 = -0.0092`
Therefore, the updated value of `u_3^((1)) = u_3^((0)) + δu^((1)) = 0.4003 - 0.0092 = 0.3911`
We can check the residual again for `F_(ext) = 3`:
`f_(i nt) = (u_3^((1)))^3+9(u_3^((1)))^2+4(u_3^((1)))`
`= 0.3911^3 + 9(0.3911)^2 +4(0.3911) = 3.0008`
Therefore, `R^(1) = f_(i nt) - F_(ext) = 3.0008 - 3 = 0.0008` which is less than 0.001
No further iterations are required
Therefore, final value of `u_3 = u_3^((1)) = 0.3911`
Results tabulated:
RESULTS PLOT
We can see that plots of these two methods, even though they are similar in nature, they drift away from each other. The implicit method, in general, provides more accurate results and can be considered to be a pseudo-yardstick here. As discussed earlier, the key to ensuring more useful results via the explicit method is to use smaller time steps. This creates a better equilibrium between the internal and external forces and as a result, a more accurate plot and results. But this does not mean the implicit is ideal in every situation. The implicit method is computationally very expensive and time-consuming, due to its iterative nature.
CONCLUSION
The given non-linear static problem was solved using implicit and explicit methods. The results were compared graphically. Between the explicit and implicit methods, the right method needs to be chosen according to the type of simulation. If the event involved in the simulation is gradual, implicit should be considered and if it's an instantaneous event, explicit would be more appropriate. The time taken to complete solving needs to be considered as well.