This site delves into the connection between mechanical engineering, physics, and mathematics using numerical methods and computer algebra systems (CAS). It presents straightforward articles that explore both analytical and numerical approaches to solving real-world engineering and scientific problems.
Showing posts with label Fluid Mechanics. Show all posts
Showing posts with label Fluid Mechanics. Show all posts
The Mach number is a dimensionless
parameter in fluid dynamics that expresses the ratio between the flow velocity
and the local speed of sound. It is named in honor of Austrian physicist and
philosopher Ernst Mach.
$$M = \frac{u}{c}$$
where: \(M\) – Mach number \(u\) – local flow velocity with respect to
the boundaries \(c\) – speed of sound in the medium
The speed of sound depends on the square
root of the thermodynamic temperature. By definition, \(\text{Mach} \; 1\) corresponds to a
flow velocity \(u\) equal to the local speed of sound. At \(\text{Mach} \; 0.5\), \(u\) is half the
speed of sound (subsonic), while at \(\text{Mach} \; 2.0\), \(u\) is twice the speed of sound
(supersonic). Flow at \(\text{Mach} \; 5.0\) or higher is classified as hypersonic.
At standard atmospheric conditions – dry air
at mean sea level and a temperature of \(15^\circ C\) – the speed of sound is approximately \(340.3 \; m/s\) or \(1,225.1 \; km/h\). The speed of sound is not constant; in gases, it
increases in proportion to the square root of the absolute temperature. Since
atmospheric temperature generally decreases with altitude above sea level, the
speed of sound correspondingly decreases as well.
Flight can be roughly classified in six speed
categories:
When an object moves faster than the speed
of sound (Mach 1), a sudden change in air pressure forms in front of it. This
change, called a shock wave,
spreads backward in a cone shape known as the Mach
cone. The shock wave creates the sonic
boom we hear when a fast aircraft passes by. As the object goes
faster, the cone becomes narrower; just above Mach 1, it looks almost flat.
For further information on the Mach number,
refer to the article “Mach number” in Wikipedia:
The Free Encyclopedia.
Understanding Mach Number Through
Visualization
Figure 1 shows a fighter jet flying at a
supersonic speed of \(Mach \; 2.0\). The cloud that forms around the aircraft at this
speed is called a vapor
cone or shock
collar. It appears as a visible cloud of condensed water vapor
caused by a sudden drop in air pressure and temperature as the jet nears the
speed of sound.
This article provides a Scilab script to
visualize the relative speed of an object and the propagation of sound waves,
as shown in Figure 2, where the object moves at Mach 2.0. The dashed
circle represents the propagation of the shock wave from its center and the
object’s previous position.
For comparison, Figure 3 illustrates an
object moving at the speed of sound along with the resulting shock wave
propagation.
Figure 3. Sonic Speed Simulation (Mach
1.0).
Scilab Code for Figure 2.
// Copyright (C) 2025 - Gani Comia// Date of creation: 19 Sep 2025// Script: Supersonic Speed Simulation and Mach Coneclear;clc;
// primary parameters
M =2; // Mach number (supersonic > 1)
S =1; // speed of sound
mu =asind(S./M); // angle of mach cone
dist =1:10; // object position by a factor of speed of sound
N =length(dist); // number of positions// secondary parameters
a =1:10// object position at x-direction
b =0// object position at y-direction// visualization of object's position and shock wave propagationclf;
g =gcf()
g.figure_size=[700,700]for i =1:N
plot(dist,0,"ko","linewidth",5)
theta =1:360
r(i)=sind(mu).*(N-i)
x(i,:)= r(i)*cosd(theta)+ a(i);
y(i,:)= r(i)*sind(theta)+ b;
plot(x(i,:),y(i,:),"r--")
a1 =legend(["Moving Object","Shock Wave"],with_box=%f)
a1.font_size=2.5endtitle("Supersonic Speed Simulation at Mach 2.0 in 1D", "fontsize",4)xlabel("Object Position / Shock Wave Propagation","fontsize",3.5)ylabel("Shock Wave Propagation","fontsize",3.5)xgrid(color("grey"),1,7)
a2 =xstring(6.5,4,"https://gani-mech-toolbox.blogspot.com")
a2.font_size=2.5// mach cone plot and visualization
x_val =[a(1)+r(1).*cosd(90-mu) a($)]
y_val =[r(1).*sind(90-mu)0]
slope =diff(y_val)./diff(x_val)disp(slope)// straight line using slope-intercept formula, y = mx + b, and// solve for y-intercept, b, at x = 0
y_int = y_val(1)- slope.*x_val(1)disp(y_int)
x_val_new =[0 a($)]
y_val_new_1 =[y_int 0]
y_val_new_2 =[-y_int 0]plot(x_val_new,y_val_new_1,"b-","linewidth",3)plot(x_val_new,y_val_new_2,"b-","linewidth",3)
a3 =xstring(x_val(1),y_val(1),"Mach Cone",mu-6)
a3.font_size=3
a4 =xstring(x_val(1),-y_val(1)-0.7,"Mach Cone",mu-56)
a4.font_size=3
ax =gca()
ax.data_bounds=[0-5 ; 126];
//disp(r)
Feel free to comment for inquiry,
clarification, correction or suggestion for improvement. Drop your email to make
a request to the author.
Disclaimer: The formulas and calculations
presented are for technical reference only. Users must verify the accuracy and
ensure compliance with applicable engineering standards, codes, and safety
requirements before practical application.
Topic: Equilibrium of Force Systems and Wind Loads
Subject: Fluid and Engineering Mechanics
Tool: QCAD & Scilab
By: Gani Comia, May 2025
Gantry Crane Structure and Equilibrium of Force Systems
Gantry cranes are travelling cranes designed for lifting and
moving heavy loads. This material handling equipment is generally used outdoors
where it is not convenient to erect an overhead runway. Its structure consists
of a girder or bridge, freestanding legs, electric motor that drives the steel
wheel, and a rail system on the ground. The
bridge or girder is carried at the ends by the legs supported by trucks with
wheels so that the crane can travel. The bridge carries a hoisting assembly
unit. The crane is driven by the motor through a gear reduction shaft. [1]. The motor is referred to as a brake
motor because it prevents movement caused by inertia when it is deenergized.
The objective of this article is to present an engineering
simulation method of analyzing the effect of a wind loads due to typhoon to an outdoor
gantry crane. The assumption is that there is a particular value of wind velocity
that will make the crane to tip over at the wind direction.
Figure 1 shows a simple illustration of an outdoor gantry
crane. On its side view is the free-body diagram (FBD) showing the
corresponding concurrent forces acting on the structure. The following
conventions will be used in the analysis: - \(W\) is the weight of the crane; \(F_{g h}\) is the horizontal induced wind load on the girder; \(F_{s h}\) and \(F_{s v}\) are the
horizontal and vertical wind load acting on the support legs respectively; and \(R_{a c}\) and \(R_{b d}\) are the reaction forces at the four wheels, \(a\), \(b\), \(c\), and \(d\).
Figure 1. Outdoor Gantry Crane and Free-Body Diagram for the
Forces System.
The solution is to determine as to whether the crane will be
tipping over at the center of rotation, \(R_{a c}\), due to moment on the direction of
the wind induced forces. Figure 2 illustrates two cases of wind velocity, \(V = 0\) and \(V > 0\). For the given domain of \(V\), there is a reaction force, \(R_{b d}\), where
in there is a resultant moment from the coplanar force system that determines
the rotation of crane from the center.
Figure 2. FDB of Gantry Crane due to Effect of Wind
Velocity.
Case #1 is a condition where the wind velocity is \(V = 0\). With
this there will be no induced forces and the force system present on FBD are
crane weight, \(W\), and the two reaction forces, \(R_{a c}\) and \(R_{b d}\), on the wheel
supporting the structure. In this condition, the equilibrium for concurrent
force systems is obtained by determining the equations that produce a zero
resultant. The resultant will be zero and equilibrium will exist when the
following equations are satisfied:
$$\sum F_x = 0 \tag{1}$$
$$\sum F_y = 0 \tag{2}$$
Where:
\(F_x\) – forces parallel to x-axis
\(F_y\) – forces parallel to y-axis
With the two conditions of equilibrium,
reaction forces can be determined using Equation (3).
$$R_{a c} = R_{b d} = \frac{W}{2} \tag{3}$$
Case #2 is the condition wherein the wind
velocity acting on the crane is inducing resultant forces that make the
structure to move by rotation. In this condition, the equilibrium in terms of
moment summations about the center on its line of action can be used to
calculate the resultant force. Hence the two other equations of equilibrium are
$$\sum M_A = 0 \tag{4}$$
$$\sum M_B = 0 \tag{5}$$
Where:
\(M_A\) – moment at center A
\(M_B\) – moment at center B
For the case #2 on Figure 2, reaction \(R_{a c}\) is used as the center of moment. For the increasing wind velocity, \(V_{wind}\),
there exist a corresponding increase in wind load represented by \(F_{g h}\), \(F_{s h}\),
and \(F_{s v}\) as they are functions of wind velocity in Equation 6.
Wind is a mass of air that moves in a mostly horizontal
direction from an area of high pressure to an area with low pressure. High
winds can be very destructive because they generate pressure against the
surface of a structure. The intensity of this pressure is the wind load. The wind
effect is dependent upon the size and shape of the structure. Calculating wind
load is necessary for the design and construction of safer, more wind-resistant
building and placement of objects such as antennas on top of the buildings.
The generic formula for wind load is
$$F = A \; P \; C_d \tag{8}$$
Where:
\(F\) – wind load or force, \(N\)
\(A\) – projected area of the object, \(m^2\)
\(P\) – wind pressure, \(N/m^2\)
\(C_d\) – drag coefficient
The generic formula is used for
estimating the wind load on a specific object and it does not meet the building
code requirements for planning a new construction.
The formula for wind pressure, \(P\), in SI
units is
$$P = \text{0.613} \; V^2 \tag{9}$$
Where:
\(V\) – wind speed, \(m/s\)
The wind speed
categorizes by the Philippine Atmospheric, Geophysical and Astronomical Services
Administration (PAGASA) will be used as basis in calculation. PAGASA is the
Philippine government agency responsible for providing weather forecasts,
typhoon warnings, and other meteorological and astronomical services. Currently
the Tropical Cyclone Wind Signal (TCWS) has been in use based on the adoption
of best practices from other tropical cyclone warning centers.
The following are the corresponding wind
threat based on the wind speed value.
TCWS #1: Wind Threat = 10.8 ~ 17.1 \(m/s\)
TCWS #2: Wind Threat = 17.2 ~ 24.4 \(m/s\)
TCWS #3: Wind Threat = 24.5 ~ 32.6 \(m/s\)
TCWS #4: Wind Threat = 32.7 ~ 51.2 \(m/s\)
TCWS #5: Wind Threat = 51.3 \(m/s\) or higher
\(C_d\) is the drag coefficient for the
object subjected under wind pressure. Drag is the force that air exerts on the
structure. Drag coefficient for the structure’s shape can be estimated roughly
using the following standard:
\(C_d\) = 0.8 for short cylinder
\(C_d\) = 1.2 for long cylinder
\(C_d\) = 1.4 for shorter flat plate
\(C_d\) = 2.0 for long flat plate
Wind Load Simulation for Gantry Crane
Figure 3 shows the effect of wind load on
the gantry crane structure in terms of the equilibrium of force
systems.
Figure 3. Wind Load Simulation using the Equilibrium
of Force System on Gantry Crane.
The graph’s domain is based on the TCWS converted
in terms of \(km/hr\). Corresponding wind speeds or velocities for the five
categories of TCWS are plotted. Wheel reaction load for \(R_{b d}\) are shown for each
TCWS. The red dash line which depicts \(R = 0\) is the critical wheel reaction load
for \(R_{b d}\) on the FBD in Figure 2. With this \(R = 0\), the critical wind velocity is
calculated to have a value of \(\text{127}\) \(km/hr\) under the TCWS #4. A horizontal wind speed which of \(V_{wind} >\) \(\text{127}\) \(km/hr\) acting perpendicularly to structure will result to tripping over the crane
to the wind direction.
With this, a conclusion can be made from
the simulation stated as follows: - a wind speed of over \(\text{127}\) \(km/hr\) which is
within the TCWS #4 will induce an equivalent critical wind load to affect the
crane stability. The use of a device in gantry crane known as “Storm Lock” is
highly recommended to increase its capacity for critical wind load.
Below is the Scilab script to generate
the simulation plot on Figure 3. Parameters can be edited for the other cases.
Scilab Script
// Copyright (C) 2025 - Gani Comia// Date of creation: 24 Apr 2025// Wind Load Effect Analysis on Outdoor Gantry Crane
clear;clc;
// gantry crane primary parameters for girders
lG =1.35+22+0.35; // m, crane girder length
hG =1; // m, crane girder height
W =14.5; // Tonne, crane weight
shapeG =2; // crane girder as long flat plate// gantry crane primary parameters for support legs
lS =5.7; // m, support leg length
dS =0.25; // m, support leg diameter
shapeS =1.2; // leg support as long cylinder tube// wind velocity domain
Vmps =0:52; // m/s, wind speed, tropical cyclone
Vkph = Vmps*3.6; // km/hr, wind speed, tropical cyclone// modules or functions// wind load analysis using generic formula// structure surface areafunctionA=area(l, h)// calculate girder area perpendicular to wind velocity// input argument:// l - surface length (m)// h - surface height (m)// ouput argument:// A - surface area (m^2)A=l.*h;
endfunction// wind pressurefunctionP=windPressure(V)// calculate the induced wind pressure// input argument:// V - wind velocity (m/s)// ouput argument:// P - wind pressure (N/m^2 or Pa)P=0.613*V.^2;
endfunction// drag coefficientfunctionCd=dragCoefficient(shape)// calculate the theortical drag coefficient// input argument: structure shape// shape = 2 for long flat plate// shape = 1.4 for shorter flat plate// shape = 1.2 for long cylinder tube// shape = 0.8 for short cylinder tube// output argument:// Cd - drag coefficient (no unit)Cd=shape;
endfunction// induced wind loadfunctionFw=windLoad(A, P, Cd)// calculate the induced wind load// input argument:// A - surface area (m^2)// P - wind pressure (N/m^2)// Cd - drag coefficient (no unit)// output argument:// Fw - wind load (N)Fw=A.*P.*Cd;
endfunction// main function// wind load calculation at the girder
A_g =area(lG,hG);
P_g =windPressure(Vmps);
Cd_g =dragCoefficient(shapeG);
Fgh =windLoad(A_g,P_g,Cd_g);
Fgh = Fgh/9806.65; // Tonne, wind load // wind load calculation at the support legs
A_s =area(lS,dS); // m^2, surface area// angle of support leg
phi = atand(((3.5/2)-0.168)/5.495); // degrees, angle
Vmps_s = Vmps./cosd(phi); // m/s, wind speed perpendicular to legs
P_s =windPressure(Vmps_s);
Cd_s =dragCoefficient(shapeS);
Fsr =windLoad(A_s,P_s,Cd_s);
// force component at legs
Fsh = Fsr*cosd(phi); Fsh =3*Fsh/9806.65; // N to tonnes
Fsv = Fsr*sind(phi); Fsv =3*Fsv/9806.65; // N to tonnes// critial load
Rbd =(W*(3.5/2)+Fsv*((3.5/2)+0.937)-Fgh*(6.25+0.375)-Fsh*(2.9+0.375))/3.5;
// Plotting calculation results
clf;
fig = gcf()
fig.figure_size=[700,700]
plot(Vkph,Rbd,"b-","linewidth",4)
note ="$\LARGE R_{wheel} = f(V_{wind})$"
legend(note,with_box=%F)
plot([0,200],[0,0],"r--","linewidth",1.8)
title("Wind Load Effect on Outdoor Gantry Crane","fontsize",4)
xlabel("Wind Velocity, Tropical Cyclone, V, (kph)","fontsize",3.5)
ylabel("Wheel Reaction Load, R, (Tonne)","fontsize",3.5)
xgrid(color("green"),1,7)
xstring(10,0,"Critical Wheel Reaction Load = 0")
xstring(120,5.0,"https://gani-mech-toolbox.blogspot.com")// Wind velocity limit of tropical cyclone
Vwind =[10.8,17.2,24.5,32.7,51.3]; // m/s, wind velocity
Vwind = Vwind*3.6; // kph, wind velocity // tropical cyclone signal number
nVwind = length(Vwind)for i =1:nVwind
R(i)= interp1(Vkph,Rbd,Vwind(i))
vLx =[Vwind(i) Vwind(i)]; vLy =[-10 R(i)]
hLx =[0 Vwind(i)]; hLy =[R(i) R(i)]
plot(vLx,vLy,"k--")
plot(hLx,hLy,"k--")
plot(Vwind(i),R(i),"marker","d","markerFaceColor","red")
xstring(Vwind(i)+8,-9.8,["Signal No. ",string(i)],-90)end// critical wind velocity
format(5)
Rcrit =0
Vcrit = interp1(Rbd,Vkph,Rcrit)
plot(Vcrit,Rcrit,"marker","S","markerFaceColor","red")
xstring(Vcrit,Rcrit,["Critical Velocity (kph) =",string(Vcrit)])
Feel free to comment for inquiry,
clarification, correction or suggestion for improvement. Drop your email to make
a request to the author.
Disclaimer: The formulas and calculations
presented are for technical reference only. Users must verify the accuracy and
ensure compliance with applicable engineering standards, codes, and safety
requirements before practical application.
References
E.A. Avallone and T. Baumeister III. Marks Standard Handbook
for Mechanical Engineers. 9th Ed. McGraw-Hill International Edition.
1987.
F.L. Singer. Engineering Mechanics. 2nd Ed.
Harper International Edition, New York, Evanston & London. 1970.
The size of the pump for the particular application is
decided base on the volumetric flow (\(Q\)) and the total dynamic head (\(TDH\)). The pump’s
motor drive and specification in terms of \(HP\) can be just determined from the pump
maker’s catalogue for the given flow rate and \(TDH\).
Volumetric flow is usually an available information based on
application. The \(TDH\) is calculated considering the following factors:
Pressure Head
Static or Elevation Head
Velocity Head
Friction Head
For this presentation, the calculation of \(TDH\) will be based
on the first three factors of Bernoulli’s Principle. The theoretical \(TDH\) can be adjusted by adding the total friction head and other factors.
Bernoulli’s principle is a fundamental concept in fluid
dynamics that describes the relationship of pressure, elevation, and velocity
of a fluid flow. It calculates the Total Dynamic Head (\(TDH\)) in pump systems. \(TDH\) is the total energy imparted by the pump to the fluid and calculated as
shown:
$$TDH\;=\;H_p\;+\;H_z\;+\;H_v \tag{1}$$
For incompressible fluids, TDH can be expressed in the simplified form of Bernoulli’s equation:
Fig. 1 shown a graph for the relationship of pressure and \(TDH\) for a given volumetric flow and different static or elevation heads. It can be a
helpful guide for a pump designer to make a decision understandable to the
stakeholders or clients. This also eliminates the impression of unjustifiable guesswork.
Fig. 1. Pump Sizing in terms of \(TDH\) based on pressure,
static, and velocity head.
Presented here is the Scilab script as a toolbox to produce a
similar graph for pump sizing using Bernoulli’s principle.
// Copyright (C) 2024 - Gani Comia// Date of creation: 16 Dec 2024// Water Pump Sizing using Bernoulli's Principle (Theory)
clear;clc;
// Pumps for fire protection// Required Pressure = 138~1034 kPa (20~150 psi)// Volumetric Flow = 22.7 m^3/hr (100 GPM)// Required: Pump's TDH// sub-routine program—total dynamic headfunctionH=tdh(P, z, V)
wDensity =9800; // N/m^3, weight density of water
g =9.8; // m/s^2, gravitational acceleration// P - kPa, pressure// z - m, elevation// V - m/s, velocityH=(1000*P)./wDensity +z+(V.^2)./(2*g)endfunction// main program// (1) velocity head, calculation
Q =23; // m^3/hr, volumetric flow rate
Q = Q/3600; // m^3/s, volumetric flow rate// pipe size at service point dia 38mm GI Pipe S40
pipeID =(38-2*1.5)/1000; // m, inside diameter of pipe
A =(%pi/4)*(pipeID.^2); // m^2, pipe ID area
V = Q./A; // m/sec, velocity// (2) elevation head, assume a range for simulation
z =0:10:50; // m, elevation range
Nz = length(z)// no. of elevation steps
disp(z)// (3) pressure head, assigned a domain for simulation// P = 300 ~ 700 kPa pressure range
P = linspace(300,700,Nz); // kPa, pressure range// (4) simulation plot
clf;
for i =1:Nz
H =tdh(P,z(i),V)
plot(P,H,"linewidth",1.75)
title("H2O Pump Sizing using the Bernoulli Principle")
xlabel("Pressure, P [ kPa ]")
ylabel("Total Dynamic Head, TDH [ m ]")end
P1 =600;
for i =1:Nz
H1 =tdh(P1,z(i),V)
note1 =["Z =",string(z(i)),"m"]
xstring(P1,H1,note1,-18)end// (5) simulation of sample parameters// for elevation = 20 m and pressure = 500 kPa
P2 =500; // kPa, sample pressure
H2 =tdh(P2,z(3),V); // m, sample tdh
x1 =[P2 P2];
y1 =[30 H2];
plot(x1,y1,"--r")// vertical line
xstring(P2,32,["P =",string(P2)],-90)
x2 =[300 P2];
y2 =[H2 H2];
plot(x2,y2,"--r")// horizontal line
format(6);
xstring(350,H2,["TDH =",string(H2)])
plot(P2,H2,"o","color","black")
note2 =["Q =",string(3600*Q),"m3 / hr"]
xstring(350,120,note2)
note3 =["https://gani-mech-toolbox.blogspot.com"]
xstring(350,115,note3)// custom markers
a = gca();
a.children.children.mark_background=9; // Fill color
This kind of pump sizing or selection can
simplify the intricacies of engineering analysis, especially for the
stakeholders with minimal or no knowledge at all in fluid mechanics.
Feel free to comment for inquiry,
clarification, or suggestion for improvement. Drop your email to request the
soft copy of the file.
Disclaimer: The formulas and calculations
presented are for technical reference only. Users must verify the accuracy and
ensure compliance with applicable engineering standards, codes, and safety
requirements before practical application.
Topic: Visual Plot of Pump's Pressure Head and Volumetric Flow;
Subject: Fluid Mechanics;
Tool: Euler Math Toolbox;
by: Gani Comia, Nov. 2024;
In one scenario, a 100 GPM volumetric flow of water is being required in a building's fire protection system. Your facility manager is asking if it can be achieved from existing < 70 GPM and if the pressure setting can be adjusted from the current 100 psi to something like 160 psi. As a technical resource person in mechanical, you need to provide an explanation with the supporting calculation and visual simulation in layman's terms.
Fig. 1. Visual representation of the pump's characteristics.
Feel free to comment for inquiry, clarification, or suggestion for improvement. Drop your email to request the softcopy of the file.
Disclaimer: The formulas and calculations
presented are for technical reference only. Users must verify the accuracy and
ensure compliance with applicable engineering standards, codes, and safety
requirements before practical application.
