CN114692468B - 一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 - Google Patents
一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 Download PDFInfo
- Publication number
- CN114692468B CN114692468B CN202210434597.5A CN202210434597A CN114692468B CN 114692468 B CN114692468 B CN 114692468B CN 202210434597 A CN202210434597 A CN 202210434597A CN 114692468 B CN114692468 B CN 114692468B
- Authority
- CN
- China
- Prior art keywords
- scale
- strain
- matrix
- calculating
- displacement
- 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.)
- Active
Links
- 239000000835 fiber Substances 0.000 title claims abstract description 154
- 239000011159 matrix material Substances 0.000 title claims abstract description 142
- 238000000034 method Methods 0.000 title claims abstract description 67
- 239000002131 composite material Substances 0.000 title claims abstract description 37
- 230000006378 damage Effects 0.000 title claims abstract description 36
- 239000011226 reinforced ceramic Substances 0.000 title claims abstract description 33
- 238000006073 displacement reaction Methods 0.000 claims abstract description 88
- 238000013179 statistical model Methods 0.000 claims abstract description 10
- 230000008569 process Effects 0.000 claims abstract description 7
- 239000011148 porous material Substances 0.000 claims description 69
- 238000004364 calculation method Methods 0.000 claims description 33
- 239000000463 material Substances 0.000 claims description 25
- 238000004422 calculation algorithm Methods 0.000 claims description 13
- 238000009826 distribution Methods 0.000 claims description 13
- 238000002591 computed tomography Methods 0.000 claims description 11
- 239000011153 ceramic matrix composite Substances 0.000 claims description 8
- 230000007547 defect Effects 0.000 claims description 7
- 230000000694 effects Effects 0.000 claims description 7
- 238000009941 weaving Methods 0.000 claims description 7
- 238000012545 processing Methods 0.000 claims description 6
- 239000006185 dispersion Substances 0.000 claims description 4
- 238000013139 quantization Methods 0.000 claims description 4
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 238000013461 design Methods 0.000 claims description 3
- 239000012783 reinforcing fiber Substances 0.000 claims description 3
- 238000012216 screening Methods 0.000 claims description 3
- 238000010586 diagram Methods 0.000 description 7
- 238000005192 partition Methods 0.000 description 4
- 238000011160 research Methods 0.000 description 4
- 239000000654 additive Substances 0.000 description 3
- 230000000996 additive effect Effects 0.000 description 3
- 238000012986 modification Methods 0.000 description 3
- 230000004048 modification Effects 0.000 description 3
- 238000000638 solvent extraction Methods 0.000 description 3
- 241001229889 Metis Species 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 238000009954 braiding Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 238000012417 linear regression Methods 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 239000000203 mixture Substances 0.000 description 2
- 238000002360 preparation method Methods 0.000 description 2
- 238000007619 statistical method Methods 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 208000027418 Wounds and injury Diseases 0.000 description 1
- 238000005299 abrasion Methods 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 238000005520 cutting process Methods 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 230000006735 deficit Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 238000011438 discrete method Methods 0.000 description 1
- 229920001971 elastomer Polymers 0.000 description 1
- 239000000806 elastomer Substances 0.000 description 1
- 239000004744 fabric Substances 0.000 description 1
- 208000014674 injury Diseases 0.000 description 1
- 230000003902 lesion Effects 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- 230000003647 oxidation Effects 0.000 description 1
- 238000007254 oxidation reaction Methods 0.000 description 1
- 238000011002 quantification Methods 0.000 description 1
- 230000002787 reinforcement Effects 0.000 description 1
- 230000003014 reinforcing effect Effects 0.000 description 1
- 238000003860 storage Methods 0.000 description 1
- 239000000126 substance Substances 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 239000011800 void material Substances 0.000 description 1
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
- G06F2113/00—Details relating to the application field
- G06F2113/26—Composites
-
- 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/02—Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
-
- 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
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)
Abstract
本发明公开了一种一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,包括建立宏观尺度构件有限元模型;建立纤维束尺度有限元模型;建立其统计模型;建立纤维丝尺度有限元模型;计算宏观尺度下所有单元的高斯点应变;计算纤维束尺度下所有单元的高斯点应变;计算其纤维束尺度位移;计算其宏观尺度位移;利用相场法求解相场并更新宏观尺度的边界条件,迭代整个过程,从而预测损伤演化。本发明可快速实现连续纤维增强陶瓷基复合材料宏观形状、位置和角度的精确建模,从底层的纤维丝尺度开始计算,全面考虑了连续纤维增强陶瓷基复合材料的纤维束尺度下纤维形状对宏观尺度损伤的影响,能够更加准确的预测连续纤维增强陶瓷基复合材料的损伤。
Description
技术领域
本发明属于连续纤维增强陶瓷基复合材料跨尺度损伤预测研究技术领域,具体是一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法。
背景技术
连续纤维增强陶瓷基复合材料自诞生以来,因其耐高温、高强度、强韧性、低密度、热膨胀小及抗磨损等优良性能,被广泛用于制造飞机、航天器、火箭等发动机和核能的关键部件。目前国内外各研究机构在材料制备、性能改良、制造工艺等相关领域开展了诸多开创性的研究工作。总体而言,我国在连续纤维增强陶瓷基复合材料基础研究、产业化和应用等方面与西方国家仍存在一定差距。
连续纤维增强陶瓷基复合材料的宏观性能及损伤机制不仅与基体材料特性密切相关,且与微结构组成、界面相组成及分布、增强相黏着/接触状态、孔隙特征和裂纹分布等有关。鉴于其制备过程影响因素众多,目前,在不同相热膨胀系数不匹配、增强相抗氧化、材料宏观参数精确预测和内部损伤跨尺度表征等方面仍存在较多亟待解决的问题。长期以来,连续纤维增强陶瓷基复合材料的材料损伤多采用实验测试获得,难以对其进行预测。鉴于其微结构缺陷存在随机性,且不同位置的编织方式存在差异,导致测试结果存在一定差异。
因此,建立连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,并系统分析各种因素的影响机制,不仅有助于指导微结构编织方式选择和材料性能改善,更可优化材料合成工艺和合成新型陶瓷基复合材料。
发明内容
本发明的目的在于提供一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,该方法结合连续纤维增强陶瓷基复合材料微结构模型参数模型,可实现陶瓷基复合材料跨尺度建模;通过相场法引入纤维束尺度计算结果对宏观尺度的损伤进行预测。
本发明解决上述技术问题所提供的技术方案是:一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,包括以下步骤:
S1、建立连续纤维增强陶瓷基复合材料的宏观尺度构件有限元模型;
S2、建立包含任意复杂编织预制体的纤维束尺度有限元模型;
S3、通过导入CT扫描得到的孔隙信息,建立其形状、位置、尺寸、分布的统计模型,并在纤维束尺度模型引入自然缺陷特征;
S4、建立纤维丝尺度有限元模型,并通过统计模型引入自然缺陷特征;
S5、引入并行计算,计算宏观尺度下所有单元的高斯点应变,作为纤维束尺度下模型的边界条件;
S6、引入并行计算,计算纤维束尺度下所有单元的高斯点应变,以作为纤维丝尺度下模型的边界条件;
S7、采用CPU+GPU并行计算算法,在纤维丝尺度下计算平均应力和平均应变返回纤维束尺度计算其纤维束尺度位移;
S8、在纤维束尺度下计算平均应力和平均应变返回宏观尺度计算其宏观尺度位移;
S9、在宏观尺度下,引入上一步计算处的位移场,利用相场法求解相场并更新宏观尺度的边界条件,迭代整个过程,从而预测损伤演化。
进一步的技术方案是,所述步骤S1具体实现方法为:建立宏观尺度模型,修改其边界从而预设一条裂纹。
进一步的技术方案是,所述步骤S2具体实现方法为:
S21、建立预制体参数化模型;
结合增强纤维和纤维束截面形状及编织特点,开展了复杂预制体参数化精确建模,实现一维、二维、2.5维、三维增韧相参数化建模,同时完成不同形状纤维束及预制体设计;
S22、预制体水平集函数描述和自动离散;
结合水平集函数φ(i)(x,y,z)的值和预制体参数化模型,建立三维四面体节点处的水平集函数值和距离值;若单元同时存在函数fR(x,y,z)<0和函数fR(x,y,z)>0的节点,则该单元穿过增韧相边界Г,被定义为相交单元;则通过单元自动离散方法对该相交单元进行细分和编号及存储更新。
进一步的技术方案是,所述步骤S3具体实现方法为:
S31、对连续纤维增韧陶瓷基复合材料进行断层CT扫描,智能识别连续纤维增韧陶瓷基复合材料内的所有自然孔隙,获得自然孔隙的图像及CT扫描孔隙数据;
S32、提取CT扫描孔隙数据,并产生有限元模型,引入如下的孔隙量化描述函数:
式中:x、y、z为孔隙笛卡尔坐标系中的坐标;Rx、Ry、Rz分别表示为x、y、z方向上的半径;
同时根据如下的孔隙随机分布参数Dv来描述其分布特征和Ev来量化孔隙逼近效果;
式中:Vp为孔隙的实际体积;Vf为拟合得到的实际体积;Vc为上述两个区域重合的体积;
S33、通过基于水平集函数的自动建模方法,建立孔隙结构;
S34、筛选出所有自然孔隙区域包含的所有四面体单元并记录;
S35、将所有自然孔隙区域的四面体单元材料中的模量设为0。
进一步的技术方案是,所述步骤S4具体实现方法为:
S41、通过基于水平集函数的自动建模方法,建立纤维丝尺度模型;
S42、通过统计模型引入纤维丝尺度的孔隙特征,通过基于水平集函数的自动建模建立孔隙模型并将孔隙模型内的单元材料模量设为0。
进一步的技术方案是,所述步骤S5具体实现方法为:
S51、通过以下公式计算宏观尺度下所有单元的高斯点应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为高斯点应变;
S52、将宏观尺度下计算得到的高斯点应变,作为纤维束尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;
S53、将宏观尺度模型的所有高斯点分成N份,采用并行命令对N份高斯点同时进行处理,设置下层纤维束尺度的边界条件。
进一步的技术方案是,所述步骤S6具体实现方法为:
S61、通过以下公式计算纤维束尺度下所有单元的高斯点应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S62、将纤维束尺度下计算得到的高斯点应变,作为纤维丝尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;
S63、将纤维束尺度模型的所有高斯点分成N份,采用并行命令对N份高斯点同时进行处理,设置下层纤维丝尺度的边界条件。
进一步的技术方案是,所述步骤S7具体实现方法为:
S71、引入并行计算用于计算位移;
S72、利用GPU加速计算上一步中拆分的刚度矩阵,使用Cude库,对拆分出刚度矩阵建立相应的线程并进行计算;
单元刚度矩阵计算公式如下:
式中:Ke为单元刚度矩阵;B为应变矩阵;D为弹性矩阵;
刚度矩阵组装计算公式如下:
Kii=∑(GTKeG)
式中:Kii为拆分的刚度矩阵;G为单元节点的自由度与结构节点的自由度之间的转化矩阵;
S73、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S74、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S75、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S76、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回纤维束尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅各比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算纤维束尺度结合纤维丝尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵。
进一步的技术方案是,所述步骤S8具体实现方法为:
S81、利用上一步计算得到的纤维束尺度下位移矩阵,通过遍历单元可得到每个单元的节点位移矩阵q;
S82、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S83、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S84、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S85、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回宏观尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅各比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算宏观尺度结合纤维束尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵。
进一步的技术方案是,所述步骤S9中的计算公式包括:
损伤场函数方程:
断裂韧性函数方程:
g(d)=(1-d)2
相场有限元控制方程:
[Kφ]{d(x)}n+1={F}n
{Fd}n=∫Ω2Hn[Nd]TdΩ
位移场控制方程:
([K1]+[K2]){u}n+1={F}n+1
式中:d(x)为损失函数;g(d)表示断裂韧性函数。
本发明的有益效果是:本发明基于水平集方法的连续纤维增强陶瓷基复合材料微结构建立纤维束尺度及宏观编织体尺度有限元模型,以水平集函数描述纤维束尺度及宏观编织体尺度各组成成分和大小,完成的连续纤维增强陶瓷基复合材料有限元数字模型;考虑了连续纤维增强陶瓷基复合材料的纤维束尺度和宏观尺度编织体的微结构特点和微裂纹与孔隙,利用旋转矩阵有效设定不同编织位置宏观横观各向同性材料属性,通过相场法计算跨尺度影响下宏观尺度所发生的损伤。
附图说明
图1是本发明的一种连续纤维增强陶瓷基复合材料微结构精确建模及宏观模量自动预测方法的流程图;
图2是本发明的宏观尺度有限元模型示意图;
图3为单元细分思路图;
图4为连续纤维增韧SiCf/SiC CMC有限元模型图;
图5为2D连续纤维增韧SiCf/SiC复合材料有限元模型图;
图6为2.5D连续纤维增强SiCf/SiC复合材料有限元模型图;
图7为3D连续纤维增韧SiCf/SiC复合材料增强相离散模型图;
图8为CMC纤维束尺度孔隙信息图;
图9为拟合部分结果图;
图10为复杂3D自然孔隙拟合及逼近效果图;
图11是孔隙拟合效果图;
图12是自然孔隙统计参数示意图;
图13是2.5D CFCMC自然孔隙参数统计;
图14是CFCMC小孔隙方向性统计分析;
图15是纤维丝尺度模型;
图16是CMC模型区域分解效果图;
图17是区域分解算法原理中网格划分及分区情况;
图18是相场法预测损伤示意图;
图19为宏观尺度上部件的损伤示意图;
图20为孔隙分布、尺寸等CMC对损伤演化影响。
具体实施方式
下面将结合附图对本发明的技术方案进行清楚、完整地描述,显然,所描述的实施例是发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
如图1所示,本发明提供了一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,包括以下步骤:
S1、建立连续纤维增强陶瓷基复合材料的宏观尺度构件有限元模型;
具体的实现步骤为:建立宏观尺度模型,修改宏观尺度模型的边界从而预设一条裂纹,即建立宏观尺度构件有限元模型;
S2、建立包含任意复杂编织预制体的纤维束尺度有限元模型;
具体包括以下步骤:
S21、建立预制体参数化模型;
结合增强纤维和纤维束截面形状及编织特点,开展了复杂预制体参数化精确建模,提出了如下的参数化模型,实现了一维(长纤维/短切纤维)、二维(编织布)、2.5维(编织穿刺)、三维(编织体)等增韧相参数化建模,同时可完成不同形状纤维束及预制体设计,其中如下第一个为一维的预制体参数化模型,第二个为二维和2.5维的预制体参数化模型,第三个为三维的预制体参数化模型;
S22、预制体水平集函数描述和自动离散;
为了建立预制体结构几何建模和材料各向异性准确幅值,上述增韧纤维外表面Γ可进一步表达为水平集函数φ(i)(x,y,z)=0;则空间任意点P(xi,yi,zi)到边界曲面的距离可描述为函数fR(x,y,z);
fR(x,y,z)=|φ(i)(x,y,z)|
因此,Г表面上的任意一点满足:
fR(x,y,z)=0
而Г表面外任意一点满足:
fR(x,y,z)>0
和Г表面内任意一点满足:
fR(x,y,z)<0
结合水平集函数φ(i)(x,y,z)的值和预制体参数化模型,建立三维四面体节点处的水平集函数值和距离值;若单元同时存在函数fR(x,y,z)<0和函数fR(x,y,z)>0的节点,则该单元穿过增韧相边界Г,被定义为相交单元;则通过单元自动离散方法对该相交单元进行细分和编号及存储更新,单元细分思路见图3;
采用水平集法描述纤维束外边界,同时结合上述3D四面体单元自动切割技术,实现了SiCf/SiC复合材料的微结构自动生成,效果见图4-图7;
S3、通过导入CT扫描得到的孔隙信息,建立其形状、位置、尺寸、分布等统计模型,并在纤维束尺度模型引入自然缺陷(孔隙、裂纹和界面)特征;
具体包括以下步骤:
S31、利用计算机断层扫描技术对连续纤维增韧陶瓷基复合材料进行断层扫描,智能识别连续纤维增韧陶瓷基复合材料内的所有自然孔隙,获得自然孔隙的图像及CT扫描孔隙数据,见图8;
S32、提取CT扫描孔隙数据,并产生有限元模型,引入如下的孔隙量化描述函数:
式中:x、y、z为孔隙笛卡尔坐标系中的坐标;Rx、Ry、Rz分别表示为x、y、z方向上的半径;
其中通过调整(Rx,Ry,Rz)和θ可调整孔隙尺寸、形状等信息;
同时根据如下的孔隙随机分布参数Dv来描述其分布特征和Ev来量化孔隙逼近效果;部分拟合图像见图9;
式中:Vp为孔隙的实际体积;Vf为拟合得到的实际体积;Vc为上述两个区域重合的体积;
具体量化过程为:1)计算孔隙域内单元中心点坐标和体积,存储于4维矩阵;2)以单元体积×(单元中心-主轴距离)为变量,通过线性回归分析寻优获得主轴1;3)在垂直于主轴1面内,再次使用线性回归分析确定主轴2,而后确定主轴3;4)确定最佳逼近半径和指数;获得的部分孔隙拟合效果见图11;
为了准确反映自然孔隙特征及其影响,同时实现该孔隙特征的参数化自动建模,本发明选定孔隙方向、形状、尺寸及旋转角度进行统计分析(见图10);针对获得的不同孔隙信息采用合适的间隔做柱状图统计,得到2.5维(编织穿刺)孔隙参数信息的分布图见图13;直方图的因变量是每个区间内孔隙数量与总数量的比值;
同时,还进一步统计了小孔隙(体积小于最大孔隙体积1/10)的方向性,如图14所示;统计小孔隙的方向和尺寸发现,其包络面近似接近球状,故可认为此类孔隙无方向性,以降低建模难度;
S33、通过基于水平集函数的自动建模方法,建立孔隙结构;
S34、筛选出所有自然孔隙区域包含的所有四面体单元并记录;
S35、将所有自然孔隙区域的四面体单元材料中的模量设为0
S4、建立纤维丝尺度有限元模型,并通过统计模型引入自然缺陷(孔隙、裂纹和界面)特征;
具体包括以下步骤:
S41、通过基于水平集函数的自动建模方法,建立纤维丝尺度模型;
S42、通过上一步的统计模型引入纤维丝尺度的孔隙特征,通过基于水平集函数的自动建模建立孔隙模型(如图15所示)并将孔隙模型内的单元材料模量设为0;
S5、采用开发的CPU+GPU并行计算算法,计算宏观尺度下所有单元的高斯点应变,作为纤维束尺度下模型的边界条件;
具体包括以下步骤:
S51、通过以下公式计算宏观尺度下所有单元的高斯点的应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S52、将宏观尺度下计算得到的高斯点的应变,作为纤维束尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;需要注意的是,x+与x-必须关于yoz平面对称,即在模型外表面上任意节点,必须在对称面上存在相同坐标的节点与之对应;而表示周期边界条件;
S53、将宏观尺度模型的所有高斯点分成N份,采用openmp中的并行命令parallelfor对N份高斯点同时进行处理,设置下层纤维束尺度的边界条件并计算其应力应变;
S6、采用开发的CPU+GPU并行计算算法,计算纤维束尺度下所有单元的高斯点应变,以作为纤维丝尺度模型的边界条件;
具体包括以下步骤:
S61、通过以下公式计算纤维束尺度下所有单元的高斯点的应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S62、将纤维束尺度下计算得到的高斯点的应变,作为纤维丝尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;
S63、将纤维束尺度模型的所有高斯点分成N份,采用OpenMP中的并行命令parallel for对N份高斯点同时进行处理,设置下层纤维丝尺度的边界条件并计算其应力应变;
S7、引入并行计算,在纤维丝尺度下计算平均应力和平均应变返回纤维束尺度计算其纤维束尺度位移;
具体包括以下步骤:
S71、引入并行计算用于计算位移:
在线性系统KU=F中,K是大小为m×m的矩阵,F是大小为m的向量,解集U是大小为m的向量,索引集{1,...,m}被划分为两个集N1={1,...,ms}和N2={ms+1,...,m};
令同样的,令/>则线性系统有以下分块矩阵形式:
其中令D为K的对角矩阵,则雅可比算法可以写作:
DUn+1=DUn+(F-KUn)
定义离散层面的残差rn=F-KUn,则上式可以写作:
Un+1=Un+D-1(F-KUn)=Un+D-1rn
对雅可比算法进行分块处理,可定义块雅可比算法:
为了使该算法形式更加紧凑,引入R1作为N到N1的限制算子,引入R2作为N到N2的限制算子;转置算子是从Ni到N的扩展算子,有/>
上式也可以写为:
即:
又有:
因此,块雅可比算法式可以写为紧凑形式:
为了与连续层面的限制加性施瓦兹算法建立联系,列出线性代数对应的参数:限制算子Ri是大小为#Ni×#N的矩形布尔矩阵;作用相当于在连续层面上从u获取ui;转置算子Ri T是扩展算子,作用相当于连续层面上的扩展算符Ei;而离散层面上的单位函数χi是对应于大小为#Ni×#Ni的对角矩阵Id,其矩阵形式是:
下面以一个一维有限元问题为例,将以最小重叠域的分区情况展示线性代数参数的具体形式。图17展示一维有限元问题的网格划分及分区情况:
节点数为5,索引集{1,2,3,4,5}被划分为两个集N1={1,2,3}和N2={3,4,5};则限制矩阵R1和R2大小都为3×5,其具体数值为:
而满足的最简单的单位矩阵则是:
而对于索引集的划分,对于简单、规则并且单元编号有序的模型可以用单元的序号或者节点坐标来区分;但该种方式不具有普遍性,对于复杂模型或单元编号无序的模型就会失效;目前对于任意索引集的区域划分操作大多是用METIS或者SCOTCH这种自动图形分割软件包进行的,本文将选用METIS作为区域划分的手段;
根据最简二分区情况作为推导过程,而实际情况往往是将划分好网格的有限元模型分为多个分区进行分区求解并迭代计算,基于下式有适用于离散层面任意分区数的限制加性施瓦兹算法:
将其中的迭代矩阵定义为限制加性施瓦兹算法的预处理器,其具体形式为:
则上式可以写作:
S72、利用GPU加速计算上一步中拆分的刚度矩阵,使用Cude库,对拆分出刚度矩阵建立相应的线程并进行计算;
单元刚度矩阵计算公式如下:
式中:Ke为单元刚度矩阵;B为应变矩阵;D为弹性矩阵;
刚度矩阵组装计算公式如下:
Kii=∑(GTKeG)
式中:Kii为拆分的刚度矩阵;G为单元节点的自由度与结构节点的自由度之间的转化矩阵;
S73、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S74、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S75、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S76、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回纤维束尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅各比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算纤维束尺度结合纤维丝尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵;
S8、在纤维束尺度下计算平均应力和平均应变返回宏观尺度计算其宏观尺度位移;
具体包括以下步骤:
S81、利用上一步计算得到的纤维束尺度下位移矩阵,通过遍历单元可得到每个单元的节点位移矩阵q;
S82、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S83、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S84、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S85、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回宏观尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅各比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算宏观尺度结合纤维束尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵;
S9、宏观尺度下,引入上一步计算处的位移场,利用相场法求解相场并更新宏观尺度的边界条件,迭代整个过程,从而预测损伤演化;
其中相场法的核心思想是用弥散区域来代表裂纹,线弹性体材料内部当损伤达到一定程度时,材料开始产生断裂,在相场法中为了进行量化,引入损伤场函数d(x)(相场函数)对损伤值进行描述,其中当d(x)=0时,认为材料内部完好无损,当d(x)=1时,认为材料内部完全损坏,即该部分产生裂纹;
损伤场描述,g(d)表示断裂韧性函数,为了满足在局部损伤区域应变能密度函数取有限值,g(1)=0;为了保证材料初始未产生损伤,g(0)=1;为了保证材料完全损伤,g(1);因此,该函数表达式可以写为:
g(d)=(1-d)2
相场有限元控制方程:
[Kφ]{d(x)}n+1={F}n
{Fd}n=∫Ω2Hn[Nd]TdΩ
位移场控制方程:
([K1]+[K2]){u}n+1={F}n+1
/>
式中:d(x)为损失函数;g(d)表示断裂韧性函数。
本发明可快速实现连续纤维增强陶瓷基复合材料宏观形状、位置和角度的精确建模,从底层的纤维束尺度开始计算,全面考虑了连续纤维增强陶瓷基复合材料的纤维束尺度下纤维形状对宏观尺度损伤的影响,能够更加准确的预测连续纤维增强陶瓷基复合材料的损伤。
以上所述,并非对本发明作任何形式上的限制,虽然本发明已通过上述实施例揭示,然而并非用以限定本发明,任何熟悉本专业的技术人员,在不脱离本发明技术方案范围内,可利用上述揭示的技术内容作出些变动或修饰为等同变化的等效实施例,但凡是未脱离本发明技术方案的内容,依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰,均仍属于本发明技术方案的范围内。
Claims (6)
1.一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,包括以下步骤:
S1、建立连续纤维增强陶瓷基复合材料的宏观尺度构件有限元模型;
建立宏观尺度模型,修改其边界从而预设一条裂纹;
S2、建立包含任意复杂编织预制体的纤维束尺度有限元模型;
S21、建立预制体参数化模型;
结合增强纤维和纤维束截面形状及编织特点,开展了复杂预制体参数化精确建模,实现一维、二维、2.5维、三维增韧相参数化建模,同时完成不同形状纤维束及预制体设计;
S22、预制体水平集函数描述和自动离散;
结合水平集函数φ(i)(x,y,z)的值和预制体参数化模型,建立三维四面体节点处的水平集函数值和距离值;若单元同时存在函数fR(x,y,z)<0和函数fR(x,y,z)>0的节点,则该单元穿过增韧相边界Г,被定义为相交单元;则通过单元自动离散方法对该相交单元进行细分和编号及存储更新;
S3、通过导入CT扫描得到的孔隙信息,建立其形状、位置、尺寸、分布的统计模型,并在纤维束尺度模型引入自然缺陷特征;
S4、建立纤维丝尺度有限元模型,并通过统计模型引入自然缺陷特征;
S41、通过基于水平集函数的自动建模方法,建立纤维丝尺度模型;
S42、通过统计模型引入纤维丝尺度的孔隙特征,通过基于水平集函数的自动建模建立孔隙模型并将孔隙模型内的单元材料模量设为0;
S5、引入并行计算,计算宏观尺度下所有单元的高斯点应变,作为纤维束尺度下模型的边界条件;
S6、引入并行计算,计算纤维束尺度下所有单元的高斯点应变,以作为纤维丝尺度下模型的边界条件;
S7、采用CPU+GPU并行计算算法,在纤维丝尺度下计算平均应力和平均应变,返回纤维束尺度计算其纤维束尺度位移;
S8、在纤维束尺度下计算平均应力和平均应变,返回宏观尺度计算其宏观尺度位移;
S9、在宏观尺度下,引入上一步计算处的位移场,利用相场法求解相场并更新宏观尺度的边界条件,迭代整个过程,从而预测损伤演化。
2.根据权利要求1所述的一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,所述步骤S3具体实现方法为:
S31、对连续纤维增韧陶瓷基复合材料进行断层CT扫描,智能识别连续纤维增韧陶瓷基复合材料内的所有自然孔隙,获得自然孔隙的图像及CT扫描孔隙数据;
S32、提取CT扫描孔隙数据,并产生有限元模型,引入如下的孔隙量化描述函数:
式中:x、y、z为孔隙笛卡尔坐标系中的坐标;Rx、Ry、Rz分别表示为x、y、z方向上的半径;
同时根据如下的孔隙随机分布参数Dv来描述其分布特征和Ev来量化孔隙逼近效果;
式中:Vp为孔隙的实际体积;Vf为拟合得到的实际体积;Vc为孔隙的实际区域和拟合得到的实际区域重合的体积;
S33、通过基于水平集函数的自动建模方法,建立孔隙结构;
S34、筛选出所有自然孔隙区域包含的所有四面体单元并记录;
S35、将所有自然孔隙区域的四面体单元材料中的模量设为0。
3.根据权利要求1所述的一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,所述步骤S5具体实现方法为:
S51、通过以下公式计算宏观尺度下所有单元的高斯点的应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S52、将宏观尺度下计算得到的高斯点应变,作为纤维束尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;
S53、将宏观尺度模型的所有高斯点分成N份,采用并行命令对N份高斯点同时进行处理,设置下层纤维束尺度的边界条件。
4.根据权利要求1所述的一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,所述步骤S6具体实现方法为:
S61、通过以下公式计算纤维束尺度下所有单元的高斯点的应变;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S62、将纤维束尺度下计算得到的高斯点应变,作为纤维丝尺度的边界条件,其对应公式为;
式中:表示在x轴正方向且与x轴垂直的有限元模型外表面上节点x+所对应的位移矢量;/>表示在x轴负方向且与x轴垂直的有限元模型外表面上节点x-所对应的位移矢量;
S63、将纤维束尺度模型的所有高斯点分成N份,采用并行命令对N份高斯点同时进行处理,设置下层纤维丝尺度的边界条件。
5.根据权利要求1所述的一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,所述步骤S7具体实现方法为:
S71、引入并行计算用于计算位移;
S72、利用GPU加速计算上一步中拆分的刚度矩阵,使用Cuda库,对拆分出刚度矩阵建立相应的线程并进行计算;
单元刚度矩阵计算公式如下:
式中:Ke为单元刚度矩阵;B为应变矩阵;D为弹性矩阵;
刚度矩阵组装计算公式如下:
Kii=∑(GTKeG)
式中:Kii为拆分的刚度矩阵;G为单元节点的自由度与结构节点的自由度之间的转化矩阵;
S73、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S74、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S75、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S76、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回纤维束尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅可比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算纤维束尺度结合纤维丝尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵。
6.根据权利要求5所述的一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法,其特征在于,所述步骤S8具体实现方法为:
S81、利用上一步计算得到的纤维束尺度下位移矩阵,通过遍历单元可得到每个单元的节点位移矩阵q;
S82、根据节点位移矩阵q计算应变ε;
ε=B·q
式中:q为节点位移矩阵;B为应变矩阵;ε为应变;
S83、根据应变ε计算应力σ;
σ=Dε
式中:σ为应力;
S84、根据应力σ和应变ε分别计算平均应力和平均应变;
平均应力矢量计算公式为:
式中:为平均应力矢量;
平均应变矢量计算公式为:
式中:为平均应变矢量;
S85、将计算得到的平均应变带回纤维束尺度,计算位移增量;
首先,将计算得到的平均应变带回宏观尺度;
式中:B为应变矩阵;为纤维丝尺度返回的平均应力;J为雅各比矩阵;h为权重系数;
再计算位移增量;
KΔU=fext-fint
式中:fext为力矩阵;ΔU为位移增量;
最后计算宏观尺度结合纤维束尺度的结合位移矩阵;
Un+1=Un+ΔU
式中:Un为位移矩阵;Un+1为结合位移矩阵。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210434597.5A CN114692468B (zh) | 2022-04-24 | 2022-04-24 | 一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210434597.5A CN114692468B (zh) | 2022-04-24 | 2022-04-24 | 一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114692468A CN114692468A (zh) | 2022-07-01 |
CN114692468B true CN114692468B (zh) | 2023-09-19 |
Family
ID=82145669
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210434597.5A Active CN114692468B (zh) | 2022-04-24 | 2022-04-24 | 一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114692468B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116230118A (zh) * | 2022-11-23 | 2023-06-06 | 西南交通大学 | 一种基于主动学习的连续纤维增强复合材料按需设计方法 |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105701312A (zh) * | 2015-12-17 | 2016-06-22 | 南京航空航天大学 | 复杂编织结构陶瓷基复合材料疲劳迟滞行为预测方法 |
CN105760605A (zh) * | 2015-12-17 | 2016-07-13 | 南京航空航天大学 | 复杂编织结构陶瓷基复合材料疲劳寿命预测方法 |
CN107273649A (zh) * | 2017-08-16 | 2017-10-20 | 中国石油大学(华东) | 一种脆性材料在高温蠕变状态下失效概率的预测方法 |
CN109583052A (zh) * | 2018-11-15 | 2019-04-05 | 西南交通大学 | 纤维增强陶瓷基复合材料微结构精确建模及模量计算方法 |
CN109657412A (zh) * | 2019-01-22 | 2019-04-19 | 北京航空航天大学 | 一种陶瓷基复合材料的疲劳寿命的多尺度预测方法 |
CN113033042A (zh) * | 2021-03-08 | 2021-06-25 | 西北工业大学 | 连续纤维增韧陶瓷基复合材料自然孔隙信息拟合方法 |
Family Cites Families (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109614755B (zh) * | 2018-12-29 | 2023-04-07 | 南京航空航天大学 | 一种通过迟滞耗散能预测编织陶瓷基复合材料高温疲劳纤维/基体界面剪应力的方法 |
CN111523237B (zh) * | 2020-04-24 | 2023-01-17 | 南京航空航天大学 | 一种编织陶瓷基复合材料蠕变断裂行为的预测方法 |
-
2022
- 2022-04-24 CN CN202210434597.5A patent/CN114692468B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105701312A (zh) * | 2015-12-17 | 2016-06-22 | 南京航空航天大学 | 复杂编织结构陶瓷基复合材料疲劳迟滞行为预测方法 |
CN105760605A (zh) * | 2015-12-17 | 2016-07-13 | 南京航空航天大学 | 复杂编织结构陶瓷基复合材料疲劳寿命预测方法 |
CN107273649A (zh) * | 2017-08-16 | 2017-10-20 | 中国石油大学(华东) | 一种脆性材料在高温蠕变状态下失效概率的预测方法 |
CN109583052A (zh) * | 2018-11-15 | 2019-04-05 | 西南交通大学 | 纤维增强陶瓷基复合材料微结构精确建模及模量计算方法 |
CN109657412A (zh) * | 2019-01-22 | 2019-04-19 | 北京航空航天大学 | 一种陶瓷基复合材料的疲劳寿命的多尺度预测方法 |
CN113033042A (zh) * | 2021-03-08 | 2021-06-25 | 西北工业大学 | 连续纤维增韧陶瓷基复合材料自然孔隙信息拟合方法 |
Non-Patent Citations (2)
Title |
---|
Numerical analysis of a ceramic matrix composites with strain induced damage using finite element method;Ratnakar Pandu等;《International Conference on Advanced Nanomaterials & Emerging Engineering Technologies》;第648-654页 * |
孔隙对陶瓷基复合材料强度影响的研究进展;鲁明远等;《材料导报》;第180-185页 * |
Also Published As
Publication number | Publication date |
---|---|
CN114692468A (zh) | 2022-07-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Scerrato et al. | On the effect of shear stiffness on the plane deformation of linear second gradient pantographic sheets | |
Gyulassy et al. | Topologically clean distance fields | |
CN112836411B (zh) | 加筋板壳结构的优化方法、装置、计算机设备和存储介质 | |
CN109344524B (zh) | 一种薄板结构加强筋分布优化方法 | |
CN114692468B (zh) | 一种连续纤维增强陶瓷基复合材料跨尺度损伤预测方法 | |
Hadi et al. | Design configurations and creation of lattice structures for metallic additive manufacturing | |
Cai et al. | A mixed cover meshless method for elasticity and fracture problems | |
CN114792037B (zh) | 一种超材料隔振器的序贯稳健性优化设计方法 | |
Tian et al. | Advanced 4‐node tetrahedrons | |
Jaworski et al. | Implementation features of composite materials effective mechanical characteristics finding method based on microlevel cellular structural models | |
CN108932385A (zh) | 一种编织复合材料内部变截面纤维束代表性体元的建模方法 | |
Kermani et al. | Free vibration optimization of a sandwich structure with a hybrid composite grid core | |
CN114741782A (zh) | 一种圆锥壳加筋结构优化方法、装置、计算机和存储介质 | |
Liu et al. | Intelligent optimization of stiffener unit cell via variational autoencoder-based feature extraction | |
Wang et al. | A compact yet flexible design space for large-scale nonperiodic 3D woven composites based on a weighted game for generating candidate tow architectures | |
Winslow et al. | Mapping two-way grids onto free-form surfaces | |
CN116386781A (zh) | 一种梯度多构型点阵结构的多尺度拓扑优化方法 | |
CN114741784B (zh) | 一种圆柱壳加筋结构优化方法、装置、计算机和存储介质 | |
Rao et al. | Multi-objective optimal design of hybrid laminate composite structures using scatter search | |
CN115547431A (zh) | 一种基于机器学习与渐进损伤的短切纤维复合材料疲劳性能预测方法 | |
Ahlquist et al. | Evolutionary structural and spatial adaptation of topologically differentiated tensile systems in architectural design | |
Wang et al. | Variable-thickness optimization method for shell structures based on a regional evolutionary control strategy | |
Karutz et al. | A quadtree data structure for the coupled finite‐element/element‐free Galerkin method | |
CN116911059B (zh) | 一种3d机织复合材料参数化建模的方法及装置 | |
Schultz | Meso-scale and multicontinuum modeling of a triaxial braided textile composite |
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 | ||
TR01 | Transfer of patent right | ||
TR01 | Transfer of patent right |
Effective date of registration: 20240311 Address after: Room 124, 1st to 2nd floors, No.124 Xingsheng Street, Huangcun Town, Daxing District, Beijing, 102627 (cluster registration) Patentee after: Beijing Huifengyuan Technology Co.,Ltd. Country or region after: China Address before: 610000 North Section of Second Ring Road, Chengdu City, Sichuan Province Patentee before: SOUTHWEST JIAOTONG University Country or region before: China |