首页> 中国专利> 一种基于液体闪烁探测器的快中子多重性测量分析方法

一种基于液体闪烁探测器的快中子多重性测量分析方法

摘要

本发明涉及一种基于液体闪烁探测器的快中子多重性测量分析方法,包括以下步骤:建立基于液体闪烁探测器的快中子多重性测量的数学模型,该数学模型为关联各重符合计数率和样品未知参数的方程组;采用多个液体闪烁探测器探测样品源中子裂变过程,得到中子计数率,进而得到各重符合计数率;将测得的各重符合计数率代入所述方程组,求解该方程组,获得样品未知参数(自发或诱发裂变率等)。采用本发明提供的方法能给出被测样品稳定可靠的分析结果。

著录项

法律信息

  • 法律状态公告日

    法律状态信息

    法律状态

  • 2017-03-22

    授权

    授权

  • 2015-07-01

    实质审查的生效 IPC(主分类):G01T3/06 申请日:20150202

    实质审查的生效

  • 2015-06-03

    公开

    公开

说明书

技术领域

本发明属于核燃料非破坏分析(Non-Destructive Assay,以下简称NDA) 技术领域,具体涉及一种基于液体闪烁探测器的快中子多重性测量分析方 法。

背景技术

核材料(主要是铀和钚)的NDA测量方法一般可分为中子测量法、γ 射线测量法、量热法等。由于中子具有较强的穿透性,使其测量质量较大 的核材料时相比其他方法具有显著的优势,某些情况下也是唯一可行的测 量技术。与γ测量不同,裂变中子能谱是连续谱,因此中子测量通常主要 是测量中子的计数率以及时间关联符合谱等确定样品相关参数。

中子测量技术主要用于含Pu/U物料的精确定量分析、常规样品分析、 库房盘存、Pu/U生产线闭合衡算等方面。另外,中子测量技术还可以用于 军控核查等相关领域。

中子测量方法主要有以下三种:总中子测量方法、符合中子测量方法 和多重性测量方法。相比总中子测量和符合中子测量方法,多重性测量适 用范围更广,并且在测量精度上有明显优势。

钚材料自发裂变中子产额较高(自发裂变率:473.5/(s*g)),而铀 材料的自发裂变中子产额很低,因此核保障领域中中子的测量始于钚样品 的测量研究。钚的中子测量一般为被动法,无需外加激发源;而铀材料的 测量一般为主动法,需要外加质询源。

多重性测量是通过测量数个中子多重性测量参数(总中子计数率,二 重符合中子计数率,三重符合中子计数率),然后求解描述这些测量参数 和样品未知参数(对钚,为自发裂变率F,泄露自增殖系数M,(α,n)反 应率α;对铀,为诱发裂变率F,泄露自增殖系数M)之间关系的方程组。 在求解出F值和M值后,即可确定240Pu的等效质量(对铀样品,则为235U 的质量)。

建立描述样品未知参数与测量参数之间关系的方程组是中子多重性测 量方法数据分析方法的核心。而该方程组是在一定的假设近似下,进行理 论推导得到的。

传统中子多重性测量技术采用3He中子探测器。3He对γ射线不敏感, 对热中子具有很大的俘获截面,而裂变中子以MeV量级的快中子为主,因 此测量装置都用很厚的聚乙烯材料包裹,以达到慢化快中子的目的。

尽管3He测量技术目前已经比较成熟,但近年来由于3He全球供应的短 缺(核武器生产大幅削减)以及价格上涨,导致设备成本较为高昂,这限 制了其大规模的推广应用,因此有必要研究其替代探测器以及相应分析技 术。而且3He测量的是慢化后的热中子,其偶然符合本底较高,这在一定 程度上影响了其测量精度。

液体闪烁体探测器价格远低于3He,并且能直接测量快中子,能最大程 度上保留时间关联中子之间的时间关联信息,同时偶然符合本底极低。另 外,近年来,与闪烁体探测器相关的n/γ甄别技术以及高时间分辨率的时 间计数器等都有较快发展。基于以上这些原因,液体闪烁体探测器成为一 种十分有发展前景的新型中子多重性计数器。

目前国内外在液体闪烁探测器中子多重性测量方面的研究还处于原型 装置的Monte Carlo模拟计算、设计和改进中,如何通过测量参数求解样 品未知参数的数据分析方法还是个亟待解决的难题。

发明内容

针对现有技术中存在的缺陷,本发明的目的是提供一种基于液体闪烁 探测器的快中子多重性测量分析方法,采用该方法能得到核燃料被测样品 未知参数的稳定可靠的分析结果。

为达到以上目的,本发明采用的技术方案是:一种基于液体闪烁探测 器的快中子多重性测量分析方法,包括以下步骤:建立基于液体闪烁探测 器的快中子多重性测量的数学模型,该数学模型为关联各重符合计数率和 样品未知参数的方程组;采用多个液体闪烁探测器探测样品源中子裂变过 程,得到中子计数率,进而得到各重符合计数率;将测得的各重符合计数 率代入所述方程组,求解该方程组,获得样品未知参数。

进一步,所述方程组的表达式为:

R1×K=F·P1×vmax·Qvmax×vmax·Tvmax×K---(1),

其中,式(1)中,

K表示液体闪烁探测器的个数,

F表示样品源的中子自/诱发裂变率,

R1×K表示一重、二重、……K重符合计数率所组成的大小为1xK的行向 量,

表示自增殖过程之后从样品中出射的中子数概率分布所组成的大 小为1xνmax'的行向量,νmax'表示向量所对应的最大出射中子数目,其中,

向量的每一项Pn的表达式为:

Pn=a″n·q2+a'n·q+an    (2),

式(2)中的q表示裂变中子引发下一次次级裂变的概率;a″n、a'n、an是 Pn与q所形成的二次函数的系数,n表示自增殖过程后的中子数, n=1,2,……,νmax';

式(1)中,是大小为νmax'×νmax'的转移矩阵,转移矩阵中 每一项Q(i,j)的表达式为:

Q(i,j)=Cijϵj(1-ϵ)(i-j),i>=j0,i<j---(3),

式(3)中的Q(i,j)表示i个中子从样品中出射之后能探测到j个中子 的概率,表示从i个中子中任选j个的可能选择数,ε为 总的探测效率;

式(1)中,是大小为νmax'×K的转移矩阵,转移矩阵的每一 项的表达式分别为:

Ti1=Σj=1K(ϵjϵ)i,i1Ti2=Σj1=1,j2=1,j1j2j1=K,j2=K(ϵj1ϵ+ϵj2ϵ)i-(K-1)Ti1,i2Ti3=Σj1=1,j2=1,j3=1,j1j2,j1j3,j2j3j1=K,j2=K,j3=K,(ϵj1ϵ+ϵj2ϵ+ϵj3ϵ)i-(K-1)Ti1-(K-2)Ti2,i3···Tis=Σj1=1,j2=1,..,js=1,j1j2,j1j3,...,j1js,j2j3,...,j1js,...,js-1jsj1=K,j2=K,..,js=K(ϵj1ϵ+ϵj2ϵ+...+ϵjsϵ)i-Σj=1j=s-1(K-j)Tij,sivmax,2sKTis=0,i<s---(4),

式(4)中的Tis代表在探测到i个中子的情况下产生s个信号的概 率,εs(s=1,2,…,K)代表第s个探测器的探测效率。

再进一步,液体闪烁探测器的个数K为4时,所述方程组的表达式为:

R1×4=F·P1×vmax·Qvmax×vmax·Tvmax×4---(5),

其中,式(5)中,的表达式为:

Ti1=(ϵ1ϵ)i+(ϵ2ϵ)i+(ϵ3ϵ)i+(ϵ4ϵ)i,i1Ti2=(ϵ1ϵ+ϵ2ϵ)i+(ϵ1ϵ+ϵ3ϵ)i+(ϵ1ϵ+ϵ4ϵ)i+(ϵ2ϵ+ϵ3ϵ)i+(ϵ2ϵ+ϵ4ϵ)i+(ϵ3ϵ+ϵ4ϵ)i-3Ti1,i2Ti3=(ϵ1ϵ+ϵ2ϵ+ϵ3ϵ)i+(ϵ1ϵ+ϵ2ϵ+ϵ4ϵ)i+(ϵ1ϵ+ϵ3ϵ+ϵ4ϵ)i+(ϵ2ϵ+ϵ3ϵ+ϵ4ϵ)i-2Ti2-3Ti1,i3Ti4=1-Ti1-Ti2-Ti3,i4Tij=0,i<j---(6),

式(6)中,ε1、ε2、ε3和ε4分别代表第1、2、3、4个液体闪烁探 测器的探测效率。

更进一步,使用测得的一重、二重、三重、四重符合计数率中的任意 三个求解方程组,通过消去式(5)中的F,得到:

f1(q,ϵ)=0f2(q,ϵ)=0---(7),

然后记f=f12+f22,令

f=f12+f22=0    (8),

用迭代法求解方程式(8),得到q、ε,进而根据式(5)求得自/诱 发裂变率F。

进一步,所述方法还包括根据式(9)计算样品源的泄漏自增殖系数M (以下简称“自增殖系数”),

M=1-q1-q·v---(9),

式(9)中,v表示每次次级裂变事件发射的裂变中子平均数。例如, 对235U,该值通常可取为2.64。

本方法在一定假设和近似的基础上,建立了一个描述液体闪烁探测器 探测裂变快中子主要过程的数学模型,该模型以矩阵乘法的形式建立了一 系列将多重性测量参数和样品未知参数联系起来的方程组,通过求解方程 组得到样品的未知参数自发裂变率或诱发裂变率(简写为“自/诱发裂变 率”)F和自增殖系数M,以确定样品中的核源的等效质量。采用本发明提 供的方法能给出被测样品稳定可靠的分析结果。

附图说明

图1近似地示出了裂变中子的测量过程;

图2是用于探测中子裂变过程的四个液体闪烁探测器所组成的探测器 系统的结构示意图。

具体实施方式

下面结合具体实施方式对本发明作进一步描述。

本发明提供的一种基于液体闪烁探测器的快中子多重性测量分析方 法,包括以下步骤:建立基于液体闪烁探测器的快中子多重性测量的数学 模型,该数学模型为关联各重符合计数率和样品未知参数的方程组;采用 多个液体闪烁探测器探测样品源中子裂变过程,得到中子计数率,进而得 到各重符合计数率;将测得的各重符合计数率代入所述方程组,求解该方 程组,获得样品未知参数。

以下具体描述本方法中数学模型的建立过程。

当一个重核发生裂变时,它以Pν概率产生ν个中子(ν最小为0,最大 值视具体核素而定,通常为6~9)。如果这些裂变中子引发其他核发生裂 变(自增殖过程),可能会产生更多中子。通常,将自发裂变或质询源中 子引起的诱发裂变被称为初级裂变事件,而裂变中子诱发的诱发裂变被称 为次级裂变事件。如果裂变中子被探测器探测到,并且其沉积的能量超过 阈值,就会产生一个中子信号。由于来自于同一个裂变事件的关联中子到 达探测器的时间间隔极短,如果其中多个中子进入同一个探测器,这些中 子将无法区分而成为一个中子信号(探测器的脉冲堆积)。本发明使用转 移矩阵来描述上述过程。

在提出本方法的数学模型时,本方法首先针对裂变中子测量过程做了 以下一些基本的假设和近似:

a、不存在裂变中子进入探测器后又被反射回样品的情形;

b、次级裂变事件发射的裂变中子具有和初级裂变事件发射的裂变中子 相同的能谱;

c、不存在俘获中子却不发生裂变的情形;

d、来自于同一个裂变链(包括初级裂变以及后续次级裂变)的中子如 果被探测到,都将处符合门宽内;

e、来自其他裂变事件的偶然符合计数可以忽略。

需要注意的是,以上假设和近似在通常的测量场景下都能满足,但不 排除在部分极端条件下,实际情景和以上假设及近似可能存在较大偏离从 而导致本方法的分析结果偏差明显变大的可能性。

现以235U源的裂变过程为例说明本方法中裂变中子的主要测量过程。 如图1所示,235U源首先在质询中子的激发下发生初级诱发裂变,产生第一 代中子,其具有存在初始裂变中子数概率分布P0,第一代中子引发自增殖 过程后产生更多中子,其具有中子数概率分布P,用转移矩阵Q描述探测 效率的影响,则探测器探测到的中子数概率分布为D=P·Q;用转移矩阵T描 述探测器脉冲堆积的影响,则每一次裂变事件对应的探测到的N重符合计 数的概率分布为:

r=P·Q·T    (1′),

上式(1′)中,r是一个1×K(K是液体闪烁探测器的个数)的行向量, 其元素分别对应一、二、三、…、K重符合计数的概率。若235U源每秒初级 裂变次数为F,则最后测量得到的各重符合计数率R满足下式:

R=F·r=F·P1×vmax·Qvmax·Tvmax×K---(1),

上式(1)中,F表示样品源的中子自/诱发裂变率,R1×K表示一重、二 重、……K重符合计数率所组成的大小为1xK的行向量,表示自增殖 过程之后从样品中出射的中子数概率分布所组成的大小为1xνmax'的行向 量,νmax'表示向量所对应的最大出射中子数目,其中,

向量的每一项Pn的表达式为:

Pn=a″n·q2+a'n·q+an    (2),

式(2)中的q表示裂变中子引发下一次次级裂变的平均概率;a″n、a'n、 an是Pn与q所形成的二次函数的系数,n表示自增殖过程后的中子数, n=1,2,……,νmax';

式(1)中,是大小为νmax'×νmax'的转移矩阵,转移矩阵中 每一项Q(i,j)的表达式为:

Q(i,j)=Cijϵj(1-ϵ)(i-j),i>=j0,i<j---(3),

式(3)中的Q(i,j)表示i个中子从样品中出射之后能探测到j个中子 的概率,表示从i个中子中任选j个的可能选择数,ε为 总的探测效率(若效率均匀性较好,ε可以预先刻度;否则需要作为未知 参数进行求解);

式(1)中,是大小为νmax'×K的转移矩阵,转移矩阵的每一 项的表达式分别为:

Ti1=Σj=1K(ϵjϵ)i,i1Ti2=Σj1=1,j2=1,j1j2j1=K,j2=K(ϵj1ϵ+ϵj2ϵ)i-(K-1)Ti1,i2Ti3=Σj1=1,j2=1,j3=1,j1j2,j1j3,j2j3j1=K,j2=K,j3=K,(ϵj1ϵ+ϵj2ϵ+ϵj3ϵ)i-(K-1)Ti1-(K-2)Ti2,i3···Tis=Σj1=1,j2=1,..,js=1,j1j2,j1j3,...,j1js,j2j3,...,j1js,...,js-1jsj1=K,j2=K,..,js=K(ϵj1ϵ+ϵj2ϵ+...+ϵjsϵ)i-Σj=1j=s-1(K-j)Tij,sivmax,2sKTis=0,i<s---(4),

式(4)中的Tis代表在探测到i个中子的情况下产生s个信号的 概率,εs(j=1,2,…,K)代表第s个探测器的探测效率。

本文中,“N重”符合计数的定义有别于传统3He测量中的定义。本文 的“N重”符合计数的定义为:如果在符合门宽(本工作为100ns,实际可 以更低)内记录到且只记录到N个中子信号,那么这N个中子信号将被视 为一个N重事件,N重符合计数增加1。

当液体闪烁探测器的个数为4时,测量到的各重符合净计数率(计数 率扣除本底后即为净计数率)可用大小为1×4的行向量R来表示,例如R(2) 表示二重符合净计数率。这样,本发明数学模型中的方程组的表达式可表 示为:

R1×4=F·P1×vmax·Qvmax×vmax·Tvmax×4---(5),

式(5)中,R1×4表示一重、二重、三重、四重符合计数率R(1)、R(2)、 R(3)、R(4)所组成的大小为1x4的行向量。

以下详细上述测量过程所涉及的其他参数的计算方法。

(ⅰ)初始裂变中子数概率分布——行向量P0

可以用行向量P0来代表初始裂变中子数概率分布:

P0=(P00,P01,P02,P03,...,P0vmax)

对于235U,上式中的νmax等于7;对于252Cf,上式中的νmax等于9。不同核 素的初始裂变中子数概率分布可从已有的文献和核数据中查询。

252Cf为例,其P0值为:

P0=[0.0021,0.0247,0.1229,0.2714,0.3076,0.1877,0.0677,0.0141, 0.0017,0.0001]。

例如,P03=0.2714,表示一次诱发裂变产生3个中子的概率是0.2714。

(ⅱ)自增殖过程之后从样品中出射的中子数概率分布——行向量P

自增殖过程在初级裂变事件之后极短的时间内发生,并且次级裂变和 初级裂变事件产生的中子能谱相近,所以初级裂变中子以及相应的次级裂 变中子可以认为来自于同一次裂变事件。

可以将自增殖过程之后从样品中出射的中子数概率分布概率记为行向 量P。P仅由初始中子数概率分布P0和参数q决定,因此可以建立一个数据 库描述不同核素的参数q和P的对应关系,见式(2)。

本工作中所建立的数据库中,最多仅考虑了第4代中子(意即增殖链 最长为4,自发裂变或质询源中子引起的诱发裂变产生的中子称为第一代 中子),与此同时,每一代中子中,最多4个中子可能引发下一级次级裂 变。

表1和表2分别给出了钚材料中239Pu和铀材料中235U的q-P对应关系 曲线的系数值。由于P1、…P14之和加上产生0个中子的概率已在0.9999 以上,因此高于14个出射中子的情形因概率极小而可以忽略不计(即vmax′取14)。由于产生0个中子的裂变事件对后续测量结果无任何影响(无 中子产生),所以在表1和表2以及相关计算中都不会将该类事件考虑在 内。

表1

表2

对于252Cf,由于其本身不存在自增殖过程,因此出射中子数概率分布 P和其初始中子数概率分布P0相同。对于其他裂变核素,也可以仿照本方 法计算出类似的q-P对应关系。

(ⅲ)探测效率与转移矩阵Q

记探测到的中子数概率分布为行向量D,D由向量P和总探测效率ε 决定:

Dm=Σn=mvmaxPnCnmem(1-e)n-m,........1mvmax---(10),

式(10)中,表示从n个中子中任选m个的可能选择数, v′max是向量P所对应的最大出射中子数目。

式(10)用下述矩阵乘法来表达比较简便:

D=P·Q    (11),

式(11)中的转移矩阵Q由总探测效率ε决定,转移矩阵Q中每一项 Q(i,j)的表达式见式(3)。

这样,大小为νmax'×νmax'的转移矩阵可表达为:

Qvmax×vmax=C11ϵ0···0C21ϵ(1-ϵ)C22ϵ2···0C31ϵ(1-ϵ)2C32ϵ2(1-ϵ)···0············Cvmax1ϵ(1-ϵ)(vmax-1)Cvmax2ϵ2(1-ϵ)(vmax-2)···Cvmaxvmaxϵvmax---(12),

(ⅳ)记录信号数分布与转移矩阵T

如前所述,时间间隔极短的多个中子进入同一个探测器时,将只能给 出一个中子信号(脉冲堆积)。类似之前的推导步骤,对于四探测器系统, 用于描述这种作用过程的转移矩阵T的表达式如下:

Ti1=(ϵ1ϵ)i+(ϵ2ϵ)i+(ϵ3ϵ)i+(ϵ4ϵ)i,i1Ti2=(ϵ1ϵ+ϵ2ϵ)i+(ϵ1ϵ+ϵ3ϵ)i+(ϵ1ϵ+ϵ4ϵ)i+(ϵ2ϵ+ϵ3ϵ)i+(ϵ2ϵ+ϵ4ϵ)i+(ϵ3ϵ+ϵ4ϵ)i-3Ti1,i2Ti3=(ϵ1ϵ+ϵ2ϵ+ϵ3ϵ)i+(ϵ1ϵ+ϵ2ϵ+ϵ4ϵ)i+(ϵ1ϵ+ϵ3ϵ+ϵ4ϵ)i+(ϵ2ϵ+ϵ3ϵ+ϵ4ϵ)i-2Ti2-3Ti1,i3Ti4=1-Ti1-Ti2-Ti3,i4Tij=0,i<j---(6),

式(6)中,ε1、ε2、ε3和ε4分别代表第1、2、3、4个液体闪烁探 测器的探测效率。

由于转移矩阵T由各个探测器的效率和总探测效率的比值决定,因此 其值可以通过各探测器的净计数计算得到,因此可视为能测量得到的已知 量,而不是需求解的未知量。

对于四探测器实验平台,由于紧凑的探测器结构设计,使得探测效率 均匀的空间有限,且探测效率和样品质量、大小等都有关,因此探测效率 难以作为可刻度的系统参数,而需要作为未知量求解。其他参数由于主要 与样品有关,所以也无法刻度,需要求解。

本发明中,可以使用测得的一重、二重、三重、四重符合计数率中的 任意三个来求解方程组,通过消去式(5)中的F,得到:

f1(q,ϵ)=0f2(q,ϵ)=0---(7),

然后记f=f12+f22,令

f=f12+f22=0    (8),

用迭代法求解方程式(8),得到q、ε,进而根据式(5)求得自/诱 发裂变率F。

此外,本发明的方法还包括计算另一样品未知参数——样品源的泄漏 自增殖系数M。泄漏自增殖系数M的表达式为:

M=1-q1-q·v---(9),

式(9)中,v表示每次次级裂变事件发射的裂变中子平均数。在q值 确定后,自增殖系数M也能确定。

实施例

本实施例采用自行搭建的测量装置对252Cf源进行了实验测量和数据分 析,以检验本方法的可行性。

(1)测量装置

测量仪器包括四个BC501A液体闪烁体探测器1(见图2)、一个n/γ 甄别模块MPD-4(4路输入和四路输出,实验时只输出中子信号)、一个多 通道的信号计时器MCS6A(6路输入和6路输出,用于记录中子信号胡到达 时间及对应探测器)和一个用于数据存储的电脑。

其中,四个BC501A液体闪烁体探测器1的水平轴线在同一平面内,四 个BC501A液体闪烁体探测器的一端围成样品腔2,样品腔2的上、下部对 应设有上、下盖3-1、3-2。

探测器的阈值都设置到约为1/4Cs(对应中子能量为0.76MeV)。为减少 γ射线的干扰,在每个探测器之前放置厚约5mm的铅片。

使用一个标称值已知的252Cf中子源(体积较小,接近点源)进行分析方 法可行性的检验。实验测量时将252Cf源放置于四个探测器围成的样品腔的 大约中心位置。

在测量数据预处理(设置符合门宽为100ns,将时间间隔在该门宽内 的N个中子信号记录为一个N重符合信号。对四个探测器组成的探测系统, N=1,2,3,4)后,可以得到一重符合计数率、二重符合计数率、三重符合计 数率以及四重符合计数率。

测量数据列于表3。

表3

数据 计数率(s-1) 1号探测器 6.646 2号探测器 7.215 3号探测器 5.393 4号探测器 5.507 一重符合计数率 20.195 二重符合计数率 2.155 三重符合计数率 0.08471 四重符合计数率 9.746e-4

(2)实验数据分析

由于252Cf源不存在自增殖过程,所以q值为0,对应M值为1,需要 求解的未知参数只有两个(F、ε)。在求解两个未知数时,只需要两个方 程。由于四重计数率较低,在约48小时的测量时间内,只有百余个计数, 统计较差,因此不予采用。在不考虑串扰时,选取一重、二重、三重计数 率三者中的两个即可求解出F和ε。

计算结果列于表4。表4中的组合“12”表示一重和二重净计数率被 用于方程组求解,其他组合意义类似。

表4

组合 F ε 12 81.20 0.0836 13 80.71 0.0843 23 80.21 0.0850

252Cf源的自发裂变率F标称值为77.7±2.7次/s(其强度按其半衰期 修正到实验时的强度)。

从表4可以看到,选取一重、二重、三重计数率中的任意两个求解, 给出的三组结果相近,但仍显示出一定差异。由于一重计数中可能包含由 于n/γ甄别时的误甄等原因带来的“错误”信号的干扰,所以采用二重和 三重计数率进行求解(组合23)给出的结果应当是更准确的。表4的结果 表明,自发裂变率F的分析值(80.21)与其标称值(77.7)十分接近。考 虑到其标称值的不确定度,可认为本方法给出的分析结果是合理可靠的。

上述实施例只是对本发明的举例说明,本发明也可以以其它的特定方 式或其它的特定形式实施,而不偏离本发明的要旨或本质特征。因此,描 述的实施方式从任何方面来看均应视为说明性而非限定性的。本发明的范 围应由附加的权利要求说明,任何与权利要求的意图和范围等效的变化也 应包含在本发明的范围内。

去获取专利,查看全文>

相似文献

  • 专利
  • 中文文献
  • 外文文献
获取专利

客服邮箱:kefu@zhangqiaokeyan.com

京公网安备:11010802029741号 ICP备案号:京ICP备15016152号-6 六维联合信息科技 (北京) 有限公司©版权所有
  • 客服微信

  • 服务号