Global Journal of Engineering Sciences (GJES)
Efficient
2nd Order Taylor Series Expansion Model to Integrate Stiff and Non-Linear ODEs:
Comparison to Rk45-Type Integrators
Authored by Ridha Djebali
Abstract
A
Taylor Series Expansion (TaSE) and three Runge-Kutta RK45 embedded pairs
numerical integrators of ordinary differential equations (ODEs) are developed.
The TaSE model is designed to solve higher order (>1) ODEs. Validation of
TaSE was carried out by help of the RK45 models as well as exact solution for
well-chosen stiff test cases. The model implementation and its local and global
order of accuracy are investigated.
Keywords and Phrases: Taylor Series Expansion
model; RK45 integrators; ODEs; Accuracy; Stiffness; Validation
Introduction
Historically,
it first started with Ludwig Prandtl in 1904 [1] when he gave a seminar at the
International Congress of Mathematics in Heidelberg, under the title “Uber
Flussigkeitsbewegungen bei sehr kleiner Rei bung”. Prandtl took the first steps
in a branch of science which turns out to be double-edged; on the one hand for
the fact of being branched into several physical areas and on the other hand
for the challenges encountered mathematically (existence and uniqueness of the
solution) then numerically (solvency, accuracy etc.). Prandtl first initiated
the theory of boundary layers problems (BLP) and laid its foundations providing
later a great research interest for the several scientific and engineering
applications such as the calculation of the drag coefficient of a lifting
airfoil or hydrofoil airplane wing, turbine blade [2]. A little bit after,
Blasius [3] gave in 1908 his famous non-linear boundary value problem (BVP)
which, among others, has been used until today as good test case to the new
numerical methods and models of integration of ordinary differential equations
(ODEs) used in the following for the resolution of partial differential
equations (PDEs).
The
investigation of differential equations may be done using the symbolic analytic
methods (exact), the geometric methods as well as the numerical methods.
Besides, Blasius problem is a BVP, posing mathematical and numerical
difficulties to its resolution that should only be numerically. Weyl [4] was
the first to prove the existence and uniqueness of the solution of the Blasius
problem where an accurate determination of the so-called shear stress f′′(0)
annoyed researchers. In fact, the Blasius BVP should be transformed to initial
value problem (IVP) to complete the explicit and missing equation to the
starting closure of the problem. Blasius gave bounds of the value between
0.3315 and 0.33175 and later the value was determined with sufficient precision
as f′′(0) =0.332057336. using, among other algorithms, the shooting method and
its well recommended variants for BVPs [5-12] as well as the finite differences
and collocation methods. Parallel to the continuous appearence of new schemes
and families of numerical methods for solving ODEs [13-16], ubiquitous stiff
problems are emerging namely in electronics, solid mechanics, plasma physics,
chemical
Foot Notes
2010
Mathematics Subject Classification. 93C10, 93C15, 93C20, 30K05
kinetics,
biology etc. The application of known methods exhibits some instabilities and
accuracy loss. To overcome problems of stability, accuracy and ODEs stiffness
due to the presence of different time scales in the handled problem it is
inevitably to develop new high order flexible and accurate integrators. Many
attempts are made to extend the notion of differential equations and the range
of available methods. The classification of solution methods to ODEs and
systems of ODEs may be performed differently based on some criteria and
properties, namely explicit embedded and implicit methods (especially for the
Runge-Kutta RK method variants), onestep, two-steps or multisteps methods
[17-19].
Over
the years, many numerical methods as listed above but of different forms,
namely Adams-Bash forth, Adams-Moulton or Multistep, are used on the same
principle. Because of its high performance at relatively little computation
cost, the RK4 finds great use in wide range of applications even originally has
been founded as fixed-step integrator. The KR4 integrator idea holds wide range
of interest and gained continuous development [20- 23] and application which
lead to several improvements such as providing adaptive step size. The latter
is a good property that should exert a good integrator over its own
time-marching; so that the step size h changes continuously to reduce the
computational cost while keeping the total error within given tolerance ε. Each
‘new’ step size is called the optimal step size subject the calculation
methods. The integrators using this idea are often called embedded pairs
method, such as the RK45 variants: RK-Cash Karp method (RK45-CK), RK-Fehlberg
method (RK45-F), RK-Dormand-Prince method (RK45-DP) etc.
However,
all numerical methods have their own potential, upsides and downsides when
solving various types of ODEs and development of new solution techniques is
important. In the present work a Taylor Series Expansion (TaSE designed) model
is developed to solve higher order (>1) ODEs and system of ODEsand was
compared to three home coded RK45 models as well, namely the RK45-CK, RK45-F
and RK45-DP. The models are compared on second and three order ODEs test cases
presenting different stiffness scales. The TaSE model is designed to solve
higher oders ODEs and system of ODEs and is expected to be at least second
order solver since it uses the explicit trapezoidal rule
The problem is very
stiff for very low values of the parameters δ called the initial radius. An
investigated method should exhibit sufficient aptitude to simulate the very
acute jump of the flame radius. In Figure 2, the two models seem to be
efficient to handle such problems.
The numerical results
of comparison of the TaSE and RK45-F models are given in Table 1. By checking
the order of the TaSE model (local error not global) gives that is of second
order accuracy. Verification can be performed on the results for the ability to
capture a peak, say for example for η = 7 and for different step sizes h =
0.0005, h = 0.005 and h = 0.05. Assuming the truncation error Δf ∼ O(ha) = bha. The plot of ln(Δf) vs. ln(h) is
given in
Figure 3. The fit
gives ln(Δf) = 1.975 ln(h) + 0.95, so the local deviation under the prescribed
conditions (ζ = 5, η = 7) is Δf ≈ 2.4647h1.966. The fit was obtained within R2
= 0.99965 confidence.
Note that the results
of Table 1 are obtained for ζ = 5, TaSE model: h = 0.5 10−3, RK45-F and
RK45-CK: fixed step h = 10−3 and ε = 10−5 and that the RK45-F integrator is
found to be slightly more fast than RK45-CK one regarding the extra number of
iterations taken sometimes for adaptive stepsize.
The obtained results
for TaSE and RK45-F models gathered with literature results are given in Table
2. The results of the proposed TaSE model are in excellent agreement with those
of the RK45-F integrator either for f and f′ even for quite coarse step-size h
= 0.01; which reduces considerably the number of iterations. The computed
quadratic (RMS) deviation is close to 1.35554E-05 for f and 6.17489E-06 for f′
which supports well on the accuracy of the proposed model.
The present example
shows how to estimate x′(t + h) to second order in a second order ODE. The
problem is solved for t∈[0,7] with step size h
= 0.001.
The results of TaSE
and RK45-DP models are gathered in Table 3 with their deviations from the exact
solution. Extra digits are given here since the comparison is performed against
exact solution. The plot of the two computational results is depicted in Figure
5. Results of Table 3 and Figure 5 reveal the powerfulness of the TaSE model.
Besides, in this test case we have reviwed the global order of accuracy of the
TaSE model. The RMS deviation over the seven computed times was calculated for
three stepsizes: h = 0.001, h = 0.005 and h = 0.01. It has been found that
ln(Δf) = 1.998 ln(h)+4.445 with a coefficient of determination R2 =1. The
expected order of accuracy is always preserved.
Conclusion
A Taylor Series Expansion
numerical integrator of ODEs is developed. The model is exerted to solve higher
order ODEs. The following concluding remarks may be drawn:
• The model has been validated
on three fourth order accurate RK models and exact solution and gave excellent
agreement for three test cases of second and three order ODEs of differents
stiffness scales.
• The local and global order
of accuracy of the model was checked and estimated to be second order accurate.
• The model is always
efficient even for first approximation of the first derivative for second order
ODEs.
• The Trapezoidal explicit
rule and the theorem of the mean value are used when needed to handle to second
order very complex and hard-to-derive terms.
For
more about Iris Publishers
please click on: http://irispublishersgroup.com/submit-manuscript.php.
To read more about this article - https://irispublishers.com/gjes/fulltext/efficient-2nd-order-taylor-series-expansion-model.ID.000647.php
Indexing List of Iris Publishers: https://medium.com/@irispublishers/what-is-the-indexing-list-of-iris-publishers-4ace353e4eee
Iris
publishers google scholar citations: https://scholar.google.co.in/scholar?hl=en&as_sdt=0%2C5&q=irispublishers&btnG=

Comments
Post a Comment