import pandas as pd
import networkx as nx
import numpy as np
import itertools
# --- Constants and Assumptions ---
# These should be clearly stated and can be modified.
VOLTAGE_KV = 10.0 # Line voltage in kV
ROOT_3 = np.sqrt(3)
BASE_DG_CAPACITY_KW = 300.0 # Initial capacity for each DG
N_DG = 8
# Failure rates from problem description
FAILURE_RATE_DG_PERCENT = 0.5 / 100.0
# FAILURE_RATE_USER_PERCENT = 0.5 / 100.0 # Not directly used in this simplified line-fault model for widespread outages
# FAILURE_RATE_SWITCH_PERCENT = 0.2 / 100.0 # Assuming switch failures manifest as line failures or inability to operate tie
FAILURE_RATE_LINE_PER_KM = 0.002 # Per km per year (assuming rates are annual)
# Costs (placeholders - these are critical for actual risk values)
# Value of Lost Load (VoLL) in monetary units per kW per hour.
# For risk = P * C, if P is annual probability, C should be impact of one event.
# Let's define C_loss as total kW unserved * a severity factor.
# Or, if we want an annual risk cost: P_annual_fault * kW_unserved * hours_outage * cost_per_kWh
# For simplicity, using $/kW of unserved load for the consequence C.
COST_VOLL_PER_KW = 10.0 # Example: $10 per kW of unserved load
AVG_OUTAGE_DURATION_H = 4 # Example: average hours for an outage, if converting to energy
# Cost of Overload (Consequence C_over)
# This can be complex: accelerated aging, tripping, damage.
# Simplified: A penalty if any line is overloaded in a given state.
COST_PENALTY_FOR_ANY_OVERLOAD = 1000.0 # Example: $1000 penalty if system is in an overloaded state
# Or, a cost per MWh of overloaded energy, or per overloaded line.
# Line and Feeder Capacities
# Main feeder rated current from problem: 220A.
# P_rated_feeder_kW = ROOT_3 * VOLTAGE_KV * FEEDER_RATED_CURRENT_A * 1.0 (pf=1)
# = 1.732 * 10 * 220 = 3810.4 kW (approx 3.8 MW, problem says 2.2MW for 220A, implies lower pf or different basis)
# Let's use current as the primary limit.
FEEDER_RATED_CURRENT_A = 220.0
# Assumption for individual line segments: For this model, we'll assume all lines
# have a rated current equal to the main feeder. This is a strong simplification.
# A more detailed model would assign ratings based on conductor types or downstream load.
LINE_RATED_CURRENT_A = 100.0 # More conservative assumption for individual segments than 220A. Needs proper engineering values.
# For lines directly from substation, perhaps 220A is more appropriate.
# Let's use a dictionary for specific line ratings if known, else default.
DEFAULT_LINE_RATED_CURRENT_A = 100.0
# Tie Line Capacity
TIE_LINE_RATED_CURRENT_A = 150.0 # Assumption, should be based on tie switch/line capacity
# DG Locations (Node IDs from 1 to 62) - Based on Figure 1 interpretation
DG_LOCATIONS_KW = {
6: BASE_DG_CAPACITY_KW, 10: BASE_DG_CAPACITY_KW, 15: BASE_DG_CAPACITY_KW,
27: BASE_DG_CAPACITY_KW, 31: BASE_DG_CAPACITY_KW, 37: BASE_DG_CAPACITY_KW,
50: BASE_DG_CAPACITY_KW, 58: BASE_DG_CAPACITY_KW
}
# Tie Switches: (node1, node2, switch_id_text) - normally open
# Interpretation based on careful review of Figure 1:
# S13-1: (13, 22) - Intra-Feeder 1 (Connects two branches of Feeder 1)
# S29-2: (29, 42) - Intra-Feeder 2 (Connects two branches of Feeder 2)
# S62-3: (62, 19) - Inter-Feeder (Connects Feeder 3 (node 62) to Feeder 1 (node 19))
TIE_SWITCHES_INFO = [
{'nodes': (13, 22), 'id': 'S13-1', 'type': 'intra-F1', 'capacity_A': TIE_LINE_RATED_CURRENT_A},
{'nodes': (29, 42), 'id': 'S29-2', 'type': 'intra-F2', 'capacity_A': TIE_LINE_RATED_CURRENT_A},
{'nodes': (62, 19), 'id': 'S62-3', 'type': 'inter-F3_F1', 'capacity_A': TIE_LINE_RATED_CURRENT_A}
]
# This interpretation means Feeder 2 cannot directly receive support from F1 or F3.
# If problem implies all feeders can support each other, TIE_SWITCHES_INFO would need redefinition.
# Substation connection points (source nodes for feeders)
# CB1 -> Node 1, CB2 -> Node 23, CB3 -> Node 43
# Node 0 will represent the main grid / infinite source.
SOURCE_NODE = 0
SUBSTATION_CONNECTIONS = {
'CB1': (SOURCE_NODE, 1),
'CB2': (SOURCE_NODE, 23),
'CB3': (SOURCE_NODE, 43)
}
# Capacity of connection from source to substation nodes (effectively feeder capacity)
SUBSTATION_LINE_CAPACITY_A = FEEDER_RATED_CURRENT_A
# --- Data Loading Functions ---
def load_load_data(filename="C题附件：有源配电网62节点系统基本参数.xlsx - 表1 有源配电网62节点系统负荷参数.csv"):
df = pd.read_csv(filename)
df.columns = ['node_id', 'load_kw']
# Convert node_id to int if it's not already
df['node_id'] = df['node_id'].astype(int)
return df.set_index('node_id')['load_kw'].to_dict()
def load_topology_data(filename="C题附件：有源配电网62节点系统基本参数.xlsx - 表2 有源配电网62节点系统拓扑参数.csv"):
df = pd.read_csv(filename)
# Rename columns for easier access (assuming standard Chinese headers)
df.columns = ['line_num', 'from_node', 'to_node', 'length_km', 'resistance_ohm', 'reactance_ohm']
# Convert relevant columns to numeric
for col in ['from_node', 'to_node', 'length_km', 'resistance_ohm', 'reactance_ohm']:
df[col] = pd.to_numeric(df[col], errors='coerce')
return df
# --- Core Power Grid Model Class ---
class PowerGridModel:
def __init__(self, load_data, topology_data, dg_locations_kw, tie_switches_info, substation_connections):
self.loads_kw = load_data
self.topology_df = topology_data
self.dg_kw = dg_locations_kw.copy() # Allow modification for different scenarios
self.tie_switches_info = tie_switches_info
self.substation_connections = substation_connections
self.graph = self._build_graph()
self.feeder_info = self._identify_feeders()
def _build_graph(self):
G = nx.Graph() # Use Graph for undirected, or DiGraph if flow direction is fixed by sources
# Add nodes with load and DG info
all_nodes = set(self.topology_df['from_node']) | set(self.topology_df['to_node'])
for node_id in all_nodes:
node_id = int(node_id) # Ensure int
G.add_node(node_id,
load_kw=self.loads_kw.get(node_id, 0),
dg_kw=self.dg_kw.get(node_id, 0))
# Add lines from topology data
for _, row in self.topology_df.iterrows():
u, v = int(row['from_node']), int(row['to_node'])
G.add_edge(u, v,
id=row['line_num'],
length_km=row['length_km'],
resistance_ohm=row['resistance_ohm'],
# reactance_ohm=row['reactance_ohm'], # Ignoring reactance as per problem
rated_current_a=DEFAULT_LINE_RATED_CURRENT_A, # Default, can be refined
failed=False)
# Add substation connections (virtual lines from a common source)
# These represent the main feeder lines from CBs
G.add_node(SOURCE_NODE, type='source')
for cb_id, (src, dest_node) in self.substation_connections.items():
G.add_edge(src, dest_node, id=cb_id, length_km=0.01, resistance_ohm=0.001, # Minimal impedance
rated_current_a=SUBSTATION_LINE_CAPACITY_A, type='substation_link', failed=False)
return G
def _get_subgraph_with_operational_lines(self, graph_to_copy, faulty_line_edge=None):
"""Creates a subgraph considering only non-failed lines and open tie switches."""
g_op = graph_to_copy.copy()
# Remove failed lines
lines_to_remove = []
if faulty_line_edge: # faulty_line_edge is (u,v)
if g_op.has_edge(*faulty_line_edge):
lines_to_remove.append(faulty_line_edge)
for u, v, data in list(g_op.edges(data=True)):
if data.get('failed', False):
lines_to_remove.append((u,v))
g_op.remove_edges_from(lines_to_remove)
# Normally, tie switches are open. For restoration, specific ones might be closed.
# This base function assumes they are open unless explicitly handled by restoration logic.
return g_op
def _identify_feeders(self):
"""Identifies nodes belonging to each feeder under normal operation (tie switches open)."""
g_normal = self._get_subgraph_with_operational_lines(self.graph)
feeder_info = {} # {'CB1': {nodes}, 'CB2': {nodes}, ...}
for cb_id, (src_node, start_node) in self.substation_connections.items():
if g_normal.has_node(start_node) and g_normal.has_node(src_node) and nx.has_path(g_normal, src_node, start_node):
# Find all nodes reachable from start_node without passing through another substation's start_node
# or the main source node again, after removing other substation links.
temp_g = g_normal.copy()
other_cb_links = []
for other_cb, (s,d) in self.substation_connections.items():
if other_cb != cb_id and temp_g.has_edge(s,d):
other_cb_links.append((s,d))
temp_g.remove_edges_from(other_cb_links)
if nx.has_path(temp_g, src_node, start_node):
# All nodes in the component connected to start_node, excluding the source itself
component_nodes = nx.node_connected_component(temp_g.subgraph(
[n for n in temp_g.nodes if n != src_node or n == start_node] # Consider start_node part of feeder
), start_node)
feeder_info[cb_id] = component_nodes
else: # Should not happen if graph is built correctly
feeder_info[cb_id] = {start_node} if start_node in g_normal else set()
else:
feeder_info[cb_id] = {start_node} if start_node in g_normal else set()
return feeder_info
def _calculate_line_current_kw(self, power_kw):
"""Calculates current (A) given power (kW) at VOLTAGE_KV (line-to-line)."""
if VOLTAGE_KV <= 0: return float('inf')
return abs(power_kw) / (ROOT_3 * VOLTAGE_KV * 1.0) # Assumed PF=1 for current calculation from P
def _get_downstream_info(self, G, line_u, line_v, source_nodes_for_feeder):
"""
Calculates total load and DG power downstream of a directed line (u,v),
assuming v is further from the source_node for this path.
G: graph to operate on (can be a faulted graph)
source_nodes_for_feeder: list of possible source nodes for the current connected component.
"""
# Temporarily make graph directed from source to loads to find downstream nodes
# This is tricky if the graph is not purely radial or has loops after closing ties.
# A simpler approach for radial sections:
# Check connectivity from sources to line_v, if line_u is removed.
# If line_v is disconnected from all sources when (u,v) is cut, then everything
# in the component of line_v is downstream.
# Create a copy of G to modify
temp_g = G.copy()
if not temp_g.has_edge(line_u, line_v):
return [], 0, 0 # Line doesn't exist
temp_g.remove_edge(line_u, line_v)
downstream_nodes = set()
# Check which side (u or v) is disconnected from the source(s)
# Assume u is closer to source, v is further.
# If v is still connected to a source, then (u,v) might be part of a loop or fed from elsewhere.
# A robust way: find path from source to v. If (u,v) is on all paths, then v is downstream of u via this line.
# For radial feeders (normal operation):
# If we consider (u,v) where u is parent of v:
component_of_v = set()
q = [line_v]
visited = {line_u, line_v} # Start by marking u as visited (as if coming from u)
# If line_v is connected to any source_node without passing through line_u
v_connected_to_source_alt_path = False
for src in source_nodes_for_feeder:
if nx.has_path(temp_g, src, line_v):
v_connected_to_source_alt_path = True
break
if v_connected_to_source_alt_path: # (u,v) is part of a loop or v is fed from elsewhere
# This simple downstream logic is insufficient for meshed networks.
# For now, assume radial for this part of flow calculation.
# A more complex flow calculation (Newton-Raphson) would be needed for meshed.
# Given problem constraints, assume feeders are normally radial.
pass # This line might not have a clear "downstream" if looped.
# If v is disconnected from source when (u,v) is cut, then its component is downstream.
# Check connectivity for v in temp_g (where (u,v) is removed)
v_still_connected = any(nx.has_path(temp_g, src, line_v) for src in source_nodes_for_feeder if src in temp_g)
if not v_still_connected: # v is now isolated from source, so its component is downstream of (u,v)
component_of_v = nx.node_connected_component(temp_g, line_v)
else: # v is still connected, means (u,v) might be redundant or complex.
# Try to determine direction based on distance from source
dist_u = float('inf')
dist_v = float('inf')
for src in source_nodes_for_feeder:
if src not in G: continue
if nx.has_path(G, src, line_u): dist_u = min(dist_u, nx.shortest_path_length(G, src, line_u))
if nx.has_path(G, src, line_v): dist_v = min(dist_v, nx.shortest_path_length(G, src, line_v))
if dist_v > dist_u : # v is downstream of u
# Find component of v if (u,v) is removed and v is not connected to source
g_temp_removed_edge = G.copy()
g_temp_removed_edge.remove_edge(line_u,line_v)
is_v_conn_to_src = False
for src_node_feeder in source_nodes_for_feeder:
if src_node_feeder in g_temp_removed_edge and nx.has_path(g_temp_removed_edge, src_node_feeder, line_v):
is_v_conn_to_src = True
break
if not is_v_conn_to_src:
component_of_v = nx.node_connected_component(g_temp_removed_edge, line_v)
else: # v is still connected, (u,v) is likely a loop closing line. Flow is complex.
# For simplicity, this function will return 0 flow for loop lines if direction is ambiguous.
return [], 0, 0
elif dist_u > dist_v: # u is downstream of v (swap them)
# similar logic for u
g_temp_removed_edge = G.copy()
g_temp_removed_edge.remove_edge(line_u,line_v)
is_u_conn_to_src = False
for src_node_feeder in source_nodes_for_feeder:
if src_node_feeder in g_temp_removed_edge and nx.has_path(g_temp_removed_edge, src_node_feeder, line_u):
is_u_conn_to_src = True
break
if not is_u_conn_to_src:
component_of_v = nx.node_connected_component(g_temp_removed_edge, line_u) # component_of_v is actually comp of u
else:
return [], 0, 0
else: # Equidistant or complex, cannot determine simple downstream for this line
return [], 0, 0
downstream_nodes = component_of_v
total_downstream_load_kw = sum(G.nodes[n]['load_kw'] for n in downstream_nodes)
total_downstream_dg_kw = sum(G.nodes[n]['dg_kw'] for n in downstream_nodes if G.nodes[n]['dg_kw'] > 0)
return list(downstream_nodes), total_downstream_load_kw, total_downstream_dg_kw
def calculate_power_flows_and_currents(self, current_graph_state, active_dgs_kw):
"""
Simplified power flow for radial networks or parts of networks.
Returns dict of line flows and currents, and substation powers.
Flows: {(u,v): power_kw} where power_kw > 0 means u to v.
Currents: {(u,v): current_A}
Substation_powers: {'CB1': power_kw_drawn_from_substation}
"""
line_flows_kw = {}
line_currents_a = {}
substation_powers_kw = {}
# Update DG outputs in the graph state
for node_id, dg_val in active_dgs_kw.items():
if node_id in current_graph_state:
current_graph_state.nodes[node_id]['dg_kw'] = dg_val
for node_id in current_graph_state.nodes(): # Reset others if not in active_dgs_kw
if node_id not in active_dgs_kw and 'dg_kw' in current_graph_state.nodes[node_id]:
if current_graph_state.nodes[node_id].get('type') != 'source': # Don't zero out if it was never a DG
current_graph_state.nodes[node_id]['dg_kw'] = 0
# Determine connected components and their sources
# This is a very simplified load flow. It assumes power flows from sources (substations)
# down to loads. DG power reduces the load seen by upstream sections.
# It does not handle loops well without iterative methods (e.g. Hardy Cross or Newton-Raphson).
# For each feeder, calculate flows assuming radial structure
# This is an approximation. A full AC or DC power flow is more accurate.
# Initialize all line flows to 0
for u, v in current_graph_state.edges():
line_flows_kw[(u,v)] = 0
line_flows_kw[(v,u)] = 0 # For bi-directional calculation needs
line_currents_a[(u,v)] = 0
# Iterate multiple times for flow distribution in case of ties or complex paths
# This is a placeholder for a proper iterative flow solution.
# For now, a topological sort based flow for radial parts.
processed_nodes_for_flow_calc = set()
for cb_id, (src_node, start_node) in self.substation_connections.items():
if not current_graph_state.has_node(start_node) or not nx.is_connected(current_graph_state.subgraph([n for n in current_graph_state.nodes() if n != SOURCE_NODE])):
# This feeder might be entirely down if start_node is disconnected from actual nodes
substation_powers_kw[cb_id] = 0
continue
# Get nodes for this feeder (dynamic based on current_graph_state)
# Nodes connected to start_node, excluding SOURCE_NODE, if start_node is connected to SOURCE_NODE
feeder_nodes_component = set()
if current_graph_state.has_edge(src_node, start_node):
temp_g_for_feeder = current_graph_state.copy()
# Remove other substation links to isolate this feeder's component
other_links_to_remove = []
for other_cb, (s,d) in self.substation_connections.items():
if other_cb != cb_id and temp_g_for_feeder.has_edge(s,d):
other_links_to_remove.append((s,d))
temp_g_for_feeder.remove_edges_from(other_links_to_remove)
if nx.has_path(temp_g_for_feeder, src_node, start_node):
try:
# Consider only the part of the graph reachable from start_node, not crossing back to SOURCE_NODE
# except via the designated start_node path.
search_nodes = [n for n in temp_g_for_feeder.nodes if n != src_node]
sub_graph_feeder = temp_g_for_feeder.subgraph(search_nodes)
if start_node in sub_graph_feeder:
feeder_nodes_component = nx.node_connected_component(sub_graph_feeder, start_node)
except nx.NetworkXError: # if start_node not in subgraph
feeder_nodes_component = set()
if not feeder_nodes_component:
substation_powers_kw[cb_id] = 0
continue
# Order nodes from furthest to closest to substation for power accumulation
# This is for radial feeders. If loops exist, this is not sufficient.
# Using BFS layers from start_node
# Net load at each node (Load - DG)
node_net_power_kw = {}
for node in feeder_nodes_component:
node_net_power_kw[node] = current_graph_state.nodes[node]['load_kw'] - current_graph_state.nodes[node]['dg_kw']
# Accumulate power up towards the substation
# This requires a tree traversal (e.g., DFS post-order traversal) from leaves to root.
# For simplicity, if the feeder is a tree rooted at start_node:
if nx.is_tree(current_graph_state.subgraph(feeder_nodes_component | {start_node})): # Check if it's a tree
# Create a directed tree towards the source for easier traversal
# This part is complex if graph is not a tree.
# For now, sum all net loads on the feeder as the substation power (approximation)
total_feeder_net_load = sum(node_net_power_kw[n] for n in feeder_nodes_component)
substation_powers_kw[cb_id] = total_feeder_net_load
# Distribute this flow down the lines (highly simplified)
# A proper method: for each line, sum net_power of all nodes in subtree rooted by that line.
# This simplified flow calculation is a major placeholder.
# For overload risk, we need per-line flows.
# Simplified: assume current_graph_state is a tree rooted at start_node for this feeder
# Use BFS to assign flow from start_node downwards
# This is not a full power flow, but an estimation for line loading.
# Build a directed graph for this feeder based on BFS from start_node
# This is just for flow assignment direction.
# Actual flow needs to sum up demands from downstream.
# For each edge (u,v) in the feeder:
# Determine parent (closer to start_node) and child
# Power on (parent, child) = sum of net loads in subtree rooted at child.
# This is complex to implement robustly here without a full flow algorithm.
# Fallback: Use the _get_downstream_info logic if possible, iterate edges
# This needs to be called carefully to avoid double counting or misdirection in non-radial.
# For now, this function will primarily return substation_powers_kw and
# leave detailed line_flows_kw and line_currents_a for a more robust implementation
# or accept its high level of approximation.
# Let's try a slightly better approximation for line flows on a tree:
# For each edge (u,v) in the feeder tree (rooted at start_node)
# Assume u is parent of v. Power(u,v) = sum of net_power for all nodes in subtree of v.
if start_node in feeder_nodes_component: # Should be
try:
dfs_edges = list(nx.dfs_edges(nx.bfs_tree(current_graph_state.subgraph(feeder_nodes_component), start_node), source=start_node))
# Calculate power for each node including its children's power
node_total_subtree_power = node_net_power_kw.copy()
for u, v in reversed(dfs_edges): # From leaves up to root
if u in node_total_subtree_power and v in node_total_subtree_power:
node_total_subtree_power[u] += node_total_subtree_power[v]
for u,v in dfs_edges: # From root down
if v in node_total_subtree_power:
flow = node_total_subtree_power[v]
line_flows_kw[(u,v)] = flow
line_flows_kw[(v,u)] = -flow # Convention for direction
line_currents_a[(u,v)] = self._calculate_line_current_kw(flow)
except Exception as e:
# print(f"Warning: Could not perform tree-based flow for {cb_id} due to {e}")
pass # Keep substation power as sum, line flows might be inaccurate
else: # Not a tree, flow calculation is more complex.
# print(f"Warning: Feeder {cb_id} is not a tree. Simplified flow may be inaccurate.")
total_feeder_net_load = sum(node_net_power_kw[n] for n in feeder_nodes_component if n in node_net_power_kw)
substation_powers_kw[cb_id] = total_feeder_net_load
# Line flows in meshed networks require iterative solvers.
# For now, this part will be very approximate for meshed sections.
# Handle reverse power flow and inter-feeder DG adjustment
# "分布式能源不得向上级电网倒送功率"
# "可以在相邻馈线间进行调节"
for cb_id in list(substation_powers_kw.keys()):
if substation_powers_kw[cb_id] < 0: # Reverse power flow
excess_dg_on_feeder = -substation_powers_kw[cb_id]
# Try to transfer to other feeders via inter-feeder tie lines
# This logic is complex and needs careful state management.
# For Q1, a simpler approach might be DG curtailment on that feeder.
# Find inter-feeder tie switches connected to this feeder
# Example: S62-3 connects Feeder 3 (node 62) to Feeder 1 (node 19)
# If Feeder 1 has excess_dg, it might try to send to Feeder 3 via (19,62)
# This is an advanced feature. For now, assume DG curtailment if倒送.
# To implement curtailment: identify DGs on this feeder, reduce their output
# proportionally until substation_powers_kw[cb_id] >= 0.
# This would require re-calculating flows.
# For now, just flag it.
# print(f"Warning: Reverse power flow on {cb_id} of {substation_powers_kw[cb_id]} kW. DG curtailment or transfer needed.")
# A simple curtailment:
dG_on_feeder_nodes = [n for n in self.feeder_info.get(cb_id, set()) if n in active_dgs_kw and active_dgs_kw[n] > 0]
total_dg_cap_on_feeder = sum(active_dgs_kw[n] for n in dG_on_feeder_nodes)
if total_dg_cap_on_feeder > 0:
curtail_ratio = min(1.0, excess_dg_on_feeder / total_dg_cap_on_feeder) if total_dg_cap_on_feeder >0 else 0
for dg_node in dG_on_feeder_nodes:
active_dgs_kw[dg_node] *= (1 - curtail_ratio)
# Flows need to be recalculated after curtailment. This suggests an iterative solution.
# For this submission, we'll assume this check is done *before* final flow calc,
# or simply note the violation.
# To avoid recursion here, this function should ideally take DGs as fixed input.
# The adjustment logic should be outside or iterative.
pass
return line_flows_kw, line_currents_a, substation_powers_kw
def calculate_overload_risk(self):
"""
Calculates overload risk for the current DG setup.
Assumes DGs are at their BASE_DG_CAPACITY_KW.
"""
# Get current operational graph (no faults, ties normally open)
g_op = self._get_subgraph_with_operational_lines(self.graph)
# Calculate power flows and currents
# Need to handle DG outputs properly.
current_dg_outputs = self.dg_kw.copy() # Use the model's current DG settings
# Iterative step for DG curtailment if reverse power flow:
# This is a simplified loop. A more robust solution uses optimization or better heuristics.
for _iter in range(3): # Max 3 iterations for adjustment
line_flows, line_currents, substation_powers = self.calculate_power_flows_and_currents(g_op, current_dg_outputs)
reverse_power_detected = False
for cb_id, power_kw in substation_powers.items():
if power_kw < -1e-3: # Small threshold for倒送
reverse_power_detected = True
# print(f"Info: Reverse power on {cb_id} ({power_kw:.2f} kW). Attempting curtailment.")
excess_dg_on_feeder = abs(power_kw)
feeder_nodes_for_cb = self.feeder_info.get(cb_id, set())
dG_on_feeder_nodes = [n for n in feeder_nodes_for_cb if n in current_dg_outputs and current_dg_outputs[n] > 0]
total_dg_cap_on_feeder = sum(current_dg_outputs[n] for n in dG_on_feeder_nodes)
if total_dg_cap_on_feeder > 1e-3 : # Avoid division by zero
curtail_amount_total = excess_dg_on_feeder
for dg_node in dG_on_feeder_nodes:
# Proportional curtailment
proportion = current_dg_outputs[dg_node] / total_dg_cap_on_feeder
curtail_this_dg = proportion * curtail_amount_total
current_dg_outputs[dg_node] = max(0, current_dg_outputs[dg_node] - curtail_this_dg)
else: # No DG to curtail, reverse power might be from other sources or model issue
pass
if not reverse_power_detected:
break
# Final flows after potential curtailment
line_flows, line_currents, substation_powers = self.calculate_power_flows_and_currents(g_op, current_dg_outputs)
overloaded_lines = []
for u, v, data in g_op.edges(data=True):
if data.get('type') == 'substation_link': continue # Don't check substation links themselves for overload here
current = line_currents.get((u,v), 0)
# If flow is from v to u, current might be stored as current_uv = -current_vu
# Take absolute value of flow for current calculation, or ensure current is always positive.
# The _calculate_line_current_kw uses abs(power_kw) so current should be positive.
rated_current = data.get('rated_current_a', DEFAULT_LINE_RATED_CURRENT_A)
if current > 1.1 * rated_current:
overloaded_lines.append({'edge': (u,v), 'current': current, 'rated': rated_current, 'over_by_%': (current/(1.1*rated_current)-1)*100 if rated_current>0 else float('inf')})
if overloaded_lines:
# print(f"System Overload Detected. Overloaded lines: {overloaded_lines}")
# P_over = 1 (for this deterministic scenario)
# C_over = fixed penalty or sum of penalties
risk_overload = 1.0 * COST_PENALTY_FOR_ANY_OVERLOAD
# Or, sum of consequences for each overloaded line, if C_over is per line.
# risk_overload = sum(some_cost_function(ol['over_by_%']) for ol in overloaded_lines)
else:
# print("System is NOT overloaded in the base case.")
risk_overload = 0.0
return risk_overload, overloaded_lines, substation_powers, current_dg_outputs
def calculate_load_loss_risk(self):
"""
Calculates total load loss risk by considering single line faults.
R_loss = sum(P_fault_i * C_loss_i)
P_fault_i = annual probability of fault i
C_loss_i = consequence of fault i (e.g., unserved_load_kw * COST_VOLL_PER_KW)
"""
total_load_loss_risk = 0.0
detailed_fault_impacts = []
# Iterate through all operational lines (excluding substation virtual links for fault simulation)
original_edges = [ (u,v,data) for u,v,data in self.graph.edges(data=True)
if data.get('type') != 'substation_link' and not data.get('is_tie', False)]
for u_fault, v_fault, line_data_faulted in original_edges:
faulty_edge = (u_fault, v_fault)
line_length_km = line_data_faulted.get('length_km', 0)
# Probability of this specific line failing (annual)
# Assuming failure rates are independent and this is the probability of this line being the one to fail.
prob_line_fault = line_length_km * FAILURE_RATE_LINE_PER_KM
if prob_line_fault == 0: continue
# --- Simulate fault ---
g_faulted = self.graph.copy()
if not g_faulted.has_edge(*faulty_edge): continue
g_faulted.edges[faulty_edge]['failed'] = True # Mark as failed
g_after_fault_isolation = self._get_subgraph_with_operational_lines(g_faulted, faulty_edge=faulty_edge)
# --- Identify initial load loss ---
initial_unserved_load_kw = 0
disconnected_load_nodes = {} # {node: load_kw}
# Check connectivity for all load nodes
for node_id, load_kw in self.loads_kw.items():
if load_kw <= 0: continue
is_connected_to_source = False
for cb_id, (src_node, start_node) in self.substation_connections.items():
if nx.has_path(g_after_fault_isolation, src_node, node_id):
is_connected_to_source = True
break
if not is_connected_to_source:
initial_unserved_load_kw += load_kw
disconnected_load_nodes[node_id] = load_kw
if initial_unserved_load_kw == 0: # Fault does not cause load loss (e.g. redundant line)
detailed_fault_impacts.append({'fault': faulty_edge, 'unserved_kw_initial': 0, 'unserved_kw_final':0, 'restored_kw':0, 'risk_contrib':0})
continue
# --- Attempt restoration via tie lines ---
# This is a complex part. Needs to:
# 1. Identify disconnected areas and loads.
# 2. Identify available tie switches that can connect these areas to healthy feeders.
# 3. Check capacity of tie lines and the supporting feeder.
# 4. Prioritize restoration (e.g., maximize load restored).
# For this model, a simplified restoration:
# Iterate over available tie switches. If closing one helps, simulate it.
# This should be greedy or more optimized.
restored_load_kw_total_for_this_fault = 0
# Create a graph state for restoration attempts
g_for_restoration = g_after_fault_isolation.copy()
# Sort disconnected loads by size (optional, for prioritization)
sorted_disconnected_loads = sorted(disconnected_load_nodes.items(), key=lambda item: item[1], reverse=True)
# Try closing tie switches one by one (if they connect a live part to a dead part)
# This is a very simplified greedy approach.
# A proper approach would evaluate all combinations or use optimization.
# Identify current live sources/feeders
live_feeder_sources = [] # (cb_id, start_node_of_live_feeder)
for cb_id, (src,start) in self.substation_connections.items():
if nx.has_path(g_for_restoration, src, start): # Check if substation itself is connected
live_feeder_sources.append(start)
for tie in self.tie_switches_info:
tie_n1, tie_n2 = tie['nodes']
tie_capacity_a = tie['capacity_A']
if not g_for_restoration.has_node(tie_n1) or not g_for_restoration.has_node(tie_n2): continue
if g_for_restoration.has_edge(tie_n1, tie_n2): continue # Already closed or part of main graph (should not be for ties)
# Check if one end is live and other is dead (or part of the disconnected component)
tie_n1_is_live = any(nx.has_path(g_for_restoration, src, tie_n1) for src in live_feeder_sources)
tie_n2_is_live = any(nx.has_path(g_for_restoration, src, tie_n2) for src in live_feeder_sources)
if tie_n1_is_live == tie_n2_is_live: continue # Both live or both dead, closing doesn't bridge outage for now
live_tie_node, dead_tie_node = (tie_n1, tie_n2) if tie_n1_is_live else (tie_n2, tie_n1)
# Check if dead_tie_node is part of the current outage we are trying to fix
# This requires knowing which component dead_tie_node belongs to.
# For now, assume if it's not live, it's part of some outage.
# Simulate closing this tie switch
g_for_restoration.add_edge(live_tie_node, dead_tie_node, id=tie['id'], type='tie_closed',
rated_current_a=tie_capacity_a, resistance_ohm=0.001, length_km=0.01)
# Check how much load can be restored through this tie without overloading tie or new path
# This requires a flow calculation on g_for_restoration.
# Simplified: Check loads now connected.
newly_restored_load_kw_this_tie = 0
temp_restored_nodes_this_tie = []
for node_id, load_val in disconnected_load_nodes.items():
if node_id not in g_for_restoration: continue # Should not happen
# Check if this node is now connected to ANY source
is_now_connected = any(nx.has_path(g_for_restoration, src, node_id) for src in live_feeder_sources)
if is_now_connected and node_id not in temp_restored_nodes_this_tie: # And not already counted as restored by previous ties
# More checks needed:
# 1. Tie line capacity: Power through (live_tie_node, dead_tie_node) <= tie_capacity_a
# 2. Path capacity on the live feeder.
# This is where the simplified flow becomes a bottleneck.
# For now, assume if connected, it can be restored up to a certain limit.
# This is a MAJOR simplification.
newly_restored_load_kw_this_tie += load_val
temp_restored_nodes_this_tie.append(node_id)
# Here, we'd need to check if adding newly_restored_load_kw_this_tie overloads the tie or feeder.
# If current_through_tie > tie_capacity_a, then not all of this load can be restored.
# This part needs a proper constrained flow allocation.
# For now, let's assume a fraction can be restored if connected, or all if small.
# This is a placeholder for a more robust restoration algorithm.
# Let's assume, for now, if connected, it's restored. If this overloads things,
# the overload risk model should capture it (but that's for normal state).
# Here, the goal is to minimize unserved load.
# If this tie leads to overload, we shouldn't use it or only partially.
# For now, naively accept all newly connected load.
if newly_restored_load_kw_this_tie > 0:
restored_load_kw_total_for_this_fault += newly_restored_load_kw_this_tie
# Update disconnected_load_nodes:
for r_node in temp_restored_nodes_this_tie:
if r_node in disconnected_load_nodes:
del disconnected_load_nodes[r_node] # No longer disconnected
else: # Closing this tie didn't help, revert
if g_for_restoration.has_edge(live_tie_node, dead_tie_node):
g_for_restoration.remove_edge(live_tie_node, dead_tie_node)
# Final unserved load for this fault scenario
final_unserved_load_kw = initial_unserved_load_kw - restored_load_kw_total_for_this_fault
final_unserved_load_kw = max(0, final_unserved_load_kw) # Cannot be negative
consequence_c_loss = final_unserved_load_kw * COST_VOLL_PER_KW
risk_contribution = prob_line_fault * consequence_c_loss
total_load_loss_risk += risk_contribution
detailed_fault_impacts.append({
'fault_type': 'line', 'component_id': faulty_edge, 'prob_fault': prob_line_fault,
'unserved_kw_initial': initial_unserved_load_kw,
'restored_kw': restored_load_kw_total_for_this_fault,
'unserved_kw_final': final_unserved_load_kw,
'consequence_c_loss': consequence_c_loss,
'risk_contribution': risk_contribution
})
# TODO: Add DG faults, Switch faults, User faults if they cause wider outages.
# For DG faults: prob_dg_fault = FAILURE_RATE_DG_PERCENT
# A DG fault primarily impacts system's ability to meet load or avoid overload.
# It doesn't directly cause load loss unless it's islanded and the DG is the only source.
# The problem implies grid-connected DGs.
return total_load_loss_risk, detailed_fault_impacts
# --- Main Execution ---
if __name__ == '__main__':
print("--- 配电网风险评估模型 Q1 ---")
# 1. Load Data
print("\n1. 加载数据...")
loads = load_load_data()
topology = load_topology_data()
# print(f"负荷数据: {len(loads)} 点")
# print(f"拓扑数据: {len(topology)} 条线路")
# 2. Initialize Power Grid Model
print("\n2. 初始化电网模型...")
grid = PowerGridModel(loads, topology, DG_LOCATIONS_KW, TIE_SWITCHES_INFO, SUBSTATION_CONNECTIONS)
# print(f"电网图: {grid.graph.number_of_nodes()} 个节点, {grid.graph.number_of_edges()} 条边")
# print(f"馈线信息: {grid.feeder_info}")
# --- 问题1: 失负荷风险和过负荷风险计算模型 ---
print("\n--- 问题1: 风险计算 ---")
# A. 过负荷风险模型 (R_over = P_over * C_over)
# For Q1, DGs are at BASE_DG_CAPACITY_KW. This is a deterministic check for this state.
# P_over = 1 if overload occurs, 0 otherwise. C_over is the penalty.
print("\nA. 计算过负荷风险...")
# Note: The calculate_power_flows_and_currents is highly simplified.
# Results for overload depend heavily on its accuracy and line ratings.
try:
risk_overload, overloaded_lines_details, substation_p, final_dg_out = grid.calculate_overload_risk()
print(f" 计算得到的过负荷风险 (R_over): ${risk_overload:.2f}")
if overloaded_lines_details:
print(f" 检测到过负荷线路 ({len(overloaded_lines_details)} 条):")
# for ol in overloaded_lines_details[:3]: # Print first 3
# print(f" - 线路 {ol['edge']}, 电流: {ol['current']:.2f}A, 额定: {ol['rated']:.2f}A, 超出: {ol['over_by_%']:.2f}%")
else:
print(" 在当前DG配置下，未检测到线路过负荷。")
# print(f" 变电站出口功率 (kW): {substation_p}")
# print(f" 最终DG出力 (kW) (可能经过削减): {final_dg_out}")
except Exception as e:
print(f" 计算过负荷风险时发生错误: {e}")
risk_overload = -1 # Indicate error
# B. 失负荷风险模型 (R_loss = sum(P_fault_i * C_loss_i))
print("\nB. 计算失负荷风险...")
# Note: Restoration logic is simplified.
try:
total_r_loss, fault_details = grid.calculate_load_loss_risk()
print(f" 计算得到的总失负荷风险 (R_loss): ${total_r_loss:.2f} (基于所选成本)")
# print("\n 部分故障场景详情:")
# for fd in fault_details[:3]: # Print first 3
# print(f" - 故障线路: {fd.get('component_id')}, "
# f"初始失负荷: {fd.get('unserved_kw_initial'):.2f} kW, "
# f"最终失负荷: {fd.get('unserved_kw_final'):.2f} kW, "
# f"风险贡献: ${fd.get('risk_contribution'):.2f}")
except Exception as e:
print(f" 计算失负荷风险时发生错误: {e}")
total_r_loss = -1 # Indicate error
print("\n--- 模型执行完毕 ---")
print("注意: 此模型包含多项简化和假设 (如线路额定电流, 成本参数, 潮流计算简化, 恢复逻辑简化).")
print("结果的准确性取决于这些假设的合理性和参数的精确性。")
根据此代码分布式能源接入配电网的风险分析
背景知识：
随着我国双碳目标的推进，可再生分布式能源在配电网中的大规模应用不可避免，这对传统配电网运行提出挑战。为了量化分析配电网中接入分布式能源的风险，需要对其进行建模与分析。
配电网发生故障后失负荷，可以通过联络线实现部分复电，供电恢复的目标是在系统拓扑结构发生变化时，将系统的经济损失降至最小，行业分类将消费者分为居民住宅、商业、政府机构和办公建筑等类别，按照这种分类，供电中断危害可依据部门客户危害度函数进行计算。
分布式能源接入配电网，由于其发电出力的波动性与不确定性，对馈线的失负荷和过负荷带来影响，为了提高配电网就地消纳能力，本竞赛题目要求分布式能源不得向上级电网倒送功率。
风险评估的通用计算公式为：
（1）
式中：代表系统风险；代表系统失负荷的发生概率；代表由系统失负荷造成的危害程度；代表系统过负荷的发生概率；代表由系统过负荷造成的危害程度。
由式（1）可知，风险评估的量化计算为该场景下各事件概率与其危害的乘积之和，为了实现风险评估，需要从事件概率计算及事件造成后果的危害度函数构建两方面入手进行分析与建模。
名称解释
失负荷：因故障导致负荷供电中断 。
过负荷：线路电流超过额定载流量10%以上。
馈线：从变电站出线开关到终端负荷的配电线路。
联络线：通过联络开关连接不同馈线的线路，正常运行时联络开关处于断开状态，根据运行方式调整的需要可以调整联络开关的状态，实现馈线间的功率转移。
配电网的运行方式: 在同一时刻, 馈线上的每个负荷与各变电站之间只有一条通路（只由一个变电站出线开关CB供电），每个分布式电源DG都可接入馈线供电，不同用户类型的停电损失会造成不同的危害度。
计算配电系统风险的约束条件如下：
（1）各个类型故障是独立发生的，同一时间同一类型只发生一个故障；
（2）不考虑无功功率和电压越限的影响，风险计算分析仅考虑有功功率和电流的影响；
（3）联络开关不考虑故障恢复的自愈系统对失负荷的影响，但是需考虑联络开关的负荷转移能力。
问题：
（1）请分别建立分布式能源接入配电网后，配电系统的失负荷风险和过负荷风险的计算模型，要求：失负荷风险模型要记及本馈线故障造成的负荷损失可以从其他相邻馈线通过联络线转供实现复电的情况；过负荷风险模型要记及本馈线的有功功率不得向上级变电站倒送的情况，但是可以在相邻馈线间进行调节；
最新发布
本文介绍了QPainterPath类的使用方法，包括如何构建复杂的图形形状、获取路径信息及路径转换等。此外，还提供了示例代码和相关函数介绍。
扫一扫
1806
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
