nrel / pras Goto Github PK

NREL's Probabilistic Resource Adequacy Suite (PRAS)

License: Other

Julia 100.00%

pras's Introduction

Probabilistic Resource Adequacy Suite

The Probabilistic Resource Adequacy Suite (PRAS) is a collection of tools for bulk power system resource adequacy analysis and capacity credit calculation. The most recent documentation report (for version 0.6) is available here.

pras's People

Contributors

Stargazers

Watchers

Forkers

yaacob3 anamileva claytonpbarrows rajatpungaliya tamas-borbath gsaenzt niklauskun ehbartolo boyuan276 sriharisundar

pras's Issues

Attempting calculate EFC for generator in rts.pras

I define two systems:

Base System:

sys = SystemModel("Documents/PRAS-master/test/PRASBase/rts.pras")
sys.generators.capacity[1,:] = 0 * sys.generators.capacity[1,:]   # remove this generator from the base system

Augmented system:

sys2 = SystemModel("Documents/PRAS-master/test/PRASBase/rts.pras")

Increase the load to get some baseline EUE:

sys2.regions.load .+= 100
sys.regions.load .+= 100

Then attempt to calculate the EFC of the sys.generators.capacity[1]

cc_result = assess(sys, sys2, EFC{EUE}(100, "1"), SequentialMonteCarlo())

Then get the warning below:

┌ Warning: Gap between upper and lower bound risk metrics is not statistically significant (p_value=0.15741688789903002), stopping bisection. The gap between capacity bounds is 100 MW, while the target stopping gap was 1 MW.
└ @ PRAS.CapacityCredit ~/.julia/packages/PRAS/jO1jk/src/CapacityCredit/EFC.jl:81
PRAS.CapacityCredit.CapacityCreditResult{EFC{EUE}, EUE{8784, 1, Hour, MWh}, MW}(EUE = 0.4±0.2 MWh/8784h, 0, 100, [0, 100], EUE{8784, 1, Hour, MWh}[EUE = 0.6±0.3 MWh/8784h, EUE = 0.3±0.2 MWh/8784h])

cc_lower, cc_upper = extrema(cc_result)
(0, 100)

Am I doing this right?

Thanks!

Workflow question: Specifying a system using HDF5 file

Hi Gord,

I was looking at the links below, but am still unsure of how to proceed with specifying a system from scratch using a hdf5 file.

https://github.com/NREL/PRAS/blob/master/SystemModel_HDF5_spec.md
https://juliaio.github.io/HDF5.jl/stable/

Could you guide me on how to specify a toy model through this workflow? I'm assuming this is the proper workflow to do larger-scale work.

Thanks,
Gui

Add more sample-level result tests

#50 fixes a really obvious error in sample-level surplus result processing that should have been caught by automated testing (in a compiled language the package wouldn't even have built). We need better test coverage for these more niche result specs.

Working with HDF5 files

I'm trying to better understand the workflow of working with HDF5 files. Would it be possible to show how you write the
toymodel.pras? Something along the lines of the Python (or Julia) code below. That way, I can adapt and start defining HDF5 files from scratch.
I realize the structure of the PRAS/HDF5 files is defined here (https://github.com/NREL/PRAS/blob/master/SystemModel_HDF5_spec.md, but still a toy model code would be very helpful.

Thanks again.

# Generate random data for recording
acc_1 = np.random.random(1000)
station_number_1 = '1'
# unix timestamp
start_time_1 = 1542000276
# time interval for recording
dt_1 = 0.04
location_1 = 'Berkeley'

acc_2 = np.random.random(500)
station_number_2 = '2'
start_time_2 = 1542000576
dt_2 = 0.01
location_2 = 'Oakland'

hf = h5py.File('station.hdf5', 'w')
hf['/acc/1/data'] = acc_1
hf['/acc/1/data'].attrs['dt'] = dt_1
hf['/acc/1/data'].attrs['start_time'] = start_time_1
hf['/acc/1/data'].attrs['location'] = location_1

hf['/acc/2/data'] = acc_2
hf['/acc/2/data'].attrs['dt'] = dt_2
hf['/acc/2/data'].attrs['start_time'] = start_time_2
hf['/acc/2/data'].attrs['location'] = location_2

hf['/gps/1/data'] = np.random.random(100)
hf['/gps/1/data'].attrs['dt'] = 60
hf['/gps/1/data'].attrs['start_time'] = 1542000000
hf['/gps/1/data'].attrs['location'] = 'San Francisco'
hf.close()

Allow choosing specific sample seeds

For more targeted analysis, allow specifying specific sample seeds to use, instead of just requiring every sample seed in 1:nsamples.

This has been implemented in PRAS-ED, but a more elegant approach would be to combine this with a broader refactoring of result collection to avoid allocating nthreads copies of each result accumulator (pass sample-level results into the result channel instead of entire result accumulators).

Implement Base.vcat for Storages, GeneratorStorages, and Lines structures

These functions will help with combining corresponding objects from multiple regions, and look similar to the following function:
function Base.vcat(gs::G...) where {N, L, T, P, G <: Generators{N,L,T,P}}

Question: Is the vision for PRAS to be exclusively a research tool or also to be used in the ISO and utility environments?

Hi, Gord
I'm an employee in the ISO-NE on the Resource Adequacy Team. Have begun trying to use PRAS for some research topics since I am also a graduate student. I am impressed that PRAS does Sequential MC. Would you say the vision/purpose of PRAS is eventually to be an open source modern scripting alternative to GE MARS that could be used in full-scale projects at ISOs?
Thanks for the great work!
Gui

assess(sys, NonSequentialMonteCarlo(), Shortfall())

Hi, Gord
Getting a bug on the NonSeqMC, I think. Not sure if it's my version of Julia, but here goes:

shortfall, = assess(sys, Convolution(), Shortfall())
eue_overall = EUE(shortfall)
lole_overall = LOLE(shortfall)
println(eue_overall)
println(lole_overall)

EUE = 0.00000 MWh/1440min
LOLE = 0.00000 event-(5min)/1440min

shortfall, = assess(sys, NonSequentialMonteCarlo(), Shortfall())
eue_overall = EUE(shortfall)
lole_overall = LOLE(shortfall)
println(eue_overall)
println(lole_overall)

UndefVarError: `NonSequentialMonteCarlo` not defined

Stacktrace:
 [1] top-level scope
   @ In[16]:1

Add SystemModel manipulation functions for common tasks

With the changes to simulation methods in PRAS 0.5, it should be easier for users to take a full-fidelity SystemModel and convert it to a simplified version. For example:

Collapse multiple regions into a single region
Extract a subset of regions from a multi-region system
Collapse all lines / generators with zero FORs into a single device per interface / region
Convert all (or specific?) generator-storage devices into decoupled generator and storage devices
Collapse all storage devices in a region into a single device (what to do for nonzero FORs?)

Metric display oddities

Looks like there are some (related?) bugs lurking in the metric show code (specifically roundresults). Most recently:

Error showing value of type ResourceAdequacy.SpatioTemporalResult{8784,1,Hour,MWh,Float64,Backcast,NonSequentialNetworkFlow}:
ERROR: ArgumentError: can't repeat a string -10 times
Stacktrace:
 [1] repeat(::String, ::Int64) at ./strings/substring.jl:176
 [2] scinote(::String) at /home/gstephen/.julia/packages/Decimals/GlFbH/src/decimal.jl:40
 [3] parse(::Type{Decimals.Decimal}, ::String) at /home/gstephen/.julia/packages/Decimals/GlFbH/src/decimal.jl:4
 [4] Type at /home/gstephen/.julia/packages/Decimals/GlFbH/src/decimal.jl:16 [inlined]
 [5] roundresults(::EUE{8784,1,Hour,MWh,Float64}) at /home/gstephen/.julia/packages/ResourceAdequacy/Rnyia/src/metrics/metrics.jl:18
 [6] show(::IOContext{REPL.Terminals.TTYTerminal}, ::EUE{8784,1,Hour,MWh,Float64}) at /home/gstephen/.julia/packages/ResourceAdequacy/Rnyia/src/metrics/EUE.jl:33
...

An error relating to invalid exponentiation of an integer has also come up. Both are likely triggered when trying to display a certain edge case value?

Trying to install PRAS, running into "Package PRAS not found", "ERROR: cannot find name corresponding to UUID 4138dd39-2aa7-5051-a626-17a0bb65d9c8 in a registry"

julia> using PRAS
│ Package PRAS not found, but a package named PRAS is available from a registry.
│ Install package?
│ (@v1.8) pkg> add PRAS
└ (y/n/o) [y]: y
Resolving package versions...
ERROR: cannot find name corresponding to UUID 4138dd39-2aa7-5051-a626-17a0bb65d9c8 in a registry
Stacktrace:
[1] pkgerror(msg::String)
@ Pkg.Types C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Types.jl:67
[2] deps_graph(env::Pkg.Types.EnvCache, registries::Vector{Pkg.Registry.RegistryInstance}, uuid_to_name::Dict{Base.UUID, String}, reqs::Dict{Base.UUID, Pkg.Versions.VersionSpec}, fixed::Dict{Base.UUID, Pkg.Resolve.Fixed}, julia_version::VersionNumber)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:488
[3] resolve_versions!(env::Pkg.Types.EnvCache, registries::Vector{Pkg.Registry.RegistryInstance}, pkgs::Vector{Pkg.Types.PackageSpec}, julia_version::VersionNumber)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:351
[4] targeted_resolve(env::Pkg.Types.EnvCache, registries::Vector{Pkg.Registry.RegistryInstance}, pkgs::Vector{Pkg.Types.PackageSpec}, preserve::Pkg.Types.PreserveLevel, julia_version::VersionNumber)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:1254
[5] tiered_resolve(env::Pkg.Types.EnvCache, registries::Vector{Pkg.Registry.RegistryInstance}, pkgs::Vector{Pkg.Types.PackageSpec}, julia_version::VersionNumber)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:1225
[6] _resolve(io::Base.TTY, env::Pkg.Types.EnvCache, registries::Vector{Pkg.Registry.RegistryInstance}, pkgs::Vector{Pkg.Types.PackageSpec}, preserve::Pkg.Types.PreserveLevel, julia_version::VersionNumber)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:1260
[7] add(ctx::Pkg.Types.Context, pkgs::Vector{Pkg.Types.PackageSpec}, new_git::Set{Base.UUID}; preserve::Pkg.Types.PreserveLevel, platform::Base.BinaryPlatforms.Platform)
@ Pkg.Operations C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\Operations.jl:1276
[8] add(ctx::Pkg.Types.Context, pkgs::Vector{Pkg.Types.PackageSpec}; preserve::Pkg.Types.PreserveLevel, platform::Base.BinaryPlatforms.Platform, kwargs::Base.Pairs{Symbol, Base.TTY, Tuple{Symbol}, NamedTuple{(:io,), Tuple{Base.TTY}}})
@ Pkg.API C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\API.jl:275
[9] add(pkgs::Vector{Pkg.Types.PackageSpec}; io::Base.TTY, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ Pkg.API C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\API.jl:156
[10] add(pkgs::Vector{Pkg.Types.PackageSpec})
@ Pkg.API C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\API.jl:145
[11] #add#27
@ C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\API.jl:144 [inlined]
[12] add
@ C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\API.jl:144 [inlined]
[13] try_prompt_pkg_add(pkgs::Vector{Symbol})
@ Pkg.REPLMode C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\Pkg\src\REPLMode\REPLMode.jl:717
[14] #invokelatest#2
@ .\essentials.jl:729 [inlined]
[15] invokelatest
@ .\essentials.jl:726 [inlined]
[16] check_for_missing_packages_and_run_hooks(ast::Any)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:176
[17] eval_user_input(ast::Any, backend::REPL.REPLBackend)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:146
[18] repl_backend_loop(backend::REPL.REPLBackend)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:247
[19] start_repl_backend(backend::REPL.REPLBackend, consumer::Any)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:232
[20] run_repl(repl::REPL.AbstractREPL, consumer::Any; backend_on_current_task::Bool)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:369
[21] run_repl(repl::REPL.AbstractREPL, consumer::Any)
@ REPL C:\Users\glara\AppData\Local\Programs\Julia-1.8.5\share\julia\stdlib\v1.8\REPL\src\REPL.jl:355
[22] (::Base.var"#967#969"{Bool, Bool, Bool})(REPL::Module)
@ Base .\client.jl:419
[23] #invokelatest#2
@ .\essentials.jl:729 [inlined]
[24] invokelatest
@ .\essentials.jl:726 [inlined]
[25] run_main_repl(interactive::Bool, quiet::Bool, banner::Bool, history_file::Bool, color_set::Bool)
@ Base .\client.jl:404
[26] exec_options(opts::Base.JLOptions)
@ Base .\client.jl:318
[27] _start()
@ Base .\client.jl:522

julia> using PRAS
│ Package PRAS not found, but a package named PRAS is available from a registry.
│ Install package?
│ (@v1.8) pkg> add PRAS
└ (y/n/o) [y]:n

(@v1.8) pkg> add PRAS
Resolving package versions...
ERROR: cannot find name corresponding to UUID 4138dd39-2aa7-5051-a626-17a0bb65d9c8 in a registry

Bump HDF5 version?

Is it okay to bump HDF5 version? Seems like we are restricted to 0.16, but Sienna ecosystem requires > 0.17. Would it be a breaking change to update to 0.17?

Tests passed with this environment:

  [abce61dc] Decimals v0.4.1
  [31c24e10] Distributions v0.25.107
  [f67ccb44] HDF5 v0.17.1
  [62286e6e] MinCostFlows v0.1.2
  [a15396b6] OnlineStats v1.6.3
  [925886fa] OnlineStatsBase v1.6.3
  [74087812] Random123 v1.7.0
  [189a3867] Reexport v1.2.2
  [2913bbd2] StatsBase v0.34.2
  [f269a46b] TimeZones v1.13.0
  [ade2ca70] Dates
  [de0858da] Printf
  [9a3f8284] Random

Prevent negative chargeable energy

If a storage device goes on outage with a positive state of charge, the calculated chargeable energy will be negative, resulting in a negative charge limit and an infeasible flow problem.

In these situations, chargeable energy should be set to zero.

Documentation of EFC arguments for capacity credit calculations

assess(EFC(UB,CI,MW range,DisreteNonParametric([list A],[list b])), EUE or LOLE, Backcast(), NonSequentialNetworkFlow(1000),Minimal(),original system,new system)

UB = the upper bound on the capacity value, often the nameplate capacity of the new unit

CI = the first of two conditions that stop the iterative calculation. Sets the confidence interval. If the reliability metrics (EUE or LOLE) calculated for both the upper and lower firm capacity bound are statistically identical within the confidence interval, the calculation stops.

MW range = the second of two conditions that stop the iterative calculation. Sets an absolute MW difference. If the reliability metrics calculated for both the upper and lower firm capacity bound are within this difference, the calculation stops.

list A = the list of regions in which the calculation places equivalent firm capacity

list B = the proportion of firm capacity the calculation places in each region. Must add up to 1.

Dramatic slowdown after upgrading to v0.6.0

To keep up with changes in v0.5.4 to v0.6.0, I made the following change to our code

      result_classic = assess(Convolution(), Temporal(), singlenode)

became

      result_classic, = assess(singlenode, Convolution(), Shortfall())

However, runtime increased by a factor of about 100x. This is an approximate timing because timing that specific line is a bit tricky -- and since it takes so long to run, I am seeing if I can avoid iterating just to measure for now.

Is this a known issue? Or did I inadvertently request additional computations alongside my syntax change? We have some comments in our code indicating that we tried to stick with "classic" because it was much faster than "modern", but I don't see anything in the NREL PRAS docs about classic vs modern. Does classic still exist?

Thank you!

Increment minimum supported Julia version to 1.6?

Julia v1.6 is now an LTS release. PRAS is not compatible with the previous 1.0 LTS release (since threaded coroutines are only available in 1.3+). Furthermore, the stack-allocated views introduced in v1.5 present some nice opportunities for cleaner code without performance penalties (we'd used views in some hot loops in the past, but moved to passing full arrays plus indices because of performance issues). So it might be nice to require at least Julia v1.6 for PRAS v0.7, whenever that may be.

Julia crashing when more than one storage is added

When running a system with 32 generators and 8 storages, I can get a result locally only when accounting for 1 storage in 3min. Otherwise with only 4threads, it is killing my julia session as soon as I had 2 storage devices (Storages{}).
Even when running in an environment with more multi-threading capability, adding storage running in Modern() makes it very long / computing heavy.

Is there any recommendation on how to run systems with multi-storages, ac-coupled, dc-coupled and unattached?
Do you have any runtime benchmark with storages? Thanks.

Suboptimal storage strategies with small times-to-go

When multiple short-duration devices have their times-to-go rounded to the nearest integer, they can result in equal TTG values that don't reflect the proper charging/discharging priorities for the devices, leading to suboptimal dispatch.

To fix this, scale up the floating-point TTGs before rounding to integer cost modifiers, to reduce quantization artifacts.

Disallow creating systems where lambda and mu are simultaneously zero

When lambda = mu = 0, long-run availability is undefined, which causes problems both for calculating convolutions as well as for initializing the system state in SequentialMonteCarlo. Resource constructors should include a check that at least one of the two values is nonzero at all times.

Sample variance correction?

Standard errors on Monte Carlo estimates are currently based on uncorrected (divide by n) sample variance. Should this be corrected instead (divide by n-1)?

[Moved over from internal NREL GitHub]

Add frequency calculation to Classic

Once we store unit-level FOR and MTTR (which we're doing for the sequential Monte Carlo work anyways), we should be able calculate LOLF as well as LOLP for non-sequential models.

[Moved over from internal NREL GitHub]

assess function failing with certain parameter combinations

The assess function fails when NonSequentialNetworkFlow() and Network() are used. ExpectedInterfaceFlow and ExpectedInterfaceUtilization work fine when SequentialCopperplate() is used, suggesting the problem is with assess().

Attached is a minimal working example.

##These analyses work:
works1 = assess(Backcast(), NonSequentialCopperplate(), Network(), mysystemmodel)
works2 = assess(REPRA(1,10), NonSequentialCopperplate(), Network(), mysystemmodel)
flow_1 = ExpectedInterfaceFlow(works1, "1", "2", DateTime(2020,4,27,13))
flow_2 = ExpectedInterfaceFlow(works2, "1", "2", DateTime(2020,4,27,13))
utilization1 = ExpectedInterfaceUtilization(works1, "1", "2", DateTime(2020,4,27,13))
utilization2 = ExpectedInterfaceUtilization(works2, "1", "2", DateTime(2020,4,27,13))

##These analyses don't work:
error1=assess(Backcast(), NonSequentialNetworkFlow(100), Network(), mysystemmodel)
error2=assess(REPRA(1, 10), NonSequentialNetworkFlow(100), Network(), mysystemmodel)

interface_errors.txt

It gives this error message:

# ERROR: NaN is not a valid interface utilization expectation
# Stacktrace:
#  [1] error(::String) at ./error.jl:33
#  [2] ExpectedInterfaceUtilization{1,1,Hour,V} where V<:Real(::Float64, ::Float64) at /Users/mschwarz/.julia/packages/ResourceAdequacy/L0wVg/src/metrics/ExpectedInterfaceUtilization.jl:8
#  [3] makemetric(::Type, ::OnlineStatsBase.Series{Number,Tuple{OnlineStatsBase.Mean{Float64,OnlineStatsBase.EqualWeight},OnlineStatsBase.Variance{Float64,OnlineStatsBase.EqualWeight}}}) at /Users/mschwarz/.julia/packages/ResourceAdequacy/L0wVg/src/utils/misc.jl:20
#  [4] _broadcast_getindex_evalf at ./broadcast.jl:578 [inlined]
#  [5] _broadcast_getindex at ./broadcast.jl:561 [inlined]
#  [6] getindex at ./broadcast.jl:511 [inlined]
#  [7] macro expansion at ./broadcast.jl:843 [inlined]
#  [8] macro expansion at ./simdloop.jl:73 [inlined]
#  [9] copyto! at ./broadcast.jl:842 [inlined]
#  [10] copyto! at ./broadcast.jl:797 [inlined]
#  [11] copy at ./broadcast.jl:773 [inlined]
#  [12] materialize at ./broadcast.jl:753 [inlined]
#  [13] finalize(::ResourceAdequacy.NonSequentialNetworkResultAccumulator{Float64,PRASBase.SystemModel{8784,1,Hour,MW,MWh,Float64},REPRA,NonSequentialNetworkFlow}) at /Users/mschwarz/.julia/packages/ResourceAdequacy/L0wVg/src/results/network_nonsequentialaccumulator.jl:155
#  [14] assess(::REPRA, ::NonSequentialNetworkFlow, ::Network, ::PRASBase.SystemModel{8784,1,Hour,MW,MWh,Float64}, ::UInt64) at /Users/mschwarz/.julia/packages/ResourceAdequacy/L0wVg/src/abstractspecs/simulation.jl:54
#  [15] assess(::REPRA, ::NonSequentialNetworkFlow, ::Network, ::PRASBase.SystemModel{8784,1,Hour,MW,MWh,Float64}) at /Users/mschwarz/.julia/packages/ResourceAdequacy/L0wVg/src/abstractspecs/simulation.jl:43
#  [16] top-level scope at none:0

Unit-level random seeds

Define unit-level random seeds to enable adding or removing resources without impacting the random outages associated with other resources. Currently, if a new unit is added, the simulation+sample seed will be used for the additional outage draws, changing which random numbers are associated with which units and thus the outage profiles of all units.

This would add a third level of seeding to a simulation:

Simulation seed: primary random seed for the entire simulation. Provided (or generated automatically) when the simulation is initiated.
Sample seed: secondary seed assigned to each sample from 1:nsamples. Generated automatically when makesamples sends simulation jobs to the threaded coroutines. Providing these directly (#36) would allow recreating specific random outage conditions
Unit seed: new, tertiary seed assigned to each unit. Fixed across runs and samples and defined in the system specification data (e.g. .pras files)

Initialize empty Generator, Storages, and GeneratorStorages objects

Useful for creating new regions which do not have a generator of each of these types.