《基于impsco和改进newmark-β算法的结构系统及激励辨识研究-麻胜兰.pdf》由会员分享，可在线阅读，更多相关《基于impsco和改进newmark-β算法的结构系统及激励辨识研究-麻胜兰.pdf（7页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、振动与冲击第36卷第15期 JOURNAL OF VIBRATION AND SHOCK V0136 No15 2017基于IMPSCO和改进Newmark-fl算法的结构系统及激励辨识研究麻胜兰，姜绍飞，陈志刚(福州大学土木工程学院，福州350108)摘 要：基于振动的土木工程结构系统识别中，常用的识别方法如基于Newmark邶算法的方法，常常需要已知激励，而事实上激励和完整的输出响应是很难测得的。基于此，在现有的Newmark邶算法中引入多项式拟合来解决结构系统识别误差带来的漂移现象，进而提出了改进的Newmark筇算法，并结合改进协同粒子群算法(IMPSCO)，给出了仅用部分结构响应的系
2、统识别和激励辨识方法。最后，通过一七层钢框架数值算例和实验室模型试验验证了方法的有效性和可行性，并探讨了噪声、输出响应完整性的影响。研究表明，所提方法能准确地实现未知激励和部分实测响应状态下的结构系统和未知激励的识别，而且具有较强的容噪性和鲁棒性。关键词：未知激励；系统识别；改进的Newmark-3算法中图分类号：TU317 文献标志码：A DOI：1013465jcnkiiVS201715004Identification of a structural system and its excitation based on IMPSCO andmodified Newmark-fl algo
3、rithmMA Shenglan，JIANG Shaofei，CHEN Zhigang(College of Civil Engineering，Fuzhou University，Fuzhou 350108，China)Abstract： In vibrationbased struetural system identification(SI)of a civil structure，the commonly used SImethod，such as，the one based on Newmark调algorithm usually needs the known excitation
4、In fact，it is difficult tomeasure the actual excitation and full output responses of a structure systemHere，a modified Newmark-13 algorithm wasproposed through introducing polynomial fitting into Newmark召algorithm to effectively solve the drift phenomenon broughtby the error of SIFurthermore，combini
5、ng the modified Newmark够algorithm with the improved multiparticle swarm COevolution optimization(IMPSCO)，an identification approach of SI and excitation adopting responses of parts of a structuresystem was presentedFinally，numerical simulation for a sevenstory steel flame and its model tests were co
6、nducted toverify the effectiveness and feasibility of the proposed methodThe effects of noise and completeness of structural responsewere also discussedThe results showed that the proposed method can be used to effectively identify a structural systemand its excitation only adopting responses of par
7、ts of the structure system；it has stronger noise tolerance and robustnessKey words：unknown excitation；system identification；modified NewmarklB algorithm结构系统参数识别是结构健康监测的核心部分，因此国内外专家提出了很多系统识别方法圳。其中，时域下基于结构振动响应的识别方法受到了广泛关注，通过直接识别结构系统的物理参数，如刚度、质量或阻尼等，可以有效确定损伤位置、损伤程度和评估结构状态。这些方法如最小二乘法(Least Square，LS)口1、
8、卡尔曼滤波(Kalman Filter，KF)1和遗传算法基金项目：国家十二五科技支撑计划(2015BAKl4802)国家自然科学基金(51278127)收稿日期：20150922修改稿收到日期：20160628第一作者麻胜兰女，博士，1986年生通信作者姜绍飞男，闽江学者特聘教授，博士，博士生导师，1969年生Email：cejsffzueduan(Genetic Algorithm，GA)o 7 J等。但分析发现，这些方法大多仅适用于输入与输出都已知的情形，而实际工程中，由于环境的复杂性，一般难以测得结构实际受到的激励大小，甚至也难以测得完整的结构输出响应。因此，研究未知激励及部分实测输出
9、响应状态下的结构系统辨识方法，意义重大。国内外专家学者开展了一系列研究，取得了一定成果。Chen等旧。利用改进的迭代最小二乘法同时识别结构系统参数和未知输入力；此后，Chen等悖。又提出利用蚁群算法同时识别结构系统参数和未知激励，并通过试验验证了方法的有效性；此外，Sun等0。将改进Newmark算法与人工蜂群算法结合，来同时识别未知激励及结构系统参数。但研究发现Chen的方法万方数据第15期 麻胜兰等：基于IMPSCO和改进Newmark-卢算法的结构系统及激励辨识研究 23并不适用于仅有部分实测响应的情形；而Sun的方法通过将加速度积分成速度和位移，进而根据动力方程求得未知激励，但实际工程
10、中，加速度积分会由于“平移项”和“趋势项”的影响1使速度、位移和识别激励发生漂移，导致识别错误。因此，本文的目的在于研究一种未知激励及部分实测响应下，能同时识别激励与结构系统参数的时域下的辨识策略。时域下的结构系统识别通常需要与各种数值积分方法相结合，其中最常用的数值方法为Newmark够算法。1。”o，但是常规Newmark-fl算法需要已知激励，为此，通常做法是通过易测得的加速度积分得到速度和位移响应，再根据结构动力方程求得激励，然后进行系统参数识别。但如上所述，加速度积分存在趋势项误差，会使速度、位移和识别激励发生漂移，并导致系统参数识别错误。鉴于此，本文借鉴多项式去趋势的思想，对常规N
11、ewmark田算法进行改进，利用多项式拟合来修正加速度积分导致的识别激励的漂移现象；在此基础上，结合改进协同粒子群算法(IMPSCO)，提出了一种适用于未知激励及部分实测响应下的结构系统参数识别与激励辨识方法，并通过一7层框架数值仿真和试验研究验证了方法的有效性和可行性。1改进的Newmark毋算法本文对Newmark-fl算法进行了改进，具体步骤如下：步骤1 初始计算加速度、速度、位移和力：确定结构的质量、阻尼和刚度矩阵M，C，K，及时间步长At；令初始加速度、速度、位移和力分别为菇。=互孑，戈。=0，圣；。=0，菇：。=0和fo=0(上标“m”和“c”分别代表实测数据和计算数据)；计算等效
12、刚度詹=K+盖c+反去矿m(当采用平均加速度法时y=12、卢=14，当采用线性加速度法时y=12、3=16)同时令口=赢+3ffC；b=扣+t(嘉。)c；步骤2由每个时间步P计算第P+1个时间步的结构响应和激励：由实测加速度反推第P+1个时问步的速度和位移分别为莺辩=互+譬(菇+菇；并)、戈；2=戈+譬(圣+立；2)计算力+。=胍；并+瓯茹+甄鼎；按照传统Newmark-fl算法计算结构响应：石p+l=戈p+AxP、五P+1=莺P+莺P和互p+l=茧P+菇P；其中A戈。：岛二L粤世、K越，2盖缸，一分，+At 1一荔)菇，互，2反云矿一麦，一南菇，(则表示函数在下一时问步的值相对上一时间步值的
13、增量)识别激励更新圻+。=胍；m1)+P+1+Kx P+1；令量；笔=量p+1和戈辩=戈叫；步骤3 对下一个时间步进行循环令P=P+1，重复步骤2中的；步骤4 输出结构响应工、圣和菇以及识别激励，。步骤5基于最小二乘原理，对厂(t)进行多项式拟合即五。=o。+ait。步骤6识别激励修正=一。步骤7将Z作为已知激励，并重复执行步骤2中的步骤直到所有时间步循环结束。步骤8输出最终结构响应戈、莺和菇以及修正的识别激励Z。总体而言，改进的目的在于解决以下两个问题：(1)激励未知下，程序的执行问题及未知激励的识别；(2)修正可能发生的识别激励的漂移现象。对于前者，本文借鉴Sun的方法，利用已知的实测加速
14、度积分成速度和位移，进而根据动力方程求得未知激励。然而，Sun在计算下一时间步的速度响应时采用的上一时间步得到的速度和位移经常规Newmark-fl算法公式(类似本文所提算法步骤步骤2步骤中的)计算得到的数值加速度，这无疑会增加计算误差，因而本文将其改用为实测加速度(步骤2步骤中的一)；此外，本文在计算中加人了识别激励的更新(步骤2步骤中的)，以进一步减少可能存在的由于积分引起的识别误差。对于后者，本文利用多项式拟合去除可能发生的激励漂移现象，得到修正的识别激励，并将其作为已知值重新计算结构响应(步骤58)。2 系统识别与激励辨识策略如上所述，本文将IMPSCO与改进Newmark-fl算法结
15、合起来，提出了一种适用于未知激励及部分实测响应状态下的结构系统识别与激励辨识策略，其流程如图1所示；IMPSCO是对MPSCO的改进，通过对局部最优粒子的替换处理，提高了算法的稳定性和鲁棒性，具体过程可详见文献14。此外，考虑到加速度相较容易获得，本文基于加速度响应进行研究。因此，在激励位置、激励部位及其相邻部位加速度已知时，即可利用万方数据24 振动与冲击 2017年第36卷本策略识别结构系统参数和未知激励。具体步骤如下：(1)对结构系统参数进行编码并确定参数搜索范围，本文考虑质量已知情形，对于阻尼，则采用Ray!eigh阻尼模型，因此需要编码的结构系统参数0为刚度ki(i=1，2，n)和阻
16、尼比亭，(J=1，2)；(2)设置IMPSCO的基本算法参数，本文设置如下：子种群数量为挝=3；子种群大小为眠=30；进化代数m。=500；学习因子c，=c：=2；最差次数限值Ie=6；埘=04，09，线性递减，且当迭代次数大于45时，取04；(3)在各参数搜索范围内随机产生若干随机粒子(即不同参数值的组合)，并按第一节，利用改进Newmark-fl算法计算每个随机粒子下的结构响应及未知激励U)；(4)以实测加速度响应和计算所得加速度响应为基础建立适应度函数F，如式(1)所示，并计算每个随机粒子的适应度值；F= kTL一 (1)f 1(血。mea_ai,7)2i ，式中，amea和ac。m分别
17、为结构实测加速度响应和数值计算的加速度响应；k是结构测点数；f为加速度响应数据的个数。适应度计算值越大，识别结果越接近实际值。(5)替换掉适应度值最差的粒子，并对每个粒子按IMPSCO的设定进行粒子进化；(6)同样按第一节，利用改进Newmark-3算法计算每个进化粒子组合下的结构响应及未知激励(Z)，并重复步骤(4)和(5)；(7)重复步骤(6)，进行粒子循环迭代进化，直到达到设定的最大迭代次数，迭代结束，输出适应度值最高的粒子，即最优系统参数和未知激励(Z蛔)。仃再n系统参数编码并确定搜索范围(III度K；阻尼c或阻尼比)产生随机粒子(初始值)或迭代(粒子初始值迭代所得)嬲踹黟悻动力响应及
18、未知激励(再dcn)lo午l适应度函数H卢)卜卜kl输入结果(K；c或；-，iden)(堕壅)图1 结构系统识别与激励辨识策略流程图Fig1 Flow chart of SI and excitation identification3数值算例31结构模型图2为一7层剪切型框架结构，刚度为k，=k，=375 kNm；质量为m1=m6=378蚝和m7=331 kg；采用Rayleigh阻尼模型，前两阶阻尼比为f。=亭：=2。顶层受随机激励，采用常规Newmark书算法计算结构动力响应，并取加速度响应用于分析，采样频率1 000 Hz，采样时间5 s。32识别过程与结果本例中，质量已知，以结构各层
19、层间刚度及阻尼比为未知，进行编码识别，即0=k。，k：，k，f。，f：。参数搜索范围取理论值的05倍2倍，即ki1875，700kNm，i=1，2，7；f001，004=1，2。按图1流程，识别结构参数p及未知激励Z，改进Newmark-fl算法中多项式拟合阶数取n=4。识别时，仅考虑第2，3及57层加速度响应以模拟响应不完整情形。无噪声下，阻尼比的识别值分别为205和209，结果较为精确；刚度识别结果如表2所示；激励识别结果如图3所示，为简便，图中仅给出005s识别结果。图2结构模型Fig2 Structural model从表2可以发现，结构刚度的识别结果与理论值非常接近，最大误差仅为02
20、0；而从图3(a)可以发现，激励的识别值也与理论值非常吻合。此外，图4为未经多项式修正时的激励识别结果(即识别过程没有加入第一节中步骤58)，可以发现即使未利用多项式进行修正时，识别激励也没有明显发生漂移，这是由于本例为数值算例，在正向计算结构数值响应与反向识别结构物理参数时，实际上都是基于Newmark邶算法的基本假设进行的，因此识别中不会出现明显的激励漂移现象，因此，识别时多项式的修正值为0，但实际工程中，结构响应是实际测得的，并不完全满足这一假万方数据第15期 麻胜兰等：基于IMPSCO和改进Newmark-3算法的结构系统及激励辨识研究 25ts(a)无噪声(b)SNR=20 dB图3
21、激励识别结果Fig3 Identified results of excitation，+S图4未修正下激励识别结果Fig4 Identified results of excitation without correction设，因而容易产生趋势误差，导致识别激励发生漂移，这将在下一节试验中得到验证。以上分析表明，无噪声下，本文所提策略可以精确识别结构系统参数和未知激励。33比较与讨论为验证本文方法的适用性和鲁棒性，从噪声水平和输出响应的完整性两个方面进行分析讨论。(1)噪声影响按式(2)向原始加速度响应中添加噪声，信噪比(SNR)分另0为25 dB和20 dB。Y。=Y：(1+8R)SNR
22、=2019(1占)=2019(AsignalA。i。) (2)式中：Yi。和Yi分别代表原始位移响应信号和噪声污染后的信号；R是均值为0，偏差为1的正态分布随机数；s为噪声水平；A为信号幅值。同样，按图l流程进行识别，刚度和激励的识别结果如表2和图3所示。可以看出，SNR=25 dB和20dB时，刚度的最大识别误差也仅为081和220，而激励的识别结果也仍然非常精确。此外，SNR=25dB时，阻尼比的识别值分别为212和217；SNR=20 dB时，则分别为224和225；实质上Rayleigh阻尼比为2。可见阻尼比的识别误差与刚度识别一致，随噪声增加，其中阻尼比误差更大一些，这是由于本例为小
23、阻尼结构，阻尼对系统响应影响较小所致，但总体上其误差尚在可接受范围内。这些表明，本文所提策略具有较好的容噪性和鲁棒性。(2)输出响应的完整性为验证输出完整性对本文方法的影响，无噪声下，分别利用第57层、第2，3，57层和第17层的加速度响应进行识别。刚度和激励识别结果如表3和图5所示。17层加速度都参与识别时，阻尼比分别为208和202，最大误差400，而仅有57层的加速度响应参与识别时，其值分别为219和225，最大误差为1250。可见参与识别的加速度响应个数对于阻尼比的识别影响较大，但从表3可以发现，其对刚度的识别影响较小。当只有57层的加速度响应参与识别时，刚度的最大识别误差也仅增大到0
24、74，激励的识别结果也仍然非常精确(图6(b)。这些表明，虽然加速度响应的完整性会影响参数及激励识别结果，但影响较小，本文策略仍具有较高识别精度。一帅一徨从一。rllLHiUHiiLlllL) 糍一乩州一M一州一娜一M一帅一篷棚一雠删州一酬测帆一一“咿呷一一r一脚一一万方数据26 振动与冲击 2017年第36卷(b)第57层图5不同测点组合下的激励识别结果Fig5 The identified results of excitation with differentmeasured data4试验验证41试验模型如图6所示，为2跨1跨的7层钢框架缩尺模型，相应的数学简化模型如图2所示。梁、柱间
25、距均为200 mm。梁为空钢管，柱为薄钢板，且均为300 W级钢材。构件特性如表4所示，据此可求得结构各层质量为m1=m2=m6=378 kg，m7=331 kg，完整状态下结构各层理论刚度均为375 kNm。利用激振器在结构顶层施加随机激励，并利用力传感器记录激振力大小，采样频率5 000 Hz。在结构各层布置传感器以测量各层加速度响应。图6实验室钢框架模型Fig6 The experimental steel flame mode表4结构构件特性Tab4 Properties of the structural members42结构系统修正利用结构刚度理论值计算出结构的理论频率并与试验中
26、测得实测频率进行对比，如表5所示。表5结构固有频率Tab5 Natural frequency of structure从表5可以看出，理论频率与实测频率偏差较大，可见理论刚度值不符合结构实际条件，这是由于试验材料自身的缺陷、焊接的影响、试验设备的测量误差、试验人员人为的因素及测量时噪声的影响等多种因素一州一徭舢百镞州fr一晰一厂；：I5O一O5一帅一徨柑百粼州f万方数据第15期 麻胜兰等：基于IMPSCO和改进NewmarklB算法的结构系统及激励辨识研究 27引起的，因此需要对其进行修正。为此，在激励已知下，利用IMPSCO进行刚度识别与修正。首先对结构各层层间刚度及阻尼比进行编码，即0=
27、k。，k：，k，孝。，孝：。参数搜索范围为壳。刚度的识别结果如表6所示。而后利用识别刚度计算结构的固有频率，如表5所示，可以看出，识别频率与实测频率更加吻合，从而表明识别的刚度值更符合实际结构状况。因此，后文分析时以此刚度识别值作为1875，700kNm，江1，2，7；Eo，1，J=1，2。 激励未知下的刚度基准值，以便对比。表6不同测点组合下的刚度识别结果Tab6 The identified results of stiffness with different measured data43识别过程与结果以第2，3及57层的加速度响应作为已知输出响应，以模拟响应不完整情形，对各层刚度及阻
28、尼比进行编码，搜索范围与42节一致，并按图1流程进行识别，多项式阶数同样取n=4。阻尼比的识别结果分别为068和202；刚度和激励的结果分别如表6和图7(b)所示。图7(a)为未经多项式修正时的激励识别结果(第一节中步骤58步)。(a)未修正(不完整响应，第2， (b)修正(不完整响应，第2，3， (c)未修正(完整响应，第173，57层加速度) 57层加速度) 层加速度)图7激励识别结果Fig7 The identified results of excitation从表6可以看出，激励未知下，冈0度的识别结果与激励已知时的识别结果较为接近，最大误差仅为984；此外，从图7(a)可以发现，识
29、别激励由于加速度积分的趋势误差发生了明显漂移；但经本文所提策略修正后，激励的识别结果明显没有了漂移现象，且与实测值较为吻合，如图7(b)所示。以上分析表明，本文策略可以有效识别结构系统参数和未知激励，并消除激励的漂移现象。44比较与讨论为了验证本文方法的鲁棒性和有效性，对同样的结构模型进行了激励未知及完整响应状态(7层加速度都已知)下的系统辨识。首先对结构各层层间刚度及阻尼比进行编码，参数搜索范围与42节一致，同样按图2流程进行识别，多项式拟合阶数同样取n=4。阻尼比的识别值分别为052和186，与43节的识别结果较为接近，且都比较小；刚度和激励的结果如表6和图7(d)所示，图7(c)同样为未
30、经多项式修正时的激励识别结果。从表6可以发现，相较完整输出响应情形，仅有部O 01 02 03 04tS(d)修正(完整响应，第17层加速度)分实测响应时，刚度的最大误差从585增加到984，可见识别误差虽有所增大，但仍有较高的识别精度；此外，从图7(b)和图7(d)对比发现，相较完整输出响应情形，部分实测响应下激励识别值与激励实测值的吻合度虽有所不如，但总体上仍然较为吻合。以上对比表明，虽然加速度响应的不完整性对参数和激励的识别有所影响，但本文所提策略仍表现出较高的识别精度和鲁棒性，能满足工程需要。5 结论本文通过研究，得出了以下结论：(1)通过改进常规Newmark-3算法并结合IMPSC
31、O，提出了系统辨识策略，有效解决了仅有部分实测输出响应下的结构系统识别和激励辨识问题。(2)响应完整性确会影响识别结果，但本文策略对于仅有部分实测输出响应下的系统参数识别仍具有较高的识别精度和较强的鲁棒性。(3)研究表明，所提改进Newmark-3算法确实能有效消除识别激励的漂移现象，保证了系统识别结果蒯一黼坐憾熬触一n1n飘崞Il叶I川W4O 薯蜘一。万方数据振动与冲击 2017年第36卷的准确性。数值模拟和试验研究表明，所提策略能有效识别结构系统参数和激励，初步验证了该方法的可行性和有效性，且具有较好的容噪性，但仍需要更多的试验和工程实践来进一步检验。此外，本文策略仅适用于激励位置已知的情
32、形，对于激励位置未知情形，还有待进一步研究和扩展。参考文献1张延年，范鹤，董锦坤基础隔震结构的系统参数优化设计J沈阳建筑大学学报，2006，22(4)：529-533ZHANG Yannian，FAN He，DONG JinkunParameteroptimization of baseisolated structure f J Journal ofShenyang Construction University，2006，22(4)：529-5332TANG H，ZHANG W，XIE L，et a1Multistage approachfor structural damage iden
33、tification using particle swalTnoptimizationJSmart Structures and Systems，2013，1 1(1)：69863WANG F G，LING X Z，XU X，et al。Structural stiffnessidentification based on the extended kalman filter researchJAbstract and Applied Analysis，2014(2)：267-2904 YANG J N，HUANG H，LIN SSequentiM nonlinear leastsquare
34、 estimation for damage identification of structuresJInternational Journal of NonLinear Mechanics，2006，41(1)：1241405LIN J W，BEqTI R，SMYTH A W，et a1On1ineidentification of nonlinear hysteretic structural systems usinga variable trace approachJEarthquake EngineeringStructural Dynamics，2001，30(9)：127913
35、036HOSHIYA M，SAITO EStructural identification by extendedKalman filterJJournal of Engineering Mechanics，1984，(上接第15页)GE D D，ZHU H P，WANG D S，et a1Seismic response analysis of damper-connected adjacent structures with stochas-tic parametersJJournal of Zhejiang Universityenee A，2010，11(6)：402-414RICHA
36、RDSON A，WALSH K K。ABDULLAH M MClosed-form equations for coupling linear structures using stiffnessand damping elementsJStructural Control&Health Monitoting，2013，20(3)：259281OK S Y，SONG J，PARK K SOptimal design of hystereticdampers connecting adjacent structures using multiobjectivegenetic algorithm
37、and stochastic linearization methodJEn-gineeting Structures，2008，30(5)：12401249BASlLI M，ANGELIS M DOptimal passive control of adjacent structures interconnected with nonlinear hysteretic devicesJJournal of SoundVibration，2007，301(12)：106一】25110(12)：175717707PERRY M J，KOH C G，CHOO Y SModified genetic
38、algorithm strategy for structural identificationJComputersStructures，2006，84(89)：529-5408CHEN J，LI JSimuhaneous identification of strueturalparameters and input time history from outputonlymeasurementsJComputational Mechanics，2004，33(5)：3653749CHEN J，CHEN X，LIU WComplete inverse method usingant colo
39、ny optimization algorithm for structural parameters andexcitation identification from output only measurementsJMathematical Problems in Engineering，2014，18(1)：24926910SUN H，BETTI RSimultaneous identification of structuralparameters and dynamic input with incomplete output-onlymeasurementsJStructural
40、 Control and Health Monitoring，2014，21：868-8891 1JUNHEE K，KIYOUNG K，HOON SAutonomous dynamicdisplacement estimation from data fusion of acceleration andintermittent displacement measurementsJ MechanicalSystems and Signal Processing，2014，42(12)：194-20512MAUNG T S，CHEN H P，AIANI A MReal time structura
41、ldamage assessment from vibration measurementsJJournalof Physics：Conference Series，201 1，305(1)：12096-1210313SHIN S B，OH S HDamage assessment of shear buildingsby synchronous estimation of stiffness and damping usingmeasured accelerationJSmart Structures and Systems，2007，3(3)：245-26114JIANG S F，Wu S
42、 Y，DONG L QA time-domain structuraldamage detection method based on improved multiparticleswarm coevolution optimization algorithmJ，MathematicalProblems in Engineering，2014，44(1)：77-851516171819NG C L，XU Y LSemiactive control of a building complexwith variable friction dampersJEngineering Structures
43、，2007，29(6)：12091225PARK K S，OK S YHybrid control approach for seismiccoupling of two similar adjacent structuresJJournal ofSound&Vibration，2015，349：117BHARTI S D，DUMNE S M，SHRIMALI M KSeismic re-sponse analysis of adjacent buildings connected with MRdampersJEngineering Structures，2010，32(8)：21222133YING Z G，NI Y Q，KO J MStochastic optimal couplingcontrol of adjacent building structuresJComputersStructures，2003，81(3031)：27752787KIM J，RYU J，CHUNG LSeismic performance of structures connected by viscoelastic dampersJEngineeringStructures，2006，28(2)：183-1951J1JJ1j屹nHrlrlrlrl万方数据