program flat_plate_convection implicit none ! Variables double precision :: L, T_s, T_inf, U_inf, rho, mu, nu, k_f, Pr, Cp double precision :: Re_L, Re_x, Re_crit, x_crit double precision :: Nu_lam, Nu_turb, Nu_mixed, h_avg double precision :: Q_total, A, W double precision :: delta_lam, delta_turb, Cf_lam, Cf_turb double precision :: delta_t_lam character(len=20) :: flow_regime integer :: i, n_points double precision :: x, dx, local_Re, local_Nu, local_h double precision :: local_Cf, local_delta ! Read inputs read(*,*) L ! Plate length [m] read(*,*) W ! Plate width [m] read(*,*) T_s ! Surface temperature [°C] read(*,*) T_inf ! Free stream temperature [°C] read(*,*) U_inf ! Free stream velocity [m/s] read(*,*) rho ! Fluid density [kg/m³] read(*,*) mu ! Dynamic viscosity [Pa·s] read(*,*) k_f ! Thermal conductivity [W/m·K] read(*,*) Pr ! Prandtl number read(*,*) Cp ! Specific heat [J/kg·K] ! Calculate properties nu = mu / rho A = L * W Re_crit = 5.0d5 ! Reynolds number at trailing edge Re_L = U_inf * L / nu ! Determine flow regime if (Re_L < Re_crit) then flow_regime = 'Laminar' else flow_regime = 'Mixed (Lam+Turb)' end if ! Critical length (transition point) x_crit = Re_crit * nu / U_inf if (x_crit > L) x_crit = L ! ============================================ ! AVERAGE NUSSELT NUMBER CALCULATIONS ! ============================================ ! Laminar correlation (entire plate laminar): Nu = 0.664 Re^0.5 Pr^(1/3) Nu_lam = 0.664d0 * Re_L**0.5d0 * Pr**(1.0d0/3.0d0) ! Turbulent correlation (entire plate turbulent): Nu = 0.037 Re^0.8 Pr^(1/3) Nu_turb = 0.037d0 * Re_L**0.8d0 * Pr**(1.0d0/3.0d0) ! Mixed correlation: Nu = (0.037 Re^0.8 - 871) Pr^(1/3) if (Re_L > Re_crit) then Nu_mixed = (0.037d0 * Re_L**0.8d0 - 871.0d0) * Pr**(1.0d0/3.0d0) else Nu_mixed = Nu_lam end if ! Average heat transfer coefficient h_avg = Nu_mixed * k_f / L ! Total heat transfer Q_total = h_avg * A * abs(T_s - T_inf) ! Boundary layer thicknesses at trailing edge if (Re_L < Re_crit) then delta_lam = 5.0d0 * L / sqrt(Re_L) delta_turb = 0.0d0 Cf_lam = 1.328d0 / sqrt(Re_L) Cf_turb = 0.0d0 else delta_lam = 5.0d0 * x_crit / sqrt(Re_crit) delta_turb = 0.37d0 * L / Re_L**0.2d0 Cf_lam = 1.328d0 / sqrt(Re_crit) Cf_turb = 0.074d0 / Re_L**0.2d0 end if ! Thermal boundary layer thickness (laminar) delta_t_lam = delta_lam / Pr**(1.0d0/3.0d0) ! ============================================ ! OUTPUT RESULTS ! ============================================ write(*,'(A)') '============================================================' write(*,'(A)') ' FORCED CONVECTION - FLOW OVER A FLAT PLATE' write(*,'(A)') '============================================================' write(*,*) write(*,'(A)') '--- INPUT PARAMETERS ----------------------------------------' write(*,'(A,F12.4,A)') ' Plate Length (L) = ', L, ' m' write(*,'(A,F12.4,A)') ' Plate Width (W) = ', W, ' m' write(*,'(A,F12.4,A)') ' Surface Area (A) = ', A, ' m²' write(*,'(A,F12.2,A)') ' Surface Temperature = ', T_s, ' °C' write(*,'(A,F12.2,A)') ' Free Stream Temperature = ', T_inf, ' °C' write(*,'(A,F12.4,A)') ' Free Stream Velocity = ', U_inf, ' m/s' write(*,*) write(*,'(A)') '--- FLUID PROPERTIES ----------------------------------------' write(*,'(A,ES12.4,A)') ' Density (ρ) = ', rho, ' kg/m³' write(*,'(A,ES12.4,A)') ' Dynamic Viscosity (μ) = ', mu, ' Pa·s' write(*,'(A,ES12.4,A)') ' Kinematic Viscosity (ν) = ', nu, ' m²/s' write(*,'(A,F12.4,A)') ' Thermal Conductivity = ', k_f, ' W/m·K' write(*,'(A,F12.4)') ' Prandtl Number (Pr) = ', Pr write(*,*) write(*,'(A)') '--- FLOW ANALYSIS -------------------------------------------' write(*,'(A,ES12.4)') ' Reynolds Number (Re_L) = ', Re_L write(*,'(A,ES12.4)') ' Critical Re = ', Re_crit write(*,'(A,A)') ' Flow Regime = ', trim(flow_regime) write(*,'(A,F12.4,A)') ' Transition Point (x_cr) = ', x_crit, ' m' write(*,*) write(*,'(A)') '--- HEAT TRANSFER RESULTS -----------------------------------' write(*,'(A,F12.4)') ' Average Nusselt Number = ', Nu_mixed write(*,'(A,F12.4,A)') ' Average h = ', h_avg, ' W/m²·K' write(*,'(A,F12.4,A)') ' Total Heat Transfer (Q) = ', Q_total, ' W' write(*,*) write(*,'(A)') '--- BOUNDARY LAYER ------------------------------------------' if (Re_L < Re_crit) then write(*,'(A,F12.6,A)') ' δ at trailing edge = ', delta_lam*1000, ' mm' write(*,'(A,F12.6,A)') ' δ_thermal at trailing = ', delta_t_lam*1000, ' mm' write(*,'(A,ES12.4)') ' Average Cf (laminar) = ', Cf_lam else write(*,'(A,F12.6,A)') ' δ_lam at transition = ', delta_lam*1000, ' mm' write(*,'(A,F12.6,A)') ' δ_turb at trailing edge = ', delta_turb*1000, ' mm' write(*,'(A,ES12.4)') ' Cf_lam (at transition) = ', Cf_lam write(*,'(A,ES12.4)') ' Cf_turb (average) = ', Cf_turb end if write(*,*) ! ============================================ ! LOCAL PROFILE DATA (for chart) ! ============================================ write(*,'(A)') '--- LOCAL PROFILE ALONG PLATE --------------------------------' write(*,'(A)') ' x [m] Re_x Nu_x h_x [W/m2.K] d [mm]' write(*,'(A)') ' ------------------------------------------------------------' n_points = 40 dx = L / dble(n_points) do i = 1, n_points x = dble(i) * dx local_Re = U_inf * x / nu if (local_Re < Re_crit) then ! Laminar local_Nu = 0.332d0 * sqrt(local_Re) * Pr**(1.0d0/3.0d0) local_delta = 5.0d0 * x / sqrt(local_Re) local_Cf = 0.664d0 / sqrt(local_Re) else ! Turbulent local_Nu = 0.0296d0 * local_Re**0.8d0 * Pr**(1.0d0/3.0d0) local_delta = 0.37d0 * x / local_Re**0.2d0 local_Cf = 0.0592d0 / local_Re**0.2d0 end if local_h = local_Nu * k_f / x write(*,'(F8.4,2X,ES12.4,2X,F10.2,4X,F10.4,4X,F10.4)') & x, local_Re, local_Nu, local_h, local_delta*1000 end do write(*,*) write(*,'(A)') '--- CORRELATIONS USED ---------------------------------------' write(*,'(A)') ' Laminar: Nu_x = 0.332 Re_x^0.5 Pr^(1/3)' write(*,'(A)') ' Turbulent: Nu_x = 0.0296 Re_x^0.8 Pr^(1/3)' write(*,'(A)') ' Average: Nu_L = (0.037 Re_L^0.8 - 871) Pr^(1/3)' write(*,'(A)') ' Valid for: 0.6 < Pr < 60, Re_crit = 5×10⁵' end program flat_plate_convection