πŸ” Search & Filter All 178+ Engineering Calculators β†’

πŸ’» 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.

Local Pressure Coefficient (Cp)

Core Numerical Engine in Fortran 90 β€’ 30 total downloads

cp_pressure.f90 
! =========================================================================
! Source File: cp_pressure.f90
! =========================================================================

!==============================================================================
! ThermoFluidCalc β€” Calculator #23 : Pressure Coefficient Cp
!==============================================================================
! Physics : Cp = (p - p_inf) / q_inf   where q_inf = 0.5 * rho * V_inf^2
!
! Modes:
!   1 = Single point evaluation
!   2 = Surface Cp distribution (N points)
!   3 = Compressible correction sweep (Prandtl-Glauert, Karman-Tsien, Laitone)
!
! Reference : Gupta, Β§5.7 and Β§7.8
!
! Build:
!   gfortran -O2 -o cp_pressure cp_pressure.f90
!==============================================================================
program cp_pressure
  implicit none

  integer, parameter :: dp = selected_real_kind(15, 307)
  integer, parameter :: MAX_PTS = 5000

  integer  :: mode, N, npts, i
  real(dp) :: p, p_inf, rho, V_inf, gam, q_inf, Cp, Mach
  real(dp) :: xc(MAX_PTS), p_local(MAX_PTS), Cp_arr(MAX_PTS)
  real(dp) :: Cp_min, Cp_max, xc_cpmin, Cp_crit
  real(dp) :: Cp0, M_max, dM, M_cur
  real(dp) :: beta, Cp_PG, Cp_KT, Cp_Lait
  real(dp) :: term_kt, term_la
  character(len=40) :: interp

  read(*,*) mode

  select case (mode)

  !=========================================================================
  ! MODE 1 : Single point
  !=========================================================================
  case (1)
    backspace(5)
    read(*,*) mode, p, p_inf, rho, V_inf, gam

    if (rho <= 0.0_dp .or. V_inf <= 0.0_dp) then
      write(*,'(A)') 'ERROR=Density and velocity must be positive.'
      stop
    end if

    q_inf = 0.5_dp * rho * V_inf**2
    Cp    = (p - p_inf) / q_inf
    Mach  = V_inf / sqrt(gam * p_inf / rho)

    ! Interpretation
    if (abs(Cp - 1.0_dp) < 0.05_dp) then
      interp = 'Stagnation point'
    else if (abs(Cp) < 0.01_dp) then
      interp = 'Freestream (neutral)'
    else if (Cp < 0.0_dp) then
      interp = 'Favorable pressure gradient (suction)'
    else
      interp = 'Adverse pressure gradient'
    end if

    write(*,'(A,I1)')       'MODE=', mode
    write(*,'(A)')          'MODE_NAME=Single Point'
    write(*,'(A,ES15.8)')   'P=', p
    write(*,'(A,ES15.8)')   'P_INF=', p_inf
    write(*,'(A,ES15.8)')   'RHO=', rho
    write(*,'(A,ES15.8)')   'V_INF=', V_inf
    write(*,'(A,ES15.8)')   'GAMMA=', gam
    write(*,'(A,ES15.8)')   'Q_INF=', q_inf
    write(*,'(A,F12.6)')    'CP=', Cp
    write(*,'(A,F12.6)')    'MACH=', Mach
    write(*,'(A,A)')        'INTERPRETATION=', trim(interp)

  !=========================================================================
  ! MODE 2 : Surface distribution
  !=========================================================================
  case (2)
    backspace(5)
    read(*,*) mode, p_inf, rho, V_inf, gam, N

    if (N < 1 .or. N > MAX_PTS) then
      write(*,'(A)') 'ERROR=Number of points must be 1-5000.'
      stop
    end if
    if (rho <= 0.0_dp .or. V_inf <= 0.0_dp) then
      write(*,'(A)') 'ERROR=Density and velocity must be positive.'
      stop
    end if

    q_inf = 0.5_dp * rho * V_inf**2
    Mach  = V_inf / sqrt(gam * p_inf / rho)

    do i = 1, N
      read(*,*) xc(i), p_local(i)
      Cp_arr(i) = (p_local(i) - p_inf) / q_inf
    end do

    ! Find min and max Cp
    Cp_min    = Cp_arr(1)
    Cp_max    = Cp_arr(1)
    xc_cpmin  = xc(1)
    do i = 2, N
      if (Cp_arr(i) < Cp_min) then
        Cp_min   = Cp_arr(i)
        xc_cpmin = xc(i)
      end if
      if (Cp_arr(i) > Cp_max) Cp_max = Cp_arr(i)
    end do

    ! Critical Cp (isentropic)
    if (Mach > 0.0_dp .and. Mach < 1.0_dp) then
      Cp_crit = (2.0_dp / (gam * Mach**2)) * &
        ( ((2.0_dp + (gam - 1.0_dp)*Mach**2) / (gam + 1.0_dp))**(gam/(gam-1.0_dp)) &
          - 1.0_dp )
    else
      Cp_crit = -999.0_dp
    end if

    write(*,'(A,I1)')       'MODE=', mode
    write(*,'(A)')          'MODE_NAME=Surface Distribution'
    write(*,'(A,ES15.8)')   'Q_INF=', q_inf
    write(*,'(A,F12.6)')    'MACH=', Mach
    write(*,'(A,I5)')       'NPTS=', N
    write(*,'(A,F12.6)')    'CP_MIN=', Cp_min
    write(*,'(A,F12.8)')    'XC_AT_CPMIN=', xc_cpmin
    write(*,'(A,F12.6)')    'CP_MAX=', Cp_max
    write(*,'(A,F12.6)')    'CP_CRIT=', Cp_crit

    write(*,'(A)') 'DATA_START'
    do i = 1, N
      write(*,'(F12.8,A,F12.6)') xc(i), ',', Cp_arr(i)
    end do
    write(*,'(A)') 'DATA_END'

  !=========================================================================
  ! MODE 3 : Compressible correction sweep
  !=========================================================================
  case (3)
    backspace(5)
    read(*,*) mode, Cp0, gam, M_max, npts

    if (npts < 2) npts = 2
    if (npts > MAX_PTS) npts = MAX_PTS
    if (M_max >= 1.0_dp) M_max = 0.99_dp
    if (M_max <= 0.0_dp) M_max = 0.9_dp

    write(*,'(A,I1)')       'MODE=', mode
    write(*,'(A)')          'MODE_NAME=Compressible Correction'
    write(*,'(A,F12.6)')    'CP0=', Cp0
    write(*,'(A,F12.6)')    'GAMMA=', gam
    write(*,'(A,F12.6)')    'M_MAX=', M_max
    write(*,'(A,I5)')       'NPTS=', npts

    dM = M_max / real(npts - 1, dp)

    write(*,'(A)') 'DATA_START'
    do i = 0, npts - 1
      M_cur = real(i, dp) * dM
      if (M_cur >= 1.0_dp) M_cur = 0.999_dp

      beta = sqrt(max(1.0_dp - M_cur**2, 1.0e-12_dp))

      ! Prandtl-Glauert
      Cp_PG = Cp0 / beta

      ! Karman-Tsien
      term_kt = beta + (M_cur**2 * Cp0) / (2.0_dp * beta * (1.0_dp + beta))
      if (abs(term_kt) > 1.0e-15_dp) then
        Cp_KT = Cp0 / term_kt
      else
        Cp_KT = Cp0 / beta
      end if

      ! Laitone
      term_la = beta + (M_cur**2 * (1.0_dp + (gam - 1.0_dp)/2.0_dp * M_cur**2) * Cp0) &
                / (2.0_dp * beta)
      if (abs(term_la) > 1.0e-15_dp) then
        Cp_Lait = Cp0 / term_la
      else
        Cp_Lait = Cp0 / beta
      end if

      write(*,'(F10.6,3(A,F12.6))') M_cur, ',', Cp_PG, ',', Cp_KT, ',', Cp_Lait
    end do
    write(*,'(A)') 'DATA_END'

  case default
    write(*,'(A)') 'ERROR=Invalid mode (must be 1-3).'
    stop
  end select

end program cp_pressure

Solver Description

Calculate local pressure coefficient from velocity/pressure differences in incompressible and compressible flows.

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 cp_pressure.f90 -o cp_pressure

Execution Command:

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

cp_pressure < input.txt

πŸ“₯ Downloads & Local Files

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

! Local pressure p (Pa)\npΓ’Λ†ΕΎ (Pa)\n(kg/m³)\nVΓ’Λ†ΕΎ (m/s)\n
102000
! Parameter 2
101325
! Parameter 3
1.225
! Parameter 4
60.0
! Parameter 5
1.4
πŸ“¦ Download Offline Bundle (ZIP) 15 dl πŸ“„ Download input.txt 7 dl πŸ“₯ Download cp_pressure.f90 8 dl