thrustCoefficientCurve.py 37 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849
  1. import pandas as pd
  2. import numpy as np
  3. import math
  4. from pathlib import Path
  5. import matplotlib.pyplot as plt
  6. import seaborn as sns
  7. from typing import Optional, Tuple, Dict, List
  8. import re
  9. class AerodynamicThrustCalculator:
  10. """气动推力系数计算器类 - 按标准机型和描述分组"""
  11. def __init__(self, csv_path: str):
  12. """
  13. 初始化计算器
  14. Args:
  15. csv_path: CSV文件路径
  16. """
  17. self.csv_path = csv_path
  18. self.data = None
  19. self.results = None
  20. self.grouped_results = None
  21. self.figure_size = (12, 8)
  22. self.total_groups = 0
  23. def load_data(self) -> bool:
  24. """从CSV文件加载数据"""
  25. try:
  26. self.data = pd.read_csv(self.csv_path)
  27. print(f"成功加载数据,共 {len(self.data)} 条记录")
  28. print(f"包含 {self.data['标准机型'].nunique()} 种标准机型")
  29. print(f"包含 {self.data['描述'].nunique()} 种不同的描述")
  30. # 计算总分组数
  31. self.total_groups = self.data.groupby(['标准机型', '描述']).ngroups
  32. print(f"共 {self.total_groups} 个分组")
  33. return True
  34. except FileNotFoundError:
  35. print(f"错误:找不到文件 {self.csv_path}")
  36. return False
  37. except Exception as e:
  38. print(f"加载数据时发生错误: {e}")
  39. return False
  40. def validate_data(self) -> bool:
  41. """验证数据完整性"""
  42. required_columns = ['标准机型', '描述', '风速', '有功功率', '空气密度']
  43. for col in required_columns:
  44. if col not in self.data.columns:
  45. print(f"错误:缺失必要列 '{col}'")
  46. return False
  47. # 检查数据类型
  48. try:
  49. self.data['风速'] = pd.to_numeric(self.data['风速'])
  50. self.data['有功功率'] = pd.to_numeric(self.data['有功功率'])
  51. self.data['空气密度'] = pd.to_numeric(self.data['空气密度'])
  52. except ValueError as e:
  53. print(f"数据类型转换错误: {e}")
  54. return False
  55. print("数据验证通过")
  56. return True
  57. def extract_rotor_diameter_from_model(self, turbine_model: str) -> float:
  58. """
  59. 从机型名称中提取叶轮直径
  60. Args:
  61. turbine_model: 机型名称,如 'G58-850', 'MY2.0se-121'
  62. Returns:
  63. 叶轮直径(米)
  64. """
  65. # 常见机型叶轮直径映射(可根据实际情况扩展)
  66. diameter_map = {
  67. # 歌美飒
  68. 'G58': 58.0, 'G58-850': 58.0,
  69. # 明阳智能
  70. 'MY2.0se-121': 121.0,
  71. 'MY3.0se-145': 145.0,
  72. 'MY4.0se-166': 166.0,
  73. 'MY6.25se-172': 172.0,
  74. # 东方电气
  75. 'DEW-G5500-183': 183.0,
  76. # 金风科技
  77. 'GW82-1500': 82.0, 'GW82/1500': 82.0,
  78. 'GW140-3300': 140.0, 'GW140/3300': 140.0,
  79. # 上海电气
  80. 'W2000-111': 111.0,
  81. # 远景能源
  82. 'EN141-2650': 141.0, 'EN-141/2.65': 141.0,
  83. 'EN182-6250': 182.0,
  84. # 国电联合动力
  85. 'CCWE3000-122.HD': 122.0,
  86. 'CCWE1500-82.DF': 82.0,
  87. # 维斯塔斯
  88. 'V52-850': 52.0,
  89. # 湘电风能
  90. 'XE72-2000': 72.0,
  91. 'XE122-2500': 122.0,
  92. # 华锐风电
  93. 'SL1500-82': 82.0, 'SL-82/1.5': 82.0,
  94. # 其他
  95. 'S48-750': 48.0
  96. }
  97. # 尝试直接匹配
  98. if turbine_model in diameter_map:
  99. return diameter_map[turbine_model]
  100. # 尝试从名称中提取数字
  101. try:
  102. # 匹配两个连字符之间的数字,如 "G58-850" 中的 58
  103. pattern = r'-(\d+\.?\d*)-'
  104. match = re.search(pattern, turbine_model)
  105. if match:
  106. return float(match.group(1))
  107. # 匹配末尾的数字,如 "MY2.0se-121" 中的 121
  108. pattern2 = r'-(\d+\.?\d*)$'
  109. match2 = re.search(pattern2, turbine_model)
  110. if match2:
  111. return float(match2.group(1))
  112. # 匹配任意位置的数字
  113. pattern3 = r'(\d+\.?\d*)'
  114. matches = re.findall(pattern3, turbine_model)
  115. if matches:
  116. # 取最大的数字作为叶轮直径(通常叶轮直径是较大的数字)
  117. diameters = [float(m) for m in matches]
  118. # 过滤掉可能的小数字(如2.0中的2)
  119. filtered_diameters = [d for d in diameters if d > 20]
  120. if filtered_diameters:
  121. return max(filtered_diameters)
  122. elif diameters:
  123. return max(diameters)
  124. print(f"警告:无法从机型名称 '{turbine_model}' 中提取叶轮直径,使用默认值 80m")
  125. return 80.0 # 默认值
  126. except Exception as e:
  127. print(f"提取叶轮直径时出错: {e}")
  128. return 80.0 # 默认值
  129. def extract_rotor_diameter_from_description(self, description: str) -> Optional[float]:
  130. """
  131. 从描述中提取叶轮直径
  132. Args:
  133. description: 描述文本
  134. Returns:
  135. 叶轮直径(米),如果无法提取则返回None
  136. """
  137. try:
  138. # 尝试从描述中提取直径信息
  139. patterns = [
  140. r'(\d+\.?\d*)[mM]\s*(直径|叶轮|叶片)',
  141. r'直径\s*(\d+\.?\d*)\s*[mM]',
  142. r'叶片长度\s*(\d+\.?\d*)\s*[mM]',
  143. r'风轮直径\s*(\d+\.?\d*)\s*[mM]',
  144. r'D\s*[=:]\s*(\d+\.?\d*)\s*[mM]'
  145. ]
  146. for pattern in patterns:
  147. match = re.search(pattern, description)
  148. if match:
  149. return float(match.group(1))
  150. return None
  151. except Exception as e:
  152. print(f"从描述中提取叶轮直径时出错: {e}")
  153. return None
  154. def get_rotor_diameter(self, turbine_model: str, description: str) -> float:
  155. """
  156. 获取叶轮直径,优先从描述中提取,然后从机型名称中提取
  157. Args:
  158. turbine_model: 机型名称
  159. description: 描述文本
  160. Returns:
  161. 叶轮直径(米)
  162. """
  163. # 首先尝试从描述中提取
  164. diameter_from_desc = self.extract_rotor_diameter_from_description(description)
  165. if diameter_from_desc is not None:
  166. print(f"从描述中提取叶轮直径: {diameter_from_desc}m (机型: {turbine_model})")
  167. return diameter_from_desc
  168. # 如果无法从描述中提取,则从机型名称中提取
  169. diameter_from_model = self.extract_rotor_diameter_from_model(turbine_model)
  170. print(f"从机型名称中提取叶轮直径: {diameter_from_model}m (机型: {turbine_model})")
  171. return diameter_from_model
  172. def calculate_thrust_coefficient_for_group(self, group_data: pd.DataFrame,
  173. turbine_model: str, description: str) -> pd.DataFrame:
  174. """
  175. 计算单个分组的气动推力系数
  176. Args:
  177. group_data: 分组数据
  178. turbine_model: 标准机型
  179. description: 描述
  180. Returns:
  181. 包含计算结果的DataFrame
  182. """
  183. # 创建副本
  184. result_df = group_data.copy()
  185. # 获取叶轮直径
  186. rotor_diameter = self.get_rotor_diameter(turbine_model, description)
  187. # 计算扫风面积 (A = π × (D/2)²)
  188. swept_area = np.pi * (rotor_diameter / 2) ** 2
  189. # 转换功率单位:kW → W
  190. power_watts = result_df['有功功率'] * 1000
  191. # 获取其他参数
  192. air_density = result_df['空气密度']
  193. wind_speed = result_df['风速']
  194. # 初始化结果列
  195. result_df['叶轮直径_D_m'] = rotor_diameter
  196. result_df['扫风面积_A_m2'] = swept_area
  197. result_df['气动推力_F_N'] = 0.0
  198. result_df['气动推力系数_Ct'] = 0.0
  199. # 过滤掉风速为0的行(避免除零错误)
  200. valid_indices = wind_speed > 0
  201. if valid_indices.any():
  202. # 计算气动推力: F = 2 * ρ * A * P / v
  203. F = (2 * air_density[valid_indices] * swept_area *
  204. power_watts[valid_indices] / wind_speed[valid_indices])
  205. # 计算气动推力系数: Ct = (2 × F) / (ρ * v² * A)
  206. Ct = (2 * F) / (air_density[valid_indices] *
  207. wind_speed[valid_indices] ** 2 * swept_area)
  208. # 赋值回结果DataFrame
  209. result_df.loc[valid_indices, '气动推力_F_N'] = F
  210. result_df.loc[valid_indices, '气动推力系数_Ct'] = Ct
  211. # 风速为0时,推力系数设为0
  212. zero_indices = wind_speed == 0
  213. if zero_indices.any():
  214. result_df.loc[zero_indices, '气动推力系数_Ct'] = 0.0
  215. # 添加分组标识(缩短描述以避免文件名过长)
  216. short_desc = description[:50].replace('/', '_').replace('\\', '_').replace(':', '_')
  217. short_desc = re.sub(r'[<>:"|?*]', '_', short_desc) # 移除Windows文件名非法字符
  218. result_df['分组标识'] = f"{turbine_model}_{short_desc}"
  219. return result_df
  220. def calculate_thrust_coefficient(self) -> bool:
  221. """按标准机型和描述分组计算气动推力系数"""
  222. if self.data is None:
  223. print("错误:请先加载数据")
  224. return False
  225. # 按标准机型和描述分组
  226. grouped_data = self.data.groupby(['标准机型', '描述'])
  227. all_results = []
  228. print(f"\n开始计算 {self.total_groups} 个分组的气动推力系数...")
  229. for (turbine_model, description), group in grouped_data:
  230. print(f" 计算: 机型={turbine_model}, 描述={description[:50]}... ({len(group)} 条数据)")
  231. # 计算当前分组
  232. group_result = self.calculate_thrust_coefficient_for_group(
  233. group, turbine_model, description
  234. )
  235. all_results.append(group_result)
  236. # 合并所有结果
  237. if all_results:
  238. self.results = pd.concat(all_results, ignore_index=True)
  239. self.grouped_results = self.results.copy()
  240. print(f"\n计算完成!共计算 {self.total_groups} 个分组")
  241. return True
  242. else:
  243. print("错误:没有数据可用于计算")
  244. return False
  245. def save_results(self, output_path: str = None) -> bool:
  246. """保存计算结果到CSV文件"""
  247. if self.results is None:
  248. print("错误:请先计算气动推力系数")
  249. return False
  250. if output_path is None:
  251. # 生成默认输出路径
  252. original_path = Path(self.csv_path)
  253. output_path = original_path.parent / f"{original_path.stem}_计算结果.csv"
  254. try:
  255. # 选择要保存的列(保持原有列)
  256. save_columns = list(self.results.columns)
  257. self.results.to_csv(output_path, index=False, encoding='utf-8-sig')
  258. print(f"完整结果已保存到: {output_path}")
  259. # 同时保存一个简化版本
  260. simplified_path = Path(output_path).parent / f"{Path(output_path).stem}_简化版.csv"
  261. simplified_cols = [
  262. '风场id', '风场名称', '标准机型', '原始机型', '风速',
  263. '有功功率', '空气密度', '叶轮直径_D_m', '扫风面积_A_m2',
  264. '气动推力_F_N', '气动推力系数_Ct', '分组标识', '描述'
  265. ]
  266. existing_cols = [col for col in simplified_cols if col in self.results.columns]
  267. self.results[existing_cols].to_csv(simplified_path, index=False, encoding='utf-8-sig')
  268. print(f"简化结果已保存到: {simplified_path}")
  269. return True
  270. except Exception as e:
  271. print(f"保存结果时发生错误: {e}")
  272. return False
  273. def save_grouped_summary(self, output_path: str = None) -> bool:
  274. """保存分组统计摘要到CSV文件"""
  275. if self.results is None:
  276. print("错误:请先计算气动推力系数")
  277. return False
  278. if output_path is None:
  279. original_path = Path(self.csv_path)
  280. output_path = original_path.parent / f"{original_path.stem}_分组摘要.csv"
  281. try:
  282. summary_data = []
  283. # 按分组标识分组
  284. for group_id, group_data in self.results.groupby('分组标识'):
  285. # 获取标准机型(从分组标识中提取)
  286. parts = group_id.split('_')
  287. turbine_model = parts[0] if parts else ''
  288. # 获取有效数据(风速>0)
  289. valid_data = group_data[group_data['风速'] > 0]
  290. if len(valid_data) > 0:
  291. summary = {
  292. '标准机型': turbine_model,
  293. '描述': group_data['描述'].iloc[0] if '描述' in group_data.columns else '',
  294. '分组标识': group_id,
  295. '数据点数': len(group_data),
  296. '有效计算点数': len(valid_data),
  297. '叶轮直径_m': group_data['叶轮直径_D_m'].iloc[0] if '叶轮直径_D_m' in group_data.columns else '未知',
  298. '扫风面积_m2': group_data['扫风面积_A_m2'].iloc[0] if '扫风面积_A_m2' in group_data.columns else '未知',
  299. '空气密度_kg_m3': group_data['空气密度'].mean(),
  300. '风速范围_m_s': f"{group_data['风速'].min():.1f}-{group_data['风速'].max():.1f}",
  301. '最大推力系数': valid_data['气动推力系数_Ct'].max(),
  302. '最小推力系数': valid_data['气动推力系数_Ct'].min(),
  303. '平均推力系数': valid_data['气动推力系数_Ct'].mean(),
  304. '推力系数标准差': valid_data['气动推力系数_Ct'].std(),
  305. '最大推力_N': valid_data['气动推力_F_N'].max(),
  306. '最小推力_N': valid_data['气动推力_F_N'].min(),
  307. '平均推力_N': valid_data['气动推力_F_N'].mean(),
  308. '额定功率_kW': group_data['有功功率'].max(),
  309. '切入风速_m_s': group_data[group_data['有功功率'] > 0]['风速'].min() if len(group_data[group_data['有功功率'] > 0]) > 0 else 0,
  310. '额定风速_m_s': group_data[group_data['有功功率'] == group_data['有功功率'].max()]['风速'].min() if len(group_data[group_data['有功功率'] == group_data['有功功率'].max()]) > 0 else 0
  311. }
  312. summary_data.append(summary)
  313. summary_df = pd.DataFrame(summary_data)
  314. # 按标准机型排序
  315. if '标准机型' in summary_df.columns:
  316. summary_df = summary_df.sort_values(['标准机型', '平均推力系数'], ascending=[True, False])
  317. summary_df.to_csv(output_path, index=False, encoding='utf-8-sig')
  318. print(f"分组摘要已保存到: {output_path}")
  319. # 同时保存每个机型的汇总统计
  320. if '标准机型' in summary_df.columns:
  321. model_summary_path = Path(output_path).parent / f"{Path(output_path).stem}_机型汇总.csv"
  322. model_summary = summary_df.groupby('标准机型').agg({
  323. '数据点数': 'sum',
  324. '有效计算点数': 'sum',
  325. '平均推力系数': 'mean',
  326. '推力系数标准差': 'mean',
  327. '平均推力_N': 'mean',
  328. '叶轮直径_m': 'first'
  329. }).reset_index()
  330. model_summary.to_csv(model_summary_path, index=False, encoding='utf-8-sig')
  331. print(f"机型汇总已保存到: {model_summary_path}")
  332. return True
  333. except Exception as e:
  334. print(f"保存分组摘要时发生错误: {e}")
  335. return False
  336. def save_thrust_curve_data(self, output_path: str = None) -> bool:
  337. """保存推力系数曲线数据(每个风速点的推力系数)"""
  338. if self.results is None:
  339. print("错误:请先计算气动推力系数")
  340. return False
  341. if output_path is None:
  342. original_path = Path(self.csv_path)
  343. output_path = original_path.parent / f"{original_path.stem}_推力曲线数据.csv"
  344. try:
  345. # 选择关键列
  346. curve_columns = [
  347. '分组标识', '标准机型', '描述', '风速',
  348. '有功功率', '气动推力系数_Ct', '气动推力_F_N',
  349. '叶轮直径_D_m', '空气密度'
  350. ]
  351. # 只保留存在的列
  352. existing_columns = [col for col in curve_columns if col in self.results.columns]
  353. curve_data = self.results[existing_columns].copy()
  354. # 按分组和风速排序
  355. if '分组标识' in curve_data.columns and '风速' in curve_data.columns:
  356. curve_data = curve_data.sort_values(['分组标识', '风速'])
  357. curve_data.to_csv(output_path, index=False, encoding='utf-8-sig')
  358. print(f"推力曲线数据已保存到: {output_path}")
  359. return True
  360. except Exception as e:
  361. print(f"保存推力曲线数据时发生错误: {e}")
  362. return False
  363. def get_summary(self) -> Dict:
  364. """获取计算结果的统计摘要"""
  365. if self.results is None:
  366. return None
  367. # 获取有效数据(风速>0)
  368. valid_data = self.results[self.results['风速'] > 0]
  369. summary = {
  370. '总记录数': len(self.results),
  371. '有效计算记录数': len(valid_data),
  372. '总分组数': self.total_groups,
  373. '标准机型数量': self.results['标准机型'].nunique(),
  374. '描述类型数量': self.results['描述'].nunique() if '描述' in self.results.columns else 0,
  375. '最大气动推力系数': valid_data['气动推力系数_Ct'].max() if len(valid_data) > 0 else 0,
  376. '最小气动推力系数': valid_data['气动推力系数_Ct'].min() if len(valid_data) > 0 else 0,
  377. '平均气动推力系数': valid_data['气动推力系数_Ct'].mean() if len(valid_data) > 0 else 0,
  378. '最大气动推力_N': f"{valid_data['气动推力_F_N'].max():.1f}" if len(valid_data) > 0 else "0",
  379. '风速范围_m_s': f"{self.results['风速'].min():.1f} - {self.results['风速'].max():.1f}",
  380. '叶轮直径范围_m': f"{self.results['叶轮直径_D_m'].min():.1f} - {self.results['叶轮直径_D_m'].max():.1f}" if '叶轮直径_D_m' in self.results.columns else "未知",
  381. '空气密度范围_kg_m3': f"{self.results['空气密度'].min():.3f} - {self.results['空气密度'].max():.3f}" if '空气密度' in self.results.columns else "未知"
  382. }
  383. return summary
  384. def print_summary(self) -> None:
  385. """打印计算摘要"""
  386. summary = self.get_summary()
  387. if summary:
  388. print("\n" + "="*60)
  389. print("气动推力系数计算摘要")
  390. print("="*60)
  391. for key, value in summary.items():
  392. print(f"{key:25}: {value}")
  393. print("="*60)
  394. def plot_all_groups_thrust_curves(self, output_dir: str = None) -> bool:
  395. """
  396. 为所有分组绘制气动推力系数曲线
  397. Args:
  398. output_dir: 输出目录
  399. """
  400. if self.results is None:
  401. print("错误:请先计算气动推力系数")
  402. return False
  403. # 设置中文字体
  404. plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans']
  405. plt.rcParams['axes.unicode_minus'] = False
  406. # 创建输出目录
  407. if output_dir is None:
  408. output_dir = "output/plots"
  409. Path(output_dir).mkdir(parents=True, exist_ok=True)
  410. # 获取所有分组
  411. groups = self.results['分组标识'].unique()
  412. total_groups = len(groups)
  413. print(f"\n开始为所有 {total_groups} 个分组绘制推力系数曲线...")
  414. # 记录绘图进度
  415. plot_count = 0
  416. failed_plots = []
  417. # 为每个分组绘制图表
  418. for i, group_id in enumerate(groups, 1):
  419. try:
  420. group_data = self.results[self.results['分组标识'] == group_id]
  421. # 创建图形
  422. fig, axes = plt.subplots(2, 2, figsize=(14, 10))
  423. # 获取机型信息
  424. turbine_model = group_data['标准机型'].iloc[0] if '标准机型' in group_data.columns else '未知'
  425. description = group_data['描述'].iloc[0] if '描述' in group_data.columns else ''
  426. short_desc = description[:30] + "..." if len(description) > 30 else description
  427. # 1. 气动推力系数曲线
  428. ax1 = axes[0, 0]
  429. valid_data = group_data[group_data['风速'] > 0]
  430. if len(valid_data) > 0:
  431. ax1.plot(valid_data['风速'], valid_data['气动推力系数_Ct'],
  432. 'b-', linewidth=2, marker='o', markersize=4)
  433. # 标记最大推力系数点
  434. max_ct_idx = valid_data['气动推力系数_Ct'].idxmax()
  435. max_ct_wind = valid_data.loc[max_ct_idx, '风速']
  436. max_ct_value = valid_data.loc[max_ct_idx, '气动推力系数_Ct']
  437. ax1.plot(max_ct_wind, max_ct_value, 'ro', markersize=8)
  438. # 添加统计信息
  439. stats_text = f'平均Ct: {valid_data["气动推力系数_Ct"].mean():.3f}\n'
  440. stats_text += f'最大Ct: {max_ct_value:.3f}\n'
  441. stats_text += f'风速范围: {group_data["风速"].min():.1f}-{group_data["风速"].max():.1f} m/s'
  442. ax1.text(0.02, 0.98, stats_text, transform=ax1.transAxes,
  443. fontsize=10, verticalalignment='top',
  444. bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
  445. ax1.set_xlabel('风速 (m/s)', fontsize=12)
  446. ax1.set_ylabel('气动推力系数 Ct', fontsize=12)
  447. ax1.set_title('气动推力系数曲线', fontsize=14, fontweight='bold')
  448. ax1.grid(True, alpha=0.3)
  449. # 2. 功率曲线
  450. ax2 = axes[0, 1]
  451. ax2.plot(group_data['风速'], group_data['有功功率'],
  452. 'g-', linewidth=2, marker='s', markersize=4)
  453. ax2.set_xlabel('风速 (m/s)', fontsize=12)
  454. ax2.set_ylabel('功率 (kW)', fontsize=12)
  455. ax2.set_title('功率曲线', fontsize=14, fontweight='bold')
  456. ax2.grid(True, alpha=0.3)
  457. # 标记额定功率
  458. rated_power = group_data['有功功率'].max()
  459. if rated_power > 0:
  460. rated_wind = group_data[group_data['有功功率'] == rated_power]['风速'].min()
  461. ax2.plot(rated_wind, rated_power, 'ro', markersize=8)
  462. ax2.annotate(f'额定: {rated_power} kW\n风速: {rated_wind:.1f} m/s',
  463. xy=(rated_wind, rated_power),
  464. xytext=(rated_wind + 2, rated_power * 0.7),
  465. arrowprops=dict(arrowstyle='->', color='red'),
  466. fontsize=10, color='red')
  467. # 3. 气动推力曲线
  468. ax3 = axes[1, 0]
  469. if len(valid_data) > 0:
  470. ax3.plot(valid_data['风速'], valid_data['气动推力_F_N'],
  471. 'r-', linewidth=2, marker='^', markersize=4)
  472. ax3.set_xlabel('风速 (m/s)', fontsize=12)
  473. ax3.set_ylabel('气动推力 (N)', fontsize=12)
  474. ax3.set_title('气动推力曲线', fontsize=14, fontweight='bold')
  475. ax3.grid(True, alpha=0.3)
  476. # 标记最大推力点
  477. max_thrust_idx = valid_data['气动推力_F_N'].idxmax()
  478. max_thrust_wind = valid_data.loc[max_thrust_idx, '风速']
  479. max_thrust_value = valid_data.loc[max_thrust_idx, '气动推力_F_N']
  480. ax3.plot(max_thrust_wind, max_thrust_value, 'bo', markersize=8)
  481. # 4. 推力系数与功率关系
  482. ax4 = axes[1, 1]
  483. if len(valid_data) > 0:
  484. scatter = ax4.scatter(valid_data['风速'], valid_data['气动推力系数_Ct'],
  485. c=valid_data['有功功率'], cmap='viridis', s=50, alpha=0.7)
  486. ax4.set_xlabel('风速 (m/s)', fontsize=12)
  487. ax4.set_ylabel('气动推力系数 Ct', fontsize=12)
  488. ax4.set_title('推力系数与风速关系(颜色表示功率)', fontsize=14, fontweight='bold')
  489. ax4.grid(True, alpha=0.3)
  490. # 添加颜色条
  491. cbar = plt.colorbar(scatter, ax=ax4)
  492. cbar.set_label('功率 (kW)', fontsize=12)
  493. # 设置主标题
  494. fig.suptitle(f'{turbine_model} - {short_desc}\n分组: {group_id}',
  495. fontsize=16, fontweight='bold')
  496. # 调整布局
  497. plt.tight_layout()
  498. # 保存图像
  499. safe_group_id = group_id.replace('/', '_').replace('\\', '_').replace(':', '_')
  500. safe_group_id = re.sub(r'[<>:"|?*]', '_', safe_group_id) # 移除Windows文件名非法字符
  501. plot_path = Path(output_dir) / f"{safe_group_id}_推力系数曲线.png"
  502. plt.savefig(plot_path, dpi=300, bbox_inches='tight')
  503. plt.close(fig)
  504. plot_count += 1
  505. if i % 10 == 0 or i == total_groups:
  506. print(f" 进度: {i}/{total_groups} ({i/total_groups*100:.1f}%) - 已保存: {plot_path}")
  507. except Exception as e:
  508. failed_plots.append((group_id, str(e)))
  509. print(f" 警告: 分组 {group_id} 绘图失败: {e}")
  510. continue
  511. # 绘制汇总图
  512. self.plot_summary_analysis(output_dir)
  513. print(f"\n图表绘制完成!")
  514. print(f" 成功绘制: {plot_count}/{total_groups} 个分组")
  515. if failed_plots:
  516. print(f" 失败分组: {len(failed_plots)} 个")
  517. # 保存失败记录
  518. failed_path = Path(output_dir) / "绘图失败记录.txt"
  519. with open(failed_path, 'w', encoding='utf-8') as f:
  520. for group_id, error in failed_plots:
  521. f.write(f"{group_id}: {error}\n")
  522. print(f" 失败记录已保存到: {failed_path}")
  523. return True
  524. def plot_summary_analysis(self, output_dir: str) -> bool:
  525. """绘制汇总分析图"""
  526. if self.results is None:
  527. return False
  528. # 设置中文字体
  529. plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans']
  530. plt.rcParams['axes.unicode_minus'] = False
  531. # 计算每个分组的统计信息
  532. group_stats = []
  533. for group_id, group_data in self.results.groupby('分组标识'):
  534. valid_data = group_data[group_data['风速'] > 0]
  535. if len(valid_data) > 0:
  536. # 提取机型名称
  537. parts = group_id.split('_')
  538. turbine_model = parts[0] if parts else group_id
  539. group_stats.append({
  540. 'group_id': group_id,
  541. 'turbine_model': turbine_model,
  542. 'avg_ct': valid_data['气动推力系数_Ct'].mean(),
  543. 'max_ct': valid_data['气动推力系数_Ct'].max(),
  544. 'avg_thrust': valid_data['气动推力_F_N'].mean(),
  545. 'max_thrust': valid_data['气动推力_F_N'].max(),
  546. 'rated_power': group_data['有功功率'].max(),
  547. 'rotor_diameter': group_data['叶轮直径_D_m'].iloc[0] if '叶轮直径_D_m' in group_data.columns else 0,
  548. 'air_density': group_data['空气密度'].mean() if '空气密度' in group_data.columns else 0
  549. })
  550. if not group_stats:
  551. return False
  552. stats_df = pd.DataFrame(group_stats)
  553. # 创建汇总分析图
  554. fig, axes = plt.subplots(2, 2, figsize=(15, 12))
  555. # 1. 各机型平均推力系数比较
  556. ax1 = axes[0, 0]
  557. if 'turbine_model' in stats_df.columns and 'avg_ct' in stats_df.columns:
  558. model_avg = stats_df.groupby('turbine_model')['avg_ct'].mean().sort_values(ascending=False)
  559. bars = ax1.barh(range(len(model_avg)), model_avg.values)
  560. ax1.set_yticks(range(len(model_avg)))
  561. ax1.set_yticklabels(model_avg.index)
  562. ax1.set_xlabel('平均气动推力系数 Ct', fontsize=12)
  563. ax1.set_title('各机型平均推力系数比较', fontsize=14, fontweight='bold')
  564. ax1.grid(True, alpha=0.3, axis='x')
  565. # 在柱状图上显示数值
  566. for i, (bar, value) in enumerate(zip(bars, model_avg.values)):
  567. ax1.text(value, i, f' {value:.3f}', va='center', fontsize=9)
  568. # 2. 推力系数分布直方图
  569. ax2 = axes[0, 1]
  570. if 'avg_ct' in stats_df.columns:
  571. ax2.hist(stats_df['avg_ct'], bins=20, edgecolor='black', alpha=0.7, color='lightcoral')
  572. ax2.set_xlabel('平均气动推力系数 Ct', fontsize=12)
  573. ax2.set_ylabel('分组数量', fontsize=12)
  574. ax2.set_title('推力系数分布直方图', fontsize=14, fontweight='bold')
  575. ax2.grid(True, alpha=0.3)
  576. # 添加统计信息
  577. stats_text = f'总分组数: {len(stats_df)}\n'
  578. stats_text += f'平均值: {stats_df["avg_ct"].mean():.3f}\n'
  579. stats_text += f'中位数: {stats_df["avg_ct"].median():.3f}\n'
  580. stats_text += f'标准差: {stats_df["avg_ct"].std():.3f}\n'
  581. stats_text += f'范围: {stats_df["avg_ct"].min():.3f} - {stats_df["avg_ct"].max():.3f}'
  582. ax2.text(0.95, 0.95, stats_text, transform=ax2.transAxes,
  583. fontsize=10, verticalalignment='top', horizontalalignment='right',
  584. bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
  585. # 3. 推力系数与叶轮直径关系
  586. ax3 = axes[1, 0]
  587. if 'rotor_diameter' in stats_df.columns and 'avg_ct' in stats_df.columns:
  588. # 过滤掉直径为0的数据
  589. valid_diameters = stats_df[stats_df['rotor_diameter'] > 0]
  590. if len(valid_diameters) > 0:
  591. scatter = ax3.scatter(valid_diameters['rotor_diameter'], valid_diameters['avg_ct'],
  592. s=100, alpha=0.6, c=valid_diameters['rated_power'], cmap='viridis')
  593. ax3.set_xlabel('叶轮直径 (m)', fontsize=12)
  594. ax3.set_ylabel('平均推力系数 Ct', fontsize=12)
  595. ax3.set_title('推力系数与叶轮直径关系(颜色表示额定功率)', fontsize=14, fontweight='bold')
  596. ax3.grid(True, alpha=0.3)
  597. # 添加颜色条
  598. cbar = plt.colorbar(scatter, ax=ax3)
  599. cbar.set_label('额定功率 (kW)', fontsize=12)
  600. # 4. 推力系数与额定功率关系
  601. ax4 = axes[1, 1]
  602. if 'rated_power' in stats_df.columns and 'avg_ct' in stats_df.columns:
  603. scatter = ax4.scatter(stats_df['rated_power'], stats_df['avg_ct'],
  604. s=100, alpha=0.6, c=stats_df['avg_thrust'], cmap='coolwarm')
  605. ax4.set_xlabel('额定功率 (kW)', fontsize=12)
  606. ax4.set_ylabel('平均推力系数 Ct', fontsize=12)
  607. ax4.set_title('推力系数与额定功率关系(颜色表示平均推力)', fontsize=14, fontweight='bold')
  608. ax4.grid(True, alpha=0.3)
  609. # 添加颜色条
  610. cbar = plt.colorbar(scatter, ax=ax4)
  611. cbar.set_label('平均推力 (N)', fontsize=12)
  612. # 设置主标题
  613. fig.suptitle(f'气动推力系数汇总分析 (共 {len(stats_df)} 个分组)',
  614. fontsize=16, fontweight='bold')
  615. # 调整布局
  616. plt.tight_layout()
  617. # 保存图像
  618. summary_path = Path(output_dir) / "汇总分析图.png"
  619. plt.savefig(summary_path, dpi=300, bbox_inches='tight')
  620. plt.close(fig)
  621. print(f" 汇总分析图已保存到: {summary_path}")
  622. return True
  623. def run(self, output_dir: str = "output", plot_results: bool = True) -> bool:
  624. """
  625. 执行完整计算流程
  626. Args:
  627. output_dir: 输出目录
  628. plot_results: 是否绘制图表(为所有分组绘制图表)
  629. """
  630. print("="*60)
  631. print("气动推力系数计算器 - 按标准机型和描述分组")
  632. print("="*60)
  633. # 创建输出目录
  634. output_path = Path(output_dir)
  635. output_path.mkdir(parents=True, exist_ok=True)
  636. if not self.load_data():
  637. return False
  638. if not self.validate_data():
  639. return False
  640. if not self.calculate_thrust_coefficient():
  641. return False
  642. self.print_summary()
  643. # 保存计算结果
  644. print("\n保存计算结果...")
  645. results_path = output_path / "气动推力系数计算结果.csv"
  646. self.save_results(results_path)
  647. # 保存分组摘要
  648. summary_path = output_path / "分组统计摘要.csv"
  649. self.save_grouped_summary(summary_path)
  650. # 保存推力曲线数据
  651. curve_path = output_path / "推力曲线数据.csv"
  652. self.save_thrust_curve_data(curve_path)
  653. # 绘制图表(为所有分组绘制)
  654. if plot_results:
  655. print("\n开始为所有分组绘制气动推力系数曲线...")
  656. plots_dir = output_path / "推力系数曲线图"
  657. # 确认是否继续绘制(如果分组很多)
  658. if self.total_groups > 50:
  659. print(f"警告: 共有 {self.total_groups} 个分组,绘制所有图表可能需要较长时间!")
  660. response = input("是否继续绘制所有图表?(y/n): ")
  661. if response.lower() != 'y':
  662. print("跳过图表绘制...")
  663. else:
  664. self.plot_all_groups_thrust_curves(plots_dir)
  665. else:
  666. self.plot_all_groups_thrust_curves(plots_dir)
  667. # 显示前几行计算结果
  668. if self.results is not None:
  669. print("\n前5行计算结果示例:")
  670. display_cols = ['标准机型', '描述', '风速', '有功功率',
  671. '气动推力系数_Ct', '气动推力_F_N', '叶轮直径_D_m']
  672. existing_cols = [col for col in display_cols if col in self.results.columns]
  673. print(self.results[existing_cols].head())
  674. print("\n" + "="*60)
  675. print("计算流程完成!")
  676. print(f"结果文件保存在: {output_path}")
  677. print("="*60)
  678. return True
  679. def main():
  680. """主函数"""
  681. # CSV文件路径
  682. csv_file_path = f"./data/全部机型功率曲线_含标准类型_解析结果.csv" # 替换为您的CSV文件路径
  683. # 创建计算器实例
  684. calculator = AerodynamicThrustCalculator(csv_file_path)
  685. # 执行计算并绘制图表(为所有分组绘制图表)
  686. success = calculator.run(
  687. output_dir="output_results", # 输出目录
  688. plot_results=True # 为所有分组绘制图表
  689. )
  690. if success:
  691. print("\n所有计算和分析已完成!")
  692. # 显示计算结果摘要
  693. summary = calculator.get_summary()
  694. if summary:
  695. print("\n最终结果摘要:")
  696. for key, value in summary.items():
  697. print(f"{key:25}: {value}")
  698. else:
  699. print("\n计算失败!请检查错误信息。")
  700. if __name__ == "__main__":
  701. # 运行主程序
  702. main()