一种中性束注入条件下核聚变产物能谱的计算方法和系统

未命名 07-27 阅读:142 评论:0


1.本发明涉及核工程与核技术技术领域,更具体的,涉及一种中性束注入条件下核聚变产物能谱的计算方法和系统。


背景技术:

2.核聚变能源作为一种非常有价值的未来能源,具有能量密度大、原料储存丰富和安全性优异等诸多优点。在能源消耗量逐渐增加的现代社会,随着科学技术的发展,人类对能源的需求越来越大。与此同时,传统化石能源对环境的破坏问题也越来越被人们重视。
3.基于以上原因,人类对核聚变能源进行相关科学研究,从而推动核聚变能源的商用化是非常有必要的。在目前的核聚变研究中,托卡马克装置作为被人们认为可行性较大的磁约束热核聚变实验装置正在被广泛用于核聚变的研究中。在托卡马克装置中,氘氚聚变会产生高能的中子和α粒子作为聚变反应的产物。两种聚变反应产物都带有巨大的能量,这些能量可以被转化成热能,随后被进一步转化成电能,实现核聚变能源发电的完整流程,因此对聚变反应产物的能谱进行计算是十分有必要的。托卡马克装置加热过程的初始阶段一般会采用欧姆加热的方法,但是欧姆加热方法具有加热效率在超过一定温度后出现下降、加热过程中对离子的加热效率偏低和不适用于聚变堆稳态运行的脉冲式工作等缺点,所以在使用欧姆加热方法之外还需要使用中性束注入作为辅助加热的手段。传统计算方法缺乏对中性束注入情况下的聚变产物能谱的研究,而且由于实验条件的限制,对于中性束注入条件下氘氚聚变能谱的研究成果较少。因此,研究一种新的适用于中性束注入条件下和聚变产物能谱的计算方法是十分有必要的。
4.关于中性束注入条件下核聚变产物能谱的计算方法具有着重要的研究价值和广阔的未来研究前景,对于后续基于托卡马克装置进行核聚变相关内容的研究有着重要的参考价值和实用价值。基于该方法得到的中性束注入条件下核聚变产物能谱的计算结果将能够作为后续托卡马克装置中核聚变的反应条件改良及稳态运行的重要参考依据。


技术实现要素:

5.本发明为了解决以上现有技术存在的不足与缺陷的问题,提供了一种中性束注入条件下核聚变产物能谱的计算方法和系统。
6.为实现上述本发明目的,采用的技术方案如下:
7.一种中性束注入条件下核聚变产物能谱的计算方法,所述的方法包括步骤如下:
8.根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;
9.根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;
10.两种反应物核聚变后得到核聚变产物,假设核聚变产物的运动方向沿任意给定的方向,根据得到的第一相对运动动能计算核聚变产物的第一动能;
11.通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;
12.按能量箱进行数据统计整理,计算得到核聚变产物的能谱。
13.优选地,所述的第一反应物、第二反应物在中性束注入条件下发生核聚变,所述的中性束注入条件认为中性束内部的中性束粒子与中性束粒子之间不发生核聚变反应,束粒子只有背景反应物之间发生核聚变反应。
14.进一步地,所述背景反应物是托卡马克装置中的背景等离子体,背景反应物遵循麦克斯韦-玻尔兹曼分布,同时假设托卡马克装置中的等离子体是各向同性。
15.进一步地,所述中性束粒子是单能粒子束、或给定速度分布的粒子束。
16.优选地,计算核聚变产物的第一动能,具体如下:
17.根据两种反应物和核聚变产物之间的质量亏损,计算得到反应后系统增加的第二动能;
18.根据动量守恒定律得到质心系下,计算两种核聚变产物组成系统的第二相对运动动能;
19.将第二相对运动动能根据两种核聚变产物的质量比分配,得到两种核聚变产物在质心系中的第一动能。
20.优选地,根据反应物粒子的分布规律得到两种反应物在实验室参考系下的归一化速度分布函数;
21.将粒子动能的概率密度函数和微观粒子运动速率的麦克斯韦-玻尔兹曼分布的标准表达式代入反应率表达式,在一定条件下进行代换计算得到最终的能谱计算式,通过对能谱计算式进行积分得到反应率计算式。
22.进一步地,在反应率计算式中通过引入狄拉克函数来避免同种粒子发生反应时的重复计数。
23.进一步地,基于r矩阵的参数化反应截面,具体通过引入伽莫夫常数bg和pad
é
多项式s(er)计算得到反应的微分截面,基于反应的微分截面,计算得到反应速率系数。
24.一种中性束注入条件下核聚变产物能谱的计算系统,包括:
25.随机抽样模块,根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;
26.第一相对运动动能计算模块,用于根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;
27.第一动能计算模块,用于根据得到的第一相对运动动能计算核聚变产物的第一动能;
28.反应率计算模块,用于通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;
29.能谱计算模块,用于按能量箱进行数据统计整理,计算得到核聚变产物的能谱。
30.优选地,还包括作图分析模块,按照得到核聚变产物的能谱分布对数据进行收集和统计,构建可视的直方图。
31.本发明的有益效果如下:
32.本发明通过对核聚变反应中反应物的动能进行计算,获得核聚变产物能谱在由质
量亏损决定的标准值处出现的展宽。通过引入狄拉克函数δ
12
,来避免同种粒子发生核聚变反应时的重复计数。通过采用基于r矩阵(r-matrix)的参数化反应截面数据,用于中性束注入条件下核聚变产物能谱的数值计算。
33.本发明所述计算方法能够准确计算出各类核聚变研究装置在不同运行模式下所获得的核聚变产物能谱,并通过建立中性束注入装置的简化计算模型,得到中性束注入情况下核聚变产物能谱的基本特征。
附图说明
34.图1是本发明所述的中性束注入条件下核聚变产物能谱的计算方法的步骤流程图。
35.图2是等离子体截面中心位置按照麦克斯韦-玻尔兹曼分布随机生成氘粒子的速度分布示意图。
36.图3是等离子体截面中心位置按照麦克斯韦-玻尔兹曼分布随机生成氚粒子的速度分布示意图。
37.图4等离子体截面中心处能量为40kev的中性束粒子注入模式下的中子归一化能谱。
38.图5等离子体截面中心处能量为60kev的中性束粒子注入模式下的中子归一化能谱
39.图6等离子体截面中心处能量为120kev的中性束粒子注入模式下的中子归一化能谱。
40.图7是中性束注入条件下核聚变产物能谱的计算系统示意图。
具体实施方式
41.下面结合附图和具体实施方式对本发明做详细描述。
42.实施例1
43.如图1所示,一种中性束注入条件下核聚变产物能谱的计算方法,所述的方法包括步骤如下:
44.根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;
45.根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;
46.两种反应物核聚变后得到核聚变产物,假设核聚变产物的运动方向沿任意给定的方向,根据得到的第一相对运动动能计算核聚变产物的第一动能;
47.通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;
48.按能量箱进行数据统计整理,计算得到核聚变产物的能谱。
49.在一个具体的实施例中,本发明考虑双体反应理论,双体反应理论适用于第一反应物和第二反应物经过反应生成第一核聚变产物和第二核聚变产物的反应类型。
50.所述的第一反应物、第二反应物在中性束注入条件下发生核聚变,所述的中性束
注入条件认为中性束内部的中性束粒子与中性束粒子之间不发生核聚变反应,束粒子只有背景反应物之间发生核聚变反应。
51.本实施例中,所述背景反应物是托卡马克装置中的背景等离子体,背景反应物遵循麦克斯韦-玻尔兹曼分布,同时假设托卡马克装置中的等离子体是各向同性。
52.本实施例中,所述中性束粒子是单能粒子束、或给定速度分布的粒子束。
53.在一个具体的实施例中,双体反应理论适用于两种反应物生成另外两种核聚变产物的反应类型。对第一反应物和第一反应物的速度进行随机抽样,生成速度和使其满足设定的速度分布函数和同时每种反应物粒子对应的权重设定为ω1和ω2。
54.根据得到每种反应物对应的每个蒙卡样本的速度及反应物的质量计算两种反应物的相对速度、与两者的质心速度,如下:
55.首先给出两种反应物的相对速度和两种反应物的质心速度
[0056][0057][0058]
其中,m1、m2为两种反应物的质量,为两种反应物在实验室参考系下的速度。
[0059]
计算在质心系下,两种反应物组成系统的第一相对运动动能:
[0060][0061]
式中,μ为约化质量,定义为:
[0062][0063]
在一个具体的实施例中,计算核聚变产物的第一动能,具体如下:
[0064]
根据两种反应物和核聚变产物之间的质量亏损,计算得到反应后系统增加的的能量;
[0065]
q=(m1+m
2-m
3-m4)c2[0066]
其中,c代表光速,m3、m4为两种核聚变产物的质量,q为能量。
[0067]
根据动量守恒定律得到质心系下,计算两种核聚变产物组成系统的第二相对运动动能:
[0068][0069]
将第二相对运动动能根据两种核聚变产物的质量比分配,得到两种核聚变产物在质心系中的第一动能。
[0070]
本实施例以第一核聚变产物为例,得到第一核聚变产物在质心系中的动能:
[0071]
[0072]
其中,
[0073][0074][0075]
式中,为第一核聚变产物在质心系下的速度,为第一核聚变产物在实验室参考系下的速度,m3表示第一核聚变产物的质量,m4表示第二核聚变产物的质量、是质心在实验室标系下的速度。
[0076]
在一个具体的实施例中,根据反应物粒子的分布规律得到两种反应物在实验室参考系下的归一化速度分布函数;
[0077]
将粒子动能的概率密度函数和微观粒子运动速率的麦克斯韦-玻尔兹曼分布的标准表达式代入反应率表达式,在一定条件下进行代换计算得到最终的能谱计算式,通过对能谱计算式进行积分得到反应率计算式。
[0078]
本实施例中,在反应率计算式中通过引入狄拉克函数来避免同种粒子发生反应时的重复计数。
[0079]
本实施例中,基于r矩阵的参数化反应截面,具体通过引入伽莫夫常数bg和pad
é
多项式s(er)计算得到反应的微分截面,基于反应的微分截面,计算得到反应速率系数。
[0080]
两种反应物在实验室参考系下的归一化速度分布函数:
[0081][0082]
通过引入狄拉克函数δ
12
来避免同种粒子发生反应时的重复计数:
[0083][0084]
其中,i、j代表反应物的类型。
[0085]
则反应率r表示单位体积、单位时间内两种反应物发生反应的次数:
[0086][0087]
式中,n1、n2是两种反应物的数密度,《σv》为反应速率系数:
[0088][0089]
其中,σ(er)为核反应截面,其与原子核自身的性质以及两个原子核相对运动速度(vr)有关。
[0090]
所以反应率r可以表示为:
[0091][0092]
通过引入一种基于r矩阵(r-matrix)的参数化反应截面来进行后续计算:
[0093][0094]
式中,bg为伽莫夫常数,s(er)为pad
é
多项式:
[0095][0096]
式中,a1、a2、a3、a4、a5、b1、b2、b3、b4均为常数,其值如下表:
[0097][0098]
在麦克斯韦-玻尔兹曼分布下温度为t的微观粒子运动速率表达式为:
[0099][0100]
式中,m是粒子质量,v为微观粒子的运动速率。
[0101]
经代换推导得到能谱的计算式:
[0102]
[0103]
反应率r为:
[0104]
r=∫r(e3)de3[0105]
在本实施例中,双体系统中质心速度为实验室参考系下第一核聚变产物的速度为:
[0106][0107]
则第一核聚变产物在实验室参考系下的动能:
[0108][0109]
对于给定速度和的反应物粒子,单位时间、单位体积及单位立体角下的反应率为:
[0110][0111]
每个蒙特卡洛事件对反应率的影响为:
[0112][0113]
在n次迭代后能得到一系列的e
3,i
和ri(i=1,2,...,n)的值,第一核聚变产物的能谱由e
3,i
的值放入能量区间中得到。
[0114]
能量区间ej至e
j+1
中的离子数rj由反应率的平均值给出:
[0115][0116]
其中,w为每个反应物粒子的权重之和、为狄拉克函数,rj的方差由ri的方差决定:
[0117][0118]
其中:
[0119][0120]
每个能量区间的相对误差可以表示为
[0121][0122]
本实施例通过对核聚变反应中反应物的动能进行计算,获得核聚变产物能谱在由质量亏损决定的标准值处出现的展宽。通过引入狄拉克函数δ
12
,来避免同种粒子发生核聚变反应时的重复计数。通过采用基于r矩阵(r-matrix)的参数化反应截面数据,用于中性束注入条件下核聚变产物能谱的数值计算。
[0123]
本实施例所述计算方法能够准确计算出各类核聚变研究装置在不同运行模式下
所获得的核聚变产物能谱,并通过建立中性束注入装置的简化计算模型,得到中性束注入情况下核聚变产物能谱的基本特征。
[0124]
本实施例通过对适用于核聚变研究装置中核聚变双体反应理论的扩展应用,能够准确计算出两种核聚变产物的能谱和反应率数据。
[0125]
本实施例采用蒙特卡罗方法在归一化分布函数下进行随机抽样建立计算方法,能够满足各种速度分布下核聚变反应产物能谱计算的需求,具有更加广泛的适用性。
[0126]
本实施例通过对不同中性束注入装置建立简化模型,使得能够对中性束注入条件下的核聚变反应产物能谱进行有效计算,准确表现出不同中性束注入装置对于核聚变产物能谱影响的差异,进一步满足核聚变相关科学研究的需要。
[0127]
本实施例通过对反应物温度的计算,可以得到明确的产物能谱展宽,使得产物能谱的能量分布更加清晰明确。
[0128]
实施例2
[0129]
基于实施例1所示的中性束注入条件下核聚变产物能谱的计算方法,本实施例给出一个应用例子如下:
[0130]
对托卡马克核聚变研究装置在中性束注入条件下(beam-thermal mode,bt模式)和装置稳态运行的模式(thermal-thermal mode,tt模式)下的氘氚核聚变反应产物能谱进行计算,并将计算结果进行对比。对所用中性束注入装置建立简化模型,该中性束注入的粒子为氘粒子。
[0131]
本实施例考虑在所有位置处的热离子与束离子的反应。第一反应物代表热氚离子,遵循麦克斯韦-玻尔兹曼分布。第二反应物代表中性束氘粒子,中性束注入装置所注入的束粒子能量分布为:能量为120kev的束粒子占63%;能量为60kev的束粒子占21%;能量为40kev的束粒子占16%。
[0132]
本实施例假设两个装置的中性束都在发射方向上有一个θ
ext
=10m rad的展宽,同时中性束在发射方向上的展宽服从的分布定义为:
[0133][0134]
本实施例为研究不同能量的中性束粒子影响下的聚变产物能谱,使用两种反应物均为麦克斯韦-玻尔兹曼分布的情况,在两种模式下对等离子体截面中心位置的中子的归一化能谱进行计算。归一化半径等离子体截面中心位置处的温度,也就是td=t
t
=36.49kev。
[0135]
(1)对第一反应物和第二反应物的速度进行随机抽样,生成和使其满足对应的分布函数和同时每个反应物粒子对应的权重为ω1和ω2。
[0136]
1.1)根据麦克斯韦-玻尔兹曼分布随机生成两种反应物粒子,即106个氘粒子和氚粒子。
[0137]
1.2)两种反应物粒子的速度分布如图2、3所示。从图2、3中可以看出随机抽样生成的粒子分布与实际的麦克斯韦-玻尔兹曼分布十分接近,因此可以认为抽样过程是准确的。
[0138]
1.3)同时令氘粒子和氚粒子的权重分别为ω1和ω2。
[0139]
(2)计算聚变产物中中子的动能e3。
[0140]
2.1)首先根据两种反应物的质量和其在实验室参考系下的速度,得出在质心系下,两种反应物组成系统的第一相对运动动能。
[0141]
2.2)根据反应物和产物之间的质量亏损,计算得到反应后系统增加的第二动能。
[0142]
2.3)根据动量守恒定律,得到质心系下两种核聚变产物组成系统的第二相对运动动能。
[0143]
2.4)根据中子和α粒子两种聚变反应产物的质量比得到中子在质心系中的第一动能。
[0144]
(3)计算反应的微分截面,进一步计算出反应率。
[0145]
3.1)引入一种基于r矩阵(r-matrix)的参数化反应截面来进行后续计算。
[0146]
3.2)通过引入伽莫夫常数bg和pad
é
多项式s(er)计算得到反应的微分截面。
[0147]
3.3)基于反应的微分截面,进一步计算得到反应速率系数。
[0148]
3.4)将反应速率系数带入反应率计算式,得到最终的反应率结果。
[0149]
(4)聚变产物能谱的计算。
[0150]
按能量箱进行数据的统计整理,计算得到聚变产物的能谱。
[0151]
(5)收集数据进行统计,做出直方图。
[0152]
在等离子体截面中心位置处核聚变得到的产物(中子)的归一化能谱如图4、5、6所示。
[0153]
实施例3
[0154]
如图4所示,一种中性束注入条件下核聚变产物能谱的计算系统,包括:
[0155]
随机抽样模块,根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;
[0156]
第一相对运动动能计算模块,用于根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;
[0157]
第一动能计算模块,用于根据得到的第一相对运动动能计算核聚变产物的第一动能;
[0158]
反应率计算模块,用于通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;
[0159]
能谱计算模块,用于按能量箱进行数据统计整理,计算得到核聚变产物的能谱。
[0160]
本实施例,还包括作图分析模块,按照得到核聚变产物的能谱分布对数据进行收集和统计,构建可视的直方图。
[0161]
显然,本发明的上述实施例仅仅是为清楚地说明本发明所作的举例,而并非是对本发明的实施方式的限定。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明权利要求的保护范围之内。

技术特征:
1.一种中性束注入条件下核聚变产物能谱的计算方法,其特征在于:所述的方法包括步骤如下:根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;两种反应物核聚变后得到核聚变产物,假设核聚变产物的运动方向沿任意给定的方向,根据得到的第一相对运动动能计算核聚变产物的第一动能;通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;按能量箱进行数据统计整理,计算得到核聚变产物的能谱。2.根据权利要求1所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:所述的第一反应物、第二反应物在中性束注入条件下发生核聚变,所述的中性束注入条件认为中性束内部的中性束粒子与中性束粒子之间不发生核聚变反应,束粒子只有背景反应物之间发生核聚变反应。3.根据权利要求2所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:所述背景反应物是托卡马克装置中的背景等离子体,背景反应物遵循麦克斯韦-玻尔兹曼分布,同时假设托卡马克装置中的等离子体是各向同性。4.根据权利要求2所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:所述中性束粒子是单能粒子束、或给定速度分布的粒子束。5.根据权利要求1所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:计算核聚变产物的第一动能,具体如下:根据两种反应物和核聚变产物之间的质量亏损,计算得到反应后系统增加的第二动能;根据动量守恒定律得到质心系下,计算两种核聚变产物组成系统的第二相对运动动能;将第二相对运动动能根据两种核聚变产物的质量比分配,得到两种核聚变产物在质心系中的第一动能。6.根据权利要求1所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:根据反应物粒子的分布规律得到两种反应物在实验室参考系下的归一化速度分布函数;将粒子动能的概率密度函数和微观粒子运动速率的麦克斯韦-玻尔兹曼分布的标准表达式代入反应率表达式,在一定条件下进行代换计算得到最终的能谱计算式,通过对能谱计算式进行积分得到反应率计算式。7.根据权利要求6所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:在反应率计算式中通过引入狄拉克函数来避免同种粒子发生反应时的重复计数。8.根据权利要求6所述的中性束注入条件下核聚变产物能谱的计算方法,其特征在于:基于r矩阵的参数化反应截面,具体通过引入伽莫夫常数b
g
和pad
é
多项式s(e
r
)计算得到反应的微分截面,基于反应的微分截面,计算得到反应速率系数。9.一种中性束注入条件下核聚变产物能谱的计算系统,其特征在于:包括:
随机抽样模块,根据第一反应物、第二反应物的速度分布函数,均随机生成n个蒙卡样本,得到每种反应物对应的每个蒙卡样本的速度、权重;第一相对运动动能计算模块,用于根据两种反应物的相对速度、与两者的质心速度,计算在质心系下,两种反应物组成系统的第一相对运动动能;第一动能计算模块,用于根据得到的第一相对运动动能计算核聚变产物的第一动能;反应率计算模块,用于通过引入狄拉克函数和基于r矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;能谱计算模块,用于按能量箱进行数据统计整理,计算得到核聚变产物的能谱。10.根据权利要求9所述的中性束注入条件下核聚变产物能谱的计算系统,其特征在于:还包括作图分析模块,按照得到核聚变产物的能谱分布对数据进行收集和统计,构建可视的直方图。

技术总结
本发明公开了一种中性束注入条件下核聚变产物能谱的计算方法和系统,如下:根据第一反应物、第二反应物的速度分布函数,均随机生成N个蒙卡样本,得到每个蒙卡样本的速度、权重;根据两种反应物的相对速度、与两者的质心速度,计算两种反应物组成系统的第一相对运动动能;根据得到的第一相对运动动能计算核聚变产物的第一动能;通过引入狄拉克函数和基于R矩阵的参数化反应截面,计算核聚变产物在第一动能下的反应速率系数,并带入反应率计算式,得到最终的反应率结果;按能量箱进行数据统计整理,计算得到核聚变产物的能谱。本发明能够准确计算出各类核聚变研究装置在不同运行模式下所获得的核聚变产物能谱。式下所获得的核聚变产物能谱。式下所获得的核聚变产物能谱。


技术研发人员:王志斌 张镇乾 陈坚 黄汉宇
受保护的技术使用者:中山大学
技术研发日:2023.04.13
技术公布日:2023/7/25
版权声明

本文仅代表作者观点,不代表航空之家立场。
本文系作者授权航家号发表,未经原创作者书面授权,任何单位或个人不得引用、复制、转载、摘编、链接或以其他任何方式复制发表。任何单位或个人在获得书面授权使用航空之家内容时,须注明作者及来源 “航空之家”。如非法使用航空之家的部分或全部内容的,航空之家将依法追究其法律责任。(航空之家官方QQ:2926969996)

飞行汽车 https://www.autovtol.com/

分享:

扫一扫在手机阅读、分享本文

相关推荐