💻 Fortran Source Code Library
We currently offer 172 open-source, production-grade Fortran codes for offline testing. Run calculations locally on your own machine, view code structure, read technical explanations, and download compilation packages including sample input files.
Fluid Statics Calculator
Core Numerical Engine in Fortran 90 • 44 total downloads
fluid_statics.f90
! =========================================================================
! Source File: fluid_statics.f90
! =========================================================================
program fluid_statics
implicit none
! Gravitational acceleration [m/s2]
double precision, parameter :: g = 9.80665d0
integer :: mode
! Mode A variables (Layered Fluids)
double precision :: P0, z_interest, P_final
integer :: N_layers
double precision, allocatable :: rho_layer(:), h_layer(:), p_interface(:)
double precision :: current_depth, accum_pressure
integer :: i, k
! Mode B variables (Submerged Surfaces)
integer :: shape
double precision :: dim1, dim2, dim3, dim4, dim5
double precision :: theta, hc, rho_fluid, P_atm
double precision :: theta_rad, area, Ixx_c, yc, ycp, hcp, FR_gauge, FR_abs
double precision :: shape_centroid_dist
! Mode C variables (Manometers)
integer :: N_left, N_right
double precision, allocatable :: rho_L(:), h_L(:)
double precision, allocatable :: rho_R(:), h_R(:)
double precision :: rho_m, dh, sum_L, sum_R, delta_P
! Read mode
read(*,*,iostat=i) mode
if (i /= 0) then
print *, "ERROR: Invalid calculation mode."
stop
end if
select case(mode)
case(1)
! ==========================================
! MODE 1: HYDROSTATIC PRESSURE (LAYERED)
! ==========================================
read(*,*,iostat=i) P0
if (i /= 0) then
print *, "ERROR: Failed to read surface pressure P0."
stop
end if
read(*,*,iostat=i) N_layers
if (i /= 0 .or. N_layers <= 0 .or. N_layers > 20) then
print *, "ERROR: Invalid number of layers."
stop
end if
allocate(rho_layer(N_layers))
allocate(h_layer(N_layers))
allocate(p_interface(0:N_layers))
do i = 1, N_layers
read(*,*,iostat=k) rho_layer(i), h_layer(i)
if (k /= 0 .or. rho_layer(i) < 0.0d0 .or. h_layer(i) <= 0.0d0) then
print *, "ERROR: Invalid layer data for layer", i
stop
end if
end do
read(*,*,iostat=i) z_interest
if (i /= 0 .or. z_interest < 0.0d0) then
print *, "ERROR: Invalid depth of interest."
stop
end if
! Compute interface pressures
p_interface(0) = P0
accum_pressure = P0
do i = 1, N_layers
accum_pressure = accum_pressure + (rho_layer(i) * g * h_layer(i)) / 1000.0d0 ! in kPa
p_interface(i) = accum_pressure
end do
! Compute pressure at z_interest
P_final = P0
current_depth = 0.0d0
do i = 1, N_layers
if (z_interest <= current_depth + h_layer(i)) then
! Target lies in current layer
P_final = p_interface(i-1) + (rho_layer(i) * g * (z_interest - current_depth)) / 1000.0d0
exit
else
current_depth = current_depth + h_layer(i)
if (i == N_layers) then
! Depth is beyond all layers, extrapolate with bottom layer density
P_final = p_interface(N_layers) + (rho_layer(N_layers) * g * (z_interest - current_depth)) / 1000.0d0
end if
end if
end do
! Print results
print *, "--- RESULTS MODE A ---"
print '(A,F12.4,A)', "Surface Pressure (P0) = ", P0, " kPa"
print '(A,F12.4,A)', "Target Depth = ", z_interest, " m"
print '(A,F12.4,A)', "Pressure at Target Depth = ", P_final, " kPa"
print *
print *, "--- INTERFACE PRESSURE DISTRIBUTION ---"
print '(A12, A15, A18)', "Boundary", "Depth [m]", "Pressure [kPa]"
current_depth = 0.0d0
print '(I12, F15.3, F18.4)', 0, current_depth, p_interface(0)
do i = 1, N_layers
current_depth = current_depth + h_layer(i)
print '(I12, F15.3, F18.4)', i, current_depth, p_interface(i)
end do
case(2)
! ==========================================
! MODE 2: SUBMERGED PLANE SURFACES
! ==========================================
read(*,*,iostat=i) shape
if (i /= 0 .or. shape < 1 .or. shape > 4) then
print *, "ERROR: Invalid shape code (1=Rect, 2=Circ, 3=Tri, 4=Trap)."
stop
end if
read(*,*,iostat=i) dim1, dim2, dim3, dim4, dim5
if (i /= 0) then
print *, "ERROR: Failed to read dimensions."
stop
end if
read(*,*,iostat=i) theta
if (i /= 0 .or. theta < 0.0d0 .or. theta > 90.0d0) then
print *, "ERROR: Inclination angle must be in [0, 90] degrees."
stop
end if
read(*,*,iostat=i) hc
if (i /= 0 .or. hc < 0.0d0) then
print *, "ERROR: Centroid depth must be non-negative."
stop
end if
read(*,*,iostat=i) rho_fluid
if (i /= 0 .or. rho_fluid <= 0.0d0) then
print *, "ERROR: Fluid density must be positive."
stop
end if
read(*,*,iostat=i) P_atm
if (i /= 0 .or. P_atm < 0.0d0) then
P_atm = 101.325d0 ! Default standard atm in kPa
end if
! Calculate shape properties
select case(shape)
case(1)
! Rectangle: dim1 = width b, dim2 = height h
if (dim1 <= 0.0d0 .or. dim2 <= 0.0d0) then
print *, "ERROR: Rectangle dimensions must be positive."
stop
end if
area = dim1 * dim2
Ixx_c = (dim1 * dim2**3) / 12.0d0
case(2)
! Circle: dim1 = Radius R
if (dim1 <= 0.0d0) then
print *, "ERROR: Circle radius must be positive."
stop
end if
area = 3.141592653589793d0 * dim1**2
Ixx_c = (3.141592653589793d0 * dim1**4) / 4.0d0
case(3)
! Triangle: dim1 = base b, dim2 = height h
if (dim1 <= 0.0d0 .or. dim2 <= 0.0d0) then
print *, "ERROR: Triangle dimensions must be positive."
stop
end if
area = 0.5d0 * dim1 * dim2
Ixx_c = (dim1 * dim2**3) / 36.0d0
case(4)
! Trapezoid: dim1 = top base a, dim2 = bottom base b, dim3 = height h
if (dim1 <= 0.0d0 .or. dim2 <= 0.0d0 .or. dim3 <= 0.0d0) then
print *, "ERROR: Trapezoid dimensions must be positive."
stop
end if
area = 0.5d0 * (dim1 + dim2) * dim3
Ixx_c = (dim3**3 * (dim1**2 + 4.0d0 * dim1 * dim2 + dim2**2)) / (36.0d0 * (dim1 + dim2))
end select
theta_rad = theta * 3.141592653589793d0 / 180.0d0
! Resultant Forces
! Force = pressure at centroid * Area
FR_gauge = (rho_fluid * g * hc * area) / 1000.0d0 ! in kN
FR_abs = ((rho_fluid * g * hc + P_atm * 1000.0d0) * area) / 1000.0d0 ! in kN
! Center of Pressure
if (theta < 0.001d0) then
! Horizontal surface
yc = hc
ycp = yc
hcp = hc
else
yc = hc / sin(theta_rad)
ycp = yc + Ixx_c / (yc * area)
hcp = ycp * sin(theta_rad)
end if
! Print results
print *, "--- RESULTS MODE B ---"
print '(A,I2)', "Shape Code = ", shape
print '(A,F12.4,A)', "Surface Area (A) = ", area, " m2"
print '(A,E14.6,A)', "Moment of Inertia (Ixx,c) = ", Ixx_c, " m4"
print '(A,F12.4,A)', "Inclination Angle (theta) = ", theta, " deg"
print '(A,F12.4,A)', "Centroid Depth (hc) = ", hc, " m"
print '(A,F12.4,A)', "Centroid distance (yc) = ", yc, " m"
print '(A,F12.4,A)', "Resultant Force (Gauge) = ", FR_gauge, " kN"
print '(A,F12.4,A)', "Resultant Force (Absolute) = ", FR_abs, " kN"
print '(A,F12.4,A)', "Center of Pressure (ycp) = ", ycp, " m"
print '(A,F12.4,A)', "Center of Pressure depth(h)= ", hcp, " m"
case(3)
! ==========================================
! MODE 3: MANOMETER DIFFERENTIAL PRESSURE
! ==========================================
read(*,*,iostat=i) N_left
if (i /= 0 .or. N_left < 0 .or. N_left > 10) then
print *, "ERROR: Invalid number of fluids in left leg."
stop
end if
allocate(rho_L(N_left))
allocate(h_L(N_left))
do i = 1, N_left
read(*,*,iostat=k) rho_L(i), h_L(i)
if (k /= 0 .or. rho_L(i) < 0.0d0 .or. h_L(i) < 0.0d0) then
print *, "ERROR: Invalid left leg fluid data at index", i
stop
end if
end do
read(*,*,iostat=i) N_right
if (i /= 0 .or. N_right < 0 .or. N_right > 10) then
print *, "ERROR: Invalid number of fluids in right leg."
stop
end if
allocate(rho_R(N_right))
allocate(h_R(N_right))
do i = 1, N_right
read(*,*,iostat=k) rho_R(i), h_R(i)
if (k /= 0 .or. rho_R(i) < 0.0d0 .or. h_R(i) < 0.0d0) then
print *, "ERROR: Invalid right leg fluid data at index", i
stop
end if
end do
read(*,*,iostat=i) rho_m, dh
if (i /= 0 .or. rho_m <= 0.0d0) then
print *, "ERROR: Invalid manometer fluid density."
stop
end if
! Calculate stages
sum_L = 0.0d0
do i = 1, N_left
sum_L = sum_L + (rho_L(i) * g * h_L(i)) / 1000.0d0 ! in kPa
end do
sum_R = 0.0d0
do i = 1, N_right
sum_R = sum_R + (rho_R(i) * g * h_R(i)) / 1000.0d0 ! in kPa
end do
! Delta P = P_left - P_right
! P_left + sum_L = P_right + sum_R + rho_m*g*dh
! Delta P = sum_R + (rho_m*g*dh)/1000 - sum_L
delta_P = sum_R + (rho_m * g * dh) / 1000.0d0 - sum_L
! Print results
print *, "--- RESULTS MODE C ---"
print '(A,F12.4,A)', "Left Leg Fluids head = ", sum_L, " kPa"
print '(A,F12.4,A)', "Right Leg Fluids head = ", sum_R, " kPa"
print '(A,F12.4,A)', "Manometer deflection head = ", (rho_m * g * dh) / 1000.0d0, " kPa"
print '(A,F12.4,A)', "Pressure Difference (DP) = ", delta_P, " kPa"
print '(A,F12.4,A)', "Pressure Difference (mH2O) = ", delta_P / (9.80665d0), " mH2O"
case default
print *, "ERROR: Calculator mode not recognized."
end select
end program fluid_statics
Solver Description
Computes hydrostatic pressure profiles in single or multi-layer liquid systems, resultant gauge and absolute forces and centers of pressure on submerged shapes (rectangular, circular, triangular, and trapezoidal), and multi-fluid U-tube manometer pressure differentials.
Key Numerical Methods & Architecture
- Input Redirection: Reads parameters sequentially from standard input (`stdin`) using Fortran sequential read (`read(*,*)`), ensuring modular integration.
- Modular Design: Formulated using pure mathematical routines, separation of equations from output formatting, and precise numerical solvers (e.g. bisection, Newton-Raphson).
- Standard Compliant: Written in clean, standards-compliant Fortran 90 to ensure cross-compiler compatibility.
🛠️ Local Compilation
To test this code on your machine, compile the source code file(s) using a standard Fortran compiler (e.g., `gfortran`).
Compilation Command:
gfortran -O3 fluid_statics.f90 -o fluid_statics_calc
Execution Command:
Execute the program by feeding the sample input file into the program using stdin redirection:
fluid_statics_calc < input.txt
📥 Downloads & Local Files
Preview of the required input file (input.txt):
! Calculation mode (1=Layered, 2=Surfaces, 3=Manometer)
1
! Surface pressure [kPa]
101.325
! Number of fluid layers
2
! Density and thickness of layers
800.0 2.0
! Depth of interest [m]
1000.0 3.0
! Parameter 6
2.5
1
! Surface pressure [kPa]
101.325
! Number of fluid layers
2
! Density and thickness of layers
800.0 2.0
! Depth of interest [m]
1000.0 3.0
! Parameter 6
2.5