diff --git a/src/core/chp.jl b/src/core/chp.jl
index c59514a9f..77a979371 100644
--- a/src/core/chp.jl
+++ b/src/core/chp.jl
@@ -22,6 +22,7 @@ conflict_res_min_allowable_fraction_of_max = 0.25
     cooling_thermal_factor::Float64 = NaN  # only needed with cooling load
     unavailability_periods::AbstractVector{Dict} = Dict[] # CHP unavailability periods for scheduled and unscheduled maintenance, list of dictionaries with keys of "['month', 'start_week_of_month', 'start_day_of_week', 'start_hour', 'duration_hours'] all values are one-indexed and start_day_of_week uses 1 for Monday, 7 for Sunday
     unavailability_hourly::AbstractVector{Int64} = Int64[] # Hourly 8760 profile of unavailability (1) and availability (0)
+    production_factor_series::Union{Nothing, Vector{<:Real}} = nothing # Optional user-provided production factor time-series (length of time_steps_per_hour * 8760). If provided, this will override the production factor calculated from unavailability.
 
     # Optional inputs:
     size_class::Union{Int, Nothing} = nothing # CHP size class for using appropriate default inputs, with size_class=0 using an average of all other size class data 
@@ -113,6 +114,7 @@ Base.@kwdef mutable struct CHP <: AbstractCHP
     can_serve_space_heating::Bool = true
     can_serve_process_heat::Bool = true
     is_electric_only::Bool = false
+    production_factor_series::Union{Nothing, Vector{<:Real}} = nothing
 
     macrs_option_years::Int = 5
     macrs_bonus_fraction::Float64 = 1.0
diff --git a/src/core/production_factor.jl b/src/core/production_factor.jl
index 6fc19e337..07749ccbf 100644
--- a/src/core/production_factor.jl
+++ b/src/core/production_factor.jl
@@ -278,9 +278,15 @@ end
 production_factor for CHP accounts for unavailability (`unavailability_periods`) of CHP due to 
 scheduled (mostly off-peak) and "unscheduled" (on-peak) maintenance. 
 Note: this same prod_factor should be applied to electric and thermal production
+
+If the user provides their own production_factor_series, that will be returned instead of generating from unavailability.
 """
 function get_production_factor(chp::AbstractCHP, year::Int=2017, outage_start_time_step::Int=0, outage_end_time_step::Int=0, ts_per_hour::Int=1)
     
+    if !(isnothing(chp.production_factor_series))
+        return chp.production_factor_series
+    end
+    
     prod_factor = [1.0 - chp.unavailability_hourly[i] for i in 1:8760 for _ in 1:ts_per_hour]
 
     # Ignore unavailability in time_step if it intersects with an outage interval(s)
diff --git a/test/test_chp_prodfactor.jl b/test/test_chp_prodfactor.jl
new file mode 100644
index 000000000..80e833811
--- /dev/null
+++ b/test/test_chp_prodfactor.jl
@@ -0,0 +1,59 @@
+using Revise
+using REopt
+using JSON
+using DelimitedFiles
+using PlotlyJS
+using Dates
+using Test
+using JuMP
+using HiGHS
+using DotEnv
+DotEnv.load!()
+
+
+###########   CHP production_factor_series test   #############
+@testset "CHP with production_factor_series" begin
+    # Create a simple scenario with electric-only CHP and a custom production factor
+    input_data = Dict(
+        "Site" => Dict(
+            "latitude" => 39.7,
+            "longitude" => -104.9
+        ),
+        "ElectricLoad" => Dict(
+            "loads_kw" => repeat([100.0], 8760),
+            "year" => 2025
+        ),
+        "ElectricTariff" => Dict(
+            "urdb_label" => "5ed6c1a15457a3367add15ae"
+        ),
+        "CHP" => Dict(
+            "is_electric_only" => true,
+            "fuel_cost_per_mmbtu" => 4.0,
+            "max_kw" => 100.0,
+            "min_kw" => 100.0,
+            "production_factor_series" => vcat(repeat([0.5], 4380), repeat([1.0], 4380))  # 50% for first half year, 100% for second half
+        )
+    )
+    
+    s = Scenario(input_data)
+    p = REoptInputs(s)
+    
+    # Verify the production factor series was properly assigned
+    @test !isnothing(s.chp.production_factor_series)
+    @test length(s.chp.production_factor_series) == 8760
+    @test s.chp.production_factor_series[1] ≈ 0.5 atol=0.001
+    @test s.chp.production_factor_series[8760] ≈ 1.0 atol=0.001
+    
+    # Run the optimization
+    m = Model(optimizer_with_attributes(HiGHS.Optimizer, "output_flag" => false, "log_to_console" => false, "mip_rel_gap" => 0.01))
+    results = run_reopt(m, p)
+    
+    # Check that production matches the pattern (lower in first half, higher in second half)
+    first_half_avg = sum(results["CHP"]["electric_production_series_kw"][1:4380]) / 4380
+    second_half_avg = sum(results["CHP"]["electric_production_series_kw"][4381:8760]) / 4380
+    
+    println("CHP production_factor_series test passed!")
+    println("CHP size: ", results["CHP"]["size_kw"], " kW")
+    println("First half avg production: ", round(first_half_avg, digits=2), " kW")
+    println("Second half avg production: ", round(second_half_avg, digits=2), " kW")
+end
\ No newline at end of file