Skip to content

Commit

Permalink
bespoke icing fixes from travis and paul
Browse files Browse the repository at this point in the history
  • Loading branch information
@grantbuster
grantbuster committed Oct 18, 2024
1 parent 5b94eb6 commit f424975
Show file tree
Hide file tree
Showing 2 changed files with 127 additions and 15 deletions.
74 changes: 59 additions & 15 deletions reV/bespoke/bespoke.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,8 @@ class BespokeMultiPlantData:
reads to a single HDF5 file.
"""

def __init__(self, res_fpath, sc_gid_to_hh, sc_gid_to_res_gid):
def __init__(self, res_fpath, sc_gid_to_hh, sc_gid_to_res_gid,
pre_load_humidity=False):
"""Initialize BespokeMultiPlantData
Parameters
Expand All @@ -75,6 +76,9 @@ def __init__(self, res_fpath, sc_gid_to_hh, sc_gid_to_res_gid):
GID values. Resource GID values should correspond to GID
values in the HDF5 file, so any GID map must be applied
before initializing :class`BespokeMultiPlantData`.
pre_load_humidity : optional, default=False
Option to pre-load relative humidity data (useful for icing
runs). If ``False``, relative humidities are not loaded.
"""
self.res_fpath = res_fpath
self.sc_gid_to_hh = sc_gid_to_hh
Expand All @@ -84,6 +88,8 @@ def __init__(self, res_fpath, sc_gid_to_hh, sc_gid_to_res_gid):
self._wind_speeds = None
self._temps = None
self._pressures = None
self._relative_humidities = None
self._pre_load_humidity = pre_load_humidity
self._time_index = None
self._pre_load_data()

Expand Down Expand Up @@ -117,6 +123,11 @@ def _pre_load_data(self):
for hh, gids in self.hh_to_res_gids.items()
}
self._time_index = res.time_index
if self._pre_load_humidity:
self._relative_humidities = {
hh: res["relativehumidity_2m", :, gids]
for hh, gids in self.hh_to_res_gids.items()
}

logger.debug(
f"Data took {(time.time() - start_time) / 60:.2f} " f"min to load"
Expand All @@ -139,13 +150,17 @@ def get_preloaded_data_for_gid(self, sc_gid):
hh = self.sc_gid_to_hh[sc_gid]
sc_point_res_gids = sorted(self.sc_gid_to_res_gid[sc_gid])
data_inds = np.searchsorted(self.hh_to_res_gids[hh], sc_point_res_gids)

rh = (None if not self._pre_load_humidity
else self._relative_humidities[hh][:, data_inds])
return BespokeSinglePlantData(
sc_point_res_gids,
self._wind_dirs[hh][:, data_inds],
self._wind_speeds[hh][:, data_inds],
self._temps[hh][:, data_inds],
self._pressures[hh][:, data_inds],
self._time_index,
rh,
)


Expand All @@ -158,7 +173,8 @@ class BespokeSinglePlantData:
"""

def __init__(
self, data_inds, wind_dirs, wind_speeds, temps, pressures, time_index
self, data_inds, wind_dirs, wind_speeds, temps, pressures, time_index,
relative_humidities=None,
):
"""Initialize BespokeSinglePlantData
Expand All @@ -169,22 +185,41 @@ def __init__(
the second dimension of `wind_dirs`, `wind_speeds`, `temps`,
and `pressures`. The GID value of data_inds[0] should
correspond to the `wind_dirs[:, 0]` data, etc.
wind_dirs, wind_speeds, temps, pressures : 2D np.array
Array of wind directions, wind speeds, temperatures, and
wind_dirs : 2D np.array
Array of wind directions. Dimensions should be correspond to
[time, location]. See documentation for `data_inds` for
required spatial mapping of GID values.
wind_speeds : 2D np.array
Array of wind speeds. Dimensions should be correspond to
[time, location]. See documentation for `data_inds` for
required spatial mapping of GID values.
temps : 2D np.array
Array oftemperatures. Dimensions should be correspond to
[time, location]. See documentation for `data_inds` for
required spatial mapping of GID values.
pressures : 2D np.array
Array of pressures. Dimensions should be correspond to
pressures, respectively. Dimensions should be correspond to
[time, location]. See documentation for `data_inds` for
required spatial mapping of GID values.
time_index : 1D np.array
Time index array corresponding to the temporal dimension of
the 2D data. Will be exposed directly to user.
relative_humidities : 2D np.array, optional
Array of relative humidities. Dimensions should be
correspond to [time, location]. See documentation for
`data_inds` for required spatial mapping of GID values.
If ``None``, relative_humidities cannot be queried.
"""

self.data_inds = data_inds
self.wind_dirs = wind_dirs
self.wind_speeds = wind_speeds
self.temps = temps
self.pressures = pressures
self.time_index = time_index
self.relative_humidities = relative_humidities
self._humidities_exist = relative_humidities is not None

def __getitem__(self, key):
dset_name, t_idx, gids = key
Expand All @@ -197,6 +232,8 @@ def __getitem__(self, key):
return self.temps[t_idx, data_inds]
if "pressure" in dset_name:
return self.pressures[t_idx, data_inds]
if self._humidities_exist and "relativehumidity" in dset_name:
return self.relative_humidities[t_idx, data_inds]
msg = f"Unknown dataset name: {dset_name!r}"
logger.error(msg)
raise ValueError(msg)
Expand Down Expand Up @@ -711,6 +748,7 @@ def get_weighted_res_ts(self, dset):
weights[i] = self.sc_point.include_mask_flat[mask].sum()

weights /= weights.sum()
data = data.astype(np.float32)
data *= weights
data = np.sum(data, axis=1)

Expand Down Expand Up @@ -896,15 +934,18 @@ def res_df(self):
if np.nanmax(pres) > 1000:
pres *= 9.86923e-6

self._res_df = pd.DataFrame(
{
"temperature": temp,
"pressure": pres,
"windspeed": ws,
"winddirection": wd,
},
index=ti,
)
data = {
"temperature": temp,
"pressure": pres,
"windspeed": ws,
"winddirection": wd,
}

if self.sam_sys_inputs.get("en_icing_cutoff"):
rh = self.get_weighted_res_ts("relativehumidity_2m")
data["relativehumidity"] = rh

self._res_df = pd.DataFrame(data, index=ti)

if "time_index_step" in self.original_sam_sys_inputs:
ti_step = self.original_sam_sys_inputs["time_index_step"]
Expand Down Expand Up @@ -2156,7 +2197,10 @@ def _pre_load_data(self, pre_load_data):

logger.info("Pre-loading resource data for Bespoke run... ")
self._pre_loaded_data = BespokeMultiPlantData(
self._res_fpath, sc_gid_to_hh, sc_gid_to_res_gid
self._res_fpath,
sc_gid_to_hh,
sc_gid_to_res_gid,
pre_load_humidity=self._project_points.sam_config_obj.icing,
)

def _hh_for_sc_gid(self, sc_gid):
Expand Down
68 changes: 68 additions & 0 deletions tests/test_bespoke.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1663,3 +1663,71 @@ def test_bespoke_5min_sample():
assert len(f["time_index-2010"]) == 8760
assert len(f["windspeed-2010"]) == 8760
assert len(f["winddirection-2010"]) == 8760


def test_bespoke_run_with_icing_cutoff():
"""Test bespoke run with icing cutoff enabled."""
output_request = ("system_capacity", "cf_mean", "cf_profile")
with tempfile.TemporaryDirectory() as td:
res_fp = os.path.join(td, "ri_100_wtk_{}.h5")
excl_fp = os.path.join(td, "ri_exclusions.h5")
shutil.copy(EXCL, excl_fp)
shutil.copy(RES.format(2012), res_fp.format(2012))
shutil.copy(RES.format(2013), res_fp.format(2013))
res_fp = res_fp.format("*")

TechMapping.run(excl_fp, RES.format(2012), dset=TM_DSET, max_workers=1)
bsp = BespokeSinglePlant(
33,
excl_fp,
res_fp,
TM_DSET,
SAM_SYS_INPUTS,
OBJECTIVE_FUNCTION,
CAP_COST_FUN,
FOC_FUN,
VOC_FUN,
BOS_FUN,
ga_kwargs={"max_time": 5},
excl_dict=EXCL_DICT,
output_request=output_request,
)

out = bsp.run_plant_optimization()
out = bsp.run_wind_plant_ts()
bsp.close()

sam_inputs_ice = copy.deepcopy(SAM_SYS_INPUTS)
sam_inputs_ice["en_icing_cutoff"] = 1
sam_inputs_ice["en_low_temp_cutoff"] = 1
sam_inputs_ice["icing_cutoff_rh"] = 90 # High values to ensure diff
sam_inputs_ice["icing_cutoff_temp"] = 10
sam_inputs_ice["low_temp_cutoff"] = 0
bsp = BespokeSinglePlant(
33,
excl_fp,
res_fp,
TM_DSET,
sam_inputs_ice,
OBJECTIVE_FUNCTION,
CAP_COST_FUN,
FOC_FUN,
VOC_FUN,
BOS_FUN,
ga_kwargs={"max_time": 5},
excl_dict=EXCL_DICT,
output_request=output_request,
)

out_ice = bsp.run_plant_optimization()
out_ice = bsp.run_wind_plant_ts()
bsp.close()

ae_dsets = [
"annual_energy-2012",
"annual_energy-2013",
"annual_energy-means",
]
for dset in ae_dsets:
assert not np.isclose(out[dset], out_ice[dset])
assert out[dset] > out_ice[dset]

0 comments on commit f424975

Please sign in to comment.