CN103303495A - 一种动力下降过程干扰力矩的估计方法 - Google Patents
一种动力下降过程干扰力矩的估计方法 Download PDFInfo
- Publication number
- CN103303495A CN103303495A CN2013101251835A CN201310125183A CN103303495A CN 103303495 A CN103303495 A CN 103303495A CN 2013101251835 A CN2013101251835 A CN 2013101251835A CN 201310125183 A CN201310125183 A CN 201310125183A CN 103303495 A CN103303495 A CN 103303495A
- Authority
- CN
- China
- Prior art keywords
- lander
- disturbance torque
- mass
- attitude
- estimation
- 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
Landscapes
- Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)
Abstract
一种动力下降过程干扰力矩的估计方法。在月球软着陆动力下降过程中,因为质心偏移产生的干扰力矩有可能会危及着陆的安全。在线估计出干扰力矩的大小是进行预警和应对的重要措施。首先利用制导和控制系统发出的发动机推力大小和脉宽指令,递推计算着陆器剩余质量,并使用地面装订的曲线(函数)计算着陆器质心高度和惯量大小;然后根据质心高度和脉宽指令计算控制力矩大小;接下来对陀螺输出进行差分获得角加速度计算值;之后,使用姿态动力学方程反算干扰力矩;最后使用滤波器降低干扰力矩估值中噪声的影响。本发明提供了一种监视着陆器干扰力矩变化的方法,有利于提高着陆过程控制系统的鲁棒性,降低着陆的风险。
Description
技术领域
本发明涉及一种月球软着陆动力下降过程中在线估计着陆器干扰力矩的方法,属于月球探测自主导航、制导和控制领域。
背景技术
月球软着陆过程主要依靠着陆器上安装的制动发动机产生制动力,在减小速度的同时降低高度来完成。对于不具备推力矢量控制能力的一类着陆器来说,由于存在制动发动机制造、安装精度的限制,以及下降过程推进剂消耗不平衡性等方面的问题,发动机推力不可能保持过质心,这样下降过程中制动推力就会产生作用于着陆器的干扰力矩。正常情况下,该干扰力矩需要由姿态控制系统进行抵消。但一旦干扰力矩超过姿态控制系统的能力,就会危及着陆过程的安全,造成姿态角速度过大(超过陀螺量程)或者姿态翻转,为此需要对下降过程干扰力矩的方向和大小进行估计,以便于在危险发生前提前采取对策。
美国的阿波罗飞船在着陆月球过程中采用的是摇摆发动机,可以通过调整推力方向来消除这种干扰力矩;而其他无人月球着陆器均采用的是固定推力方向的制动发动机,会在着陆过程中产生这种干扰力矩,但一般依靠具有足够能力的姿控系统来进行抵消。目前通过已公开的资料,尚未见到使用干扰力矩估计方法来进行着陆过程故障检测和处理的报道。
发明内容
本发明所解决的技术问题:针对月球软着陆动力下降过程由于制动发动机推力不过质心可能产生较大干扰力矩,危及着陆安全的问题,提供一种动力下降过程中着陆器干扰力矩的估计方法,提高着陆器着陆的安全性。
本发明所采用的技术方案为:一种动力下降过程中着陆器干扰力矩的估计方法,利用GNC姿控指令计算作用于着陆器的控制力矩大小,利用陀螺测量差分估计星体角加速度,依靠动力学模型解算干扰力矩大小,并引入滤波器消除角速度差分和推力器偏差造成的噪声,对干扰力矩估计进行平滑。实现步骤如下:
(1)着陆器质量特性估计
根据发动机比冲和推进剂消耗的关系,估算推进剂剩余质量;然后根据事先存储在星载计算机上的质心-剩余质量变化曲线以及惯量-剩余质量变化曲线,计算出质心的位置和惯量的大小。
(2)姿控力矩计算
上一个周期,控制计算机根据姿态和姿态角速度误差计算出了滚动、俯仰和偏航三个轴上的喷气脉宽;根据(1)中质心位置的计算结果,并结合着陆器上所配备的姿控推力器的安装位置和指向,计算出一个周期内着陆器姿控喷气所产生的平均控制力矩。
(3)角加速度计算
对上两个周期陀螺姿态角速度的测量结果进行差分,获得上一周期平均角加速度的测量值。
(4)利用动力学对干扰力矩的解算
利用(1)~(3)中计算出的探测器惯量、控制力矩、角速度、角加速度,通过姿态动力学方程反算干扰力矩的大小。
(5)干扰力矩估值的平滑
将干扰力矩在着陆器三个体轴的分量分别通过一个低通滤波器进行滤波,获得干扰力矩的估计值。
所述步骤(1)中估算推进剂剩余质量利用火箭公式估算着,用事先装订的剩余质量-质心高度和剩余质量-惯量曲线(函数)估算着陆器质量特性,避免从力学原理出发的复杂的质量特性计算;利用姿态动力学方程,使用陀螺测量和姿态控制脉宽指令作为输入,反算干扰力矩大小,并使用滤波技术降低干扰力矩估计中的噪声。
本发明与现有技术相比的优点在于:本发明利用着陆器装备的常用姿态敏感器,使用动力学方程,实时在线估计干扰力矩大小,提高着陆过程控制系统的鲁棒性,实现了对自身所受干扰力矩的估计,从而为监视着陆器受力情况,提前发现姿控系统可能面临的危险提供了技术手段,有利于降低月球着陆动力下降过程的风险,提高着陆的安全性。
附图说明
图1是本发明方法的实现流程图;
图2是着陆器质量特性计算的原理图;
图3是姿控喷气脉宽的示意图;
图4是着陆器质心高度随剩余质量的变化曲线;
图5是着陆器三个主惯量随剩余质量的变化曲线;
图6是动力下降过程实际干扰力矩的变化曲线;
图7是估计出的动力下降过程干扰力矩变化曲线。
具体实施方式
着陆器处于月球软着陆动力下降过程。着陆器上装备有陀螺,下降和减速由制动发动机完成,姿态采用喷气控制。
如图1所示,本发明具体计算过程如下:
(1)着陆器质量特性估计
着陆过程推进剂的消耗速度取决于制动发动机以及姿控发动机推力的大小和它们的比冲。假定当前时间为tk,上两周期(用时间段[tk-2,tk-1]表示)制导向制动发动机发出的推力大小指令为Fcmd,k-2;姿态控制向姿控发动机发出的各轴喷气脉宽为Ton,k-2,i(i=x,y,z表示着陆器的三个坐标轴),而各轴姿控发动机输出的推力大小为FRCS,i(i=x,y,z);另外设制动发动机的比冲为Isp,main,各轴姿控发动机的比冲为Isp,RCS,i,着陆器的质量为m,则质量变化率为
式(1)就是所谓的火箭公式。
设控制周期为T,上两控制周期初始时刻tk-2时的质量为mk-2,那么在上一周期tk-1时刻着陆器质量变为
根据式(2)就可以递推获得当前探测器的质量估计。
之后,根据事先存储在星载计算机上的剩余质量-质心高度曲线(用函数ro=f(m)表示),以及剩余质量-整器惯量曲线(用函数I=h(m)表示),求出tk-1时刻的质心位置和惯量
剩余质量-质心高度曲线和剩余质量-整器惯量曲线,应根据着陆器的结构设计中的若干参数,包括着陆器干器(不含推进剂)质量mb,干器质量中心相对于着陆器本体的位置rb,干器惯量Ib,储箱个数n和形状,以及不同贮箱之间推进剂消耗的关系,在地面事先计算出来。
下面给出了一种平铺四球形贮箱着陆器结构下函数f()和h()的计算方法。如图2所示,着陆器在同一高度上安装有尺寸相同的四个球形贮箱。贮箱半径为R,贮箱几何中心在着陆器本体坐标系(O-XYZ)的坐标分别用矢量rt,i(i=1,…,4)表示。设四个贮箱内推进剂的密度分别为ρi,且四个贮箱内推进剂消耗的体积一样,那么第i个贮箱内剩余推进剂的质量mi可如下计算
mi=ρiV
其中,V表示贮箱内推进剂的体积。
对于球形贮箱,推进剂沉于贮箱下部,构成球冠。设球冠的高度为h,则球冠的体积公式为
这是一个关于h的3次方程,可以用数值方法(例如牛顿迭代法)解出h。接下来根据球冠的重心计算公式,可以求出球冠的质心距离贮箱球心的距离为
那么,可算出第i个贮箱内推进剂质心相对着陆器本体坐标系原点的矢径为
这样就可以求出整器质心(包括干器和推进剂)在着陆器本体坐标系下的位置ro
着陆器整器的惯量可如下计算
其中,E为3×3单位矩阵。
式(4)-(9)实际构成了函数ro=f(m)和I=h(m)。从原理上说(4)-(9)式可直接用于星上计算;但通常星载计算机的能力有限,那么则可以在地面利用式(4)-(9)事先求出不同m下ro和I取值,并用曲线进行拟合。星载计算机根据拟合后的曲线计算ro和I。
(2)姿控力矩计算
姿态控制系统每个周期都会发出三轴喷气脉宽指令,根据它就可以计算出一个控制周期内的等效控制力矩。但是,由于从控制指令计算到推力器响应存在延时τ,因此一个周期内喷气产生的实际控制力矩与指令并不匹配,需要补偿延时的效果。假设上一周期和上两周期,姿态控制器输出的喷气脉宽分别为[Ton,x,k-1,Ton,y,k-1,Ton,z,k-1]和[Ton,x,k-2,Ton,y,k-2,Ton,z,k-2],下标x,y,z分别表示三个坐标轴。以某一个轴为例,它的喷气效果如图3所示,为简单起见图3中的喷气脉宽省略了代表具体某一轴的标号。
对于上一周期来说,实际产生的喷气是由再上一周期残余的喷气时间和该周期内实际产生的喷气时间的累加。因此,该周期内的实际喷气时间为
Ton,i=max(Ton,i,k-2+τ-T,0)+min(T-τ,Ton,i,k-1),(i=x,y,z) (10)
那么一个周期内等效的控制力矩矢量为
MRCS,i是由负责第i轴姿态控制的推力器所产生的矩,它可以如下计算:
假设该推力器在着陆器本体坐标系下的安装位置为rRCS,i,安装方向为pRCS,i,那么它相对整器质心产生的平均力矩为
MRCS,i=-(rRCS,i-ro,k-1)×pRCS,i (12)
(3)角加速度计算
角加速度采用陀螺差分计算而来。假设陀螺测量的是角度增量,即获得了上两周期内的角度增量分别为Δgk-1和Δgk-2,那么这两个周期内的平均角速度为
利用差分可以得到上一周期平均角加速度的近似值为
(4)利用动力学对干扰力矩的解算
将着陆器视为刚体,则根据刚体动力学方程
可以直接反算出干扰力矩的大小,即
(5)干扰力矩估值的平滑
在式(16)的计算过程中,由于角速度差分、姿控发动机推力大小的扰动,使得计算出的干扰力矩存在噪声,为此需要引入适当的滤波器进行平滑。滤波器的选择可以有很多种,这里使用一个最简单的一阶低通滤波进行平滑。
其中,K为0-1之间的常数。
随着时间的递推,重复上述步骤,就完成了干扰力矩的实时估计。
仿真分析。以月球着陆器为例对本发明提出的动力下降过程干扰力矩估计方法进行验证。着陆器剩余质量-质心高度、剩余质量-惯量变化曲线分别如图4和图5所示。
采用本发明的方法对干扰力矩进行估计,滤波系数K取为0.02。由质心变化产生的干扰力矩及其估计值分别如图6和7所示。可以看到,本发明能够较为准确的估计出干扰力矩的变化情况。
本发明未详细阐述部分属于本领域公知技术。
Claims (2)
1.一种动力下降过程干扰力矩的估计方法,其特征在于实现步骤如下:
(1)着陆器质量特性估计
根据发动机比冲和推进剂消耗的关系,估算推进剂剩余质量;然后根据事先存储在星载计算机上的质心-剩余质量变化曲线以及惯量-剩余质量变化曲线,计算出着陆器质心的位置和着陆器惯量的大小;
(2)姿控力矩计算
上一个周期,控制计算机根据姿态和姿态角速度误差计算出着陆器滚动、俯仰和偏航三个轴上的喷气脉宽;根据步骤(1)中着陆器的质心位置计算结果,并结合着陆器上所配备的姿控推力器的安装位置和指向,计算出一个周期内着陆器姿控喷气所产生的平均控制力矩;
(3)角加速度计算
对上两个周期着陆器中陀螺姿态角速度的测量结果进行差分,获得上一周期角加速度的测量值;
(4)利用动力学对干扰力矩的解算
利用步骤(1)~(3)中计算出的着陆器惯量、平均控制力矩和角加速度,通过姿态动力学方程反算出着陆器干扰力矩的大小。
(5)干扰力矩估值的平滑
将干扰力矩在着陆器滚动、俯仰和偏航三个轴的分量分别通过一个低通滤波器进行滤波,获得着陆器干扰力矩的估计值。
2.根据权利要求1所述的动力下降过程干扰力矩的估计方法,其特征在于:所述步骤(1)中估算推进剂剩余质量利用火箭公式估算着,用事先装订的剩余质量-质心高度和剩余质量-惯量曲线估算着陆器质量特性。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310125183.5A CN103303495B (zh) | 2013-04-11 | 2013-04-11 | 一种动力下降过程干扰力矩的估计方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201310125183.5A CN103303495B (zh) | 2013-04-11 | 2013-04-11 | 一种动力下降过程干扰力矩的估计方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN103303495A true CN103303495A (zh) | 2013-09-18 |
CN103303495B CN103303495B (zh) | 2015-07-08 |
Family
ID=49129301
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201310125183.5A Active CN103303495B (zh) | 2013-04-11 | 2013-04-11 | 一种动力下降过程干扰力矩的估计方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN103303495B (zh) |
Cited By (13)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103662090A (zh) * | 2013-12-13 | 2014-03-26 | 北京控制工程研究所 | 一种智能动力下降轨迹在线规划方法 |
CN103926835A (zh) * | 2014-04-04 | 2014-07-16 | 北京航空航天大学 | 一种基于干扰观测器的着陆器动力下降段优化控制方法 |
CN105116906A (zh) * | 2015-07-17 | 2015-12-02 | 中国空间技术研究院 | 基于向量理论的航天器变轨发动机干扰力矩计算方法 |
CN105629987A (zh) * | 2016-03-31 | 2016-06-01 | 北京航空航天大学 | 一种火星着陆器抗干扰容错控制方法 |
CN107792396A (zh) * | 2017-09-25 | 2018-03-13 | 上海卫星工程研究所 | 发动机工作时干扰力矩实时估计方法 |
CN107839903A (zh) * | 2017-09-25 | 2018-03-27 | 上海卫星工程研究所 | 单翼卫星转移轨道段帆板放气时间估计方法 |
CN108279695A (zh) * | 2018-01-25 | 2018-07-13 | 北京控制工程研究所 | 一种航天器干扰力矩的快速在轨闭环辨识方法、系统和介质 |
CN108459610A (zh) * | 2018-02-27 | 2018-08-28 | 北京控制工程研究所 | 一种着陆器动力下降液体晃动抑制方法 |
CN109018440A (zh) * | 2018-07-23 | 2018-12-18 | 哈尔滨工业大学 | 垂直起降运载火箭子级着陆段精确垂直软着陆四次多项式制导方法 |
CN109634306A (zh) * | 2018-12-28 | 2019-04-16 | 北京星际荣耀空间科技有限公司 | 飞行器控制参数的确定方法及装置 |
CN110510155A (zh) * | 2019-09-02 | 2019-11-29 | 北京深蓝航天科技有限公司 | 在线迭代的着陆制导方法及包含该方法的系统、存储器 |
CN110955974A (zh) * | 2019-11-29 | 2020-04-03 | 清华大学 | 一种火箭回收仿真平台及实现方法 |
CN114019936A (zh) * | 2021-09-30 | 2022-02-08 | 北京控制工程研究所 | 一种适用于深空着陆的平移推力器故障诊断与重组方法 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5667167A (en) * | 1994-09-02 | 1997-09-16 | Kistler Aerospace Corporation | Methods and apparatus for reusable launch platform and reusable spacecraft |
CN1983098A (zh) * | 2005-12-14 | 2007-06-20 | 上海微小卫星工程中心 | 主动磁控为主的微小卫星姿态控制方法及系统 |
CN101145027A (zh) * | 2007-10-17 | 2008-03-19 | 北京航空航天大学 | 一种高精度磁悬浮控制力矩陀螺框架伺服控制系统 |
CN101226561A (zh) * | 2007-12-28 | 2008-07-23 | 南京航空航天大学 | 用于微型航天器姿态轨道控制系统的微型仿真支持系统及工作方法 |
CN102620892A (zh) * | 2011-12-15 | 2012-08-01 | 上海卫星工程研究所 | 一种转动部件的动平衡测试方法 |
-
2013
- 2013-04-11 CN CN201310125183.5A patent/CN103303495B/zh active Active
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5667167A (en) * | 1994-09-02 | 1997-09-16 | Kistler Aerospace Corporation | Methods and apparatus for reusable launch platform and reusable spacecraft |
CN1983098A (zh) * | 2005-12-14 | 2007-06-20 | 上海微小卫星工程中心 | 主动磁控为主的微小卫星姿态控制方法及系统 |
CN101145027A (zh) * | 2007-10-17 | 2008-03-19 | 北京航空航天大学 | 一种高精度磁悬浮控制力矩陀螺框架伺服控制系统 |
CN101226561A (zh) * | 2007-12-28 | 2008-07-23 | 南京航空航天大学 | 用于微型航天器姿态轨道控制系统的微型仿真支持系统及工作方法 |
CN102620892A (zh) * | 2011-12-15 | 2012-08-01 | 上海卫星工程研究所 | 一种转动部件的动平衡测试方法 |
Cited By (19)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103662090B (zh) * | 2013-12-13 | 2015-04-22 | 北京控制工程研究所 | 一种智能动力下降轨迹在线规划方法 |
CN103662090A (zh) * | 2013-12-13 | 2014-03-26 | 北京控制工程研究所 | 一种智能动力下降轨迹在线规划方法 |
CN103926835A (zh) * | 2014-04-04 | 2014-07-16 | 北京航空航天大学 | 一种基于干扰观测器的着陆器动力下降段优化控制方法 |
CN103926835B (zh) * | 2014-04-04 | 2016-06-29 | 北京航空航天大学 | 一种基于干扰观测器的着陆器动力下降段优化控制方法 |
CN105116906A (zh) * | 2015-07-17 | 2015-12-02 | 中国空间技术研究院 | 基于向量理论的航天器变轨发动机干扰力矩计算方法 |
CN105629987A (zh) * | 2016-03-31 | 2016-06-01 | 北京航空航天大学 | 一种火星着陆器抗干扰容错控制方法 |
CN107839903B (zh) * | 2017-09-25 | 2020-07-14 | 上海卫星工程研究所 | 单翼卫星转移轨道段帆板放气时间估计方法 |
CN107792396A (zh) * | 2017-09-25 | 2018-03-13 | 上海卫星工程研究所 | 发动机工作时干扰力矩实时估计方法 |
CN107839903A (zh) * | 2017-09-25 | 2018-03-27 | 上海卫星工程研究所 | 单翼卫星转移轨道段帆板放气时间估计方法 |
CN108279695A (zh) * | 2018-01-25 | 2018-07-13 | 北京控制工程研究所 | 一种航天器干扰力矩的快速在轨闭环辨识方法、系统和介质 |
CN108279695B (zh) * | 2018-01-25 | 2021-03-26 | 北京控制工程研究所 | 一种航天器干扰力矩的快速在轨闭环辨识方法、系统和介质 |
CN108459610A (zh) * | 2018-02-27 | 2018-08-28 | 北京控制工程研究所 | 一种着陆器动力下降液体晃动抑制方法 |
CN108459610B (zh) * | 2018-02-27 | 2021-02-05 | 北京控制工程研究所 | 一种着陆器动力下降液体晃动抑制方法 |
CN109018440A (zh) * | 2018-07-23 | 2018-12-18 | 哈尔滨工业大学 | 垂直起降运载火箭子级着陆段精确垂直软着陆四次多项式制导方法 |
CN109634306A (zh) * | 2018-12-28 | 2019-04-16 | 北京星际荣耀空间科技有限公司 | 飞行器控制参数的确定方法及装置 |
CN110510155A (zh) * | 2019-09-02 | 2019-11-29 | 北京深蓝航天科技有限公司 | 在线迭代的着陆制导方法及包含该方法的系统、存储器 |
CN110955974A (zh) * | 2019-11-29 | 2020-04-03 | 清华大学 | 一种火箭回收仿真平台及实现方法 |
CN114019936A (zh) * | 2021-09-30 | 2022-02-08 | 北京控制工程研究所 | 一种适用于深空着陆的平移推力器故障诊断与重组方法 |
CN114019936B (zh) * | 2021-09-30 | 2023-12-29 | 北京控制工程研究所 | 一种适用于深空着陆的平移推力器故障诊断与重组方法 |
Also Published As
Publication number | Publication date |
---|---|
CN103303495B (zh) | 2015-07-08 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN103303495B (zh) | 一种动力下降过程干扰力矩的估计方法 | |
CN103926835B (zh) | 一种基于干扰观测器的着陆器动力下降段优化控制方法 | |
CN105700536B (zh) | 基于绳系拖曳系统的主动星姿态和系绳摆振联合控制方法 | |
CN102163057B (zh) | 一种自治飞艇平面路径跟踪控制方法 | |
CN100495261C (zh) | 月球软着陆制导、导航与控制半物理仿真试验系统 | |
CN102981507B (zh) | 一种软着陆自主障碍规避常推力器控制方法 | |
CN104317300B (zh) | 一种基于模型预测控制的平流层飞艇平面路径跟踪控制方法 | |
CN104142686A (zh) | 一种卫星自主编队飞行控制方法 | |
CN103324083B (zh) | 基于鲁棒观测器的非线性船舶运动控制方法 | |
CN103955223B (zh) | 一种深空探测软着陆过程的姿轨耦合控制方法 | |
CN106773713A (zh) | 针对欠驱动海洋航行器的高精度非线性路径跟踪控制方法 | |
CN103984237A (zh) | 基于运动状态综合识别的轴对称飞行器三通道自适应控制系统设计方法 | |
CN105955028A (zh) | 一种航天器在轨规避制导控制一体化算法 | |
CN105182989A (zh) | 一种风场影响下的飞机姿态控制方法 | |
CN103064423B (zh) | 多约束多航天器飞行间距预示及碰撞规避方法 | |
CN106767797A (zh) | 一种基于对偶四元数的惯性/gps组合导航方法 | |
CN104345735A (zh) | 一种基于落脚点补偿器的机器人行走控制方法 | |
Chung et al. | Posture stabilization strategy for a trotting point-foot quadruped robot | |
CN113867143B (zh) | 地外天体安全软着陆解析避障制导方法 | |
CN106354147A (zh) | 一种基于发动机预摆的大质心横移下的姿轨控耦合控制系统与方法 | |
CN104932266B (zh) | 一种基于前馈补偿的着陆器进入段精确控制方法 | |
CN103708045A (zh) | 一种探月飞船跳跃式再入的在线参数辨识方法 | |
CN103662092A (zh) | 一种衔接主减速和接近段的预测校正方法 | |
EP4033205A1 (en) | Systems and methods for model based vehicle navigation | |
CN103863578B (zh) | 火星着陆器喷气推力器和控制力矩陀螺复合控制系统 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant |