All Courses
All Courses
Courses by Software
Courses by Semester
Courses by Domain
Tool-focused Courses
Machine learning
POPULAR COURSES
Success Stories
Finite volume method The Finite Volume Method (FVM) is a method of representing and evaluating partial differential equations in the form of algebraic equations. In the finite volume method volume integrals in partial differential equations that contain a divergence term are converted to surface integrals using the divergence…
Faizan Akhtar
updated on 19 Mar 2021
Finite volume method
The Finite Volume Method (FVM) is a method of representing and evaluating partial differential equations in the form of algebraic equations. In the finite volume method volume integrals in partial differential equations that contain a divergence term are converted to surface integrals using the divergence theorem. These terms are then evaluated as fluxes at the surfaces of each finite volume because the flux leaving the particular volume is identical to the flux entering an adjacent volume, these methods are conservative.
It can be easily formulated for unstructured mesh where the axes are not mutually orthogonal to each other, which makes this method highly popular in industries and academics.
Differences between FDM and FVM
S No | Finite Difference Method | Finite Volume Method | ||||
1 | It is evaluated at the nodes | It is evaluated at the volume surrounding the nodes | ||||
![]() |
![]() |
|||||
2 | It is well suited for structured mesh | The process is much more efficient for the unstructured mesh (complex geometry, arbitrary shape). | ||||
3 | It can be applied to the unstructured mesh but the process becomes tedious. | The discretization process looks much simpler and easier when FVM is applied to complex geometry shapes. | ||||
4 | It is computationally more expensive for the unstructured mesh. | It is computationally less expensive for the unstructured mesh | ||||
5 | It has failed to preserve the concept of conservation. | It preserves the concept of conservation. | ||||
6 | It is highly undesirable in the case of moving meshes. | It is desirable in the case of moving meshes. |
Application of finite volume method to 1-D steady-state heat conduction problem
The governing equation for one-dimensional steady-state heat conduction equation with source term is given as
ddx(α∗dTdx)+S=0
where 'T' is the temperature of the rod. The boundary values of temperature at A and B are prescribed.
α is thermal diffusivity.
S is the source term.
The first step in the finite control volume is to discretize the domain into discrete control volumes as shown in figure
The number of nodal points is placed between A & B. The boundary of the control volume is positioned midway between the two adjacent nodes, thus each node is surrounded by the control volume. The control volume should be set up near the edges of the domain such that the physical boundary and the control volume boundary should coincide with each other.
A general nodal point is identified as P and its neighboring nodes to the west and east are identified as W and E respectively. The west side of the control volume is identified as w whereas the east side is identified as e. The distance between the nodes W and P, and between nodes P and E are identified as δxWPandδxPE respectively. Similarly between w and P, and between P and e are identified as δxwPandδxPe. The control volume width △x=δxwe.
The key step in the finite volume method is to integrate the governing equation over a control volume to yield a discretize equation at its nodal point P. For the control volume defined above this gives
∫△Vddx(α∗dTdx)∗dV+∫△VSdV=0
dV=A∗dx
The equation becomes
(α∗A∗dTdx)e−(α∗A∗dTdx)w+S△V=0
A is the cross sectional area of the control volume
â–³V is the control volume.
S is the average value of the source S over the control volume.
The above equation states that the diffusive flux of temperature leaving the east face of the control volume minus the diffusive flux of temperature entering the west face is equal to the generation of temperature. It constitutes a balance equation of temperature over a control volume.
In order to derive useful forms of the discretized equation, the interface diffusion coefficient α and temperature gradient dTdx at east and west face are required`. Linear approximations seem to be the obvious and simplest way of calculating interface values and gradients. This practice is called a central differencing scheme.
αw=αW+αP2
αe=αE+αP2
And the diffusive flux terms are evaluated as
(α∗A∗dTdx)e=αe∗Ae∗TE−TPδxPE
(α∗A∗dTdx)w=αw∗Aw∗TP−TWδxWP
The finite volume method approximates the source term as
S∗△V=Su+Sp∗TP
The final equation can be rearranged as
(αe∗Ae∗TE−TPδxPE)−(αw∗Aw∗TP−TWδxWP)+(Su+Sp∗TP)=0
This can be arranged as
(αe∗AeδxPE+αw∗AwδxWP−Sp)∗TP=(αw∗AwδxWP)∗TW+(αe∗AeδxPE)∗TE+Su
Identifying the coefficient of TP,TW,TEas aW,aP,aE and rearranging the equation as under
aP∗TP=aW∗TW+aE∗TE+Su
where
aW | aP | aE |
αw∗AwδxWP | αe∗AeδxPE | aW+aP−Sp |
Interpolation schemes
The approximation of the surface and volume integrals may require values of the variable at locations other than the computational nodes of the control volume (CV). Values at this location are obtained using interpolation formulae.
The interpolation schemes are listed below
Upwind biased interface values
v>0
ui−0.5≈ui−1
ui+0.5≈ui
Backward differencing
Ic=v∗ui−ui−1dx
Forward differencing
v<0
ui−0.5≈ui
ui+0.5≈ui+1
Ic=v∗ui+1−uidx
Taylor series expansion
vui+0.5=vui+v△x2⋅(∂u∂x)i+0.5+v(△x)28⋅(∂2u∂x2)i+0.5+....
a first order accurate flux approximation, the leading truncation error resembles a diffusive flux d∂u∂x with d=v△x2 being the numerical diffusion coefficient.
Interpolation polynomial
p1(x)=uLxR−xxR−xL+uRx−xLxR−xL
Average interface values
ui−0.5≈ui−1+ui2,ui+0.5≈ui+ui+12
Ic=vui+1−ui−12△x
ui+1=ui+0.5+△x2∗(∂u∂x)i+0.5+(△x)28(∂2u∂x2)i+0.5+......
ui=ui+0.5−△x2∗(∂u∂x)i+0.5+(△x)28(∂2u∂x2)i+0.5+......
ui+0.5=ui+ui+12−(△x)28⋅(∂2u∂x2)i+0.5 (second-order accuracy)
Quadratic upwind interpolation for convective kinematics
p2(x)=uLx−xMxL−xMxR−xxR−xL+uMx−xLxM−xLxR−xxR−xM+uRx−xLxR−xLx−xMxR−xM
Upwind biased interface values
v>0
ui−0.5≈3ui+6ui−1−ui−28
ui+0.5≈3ui+1+6ui−ui−18
Ic≈v3ui+1−7ui−1+3ui+ui−28△x
v<0
ui−0.5≈3ui−1+6ui−ui+18
ui+0.5≈3ui+6ui+1−ui+28
Ic≈−v3ui−1−7ui+1+3ui+ui+28△x
Flux limiters
They are used in high-resolution schemes to avoid spurious oscillations (wiggles) that would otherwise occur with high-order spatial discretization schemes due to shocks, discontinuities, or sharp changes in the solution domain. For smoothly changing solution flux limiters are not used and the spatial derivatives can be implemented with high-resolution schemes. If there are discontinuities in the domain the fux limiters change the solution to low and high-resolution schemes depending upon the gradients in the nearby surroundings.
Matlab representation of the equation
Case-1
Absence of source term in the heat conduction equation
ddx(K∗dTdx)=0
Integrating above the control volume dVand discretization of the above equation will yield
K⋅A⋅(TE−TPdx)−K⋅A⋅(TP−TWdx=0
2∗K∗Adx∗TP=K∗Adx∗TE+K⋅Adx∗TW
The above equation represents the solution of temperature at the inner nodes.
The node associated with the left boundary
2∗K∗Adx∗TP=K∗Adx∗TA+K⋅Adx∗TE
The node associated with the right boundary
2∗K∗Adx∗TP=K∗Adx∗TB+K⋅Adx∗TW
Matlab Code
% one dimensional heat conduction steady state
% finite volume method
% no source term in the governing equation
%governing equation=((d/dx(K*dT/dx))=0;
clear
close all
clc
n=25;
L=1;
x=linspace(0,L,n);
dx=x(2)-x(1);
% defining boundaries of temperature
T(1)=100;
T(n)=500;
T_old=T;
% exact solution
T_exact=400.*x+100;
% thermal conductivity
K=1000;
Area=1e-2;
q=100;
source_term=q*Area*dx;
error=1;
tol=1e-7;
iter_number=0;
while (error>tol)
for i=1:n
if (i>2) && (i<n-1)
% calculating the coefficient of temperature
a_p=2*K*Area/dx;
a_e=K*Area/dx;
a_w=K*Area/dx;
T(i)=(a_w*T(i-1)+a_e*T(i+1))/a_p;
end
if i==2
a_p=2*K*Area/dx;
a_e=K*Area/dx;
a_w=0;
T(i)=(a_e*T(i+1)+(K*Area*100)/dx)/a_p;
end
if i==n-1
a_p=2*K*Area/dx;
a_e=0;
a_w=K*Area/dx;
T(i)=(a_w*T(i-1)+(K*Area*500)/dx)/a_p;
end
end
% calculating error
iter_number=iter_number+1;
error=max(abs(T_old-T));
T_old=T;
end
plot(x,T)
hold on
plot(x,T_exact,'o')
title('Comparison of numerical solution with exact solution in the abscence of source term in the governing equation')
legend('Numerical solution','Exact solution')
ylabel('Temperature')
xlabel('Length of rod')
pause(0.03)
Output graph
The number of grid points=25
Case-2
Presence of source term in the heat conduction equation
ddx(K∗dTdx)+q=0
2∗K∗Adx∗TP=K∗Adx∗TE+K⋅Adx∗TW+q∗A∗dx
The above equation represents the solution of temperature at the inner nodes.
The node associated with the left boundary
2∗K∗Adx∗TP=K∗Adx∗TA+K⋅Adx∗TE+q∗A∗dx
The node associated with the right boundary
2∗K∗Adx∗TP=K∗Adx∗TB+K⋅Adx∗TW+q∗A∗dx
Matlab Code
% one dimensional heat conduction steady state
% finite volume method
% source term in the governing equation
%governing equation=((d/dx(K*dT/dx))+q=0;
clear
close all
clc
n=25;
L=2;
x=linspace(0,L,n);
dx=x(2)-x(1);
% defining boundaries of temperature
T(1)=100;
T(n)=500;
T_old=T;
% thermal conductivity
K=0.5;
Area=1;
q=1000;
source_term=q*Area*dx;
% exact solution
T_exact=((400)/(L)+(500/K)*(L-x)).*x+100;
error=1;
tol=1e-7;
iter_number=0;
while (error>tol)
for i=1:n
if (i>2) && (i<n-1)
% calculating the coefficient of temperature
a_p=2*K*Area/dx;
a_e=K*Area/dx;
a_w=K*Area/dx;
T(i)=(a_w*T(i-1)+a_e*T(i+1)+source_term)/a_p;
end
if i==2
a_p=2*K*Area/dx;
a_e=K*Area/dx;
a_w=0;
T(i)=(a_e*T(i+1)+(K*Area*100)/dx+source_term)/a_p;
end
if i==n-1
a_p=2*K*Area/dx;
a_e=0;
a_w=K*Area/dx;
T(i)=(a_w*T(i-1)+(K*Area*500)/dx+source_term)/a_p;
end
end
% calculating error
iter_number=iter_number+1;
error=max(abs(T_old-T));
T_old=T;
end
figure(2);
plot(x,T)
hold on
plot(x,T_exact,'o')
title('Comparison of numerical solution with exact solution in the presence of source term in the governing equation')
legend('Numerical solution','Exact solution')
ylabel('Temperature')
xlabel('Length of rod')
pause(0.03)
Output graph
The number of grid points=5
The number of grid points=10
The number of grid points=25
Conclusion
It can be inferred from the graph that the accuracy of the numerical solution increases with an increase in the number of nodes.
Leave a comment
Thanks for choosing to leave a comment. Please keep in mind that all the comments are moderated as per our comment policy, and your email will not be published for privacy reasons. Please leave a personal & meaningful conversation.
Other comments...
Week-4 : Basic Calibration of Single cylinder CI-Engine
Aim: Basic Calibration of Single cylinder CI-Engine Objective : Explore Tutorial No 1- CI final 1.Compare SI vs CI and list down differences (assignment no 2-SI) 2. Comments on the following parameters BSFC Exhaust Temperature A/F ratios 3.Change MFB 50 and observe impact on performance Introduction Difference…
11 Nov 2021 05:26 AM IST
Week 2 : Basic Calibration of Single cylinder SI-Engine
Aim: Basic Calibration of Single-cylinder SI-Engine Objective: 1. Run the case at 1800 rpm and list down important parameters (20 Marks) air flow rate BMEP BSFC In-cylinder pressure 2. Increase the power output at 3600 rpm by 10% (30 Marks) Introduction A spark-ignition engine (SI engine) is…
22 Oct 2021 07:11 AM IST
Week 1 : Exploring the GUI of GT-POWER
Aim: Exploring the GUI of GT-suite GT-suite can be used in a wide range of applications. There are several industries that are using GT-suite for various applications On-highway and off-highway vehicle Marine and rail Industrial machinery Aerospace Power generation Multiphysics platform GT-Suite has a versatile multiphysics…
12 Oct 2021 02:52 PM IST
Week 8: Literature review - RANS derivation and analysis
Aim: Literature review - RANS derivation and analysis Objective: Apply Reynolds decomposition to the NS equations and come up with the expression for Reynolds stress. Explain your understanding of the terms Reynolds stress What is turbulent viscosity? How is it different from molecular viscosity? Introduction…
01 Sep 2021 07:52 AM IST
Related Courses
0 Hours of Content
Skill-Lync offers industry relevant advanced engineering courses for engineering students by partnering with industry experts.
© 2025 Skill-Lync Inc. All Rights Reserved.