diff --git a/cases.csv b/cases.csv index b057975e..512ca9b2 100644 --- a/cases.csv +++ b/cases.csv @@ -270,9 +270,12 @@ GSw_PRM_StressOutages,Whether to apply the availability factor (forced + schedul GSw_PRM_StressSeedLoadLevel,Region hierarchy level at which to include peak coincident load days as seeded stress periods (or False to ignore peaks),false; False; FALSE; r; nercr; transreg; transgrp; cendiv; st; interconnect; country; usda_region; ccreg,transgrp, GSw_PRM_StressSeedMinRElevel,Region hierarchy level at which to include minimum wind and solar capacity factor days as seeded stress periods (or False to ignore min-wind/solar CF days),false; False; FALSE; r; nercr; transreg; transgrp; cendiv; st; interconnect; country; usda_region; ccreg,interconnect, GSw_PRM_StressStorageCutoff,"How to select ""shoulder"" stress periods giving storage time to recharge before/after high-unserved-energy periods. Two-part switch separated by _. The first argument is ""EUE"" or ""capacity"" or ""absolute"". If ""EUE"" the second argument specifies a PRAS storage headspace / EUE threshold [fraction]; if ""cap"" it specifies a headspace / storage capacity threshold [fraction]; if ""abs"" it specifies absolute number of periods before/after [integer]. Turned off if set to ""off"".",N/A,EUE_0.1, -GSw_PRM_StressThresholdMetrics,"/-delimited list of metrics for identifying stress periods (supported options are: NEUE, LOLH)",N/A,NEUE, -GSw_PRM_StressThresholdLOLH,LOLH threshold [hours/year] above which to re-solve the latest model year with new stress periods; formulated as HierarchyLevel_LOLH_PeriodAggMethod where HierarchyLevel is a column in hierarchy.csv; LOLH is loss of load event-hours per year [event-h/year]; PeriodAggMethod is 'sum' or 'max' over the hours in each period (only used in period selection) (see README.md for detailed notes),N/A,transgrp_2.4_sum, -GSw_PRM_StressThresholdNEUE,Annual NEUE level [ppm] threshold above which to re-solve the latest model year with new stress periods; formulated as HierarchyLevel_NEUEppm_PeriodAggMethod where HierarchyLevel is a column in hierarchy.csv; NEUEppm is normalized expected unserved energy in parts per million [ppm]; PeriodAggMethod is 'sum' or 'max' over the hours in each period (only used in period selection) (see README.md for detailed notes),N/A,transgrp_1_sum, +GSw_PRM_StressThresholdMetrics,"/-delimited list of resource adequacy metrics. NEUE and LOLH identify stress periods; CVAR and NCVAR are check-only metrics and do not add stress periods or update PRM.",N/A,NEUE, +GSw_PRM_StressThresholdLOLH,LOLH threshold [hours/year] above which to re-solve the latest model year with new stress periods; formulated as HierarchyLevel_LOLH_PeriodAggMethod where HierarchyLevel is a column in hierarchy.csv; LOLH is loss of load event-hours per year [event-h/year]; PeriodAggMethod is 'sum' or 'max' over the hours in each period (only used in period selection) (see README.md for detailed notes),N/A,transgrp_2.4_sum, +GSw_PRM_StressThresholdNEUE,Annual NEUE level [ppm] threshold above which to re-solve the latest model year with new stress periods; formulated as HierarchyLevel_NEUEppm_PeriodAggMethod where HierarchyLevel is a column in hierarchy.csv; NEUEppm is normalized expected unserved energy in parts per million [ppm]; PeriodAggMethod is 'sum' or 'max' over the hours in each period (only used in period selection) (see README.md for detailed notes),N/A,transgrp_1_sum, +GSw_PRM_CVARAlpha,Alpha value used for CVAR/NCVAR calculations,float,0.95, +GSw_PRM_StressThresholdNCVAR,Annual NCVAR level [ppm] used only for NCVAR target checks; formulated as HierarchyLevel_NCVARppm_NCVAR_cvar where HierarchyLevel is a column in hierarchy.csv,N/A,transgrp_10_NCVAR_cvar, +GSw_PRM_StressThresholdCVAR,Annual CVAR level [MWh/year] used only for CVAR target checks; formulated as HierarchyLevel_CVARmwh_CVAR_cvar where HierarchyLevel is a column in hierarchy.csv,N/A,country_50000_CVAR_cvar, GSw_PRM_UpdateFraction,Fraction to add to the PRM if a region fails RA threshold (only used if GSw_PRM_UpdateMethod = 1),float,0.02, GSw_PRM_UpdateMethod,Option to update PRM: (0) no update; (1) static update set by GSw_PRM_UpdateFraction; (2) dynamic update informed by PRAS; (3) dynamic update but only after all new stress periods have been added,0; 1; 2; 3,0, GSw_PRMTRADE_level,hierarchy level within which to allow PRM trading,r; nercr; transreg; transgrp; cendiv; st; interconnect; country; usda_region,country, diff --git a/reeds/resource_adequacy/ra_calcs.py b/reeds/resource_adequacy/ra_calcs.py index 214f26ab..33452980 100644 --- a/reeds/resource_adequacy/ra_calcs.py +++ b/reeds/resource_adequacy/ra_calcs.py @@ -24,6 +24,7 @@ def run_pras( write_surplus=False, write_energy=False, write_shortfall_samples=False, + write_shortfall_samples_totals=False, write_availability_samples=False, **kwargs, ): @@ -68,7 +69,9 @@ def run_pras( f"--write_surplus={int(write_surplus)}", f"--write_energy={int(write_energy)}", f"--write_shortfall_samples={int(write_shortfall_samples)}", + f"--write_shortfall_samples_totals={int(write_shortfall_samples_totals)}", f"--write_availability_samples={int(write_availability_samples)}", + f"--cvar_alpha={float(sw['GSw_PRM_CVARAlpha'])}", f"--iteration={iteration}", f"--samples={sw['pras_samples']}", f"--overwrite={int(overwrite)}", @@ -153,13 +156,25 @@ def main(t, tnext, casedir, iteration=0): 2: True, }[int(sw['pras'])] if pras_this_solve_year or int(sw.GSw_PRM_StressIterateMax): + stress_metrics = [ + m.strip().upper() + for m in sw.GSw_PRM_StressThresholdMetrics.split('/') + if m.strip() + ] + write_shortfall_samples_totals = ( + any(m in {'CVAR', 'NCVAR'} for m in stress_metrics) + or int(sw.GSw_PRM_UpdateMethod) > 1 + ) + result = run_pras( casedir, t, iteration=iteration, write_flow=(True if t == max(solveyears) else False), write_energy=True, write_shortfall_samples=(True if int(sw.GSw_PRM_UpdateMethod) > 1 else False), + write_shortfall_samples_totals=write_shortfall_samples_totals, ) + if result.returncode: raise Exception( f"run_pras.jl returned code {result.returncode}. Check gamslog.txt for error trace." diff --git a/reeds/resource_adequacy/run_pras.jl b/reeds/resource_adequacy/run_pras.jl index ec516737..80ac0387 100644 --- a/reeds/resource_adequacy/run_pras.jl +++ b/reeds/resource_adequacy/run_pras.jl @@ -1,5 +1,6 @@ #%% Imports import ArgParse +import CSV import DataFrames import Logging import LoggingExtras @@ -70,10 +71,20 @@ function parse_commandline() default = 0 required = false "--write_shortfall_samples" - help = "Write the sample-level shortfall" + help = "Write per-sample hourly shortfall by region" arg_type = Int default = 0 required = false + "--write_shortfall_samples_totals" + help = "Write per-sample total shortfall by region over the full PRAS time period" + arg_type = Int + default = 0 + required = false + "--cvar_alpha" + help = "Alpha for CVaR (e.g., 0.95)" + arg_type = Float64 + default = 0.95 + required = false "--write_availability_samples" help = "Write the sample-level generator and storage availability" arg_type = Int @@ -182,7 +193,7 @@ function run_pras(pras_system_path::String, args::Dict) if args["write_energy"] == 1 resultspec["energy"] = PRAS.StorageEnergy() end - if args["write_shortfall_samples"] == 1 + if args["write_shortfall_samples"] == 1 || args["write_shortfall_samples_totals"] == 1 resultspec["short_samples"] = PRAS.ShortfallSamples() end if args["write_availability_samples"] == 1 @@ -208,6 +219,41 @@ function run_pras(pras_system_path::String, args::Dict) @info "$(PRAS.EUE(results["short"]))" @info "$(PRAS.NEUE(results["short"]))" + # CVAR results + if haskey(results, "short_samples") + _,_,_,_,energyunit = PRAS.get_params(sys) + alpha = Float64(args["cvar_alpha"]) + cvar_obj = PRAS.CVAR(energyunit, results["short_samples"], alpha) + + cvar_value = PRAS.val(cvar_obj.cvar) + cvar_stderr = PRAS.stderror(cvar_obj.cvar) + cvar_var = cvar_obj.var + + ### Normalize CVaR by total load over the full PRAS time period and report NCVAR in ppm. + total_load = sum(sys.regions.load) + ncvar_value = cvar_value / total_load * 1e6 + ncvar_stderr = cvar_stderr / total_load * 1e6 + ncvar_var = cvar_var / total_load * 1e6 + + @info "CVAR = $(cvar_value)±$(cvar_stderr) MWh; VaR = $(cvar_var) MWh; alpha = $(alpha)" + @info "NCVAR = $(ncvar_value)±$(ncvar_stderr) ppm; VaR = $(ncvar_var) ppm; alpha = $(alpha)" + + ### Write risk metrics to CSV + base = replace(pras_system_path, ".pras"=>"") + riskfile = base * "-risk_metrics.csv" + dfrisk = DF.DataFrame( + metric = ["CVAR", "NCVAR"], + alpha = [alpha, alpha], + region = ["USA", "USA"], + unit = ["MWh", "ppm"], + value = [cvar_value, ncvar_value], + stderr = [cvar_stderr, ncvar_stderr], + var = [cvar_var, ncvar_var], + ) + CSV.write(riskfile, dfrisk) + @info("Wrote risk metrics to $(riskfile)") + end + ## Filter out DC regions used for VSC HVDC transmission regions = [r for r in sys.regions.names if !(occursin("|", r))] @@ -285,6 +331,7 @@ function run_pras(pras_system_path::String, args::Dict) end @info("Wrote PRAS surplus to $(surplusfile)") end + ### Storage energy if args["write_energy"] == 1 dfenergy = DF.DataFrame() @@ -302,31 +349,41 @@ function run_pras(pras_system_path::String, args::Dict) @info("Wrote PRAS storage energy to $(energyfile)") end - ### Sample-level shortfall + ### Per-sample hourly shortfall by region if args["write_shortfall_samples"] == 1 - dictshort = Dict(s => DF.DataFrame() for s = 1:args["samples"]) - for s in range(1, args["samples"]) - dictshort[s] = DF.DataFrame( - transpose(getindex.(results["short_samples"][:, :], s)), - sys.regions.names - ) - # subset to regions (filter out DC regions) - dictshort[s] = dictshort[s][:,findall(regions .∈ Ref(sys.regions.names))] - end - ## Write it + sf = results["short_samples"] + ## Use filtered system regions to avoid DC converter pseudo-regions without load. + region_names = sf.regions.names + idx = [findfirst(==(r), region_names) for r in regions] + shortfile = replace(outfile, ".h5"=>"-shortfall_samples.h5") HDF5.h5open(shortfile, "w") do f - ## Create a group for each sample. Within each group, write an array for each region. - for s in range(1, args["samples"]) + for s in 1:args["samples"] HDF5.create_group(f, "$s") - for column in DF._names(dictshort[s]) - f["$s"]["$column", compress=4] = convert(Array, dictshort[s][!, column]) + for (r, i_r) in zip(regions, idx) + arr = Float64.(sf.shortfall[i_r, :, s]) + f["$s"]["$r", compress=4] = arr end end end @info("Wrote PRAS shortfall by sample to $(shortfile)") end + ### Per-sample total shortfall by region over the full PRAS time period, needed for CVaR + if args["write_shortfall_samples_totals"] == 1 + sf = results["short_samples"] + ## Use filtered system regions to avoid DC converter pseudo-regions without load. + totalsfile = replace(outfile, ".h5"=>"-shortfall_totals_by_sample.h5") + HDF5.h5open(totalsfile, "w") do f + f["sample", compress=4] = collect(1:args["samples"]) + f["USA", compress=4] = Float64.(sf[]) + for r in regions + f["$r", compress=4] = Float64.(sf[r]) + end + end + @info("Wrote PRAS shortfall totals by sample to $(totalsfile)") + end + ### Sample-level generator and storage availability if args["write_availability_samples"] == 1 dictavail = Dict(s => DF.DataFrame() for s = 1:args["samples"]) @@ -457,6 +514,7 @@ if abspath(PROGRAM_FILE) == @__FILE__ # "write_surplus" => 0, # "write_energy" => 0, # "write_shortfall_samples" => 1, + # "write_shortfall_samples_totals" => 1, # "write_availability_samples" => 0, # "overwrite" => 1, # "debug" => 0, @@ -466,6 +524,7 @@ if abspath(PROGRAM_FILE) == @__FILE__ # "pras_existing_unit_size" => 1, # "pras_max_unitsize_prm" => 1, # "pras_seed" => 1, + # "cvar_alpha" => 0.95, # ) # reedscase = args["reedscase"] # solve_year = args["solve_year"] @@ -478,11 +537,11 @@ if abspath(PROGRAM_FILE) == @__FILE__ #%% Include ReEDS2PRAS include(joinpath( - args["reedscase"], "reeds", "resource_adequacy", "reeds2pras", "src", "ReEDS2PRAS.jl" + args["reeds_path"], "reeds", "resource_adequacy", "reeds2pras", "src", "ReEDS2PRAS.jl" )) #%% Run it main(args) #%% -end +end \ No newline at end of file diff --git a/reeds/resource_adequacy/stress_periods.py b/reeds/resource_adequacy/stress_periods.py index 804c2da7..1e12b853 100644 --- a/reeds/resource_adequacy/stress_periods.py +++ b/reeds/resource_adequacy/stress_periods.py @@ -17,6 +17,8 @@ # import importlib # importlib.reload(functions) +CVAR_METRICS = {'CVAR', 'NCVAR'} + #%%### Functions def get_pras_stress_metric(case, t, iteration=0, stress_metric:Literal['EUE','LOLE']='EUE'): @@ -180,6 +182,127 @@ def get_and_write_neue(sw, write=True): eue.to_csv(os.path.join(sw['casedir'],'outputs','eue.csv')) return neue +def get_cvar_alpha(sw): + alpha = float(sw.GSw_PRM_CVARAlpha) + if not (0 <= alpha < 1): + raise ValueError(f"GSw_PRM_CVARAlpha must be in [0, 1). Got {alpha}") + return alpha + + +def get_shortfall_totals_by_sample(case, t, iteration=0): + filepath = os.path.join(case, 'handoff', 'PRAS', f'PRAS_{t}i{iteration}-shortfall_totals_by_sample.h5') + if not os.path.isfile(filepath): + raise FileNotFoundError(f"{filepath} not found. Re-run PRAS with --write_shortfall_samples_totals 1.") + df = reeds.io.read_pras_results(filepath) + df.columns = df.columns.astype(str) + if 'sample' in df.columns: + df = df.set_index('sample') + elif df.index.name != 'sample': + df.index = pd.RangeIndex(1, len(df) + 1, name='sample') + return df.apply(pd.to_numeric, errors='coerce').clip(lower=0) + + +def _sample_cvar(samples, alpha=0.95): + x = pd.Series(samples).dropna().astype(float) + if x.empty: + return np.nan + n_tail = max(1, int(np.ceil((1 - alpha) * len(x)))) + return x.sort_values(ascending=False).iloc[:n_tail].mean() + + +def get_annual_cvar_stress_metric(case, t, stress_metric='NCVAR', iteration=0, alpha=0.95): + stress_metric = stress_metric.upper() + if stress_metric not in CVAR_METRICS: + raise NotImplementedError(f"get_annual_cvar_stress_metric only supports {CVAR_METRICS}. Got {stress_metric}.") + + ### Read per-sample total shortfall by region from PRAS. + shortfall_samples = get_shortfall_totals_by_sample(case=case, t=t, iteration=iteration).drop(columns=['USA'], errors='ignore') + + ### Read hourly PRAS load. This is used to normalize CVAR into NCVAR. + dfload = reeds.io.read_h5py_file(os.path.join(case, 'handoff', 'reeds_data', f'pras_load_{t}.h5')) + + ### Calculate CVAR/NCVAR at each supported hierarchy level. + levels = ['country','interconnect','nercr','transreg','transgrp','st','r'] + _metric = {} + for hierarchy_level in levels: + ### Aggregate sample-level regional shortfall to the requested hierarchy level. + rmap = reeds.io.get_rmap(case=case, hierarchy_level=hierarchy_level) + regions = [c for c in shortfall_samples.columns if c in rmap.index] + if not regions: + continue + + shortfall_agg = shortfall_samples[regions].rename(columns=rmap).groupby(axis=1, level=0).sum() + + ### CVAR is the average of the highest-shortfall samples in the alpha tail. + cvar = shortfall_agg.apply(lambda s: _sample_cvar(s, alpha=alpha), axis=0) + + if stress_metric == 'NCVAR': + ### NCVAR normalizes CVAR by total load over the PRAS time period and reports ppm. + load_regions = [c for c in dfload.columns if c in rmap.index] + load_agg = dfload[load_regions].rename(columns=rmap).groupby(axis=1, level=0).sum().sum() + cvar = cvar / load_agg.reindex(cvar.index) * 1e6 + + _metric[hierarchy_level,'cvar'] = cvar + + return pd.concat(_metric, names=['level','metric','region']).rename(stress_metric) + +def evaluate_cvar_target_check(sw, t, iteration=0, stress_metrics=None): + if stress_metrics is None: + stress_metrics = sw.GSw_PRM_StressThresholdMetrics.split('/') + + ### Only evaluate CVAR/NCVAR here. Standard metrics are handled by _evaluate_stress_threshold_criterion(), which can add stress periods. + stress_metrics = [m.upper() for m in stress_metrics if str(m).strip()] + cvar_metrics = [m for m in stress_metrics if m in CVAR_METRICS] + if not cvar_metrics: + return pd.DataFrame() + + records = [] + for metric in cvar_metrics: + switch_name = f'GSw_PRM_StressThreshold{metric}' + if switch_name not in sw.index: + print(f"Warning: {metric} is listed in GSw_PRM_StressThresholdMetrics, but {switch_name} is not defined. Skipping {metric} target check.") + continue + + for criterion in str(sw[switch_name]).split('/'): + criterion = criterion.strip() + if not criterion: + continue + + ### CVAR/NCVAR criteria are check-only and use: HierarchyLevel_Threshold_Metric_cvar + hierarchy_level, stress_level, stress_metric, metric_name = criterion.split('_') + stress_metric = stress_metric.upper() + metric_name = metric_name.lower() + + if stress_metric not in CVAR_METRICS: + raise ValueError(f"Invalid CVAR criterion: {criterion}. Metric must be one of {CVAR_METRICS}.") + if metric_name != 'cvar': + raise ValueError(f"Invalid CVAR criterion: {criterion}. The fourth field must be 'cvar'.") + + ### Read the annual CVAR/NCVAR metric written earlier in main(). + stress_vals = pd.read_csv( + os.path.join(sw.casedir, 'outputs', f'{stress_metric.lower()}_{t}i{iteration}.csv'), + index_col=['level','metric','region'], + ).squeeze(1) + + ### Compare each region against the check-only threshold. + this_test = stress_vals.xs((hierarchy_level, 'cvar'), level=['level','metric']) + threshold = float(stress_level) + + for region, value in this_test.items(): + records.append({'criterion':criterion, 'level':hierarchy_level, 'metric':stress_metric, 'region':region, 'alpha':get_cvar_alpha(sw), 'threshold':threshold, 'value':value, 'passed':value <= threshold}) + + failed = this_test.loc[this_test > threshold] + if len(failed): + print(f"GSw_PRM_StressThreshold = {criterion} failed for:") + print(failed) + print(f"{stress_metric} is check-only: no stress periods and no PRM increment will be added.") + else: + print(f"GSw_PRM_StressThreshold = {criterion} passed") + + dfcheck = pd.DataFrame(records) + if not dfcheck.empty: + dfcheck.to_csv(os.path.join(sw.casedir, 'outputs', f'cvar_target_check_{t}i{iteration}.csv'), index=False) + return dfcheck def get_annual_stress_metric(case, t, stress_metric, iteration=0): """ @@ -397,6 +520,10 @@ def _evaluate_stress_threshold_criterion( def get_stress_metrics_sorted_periods(sw, t, iteration): + stress_metric_switches = [m.upper() for m in sw.GSw_PRM_StressThresholdMetrics.split('/') if str(m).strip() and m.upper() not in CVAR_METRICS] + if not stress_metric_switches: + return {}, None, None + ### Get storage state of charge (SOC) to use in selection of "shoulder" stress periods dfenergy = reeds.io.read_pras_results( os.path.join(sw['casedir'], 'handoff', 'PRAS', f"PRAS_{t}i{iteration}-energy.h5") @@ -418,11 +545,7 @@ def get_stress_metrics_sorted_periods(sw, t, iteration): ### Check all stress criteria; for regions that fail, add new stress periods stress_criteria = {'stress_sorted_periods': {}, 'failed': {}, 'high_stress_periods': {}, 'shoulder_periods': {}} - - # Validation check: Display any GSw_PRM_StressThreshold{metric} - # that is not specified in GSw_PRM_StressThresholdMetrics - stress_metric_switches = sw.GSw_PRM_StressThresholdMetrics.split('/') - + # stress periods column names for writing outputs stress_metrics_units = {'EUE':'MWh/year', 'NEUE':'ppm', 'LOLH':'event-h/year', 'LOLE':'event-day/year'} stress_metrics_col_names = {m:f'{m}_{stress_metrics_units[m]}' for m in @@ -517,27 +640,15 @@ def prm_increment_pras(sw, t, iteration, combined_periods_write, failed_regions) ) ## shortfall data - # read the net shortfall (positive) and net surplus (negative) results - # by sample from PRAS run (MWh) - filepath = os.path.join(sw['casedir'], 'handoff', 'PRAS', - f'PRAS_{sw["t"]}i{iteration}-shortfall_samples.h5') - net_short = reeds.io.read_pras_results(filepath) - # get number of samples - n_samples = len(net_short) - # collapse dict of dataframes by sample in 1 dataframe (keep index to preserve hours) - net_short = pd.concat( - (df.assign(**{"sample": k}) for k, df in net_short.items()), ignore_index=False) - # convert to long format with shortfall by sample, hour, and r - net_short.index.names=['hour'] - net_short = net_short.reset_index().set_index(['sample','hour']) - net_short = net_short.sort_index(level=['sample', 'hour'], ascending=[True, True]) - net_short = net_short.melt( - ignore_index=False, var_name='r', value_name='net_short_mwh').reset_index() - - # zero-out negative values (net surplus) for determining regional unserved energy totals - net_short['net_short_mwh'] = net_short['net_short_mwh'].clip(lower=0) - # calaculate total regional net shortfall for all hours by sample - net_short_crit = net_short.groupby(['r','sample'], as_index=False)['net_short_mwh'].sum() + net_short_totals = get_shortfall_totals_by_sample(case=sw['casedir'], t=t, iteration=iteration) + n_samples = len(net_short_totals) + net_short_crit = ( + net_short_totals + .drop(columns=['USA'], errors='ignore') + .reset_index() + .melt(id_vars='sample', var_name='r', value_name='net_short_mwh') + ) + net_short_crit['net_short_mwh'] = net_short_crit['net_short_mwh'].clip(lower=0) ## get load data dfload = reeds.io.read_file( @@ -697,7 +808,10 @@ def main(sw, t, iteration=0, logging=True): """ #%% Write consolidated stress metrics so far try: - ## TODO: Check if get_and_write_neue() is still needed + ## TODO: Check if get_and_write_neue() is still needed. + # get_and_write_neue() writes outputs/neue.csv and outputs/eue.csv across + # completed PRAS years. The annual metric files below are still written per + # solve year / iteration. _neue_simple = get_and_write_neue(sw, write=True) ## TODO: check if need to refactor or remove @@ -714,6 +828,24 @@ def main(sw, t, iteration=0, logging=True): os.path.join(sw.casedir, 'outputs', f"{stress_metric.lower()}_{t}i{iteration}.csv") ) + # CVAR / NCVAR are additional check-only metrics. They do not replace or modify the standard EUE / NEUE / LOLH calculations above. + cvar_metrics = [m.upper() for m in sw.GSw_PRM_StressThresholdMetrics.split('/') if str(m).strip() and m.upper() in CVAR_METRICS] + + for stress_metric in cvar_metrics: + print(f"Calculating and writing annual {stress_metric} for iteration {iteration}") + + dfmetric = get_annual_cvar_stress_metric(case=sw.casedir, t=t, stress_metric=stress_metric, iteration=iteration, alpha=get_cvar_alpha(sw),) + + # CVAR is an absolute MWh metric calculated across all resource-adequacy years, so convert it to annual average. + # NCVAR is already normalized by total load and remains in ppm. + if stress_metric == 'CVAR': + dfmetric /= len(sw['resource_adequacy_years']) + + dfmetric.round(2).to_csv(os.path.join(sw.casedir,'outputs', f"{stress_metric.lower()}_{t}i{iteration}.csv",)) + + # CVAR / NCVAR are check-only: no stress periods and no PRM increment. + evaluate_cvar_target_check(sw=sw, t=t, iteration=iteration, stress_metrics=cvar_metrics,) + except Exception as err: if int(sw['pras']) == 2: print(traceback.format_exc()) @@ -767,7 +899,6 @@ def main(sw, t, iteration=0, logging=True): #%% Done return - # if __name__ == '__main__': # #%%### option to run script directly for debugging # casedir = "/path/to/ReEDS/runs/runname"