Explicit vs Implicit Method
EXPLICIT AND IMPLICIT ANALYSIS
AIM:
The aim of this project is to compare the difference between Explicit and Implicit solver methods and to determine which provides the most feasible solution by solving them numerically.
To achieve this, we perform a Non-linear static analysis on a bar element constrained on one side while subjected to an external load on the other end using an Incremental load scheme using both Implicit and Explicit method and then solve them numerically. These results are then compared with the exact solution using MATLAB to determine the most feasible method.
IMPLICIT ANALYSIS
In a Non-linear analysis, the Implicit method is a solver method used to find the solution of a time-dependent problem by considering the state or solution of the system only at the current time step. Moreover, when using the implicit method as a solver technique, it has a convergence criterion which aims at converging the resultant values at each time step before proceeding to the next time step.
Implicit analysis solves for the domain variables using the below equation,
`u(n+1)=u(n)+ ∆T ×F(u(n+1))`
Here, to compute the required parameter, in this case, the displacement, at the current time step i.e u(n+1), the implicit solver method makes use of the solution obtained from that specific time step i.e F(u(n+1)).
EXPLICIT ANALYSIS
In a Non-linear analysis, the Explicit method is a solver method that is much similar to that of an implicit method. However, unlike the implicit method, this solver technique computes the solution of a time-dependent problem at the current time step by considering the state or solution of the system at the previous step. Moreover, it does not follow a convergence criterion at each increment but rather computes the local variables at that step and moves to the next step.
Explicit analysis solves for the local domain variables using the below equation,
`u(n+1)=u(n)+ ∆T×F(u(n))`
Here, to compute the required parameter, in this case, the displacement, at the current time step i.e u(n+1), the explicit solver method makes use of the solution obtained from the previous time step i.e F(u(n))
EXPLICIT VS IMPLICIT ANALYSIS - A PRACTICAL COMPARISON:
OBJECTIVE:
To perform a Non-linear static analysis on a bar element constrained on one side while subjected to an external load on the other end using an Incremental load scheme using both Implicit and Explicit method and then solve them numerically. These results are then compared with the exact ideal solution using MATLAB to determine the most feasible method.
EXPLANATION:
Let us consider a bar element constrained on all DOFs on one side whereas the opposite side is un-constrained. Consider a force acting on the bar in the axial direction causing Tension in the bar.
The internal force f(u) experienced on the bar is a non-linear (higher-order) function of displacement.
`f(u)=u^3+9u^2+4 u`
We will solve this non-linear static problem numerically using both Implicit and Explicit method to determine the most feasible and approximate solution.
NUMERICAL SOLUTION:
Given, Internal force, `f(u)=u^3+9u^2+4 u` -----(1)
We can compute the stiffness from the above equation as,
`k= (df)/(du)=3u^2+18u+4` -----(2)
Wkt, `∆F=k×∆u` -----(3)
Where, F -> External load ΔF-> Load Increment
f(u) -> Internal force k -> Stiffness
u -> Displacement Δu -> Displacement Increment
EXPLICIT METHOD:
In the explicit method, we implement the “Incremental Explicit Load Scheme” where three equal steps of load increment of ΔF = 1 are used to load the bar in tension.
Our objective is to compute or solve for the displacement caused due to each individual load step and the internal force exerted by the bar to attain equilibrium. Hence, for this case of Non-linear static problem, we compare the internal and external forces which when are equal cause the system to be in static condition i.e in equilibrium.
Step 1: At first load step, ΔF = 1, u0 = 0
From (2), `k(u_0 )=3(0)+18(0)+4 = 4`
From (3), `∆u_1= (∆F)/k= 1/4=0.25`
Therefore,`u_1= u_0+ ∆u_1`
Where, u0 -> Initial displacement at no load
u1 -> Final displacement after 1st step load
Δu1 -> Change in displacement caused due to step load, ΔF=1
`u1 = 0+0.25=0.25`
Step 2: At second load step, ΔF = 1, u1 = 0.25
From (2),`k(u_1 )=3(0.25)^2+18(0.25)+4 = 8.687`
From (3), `∆u_2= (∆F)/(k(u_1))= 1/8.687=0.115`
Therefore,`u_2= u_1+ ∆u_2`
`u1 = 0.25+0.115=0.365`
Step 3: At third load step, ΔF = 1, u1 = 0.25
From (2),`k(u_2 )=3(0.365)^2+18(0.365)+4 = 10.97`
From (3), `∆u_3= (∆F)/(k(u_2))= 1/10.97=0.091`
`u_3= u_2+ ∆u_3`
`u1 = 0.365+0.091=0.456`
Hence, Total external force, F = ΔF+ ΔF+ ΔF = 1+1+1
F = 3
Internal force,`f=f(u_3 )=f(0.456)`
`f(u_3 )= 〖(0.456)〗^3+9×〖(0.456)〗^2+4×(0.456)`
`f(u_3 )=3.79`
Therefore, `F ≠f`
Hence, the system is not in equilibrium.
|
Step (i) |
ΔFi |
ui |
(F external)i |
(F internal)i |
R = f - F |
|
1 |
1 |
0.25 |
1 |
1.578 |
0.578 |
|
2 |
1 |
0.365 |
2 |
2.707 |
0.707 |
|
3 |
1 |
0.456 |
3 |
3.79 |
0.79 |
RESULT:
This deviation of the value of Internal force from the External force is called the Residue or error value. This high deviation in result implies Explicit analysis is a poor solver method.
However, we only computed the values from 3 iterations, this could be the reason why the values did not converge. However, upon performing the same analysis with Implicit solver method for about, say 20 Iterations, we can see a significant change in values with the residues decreasing as the iteration progresses.
IMPLICIT METHOD:
In the implicit method, we implement the “Incremental Explicit Load Scheme” where three equal steps of load increment of ΔF = 1 are used to load the bar in tension. However, unlike the explicit method, here we perform a Newton-Raphson iteration at the end of each time step to minimize the Residue and enforce equilibrium at the system before proceeding the next load step.
Wkt, Residue, `R(u)=f-F`
Upon expanding the above equation as a Taylor series, we get,
`R_(j+1)= R_j (u_j )+ (dR (u_j))/du 〖δu〗_(j+1)+⋯`
Hence ignoring higher-order terms, we get,
`R_(j+1)= R_j (u_j )+ (dR (u_j))/du 〖δu〗_(j+1) = R_j+k(u_j ) 〖δu〗_(j+1)`
To attain equilibrium,`R_(j+1)=0`
Hence solving for `〖δu〗_(j+1)`, we get,
`〖δu〗_(j+1)= 〖-[K(u_j )]〗^(-1)×R_j`
NOTE: In this Implicit method, we need to reduce the Residue below a tolerance level, because we cannot eliminate the residue completely using few Newton-Raphson iterations done numerically.
Hence, in this case, let us assume the tolerance limit as`10^(-2)`.
Step 1: At first load step, ΔF = 1, u0 = 0
From (2),`k(u_0 )=3(0)+18(0)+4 = 4`
From (3), `∆u_1= (∆F)/k= 1/4=0.25`
Therefore, `u_1= u_0+ ∆u_1`
Where, u0 -> Initial displacement at no load
u1 -> Final displacement after 1st step load
Δu1 -> Change in displacement caused due to step load, ΔF=1
`u1 = 0+0.25=0.25`
Internal force, `f(u_1 )= 〖(0.25)〗^3+9×〖(0.25)〗^2+4×(0.25) = 1.578`
Check Residue: `R(u)=f-F=1.578-1=0.578`
`R^((0))> 10^(-2)` -> Perform Newton-Raphson Iteration
Newton-Raphson:
Iteration 1: `u_1= 〖u_1〗^((0))=0.25`
`〖δu〗^((1) )= 〖-[K(u_j )]〗^(-1)×R_j`
`〖δu〗^((1) )= 〖-[8.68]〗^(-1)×0.578= -0.0665`
`〖u_1〗^((1))=〖u_1〗^((0))+ 〖δu〗^((1) )= 0.25+(-0.0665)=0.1835`
`f (〖u_1〗^((1)))= 1.0432`
Check Residue: `R(u)=f-F=1.0432-1=0.0432`
`R^((1))> 10^(-2)` -> Perform Newton-Raphson Iteration
Iteration 2: `〖u_1〗^((1))=1.0432`
`〖δu〗^((2) )= 〖δu〗^((1) ) 〖-[K(〖u_1〗^((1)) )]〗^(-1)×R^((1))`
`〖δu〗^((2) )=-0.0665 〖-[26.042]〗^(-1)×0.0432= -0.068159`
`〖u_1〗^((2))=〖u_1〗^((0))+ 〖δu〗^((2) )= 0.25+(-0.068159)=0.18184`
`f (〖u_1〗^((2)))= 1.0310`
Check Residue: `R(u)=f-F=1.0310-1=0.0310`
`R^((1))> 10^(-2)` -> Perform Newton-Raphson Iteration
Iteration 3: `〖u_1〗^((2))=0.18184`
`〖δu〗^((3) )= 〖δu〗^((2) ) 〖-[K(〖u_1〗^((2)) )]〗^(-1)×R^((2))`
`〖δu〗^((3) )=-0.068159 〖-[7.3723]〗^(-1)×0.0310= -0.072364`
`〖u_1〗^((3))=〖u_1〗^((0))+ 〖δu〗^((3) )= 0.25+(-0.072364)=0.17764`
`f (〖u_1〗^((3) ))= 1.0002`
Check Residue: `R(u)=f-F=1.0002-1=0.0002`
`R^((1))< 10^(-2)`
Therefore, `u_1=〖 u_1〗^((3))=0.17764`
Step 2: At second load step, ΔF = 1, u1 = 0.17764
From (2),`k(u_1 )=3(0.17764)+18(0.17764)+4 = 7.28`
From (3), `∆u_2= (∆F)/k= 1/7.28=0.13736`
`u_2= u_1+ ∆u_2`
`u_2=0.177+0.13736=0.31436`
Internal force, `f(u_2 )= 〖(0.31)〗^3+9×〖(0.31)〗^2+4×(0.31) = 2.1779`
Check Residue: `R(u)=f-F=2.1779-2=0.1779`
`R^((0))> 10^(-2)` -> Perform Newton-Raphson Iteration
Newton-Raphson:
Iteration 1: `u_2= 〖u_2〗^((0))=0.31436`
`〖δu〗^((1) )= 〖-[K(u_j )]〗^(-1)×R_j`
`〖δu〗^((1) )= 〖-[9.9549]〗^(-1)×0.1779= -0.017871`
`〖u_2〗^((1))=〖u_2〗^((0))+ 〖δu〗^((1) )= 0.31436+(-0.017871)=0.29649`
`f (〖u_2〗^((1) ))= 2.0032`
Check Residue: `R(u)=f-F=2.0032-1=0.0032`
`R^((1))< 10^(-2)`
Therefore, `u_2=0.29649`
Step 3: At third load step, ΔF = 1, u2 = 0.29649
From (2), `k(u_2 )=3( 0.29649)+18( 0.29649)+4 = 9.6005`
From (3), `∆u_3= (∆F)/k= 1/9.6005=0.10416`
`u_3= u_2+ ∆u_3`
`u3 = 0.29649+0.10416=0.40065`
Internal force, `f(u_3 )= 〖(0.40065)〗^3+9×〖(0.40065)〗^2+4×(0.40065) = 3.1116`
Check Residue: `R(u)=f-F=3.1116-3=0.1116`
`R^((0))> 10^(-2)` -> Perform Newton-Raphson Iteration
Newton-Raphson:
Iteration 1: `u_3= 〖u_3〗^((0))=0.40065`
`〖δu〗^((1) )= 〖-[K(u_j )]〗^(-1)×R_j`
`〖δu〗^((1) )= 〖-[11.693]〗^(-1)×0.1116= -0.009544`
`〖u_3〗^((1))=〖u_3〗^((0))+ 〖δu〗^((1) )= 0.40065+(-0.0095440)=0.39111`
`f (〖u_3〗^((1) ))= 3.0009`
Check Residue: `R(u)=f-F=3.0009-3=0.0009`
`R^((1))< 10^(-2)`
Therefore, `u_3=0.39111`
|
Step (i) |
ΔFi |
ui |
(F external)i |
(F internal)i |
R = f - F |
|
1 |
1 |
0.17764 |
1 |
1.0002 |
0.0002 |
|
2 |
1 |
0.29649 |
2 |
2.0032 |
0.0032 |
|
3 |
1 |
0.39111 |
3 |
3.0009 |
0.0009 |
RESULT:
Hence in all 3 load steps by performing Newton-Raphson iterations we were able to achieve a static equilibrium condition.
COMPARISON:
Let us compare our obtained results using MATLAB to get a clear understanding of how each of these methods correlates with each other and with the exact solution.
The following MATLAB command can either be executed in a .m script or in the Command window.
>> u1 = [0,0.25,0.365,0.456]; %Explicit method
>> u2 = [0,0.17764,0.29649,0.39111] %Implicit method
>> F = [0,1,2,3];
>>subplot(2,1,1),plot(u1,F,'ro',u,F,'r-.'),title('Explicit Method'),xlabel('Displacement (u)'),ylabel('Force (F)'),...
subplot(2,1,2),plot(u2,F,'g^',u1,F,'g--'),title('Implicit Method'),xlabel('Displacement (u)'),ylabel('Force (F)')
The above code will generate 2 subplots with the values obtained using numerical methods for both the explicit and implicit analysis and plots it with the external force applied to the bar in increments.
Hence, the resultant plot will be,
Similarly, on plotting them together, we get,
>> plot(u1,F,'r-.',u,F,'ko','linewidth',2),hold on,plot(u2,F,'g-',u1,F,'k^','linewidth',2),...
title('IMPLICIT VS EXPLICIT ANALYSIS'),xlabel('Displacement (u)'),ylabel('Force (F)')
However, these plots don't tell us the entire picture. Therefore, to get a clear comparison on which method provides the most approximate solution, we can plot the obtained values and plot along with the exact solution.
But before plotting we need to find the solution of the given equation,
>> syms u
>> f(u) = [u^3+9*u^2+4*u];
>> pretty(f)
3 2
u + 9 u + 4 u
Now we can substitute the obtained displacement values in the above function to get the resultant internal forces. Since the system is in equilibrium the net force i.e external and internal must be equal. This then can be plotted to get the exact solution.
>>d = [u1 u2];
>>d = sort(d,2)
d =
0 0 0.1776 0.2500 0.2965 0.3650 0.3911 0.4560
>>f(u) = subs(f,d);
To plot all three graphs,
>> plot(u1,F,'r-.',u1,F,'ko','linewidth',1.2),hold on,plot(u2,F,'g--',u2,F,'k^','linewidth',1.2),hold on,plot(ans,f,'k-','linewidth',1.2),...
legend('Explicit method','','Implicit method','','Exact solution'),title('IMPLICIT VS EXPLICIT ANALYSIS'),xlabel('Displacement (u)'),ylabel('Force (F)')
CONCLUSION:
From the above graph, we can clearly see that the results obtained from Implicit method almost overlaps with the exact solution. The small deviation that we can observe is due to our criterion on a tolerance limit of 0.01. We can also observe that the explicit method results are far off from the desired solution. However, as stated before, if we increase the number of iterations performed in the explicit analysis we can see the results converge to the exact solution.
Hence, we can conclude that the Implicit method provides a much more approximate solution compared to the solution obtained by the explicit method. It suffices to say that this is because the implicit method checks for a convergence criterion at each load increments and hence minimizes the persisting residue using the Newton-Raphson method before moving onto the next step.