💻 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.
Plate Heat Exchanger
Core Numerical Engine in Fortran 90 • 49 total downloads
hx_plate.f90
! =========================================================================
! Source File: hx_plate.f90
! =========================================================================
program hx_plate
implicit none
! Inputs
integer :: mode, Nplates
double precision :: T_hi, T_ho, T_ci, T_co
double precision :: m_h, m_c
double precision :: rho_h, Cp_h, mu_h, k_h, mu_wh
double precision :: rho_c, Cp_c, mu_c, k_c, mu_wc
double precision :: W, Lp, t_plate, b_depth, pc_pitch, beta_deg, d_port_mm
double precision :: dP_max_kpa
double precision :: Rfh, Rfc, kw
! Calculated properties
double precision :: X, phi, Dh, Q_h, Q_c, Q, dT1, dT2, LMTD_cf
double precision :: hs_h, hs_c, dP_h, dP_c, dP_port_h, dP_port_c, U_f, U_c, A_prov, A_req
double precision :: Re_h, Re_c, vch_h, vch_c
integer :: Nplates_curr, iostat_val
logical :: has_cross_error
double precision, parameter :: pi = 3.141592653589793d0
double precision :: dP_max
! Read all inputs from stdin
read(*,*,iostat=iostat_val) mode
if (iostat_val /= 0) then
write(*,*) 'ERROR: Invalid inputs'
stop
end if
read(*,*,iostat=iostat_val) Nplates
read(*,*,iostat=iostat_val) T_hi
read(*,*,iostat=iostat_val) T_ho
read(*,*,iostat=iostat_val) T_ci
read(*,*,iostat=iostat_val) T_co
read(*,*,iostat=iostat_val) m_h
read(*,*,iostat=iostat_val) m_c
read(*,*,iostat=iostat_val) rho_h
read(*,*,iostat=iostat_val) Cp_h
read(*,*,iostat=iostat_val) mu_h
read(*,*,iostat=iostat_val) k_h
read(*,*,iostat=iostat_val) mu_wh
read(*,*,iostat=iostat_val) rho_c
read(*,*,iostat=iostat_val) Cp_c
read(*,*,iostat=iostat_val) mu_c
read(*,*,iostat=iostat_val) k_c
read(*,*,iostat=iostat_val) mu_wc
read(*,*,iostat=iostat_val) W
read(*,*,iostat=iostat_val) Lp
read(*,*,iostat=iostat_val) t_plate
read(*,*,iostat=iostat_val) b_depth
read(*,*,iostat=iostat_val) pc_pitch
read(*,*,iostat=iostat_val) beta_deg
read(*,*,iostat=iostat_val) d_port_mm
read(*,*,iostat=iostat_val) dP_max_kpa
read(*,*,iostat=iostat_val) Rfh
read(*,*,iostat=iostat_val) Rfc
read(*,*,iostat=iostat_val) kw
! Defaults handling
if (mu_wh <= 0.0d0) mu_wh = mu_h
if (mu_wc <= 0.0d0) mu_wc = mu_c
! Basic sanity checks
if (W <= 0.0d0 .or. Lp <= 0.0d0 .or. b_depth <= 0.0d0 .or. pc_pitch <= 0.0d0 .or. beta_deg <= 0.0d0) then
write(*,*) 'ERROR: Geometrical plate values must be positive and non-zero.'
stop
end if
! Compute heat duty (average of hot and cold sides)
Q_h = m_h * Cp_h * abs(T_hi - T_ho)
Q_c = m_c * Cp_c * abs(T_co - T_ci)
Q = (Q_h + Q_c) / 2.0d0
! Counter-flow LMTD
dT1 = T_hi - T_co
dT2 = T_ho - T_ci
has_cross_error = .false.
if (dT1 <= 0.0d0 .or. dT2 <= 0.0d0) then
LMTD_cf = 0.0d0
has_cross_error = .true.
else if (abs(dT1 - dT2) < 1.0d-6) then
LMTD_cf = dT1
else
LMTD_cf = (dT1 - dT2) / log(dT1 / dT2)
end if
if (has_cross_error) then
write(*,*) 'ERROR: Temperature cross error detected. LMTD is undefined.'
stop
end if
! Convert dP_max from kPa to Pa
dP_max = dP_max_kpa * 1000.0d0
! Initial geometrical constants
X = (pi * b_depth) / pc_pitch
phi = (1.0d0 / 6.0d0) * (1.0d0 + sqrt(1.0d0 + X**2) + 4.0d0 * sqrt(1.0d0 + X**2 / 2.0d0))
Dh = 2.0d0 * b_depth / phi
if (mode == 1) then
! Verification mode: evaluate at input plate count
call calculate_phe(Nplates, hs_h, hs_c, dP_h, dP_c, dP_port_h, dP_port_c, U_f, U_c, A_prov, A_req, Re_h, Re_c, vch_h, vch_c)
else
! Sizing mode: iterate plate count to meet area and pressure limits
Nplates_curr = Nplates
if (mod(Nplates_curr, 2) == 0) Nplates_curr = Nplates_curr + 1
if (Nplates_curr < 5) Nplates_curr = 5
do while (Nplates_curr <= 50000)
call calculate_phe(Nplates_curr, hs_h, hs_c, dP_h, dP_c, dP_port_h, dP_port_c, U_f, U_c, A_prov, A_req, Re_h, Re_c, vch_h, vch_c)
if (A_prov >= A_req .and. dP_h <= dP_max .and. dP_c <= dP_max) then
Nplates = Nplates_curr
exit
end if
Nplates_curr = Nplates_curr + 2
end do
if (Nplates_curr > 50000) then
write(*,*) 'ERROR: Plate sizing loop did not converge within 50,000 plates. Pressure drop limit may be too restrictive.'
stop
end if
end if
! Print outputs to stdout
write(*,'(A,I5)') 'N_PLATES=', Nplates
write(*,'(A,F18.4)') 'ENLARGEMENT_PHI=', phi
write(*,'(A,F18.6)') 'DH=', Dh
write(*,'(A,F18.4)') 'Q=', Q
write(*,'(A,F18.4)') 'LMTD=', LMTD_cf
write(*,'(A,F18.4)') 'A_PROV=', A_prov
write(*,'(A,F18.4)') 'A_REQ=', A_req
write(*,'(A,F18.4)') 'U_CLEAN=', U_c
write(*,'(A,F18.4)') 'U_FOULED=', U_f
write(*,'(A,F18.4)') 'H_HOT=', hs_h
write(*,'(A,F18.4)') 'H_COLD=', hs_c
write(*,'(A,F18.4)') 'RE_HOT=', Re_h
write(*,'(A,F18.4)') 'RE_COLD=', Re_c
write(*,'(A,F18.4)') 'V_HOT=', vch_h
write(*,'(A,F18.4)') 'V_COLD=', vch_c
write(*,'(A,F18.4)') 'DP_HOT_TOTAL=', dP_h
write(*,'(A,F18.4)') 'DP_COLD_TOTAL=', dP_c
write(*,'(A,F18.4)') 'DP_HOT_PORT=', dP_port_h
write(*,'(A,F18.4)') 'DP_COLD_PORT=', dP_port_c
write(*,'(A,F18.4)') 'DP_HOT_CHAN=', dP_h - dP_port_h
write(*,'(A,F18.4)') 'DP_COLD_CHAN=', dP_c - dP_port_c
contains
subroutine calculate_phe(Np_val, hs_h, hs_c, dP_h, dP_c, dP_port_h, dP_port_c, U_f, U_c, A_prov_val, A_req_val, Re_h_val, Re_c_val, vch_h_val, vch_c_val)
integer, intent(in) :: Np_val
double precision, intent(out) :: hs_h, hs_c, dP_h, dP_c, dP_port_h, dP_port_c, U_f, U_c, A_prov_val, A_req_val
double precision, intent(out) :: Re_h_val, Re_c_val, vch_h_val, vch_c_val
integer :: Nch
double precision :: Ach, Pr_h, Pr_c, Kf_h, p_h, f_turb_h, f_lam_h, f_h
double precision :: Kf_c, p_c, f_turb_c, f_lam_c, f_c
double precision :: beta_rad, Nu_h, Nu_c, Rw
double precision :: U_clean_inv, U_fouled_inv, d_port, v_port_h, v_port_c
double precision :: dP_channel_h, dP_channel_c
Nch = (Np_val - 1) / 2
if (Nch < 1) Nch = 1
Ach = b_depth * W
vch_h_val = m_h / (rho_h * Ach * dble(Nch))
vch_c_val = m_c / (rho_c * Ach * dble(Nch))
Re_h_val = rho_h * vch_h_val * Dh / mu_h
Re_c_val = rho_c * vch_c_val * Dh / mu_c
Pr_h = mu_h * Cp_h / k_h
Pr_c = mu_c * Cp_c / k_c
! Chevron friction factors smooth continuous fit
Kf_h = 0.15d0 + 0.08d0 * exp(0.06d0 * beta_deg)
p_h = 0.30d0 - 0.002d0 * beta_deg
f_turb_h = Kf_h * Re_h_val**(-p_h)
f_lam_h = (10.0d0 + 0.7d0 * beta_deg) / Re_h_val
f_h = f_turb_h
if (f_lam_h > f_h) f_h = f_lam_h
Kf_c = 0.15d0 + 0.08d0 * exp(0.06d0 * beta_deg)
p_c = 0.30d0 - 0.002d0 * beta_deg
f_turb_c = Kf_c * Re_c_val**(-p_c)
f_lam_c = (10.0d0 + 0.7d0 * beta_deg) / Re_c_val
f_c = f_turb_c
if (f_lam_c > f_c) f_c = f_lam_c
! Nusselt numbers
beta_rad = beta_deg * pi / 180.0d0
Nu_h = 0.205d0 * Pr_h**(1.0d0/3.0d0) * (mu_h / mu_wh)**0.14d0 * (f_h * Re_h_val**2 * sin(2.0d0 * beta_rad))**(1.0d0/3.0d0)
Nu_c = 0.205d0 * Pr_c**(1.0d0/3.0d0) * (mu_c / mu_wc)**0.14d0 * (f_c * Re_c_val**2 * sin(2.0d0 * beta_rad))**(1.0d0/3.0d0)
hs_h = Nu_h * k_h / Dh
hs_c = Nu_c * k_c / Dh
Rw = t_plate / kw
U_clean_inv = 1.0d0 / hs_h + Rw + 1.0d0 / hs_c
U_c = 1.0d0 / U_clean_inv
U_fouled_inv = 1.0d0 / hs_h + Rfh + Rw + Rfc + 1.0d0 / hs_c
U_f = 1.0d0 / U_fouled_inv
A_req_val = Q / (U_f * LMTD_cf)
A_prov_val = dble(Np_val - 2) * W * Lp * phi
! Port velocity
d_port = d_port_mm / 1000.0d0
v_port_h = m_h / (rho_h * pi * d_port**2 / 4.0d0)
v_port_c = m_c / (rho_c * pi * d_port**2 / 4.0d0)
! Pressure drops
dP_channel_h = f_h * (Lp / Dh) * (rho_h * vch_h_val**2 / 2.0d0) * (mu_h / mu_wh)**(-0.14d0)
dP_channel_c = f_c * (Lp / Dh) * (rho_c * vch_c_val**2 / 2.0d0) * (mu_c / mu_wc)**(-0.14d0)
dP_port_h = 1.4d0 * 1.0d0 * (rho_h * v_port_h**2 / 2.0d0)
dP_port_c = 1.4d0 * 1.0d0 * (rho_c * v_port_c**2 / 2.0d0)
dP_h = dP_channel_h + dP_port_h
dP_c = dP_channel_c + dP_port_c
end subroutine calculate_phe
end program hx_plate
Solver Description
Chevron Plate Heat Exchanger (PHE) rating and sizing. Evaluates enlargement factor $\phi$, hydraulic diameter $D_h$, convection coefficients using Shah & Martin correlations, and port/channel pressure drops.
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 hx_plate.f90 -o hx_plate
Execution Command:
Execute the program by feeding the sample input file into the program using stdin redirection:
hx_plate < input.txt
📥 Downloads & Local Files
Preview of the required input file (input.txt):
! Mode (1=Verify / Rating, 2=Sizing)
1
! Initial plate count Nplates
29
! Hot fluid inlet temp [°C]
90.0
! Hot fluid outlet temp [°C]
40.0
! Cold fluid inlet temp [°C]
20.0
! Cold fluid outlet temp [°C]
35.0
! Hot fluid mass flow rate [kg/s]
2.0
! Cold fluid mass flow rate [kg/s]
3.0
! Hot fluid density [kg/m3]
990.0
! Hot fluid specific heat [J/kg-K]
4180.0
! Hot fluid viscosity [Pa-s]
0.0008
! Hot fluid conductivity [W/m-K]
0.6
! Hot wall viscosity [Pa-s]
0.0008
! Cold fluid density [kg/m3]
998.0
! Cold fluid specific heat [J/kg-K]
4180.0
! Cold fluid viscosity [Pa-s]
0.001
! Cold fluid conductivity [W/m-K]
0.6
! Cold wall viscosity [Pa-s]
0.001
! Plate width W [m]
0.35
! Thermal height Lp [m]
0.85
! Plate thickness t_plate [m]
0.0005
! Pressing depth b [m]
0.003
! Corrugation wavelength pc [m]
0.008
! Chevron angle beta [deg]
45.0
! Port diameter [mm]
50.0
! Max allowable pressure drop [kPa]
150.0
! Hot side fouling factor [m2-K/W]
0.0001
! Cold side fouling factor [m2-K/W]
0.0001
! Plate wall thermal conductivity [W/m-K]
50.0
1
! Initial plate count Nplates
29
! Hot fluid inlet temp [°C]
90.0
! Hot fluid outlet temp [°C]
40.0
! Cold fluid inlet temp [°C]
20.0
! Cold fluid outlet temp [°C]
35.0
! Hot fluid mass flow rate [kg/s]
2.0
! Cold fluid mass flow rate [kg/s]
3.0
! Hot fluid density [kg/m3]
990.0
! Hot fluid specific heat [J/kg-K]
4180.0
! Hot fluid viscosity [Pa-s]
0.0008
! Hot fluid conductivity [W/m-K]
0.6
! Hot wall viscosity [Pa-s]
0.0008
! Cold fluid density [kg/m3]
998.0
! Cold fluid specific heat [J/kg-K]
4180.0
! Cold fluid viscosity [Pa-s]
0.001
! Cold fluid conductivity [W/m-K]
0.6
! Cold wall viscosity [Pa-s]
0.001
! Plate width W [m]
0.35
! Thermal height Lp [m]
0.85
! Plate thickness t_plate [m]
0.0005
! Pressing depth b [m]
0.003
! Corrugation wavelength pc [m]
0.008
! Chevron angle beta [deg]
45.0
! Port diameter [mm]
50.0
! Max allowable pressure drop [kPa]
150.0
! Hot side fouling factor [m2-K/W]
0.0001
! Cold side fouling factor [m2-K/W]
0.0001
! Plate wall thermal conductivity [W/m-K]
50.0