diff --git a/README.md b/README.md index d0d6690..f5a26c0 100644 --- a/README.md +++ b/README.md @@ -1,5 +1,35 @@ # Functional-Recovery---Python This is translation of Matlab codebase into Python for quantifying building-specific functional recovery and reoccupancy based on a probabilistic performance-based earthquake engineering framework. + +## Requirements +- The `requirements.txt` file defines the Python package dependencies required to run this codebase. Follow the instructions below to install all required depenedancies listed in the 'requirements.txt' file. +- Recommended Python version: `3.9` (the codebase was developed and tested with Python 3.9). + +Installation (using a virtual environment is recommended): + +```powershell +# create a virtual environment +python -m venv .venv + +# activate the virtual environment (PowerShell) +.\.venv\Scripts\Activate.ps1 + +# upgrade pip (optional but recommended) +python -m pip install --upgrade pip + +# install dependencies from requirements.txt +pip install -r requirements.txt +``` + +If you prefer conda: + +```bash +conda create -n frec python=3.9 +conda activate frec +pip install -r requirements.txt +``` + +If you run into platform-specific dependency issues, please refer to the package error messages and install any missing system libraries before re-running `pip install -r requirements.txt`. Original Matlab code is from Dr. Dustin Cook's Github directory https://github.com/OpenPBEE/PBEE-Recovery. ### Method Description diff --git a/fn_red_tag.py b/fn_red_tag.py index 939ea87..18603f5 100644 --- a/fn_red_tag.py +++ b/fn_red_tag.py @@ -1,4 +1,4 @@ -def fn_red_tag( calculate_red_tag, damage, comps, simulated_replacement_time): +def fn_red_tag( calculate_red_tag, damage, comps, simulated_replacement_time, red_tag_options ): '''Perform the ATC-138 functional recovery time assessement given similation of component damage for a single shaking intensity @@ -27,9 +27,10 @@ def fn_red_tag( calculate_red_tag, damage, comps, simulated_replacement_time): import numpy as np - def simulate_tagging(damage, comps, sc_ids, sc_thresholds): + def simulate_tagging(damage, comps, sc_ids, sc_thresholds, red_tag_options ): red_tag_impact = np.zeros(np.shape(damage['tenant_units'][0]['qnt_damaged'])) # num reals by num comp_ds + red_tag_impact_cb = np.zeros(np.shape(damage['tenant_units'][0]['qnt_damaged'])) # num reals by num comp_ds, specific for counting red tags for coupling beams num_reals = len(red_tag_impact) # Go through each structural system and calc the realizations where the @@ -37,8 +38,12 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): sc_tag = np.zeros([num_reals,len(sc_ids)]) for sc in range(len(sc_ids)): - story_tag = np.zeros([num_reals,len(damage['story'])]) + + # resets the coupling beam counter for every safety class + coupling_beam_dmg = np.zeros([num_reals,3]) # three directions + coupling_beam_qty = np.zeros([1,3]) # three directions + red_tag_impact_cb_sc = np.zeros(np.shape(damage['tenant_units'][0]['qnt_damaged'])) # num reals by num comp_ds for s in range(len(damage['story'])): sc_filt = np.array(damage['comp_ds_table']['safety_class']) >= sc_ids[sc] @@ -54,6 +59,10 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): sys_tag = np.zeros([num_reals, len(structural_systems)]).astype(bool) + if red_tag_options['ignore_coupling_beam_for_red_tag']: + # just ignore system 12 (coupling beams) + structural_systems = structural_systems[structural_systems != 12] + for sys in range(len(structural_systems)): ss_filt_ds = np.logical_or(np.array(damage['comp_ds_table']['structural_system']) == structural_systems[sys], np.array(damage['comp_ds_table']['structural_system_alt']) == structural_systems[sys]) ss_filt_comp = np.logical_or(np.array(comps['comp_table']['structural_system']) == structural_systems[sys], np.array(comps['comp_table']['structural_system_alt']) == structural_systems[sys]) @@ -61,7 +70,7 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): # Check damage among each series within this structural system series = np.unique(np.array(damage['comp_ds_table']['structural_series_id'])[ss_filt_ds]) ser_dmg = np.zeros([num_reals,len(series)]) - ser_qty = np.zeros([num_reals,len(series)]) + ser_qty = np.zeros([1,len(series)]) for ser in range(len(series)): ser_filt_ds = np.array(damage['comp_ds_table']['structural_series_id']) == series[ser] ser_filt_comp = np.array(comps['comp_table']['structural_series_id']) == series[ser] @@ -73,17 +82,24 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): ser_qty[:,ser] = np.sum(num_comps[ser_filt_comp & ss_filt_comp]) # Check if this system is causing a red tag - sys_dmg = np.nanmax(ser_dmg, axis = 1) - sys_qty = np.nanmax(ser_qty, axis = 1) - sys_ratio = sys_dmg / sys_qty - sys_tag[:,sys] = sys_ratio > sc_thresholds[sc] - - '''Calculate the impact that each component has on red tag - (boolean, 1 = affects red tag, 0 = does not affect) - Take all damage that is part of this system at this story - in this direction that is damaged to this safety class - level, only where damage exceeds tagging threshold''' - red_tag_impact = np.fmax(red_tag_impact, 1*sys_tag[:,sys].reshape(len(sys_tag),1) * ss_filt_ds.reshape(1,len(ss_filt_ds)) * sc_filt.reshape(1,len(ss_filt_ds)) * (sc_dmg>0)) + if structural_systems[sys] == 12 and bool(red_tag_options.get('tag_coupling_beams_over_height', False)): + # HARED CODED CHECK FOR COUPLING BEAMS + # just do counts for now instead of adding to red tag per story check + coupling_beam_dmg[:, direc-1] = coupling_beam_dmg[:,direc-1] + np.nanmax(ser_dmg, axis = 1) + coupling_beam_qty[0, direc-1] = coupling_beam_qty[0,direc-1] + np.nanmax(ser_qty, axis = 1) + red_tag_impact_cb_sc = np.fmax(red_tag_impact_cb_sc, ss_filt_ds.reshape(1,-1) * sc_filt.reshape(1,-1) * (sc_dmg>0)) + else: + sys_dmg = np.nanmax(ser_dmg, axis = 1) + sys_qty = np.nanmax(ser_qty, axis = 1) + sys_ratio = sys_dmg / sys_qty + sys_tag[:,sys] = sys_ratio > sc_thresholds[sc] + + '''Calculate the impact that each component has on red tag + (boolean, 1 = affects red tag, 0 = does not affect) + Take all damage that is part of this system at this story + in this direction that is damaged to this safety class + level, only where damage exceeds tagging threshold''' + red_tag_impact = np.fmax(red_tag_impact, 1*sys_tag[:,sys].reshape(-1,1) * ss_filt_ds.reshape(1,-1) * sc_filt.reshape(1,-1) * (sc_dmg>0)) # Combine across all systems in this direction @@ -96,7 +112,19 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): # Combine across all stories for this safety class sc_tag[:,sc] = np.nanmax(story_tag, axis = 1) + + # Additional Hard Coded Check for Coupling Beams + if np.any(np.array(damage['comp_ds_table']['structural_system']) == 12): + CB_damage_ratio = coupling_beam_dmg / coupling_beam_qty + CB_tag_sc = np.nanmax(CB_damage_ratio, axis = 1) > sc_thresholds[sc] + sc_tag[:,sc] = np.fmax(sc_tag[:,sc], CB_tag_sc) # add it back to sum + red_tag_impact_cb = np.fmax( + red_tag_impact_cb, + red_tag_impact_cb_sc * sc_tag[:, sc].reshape(-1, 1) + ) + red_tag_impact = np.fmax(red_tag_impact, red_tag_impact_cb) + # Combine all safety class checks into one simulated red tag red_tag = np.nanmax(sc_tag, axis = 1) @@ -116,10 +144,11 @@ def simulate_tagging(damage, comps, sc_ids, sc_thresholds): # Simulate Red Tags sc_ids = np.array([1, 2, 3, 4]) sc_thresholds = np.array([0.5, 0.25, 0.1, 0]) - red_tag, red_tag_impact = simulate_tagging(damage, comps, sc_ids, sc_thresholds) + red_tag, red_tag_impact = simulate_tagging(damage, comps, sc_ids, sc_thresholds, red_tag_options ) # Inspection is flagged for 50% of the red tag thresholds - inspection_tag = simulate_tagging(damage, comps, sc_ids, 0.5*sc_thresholds)[0] + inspection_thresholds = 0.5*sc_thresholds + inspection_tag, _ = simulate_tagging(damage, comps, sc_ids, inspection_thresholds, red_tag_options ) else: diff --git a/functionality/fn_calculate_reoccupancy.py b/functionality/fn_calculate_reoccupancy.py index 6aa1b33..6af1dff 100644 --- a/functionality/fn_calculate_reoccupancy.py +++ b/functionality/fn_calculate_reoccupancy.py @@ -67,18 +67,24 @@ def fn_calculate_reoccupancy(damage, damage_consequences, utilities, recovery_day['tenant_safety'], comp_breakdowns['tenant_safety'] = other_functionality_functions.fn_tenant_safety( damage, building_model, functionality_options, tenant_units) ## Combine Check to determine the day the each tenant unit is reoccupiable - # Check the day the building is safe - day_building_safe = np.fmax(recovery_day['building_safety']['red_tag'], - np.fmax(recovery_day['building_safety']['shoring'], - np.fmax(recovery_day['building_safety']['hazardous_material'], - np.fmax(recovery_day['building_safety']['entry_door_access'], - recovery_day['building_safety']['fire_suppression'])))) - # Check the day each story is accessible - day_story_accessible = np.fmax(recovery_day['story_access']['stairs'], np.fmax(recovery_day['story_access']['stair_doors'], np.fmax(recovery_day['story_access']['flooding'], recovery_day['story_access']['horizontal_egress']))) - - # Check the day each tenant unit is safe - day_tenant_unit_safe = np.fmax(recovery_day['tenant_safety']['interior'], np.fmax(recovery_day['tenant_safety']['exterior'], recovery_day['tenant_safety']['hazardous_material'])) + # Go through each of the building safety checks and combine them to check the day the building is safe (max of all checks) + fault_tree_events_building_safety = list(recovery_day['building_safety'].keys()) + day_building_safe = recovery_day['building_safety'][fault_tree_events_building_safety[0]] + for i in range(1, len(fault_tree_events_building_safety)): + day_building_safe = np.fmax(day_building_safe, recovery_day['building_safety'][fault_tree_events_building_safety[i]]) + + # Go through each of the story access checks and combine them to check the day each story is accessible (max of all checks) + fault_tree_events_story_access = list(recovery_day['story_access'].keys()) + day_story_accessible = recovery_day['story_access'][fault_tree_events_story_access[0]] + for i in range(1, len(fault_tree_events_story_access)): + day_story_accessible = np.fmax(day_story_accessible, recovery_day['story_access'][fault_tree_events_story_access[i]]) + # Go through each of the tenant unit safety checks and combine them to check the day each tenant unit is safe (max of all checks) + fault_tree_events_tenant_unit_safe = list(recovery_day['tenant_safety'].keys()) + day_tenant_unit_safe = recovery_day['tenant_safety'][fault_tree_events_tenant_unit_safe[0]] + for i in range(1, len(fault_tree_events_tenant_unit_safe)): + day_tenant_unit_safe = np.fmax(day_tenant_unit_safe, recovery_day['tenant_safety'][fault_tree_events_tenant_unit_safe[i]]) + # Combine checks to determine when each tenant unit is re-occupiable day_tenant_unit_reoccupiable = np.fmax(np.fmax(day_building_safe.reshape(len(day_building_safe),1), day_story_accessible), day_tenant_unit_safe) diff --git a/functionality/fn_check_habitability.py b/functionality/fn_check_habitability.py index d4e844a..3372ee2 100644 --- a/functionality/fn_check_habitability.py +++ b/functionality/fn_check_habitability.py @@ -41,9 +41,9 @@ def fn_check_habitability( damage, damage_consequences, reoc_meta, func_meta, # Go through each of the tenant function branches and combines checks day_tenant_unit_reoccupiable = np.zeros([num_reals,1]) - fault_tree_events_LV1 = list(comp_breakdowns.keys()) + fault_tree_events_LV1 = list(recovery_day.keys()) for i in range(len(fault_tree_events_LV1)): - fault_tree_events_LV2 = list(comp_breakdowns[fault_tree_events_LV1[i]].keys()) + fault_tree_events_LV2 = list(recovery_day[fault_tree_events_LV1[i]].keys()) for j in range(len(fault_tree_events_LV2)): day_tenant_unit_reoccupiable = np.fmax(day_tenant_unit_reoccupiable, recovery_day[fault_tree_events_LV1[i]][fault_tree_events_LV2[j]].reshape(num_reals, int(np.size(recovery_day[fault_tree_events_LV1[i]][fault_tree_events_LV2[j]])/num_reals))) diff --git a/inputs/Inputs2Copy/optional_inputs.py b/inputs/Inputs2Copy/optional_inputs.py index f6be343..34cbebe 100644 --- a/inputs/Inputs2Copy/optional_inputs.py +++ b/inputs/Inputs2Copy/optional_inputs.py @@ -71,6 +71,10 @@ "calculate_red_tag" : 1, "red_tag_clear_time" : 7, "red_tag_clear_beta" : 0.6, +"red_tag_options" : { + "tag_coupling_beams_over_height" : True, + "ignore_coupling_beam_for_red_tag" : False + }, "include_local_stability_impact" : 1, "include_flooding_impact": 1, "egress_threshold" : 0.5, @@ -80,14 +84,15 @@ "exterior_safety_threshold" : 0.1, "interior_safety_threshold" : 0.25, "door_access_width_ft" : 9, -"habitability_requirements": {"electrical" : 0, - "water_potable" : 0, - "water_sanitary" : 0, - "hvac_ventilation" : 0, - "hvac_heating" : 0, - "hvac_cooling" : 0, - "hvac_exhaust" : 0 - }, +"habitability_requirements": { + "electrical" : 0, + "water_potable" : 0, + "water_sanitary" : 0, + "hvac_ventilation" : 0, + "hvac_heating" : 0, + "hvac_cooling" : 0, + "hvac_exhaust" : 0 + }, "water_pressure_max_story" : 4, "heat_utility" : 'gas', } diff --git a/main_PBEE_recovery.py b/main_PBEE_recovery.py index 30cb66f..dce3288 100644 --- a/main_PBEE_recovery.py +++ b/main_PBEE_recovery.py @@ -71,7 +71,8 @@ def main_PBEE_recovery(damage, damage_consequences, building_model, ## Calculate Red Tags RT, RTI, IT = fn_red_tag(functionality_options['calculate_red_tag'], damage, building_model['comps'], - np.array(damage_consequences['simulated_replacement_time'])) + np.array(damage_consequences['simulated_replacement_time']), + functionality_options['red_tag_options']) damage_consequences['red_tag'] = RT damage_consequences['red_tag_impact'] = RTI diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..af3fc8b --- /dev/null +++ b/requirements.txt @@ -0,0 +1,5 @@ +numpy +pandas +scipy +matplotlib +seaborn \ No newline at end of file