Fortran Source Code & Offline Testing

Pipe Flow & Friction Head Loss

Core Numerical Engine in Fortran 90 • 0 total downloads

pipe_flow_calculator.f90

! =========================================================================
! Source File: pipe_flow_calculator.f90
! =========================================================================

program pipe_flow_calculator
    implicit none
    
    double precision :: D, L, roughness, rho, mu, nu
    double precision :: Q_flow, U_m, Re, f, f_old, dP, dP_kPa
    double precision :: A_c, pi, V_dot
    double precision :: head_loss, pump_power
    double precision :: f_lam, f_turb_smooth, f_turb_rough
    double precision :: rel_rough
    integer :: i, n_iter, converged
    character(len=30) :: flow_regime
    
    ! Minor loss variables
    integer :: n_fittings
    double precision :: K_total, dP_minor, head_minor
    
    pi = 3.14159265358979d0
    
    ! Read inputs
    read(*,*) D          ! Pipe inner diameter [m]
    read(*,*) L          ! Pipe length [m]
    read(*,*) roughness  ! Surface roughness eps [m]
    read(*,*) rho        ! Fluid density [kg/m3]
    read(*,*) mu         ! Dynamic viscosity [Pa.s]
    read(*,*) V_dot      ! Volume flow rate [m3/s]
    read(*,*) K_total    ! Sum of minor loss coefficients
    
    ! Computed properties
    nu = mu / rho
    A_c = pi * D**2 / 4.0d0
    U_m = V_dot / A_c
    Re = rho * U_m * D / mu
    rel_rough = roughness / D
    
    ! ============================================
    ! FLOW REGIME DETERMINATION
    ! ============================================
    if (Re < 2300.0d0) then
        flow_regime = 'Laminar'
    else if (Re < 4000.0d0) then
        flow_regime = 'Transitional'
    else
        flow_regime = 'Turbulent'
    end if
    
    ! ============================================
    ! FRICTION FACTOR CALCULATION
    ! ============================================
    
    ! Laminar friction factor
    f_lam = 64.0d0 / Re
    
    ! Turbulent: Colebrook-White equation (iterative)
    ! 1/√f = -2.0 log1₀(eps/D/3.7 + 2.51/(Re√f))
    
    if (Re >= 2300.0d0) then
        ! Initial guess: Swamee-Jain approximation
        f = 0.25d0 / (log10(rel_rough/3.7d0 + 5.74d0/Re**0.9d0))**2
        
        converged = 0
        n_iter = 0
        do i = 1, 100
            f_old = f
            f = 1.0d0 / (-2.0d0 * log10(rel_rough/3.7d0 + 2.51d0/(Re*sqrt(f_old))))**2
            n_iter = i
            if (abs(f - f_old) / f < 1.0d-8) then
                converged = 1
                exit
            end if
        end do
        f_turb_rough = f
    else
        f = f_lam
        f_turb_rough = 0.0d0
        n_iter = 0
        converged = 1
    end if
    
    ! Smooth pipe (Petukhov)
    if (Re >= 2300.0d0) then
        f_turb_smooth = (0.790d0 * log(Re) - 1.64d0)**(-2.0d0)
    else
        f_turb_smooth = 0.0d0
    end if
    
    ! ============================================
    ! PRESSURE DROP (Darcy-Weisbach)
    ! ============================================
    ! deltaP = f x (L/D) x (rhoU2/2)
    dP = f * (L / D) * rho * U_m**2 / 2.0d0
    dP_kPa = dP / 1000.0d0
    head_loss = dP / (rho * 9.81d0)
    
    ! Minor losses
    dP_minor = K_total * rho * U_m**2 / 2.0d0
    head_minor = dP_minor / (rho * 9.81d0)
    
    ! Total
    pump_power = (dP + dP_minor) * V_dot
    
    ! ============================================
    ! OUTPUT
    ! ============================================
    write(*,'(A)') '============================================================'
    write(*,'(A)') '   PIPE FLOW CALCULATOR — DARCY-WEISBACH METHOD'
    write(*,'(A)') '============================================================'
    write(*,*)
    write(*,'(A)') '--- PIPE GEOMETRY ---'
    write(*,'(A,F12.2,A)')  '  Inner Diameter (D)      = ', D*1000, ' mm'
    write(*,'(A,F12.4,A)')  '  Pipe Length (L)         = ', L, ' m'
    write(*,'(A,ES12.4,A)') '  Surface Roughness (eps)   = ', roughness, ' m'
    write(*,'(A,ES12.4)')   '  Relative Roughness eps/D  = ', rel_rough
    write(*,'(A,ES12.4,A)') '  Cross-Section Area      = ', A_c, ' m2'
    write(*,*)
    write(*,'(A)') '--- FLUID PROPERTIES ---'
    write(*,'(A,F12.2,A)')  '  Density (rho)             = ', rho, ' kg/m3'
    write(*,'(A,ES12.4,A)') '  Dynamic Viscosity (mu)   = ', mu, ' Pa.s'
    write(*,'(A,ES12.4,A)') '  Kinematic Viscosity (nu) = ', nu, ' m2/s'
    write(*,*)
    write(*,'(A)') '--- FLOW CONDITIONS ---'
    write(*,'(A,ES12.4,A)') '  Volume Flow Rate (Q)    = ', V_dot, ' m3/s'
    write(*,'(A,F12.4,A)')  '  Volume Flow Rate        = ', V_dot*1000*60, ' L/min'
    write(*,'(A,F12.4,A)')  '  Mean Velocity (U)       = ', U_m, ' m/s'
    write(*,'(A,ES12.4)')   '  Reynolds Number (Re)    = ', Re
    write(*,'(A,A)')        '  Flow Regime             = ', trim(flow_regime)
    write(*,*)
    write(*,'(A)') '--- FRICTION FACTOR ---'
    write(*,'(A,ES12.6)')   '  f (Darcy)               = ', f
    write(*,'(A,ES12.6)')   '  f_laminar (64/Re)       = ', f_lam
    if (Re >= 2300.0d0) then
        write(*,'(A,ES12.6)')   '  f_smooth (Petukhov)     = ', f_turb_smooth
        write(*,'(A,ES12.6)')   '  f_rough (Colebrook)     = ', f_turb_rough
        write(*,'(A,I4,A)')     '  Colebrook iterations    = ', n_iter, ' (converged)'
    end if
    write(*,*)
    write(*,'(A)') '--- MAJOR LOSSES (friction) ---'
    write(*,'(A,F12.2,A)')  '  Pressure Drop (deltaP)     = ', dP, ' Pa'
    write(*,'(A,F12.4,A)')  '  Pressure Drop           = ', dP_kPa, ' kPa'
    write(*,'(A,F12.4,A)')  '  Head Loss (h_f)         = ', head_loss, ' m'
    write(*,*)
    write(*,'(A)') '--- MINOR LOSSES (fittings) ---'
    write(*,'(A,F12.4)')    '  Total K coefficient     = ', K_total
    write(*,'(A,F12.2,A)')  '  Minor Pressure Drop     = ', dP_minor, ' Pa'
    write(*,'(A,F12.4,A)')  '  Minor Head Loss         = ', head_minor, ' m'
    write(*,*)
    write(*,'(A)') '--- TOTAL RESULTS ---'
    write(*,'(A,F12.2,A)')  '  Total Pressure Drop     = ', dP + dP_minor, ' Pa'
    write(*,'(A,F12.4,A)')  '  Total Pressure Drop     = ', (dP + dP_minor)/1000, ' kPa'
    write(*,'(A,F12.4,A)')  '  Total Head Loss         = ', head_loss + head_minor, ' m'
    write(*,'(A,F12.4,A)')  '  Required Pump Power     = ', pump_power, ' W'
    write(*,*)
    
    ! Moody chart data points
    write(*,'(A)') '--- FRICTION FACTOR vs Re (for Moody chart) ---'
    write(*,'(A)') '  Re            f_smooth      f_rough       f_laminar'
    write(*,'(A)') '  ---'
    
    do i = 1, 30
        ! Logarithmic spacing from Re=500 to Re=1e7
        Re = 10.0d0**(2.7d0 + dble(i-1) * 4.3d0 / 29.0d0)
        
        f_lam = 64.0d0 / Re
        
        ! Colebrook for rough pipe
        f = 0.25d0 / (log10(rel_rough/3.7d0 + 5.74d0/Re**0.9d0))**2
        do n_iter = 1, 50
            f_old = f
            f = 1.0d0 / (-2.0d0 * log10(rel_rough/3.7d0 + 2.51d0/(Re*sqrt(f_old))))**2
            if (abs(f - f_old) / f < 1.0d-8) exit
        end do
        f_turb_rough = f
        
        ! Smooth
        f_turb_smooth = (0.790d0 * log(Re) - 1.64d0)**(-2.0d0)
        
        if (Re < 2300.0d0) then
            write(*,'(ES12.4,2X,F12.8,2X,F12.8,2X,F12.8)') Re, f_turb_smooth, f_turb_rough, f_lam
        else
            write(*,'(ES12.4,2X,F12.8,2X,F12.8,2X,A)') Re, f_turb_smooth, f_turb_rough, '    ---'
        end if
    end do
    
    write(*,*)
    write(*,'(A)') '--- EQUATIONS USED ---'
    write(*,'(A)') '  Darcy-Weisbach: deltaP = f(L/D)(rhoU2/2)'
    write(*,'(A)') '  Colebrook:      1/√f = -2log1₀(eps/D/3.7 + 2.51/(Re√f))'
    write(*,'(A)') '  Minor losses:   deltaP_minor = ΣK(rhoU2/2)'
    write(*,'(A)') '  Head loss:      h_f = deltaP/(rhog)'
    
end program pipe_flow_calculator

Solver Description

Calculates pressure drop and friction head loss in circular pipes. Evaluates velocity ($V$) and Reynolds number ($Re_D$). Iteratively solves the implicit Colebrook-White equation for friction factor ($f$) using a Newton-Raphson numerical scheme, and computes head loss ($h_f$) via the Darcy-Weisbach equation.

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 -ffree-line-length-none pipe_flow_calculator.f90 -o pipe_flow_calc

Execution Command:

Execute the program by feeding the sample input file into the program using stdin redirection:

pipe_flow_calc < input.txt

📥 Downloads & Local Files

Preview of the required input file (input.txt):

! Pipe Inner Diameter ($D$) [mm]
100.0
! Pipe Length ($L$) [m]
100.0
! Surface Roughness ($\varepsilon$) [mm]
0.15
! Fluid Density ($\rho$) [kg/m³]
1000.0
! Fluid Viscosity ($\mu$) [Pa-s]
0.001
! Volumetric Flow Rate ($Q$) [m³/s]
0.02

📦 Download Offline Bundle (ZIP) 0 dl 📄 Download input.txt 0 dl 📥 Download pipe_flow_calculator.f90 0 dl

💻 Fortran Source Code Library

Pipe Flow & Friction Head Loss

Solver Description

Key Numerical Methods & Architecture

🛠️ Local Compilation

📥 Downloads & Local Files