🔬
Solver Purpose & Physical Scope
Models turbulent submerged free jet hydrodynamics using Schlichting self-similarity and Tollmien potential core decay. Computes centerline velocity decay, radial velocity profiles, fluid entrainment factors, centerline stagnation dynamic pressures, nozzle reaction thrust, and normal impact forces on target surfaces.
📂 Discipline: Fluid
⚡ Precision: IEEE-754 64-bit Real(`real(8)`)
📥 Total Downloads: 254 times
📄 Source File:
free_jet_impact.f90
📁 calcul/FluidMechanics /
free_jet_impact.f90
program free_jet_impact
implicit none
integer :: i, iostat_val, n_points, vane_code
double precision :: D0_mm, U0, rho, mu, X_mm, D0_m, X_m
double precision :: A0, Q0, mdot0, Re0, F_thrust
double precision :: L_core_m, L_core_mm, u_max_target, decay_factor
double precision :: entrainment_ratio, Q_target, r_half_m, r_half_mm
double precision :: q_dyn_Pa, P_stag_kPa, impact_force_N
double precision :: x_pos_m, u_pos, r_curr_mm, eta, u_radial
double precision, parameter :: PI = 3.141592653589793d0
! Read inputs
read(*,*,iostat=iostat_val) D0_mm
if (iostat_val /= 0) then
write(*,*) 'ERROR: Invalid nozzle diameter.'
stop
end if
read(*,*,iostat=iostat_val) U0
read(*,*,iostat=iostat_val) rho
read(*,*,iostat=iostat_val) mu
read(*,*,iostat=iostat_val) X_mm
read(*,*,iostat=iostat_val) vane_code ! 1 = flat 90 deg, 2 = inclined 45 deg, 3 = Pelton 180 deg
if (iostat_val /= 0) then
write(*,*) 'ERROR: Failed to read all free jet inputs.'
stop
end if
! Validation
if (D0_mm <= 0.0d0 .or. U0 <= 0.0d0 .or. rho <= 0.0d0 .or. mu <= 0.0d0 .or. X_mm <= 0.0d0) then
write(*,*) 'ERROR: Geometry, velocity, and properties must be positive.'
stop
end if
! Conversions
D0_m = D0_mm * 1.0d-3
X_m = X_mm * 1.0d-3
A0 = (PI / 4.0d0) * (D0_m**2)
Q0 = A0 * U0
mdot0 = rho * Q0
Re0 = (rho * U0 * D0_m) / mu
F_thrust = mdot0 * U0
! Potential Core
L_core_m = 5.5d0 * D0_m
L_core_mm = L_core_m * 1.0d3
! Centerline Velocity at Target Distance X
if (X_m <= L_core_m) then
u_max_target = U0
entrainment_ratio = 1.0d0
else
decay_factor = 5.8d0 / (X_m / D0_m)
u_max_target = min(U0, U0 * decay_factor)
entrainment_ratio = max(1.0d0, 0.32d0 * (X_m / D0_m))
end if
Q_target = Q0 * entrainment_ratio
r_half_m = 0.097d0 * X_m
r_half_mm = max(D0_m / 2.0d0, r_half_m) * 1.0d3
q_dyn_Pa = 0.5d0 * rho * (u_max_target**2)
P_stag_kPa = q_dyn_Pa / 1.0d3
! Reaction Impact Force
if (vane_code == 3) then
! Pelton 180 deg bucket
impact_force_N = 2.0d0 * mdot0 * min(U0, u_max_target * 1.15d0) * 0.95d0
else if (vane_code == 2) then
! 45 deg plate
impact_force_N = mdot0 * min(U0, u_max_target * 1.15d0) * sin(45.0d0 * PI / 180.0d0)
else
! Flat 90 deg plate
impact_force_N = mdot0 * min(U0, u_max_target * 1.10d0)
end if
! Output Results
write(*,'(A)') '============================================================'
write(*,'(A)') ' THERMOFLUIDCALC — TURBULENT FREE JET & IMPINGEMENT ENGINE'
write(*,'(A)') '============================================================'
write(*,'(A,ES12.4,A)') 'Nozzle Exit Diameter D0 = ', D0_mm, ' mm'
write(*,'(A,ES12.4,A)') 'Initial Jet Velocity U0 = ', U0, ' m/s'
write(*,'(A,ES12.4)') 'Reynolds Number Re0 = ', Re0
write(*,'(A,ES12.4,A)') 'Mass Flow Rate (mdot) = ', mdot0, ' kg/s'
write(*,'(A,ES12.4,A)') 'Volumetric Flow Rate (Q0) = ', Q0, ' m3/s'
write(*,'(A)') '------------------------------------------------------------'
write(*,'(A,ES12.4,A)') 'Potential Core Length = ', L_core_mm, ' mm'
write(*,'(A,ES12.4,A)') 'Target Distance X = ', X_mm, ' mm'
write(*,'(A,ES12.4,A)') 'Centerline Velocity u_max = ', u_max_target, ' m/s'
write(*,'(A,ES12.4,A)') 'Jet Half-Width Radius r1/2= ', r_half_mm, ' mm'
write(*,'(A,ES12.4)') 'Fluid Entrainment Q(X)/Q0 = ', entrainment_ratio
write(*,'(A)') '------------------------------------------------------------'
write(*,'(A,ES12.4,A)') 'Stagnation Pressure = ', P_stag_kPa, ' kPa'
write(*,'(A,ES12.4,A)') 'Nozzle Thrust Reaction F0 = ', F_thrust, ' N'
write(*,'(A,ES12.4,A)') 'Target Impingement Force = ', impact_force_N, ' N'
write(*,'(A)') '============================================================'
write(*,*)
write(*,'(A)') '--- AXIAL CENTERLINE VELOCITY DECAY U_MAX(X) ---'
write(*,'(A)') ' x/D0 x_mm u_max(m/s)'
n_points = 20
do i = 1, n_points
x_pos_m = (D0_m * 0.5d0) + (35.0d0 * D0_m) * (dble(i - 1) / dble(n_points - 1))
if (x_pos_m <= L_core_m) then
u_pos = U0
else
u_pos = min(U0, U0 * (5.8d0 / (x_pos_m / D0_m)))
end if
write(*,'(2X,F10.2,2X,F10.2,2X,F10.2)') x_pos_m / D0_m, x_pos_m * 1.0d3, u_pos
end do
write(*,*)
write(*,'(A)') '--- SCHLICHTING RADIAL VELOCITY PROFILE AT TARGET X ---'
write(*,'(A)') ' r_mm u(r)(m/s)'
do i = -15, 15
r_curr_mm = (r_half_mm * 2.2d0) * (dble(i) / 15.0d0)
eta = abs(r_curr_mm) / max(0.1d0, r_half_mm)
u_radial = u_max_target * ((1.0d0 + 0.414d0 * (eta**2))**(-2.0d0))
write(*,'(2X,F10.2,2X,F10.2)') r_curr_mm, u_radial
end do
end program free_jet_impact
💻 How to Compile & Run Locally
1. Compilation (GNU Fortran / Intel oneAPI):
gfortran -O3 free_jet_impact.f90 -o free_jet_impact
2. Execution with input.txt redirection:
free_jet_impact < input.txt
📄 Sample input.txt File Structure
Sample Data:
25.0 30.0 998.0 1.002e-3 300.0 1
Parameter Description:
Nozzle Exit Diameter D0 [mm] Initial Jet Velocity U0 [m/s] Fluid Density [kg/m3] Fluid Viscosity [Pa-s] Target Distance X [mm] Target Geometry (1=Flat 90 deg, 2=Inclined 45 deg, 3=Pelton 180 deg)