940A. Points on the line(思维逆向)

本文探讨了一道关于在坐标轴上删除部分点以确保剩余点集直径不超过特定值的问题。介绍了解题思路,并提供了一个使用C++实现的具体算法示例。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

We've got no test cases. A big olympiad is coming up. But the problemsetters' number one priority should be adding another problem to the round.

The diameter of a multiset of points on the line is the largest distance between two points from this set. For example, the diameter of the multiset {1, 3, 2, 1} is 2.

Diameter of multiset consisting of one point is 0.

You are given n points on the line. What is the minimum number of points you have to remove, so that the diameter of the multiset of the remaining points will not exceed d?

Input

The first line contains two integers n and d (1 ≤ n ≤ 100, 0 ≤ d ≤ 100) — the amount of points and the maximum allowed diameter respectively.

The second line contains n space separated integers (1 ≤ xi ≤ 100) — the coordinates of the points.

Output

Output a single integer — the minimum number of points you have to remove.

Examples
Input
Copy
3 1
2 1 4
Output
1
Input
Copy
3 0
7 7 7
Output
0
Input
Copy
6 3
1 3 4 6 9 10
Output
3
Note

In the first test case the optimal strategy is to remove the point with coordinate 4. The remaining points will have coordinates 1 and 2, so the diameter will be equal to 2 - 1 = 1.

In the second test case the diameter is equal to 0, so its is unnecessary to remove any points.

In the third test case the optimal strategy is to remove points with coordinates 1, 9 and 10. The remaining points will have coordinates 3, 4 and 6, so the diameter will be equal to 6 - 3 = 3.

题意:问最少删除几个数字,使得剩余数字的 最大值   -   最小值 < = m;

思路:思维逆向,求符合题意的最多数字

 #include<bits/stdc++.h>
using namespace std;
const int N = 105;
int arr[N];
int main(){
	int n,m;
	scanf("%d%d",&n,&m);
	for(int i=1;i<=n;i++)
	   scanf("%d",&arr[i]);
	sort(arr+1,arr+n+1);
	int ans=0;
	for(int i=1;i<=n;i++){
		for(int j=i;j<=n;j++)
		if(arr[j]-arr[i]<=m)
		   ans=max(ans,j-i+1);
	}
	printf("%d\n",n-ans);
	return 0;
}

# -*- coding: utf-8 -*- """ Traffic Monitoring System v2.5 (Camera Version) Features: - Real-time camera object detection - YOLOv8 object detection - ByteTrack tracking - Matching colors for boxes and trajectories - Huawei Cloud IoT data reporting - Movable table - Trajectory display """ import os import sys import cv2 import numpy as np import torch import supervision as sv from ultralytics import YOLO from typing import Dict, Optional, Tuple, List import warnings import time import logging import socket import json from pathlib import Path from collections import defaultdict, deque # Import Huawei IoT SDK with encoding handling try: from iot_device_sdk_python.client.client_conf import ClientConf from iot_device_sdk_python.client.connect_auth_info import ConnectAuthInfo from iot_device_sdk_python.iot_device import IotDevice from iot_device_sdk_python.client.request.service_property import ServiceProperty except ImportError as e: logging.error(f"Huawei IoT SDK import error: {str(e)}") sys.exit(1) # Configure global logging with UTF-8 encoding logging.basicConfig( level=logging.INFO, format='%(asctime)s - %(name)s - %(levelname)s - %(message)s', handlers=[ logging.StreamHandler(), logging.FileHandler('traffic_monitor.log', encoding='utf-8', mode='a') ] ) logger = logging.getLogger("TrafficMonitor") class Config: """System configuration parameters""" MODEL_PATH = r"C:\Users\yuanb\Desktop\物联网设计大赛无pyqt同时MQTT传输\runs\detect\train_v8\weights\best.pt" CLASS_NAMES = { 0: 'Pedestrian', 1: 'Person', 2: 'Bicycle', 3: 'Car', 4: 'Van', 5: 'Truck', 6: 'Tricycle', 7: 'CoveredTricycle', 8: 'Bus', 9: 'Motorcycle' } # Color palette for different classes CLASS_COLORS = { 0: (255, 0, 0), # Red - Pedestrian 1: (0, 255, 0), # Green - Person 2: (0, 0, 255), # Blue - Bicycle 3: (255, 255, 0), # Yellow - Car 4: (255, 0, 255), # Magenta - Van 5: (0, 255, 255), # Cyan - Truck 6: (255, 128, 0), # Orange - Tricycle 7: (128, 0, 255), # Purple - CoveredTricycle 8: (0, 128, 255), # Light Blue - Bus 9: (255, 0, 128) # Pink - Motorcycle } IOT_CONFIG = { "server_uri": "117.78.5.125", "port": 1883, "device_id": "681eeebc9314d118511a5690_DEV4", "secret": "Wjy15640916538", "service_id": "鲁班猫4", "use_tls": False, "keep_alive": 300, "clean_session": False } REPORT_INTERVAL = 10 MIN_REPORT_INTERVAL = 30 LINE_COLOR = sv.Color(r=0, g=255, b=0) TEXT_COLOR = (255, 255, 255) # White text for better visibility FONT_SCALE = 0.7 TABLE_POSITION = (10, 400) # Initial position, will be updated by mouse TABLE_ROW_HEIGHT = 30 TABLE_COL_WIDTH = 180 # Increased width TABLE_COLORS = { 'header': (0, 100, 255), 'row_even': (50, 50, 50), 'row_odd': (80, 80, 80), 'text': (255, 255, 255) } TRACK_HISTORY_LENGTH = 30 # Number of trajectory points to keep class LineCounter: """Improved counting logic with bidirectional counting and classification""" def __init__(self, start_point: sv.Point, end_point: sv.Point): self.start = start_point self.end = end_point self.tracker_state = {} self.counts = { 'Inbound_Vehicle': 0, 'Outbound_Vehicle': 0, 'Inbound_Pedestrian': 0, 'Outbound_Pedestrian': 0, 'Total_Vehicles': 0, 'Total_Pedestrians': 0 } self.class_mapping = { 'Vehicle': {3, 4, 5, 6, 7, 8, 9}, 'Pedestrian': {0, 1} } self.history = [] def _get_side(self, point: sv.Point) -> bool: """Calculate which side of the line the point is on (cross product method)""" vector_line = (self.end.x - self.start.x, self.end.y - self.start.y) vector_point = (point.x - self.start.x, point.y - self.start.y) cross_product = vector_line[0] * vector_point[1] - vector_line[1] * vector_point[0] return cross_product > 0 def update(self, detections: sv.Detections): """Update counter state""" current_counts = { 'Inbound_Vehicle': 0, 'Outbound_Vehicle': 0, 'Inbound_Pedestrian': 0, 'Outbound_Pedestrian': 0, } for xyxy, _, confidence, class_id, tracker_id in detections: if tracker_id is None or confidence < 0.5: continue # Calculate center point x1, y1, x2, y2 = xyxy center = sv.Point(x=int((x1 + x2) / 2), y=int((y1 + y2) / 2)) current_side = self._get_side(center) # Initialize new object if tracker_id not in self.tracker_state: self.tracker_state[tracker_id] = current_side continue # Check if crossed the line if self.tracker_state[tracker_id] != current_side: if class_id in self.class_mapping['Vehicle']: direction = 'Inbound' if current_side else 'Outbound' current_counts[f'{direction}_Vehicle'] += 1 self.counts[f'{direction}_Vehicle'] += 1 self.counts['Total_Vehicles'] += 1 elif class_id in self.class_mapping['Pedestrian']: direction = 'Inbound' if current_side else 'Outbound' current_counts[f'{direction}_Pedestrian'] += 1 self.counts[f'{direction}_Pedestrian'] += 1 self.counts['Total_Pedestrians'] += 1 self.tracker_state[tracker_id] = current_side if any(current_counts.values()): self.history.append({ 'timestamp': time.time(), 'counts': current_counts.copy() }) def get_counts(self) -> Dict[str, int]: """Get current statistics""" return self.counts.copy() def get_recent_counts(self, seconds=10) -> Dict[str, int]: """Get counts from recent specified time period""" recent_counts = { 'Inbound_Vehicle': 0, 'Outbound_Vehicle': 0, 'Inbound_Pedestrian': 0, 'Outbound_Pedestrian': 0 } cutoff = time.time() - seconds for record in reversed(self.history): if record['timestamp'] < cutoff: break for key in recent_counts: recent_counts[key] += record['counts'].get(key, 0) return recent_counts class HuaweiIoTClient: """Huawei Cloud IoT Client""" def __init__(self, config: dict): self.config = config self.device = None self.max_retries = 5 self.base_reconnect_delay = 5 self.last_report_time = 0 self._initialize_client() def _initialize_client(self): """Initialize client connection""" try: auth_info = ConnectAuthInfo() auth_info.server_uri = self.config["server_uri"] auth_info.port = self.config["port"] auth_info.id = self.config["device_id"] auth_info.secret = self.config["secret"] self.config["use_tls"] = False client_conf = ClientConf(auth_info) client_conf.clean_session = self.config["clean_session"] client_conf.keep_alive = self.config["keep_alive"] client_conf.will_topic = f"$oc/devices/{self.config['device_id']}/sys/events/down" client_conf.will_message = json.dumps({ "type": "ABNORMAL_DISCONNECT", "timestamp": time.strftime("%Y%m%dT%H%M%SZ", time.gmtime()) }) self.device = IotDevice(client_conf) self._connect_with_retry() except Exception as e: logger.error(f"Client initialization failed: {str(e)}") raise def _check_network(self) -> bool: """Check network connectivity""" try: with socket.create_connection( (self.config["server_uri"], self.config["port"]), timeout=10 ): return True except Exception as e: logger.warning(f"Network check failed: {str(e)}") return False def _connect_with_retry(self): """Exponential backoff reconnection mechanism""" retry_count = 0 while retry_count < self.max_retries: try: if not self._check_network(): raise ConnectionError("Network unavailable") result = self.device.connect() if result == 0: logger.info(f"Connection successful (attempt {retry_count + 1}/{self.max_retries})") return True wait_time = min(self.base_reconnect_delay * (2 ** retry_count), 60) logger.warning(f"Connection failed, retrying in {wait_time} seconds... (error code: {result})") time.sleep(wait_time) retry_count += 1 except Exception as e: wait_time = min(self.base_reconnect_delay * (2 ** retry_count), 60) logger.error(f"Connection exception: {str(e)}, retrying in {wait_time} seconds") time.sleep(wait_time) retry_count += 1 logger.error("Max retries reached, giving up") return False def is_connected(self) -> bool: """Check connection status""" try: return (self.device is not None and hasattr(self.device, 'get_client') and self.device.get_client().is_connected()) except Exception as e: logger.warning(f"Connection status check exception: {str(e)}") return False def report_traffic_data(self, counts: Dict[str, int]) -> bool: """Data reporting method""" current_time = time.time() if current_time - self.last_report_time < Config.MIN_REPORT_INTERVAL: return False if not self.is_connected(): logger.warning("Connection lost, attempting to reconnect...") if not self._connect_with_retry(): return False try: service = ServiceProperty() service.service_id = self.config["service_id"] service.properties = { "vehicle_in": counts.get('Inbound_Vehicle', 0), "vehicle_out": counts.get('Outbound_Vehicle', 0), "vehicle_total": counts.get('Total_Vehicles', 0), "person_in": counts.get('Inbound_Pedestrian', 0), "person_out": counts.get('Outbound_Pedestrian', 0), "person_total": counts.get('Total_Pedestrians', 0), "timestamp": int(time.time() * 1000) } service.event_time = time.strftime("%Y%m%dT%H%M%SZ", time.gmtime()) result = self.device.get_client().report_properties([service]) if result == 0: self.last_report_time = current_time logger.info(f"Data reported successfully: {service.properties}") return True logger.error(f"Reporting failed, return code: {result}") return False except Exception as e: logger.error(f"Reporting exception: {str(e)}") return False class TrafficMonitor: """Main traffic monitoring program""" def __init__(self): # Hardware configuration self.device = 'cuda' if torch.cuda.is_available() else 'cpu' logger.info(f"Running on: {self.device.upper()}") # Model loading self.model = self._load_model() self.tracker = sv.ByteTrack() self.line_counter = None self.frame_count = 0 self.current_tracks = {} self.track_history = defaultdict(lambda: deque(maxlen=Config.TRACK_HISTORY_LENGTH)) # Visualization tools with custom color settings self.box_annotator = sv.BoundingBoxAnnotator(color_lookup=sv.ColorLookup.TRACK) self.label_annotator = sv.LabelAnnotator(color_lookup=sv.ColorLookup.TRACK) self.line_annotator = sv.LineZoneAnnotator( thickness=2, text_thickness=1, text_scale=0.5 ) # Table properties self.table_position = list(Config.TABLE_POSITION) self.is_dragging = False self.drag_offset = (0, 0) # Huawei Cloud connection try: self.iot_client = HuaweiIoTClient(Config.IOT_CONFIG) except Exception as e: logger.error(f"Huawei IoT client initialization failed: {str(e)}") self.iot_client = None def _load_model(self) -> YOLO: """Load YOLO model""" try: model = YOLO(Config.MODEL_PATH) model.fuse() logger.info(f"Model loaded successfully, detection classes: {model.names}") return model except Exception as e: logger.error(f"Model loading failed: {str(e)}") sys.exit(1) def set_counting_line(self, start: Tuple[int, int], end: Tuple[int, int]): """Set counting line coordinates""" self.line_counter = LineCounter( start_point=sv.Point(*start), end_point=sv.Point(*end) ) logger.info(f"Counting line set: {start} -> {end}") def _draw_stats(self, frame: np.ndarray, counts: Dict[str, int]) -> np.ndarray: """Draw statistics on frame""" stats = [ f"Inbound Vehicle: {counts['Inbound_Vehicle']}", f"Outbound Vehicle: {counts['Outbound_Vehicle']}", f"Total Vehicles: {counts['Total_Vehicles']}", f"Inbound Pedestrian: {counts['Inbound_Pedestrian']}", f"Outbound Pedestrian: {counts['Outbound_Pedestrian']}", f"Total Pedestrians: {counts['Total_Pedestrians']}" ] for i, text in enumerate(stats): cv2.putText( frame, text, (10, 30 + i * 30), cv2.FONT_HERSHEY_SIMPLEX, Config.FONT_SCALE, Config.TEXT_COLOR, 2 ) return frame def _draw_data_table(self, frame: np.ndarray, counts: Dict[str, int], recent_counts: Dict[str, int]) -> np.ndarray: """Draw movable data table""" x, y = self.table_position row_height = Config.TABLE_ROW_HEIGHT col_width = Config.TABLE_COL_WIDTH # Calculate table dimensions num_cols = 3 num_rows = 5 # Header + 4 data rows table_width = num_cols * col_width table_height = num_rows * row_height # Draw table background cv2.rectangle(frame, (x, y), (x + table_width, y + table_height), (40, 40, 40), -1) # Table header headers = ["Metric", "Total", "Last 10s"] for i, header in enumerate(headers): cv2.rectangle(frame, (x + i * col_width, y), (x + (i + 1) * col_width, y + row_height), Config.TABLE_COLORS['header'], -1) cv2.putText(frame, header, (x + i * col_width + 10, y + int(row_height * 0.7)), cv2.FONT_HERSHEY_SIMPLEX, 0.6, Config.TABLE_COLORS['text'], 1) # Data rows rows = [ ("Inbound Vehicle", counts['Inbound_Vehicle'], recent_counts['Inbound_Vehicle']), ("Outbound Vehicle", counts['Outbound_Vehicle'], recent_counts['Outbound_Vehicle']), ("Inbound Pedestrian", counts['Inbound_Pedestrian'], recent_counts['Inbound_Pedestrian']), ("Outbound Pedestrian", counts['Outbound_Pedestrian'], recent_counts['Outbound_Pedestrian']) ] for row_idx, (label, total, recent) in enumerate(rows): bg_color = Config.TABLE_COLORS['row_even'] if row_idx % 2 == 0 else Config.TABLE_COLORS['row_odd'] # Draw row background cv2.rectangle(frame, (x, y + (row_idx + 1) * row_height), (x + table_width, y + (row_idx + 2) * row_height), bg_color, -1) # Draw cell content cv2.putText(frame, label, (x + 10, y + (row_idx + 1) * row_height + int(row_height * 0.7)), cv2.FONT_HERSHEY_SIMPLEX, 0.5, Config.TABLE_COLORS['text'], 1) cv2.putText(frame, str(total), (x + col_width + 10, y + (row_idx + 1) * row_height + int(row_height * 0.7)), cv2.FONT_HERSHEY_SIMPLEX, 0.5, Config.TABLE_COLORS['text'], 1) cv2.putText(frame, str(recent), (x + 2 * col_width + 10, y + (row_idx + 1) * row_height + int(row_height * 0.7)), cv2.FONT_HERSHEY_SIMPLEX, 0.5, Config.TABLE_COLORS['text'], 1) # Draw drag handle (small rectangle at top-right) handle_size = 20 cv2.rectangle(frame, (x + table_width - handle_size, y), (x + table_width, y + handle_size), (0, 255, 255), -1) return frame def _handle_table_drag(self, event, x, y, flags, param): """Mouse callback for table dragging""" table_x, table_y = self.table_position table_width = 3 * Config.TABLE_COL_WIDTH table_height = 5 * Config.TABLE_ROW_HEIGHT handle_size = 20 # Check if mouse is in drag handle area handle_rect = (table_x + table_width - handle_size, table_y, table_x + table_width, table_y + handle_size) if event == cv2.EVENT_LBUTTONDOWN: if (handle_rect[0] <= x <= handle_rect[2] and handle_rect[1] <= y <= handle_rect[3]): self.is_dragging = True self.drag_offset = (x - table_x, y - table_y) elif event == cv2.EVENT_LBUTTONUP: self.is_dragging = False elif event == cv2.EVENT_MOUSEMOVE: if self.is_dragging: # Update table position with offset new_x = x - self.drag_offset[0] new_y = y - self.drag_offset[1] # Keep table within frame bounds frame_height, frame_width = param.shape[:2] new_x = max(0, min(new_x, frame_width - table_width)) new_y = max(0, min(new_y, frame_height - table_height)) self.table_position = [new_x, new_y] def _draw_tracks(self, frame: np.ndarray, detections: sv.Detections) -> np.ndarray: """Draw object trajectories with matching colors""" for xyxy, _, confidence, class_id, tracker_id in detections: if tracker_id is None or confidence < 0.5: continue # Get color based on class color = Config.CLASS_COLORS.get(class_id, (255, 255, 255)) # Default to white if class not found # Calculate center point x1, y1, x2, y2 = xyxy center = (int((x1 + x2) / 2), int((y1 + y2) / 2)) # Update track history self.track_history[tracker_id].append(center) # Draw trajectory line points = list(self.track_history[tracker_id]) for i in range(1, len(points)): cv2.line(frame, points[i - 1], points[i], color, thickness=2) # Draw bounding box with same color cv2.rectangle(frame, (int(x1), int(y1)), (int(x2), int(y2)), color, 2) # Display ID and class label = f"{Config.CLASS_NAMES.get(class_id, 'Unknown')} ID:{tracker_id}" cv2.putText(frame, label, (int(x1), int(y1) - 10), cv2.FONT_HERSHEY_SIMPLEX, 0.5, color, 2) return frame def _update_track_info(self, detections: sv.Detections): """Update tracking object information""" current_tracks = {} for xyxy, _, confidence, class_id, tracker_id in detections: if tracker_id is not None: current_tracks[tracker_id] = (class_id, confidence) self.current_tracks = current_tracks def process_frame(self, frame: np.ndarray) -> np.ndarray: """Process single frame""" try: # Object detection results = self.model(frame, imgsz=640, verbose=False)[0] detections = sv.Detections.from_ultralytics(results) # Object tracking with ByteTrack detections = self.tracker.update_with_detections(detections) self._update_track_info(detections) # Initialize count results counts = { 'Inbound_Vehicle': 0, 'Outbound_Vehicle': 0, 'Inbound_Pedestrian': 0, 'Outbound_Pedestrian': 0, 'Total_Vehicles': 0, 'Total_Pedestrians': 0 } recent_counts = { 'Inbound_Vehicle': 0, 'Outbound_Vehicle': 0, 'Inbound_Pedestrian': 0, 'Outbound_Pedestrian': 0 } if self.line_counter: # Update counter self.line_counter.update(detections) counts = self.line_counter.get_counts() recent_counts = self.line_counter.get_recent_counts(seconds=10) # Draw counting line cv2.line( frame, (self.line_counter.start.x, self.line_counter.start.y), (self.line_counter.end.x, self.line_counter.end.y), Config.LINE_COLOR.as_bgr(), 2 ) # Periodic data reporting if self.frame_count % Config.REPORT_INTERVAL == 0 and self.iot_client: self.iot_client.report_traffic_data(counts) # Draw object trajectories and boxes (now handled together in _draw_tracks) frame = self._draw_tracks(frame, detections) # Display statistics and table frame = self._draw_stats(frame, counts) frame = self._draw_data_table(frame, counts, recent_counts) # Set mouse callback for table dragging cv2.setMouseCallback("Traffic Monitoring", self._handle_table_drag, frame) self.frame_count += 1 return frame except Exception as e: logger.error(f"Frame processing exception: {str(e)}") return frame def interactive_line_setup(frame: np.ndarray) -> Optional[Tuple[Tuple[int, int], Tuple[int, int]]]: """Interactive counting line setup""" points = [] instructions = [ "1. Click to set line start point", "2. Click to set line end point", "3. Press ESC to confirm" ] def mouse_callback(event, x, y, flags, param): if event == cv2.EVENT_LBUTTONDOWN: points.append((x, y)) cv2.circle(frame, (x, y), 5, (0, 0, 255), -1) if len(points) > 1: cv2.line(frame, points[0], points[1], (0, 255, 0), 2) cv2.namedWindow("Set Counting Line") cv2.setMouseCallback("Set Counting Line", mouse_callback) for i, text in enumerate(instructions): cv2.putText(frame, text, (10, 30 + i * 30), cv2.FONT_HERSHEY_SIMPLEX, 0.7, (0, 255, 255), 2) while True: cv2.imshow("Set Counting Line", frame) key = cv2.waitKey(1) if key == 27 or len(points) >= 2: # ESC key or two points selected break cv2.destroyAllWindows() return (points[0], points[1]) if len(points) >= 2 else None def main(): """Main program entry""" try: # Open default camera (index 0) cap = cv2.VideoCapture(0) if not cap.isOpened(): logger.error("Unable to open camera") return # Set camera resolution (optional, adjust as needed) cap.set(cv2.CAP_PROP_FRAME_WIDTH, 1280) cap.set(cv2.CAP_PROP_FRAME_HEIGHT, 720) # Initialize monitor monitor = TrafficMonitor() # Read a few frames to allow camera to adjust for _ in range(5): ret, frame = cap.read() if not ret: logger.error("Failed to capture initial frames from camera") return # Set counting line using a frame from the camera line_points = interactive_line_setup(frame.copy()) if line_points: monitor.set_counting_line(*line_points) else: logger.info("No counting line set, will only perform object detection") # Main loop while cap.isOpened(): ret, frame = cap.read() if not ret: logger.warning("Camera frame read failed") break # Process frame processed_frame = monitor.process_frame(frame) # Display results cv2.imshow("Traffic Monitoring", processed_frame) key = cv2.waitKey(1) if key == ord('q') or key == 27: # 'q' or ESC to quit logger.info("Program terminated by user") break elif key == ord('r'): # 'r' to reset counting line line_points = interactive_line_setup(frame.copy()) if line_points: monitor.set_counting_line(*line_points) except KeyboardInterrupt: logger.info("Program interrupted by user") except Exception as e: logger.critical(f"Program exception: {str(e)}", exc_info=True) finally: # Release resources if 'cap' in locals() and cap.isOpened(): cap.release() cv2.destroyAllWindows() logger.info("Program exited") if __name__ == "__main__": # Check required dependencies try: import iot_device_sdk_python except ImportError: logger.error("Huawei IoT SDK not found, please install: pip install iot_device_sdk_python") sys.exit(1) # Ignore warning messages warnings.filterwarnings("ignore", category=UserWarning) main()为我修改这个代码改为画两条线,第一个点和第二个点是第一条线,第三个点和第四个点是第二条线,同时计算第一个点的正向行走人数-第二个点正向行走的绝对值,第二个点-第一个点正向行走的绝对值,第一条线正、逆向行走人数,总人数。第二条线正、逆向行走人数,总人数,尽量减少源代码的修改,MQTT协议传送部分不要过大改动
最新发布
08-16
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("结果的准确性取决于这些假设的合理性和参数的精确性。") 对此代码进行分析
05-12
<think>嗯,用户提到了关于LINE_LIST标记的两个错误信息:“Number of points should be a multiple of 2 for LINE_LIST marker”和“Points should not be empty for specified marker type”。我需要解释这两个错误的原因和解决方法。 首先,LINE_LIST是ROS中RViz的一个标记类型,用于绘制线段列表。每个线段需要两个点,所以总点数必须是偶数。如果用户提供的点数不是2的倍数,就会触发第一个错误。比如,如果有3个点,只能组成1条线段,剩下的第3个点没有配对,这时候就会报错。 第二个错误是说点列表不能为空。也就是说,用户在使用LINE_LIST时,必须至少提供两个点来构成一条线段。如果点列表为空或没有正确初始化,RViz就无法渲染任何内容,导致这个错误。 接下来,我需要给出解决方案。首先检查点数量是否是2的倍数,然后确认点列表是否为空。可能需要添加或删除点来满足条件。同时,还要确保正确初始化点列表,比如在ROS中正确使用geometry_msgs::Point来填充数据。 另外,用户可能是在开发ROS节点时遇到这个问题,特别是使用visualization_msgs::Marker来发布标记。需要提醒他们在发布前验证数据,避免运行时错误。可能用户对RViz的标记类型不太熟悉,或者数据处理时出现了疏漏,比如动态生成点时没有正确配对。 最后,举一个例子会更清楚。比如用四个点组成两条线段,这样点数就是2的倍数,不会出错。同时,确保每个点的坐标正确填充,避免空值或未初始化的数据。</think>好的,我将分步解释您提到的两个与 `LINE_LIST` 标记相关的错误信息: --- ### **1. "Number of points should be a multiple of 2 for LINE_LIST marker"** #### **原因分析** - **`LINE_LIST` 的定义**:在ROS的RViz可视化工具中,`LINE_LIST` 是一种标记类型,用于绘制**离散的线段集合**。每条线段需要两个点(起点和终点),因此点列表的总数必须是偶数。 - **数学本质**:若点数为 $N$,则线段数量为 $\frac{N}{2}$。若 $N$ 不是偶数,最后一对点会不完整,导致无法绘制线段。 #### **示例错误** - 错误输入:3个点 → 只能组成1条线段(使用前2个点),第3个点无法配对。 - 正确输入:4个点 → 可组成2条线段:$P_0 \rightarrow P_1$ 和 $P_2 \rightarrow P_3$。 #### **解决方法** 1. 检查点列表长度是否为偶数。 2. 若动态生成点,确保每次添加点成对出现。 3. 删除多余的点或补充缺失的点。 --- ### **2. "Points should not be empty for specified marker type"** #### **原因分析** - **空点列表**:当尝试发布 `LINE_LIST` 标记时,点列表(`points` 字段)未初始化或为空。 - **逻辑矛盾**:`LINE_LIST` 必须至少包含2个点(1条线段),否则无法绘制任何内容。 #### **示例错误** - 错误输入:发布空点列表 `[]` → RViz无法渲染。 - 正确输入:至少包含2个点,例如 `[Point(x=0,y=0), Point(x=1,y=1)]`。 #### **解决方法** 1. 初始化点列表时至少添加2个点。 2. 检查数据生成逻辑,确保点未被意外清空。 3. 添加条件判断:若点列表为空,跳过发布或填充默认值。 --- ### **完整验证流程** 1. **数据生成阶段**: ```python points = [] for i in range(num_lines): start_point = Point(x=i, y=0) # 线段起点 end_point = Point(x=i, y=1) # 线段终点 points.append(start_point) points.append(end_point) # 成对添加 ``` 2. **发布前检查**: ```python if len(points) % 2 != 0: rospy.logwarn("点数需为偶数!当前点数:%d", len(points)) elif not points: rospy.logerr("点列表为空!") else: marker.points = points publisher.publish(marker) ``` --- ### **技术背景** - **RViz标记类型**:`LINE_LIST` 适用于离散线段(如传感器检测到的边缘),而 `LINE_STRIP` 适用于连续折线。 - **数学验证**:设点集为 $P = \{P_0, P_1, \dots, P_{n-1}\}$,当 $n \mod 2 = 0$ 时,线段数量为 $\frac{n}{2}$。 --- 通过以上步骤,可确保 `LINE_LIST` 正确渲染。请检查代码中点的生成和添加逻辑是否符合这些规则。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值