import numpy as np
import pypower.api as pp
import gurobipy as gp
import matplotlib.pyplot as plt
import seaborn as sns
import pandas as pd
import time
import warnings
import matplotlib
import math
import traceback
# 设置支持中文的字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 忽略特定警告
warnings.filterwarnings("ignore", category=RuntimeWarning)
# =================== 数值稳定性增强函数 ===================
def safe_divide(a, b, default=1e-6):
"""安全除法函数,避免除以零"""
if np.isclose(b, 0, atol=1e-10):
return default
return a / b
def clean_ppc_data(ppc):
"""
清洗电力系统数据,处理无效值 (增强版)
"""
branch = ppc['branch'].copy()
bus = ppc['bus'].copy()
gen = ppc['gen'].copy()
# 1. 处理分支数据
branch[:, 2] = np.where(np.abs(branch[:, 2]) < 1e-5, 1e-5, branch[:, 2]) # 电阻
branch[:, 3] = np.where(np.abs(branch[:, 3]) < 1e-5, 0.01, branch[:, 3]) # 电抗
branch[:, 4] = np.where(np.abs(branch[:, 4]) < 1e-5, 1e-5, branch[:, 4]) # 充电电纳
branch[:, 8] = np.where(np.abs(branch[:, 8]) < 1e-5, 1.0, branch[:, 8]) # 变比
branch[:, 9] = np.nan_to_num(branch[:, 9], nan=0.0) # 角度
branch[:, 5] = np.where(branch[:, 5] <= 0, 1e6, branch[:, 5]) # 线路容量
# 2. 处理母线数据
bus[:, 11] = np.clip(bus[:, 11], 0.9, 1.1) # VMAX
bus[:, 12] = np.clip(bus[:, 12], 0.9, 1.1) # VMIN
# 3. 处理发电机数据
gen[:, 9] = np.maximum(gen[:, 9], 0) # PMIN
gen[:, 8] = np.maximum(gen[:, 8], gen[:, 9] + 1e-3) # PMAX
gen[:, 5] = np.maximum(gen[:, 5], 0.9) # VMIN
gen[:, 4] = np.minimum(np.maximum(gen[:, 4], gen[:, 5] + 0.01), 1.1) # VMAX
return {
'baseMVA': ppc['baseMVA'],
'bus': bus,
'gen': gen,
'branch': branch,
'gencost': ppc.get('gencost', None)
}
# =================== 导纳矩阵计算函数 (优化后) ===================
def makeYbus(ppc):
"""构建导纳矩阵 (优化版本)"""
baseMVA = ppc['baseMVA']
bus = ppc['bus']
branch = ppc['branch']
# 数据清洗
ppc = clean_ppc_data(ppc)
bus = ppc['bus']
branch = ppc['branch']
n_bus = bus.shape[0]
Ybus = np.zeros((n_bus, n_bus), dtype=complex)
# 处理变压器和线路参数
for i in range(branch.shape[0]):
f = int(branch[i, 0] - 1)
t = int(branch[i, 1] - 1)
# 获取线路参数 (使用安全值)
r = max(branch[i, 2], 1e-5) # 电阻最小值
x = branch[i, 3]
b = branch[i, 4] # 线路充电电纳
tap = branch[i, 8] if branch[i, 8] != 0 else 1.0
shift = np.deg2rad(branch[i, 9]) # 弧度
# 处理零阻抗
if np.isclose(r, 0, atol=1e-10) and np.isclose(x, 0, atol=1e-10):
y_series = 1e5 # 大导纳值
else:
y_series = 1 / complex(r, x)
# 移相器处理
if shift != 0:
phase_shift = np.exp(1j * shift)
else:
phase_shift = 1.0
# 标准导纳矩阵元素计算
Ytt = y_series + 1j * b / 2
Yff = Ytt / (tap * np.conj(phase_shift)) / tap
Yft = -y_series / (np.conj(phase_shift) * tap)
Ytf = -y_series / (phase_shift * tap)
# 更新导纳矩阵
Ybus[f, f] += Yff
Ybus[t, t] += Ytt
Ybus[f, t] += Yft
Ybus[t, f] += Ytf
return Ybus
# =================== Gurobi求解AC OPF (优化后 - 极坐标形式) ===================
def run_ac_opf_gurobi(ppc):
try:
# 获取系统数据
baseMVA = ppc['baseMVA']
bus = ppc['bus']
gen = ppc['gen']
branch = ppc['branch']
gencost = ppc.get('gencost', None)
# 确保数据已清洗
ppc = clean_ppc_data(ppc)
bus = ppc['bus']
gen = ppc['gen']
branch = ppc['branch']
# 系统参数
n_bus = bus.shape[0]
n_gen = gen.shape[0]
n_branch = branch.shape[0]
# 默认成本函数
if gencost is None:
gencost = np.zeros((n_gen, 6))
for i in range(n_gen):
gencost[i] = [1, 0, 0, 3, 0.01, 0.1] # 二次成本函数
# 构建导纳矩阵
Ybus = makeYbus(ppc)
G = np.real(Ybus)
B = np.imag(Ybus)
# 创建Gurobi模型
model = gp.Model('AC_OPF_Polar')
model.Params.NumericFocus = 3
model.Params.FeasibilityTol = 1e-6
model.Params.OptimalityTol = 1e-6
model.Params.NonConvex = 2
model.Params.BarConvTol = 1e-8
model.Params.TimeLimit = 600 # 600秒时间限制
# =================== 创建变量 (极坐标形式) ===================
# 电压幅值和角度
v_mag = model.addVars(n_bus, lb=bus[:, 12], ub=bus[:, 11], name='v_mag') # |V|
theta = model.addVars(n_bus, lb=-2*np.pi, ub=2*np.pi, name='theta') # θ
# 发电机出力
Pg = model.addVars(n_gen, name='Pg') # 发电机有功
Qg = model.addVars(n_gen, name='Qg') # 发电机无功
# =================== Jabr变量 (Page 32) ===================
v_sq = model.addVars(n_bus, name='v_sq') # v² = |V|²
c = model.addVars(n_branch, lb=-1e3, ub=1e3, name='c') # |Vₖ||Vₘ|cosθ
s = model.addVars(n_branch, lb=-1e3, ub=1e3, name='s') # |Vₖ||Vₘ|sinθ
# =================== 设置参考母线 ===================
ref_bus = 0
model.addConstr(theta[ref_bus] == 0.0, name='theta_ref')
# =================== 变量关系约束 ===================
# 电压幅值平方关系
for i in range(n_bus):
model.addConstr(v_sq[i] == v_mag[i] * v_mag[i], name=f'v_sq_{i}')
# Jabr变量关系
for idx in range(n_branch):
f = int(branch[idx, 0] - 1)
t = int(branch[idx, 1] - 1)
# |Vₖ||Vₘ|cos(θₖ - θₘ)
model.addConstr(c[idx] == v_mag[f] * v_mag[t] * gp.cos(theta[f] - theta[t]),
name=f'c_{f}_{t}')
# |Vₖ||Vₘ|sin(θₖ - θₘ)
model.addConstr(s[idx] == v_mag[f] * v_mag[t] * gp.sin(theta[f] - theta[t]),
name=f's_{f}_{t}')
# Jabr不等式 (Page 32)
model.addConstr(c[idx]*c[idx] + s[idx]*s[idx] <= v_sq[f] * v_sq[t],
name=f'Jabr_{f}_{t}')
# =================== 发电机约束 ===================
gen_bus = gen[:, 0].astype(int) - 1
for i in range(n_gen):
bus_idx = int(gen_bus[i])
# 有功约束
model.addConstr(Pg[i] >= gen[i, 9] / baseMVA, name=f'Pgmin_{i}')
model.addConstr(Pg[i] <= gen[i, 8] / baseMVA, name=f'Pgmax_{i}')
# 无功约束
model.addConstr(Qg[i] >= gen[i, 4] / baseMVA, name=f'Qgmin_{i}')
model.addConstr(Qg[i] <= gen[i, 3] / baseMVA, name=f'Qgmax_{i}')
# =================== 功率平衡约束 ===================
Pd = bus[:, 2] / baseMVA
Qd = bus[:, 3] / baseMVA
# 预计算发电机位置
gen_at_bus = {}
for i in range(n_bus):
gen_at_bus[i] = []
for k in range(n_gen):
bus_idx = int(gen_bus[k])
gen_at_bus[bus_idx].append(k)
# 节点注入功率计算
for i in range(n_bus):
# 有功平衡
Pi_expr = 0.0
for j in range(n_bus):
if abs(G[i, j]) > 1e-10 or abs(B[i, j]) > 1e-10:
if i == j:
Pi_expr += G[i, j] * v_sq[i]
else:
# 找到连接i和j的支路
branch_found = False
for idx in range(n_branch):
f = int(branch[idx, 0] - 1)
t = int(branch[idx, 1] - 1)
if (f == i and t == j) or (f == j and t == i):
if f == i and t == j:
Pi_expr += G[i, j] * v_sq[i] + G[i, j] * c[idx] + B[i, j] * s[idx]
else:
Pi_expr += G[i, j] * v_sq[i] + G[i, j] * c[idx] - B[i, j] * s[idx]
branch_found = True
break
if not branch_found:
# 如果没有直接支路,使用近似
theta_diff = theta[i] - theta[j]
Pi_expr += v_mag[i] * v_mag[j] * (G[i, j] * gp.cos(theta_diff) +
B[i, j] * gp.sin(theta_diff))
gen_P = gp.quicksum(Pg[k] for k in gen_at_bus.get(i, []))
model.addConstr(Pi_expr == gen_P - Pd[i], name=f'P_balance_{i}')
# 无功平衡
Qi_expr = 0.0
for j in range(n_bus):
if abs(G[i, j]) > 1e-10 or abs(B[i, j]) > 1e-10:
if i == j:
Qi_expr += -B[i, j] * v_sq[i]
else:
# 找到连接i和j的支路
branch_found = False
for idx in range(n_branch):
f = int(branch[idx, 0] - 1)
t = int(branch[idx, 1] - 1)
if (f == i and t == j) or (f == j and t == i):
if f == i and t == j:
Qi_expr += -B[i, j] * v_sq[i] - B[i, j] * c[idx] + G[i, j] * s[idx]
else:
Qi_expr += -B[i, j] * v_sq[i] - B[i, j] * c[idx] - G[i, j] * s[idx]
branch_found = True
break
if not branch_found:
# 如果没有直接支路,使用近似
theta_diff = theta[i] - theta[j]
Qi_expr += v_mag[i] * v_mag[j] * (G[i, j] * gp.sin(theta_diff) -
B[i, j] * gp.cos(theta_diff))
gen_Q = gp.quicksum(Qg[k] for k in gen_at_bus.get(i, []))
model.addConstr(Qi_expr == gen_Q - Qd[i], name=f'Q_balance_{i}')
# =================== 线路潮流约束 ===================
for idx in range(n_branch):
f = int(branch[idx, 0] - 1)
t = int(branch[idx, 1] - 1)
rateA = branch[idx, 5] / baseMVA # 线路容量 (标幺值)
# 计算线路潮流 (Page 32)
P_ft = G[f, f] * v_sq[f] + G[f, t] * c[idx] + B[f, t] * s[idx]
Q_ft = -B[f, f] * v_sq[f] - B[f, t] * c[idx] + G[f, t] * s[idx]
# 损耗不等式 (Page 37)
P_tf = G[t, t] * v_sq[t] + G[t, f] * c[idx] - B[t, f] * s[idx]
model.addConstr(P_ft + P_tf >= 0, name=f'loss_{f}_{t}')
# 线路容量约束
if rateA > 0:
model.addConstr(P_ft*P_ft + Q_ft*Q_ft <= rateA*rateA, name=f'line_cap_{idx}')
# =================== 目标函数 ===================
cost = gp.QuadExpr()
for i in range(n_gen):
model_type = gencost[i, 0]
ncost = int(gencost[i, 3])
if model_type == 1 and ncost == 3: # 二次成本
c2 = gencost[i, 4]
c1 = gencost[i, 5]
cost += c2 * (baseMVA ** 2) * Pg[i] * Pg[i] + c1 * baseMVA * Pg[i]
elif model_type == 1 and ncost == 2: # 线性成本
c1 = gencost[i, 4]
cost += c1 * baseMVA * Pg[i]
else: # 默认线性
c1 = gen[i, 6] if gen.shape[1] > 6 else 10.0
cost += c1 * baseMVA * Pg[i]
model.setObjective(cost, gp.GRB.MINIMIZE)
# =================== 求解模型 ===================
start_time = time.time()
model.optimize()
solve_time = time.time() - start_time
# 结果状态处理
solution_status = "最优解"
if model.status == gp.GRB.OPTIMAL:
success = True
elif model.status == gp.GRB.TIME_LIMIT and model.SolCount > 0:
success = True
solution_status = f"时间限制内找到可行解 (目标值: {model.ObjVal:.4f})"
else:
success = False
solution_status = f"求解失败 (状态: {model.status})"
# 处理结果
results = {
'success': success,
'status': solution_status,
'f': model.ObjVal if hasattr(model, 'ObjVal') else float('nan'),
'gen': np.zeros((n_gen, 4)),
'bus': np.zeros((n_bus, 4)),
'solve_time': solve_time
}
# 提取结果
if model.SolCount > 0:
# 发电机结果
for i in range(n_gen):
results['gen'][i, 0] = gen_bus[i] + 1
results['gen'][i, 1] = Pg[i].X * baseMVA
results['gen'][i, 2] = Qg[i].X * baseMVA
# 母线结果
for i in range(n_bus):
results['bus'][i, 0] = i + 1
results['bus'][i, 1] = v_mag[i].X
results['bus'][i, 2] = np.rad2deg(theta[i].X)
else:
# 没有找到可行解,使用初始值
for i in range(n_gen):
results['gen'][i, 0] = gen_bus[i] + 1
results['gen'][i, 1] = gen[i, 9]
results['gen'][i, 2] = 0.0
for i in range(n_bus):
results['bus'][i, 0] = i + 1
results['bus'][i, 1] = 1.0
results['bus'][i, 2] = 0.0
return results
except Exception as e:
traceback.print_exc()
n_bus = ppc['bus'].shape[0] if 'bus' in ppc else 1
n_gen = ppc['gen'].shape[0] if 'gen' in ppc else 1
return {
'success': False,
'status': f"异常: {str(e)}",
'f': float('nan'),
'gen': np.zeros((n_gen, 4)) * float('nan'),
'bus': np.zeros((n_bus, 4)) * float('nan'),
'solve_time': 0
}
# =================== PyPower求解AC OPF ===================
def run_ac_opf_pypower(ppc):
"""使用PyPower求解AC最优潮流(增强错误处理)"""
try:
# 确保数据已清洗
ppc = clean_ppc_data(ppc)
# 设置求解选项
ppopt = pp.ppoption(VERBOSE=0, OUT_ALL=0, OPF_FLOW_LIM=1)
# 求解AC OPF
start_time = time.time()
results = pp.runopf(ppc, ppopt)
solve_time = time.time() - start_time
if results['success']:
return {
'success': True,
'status': "最优解",
'f': results['f'],
'gen': results['gen'],
'bus': results['bus'],
'solve_time': solve_time
}
else:
# 尝试收集部分结果
partial_results = {
'success': False,
'status': "求解失败",
'f': float('nan'),
'solve_time': solve_time
}
# 尝试收集可用的结果
try:
if 'gen' in results:
partial_results['gen'] = results['gen']
else:
partial_results['gen'] = ppc['gen'].copy()
partial_results['gen'][:, 1] = ppc['gen'][:, 9]
partial_results['gen'][:, 2] = 0.0
except:
partial_results['gen'] = ppc['gen'].copy()
partial_results['gen'][:, 1] = ppc['gen'][:, 9]
partial_results['gen'][:, 2] = 0.0
try:
if 'bus' in results:
partial_results['bus'] = results['bus']
else:
partial_results['bus'] = ppc['bus'].copy()
partial_results['bus'][:, 7] = 1.0
partial_results['bus'][:, 8] = 0.0
except:
partial_results['bus'] = ppc['bus'].copy()
partial_results['bus'][:, 7] = 1.0
partial_results['bus'][:, 8] = 0.0
return partial_results
except Exception as e:
import traceback
traceback.print_exc()
# 创建近似结果
n_bus = ppc['bus'].shape[0] if 'bus' in ppc else 1
n_gen = ppc['gen'].shape[0] if 'gen' in ppc else 1
# 创建近似解
gen_results = np.zeros((n_gen, 4))
if n_gen > 0:
gen_results[:, 0] = ppc['gen'][:, 0]
gen_results[:, 1] = ppc['gen'][:, 9]
gen_results[:, 2] = 0.0
bus_results = np.zeros((n_bus, 4))
if n_bus > 0:
bus_results[:, 0] = ppc['bus'][:, 0]
bus_results[:, 1] = 1.0
bus_results[:, 2] = 0.0
return {
'success': False,
'status': f"异常: {str(e)}",
'f': float('nan'),
'gen': gen_results,
'bus': bus_results,
'solve_time': 0
}
# =================== 比较分析和可视化函数 ===================
def run_comparison_analysis(case_name):
"""运行两种求解器的比较分析(增强鲁棒性)"""
try:
# 加载测试案例
if case_name == 'case9':
ppc_data = pp.case9()
elif case_name == 'case14':
ppc_data = pp.case14()
elif case_name == 'case30':
ppc_data = pp.case30()
elif case_name == 'case118':
ppc_data = pp.case118()
else:
raise ValueError(f"未知测试案例: {case_name}")
baseMVA = ppc_data['baseMVA']
n_bus = ppc_data['bus'].shape[0]
n_gen = ppc_data['gen'].shape[0]
# 使用Gurobi求解
print(f"\n使用Gurobi求解{case_name}...")
gurobi_results = run_ac_opf_gurobi(ppc_data)
# 使用PyPower求解
print(f"\n使用PyPower求解{case_name}...")
pypower_results = run_ac_opf_pypower(ppc_data)
# 创建结果比较数据框
# 1. 母线结果比较
bus_comparison = pd.DataFrame({
'Bus': range(1, n_bus + 1)
})
# 添加Gurobi结果
bus_comparison['V_gurobi'] = gurobi_results['bus'][:, 1]
bus_comparison['Angle_gurobi'] = gurobi_results['bus'][:, 2]
# 添加PyPower结果
if pypower_results['bus'].shape[1] > 8:
bus_comparison['V_pypower'] = pypower_results['bus'][:, 7]
bus_comparison['Angle_pypower'] = pypower_results['bus'][:, 8]
else:
bus_comparison['V_pypower'] = pypower_results['bus'][:, 1]
bus_comparison['Angle_pypower'] = np.zeros(n_bus)
# 计算差异
bus_comparison['V_diff'] = np.abs(bus_comparison['V_gurobi'] - bus_comparison['V_pypower'])
bus_comparison['Angle_diff'] = np.abs(bus_comparison['Angle_gurobi'] - bus_comparison['Angle_pypower'])
# 2. 发电机结果比较
gen_comparison = pd.DataFrame({
'Gen': ppc_data['gen'][:, 0]
})
# 添加Gurobi结果
gen_comparison['P_gurobi'] = gurobi_results['gen'][:, 1]
gen_comparison['Q_gurobi'] = gurobi_results['gen'][:, 2]
# 添加PyPower结果
gen_comparison['P_pypower'] = pypower_results['gen'][:, 1]
gen_comparison['Q_pypower'] = pypower_results['gen'][:, 2]
# 计算差异
gen_comparison['P_diff'] = np.abs(gen_comparison['P_gurobi'] - gen_comparison['P_pypower'])
gen_comparison['Q_diff'] = np.abs(gen_comparison['Q_gurobi'] - gen_comparison['Q_pypower'])
# 打印求解状态
print(f"\n求解状态 ({case_name}):")
print(f" Gurobi: {gurobi_results['status']}")
print(f" PyPower: {pypower_results['status']}")
# 打印求解时间
print(f"\n求解时间 ({case_name}):")
print(f" Gurobi: {gurobi_results['solve_time']:.4f} 秒")
print(f" PyPower: {pypower_results['solve_time']:.4f} 秒")
# 打印目标函数值
if not np.isnan(gurobi_results['f']):
print(f"\nGurobi目标函数值: {gurobi_results['f']:.2f}")
if not np.isnan(pypower_results['f']):
print(f"PyPower目标函数值: {pypower_results['f']:.2f}")
return bus_comparison, gen_comparison, gurobi_results, pypower_results
except Exception as e:
import traceback
traceback.print_exc()
print(f"在分析案例 {case_name} 时发生错误: {str(e)}")
# 创建空结果
bus_comp = pd.DataFrame(
columns=['Bus', 'V_gurobi', 'Angle_gurobi', 'V_pypower', 'Angle_pypower', 'V_diff', 'Angle_diff'])
gen_comp = pd.DataFrame(columns=['Gen', 'P_gurobi', 'Q_gurobi', 'P_pypower', 'Q_pypower', 'P_diff', 'Q_diff'])
return bus_comp, gen_comp, None, None
def plot_comparison_results(bus_comp, gen_comp, case_name):
"""可视化比较结果"""
plt.figure(figsize=(15, 10))
# 电压比较图
plt.subplot(2, 2, 1)
plt.plot(bus_comp['Bus'], bus_comp['V_gurobi'], 'bo-', label='Gurobi')
plt.plot(bus_comp['Bus'], bus_comp['V_pypower'], 'r^--', label='PyPower')
plt.xlabel('母线编号')
plt.ylabel('电压幅值 (p.u.)')
plt.title(f'{case_name} - 电压幅值比较')
plt.legend()
plt.grid(True)
# 电压差异图
plt.subplot(2, 2, 2)
plt.bar(bus_comp['Bus'], bus_comp['V_diff'])
plt.xlabel('母线编号')
plt.ylabel('电压差异 (p.u.)')
plt.title(f'{case_name} - 电压幅值差异')
plt.grid(True)
# 有功功率比较图
plt.subplot(2, 2, 3)
plt.plot(gen_comp['Gen'], gen_comp['P_gurobi'], 'bo-', label='Gurobi')
plt.plot(gen_comp['Gen'], gen_comp['P_pypower'], 'r^--', label='PyPower')
plt.xlabel('发电机母线')
plt.ylabel('有功出力 (MW)')
plt.title(f'{case_name} - 发电机有功比较')
plt.legend()
plt.grid(True)
# 有功功率差异图
plt.subplot(2, 2, 4)
plt.bar(gen_comp['Gen'], gen_comp['P_diff'])
plt.xlabel('发电机母线')
plt.ylabel('有功差异 (MW)')
plt.title(f'{case_name} - 发电机有功差异')
plt.grid(True)
plt.tight_layout()
plt.savefig(f'{case_name}_comparison.png', dpi=300)
plt.close()
# =================== 主程序 ===================
if __name__ == "__main__":
# 测试不同案例
test_cases = ['case9', 'case14', 'case30', 'case118']
for case in test_cases:
print(f"\n{'=' * 50}")
print(f"开始分析案例: {case}")
print(f"{'=' * 50}")
# 运行比较分析
bus_comp, gen_comp, gurobi_res, pypower_res = run_comparison_analysis(case)
# 打印节点结果
if not bus_comp.empty:
print("\n节点结果:")
print(
f"{'母线':<5} | {'电压(Gurobi)':<12} | {'电压(PyPower)':<12} | {'角度(Gurobi)':<12} | {'角度(PyPower)':<12}")
for _, row in bus_comp.iterrows():
print(f"{int(row['Bus']):<5} | {row['V_gurobi']:.6f} | {row['V_pypower']:.6f} | "
f"{row['Angle_gurobi']:.6f} | {row['Angle_pypower']:.6f}")
# 打印发电机结果
if not gen_comp.empty:
print("\n发电机结果:")
print(f"{'母线':<5} | {'Pg(Gurobi)':<12} | {'Pg(PyPower)':<12} | {'Qg(Gurobi)':<12} | {'Qg(PyPower)':<12}")
for _, row in gen_comp.iterrows():
print(f"{int(row['Gen']):<5} | {row['P_gurobi']:.6f} | {row['P_pypower']:.6f} | "
f"{row['Q_gurobi']:.6f} | {row['Q_pypower']:.6f}")
# 保存结果到CSV
if not bus_comp.empty:
bus_comp.to_csv(f'{case}_bus_comparison.csv', index=False)
if not gen_comp.empty:
gen_comp.to_csv(f'{case}_gen_comparison.csv', index=False)
# 可视化结果
if not bus_comp.empty and not gen_comp.empty:
try:
plot_comparison_results(bus_comp, gen_comp, case)
# 打印最大差异
print(f"\n最大电压差异: {bus_comp['V_diff'].max():.6f} p.u.")
print(f"最大角度差异: {bus_comp['Angle_diff'].max():.6f} 度")
print(f"最大有功差异: {gen_comp['P_diff'].max():.6f} MW")
print(f"最大无功差异: {gen_comp['Q_diff'].max():.6f} MVar")
except Exception as e:
print(f"可视化失败: {str(e)}") 这是我的代码,出现了不能使用cos的报错
最新发布