基于H_理论的结构鲁棒性分析.pdf

返回 相似 举报
基于H_理论的结构鲁棒性分析.pdf_第1页
第1页 / 共6页
基于H_理论的结构鲁棒性分析.pdf_第2页
第2页 / 共6页
基于H_理论的结构鲁棒性分析.pdf_第3页
第3页 / 共6页
亲,该文档总共6页,到这儿已超出免费预览范围,如果喜欢就下载吧!
资源描述:
建筑结构学报Journal of Building Structures 第 33 卷 第 5 期 2012 年 5 月 Vol. 33No. 5May 2012 012 文章编号 1000-6869 2012 05-0087-06 基于 H∞理论的结构鲁棒性分析 张成 1,吴 慧 2,高博青1,董石麟1 1. 浙江大学 建筑工程学院,浙江杭州 310058; 2. 浙江财经学院 工商管理学院,浙江杭州 310018 摘要 基于 H∞理论, 提出一种定量评价结构鲁棒性的新方法。采用状态空间模型描述结构系统, 基于 H ∞最优, 采用系统传 递函数的 H∞范数作为结构鲁棒性的定量评价指标。对线性系统, 分离刚度矩阵的不确定性, 给出了 H ∞结构鲁棒性指标的 计算方法; 对非线性系统, 引入 L2性能准则表达鲁棒性。通过单自由度体系和桁架结构, 明确了鲁棒性指标的物理意义, 分 析影响结构鲁棒性的因素。结果表明 基于 H∞理论的结构鲁棒性指标代表了结构稳态振动反应的最大振幅与干扰幅值之 比, 可以反映外部干扰和结构内部的不确定性与结构的响应是否成比例; 提高结构整体承载力储备、 耗能能力以及关键构 件的冗余度可以增强结构的鲁棒性; 且 H∞鲁棒性指标对结构参数变化较为敏感。 关键词 鲁棒性; H∞控制;承载力储备;冗余度;损伤 中图分类号 TU323. 4TU352文献标志码 A H∞-based structural robustness analysis ZHANG Cheng1, WU Hui2,GAO Boqing1, DONG Shilin1 1. College of Civil Engineering and Architecture,Zhejiang University,Hangzhou 310058,China; 2. Business Adiministration College,Zhejiang University of Finance & Economics,Hangzhou 310018,China AbstractBased on the robust H∞control theory,a framework for the quantitative assessment of structural robustness is proposed. Structures are described by state space model first. Then following the idea of the H∞optimization, structural robustness is assessed by H∞norm of the system transfer function. For linear systems,the uncertainty of stiffness matrix is separated,and the calculation of the H∞structural robustness index is given. For nonlinear systems,the L2perance criterion is introduced to express structural robustness. Through analysis of a single degree freedom system and a truss structure,properties of the H∞structural robustness index are clarified,and factors affecting structural robustness are studied. The results show that,the H∞structural robustness index is the ratio of the maximum steady- state vibration amplitude and the interference amplitude,which indicates whether the structural response is proportional to the internal and external uncertainty interference. Improvements of both the overall structure perance like load- carrying capacity and energy dissipation capacity reservations,and the local redundancy of key components can enhance structural robustness. Meanwhile the H∞structural robustness index is sensitive to structural parameters. Keywordsrobustness; H∞control;load- carrying capacity reservation;redundancy;damage 基金项目 国家自然科学基金项目 51178414 , 浙江省自然科学基金项目 Y1110438 。 作者简介 张成 1986 , 男, 河南郑州人, 博士研究生。E- mail c- zhang zju. edu. cn 通讯作者 高博青 1963 , 男, 浙江慈溪人, 教授。E- mail bqgao zju. edu. cn 收稿日期 2011 年 6 月 78 0引言 鲁棒性 robustness 概念于上世纪 60 年代提出 后, 已在系统控制等领域得到广泛应用。2001 年美 国世贸中心大楼发生连续倒塌后, 建筑结构的鲁棒 性开始受到重视, 并且成为保证结构在不确定环境 下正常使用的重要性能指标 [1- 2 ], 我国结构设计规范 中的 “强节点弱构件” 、 “设置多道防线” 等概念也蕴 含了鲁棒性的思想。目前, 将结构鲁棒性定义为结 构不产生与起因不相称 disproportionate破坏能 力 [3 ], 并从抗倒塌性能[4 ]、 抵抗意外局部破坏能 力 [5 ]、 参数敏感性[6 ]等多种角度, 对结构鲁棒性做出 解释。但是这些定性的描述还无法为工程应用提供 统一的标准, 需要明确的定量评价指标, 评估不同结 构设计方案的优劣, 进而进行鲁棒性设计。 根据已有的结构性能指标可以导出结构的鲁棒 性, 例如文献[ 7- 8] 分别从承载力储备和能量吸收的 角度计算鲁棒性。考虑不确定性因素, 文献[ 9] 采用 直接风险与总风险的比值计算鲁棒性; 文献[ 10]则 将参数的不确定性归为模糊性, 用结构输入和输出 的香农信息熵的比值表示结构的鲁棒性, 但上述指 标不能完全揭示结构鲁棒性的本质。 20 世纪 80 年代初提出的以 H∞范数为性能指标 的 H∞控制理论是目前控制领域中比较成功且比较 完善的理论体系, 许多鲁棒性能准则均可以用 H∞范 数约束条件来描述 [11- 12 ]。本文基于 H ∞ 理论, 从结构 鲁棒性的定义出发, 提出一种定量评估结构鲁棒性 的新方法, 采用结构系统传递函数的 H∞范数作为结 构鲁棒性的定量评价指标, 依据单自由度体系说明 其物理意义, 通过桁架结构分析鲁棒性与结构承载 力储备、 冗余度等性能的关系, 验证该评价方法的有 效性。 1H∞基本理论 1. 1线性不确定系统的数学描述 线性时不变结构系统的状态和输出方程为 η t A0η t B0u t y t C0η t D0u t { 1 其中η t 、 u t 、y t分别表示结构系统的状态向 量、 输入向量和输出向量; A0、 B0、 C0、 D0是具有适当 维数的常数系统矩阵。 对式 1 进行拉普拉斯变换, 可以得到传递函数 矩阵 G s C0 sI - A0 -1B 0 D0 2 其中s 为 复 频 率。根 据 叠 加 原 理, 有 Y s G s U s , 其中 U s 、 Y s分别是 u t 、y t的 拉普拉斯变换。从数学角度, 结构是输入函数空间 U s到输出函数空间 Y s的一个映射。 结构往往受到外部不确定荷载 w t的干扰, 实 际 状 态 方 程 应 为 η t A0η t B1w t B2u t ; 并存在由几何误差、 材料参数误差等引起的 结构内部不确定性, 即实际系统矩阵存在不确定性 ΔA、 ΔB2、 ΔC、 ΔD。 因此实际结构系统可以描述为 η t A0 ΔAη t B1w t B2 ΔB2 u t y t C0 ΔCη t D0 ΔD u t { 3 1. 2H∞最优问题 常规设计只有当结构系统由式 1 精确描述时 才是最优设计, 而 H∞理论的基本思想是将结构系统 看作一个包含不确定性的结构族, 通过设计使结构 族中的所有对象 包括实际结构 均满足期望的性能 指标 [12 ]。 设从外部干扰 w t到由于干扰产生的系统输 出 Δy t的传递函数为 GwΔy, H∞最优问题是使 GwΔy 的 H∞范数极小, 即 min‖G wΔy‖∞, 从而保证结构族 中最劣的性能指标为最优, 使系统的干扰抑制性能 满足要求。其中 H∞范数定义为 ‖G s ‖∞sup ω∈ [0, ∞ σmax[ G jω ] 4 其中 sup 表示取上确界; σmax为矩阵的最大奇异值, j 为虚数单位; ω 为实变量。 2结构鲁棒性的定量评价方法 结构鲁棒性可以通过计算 “破坏后果” 与 “起因” 的比例得到, 即由不确定性干扰引起的结构输出响 应与输入干扰的比值, 考虑其中的最不利情况, 即为 结构系统传递函数的 H∞范数。 2. 1线性结构系统的鲁棒性 在理想弹性材料和小变形的假设下, 具有 n 个自 由度的结构运动方程为 M x t C x t Kx t u t 5 其中 M、 C、 K 分别为结构的质量、 阻尼和刚度矩阵; x t 、 u t分别为节点的位移和荷载向量。令状态 向量η t x t  x t[], 输入向量 u u t , 输出向量 y x t , 则系统的状态方程为 η t A0η t B0u t 6 其中 A0 0I - M -1K - M -1 [] C , B0 0 M - [] 1。 88 考虑外部不确定荷载 w t的干扰, 设输入 u t 到输出 y t的传递函数为 Guy, 干扰 w t到干扰产 生的输出 Δyd t的传递函数为 GwΔy d,即 y t Δy d t[] Guy0 0GwΔy [] d u t w t[]。 对线性系统, Guy 与 GwΔy d是 等价的。因此结构的鲁棒性可以表示成 IRd ‖GwΔy d‖∞ ‖Guy‖∞ 7 目前, 没有统一的方法可以同时考虑质量、 刚度 和阻尼的不确定性矩阵进行合理分解 [11 ], 因此假设 不考虑结构的质量和阻尼矩阵的不确定性, 而刚度 矩阵 K 存在范数有界的不确定性, 其可能的最大变 异为 δK δ ≥ 0 , 则式 5 相应地修改为 M x t C x t K ΔK x t u t 8 其中ΔK ΣKδK,ΣK为 不 确 定 实 标 量, 满 足 ‖ΣK‖ ≤1。 将式 8 写成状态方程的形式为 η t A0 ΔAη t B0 ΔB u t 9 其中 ΔA 00 - δΣKM -1K [] 0 ,ΔB 0。为将 ΔA、 ΔB 中的不确定性分离, 可以做变换 [ ΔAΔB] EΣ[ FaFb] 10 其中 E、 Fa、 Fb是具有适当维数的已知矩阵,Σ是未 知函数阵, 且满足Σ∈ Ω {Σ Σ TΣ ≤ I } 。令 E 00 - M -1K [] 0 [13 ],Σ ΣKInn0 [] 00 , 则 Fa δI nn 0 [] 00 ,Fb 0。通过定义增广系统 GuΔy m A0EB0 Fa0Fb I           00 [14 ] , 可以计算参数不确定性对结构鲁 棒性的影响为 IRm ‖GuΔy m‖∞ δ ‖Guy‖∞ 11 则结构系统的鲁棒性为 IR IRd IRm 12 2. 2非线性结构系统的鲁棒性 非线性结构系统的传递函数矩阵与输入有关, 因此不能根据式 12 直接计算鲁棒性指标。引入 L2 性能准则 [12 ], 则鲁棒性可以用结构系统的 L 2诱导范 数表示 IR ‖GwΔy‖∞sup ‖w‖2 ≠0 ‖Δy t ‖2 ‖w t ‖2 13 其中 ‖ ‖2 ∫ ∞ 0 ‖ ‖ 2d t 1/2 ∫ ∞ 0 TQ dt 1/2, Q 为加权矩阵。 2. 3H∞鲁棒性评价指标的性质 将外部荷载和内部参数的不确定性因素统一用 珚w 表示, 则 H∞鲁棒性评价指标 IR有以下性质 1 鲁棒性的非负性, 即 IR≥ 0珚w 2 在等量的不确定性干扰下, 结构性能变化越 小, IR越小, 结构的鲁棒性越好, 即 Δy 1< Δy2 I R1< IR2珚 w1 珚w2 3 如果在微小干扰下, 结构仍然存在非 0 的输 出响应, 则 IR趋于∞, 结构具有最弱的鲁棒性, 即 Δy> 0IR→ ∞ 珚w→0 4 如果在非 0 的干扰下, 结构性能只有微小变 化, 则 IR趋于 0, 结构具有最强的鲁棒性, 即 Δy →0IR→0珚w > 0 3H∞结构鲁棒性评价指标的物理 意义 为明晰结构鲁棒性的物理意义, 以单自由度体 系为例, 其运动方程为 m x t c x t k0 x t u t 14 其中 m、 c、 k0分别为结构的质量、 阻尼和刚度, x t 、 u t分别为节点位移和节点荷载。 令状态向量 x t [ x t  x t ] T,输出向量 y t x t , 输入向量 u t u t , 则式 14 表达 的单自由度体系的传递函数为 G s 1 ms2 cs k0 15 代入式 7 , 有 IRd ‖GwΔy d‖∞ sup ω∈ [0, ∞ 1 k0- mω2 jωc 16 由于实际结构的阻尼比一般较小, 假设 c2< 2mk 即 ξ < 1/槡2 , 则当 ω 2mk0- c2 2m 槡 2 时, 有 IRd 1 c k0 m - c2 4m 槡 2 17 进一步考虑结构刚度的不确定性, 设实际刚度 k 1 Σk δ k 0, 此时有 IRm ‖GuΔy m‖∞ sup ω∈ [0, ∞ δ k0- mω2 jωc 18 同 理,当 ω 2mk0- c2 2m 槡 2 时, IRm δ/ c k0 m - c2 4m 槡 2 。 98 因此, 单自由度体系的鲁棒性为 IR IRd IRm 1 δ c k0 m - c2 4m 槡 2 19 当外部激励为 u u0sinωt 的简谐荷载时, ω 2mk - c2 / 2m2 槡 为共振频率; 由于对传递函数 取上确界, 考虑存在不确定性因素的最不利情况, 因 此, 鲁棒性指标 IRd数值上等于干扰产生的稳态振动 反应的最大振幅与干扰幅值之比, 反映了结构的响 应与干扰的关系; IRm则是考虑刚度不确定性的最不 利情况, 即刚度下降 δ 时系统的振幅。如果不考虑阻 尼, 即取 c 0, 则式中 IR→∞, 意味着即使荷载幅值 很小, 如果其频率与结构自振频率相同, 也会导致振 动发散即结构失效, 此时其鲁棒性很弱。 4影响结构鲁棒性因素 考虑结构的承载力储备、 冗余度、 耗能能力及构 件损伤等影响结构的鲁棒性因素, 以如图 1 所示的 平面十杆桁架结构为例, 定量分析其关联性。 以如图 1 所示的平面桁架结构为例, 令构件弹 性模量 E 2. 06 105MPa, 截面积 A 800 mm2 , 长 度 l 2 m, 密度 ρ 7 850 kg/m3, 采用 Rayleigh 阻尼, 阻尼比 ξ 0. 02, 取前两阶自振频率计算结构阻尼。 计算结构的鲁棒性时, 首先在 ANSYS 中建立结构的 有限元模型, 并提取刚度和质量矩阵, 进而利用 Riccati 方程, 在 MATLAB 中编程计算系统的 H∞范 数, 可以得到桁架结构的鲁棒性 IR06. 428 mm/kN。 图 1平面十杆桁架结构 Fig. 1Planar 10- bar truss 考察刚度不确定性对鲁棒性的影响, 采用式 12 计算结构的鲁棒性 IR,结果如图 2 所示。由图 可见, IR随 δ 增大呈线性变化, 当 δ 0 时, I R IR0。 假设弹性模量在[ E 1 - δ , E 1 δ ]内均匀 分布, 随机生成 E 后, 结构只存在外部不确定性, 通 过式 7 计算鲁棒性指标。取 δ 0. 1、 0. 2、 0. 3、 0. 4, 分别进行 1 000 次模拟, 结果如图 2 所示。可见采用 式 12 计算刚度不确定结构的鲁棒性, 能够实现对 最不利情况的包络, 具有一定安全储备。 图 2刚度不确定性对鲁棒性的影响 Fig. 2Effect of stiffness uncertainty on structural robustness 4. 1承载力储备 通常认为提高结构的承载力储备可以使结构在 不确定性干扰下具有一定的裕度, 从而避免结构整 体或连续倒塌 [4, 7 ]。在特定分布形式的荷载作用下, 结构的鲁棒性可采用式 13 计算。将输入节点荷载 写做 u umu0, 其中 um为荷载幅值, u0为归一化的 最不利荷载分布向量, 则输出节点位移 y umLu0, 其 中 L 为结构的柔度矩阵。如果输入向量和输出向量 的权重矩阵均取单位矩阵, 则结构的鲁棒性为 IR 1 δ sup ‖u‖2 ≠0 ‖y t ‖2 ‖u t ‖2 1 δ uT 0L TLu 0 uT 0u 0 1/2 20 杆件内力向量 f 和节点荷载向量 u 满足关系 u Aeqf 21 其中Aeq是与节点坐标和杆件长度有关的平衡矩 阵。设杆件破坏应力为 σcr, 结构到达弹性承载力 ucr, 此时 u ucru0, f σcrAf0, 其中 f0为归一化的杆 件内力分布向量, 代入式 21 有 ucru0 AeqσcrAf0 22 进一步整理, 有 ucr uT 0Aeqf0 uT 0u0 σcrA C1σcrA 23 式中 C1 uT 0Aeqf0 uT 0u0 。 结构设计承载力可以表达为 ud C1σdA 24 其中 σd为杆件在设计荷载作用下的最大应力, σd≤σcr。 所以, 结构承载力储备为 uR C1 σ cr - σ d A 25 将式 25 代入式 20 , 可得鲁棒性指标 IR C1C2 1 δ σ cr - σ d l E 1 uR 26 其中 C1、 C2是与结构刚度和荷载分布有关的常数。 09 因此, 在相同的刚度和荷载作用下, 结构承载力储备 越大, 则 IR越小, 结构鲁棒性越好。 4. 2耗能能力 结构安全除了结构破坏时的承载能力, 还包括 极限变形能力, 两者之积即为吸收的能量 [4, 8 ]。结构 达到弹性极限状态时, 吸收能量为 W 1 2 uTy 1 2 u2 cru T 0Lu0 C2 1C3 σ 2 crAl 2E 27 其中 C3为与结构刚度和荷载分布有关的常数。将 式 27 代入式 20 , 则 IR的表达式为 IR C2 1C2C3 1 δ 2 σcrl E 2 1 W 28 因此, 结构耗能能力越强, 则 IR越小, 结构的鲁 棒性越好。 4. 3冗余度 图 1 中 的 桁 架 为 二 次 超 静 定 结 构 IR0 5. 672 mm/kN , 可以移除 1 至 2 根杆件, 以考察冗余 度对结构鲁棒性的影响, 结果见表 1。 表 1冗余度与结构鲁棒性关系 Table 1Effect of redundancy on structural robustness 模型 IR/ mm kN -1 模型 IR/ mm kN -1 ①移除 9 杆或 10 杆 8. 862 ②移除 7 杆或 8 杆 8. 732 ③移除 8、 9 杆或 7、 10 杆 11. 204 ④移除 7、 9 杆或 8、 10 杆 11. 692 ⑤移除 2 杆或 4 杆 9. 351 ⑥移除 1 杆或 3 杆 23. 751 ⑦移除 1、 4 杆或 2、 3 杆 26. 097 ⑧移除 1、 2 杆或 3、 4 杆 29. 580 可见移除构件减少结构的冗余度后, 结构的鲁 棒性 IR会有不同程度的削弱。移除 1 根斜腹杆后, IR较完整结构的 IR0分别增大 36 ~38; 移除 2 根 斜腹杆后, IR增大 74 ~82, 而成为静定结构后, 2 根斜腹杆采用对称布置 模型③ 或平行布置 模型 ④ , 对结构鲁棒性的影响不大。移除非支座处的 1 根弦杆后 模型⑤ , IR较 IR0增大 45, 表明弦杆较 斜腹杆对结构鲁棒性更为重要; 移除支座处的弦杆, IR会迅速增大至 IR0的 3. 7 倍 模型⑥ , 因为此时如 果在支座处再次发生单杆局部破坏, 将导致结构整 体倒塌, 即不相称的破坏。移除 2 根弦杆后 模型 ⑦、 ⑧ , 结构的鲁棒性将进一步下降。以上结果表 明, IR对结构冗余度的变化较敏感。 4. 4构件损伤 意外事件如爆炸、 设计施工失误等可能导致构 件产生损伤, 进而导致结构的不相称破坏。用构件 截面积减小表示构件损伤程度 D 1 - A/A, 其中 A 为损伤后构件截面积。计算构件损伤对结构鲁棒性 的影响, 结果如图 3 所示, 其中纵坐标表示鲁棒性指 标的相对变化。可见随着损伤程度的增大, 结构的 鲁棒性受到相应的削弱; 相比 9 杆单独损伤的情况, 8、 9 两根杆同时受到损伤时结构鲁棒性下降更为明 显; 鲁棒性对弦杆损伤更为敏感。 图 3构件损伤对结构鲁棒性影响 Fig. 3Effect of components damage on structural robustness 5结论 1 基于鲁棒 H∞控制理论, 从定义出发提出一 种定量评估结构鲁棒性的新方法。采用系统传递函 数的 H∞范数作为结构鲁棒性的定量评价指标, 可以 反映外部的干扰和结构内部的不确定性与结构的响 应是否成比例。 2H∞范数作为结构鲁棒性的定量评价指标, 可以揭示结构的承载力、 耗能、 冗余度等物理指标的 变化, 而且均很好地反映这些物理量与结构鲁棒性 的关联性。 3 研究表明, 冗余度本质上反映的是结构局部 性能, 如对算例中桁架, 增加腹杆对提高结构鲁棒性 的作用不大。 4 意外事件导致构件损伤, 会削弱结构的鲁棒 性。构件损伤的位置及程度不同, 鲁棒性指标会有 相应的变化。 参考文献 [ 1] EN 1991- 1- 7Eurocode 1actions on structures[ S] . Brussels EuropeanCommitteeforStandardization, 2006. 19 [ 2] ASCE/SE I 7- 10Minimum design loads for buildings and other structures[ S] . Reston ASCE, 2010. [ 3] Starossek U,Haberland M. Disproportionate collapse terminology and procedures[ J] . Journal of Perance of Constructed Facilities, 2010, 24 6 519- 528. [ 4] 叶列平,程光煜,陆新征,等. 论结构抗震的鲁棒性 [ J] . 建筑结构, 2008, 38 6 11- 15. Ye Lieping, Cheng Guangyu,Lu Xinzheng,et al. Introduction of robustnessforseismicstructures[ J] .Building Structure, 2008, 38 6 11- 15. in Chinese [ 5] England J,Agarwal J,Blockley D. The vulnerability of structurestounforeseenevents[ J] .Computers& Structures, 2008, 86 10 1042- 1051. [ 6] Lee M,Kelly D W,Degenhardt R,et al. A study on the robustness of two stiffened composite fuselage panels [ J] . Composite Structures, 2010, 92 2 223- 232. [ 7] Yan D,Chang C C. Vulnerability assessment of single- pylon cable- stayed bridges using plastic limit analysis [ J] .Engineering Structures,2010,32 8 2049- 2056. [ 8] 方召欣,李惠强. 基于能量观点的结构安全性与鲁 棒性[J] . 建筑结构学报,2007,28 增刊 269- 273. FANGZhaoxin, LIHuiqiang.Safetyand robustness of structures from the viewpoint of energy [J] .JournalofBuildingStructures, 2007, 28 Suppl. 269- 273. in Chinese [ 9] Baker J W, Schubert M, FaberMH.Onthe assessment of robustness[J] . Structural Safety,2008, 30 3 253- 267. [ 10] Beer M,Liebscher M. Designing robust structuresa nonlinear simulation based approach[J] .Computers and Structures, 2008, 86 10 1102- 1122. [ 11] 宁响亮,刘红军,谭平,等. 基于 LMI 的结构振动多 目标鲁棒 H2/H∞控制[ J] . 振动工程学报, 2010, 23 2 167- 172. NING Xiangliang,LIU Hongjun, TAN Ping,et al. Multi- objective robust H2/H∞control for structural vibration based on LMI[J] . Journal of Vibration Engineering,2010,23 2 167- 172. in Chinese [ 12] 梅生伟,申铁龙,刘康志. 现代鲁棒控制理论与应 用[M] . 北京清华大学出版社,200861-63. MEI Shengwei,SHEN Tielong,LIU Kangzhi. Modern robust control theory and application[M] .Beijing Tsinghua University Press, 2008 61-63. in Chinese [ 13] 宋刚,吴志刚,林家浩. 考虑参数不确定性的结构 抗震鲁棒 H∞控制[J] . 地震工程与工程振动, 2008, 28 6 226- 232. SONG Gang,WU Zhigang, LIN Jiahao. Robust H∞control for earthquake- excited structures with parameter uncertainties[J] . Journal of Earthquake Engineering and Engineering Vibration, 2008, 28 6 226- 232. in Chinese [ 14] 王德进. H2和 H∞优化控制理论[ M] . 哈尔滨 哈尔 滨工业大学出版社,2001128- 132 WANG Dejin. H2and H∞optimal control theory[ M] . HarbinHarbin Institute of Technology Press,2001128- 132. in Chinese 29
展开阅读全文

资源标签

最新标签

长按识别或保存二维码,关注学链未来公众号

copyright@ 2019-2020“矿业文库”网

矿业文库合伙人QQ群 30735420