CN112949149A - 配气站设备工艺区管道地震易损性分析方法及分析系统 - Google Patents
配气站设备工艺区管道地震易损性分析方法及分析系统 Download PDFInfo
- Publication number
- CN112949149A CN112949149A CN202110533848.0A CN202110533848A CN112949149A CN 112949149 A CN112949149 A CN 112949149A CN 202110533848 A CN202110533848 A CN 202110533848A CN 112949149 A CN112949149 A CN 112949149A
- Authority
- CN
- China
- Prior art keywords
- pipeline
- earthquake
- seismic
- vulnerability
- curve
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/18—Network design, e.g. design based on topological or interconnect aspects of utility systems, piping, heating ventilation air conditioning [HVAC] or cabling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/04—Constraint-based CAD
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/14—Pipes
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Geometry (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Computer Networks & Wireless Communication (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明公开了配气站设备工艺区管道地震易损性分析方法及分析系统,建立配气站设备工艺区管道的有限元模型;对地震动进行调幅;进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值;绘制IDA曲线;定义管道在遭受地震作用下的各破坏状态及对应区间,并基于IDA曲线得到极限状态下地震动的强度指标;建立管道地震概率需求模型并求解,得到地震易损性函数,绘制地震易损性曲线。本发明用以解决现有技术中没有专门针对配气站设备工艺区管道的地震易损性分析方法的问题,实现为配气站设备工艺区管道提供针对性的地震易损性分析方法,为管理者提供科学制定震前防灾维护、震后应急管理策略方法的目的。
Description
技术领域
本发明涉及天然气配气站领域,具体涉及配气站设备工艺区管道地震易损性分析方法及分析系统。
背景技术
历年来,燃气供应系统因地震破坏而造成安全功能丧失的事例比比皆是:1994年的美国洛杉矶地震,致使燃气系统发生了高达15×104处漏气,导致数起火灾;1995年日本阪神地区发生大地震,地震造成煤气管道破裂致使煤气泄漏,共有459处起火,燃烧面积达数万平方米,造成大量的人员伤亡。2008年5月12日四川汶川发生里氏8.0级特大地震,都江堰市燃气管道破坏非常严重,地下管道破裂10余处,需重建的地下管道达50千米,全市燃气输配系统遭受经济损失约6700万元。因此,因地震导致的燃气供应系统安全事故会给人民的生命财产安全带来严重损害。在这种情况下,对燃气供应系统的震害风险进行科学、准确的评估,为燃气企业决策层和管理人员提供科学的评价方法是十分必要的。
设备工艺区在天然气配气站内承担着储存、调压、计量、加臭的主要任务。设备工艺区管道是指配气站站场内部连接各种设备、进行天然气输送处理的管道,与常规埋地管道不同的是,配气站设备工艺区管道在短距离(小范围)内连接了多种设备,如分离、调压、计量和各种阀门等,管道设置相较常规埋地管道而言十分复杂,因此其在遭受地震后的损害风险极大。然而现有技术中还没有专门针对配气站设备工艺区的地震易损性分析方法,无法为决策层和管理人员提供科学合理的设备工艺区防震害依据。
发明内容
本发明提供配气站设备工艺区管道地震易损性分析方法及分析系统,以解决现有技术中没有专门针对配气站设备工艺区管道的地震易损性分析方法、无法为决策层和管理人员提供科学合理的设备工艺区防震害依据的问题,实现为配气站设备工艺区管道提供具有针对性的地震易损性分析方法,为管理者提供科学制定震前防灾维护、震后应急管理策略方法的目的。
本发明通过下述技术方案实现:
配气站设备工艺区管道地震易损性分析方法,包括以下步骤:
S1、获取配气站设备工艺区管道的材料参数,建立管道几何模型;对几何模型划分网格、设置约束条件,得到有限元模型;
S2、对地震动进行调幅;
S3、将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值;以最大范式等效应力作为X轴,以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
S4、定义管道在遭受地震作用下的各破坏状态及对应区间,根据各破坏状态的对应区间得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
S5、建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
S6、根据得到的地震易损性函数,绘制配气站设备工艺区管道的地震易损性曲线。
针对现有技术中没有专门针对配气站设备工艺区的地震易损性分析方法、无法为决策层和管理人员提供科学合理的设备工艺区防震害依据的问题,本发明首先提出一种配气站设备工艺区管道地震易损性分析方法,本方法首先获取配气站设备工艺区管道的材料参数,材料参数可以根据已知数据直接获得,也可以采用力学试验方式实验得到,然后建立几何模型,并对几何模型划分网格、设置约束条件,得到有限元模型。其中网格划分方式需要综合考虑划分精度、CPU计算时间和数据存储空间的实际情况,约束条件根据配气站实际情况进行设置。本申请中,考虑到在对管道结构进行增量动力分析时,需要计算到管道结构的定义破坏点,但是当地震震级过小、释放的能量有限时,计算结果将很难使管道结构达到破坏,达不到预期的易损性分析结果,因此本申请需要对地震动进行调幅,再将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值,以最大范式等效应力作为管道损伤指标进行后续分析,并以最大范式等效应力和对应的地面峰值加速度PGA建立IDA曲线。然后定义管道在遭受地震作用下的状态,称之为破坏状态,并确定各破坏状态对应的区间,基于确定的区间,找出相邻两个区间的临界点,即为极限状态点,再根据IDA曲线得到极限状态下地震动的强度指标。最后建立管道地震概率需求模型,根据前述极限状态的计算结果,求解管道地震概率需求模型,得到地震易损性函数,绘制配气站设备工艺区管道的地震易损性曲线。本申请最终得到的地震易损性曲线,能够作为配气站设备工艺区震害分析的标准,克服了现有技术中没有专门针对配气站设备工艺区管道的地震易损性分析方法的缺陷,为配气站设备工艺区管道提供了具有针对性的地震易损性分析方法,为管理者制定震前防灾维护、震后应急管理策略提供了充分依据。
进一步的,所述材料参数包括弹性模量、拉伸屈服强度、塑性应变、泊松比和密度;所述几何模型包括管道本体、管道支架和阀门设备;所述管道本体、管道支架采用六面体网格进行划分,所述阀门设备采用四面体网格进行划分;所述约束条件包括:所述管道本体的两端设置为对称约束,所述管道支架设置为六个自由度完全约束。本方案中将管道材料参数分为弹性阶段的弹性模量、拉伸屈服强度,塑性阶段的塑性应变,以及基础参数泊松比和密度三部分,充分考虑了管道材料的物理特性,为后续建模过程做好充分准备。本方案中的几何模型包括管道本体、管道支架和阀门设备,虽然本申请的研究对象只是管道,但是对于配气站设备工艺区而言,其在小范围内设备十分密集,仅仅考虑管道本体难以满足基本的震害分析需求,故将与管道本体配套的管道支架和阀门设备纳入几何模型进行综合考量;其中,对于管道本体、管道支架,由于其结构相对简单,选择单元节点数适中、计算量中等的六面体网格进行划分,而对于结构相对复杂的阀门设备,若采用六面体网格划分容易导致质量较差,故采用四面体单元网格进行划分。此外,本方案在设置约束条件时,将管道支架四个底面的六个自由度完全约束,同时,将管道本体两端设置为对称约束,即约束管道两端x轴向的平动自由度和yz轴的转动自由度,从而准确模拟被分析管道在设备工艺区内的实际约束状态,有利于提高建模准确性。
进一步的,步骤S1对几何模型划分网格的过程中,设计不同数量的网格划分方案,对模型进行网格敏感性分析,基于网格敏感性分析结果,得出网格划分最佳方案。
网格划分的好坏将直接影响后续瞬态时程分析的精度,错误的网格划分方案甚至会导致模拟计算无法进行。作为本领域公知常识,网格划分越精细,网格数量越大,计算结果越接近实际情况。但实际操作时,却并非网格越多越好,网格数量越大,计算所需要的内存、CPU、硬盘等将会越大,过量的网格可能超出计算机性能极限,导致计算机死机,计算无法运行;同时,本方案发明人在大量研究过程中发现,由于计算误差始终存在,若非敏感区域内网格过多,反而可能导致误差累积过多而降低了计算精度。因此,为了确保所选择的网格数目及划分方案是最优或较优的,本方案在对几何模型划分网格的过程中,设计不同数量的网格划分方案,对模型进行网格敏感性分析,基于网格敏感性分析结果选择网格划分最佳方案。具体选择方法可包括:首先设计不同数量的网格划分方案;然后得出在各方案下模型的最大范式等效应力随时间变化曲线、最大等效应变随时间变化曲线、最大总变形随时间变化曲线;最后根据各曲线进行判断:首先排除在任一曲线内具有突变或明显异常的方案;然后对比剩余方案在各曲线内峰值出现的时间点,排除峰值出现时间点与其余方案具有明显差异的方案;最后在剩余的方案内,挑选网格数量更少、计算要求更低的网格划分方案作为最佳方案。
进一步的,对步骤S2地震动进行调幅的方法包括以下步骤:
S201、将已知的时程曲线加载到有限元模型中,得到第一次计算的最大范式等效应力、最大等效应变;并对配气站设备工艺区管道进行模态分析,得到一阶自振频率,将所述一阶自振频率转化为基本周期T1;
S202、将地震动记录转化为对应的反应谱曲线,找到在T1处的反应谱曲线加速度值,计算得到地震记录在基本周期T1处的所有反应谱的几何平均数Sa(T1) geomean :
式中,Sa(T1) i 为第i个反应谱在基本周期T1时的反应谱值,i取1,2,...n,n为反应谱总数;
S203、计算第i个反应谱的调幅因子SF i :SF i =Sa(T1) geomean /Sa(T1) i ;
S204、计算得到第i个反应谱调幅后的值Sa(T1) i+1:Sa(T1) i+1=SF i ·Sa(T1) i ;
S205、计算第i+1个反应谱的调幅因子SF i+1:SF i+1=ϛ i ·SF i ;式中,ϛ i 代表调幅放大系数,取值为ϛ i= ϛ1·(σ s / a'·σ i );σ s 为材料屈服强度,σ i 为第一次计算得到的最大范式等效应力,a'为放大精度系数;
S206、循环步骤S204~S205,直至完成对所有反应谱的调幅。
本申请中由于已知参数较少,传统调幅方法难以适用,为此本方案对地震动的调幅过程进行具体限定,其中首先将已知的时程曲线加载到有限元模型中进行第一次的时程分析,得到第一次计算的最大范式等效应力、最大等效应变;同时对配气站设备工艺区管道进行模态分析,得到一阶自振频率,将所述一阶自振频率转化为基本周期T1。配气站的设备工艺区管道可以看做是具有确定性周期的结构,因此本方案的调幅方法基于其基本周期T1实现,将所有的反应谱曲线作为一组,编号为1,2,3...n,i从1开始取值,依次计算每个反应谱的调幅因子,并计算得到每个反应谱调幅后的值后,以第1次反应谱的调幅因子作为依据,计算第2次反应谱的调幅因子,以此类推,直至对管道结构进行增量动力分析,计算到管道结构的定义破坏点时,完成对所有反应谱曲线的调幅即可。其中放大精度系数a'可以通过经验赋值或专家赋值方式得到。
进一步的,所述步骤S4包括以下步骤:
S401、将管道在遭受地震作用下的破坏状态定义为以下几种:基本完好、轻微破坏、中等破坏和严重破坏;各自对应的区间为:
当σ M <[σ r ],所述破坏状态为基本完好;当[σ r ]≤σ M ≤[σ e ],所述破坏状态为轻微破坏;当[σ e ]<σ M <[σ s ],所述破坏状态为中等破坏;当σ M >[σ s ],所述破坏状态为严重破坏;其中,[σ r ]=α1α2[σ s ],[σ e ]= α2[σ s ];其中,σ M 为管道在遭受地震作用下的范式等效应力;σ s 为材料屈服强度;α1、α2均为管道设计安全系数;
S402、取σ M =σ r ,σ M =σ e ,σ M =σ s 时,为极限状态点;
S403、将各极限状态点利用插值法代入IDA曲线中,得到各极限状态点对应的地面峰值加速度PGA,即为各极限状态下地震动的强度指标。
本方案限定了步骤S4的具体细化过程,其中首先给出了设备工艺区管道在遭受地震作用下的破坏状态的定义标准,其中以范式等效应力作为管道损伤指标,模拟管道在遭受调幅后各等级的地震动作用下所受到的范式等效应力,此时的范式等效应力处于不同区间,则对应管道处于不同的破坏状态。当管道处于相邻两个破坏状态的临界点时,即σ M =σ r ,σ M =σ e ,σ M =σ s 时,则为极限状态,利用插值法代入IDA曲线中即可得到各极限状态点对应的地面峰值加速度PGA,以地面峰值加速度PGA作为极限状态下地震动的强度指标。其中,管道设计安全系数可采用经验赋值或专家赋值或查阅相关规范方式得到,此不赘述。
式中,σ c ——定义的极限状态点所对应的管道损伤指标;σmax——地震作用下管道损伤指标;P f ——地震作用下管道损伤指标σmax超过定义的极限状态点所对应的管道损伤指标σ c 的条件概率;Ф( )——表示正态分布函数;α——待定系数;β——待定系数;PGA——地面峰值加速度;θ σc ——σ c 的对数标准差;θ σmax——σmax的对数标准差。
本方案所提供的管道地震概率需求模型为条件概率模型,其设计原理在于:管道的地震易损性曲线反应的是管道遭受不同强度地震作用时,管道损伤超过定义的管道损伤指标的条件概率。本方案中σ c 表示定义的极限状态点所对应的管道损伤指标,即各极限状态点对应的范式等效应力;σmax表示地震作用下管道损伤指标,即管道在遭受地震时所受的最大范式等效应力。此外,本模型中待定系数α、β可通过任意现有数学方法进行计算求解,凡是本领域能够通过任意现有技术可以实现的求解方法均可适用于本方案中,均属于本方案的保护范围之内。
进一步的,步骤S5中,求解管道地震概率需求模型的方法包括以下步骤:
确定待定系数α、β:
S501、设地面峰值加速度PGA和地震作用下管道损伤指标σmax服从指数关系:
σmax= α(PGA) β ;取a=lnα,b=β,得到:lnσmax=a+ b·ln(PGA);
S502、将地面峰值加速度PGA、地震作用下管道损伤指标σmax分别取对数,得到ln(PGA)和ln(σmax);将得到的ln(PGA)和ln(σmax)代入公式lnσmax=a+ b·ln(PGA)中进行回归分析,得到管道地震概率需求模型的回归曲线图;
S503、基于所述回归曲线图,得到a、b的取值;通过a、b的取值计算出待定系数α、β。
本方案给出了对待定系数α、β的具体计算方法之一,其中首先设定σmax= α(PGA) β ,其次引入了参数a、b,基于对数正态分布,将前述公式表达为如下指数关系:lnσmax=a+ b·ln(PGA),最后将PGA、σmax代入式中进行回归分析,基于回归曲线图即可求得a、b的取值,进而得到待定系数α、β的取值。
进一步的,步骤S5中,求解管道地震概率需求模型的方法还包括以下步骤:
配气站设备工艺区管道地震易损性分析系统,包括:
第一建模模块,用于获取配气站设备工艺区管道的材料参数,建立管道几何模型;并用于对所述几何模型划分网格、设置约束条件,得到有限元模型;
调幅模块,用于对地震动进行调幅处理;
IDA曲线模块,用于将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值,并以最大范式等效应力作为X轴、以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
极限状态计算模块,用于得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
第二建模模块,用于建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
输出模块,用于根据地震易损性函数,输出配气站设备工艺区管道的地震易损性曲线。
本发明与现有技术相比,具有如下的优点和有益效果:
1、本发明配气站设备工艺区管道地震易损性分析方法及分析系统,得到的地震易损性曲线能够作为配气站设备工艺区震害分析的标准,克服了现有技术中没有专门针对配气站设备工艺区管道的地震易损性分析方法的缺陷,为配气站设备工艺区管道提供了具有针对性的地震易损性分析方法,为管理者制定针对配气站设备工艺区内的管道的震前防灾维护、震后应急管理策略提供了充分依据。
2、本发明配气站设备工艺区管道地震易损性分析方法及分析系统,针对性的设置了配气站设备工艺区管道材料参数、几何模型、网格划分、约束条件等,显著提高了对设备工艺区管道建模的针对性,保证了后续对模型计算的准确可靠。
3、本发明配气站设备工艺区管道地震易损性分析方法及分析系统,给出了在已知参数较少的前提下,针对设备工艺区管道震害分析的地震动调幅方法,以及针对设备工艺区管道的管道地震概率需求模型及其求解方法,最终得到设备工艺区管道的地震易损性曲线。
附图说明
此处所说明的附图用来提供对本发明实施例的进一步理解,构成本申请的一部分,并不构成对本发明实施例的限定。在附图中:
图1为本发明具体实施例的流程示意图;
图2为本发明具体实施例的系统示意图;
图3为本发明具体实施例中的几何模型示意图;
图4为本发明具体实施例中得到的范式等效应力云图;
图5为本发明具体实施例中得到的等效应变云图;
图6为本发明具体实施例中得到的IDA曲线示意图;
图7为本发明具体实施例中地震易损性曲线示意图;
图8为本发明具体实施例中地震概率需求模型回归曲线图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明白,下面结合实施例和附图,对本发明作进一步的详细说明,本发明的示意性实施方式及其说明仅用于解释本发明,并不作为对本发明的限定。
实施例1:
如图1所示的配气站设备工艺区管道地震易损性分析方法,包括:
S1、获取配气站设备工艺区管道的材料参数,建立管道几何模型;对几何模型划分网格、设置约束条件,得到有限元模型;
S2、对地震动进行调幅;
S3、将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值;以最大范式等效应力作为X轴,以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
S4、定义管道在遭受地震作用下的各破坏状态及对应区间,根据各破坏状态的对应区间得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
S5、建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
S6、根据得到的地震易损性函数,绘制配气站设备工艺区管道的地震易损性曲线。
实施例2:
在实施例1的基础上,本实施例选择某天然气公司某配气站设备工艺区为研究对象,该配气站位于四川盆地北部,设计规模为4.0×104Nm3/d,调压前设计压力为4.0MPa,调压后设计压力为1.6MPa。该配气站于2014年建设完成,接上游配气站来气,进站压力为1.5MPa。工艺经过滤、调压、计量后供城镇居民用气;另一路输送至下游配气站。站内主要设备分为撬装设备以及放空系统。撬装设备包含:清管球接收装置撬、发球筒接收装置撬,调压计量撬、加臭装置;放空系统包含:回流阻火器、放空立管。
本实施例以该配气站设备工艺区内最具代表性的调压系统管道作为示例进行计算分析,其余管道同理:
第一步:
获取管道的材料参数:本实施例的材料主要为结构钢(支架用)和20号钢(管道本体用)两种,定义材料参数时分两个进行,首先是线弹性阶段,依据对20号钢管进行的力学试验,得到20号钢管的力学性能参数和本构关系,其管材的弹性模量为201.679GPa,拉伸屈服强度为392.6MPa。详细的线弹性阶段建模参数如表1和表2所示。
表1 20号钢线弹性阶段建模参数
材料名称 | 弹性模量(GPa) | 泊松比 | 密度(kg/m<sup>3</sup>) | σ<sub><i>s</i></sub>(MPa) |
20号钢 | 201.679 | 0.3 | 7850 | 392.6 |
表2 支架结构钢建模参数
材料名称 | 弹性模量(GPa) | 泊松比 | 密度(kg/m<sup>3</sup>) | σ<sub><i>s</i></sub>(MPa) |
结构钢 | 200 | 0.3 | 7850 | 250 |
各材料的塑性阶段的材料参数赋值根据力学试验得到的真应力-应变曲线进行转化,得到各向同性强化曲线,即可得到塑性应变。
建立管道几何模型:本实施例设备工艺区中调压系统管道采用“一用一备”设计,且管道布置、管径、壁厚、设备型号等条件都基本相同,故选择其中一根管路进行建模。其中,模型全长2412mm,变径前管道尺寸为89*4mm(外径*壁厚),变径后管道尺寸为57*4mm;管道支架高60mm,外壁与管道外壁相切。本实施例对模型进行了一定的简化假设,不考虑调压器、紧急切断阀和球阀的内部结构,将其简化为有一定质量的体积块,但体积块仍为实体单元模型。最终得到的几何模型如图3所示,其中图3为对建模结果经黑白处理后的示意图。
对几何模型划分网格:本模型采用实体单元solid186,它是一种高阶3D-20节点实体单元,由20个节点定义,每个节点有3个自由度,该单元支持塑性、超弹性、蠕变、大挠度和大应变运算。相较于六面体网格,在获得同等结果精度条件下,四面体网格需要更多的单元节点数,因而将耗费更长的CPU计算时间和更多的数据存储空间。与此同时,瞬态动力学分析需要均匀尺寸的网格,六面体网格仍然是首选。但是在本模型中,四个阀门设备由于结构形状较为复杂,六面体网格划分的质量较差,故采用四面体单元网格划分;模型其余部分均为六面体单元网格。
设置约束条件:将管道支架四个底面的六个自由度完全约束,同时,将管道本体两端设置为对称约束,即约束管道两端x轴向的平动自由度和yz轴的转动自由度。
第二步:
对地震动进行调幅:
首先将已知的时程曲线(1971年的San Fernando时程曲线)加载到有限元模型中,得到第一次计算的最大范式等效应力、最大等效应变;图4为得到的范式等效应力云图,从图4中可以看出,该调压系统的最大范式等效应力出现在管道内表面,最大值为20.3MPa,而地震动的峰值加速度为0.15g;图5为得到的等效应变云图,从图5中可以看出,最大等效应变位置与应变云图出现最大值的位置相同,均为管道内表面,最大应变为9.98×10-5。
需要说明的是,图4为本实施例通过软件计算得到的范式等效应力云图经黑白处理后的示意图;图5为本实施例通过软件计算得到的等效应变云图经黑白处理后的示意图。
调幅前,采用ANSYS workbench 对调压系统模型进行模态分析,得到调压系统的一阶自振频率为29.792Hz,即基本周期T1为0.034s。
将地震动记录转化为其对应的反应谱曲线,接下来找到在T1处的反应谱加速度值,计算得到20条反应谱在T1时刻的平均反应谱值;然后将每一条不同的地震动记录与相对应的调幅因子相乘,使每一条地震动记录在T1处的反应谱值与该组地震记录在T1处的平均值相等,便得到了调幅后的反应谱曲线。具体计算公式如下:
SF i =Sa(T1) geomean /Sa(T1) i ;
式中,Sa(T1) geomean 为一组地震记录在基本周期T1处的所有反应谱的几何平均数,SF i 为第i个反应谱的调幅因子;Sa(T1) i 为第i个反应谱在基本周期T1时的反应谱值,n为反应谱总数。
本实施例通过计算得到Sa(T1) geomean = 0.31g,同时得到本研究中一组地震动的调幅因子,由于地震动一共有三个方向,故每一个方向均有一个对应的调幅因子,将20组地震动的所有调幅因子总结如表3所示。
表3 首次地震动调幅对应的SF i 值
需要说明的是,表3中“时间”表示地震动发生年份,“事件”表示发生的哪次地震,“地震台站”表示检测到对应地震的地震台站的名称。其中,“事件”下方对应的英文为不同地震事件的命名,本领域内均用英文表示;“地震台站”下方对应的英文为监测到对应地震事件的地震动的地震台的名称,本领域内均用英文表示。
同时,本实施例选择变步长法作为地震动强度调幅的计算准则,调幅过程如下所示:
Sa(T1) i+1=SF i ·Sa(T1) i ;SF i+1=ϛ i ·SF;式中,ϛ i 代表调幅放大系数,取值为ϛ i= ϛ1·(σ s / a'·σ i );为材料屈服强度,σ s 为材料屈服强度,σ i 为第一次计算得到的最大范式等效应力,a'为放大精度系数。
以此类推,即可得到最终调幅结果。
第三步:将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值;以最大范式等效应力作为X轴,以地面峰值加速度PGA作为Y轴,绘制IDA曲线。
本实施例通过前述分析计算可得到在不同地震动下调压系统模型的IDA曲线。IDA曲线是将如表3中选择的7条地震动记录,按照上述方法进行调幅;然后将调幅后的地震动时程加载到管道结构上,进行瞬态时程分析,得到在任意的地震动下管道结构的最大范式等效应力响应值;最后将最大范式等效应力作为曲线X轴,地震动PGA作为曲线Y轴绘制IDA曲线。曲线各点的连接方式为B样条插值。绘制得到的IDA曲线如图6所示。由于IDA曲线离散型较大,故选择中位曲线(即图6中加粗的线条)进行后续计算。
需要说明的是,图6中的英文部分表示的是每条IDA曲线所对应的时间和事件和/或地震台站,对应至表3中的记载。其中Kern County、San Fernando、ImperialValley-06、Irpinia, Italy-01、Irpinia, Italy-02、Morgan Hill为不同地震事件的命名,本领域内均用英文表示;其中Bisaccia、Calitri、Gilroy Array #3为对应地震台的名称,本领域内均用英文表示。
第四步:
定义管道在遭受地震作用下的四种破坏状态及对应区间:
当σ M <[σ r ],所述破坏状态为基本完好;
当[σ r ]≤σ M ≤[σ e ],所述破坏状态为轻微破坏;
当[σ e ]<σ M <[σ s ],所述破坏状态为中等破坏;
当σ M >[σ s ],所述破坏状态为严重破坏;
其中,[σ r ]=α1α2[σ s ],[σ e ]= α2[σ s ];
其中,σ M 为管道在遭受地震作用下的范式等效应力;σ s 为材料屈服强度;α1、α2均为管道设计安全系数,根据地区类别进行取值;本实施例中参考《输气管道工程设计规范》(GB50251-2015),该配气站所处地区为Ⅱ类地区,所以α1、α2均取0.6。
取σ M =σ r ,σ M =σ e ,σ M =σ s 时,为IDA曲线中的极限状态点;在得到的IDA曲线中,需要找到三种极限状态对应的PGA的值,即确定地震动的强度指标。本实施例的极限状态点总结在表4中。根据IDA曲线利用插值法可以得到调压系统结构的极限状态下地震动强度指标如表5所示。
表4 极限状态描述
破坏状态 | 管道损坏程度 | 判断指标(MPa) |
基本完好 | 管道未损坏 | σ<sub><i>M</i></sub><141.336 |
轻微破坏 | 管道可能出现部分损坏,未发生泄漏 | 141.336≤σ<sub><i>M</i></sub>≤235.56 |
中等破坏 | 管道部分损坏,轻微泄漏,修理后仍可使用 | 235.56<σ<sub><i>M</i></sub>≤392.6 |
严重破坏 | 管道严重损坏,大量泄漏,需要及时更换 | σ<sub><i>M</i></sub>>392.6 |
表5 极限状态地震动强度指标
第五步:
建立管道地震概率需求模型:
式中,σ c ——定义的极限状态点所对应的管道损伤指标;σmax——地震作用下管道损伤指标;P f ——地震作用下管道损伤指标σmax超过定义的极限状态点所对应的管道损伤指标σ c 的条件概率;Ф( )——表示正态分布函数;α——待定系数;β——待定系数;PGA——地面峰值加速度;θ σc ——σ c 的对数标准差;θ σmax——σmax的对数标准差。
求解管道地震概率需求模型:
在上述管道地震概率需求模型中,存在α、β、三个待定系数。通过前述计算,我们已经得到了本调压系统在7条地震动下的IDA曲线,将计算得到的地震动强度指标(地面峰值加速度PGA)和调压系统管道的损伤指标σmax分别取对数,得到ln(PGA)和ln(σmax),然后将其代入公式lnσmax=a+ b·ln(PGA)中进行回归分析,其中a=lnα,b=β。公式lnσmax=a+ b·ln(PGA)是通过地面峰值加速度PGA和地震作用下管道损伤指标σmax的指数关系σmax= α(PGA) β 公式服从对数正态分布转换而来。
回归分析代入的数据如表6所示,表6中ln(PGA)为对地面峰值加速度PGA取对数的值,ln(σmax)为对损伤指标σmax取对数的值。通过回归分析,可以得到调压系统管道的地震概率需求模型回归曲线图如图8所示;由回归分析可得a=2.70891,b=0.52087,即可得到α=e a =15.0129,β= b=0.52087。
表6 回归分析参数
从图7中可以看出,作为本实施例的研究对象,该设备工艺区的调压系统管道在遭受不同峰值的地震动时,在三种极限状态下的超越概率:
(1)在PGA为0.5g时,调压系统管道超越基本完好的概率为0.68,超越中等破坏的概率为0.53,超越严重破坏的概率为0.41;
(2)在PGA为1.0g时,调压系统管道超越基本完好的概率为0.8,超越中等破坏的概率为0.67,超越严重破坏的概率为0.55;
(3)在PGA为2.0g时,调压系统管道超越基本完好的概率为0.89,超越中等破坏的概率为0.78,超越严重破坏的概率为0.7。
实施例3:
在实施例2的基础上,对几何模型划分网格的过程中,设计不同数量的网格划分方案,对模型进行网格敏感性分析,基于网格敏感性分析结果,得出网格划分最佳方案。具体过程示例如下:
根据本实施例所建立的有限元模型的实际情况和试算的结果设计出如表7所示五种网格划分方案,充分考虑管道和支架的网格精度,计算得到调压系统管道在各方案下的最大Von-mises应力随时间变化曲线、最大等效应变随时间变化曲线、最大总变形随时间变化曲线。
表7 网格划分方案
首先根据最大总变形随时间变化曲线排除最大总变形随时间的变化关系有突变的方案,再基于最大Von-mises应力随时间变化曲线、最大等效应变随时间变化曲线的峰值出现时间点是否相同判断模型计算结果的准确性,最后选择网格数量相对更少的方案作为网格划分理想最佳方案。
实施例4:
配气站设备工艺区管道地震易损性分析系统,如图2所示,包括:
第一建模模块,用于获取配气站设备工艺区管道的材料参数,建立管道几何模型;并用于对所述几何模型划分网格、设置约束条件,得到有限元模型;
调幅模块,用于对地震动进行调幅处理;
IDA曲线模块,用于将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值,并以最大范式等效应力作为X轴、以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
极限状态计算模块,用于得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
第二建模模块,用于建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
输出模块,用于根据地震易损性函数,输出配气站设备工艺区管道的地震易损性曲线。
实施例5:
一种计算机可读存储介质,包括计算机程序,所述计算机程序运行时,控制所述存储介质所在的设备执行如实施例1至3中任一记载的设备工艺区管道地震易损性分析方法。
以上所述的具体实施方式,对本发明的目的、技术方案和有益效果进行了进一步详细说明,所应理解的是,以上所述仅为本发明的具体实施方式而已,并不用于限定本发明的保护范围,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。
需要说明的是,在本文中,诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来,而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且,术语“包括”、“包含”或者其任何其它变体,意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、物品或者设备所固有的要素。
Claims (10)
1.配气站设备工艺区管道地震易损性分析方法,其特征在于,包括以下步骤:
S1、获取配气站设备工艺区管道的材料参数,建立管道几何模型;对几何模型划分网格、设置约束条件,得到有限元模型;
S2、对地震动进行调幅;
S3、将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值;以最大范式等效应力作为X轴,以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
S4、定义管道在遭受地震作用下的各破坏状态及对应区间,根据各破坏状态的对应区间得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
S5、建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
S6、根据得到的地震易损性函数,绘制配气站设备工艺区管道的地震易损性曲线。
2.根据权利要求1所述的配气站设备工艺区管道地震易损性分析方法,其特征在于,
所述材料参数包括弹性模量、拉伸屈服强度、塑性应变、泊松比和密度;
所述几何模型包括管道本体、管道支架和阀门设备;
所述管道本体、管道支架采用六面体网格进行划分,所述阀门设备采用四面体网格进行划分;
所述约束条件包括:所述管道本体的两端设置为对称约束,所述管道支架设置为六个自由度完全约束。
3.根据权利要求1所述的配气站设备工艺区管道地震易损性分析方法,其特征在于,步骤S1对几何模型划分网格的过程中,设计不同数量的网格划分方案,对模型进行网格敏感性分析,基于网格敏感性分析结果,得出网格划分最佳方案。
4.根据权利要求1所述的配气站设备工艺区管道地震易损性分析方法,其特征在于,步骤S2对地震动进行调幅的方法包括以下步骤:
S201、将已知的时程曲线加载到有限元模型中,得到第一次计算的最大范式等效应力、最大等效应变;并对配气站设备工艺区管道进行模态分析,得到一阶自振频率,将所述一阶自振频率转化为基本周期T1;
S202、将地震动记录转化为对应的反应谱曲线,找到在T1处的反应谱曲线加速度值,计算得到地震记录在基本周期T1处的所有反应谱的几何平均数Sa(T1) geomean :
式中,Sa(T1) i 为第i个反应谱在基本周期T1时的反应谱值,i取1,2,...n,n为反应谱总数;
S203、计算第i个反应谱的调幅因子SF i :SF i =Sa(T1) geomean /Sa(T1) i ;
S204、计算得到第i个反应谱调幅后的值Sa(T1) i+1:Sa(T1) i+1=SF i ·Sa(T1) i ;
S205、计算第i+1个反应谱的调幅因子SF i+1:SF i+1=ϛ i ·SF i ;式中,ϛ i 代表调幅放大系数,取值为ϛ i= ϛ1·(σ s / a'·σ i );σ s 为材料屈服强度,σ i 为第一次计算得到的最大范式等效应力,a'为放大精度系数;
S206、循环步骤S204~S205,直至完成对所有反应谱的调幅。
5.根据权利要求1所述的配气站设备工艺区管道地震易损性分析方法,其特征在于,所述步骤S4包括以下步骤:
S401、将管道在遭受地震作用下的破坏状态定义为以下几种:基本完好、轻微破坏、中等破坏和严重破坏;
当σ M <[σ r ],所述破坏状态为基本完好;
当[σ r ]≤σ M ≤[σ e ],所述破坏状态为轻微破坏;
当[σ e ]<σ M <[σ s ],所述破坏状态为中等破坏;
当σ M >[σ s ],所述破坏状态为严重破坏;
其中,[σ r ]=α1α2[σ s ],[σ e ]= α2[σ s ];
其中,σ M 为管道在遭受地震作用下的范式等效应力;σ s 为材料屈服强度;α1、α2均为管道设计安全系数;
S402、取σ M =σ r ,σ M =σ e ,σ M =σ s 时,为极限状态点;
S403、将各极限状态点利用插值法代入IDA曲线中,得到各极限状态点对应的地面峰值加速度PGA,即为各极限状态下地震动的强度指标。
7.根据权利要求6所述的配气站设备工艺区管道地震易损性分析方法,其特征在于,步骤S5中,求解管道地震概率需求模型的方法包括以下步骤:
确定待定系数α、β:
S501、设地面峰值加速度PGA和地震作用下管道损伤指标σmax服从指数关系:
σmax= α(PGA) β ;取a=lnα,b=β,得到:lnσmax=a+ b·ln(PGA);
S502、将地面峰值加速度PGA、地震作用下管道损伤指标σmax分别取对数,得到ln(PGA)和ln(σmax);将得到的ln(PGA)和ln(σmax)代入公式lnσmax=a+ b·ln(PGA)中进行回归分析,得到管道地震概率需求模型的回归曲线图;
S503、基于所述回归曲线图,得到a、b的取值;通过a、b的取值计算出待定系数α、β。
10.配气站设备工艺区管道地震易损性分析系统,其特征在于,包括:
第一建模模块,用于获取配气站设备工艺区管道的材料参数,建立管道几何模型;并用于对所述几何模型划分网格、设置约束条件,得到有限元模型;
调幅模块,用于对地震动进行调幅处理;
IDA曲线模块,用于将调幅后的地震动时程加载到管道的有限元模型上,进行瞬态时程分析,得到在任意地震动下管道的最大范式等效应力响应值,并以最大范式等效应力作为X轴、以地面峰值加速度PGA作为Y轴,绘制IDA曲线;
极限状态计算模块,用于得到极限状态点,基于IDA曲线得到极限状态下地震动的强度指标;
第二建模模块,用于建立管道地震概率需求模型,求解管道地震概率需求模型,得到地震易损性函数;
输出模块,用于根据地震易损性函数,输出配气站设备工艺区管道的地震易损性曲线。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110533848.0A CN112949149B (zh) | 2021-05-17 | 2021-05-17 | 配气站设备工艺区管道地震易损性分析方法及分析系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110533848.0A CN112949149B (zh) | 2021-05-17 | 2021-05-17 | 配气站设备工艺区管道地震易损性分析方法及分析系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112949149A true CN112949149A (zh) | 2021-06-11 |
CN112949149B CN112949149B (zh) | 2021-07-23 |
Family
ID=76233897
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110533848.0A Active CN112949149B (zh) | 2021-05-17 | 2021-05-17 | 配气站设备工艺区管道地震易损性分析方法及分析系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112949149B (zh) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116432335A (zh) * | 2023-03-03 | 2023-07-14 | 江阴市南方管件制造有限公司 | 一种基于有限元分析的化工用非标管件的制备工艺 |
CN117951921A (zh) * | 2024-03-26 | 2024-04-30 | 云南农业大学 | 一种边坡地震易损性评估方法 |
Citations (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106649954A (zh) * | 2016-10-08 | 2017-05-10 | 中冶华天工程技术有限公司 | 一种基于扩展pbee2理论框架的地震易损性分析方法 |
KR20180089598A (ko) * | 2017-01-31 | 2018-08-09 | 연세대학교 산학협력단 | 구조재 및 비구조재의 수리비용을 고려한 건축물 내진성능평가 구성 방법 |
CN110321653A (zh) * | 2019-07-11 | 2019-10-11 | 东北林业大学 | 一种考虑初始损伤状态的地震序列下结构易损性分析方法 |
CN111680349A (zh) * | 2020-06-03 | 2020-09-18 | 上海市建筑科学研究院有限公司 | 一种砌体结构平面外破坏地震易损性分析方法 |
CN112016857A (zh) * | 2020-11-02 | 2020-12-01 | 西南石油大学 | 基于云理论的聚乙烯管道地震易损性评估方法 |
CN112270125A (zh) * | 2020-10-23 | 2021-01-26 | 同济大学 | 一种基于有向图逻辑模型和蒙特卡洛模拟的系统地震易损性分析方法 |
CN112381475A (zh) * | 2021-01-15 | 2021-02-19 | 西南石油大学 | 一种配气站抗震安全评价方法及评价系统 |
CN112541286A (zh) * | 2020-11-23 | 2021-03-23 | 同济大学 | 一种基于增量动力分析法的隧道地震易损性分析方法 |
-
2021
- 2021-05-17 CN CN202110533848.0A patent/CN112949149B/zh active Active
Patent Citations (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106649954A (zh) * | 2016-10-08 | 2017-05-10 | 中冶华天工程技术有限公司 | 一种基于扩展pbee2理论框架的地震易损性分析方法 |
KR20180089598A (ko) * | 2017-01-31 | 2018-08-09 | 연세대학교 산학협력단 | 구조재 및 비구조재의 수리비용을 고려한 건축물 내진성능평가 구성 방법 |
CN110321653A (zh) * | 2019-07-11 | 2019-10-11 | 东北林业大学 | 一种考虑初始损伤状态的地震序列下结构易损性分析方法 |
CN111680349A (zh) * | 2020-06-03 | 2020-09-18 | 上海市建筑科学研究院有限公司 | 一种砌体结构平面外破坏地震易损性分析方法 |
CN112270125A (zh) * | 2020-10-23 | 2021-01-26 | 同济大学 | 一种基于有向图逻辑模型和蒙特卡洛模拟的系统地震易损性分析方法 |
CN112016857A (zh) * | 2020-11-02 | 2020-12-01 | 西南石油大学 | 基于云理论的聚乙烯管道地震易损性评估方法 |
CN112541286A (zh) * | 2020-11-23 | 2021-03-23 | 同济大学 | 一种基于增量动力分析法的隧道地震易损性分析方法 |
CN112381475A (zh) * | 2021-01-15 | 2021-02-19 | 西南石油大学 | 一种配气站抗震安全评价方法及评价系统 |
Non-Patent Citations (4)
Title |
---|
YING WU 等: "Seismic vulnerability analysis of buried polyethylene pipeline based on finite element method", 《INTERNATIONAL JOURNAL OF PRESSURE VESSELS AND PIPING》 * |
刘洪波 等: "框架结构的输入地震动记录调幅方法试验研究", 《地震工程与工程振动》 * |
宁成千 等: "埋地燃气管道易损性的模糊综合评价", 《内蒙古石油化工》 * |
李孟达: "框架—剪力墙结构地震易损性分析与抗震性能评估", 《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》 * |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116432335A (zh) * | 2023-03-03 | 2023-07-14 | 江阴市南方管件制造有限公司 | 一种基于有限元分析的化工用非标管件的制备工艺 |
CN117951921A (zh) * | 2024-03-26 | 2024-04-30 | 云南农业大学 | 一种边坡地震易损性评估方法 |
CN117951921B (zh) * | 2024-03-26 | 2024-06-04 | 云南农业大学 | 一种边坡地震易损性评估方法 |
Also Published As
Publication number | Publication date |
---|---|
CN112949149B (zh) | 2021-07-23 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN112949149B (zh) | 配气站设备工艺区管道地震易损性分析方法及分析系统 | |
CN103292762B (zh) | 判别大坝稳定性的位移监测方法 | |
CN101071098A (zh) | 地下钢质燃气管网管道腐蚀预测系统 | |
CN111914457B (zh) | 一种输电线路塔基边坡稳定性判别方法、装置和存储介质 | |
CN111027240A (zh) | 埋地管道安全评估方法及相关设备 | |
CN109635325A (zh) | 基于复合水动力与位移监测的水库型滑坡稳定性预测方法 | |
CN112989535A (zh) | 一种基于爆管检测效益的供水管网压力监测点优化布局方法 | |
CN113536646B (zh) | 一种单层球壳地震失效荷载计算方法 | |
CN106777488B (zh) | 一种桥梁安全性评估方法和系统 | |
CN111414692B (zh) | 基于贝叶斯修正模型的压力表校验台可靠性评估方法 | |
Thöns | Monitoring based condition assessment of offshore wind turbine support structures | |
CN106342313B (zh) | 一种基于层次分析法的测试性指标分配方法 | |
CN117494950B (zh) | 一种光储充检微电网一体站运行安全评价方法 | |
CN113836625A (zh) | 一种基于能力谱法的输电塔抗震性能等级划分方法及系统 | |
CN108805192B (zh) | 基于分层网络结构的监测数据分析方法 | |
CN103064999A (zh) | 一种用于抽水蓄能电站地下厂房结构的模型修正方法 | |
CN113609721A (zh) | 多类型极端灾害电气互联系统韧性计算方法及装置 | |
CN116861762A (zh) | 一种基于贝叶斯更新和水致劣化模型的岸坡变形预警方法 | |
CN116776435A (zh) | 一种盾构隧道纵向多级损伤指标阈值确定方法 | |
CN111341396A (zh) | 一种大气环境对材料腐蚀安全评估方法及系统 | |
Li et al. | Reliability and sensitivity analyses of monopile supported offshore wind turbines based on probability density evolution method with pre-screening of controlling parameters | |
KR102521589B1 (ko) | 상수도 관망 및 관로 설계 방법 및 장치 | |
CN114565257A (zh) | 油气并行管道风险评价方法和管控方法 | |
CN114880904A (zh) | 一种橡胶衬套的大规模有限元分析方法 | |
Xu et al. | Incremental dynamic analysis of SRC frame-bent structures in CAP1400 NPP |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |