๐ป 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.
Aerodynamic Lift & Drag Solver
Core Numerical Engine in Fortran 90 โข 25 total downloads
lift_drag.f90
! =========================================================================
! Source File: lift_drag.f90
! =========================================================================
!==============================================================================
! ThermoFluidCalc โ Calculator #24 : Lift CL & Induced Drag CDi
!==============================================================================
! Physics :
! CL = CL_alpha * (alpha - alpha_0) (linear lift model)
! CDi = CL^2 / (pi * e * AR) (induced drag)
! CD = CD0 + CDi (total drag, parabolic polar)
! AR = b^2 / S (aspect ratio)
! L = CL * q * S , D = CD * q * S
! L/D = CL / CD
!
! Modes:
! 1 = Single point (one angle of attack)
! 2 = Alpha sweep (range of angles)
!
! Reference : Gupta, ยง7.8, Eq. 7.87
!
! Build:
! gfortran -O2 -o lift_drag lift_drag.f90
!==============================================================================
program lift_drag
implicit none
integer, parameter :: dp = selected_real_kind(15, 307)
real(dp), parameter :: PI = 3.141592653589793238_dp
integer, parameter :: MAX_PTS = 5000
integer :: mode, npts, i
real(dp) :: rho, V, S, b, e, CD0, CL_alpha, alpha0_deg, alpha_deg
real(dp) :: alpha_min, alpha_max, dalpha
real(dp) :: AR, q, CL, CDi, CD, L_force, D_force, LD
real(dp) :: alpha_rad, alpha0_rad
real(dp) :: best_LD, best_alpha, best_CL
read(*,*) mode
select case (mode)
!==========================================================================
! MODE 1 : Single point
!==========================================================================
case (1)
backspace(5)
read(*,*) mode, rho, V, S, b, e, CD0, CL_alpha, alpha_deg, alpha0_deg
call validate_inputs(rho, V, S, b, e)
AR = b**2 / S
q = 0.5_dp * rho * V**2
alpha_rad = alpha_deg * PI / 180.0_dp
alpha0_rad= alpha0_deg * PI / 180.0_dp
CL = CL_alpha * (alpha_rad - alpha0_rad)
CDi = CL**2 / (PI * e * AR)
CD = CD0 + CDi
L_force = CL * q * S
D_force = CD * q * S
if (CD > 0.0_dp) then
LD = CL / CD
else
LD = 0.0_dp
end if
write(*,'(A,I1)') 'MODE=', mode
write(*,'(A)') 'MODE_NAME=Single Point'
write(*,'(A,F12.6)') 'AR=', AR
write(*,'(A,ES15.8)') 'Q=', q
write(*,'(A,F12.6)') 'ALPHA_DEG=', alpha_deg
write(*,'(A,F12.6)') 'ALPHA0_DEG=', alpha0_deg
write(*,'(A,F12.6)') 'CL=', CL
write(*,'(A,F12.6)') 'CDI=', CDi
write(*,'(A,F12.6)') 'CD0=', CD0
write(*,'(A,F12.6)') 'CD=', CD
write(*,'(A,ES15.8)') 'LIFT=', L_force
write(*,'(A,ES15.8)') 'DRAG=', D_force
write(*,'(A,F12.4)') 'LD=', LD
!==========================================================================
! MODE 2 : Alpha sweep
!==========================================================================
case (2)
backspace(5)
read(*,*) mode, rho, V, S, b, e, CD0, CL_alpha, alpha0_deg, &
alpha_min, alpha_max, npts
call validate_inputs(rho, V, S, b, e)
if (npts < 2) npts = 2
if (npts > MAX_PTS) npts = MAX_PTS
AR = b**2 / S
q = 0.5_dp * rho * V**2
alpha0_rad = alpha0_deg * PI / 180.0_dp
dalpha = (alpha_max - alpha_min) / real(npts - 1, dp)
best_LD = -1.0e30_dp
best_alpha = alpha_min
best_CL = 0.0_dp
write(*,'(A,I1)') 'MODE=', mode
write(*,'(A)') 'MODE_NAME=Alpha Sweep'
write(*,'(A,F12.6)') 'AR=', AR
write(*,'(A,ES15.8)') 'Q=', q
write(*,'(A,F12.6)') 'CD0=', CD0
write(*,'(A,I5)') 'NPTS=', npts
write(*,'(A)') 'DATA_START'
do i = 0, npts - 1
alpha_deg = alpha_min + real(i, dp) * dalpha
alpha_rad = alpha_deg * PI / 180.0_dp
CL = CL_alpha * (alpha_rad - alpha0_rad)
CDi = CL**2 / (PI * e * AR)
CD = CD0 + CDi
if (CD > 0.0_dp) then
LD = CL / CD
else
LD = 0.0_dp
end if
if (LD > best_LD) then
best_LD = LD
best_alpha = alpha_deg
best_CL = CL
end if
write(*,'(F8.3,4(A,F12.6))') alpha_deg, ',', CL, ',', CDi, ',', CD, ',', LD
end do
write(*,'(A)') 'DATA_END'
write(*,'(A,F12.4)') 'MAX_LD=', best_LD
write(*,'(A,F8.3)') 'ALPHA_MAXLD=', best_alpha
write(*,'(A,F12.6)') 'CL_MAXLD=', best_CL
case default
write(*,'(A)') 'ERROR=Invalid mode (must be 1 or 2).'
stop
end select
contains
subroutine validate_inputs(rho, V, S, b, e)
real(dp), intent(in) :: rho, V, S, b, e
if (rho <= 0.0_dp) then
write(*,'(A)') 'ERROR=Density must be positive.'; stop
end if
if (V <= 0.0_dp) then
write(*,'(A)') 'ERROR=Velocity must be positive.'; stop
end if
if (S <= 0.0_dp) then
write(*,'(A)') 'ERROR=Wing area must be positive.'; stop
end if
if (b <= 0.0_dp) then
write(*,'(A)') 'ERROR=Wingspan must be positive.'; stop
end if
if (e <= 0.0_dp .or. e > 1.0_dp) then
write(*,'(A)') 'ERROR=Oswald efficiency must be in (0, 1].'; stop
end if
end subroutine
end program lift_drag
Solver Description
Compute lift force, drag force, and respective coefficients for 2D airfoils and 3D wings.
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 lift_drag.f90 -o lift_drag
Execution Command:
Execute the program by feeding the sample input file into the program using stdin redirection:
lift_drag < input.txt
๐ฅ Downloads & Local Files
Preview of the required input file (input.txt):
! (kg/mรยณ)\nV (m/s)\nS wing area (mรยฒ)\nb wingspan (m)\ne Oswald eff.\nCD0\nCL (per rad)\n0 (รยฐ)\n0 (รยฐ)
1.225
! Parameter 2
60.0
! Parameter 3
16.2
! Parameter 4
11.0
! Parameter 5
0.8
! Parameter 6
0.034
! Parameter 7
5.7
! Parameter 8
5.0
! Parameter 9
-2.0
1.225
! Parameter 2
60.0
! Parameter 3
16.2
! Parameter 4
11.0
! Parameter 5
0.8
! Parameter 6
0.034
! Parameter 7
5.7
! Parameter 8
5.0
! Parameter 9
-2.0