💻 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.

Absorption Refrigeration

Core Numerical Engine in Fortran 90 • 42 total downloads

absorption_refrig.f90
! =========================================================================
! Source File: absorption_refrig.f90
! =========================================================================

program absorption_refrig
    implicit none
    integer :: system_type, i, iostat_val
    double precision :: T_gen, T_abs, T_cond, T_evap
    double precision :: Q_gen, Q_abs, Q_cond, Q_evap, W_pump
    double precision :: COP, COP_carnot, eta_exergy
    double precision :: x_strong, x_weak, f_circ, m_ref, m_sol
    double precision :: h_ref_vap, h_ref_liq, h_sol_strong, h_sol_weak
    double precision :: P_high, P_low, dT_overlap
    double precision :: T_i, COP_i, f_i
    character(len=40) :: system_name

    read(*,*,iostat=iostat_val) system_type
    if (iostat_val /= 0) then
        write(*,*) 'ERROR: Invalid system type input.'
        stop
    end if
    read(*,*,iostat=iostat_val) T_gen
    read(*,*,iostat=iostat_val) T_abs
    read(*,*,iostat=iostat_val) T_cond
    read(*,*,iostat=iostat_val) T_evap
    read(*,*,iostat=iostat_val) m_ref
    if (iostat_val /= 0) then
        write(*,*) 'ERROR: Failed to read all absorption cycle inputs.'
        stop
    end if
    if (T_gen<=0.0d0.or.T_abs<=0.0d0.or.T_cond<=0.0d0.or.T_evap<=0.0d0) then
        write(*,*) 'ERROR: All temperatures must be positive (K).'
        stop
    end if
    if (T_gen <= T_cond) then
        write(*,*) 'ERROR: Generator temp must exceed condenser temp.'
        stop
    end if
    if (T_evap >= T_abs) then
        write(*,*) 'ERROR: Evaporator temp must be below absorber temp.'
        stop
    end if
    if (m_ref <= 0.0d0) m_ref = 1.0d0

    select case(system_type)
    case(1)
        system_name = 'Aqua-Ammonia (NH3-H2O)'
        call aqua_ammonia_cycle(T_gen, T_abs, T_cond, T_evap, m_ref, &
            x_strong, x_weak, f_circ, Q_gen, Q_abs, Q_cond, Q_evap, &
            W_pump, COP, P_high, P_low, &
            h_ref_vap, h_ref_liq, h_sol_strong, h_sol_weak)
    case(2)
        system_name = 'LiBr-H2O (water is refrigerant)'
        call libr_h2o_cycle(T_gen, T_abs, T_cond, T_evap, m_ref, &
            x_strong, x_weak, f_circ, Q_gen, Q_abs, Q_cond, Q_evap, &
            W_pump, COP, P_high, P_low, &
            h_ref_vap, h_ref_liq, h_sol_strong, h_sol_weak)
    case default
        write(*,*) 'ERROR: System type must be 1 (NH3-H2O) or 2 (LiBr-H2O).'
        stop
    end select

    ! Carnot COP for absorption
    ! COP_Carnot = (T_evap/(T_cond-T_evap)) * ((T_gen-T_abs)/T_gen)
    COP_carnot = (T_evap/(T_cond - T_evap)) * ((T_gen - T_abs)/T_gen)
    if (COP_carnot < 0.0d0) COP_carnot = 0.0d0
    eta_exergy = COP / max(COP_carnot, 1.0d-10)

    dT_overlap = T_abs - T_evap   ! temperature lift

    write(*,'(A)') '============================================================'
    write(*,'(A)') '   ABSORPTION REFRIGERATION CYCLE ENGINE'
    write(*,'(A)') '============================================================'
    write(*,*)
    write(*,'(A)') '--- INPUTS --------------------------------------------------'
    write(*,'(A,A)')        '  System Type               = ', trim(system_name)
    write(*,'(A,F12.2,A)')  '  Generator Temperature     = ', T_gen, ' K'
    write(*,'(A,F12.2,A)')  '  Absorber Temperature      = ', T_abs, ' K'
    write(*,'(A,F12.2,A)')  '  Condenser Temperature     = ', T_cond, ' K'
    write(*,'(A,F12.2,A)')  '  Evaporator Temperature    = ', T_evap, ' K'
    write(*,'(A,ES12.4,A)') '  Refrigerant Flow          = ', m_ref, ' kg/s'
    write(*,*)
    write(*,'(A)') '--- SOLUTION PROPERTIES -------------------------------------'
    write(*,'(A,F10.4)')    '  Strong Solution x_s       = ', x_strong
    write(*,'(A,F10.4)')    '  Weak Solution x_w         = ', x_weak
    write(*,'(A,F10.4)')    '  Circulation Ratio f       = ', f_circ
    write(*,'(A,ES12.4,A)') '  Solution Flow (strong)    = ', m_ref*f_circ, ' kg/s'
    write(*,'(A,ES12.4,A)') '  High Pressure             = ', P_high, ' Pa'
    write(*,'(A,ES12.4,A)') '  Low Pressure              = ', P_low, ' Pa'
    write(*,*)
    write(*,'(A)') '--- HEAT DUTIES ---------------------------------------------'
    write(*,'(A,ES12.4,A)') '  Generator Heat Qgen       = ', Q_gen, ' W'
    write(*,'(A,ES12.4,A)') '  Absorber Heat Qabs        = ', Q_abs, ' W'
    write(*,'(A,ES12.4,A)') '  Condenser Heat Qcond      = ', Q_cond, ' W'
    write(*,'(A,ES12.4,A)') '  Evaporator Heat Qevap     = ', Q_evap, ' W'
    write(*,'(A,ES12.4,A)') '  Pump Work Wpump           = ', W_pump, ' W'
    write(*,*)
    write(*,'(A)') '--- PERFORMANCE ---------------------------------------------'
    write(*,'(A,F10.4)')    '  COP (cooling)             = ', COP
    write(*,'(A,F10.4)')    '  COP Carnot (absorption)   = ', COP_carnot
    write(*,'(A,F10.4)')    '  Exergetic Efficiency      = ', eta_exergy
    write(*,'(A,F10.2,A)')  '  Temperature Lift          = ', dT_overlap, ' K'
    write(*,*)
    write(*,'(A)') '  Energy Balance Check:'
    write(*,'(A,ES12.4,A)') '    Qgen + Qevap + Wpump   = ', Q_gen+Q_evap+W_pump, ' W'
    write(*,'(A,ES12.4,A)') '    Qcond + Qabs            = ', Q_cond+Q_abs, ' W'
    write(*,*)

    ! COP vs Generator Temperature sweep
    write(*,'(A)') '--- COP VS GENERATOR TEMPERATURE SWEEP ----------------------'
    write(*,'(A)') '  T_gen[K]      COP           f_circ        COP_Carnot'
    write(*,'(A)') '  -----------------------------------------------------------'
    do i = 1, 50
        T_i = T_cond + 5.0d0 + (T_gen + 50.0d0 - T_cond - 5.0d0)*dble(i-1)/49.0d0
        select case(system_type)
        case(1)
            call aqua_ammonia_cop(T_i, T_abs, T_cond, T_evap, COP_i, f_i)
        case(2)
            call libr_h2o_cop(T_i, T_abs, T_cond, T_evap, COP_i, f_i)
        end select
        write(*,'(F12.2,2X,F10.4,2X,F10.4,2X,F10.4)') T_i, COP_i, f_i, &
            (T_evap/(T_cond-T_evap))*((T_i-T_abs)/T_i)
    end do
    write(*,*)

    ! COP vs Evaporator Temperature sweep
    write(*,'(A)') '--- COP VS EVAPORATOR TEMPERATURE SWEEP ---------------------'
    write(*,'(A)') '  T_evap[K]     COP           COP_Carnot'
    write(*,'(A)') '  -------------------------------------------'
    do i = 1, 40
        T_i = 240.0d0 + (T_abs - 5.0d0 - 240.0d0)*dble(i-1)/39.0d0
        select case(system_type)
        case(1)
            call aqua_ammonia_cop(T_gen, T_abs, T_cond, T_i, COP_i, f_i)
        case(2)
            call libr_h2o_cop(T_gen, T_abs, T_cond, T_i, COP_i, f_i)
        end select
        write(*,'(F12.2,2X,F10.4,2X,F10.4)') T_i, COP_i, &
            (T_i/max(T_cond-T_i,1.0d0))*((T_gen-T_abs)/T_gen)
    end do
    write(*,*)
    write(*,'(A)') '--- CORRELATIONS USED ---------------------------------------'
    write(*,'(A)') '  COP = Qevap / (Qgen + Wpump).'
    write(*,'(A)') '  COP_Carnot = (Te/(Tc-Te))*(Tg-Ta)/Tg.'
    write(*,'(A)') '  f = x_strong / (x_strong - x_weak).'
    write(*,'(A)') '  Simplified property correlations for educational use.'

contains

    subroutine aqua_ammonia_cycle(Tg, Ta, Tc, Te, mref, xs, xw, f, &
            Qg, Qa, Qc, Qe, Wp, cop, Ph, Pl, hrv, hrl, hss, hsw)
        implicit none
        double precision, intent(in) :: Tg, Ta, Tc, Te, mref
        double precision, intent(out) :: xs, xw, f, Qg, Qa, Qc, Qe, Wp, cop
        double precision, intent(out) :: Ph, Pl, hrv, hrl, hss, hsw
        double precision :: h_evap_nh3, Tg_C, Ta_C, Tc_C, Te_C

        Tg_C = Tg - 273.15d0
        Ta_C = Ta - 273.15d0
        Tc_C = Tc - 273.15d0
        Te_C = Te - 273.15d0

        ! Simplified NH3 saturation pressures
        Ph = 1.0d6 * exp(11.0d0 - 3100.0d0/Tc)   ! condenser
        Pl = 1.0d6 * exp(11.0d0 - 3100.0d0/Te)   ! evaporator
        if (Ph < Pl) Ph = Pl * 2.0d0

        ! Simplified concentration correlations
        xs = 0.60d0 - 0.002d0*(Ta_C - 30.0d0)    ! strong solution (rich in NH3)
        xw = 0.30d0 + 0.002d0*(Tg_C - 100.0d0)   ! weak solution (lean)
        if (xs > 0.70d0) xs = 0.70d0
        if (xs < 0.30d0) xs = 0.30d0
        if (xw > xs - 0.02d0) xw = xs - 0.02d0
        if (xw < 0.05d0) xw = 0.05d0

        f = xs / max(xs - xw, 1.0d-10)

        ! Simplified enthalpies (kJ/kg → J/kg)
        h_evap_nh3 = 1370.0d3   ! latent heat NH3
        hrv = 1450.0d3 + 2100.0d0*(Te_C - (-33.0d0))   ! vapor from evaporator
        hrl = 300.0d3 + 4500.0d0*(Tc_C - 25.0d0)       ! liquid from condenser
        hss = -100.0d3 + 3800.0d0*(Ta_C - 25.0d0)       ! strong solution
        hsw = 100.0d3 + 4200.0d0*(Tg_C - 100.0d0)       ! weak solution

        Qe = mref * h_evap_nh3
        Qc = mref * (hrv - hrl)
        Qg = mref * (hrv - hrl) + mref * f * (hsw - hss)
        if (Qg < Qe) Qg = Qe * 1.5d0
        Qa = Qg + Qe - Qc
        Wp = mref * f * 0.001d0 * (Ph - Pl) / 0.80d0   ! pump work
        cop = Qe / max(Qg + Wp, 1.0d-10)
    end subroutine aqua_ammonia_cycle

    subroutine libr_h2o_cycle(Tg, Ta, Tc, Te, mref, xs, xw, f, &
            Qg, Qa, Qc, Qe, Wp, cop, Ph, Pl, hrv, hrl, hss, hsw)
        implicit none
        double precision, intent(in) :: Tg, Ta, Tc, Te, mref
        double precision, intent(out) :: xs, xw, f, Qg, Qa, Qc, Qe, Wp, cop
        double precision, intent(out) :: Ph, Pl, hrv, hrl, hss, hsw
        double precision :: Tg_C, Ta_C, Tc_C, Te_C

        Tg_C = Tg - 273.15d0
        Ta_C = Ta - 273.15d0
        Tc_C = Tc - 273.15d0
        Te_C = Te - 273.15d0

        ! Water saturation pressures (simplified Antoine)
        Ph = 610.78d0 * exp(17.27d0*Tc_C/(Tc_C+237.3d0))
        Pl = 610.78d0 * exp(17.27d0*Te_C/(Te_C+237.3d0))
        if (Ph < Pl) Ph = Pl * 2.0d0

        ! LiBr concentrations (mass fraction LiBr)
        ! Strong = high LiBr (leaves generator), Weak = low LiBr (leaves absorber)
        xw = 0.55d0 + 0.001d0*(Ta_C - 30.0d0)
        xs = 0.62d0 + 0.001d0*(Tg_C - 80.0d0)
        if (xs > 0.70d0) xs = 0.70d0
        if (xs < xw + 0.02d0) xs = xw + 0.02d0
        if (xw < 0.40d0) xw = 0.40d0

        ! f = xs/(xs-xw) for LiBr system (based on LiBr mass balance)
        f = xs / max(xs - xw, 1.0d-10)

        ! Simplified water enthalpies
        hrv = 2500.0d3 + 1900.0d0*(Te_C - 5.0d0)
        hrl = 105.0d3 + 4186.0d0*(Tc_C - 25.0d0)
        hss = 150.0d3 + 2000.0d0*(Tg_C - 80.0d0)
        hsw = 80.0d3 + 2200.0d0*(Ta_C - 30.0d0)

        Qe = mref * (hrv - hrl)
        Qc = mref * (hrv - hrl)
        Qg = mref * hrv + mref*(f-1.0d0)*hss - mref*f*hsw
        if (Qg < 0.0d0) Qg = Qe * 1.4d0
        Qa = Qg + Qe - Qc
        Wp = mref * f * 0.001d0 * (Ph - Pl) / 0.75d0
        cop = Qe / max(Qg + Wp, 1.0d-10)
    end subroutine libr_h2o_cycle

    subroutine aqua_ammonia_cop(Tg, Ta, Tc, Te, cop_out, f_out)
        implicit none
        double precision, intent(in) :: Tg, Ta, Tc, Te
        double precision, intent(out) :: cop_out, f_out
        double precision :: xs, xw, Qg, Qa, Qc, Qe, Wp, Ph, Pl
        double precision :: hrv, hrl, hss, hsw
        call aqua_ammonia_cycle(Tg, Ta, Tc, Te, 1.0d0, xs, xw, f_out, &
            Qg, Qa, Qc, Qe, Wp, cop_out, Ph, Pl, hrv, hrl, hss, hsw)
    end subroutine aqua_ammonia_cop

    subroutine libr_h2o_cop(Tg, Ta, Tc, Te, cop_out, f_out)
        implicit none
        double precision, intent(in) :: Tg, Ta, Tc, Te
        double precision, intent(out) :: cop_out, f_out
        double precision :: xs, xw, Qg, Qa, Qc, Qe, Wp, Ph, Pl
        double precision :: hrv, hrl, hss, hsw
        call libr_h2o_cycle(Tg, Ta, Tc, Te, 1.0d0, xs, xw, f_out, &
            Qg, Qa, Qc, Qe, Wp, cop_out, Ph, Pl, hrv, hrl, hss, hsw)
    end subroutine libr_h2o_cop

end program absorption_refrig

Solver Description

Model aqua-ammonia (NH3-H2O) and lithium bromide (LiBr-H2O) absorption refrigeration systems. Evaluates cycle thermodynamic state properties, circulation ratio, heat duty at each component (generator, condenser, evaporator, absorber), pump work, cooling Coefficient of Performance (COP), Carnot COP limit, and second-law exergy efficiency.

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

Execution Command:

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

absorption_refrig < input.txt

📥 Downloads & Local Files

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

! Working pair (1=NH3-H2O, 2=LiBr-H2O)
1
! Generator temperature T_gen [K]
403.15
! Absorber temperature T_abs [K]
308.15
! Condenser temperature T_cond [K]
308.15
! Evaporator temperature T_evap [K]
263.15
! Refrigerant mass flow rate [kg/s]
1.0