3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
79
ERROR ESTIMATION OF BILINEAR GALERKIN FINITE
ELEMENT METHOD FOR 2D THERMAL PROBLEMS
S. M. Afzal Hoq
Department of Mechanical Engineering
International Islamic University Malaysia, Kualalumpur, (Malaysia).
E-mail: afzalhoqsu@gmail.com ORCID: https://orcid.org/0000-0002-5917-2507
Abdurahim Okhunov
Department of Science and Engineering
International Islamic University Malaysia, Kualalumpur, (Malaysia).
E-mail: abdurahimokhun@iium.edu.my ORCID: https://orcid.org/0000-0001-7092-5699
C. P. Tso
Faculty of Engineering and Technology
Multimedia University, (Malaysia).
E-mail: cptso@mmu.edu.my ORCID: https://orcid.org/0000-0003-3884-5992
Recepción:
03/02/2020
Aceptación:
14/02/2020
Publicación:
13/03/2020
Citación sugerida:
Hoq, S. M. A., Okhunov, A., y Tso, C. P. (2020). Error estimation of bilinear Galerkin nite element
method for 2D thermal problems. 3C Tecnología. Glosas de innovación aplicadas a la pyme, 9(1), 79-93. http://doi.
org/10.17993/3ctecno/2020.v9n1e33.79-93
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
80
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
ABSTRACT
This study demonstrates a two-dimensional steady state heat conduction Laplace partial dierential
equation solution using the bilinear Galerkin nite element method. Heat transfer analysis is of vital
importance in many engineering applications and devising computationally inexpensive numerical
methods while maintaining accuracy is one of the primary concerns. The method uses structured mesh
grid over a two-dimensional rectangular domain and solved using a stiness matrix for the bilinear
elements, calculated using the proposed modied numerical scheme. Several numerical experiments are
conducted by controlling the number of nodes and changing element sizes of the presented scheme, and
comparison made between analytical solution and software generated solution.
KEYWORDS
Galerkin method, Bilinear element, Heat conduction, Error analysis, Partial dierential equation.
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
81
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
1. INTRODUCTION
Finite element Analysis is a robust numerical approach for solving partial dierential equations (PDEs)
used in various branches of science such as solid mechanics, uid mechanics, electromagnetics,
thermodynamics etc. (Kreyszig, 2011; Deb, Babuška, & Oden, 2001). One of the most powerful
techniques for solving PDEs with weak formulation is using a weighed residual method called the Galerkin
nite element method (GFEM). The formulation requires generating a basis function (Ainsworth &
Oden, 1997) based on the elemental boundary conditions. This trial function is substituted in the
partial dierential equation and the rst derivative of the trial function is taken for each nodal variable
(Burden & Faires, 2001) to construct the residual function. The weighed form of the residual function
for the whole domain is integrated by setting it equal to zero. Green’s theorem can be applied over the
boundary if necessary (Afzal, Sulaeman, & Okhunov, 2016). The corresponding numerical model is set
up by discretizing the rectangular domain into smaller elements where each element consists of nodal
coordinates and nodal variables, which are used to perform the Galerkin approximation of the PDEs.
This is followed by the generation of an element matrix and vector matrix of the boundary by integrating
the total number of elements. The set of linear equations represented by the matrices are consecutively
solved using the Galerkin approach. For comparison purposes, an exact solution already available for
the non-linear PDE (time independent and no heat source) used in 2D heat conduction rectangular
domain with both Dirichlet and Neumann boundary condition, is presented. Finally, a stiness matrix
applicable for homogenous rectangular domain consisting of structured mesh grid elements is presented,
the solution scheme of which signicantly reduces the CPU performance cost.
2. MATHEMATICAL CONSTRUCTION
A homogeneous domain represented by a time independent heat conduction problem with zero heat
extraction is considered and can be mathematically formulated by Laplace equation with Neumann
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
82
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
or Dirichlet boundary conditions as shown in Figure 1 (Gerald & Wheatly, 2003; Gockenbach, 2002;
Burden & Faires, 2001; Kreyszig, 2011).
Figure 1. Physical domain of Ω bounded by Γ.
Two dimensional time independent heat-conduction problem can be represented by the following partial
dierential equation.
(1)
where,
and
The boundary conditions are given by:
on
(2)
where u
b
and f
b
are Dirichlet and Neumann boundary conditions respectively.
The weighed residual w is applied on Eqs (1) and (2) to generate its strong formulation as given:
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
83
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
(3)
Consequently, the weak formulation is generated using integration by part over eq. (3) as given below:
(4)
The stiness matrix is generated by discretising the solution domain into smaller elements and performing
integration for each element.
3. BASIS FUNCTION AND STIFFNESS MATRIX
The 2D computational domain is represented by bilinear rectangular elements where each element is
constructed with a node at each of its four edges as shown in Figure 2. The nodal points are represented
by (x
1
, y
1
), (x
2
, y
2
), (x
3
, y
3
) and (x
4
, y
4
), while the nodal variables are represented by u
1
, u
2
, u
3
and u
4
. The
nodal variable u at any given location (x,y) for a bilinear element is approximated by the basis function,
which is written as:
u=a
i
+a
i
x+a
i
y+ a
i
xy (5)
Figure 2. Bilinear element.
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
84
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
The above-mentioned basis function can be represented in matrix form:
(6)
Considering each bilinear element composed of four nodal points and nodal variables, the corresponding
basis function in matrix form is represented as:
(7)
The corresponding shape function is formulated as:
(8)
The individual shape function is based on the four geometric bilinear coordinates for each element
located in a structured mesh grid:
(9)
(10)
(11)
(12)
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
85
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
Where the rectangle Area, A is given by:
(13)
As it can be observed from Eq. (5) to (11), the trial function u is dependent on the nodal variables at the
corners of the bilinear nite elements and the shape functions. The rst derivative of this trial function
gives the weighed function w in Eq. (4)
(14)
In order to determine the elemental stiness matrix, the rst integral of the weak formulation in Eq. (4)
is carried out sequentially given below:
(15)
(16)
The matrix generated is a 4x4 matrix as given by:
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
86
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
(17)
The diagonal values as follows:
(18)
(19)
(20)
(21)
Since:
k
22
=k
11
, k
23
=k
14
, k
24
=k
13
, k
33
=k
11
, k
34
=k
12
and k
44
=k
11
the matrix is symmetric about its diagonal.
(22)
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
87
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
4. RESULTS AND DISCUSSION
In this section a set of numerical experiments are carried out for the evaluation and analysis of the
results. A rectangular 2D domain is developed and steady state heat conduction with no heat source is
applied. The domain measures 5 units and 10 units in the x and y direction respectively. GFEM is used
by implementing the weighed residual approximation approach on each bilinear element for solving
the Laplace PDE governing the 2D heat conduction problem. Both Dirichlet and Neumann boundary
conditions are specied for the solution domain as given in Figure 3. The boundary conditions are below:
Figure 3. Steady state heat conduction problem of rectangular domain.
Solving the above mentioned problem required calculating the stiness matrix of the bilinear elements
using the Galerkin nite element approach. The accuracy of the approach largely depends on the element
size and number of elements considered. It has been observed increasing the number of elements and
decreasing the element size increases the accuracy of the solution by achieving higher convergence rates.
3C Tecnología. Glosas de innovación aplicadas a la pyme. ISSN: 2254 – 4143 Ed. 33 Vol. 9 N.º 1 Marzo - Junio
88
http://doi.org/10.17993/3ctecno/2020.v9n1e33.79-93
Comparison is made with the exact solution of the Eq. (23) below and Figure 4 shows the graphical
analysis.
(23)
Figure 4. Comparison between Galerkin FE solution and exact solution.
The Figure 5 represents the solution domain discretized into a numbers of smaller bilinear elements
each consisting of nodes at the corners.
Figure 5. Physical domain with bilinear mesh.