基于自适应多目标粒子群的含油污水处理过程优化控制方法与流程

文档序号:17497524发布日期:2019-04-23 21:44阅读:139来源:国知局
基于自适应多目标粒子群的含油污水处理过程优化控制方法与流程

本发明涉及一种基于自适应多目标粒子群的含油污水处理过程优化控制方法。本发明利用基于自适应多目标粒子群的优化控制方法实现含油污水处理过程中溶解氧do和硝态氮sno浓度的控制,溶解氧do和硝态氮sno的浓度是含油污水处理过程中的关键控制参数,对含油污水处理的效果、出水水质,含油污水处理过程的能耗都有着重要影响。将基于自适应多目标粒子群的优化控制方法应用在含油污水处理过程中,实现溶解氧do和硝态氮sno浓度的优化控制,在降低含油污水处理过程能耗的同时节约投资和运行成本,保证污水处理厂的平稳高效运行,既属于水处理领域,又属于智能控制领域。



背景技术:

随着我国经济社会的高速发展,石油工业进程的加快,人类活动对水环境造成的污染不断加剧,含油污水对人类、动物和植物乃至整个生态系统都产生了极其严重的影响。与此同时,国民经济的增长和公众环保意识的增强,使含油污水处理自动化技术迎来了前所未有的发展机遇;怎样防治水体污染,如何及时有效的处理含油污水并再次利用,成为我国迫在眉睫的问题;然而含油污水处理过程电能消耗大、运行成本高,研究污水处理过程优化控制实现节能降耗的意义重大,是未来污水处理行业必然的发展趋势。

污水生化处理过程的实质是利用污泥中的微生物的生命活动分解污水中的有机污染物,使污水得到净化。含油污水处理过程中,主要控制变量是溶解氧do和硝态氮sno浓度。溶解氧do和硝态氮sno浓度的高低变化会直接影响硝化过程和反硝化过程的进行,进而对含油污水处理过程中的能耗产生影响。硝化反应主要是在有氧条件下进行的,当溶解氧do浓度变大时,出水氨氮和总氮的浓度呈现下降趋势,当溶解氧do浓度增长到一定范围时,出水氨氮的变化幅度开始减弱,总氮同时也受到硝态氮的影响,硝态氮浓度增加的同时总氮浓度也会升高。然而,含油污水处理过程中反硝化反应主要是在缺氧环境下进行的,缺氧区的硝态氮sno浓度是衡量脱氮效果的重要指标,它反映了反硝化反应过程的进程,将硝态氮sno浓度控制在一个合适的范围内,能够提高反硝化反应的潜力。因此,实时动态优化控制溶解氧do和硝态氮sno浓度,对于保证出水水质,节能降耗是非常有必要的。由于含油污水的悬浮物,色度含量高,有机物的成分多样化,含油污水处理过程的非线性,时变性,动态不确定等特点,增加了溶解氧do和硝态氮sno的控制难度;近年来部分学者又以bsm1为基础,采用多目标遗传算法(moga)对含油污水处理过程进行优化,使曝气能耗和泵送能耗达到最小。但这些优化方案大多属于静态优化、稳态优化,难以根据进水水质水量的变化进行动态实时调整,而且大多采用遗传算法实现。与遗传算法相比,粒子群算法收敛速度快,不需要复杂的交叉变异操作,算法简单、参数少、易于实现,但是标准的粒子群算法容易出现早熟收敛、搜索精度较低的问题,因此设计一种自适应多目标粒子群算法,提高算法的收敛精度,从而能够很好地实现溶解氧do和硝态氮sno浓度的优化控制,降低运行成本,具有很好的实际应用价值。



技术实现要素:

本发明所要解决的技术问题是现有技术中标准的粒子群算法容易出现早熟收敛、搜索精度较低的问题,提供一种新的基于自适应多目标粒子群的含油污水处理过程优化控制方法。该方法具有不会出现早熟收敛、搜索精度较高的优点。

为解决上述问题,本发明采用的技术方案如下:一种基于自适应多目标粒子群的含油污水处理过程优化控制方法,首先通过模糊神经网络建立含油污水处理过程的溶解氧和硝态氮的优化设定值与曝气能耗和泵送能耗的目标函数;其次,针对标准粒子群算法容易早熟收敛和收敛精度低的缺点,采用自适应多目标粒子群优化方法,实现对含油污水处理目标函数的优化,同时获得溶解氧和硝态氮的优化设定值;最后,利用控制器对溶解氧和硝态氮的优化设定值进行跟踪控制,完成含油污水处理过程的多目标优化控制;

具体包括以下步骤:

1)设计用于含油污水处理过程优化控制的目标函数;

(2)利用模糊神经网络建立溶解氧和硝态氮的优化设定值与曝气能耗和泵送能耗间的关系表达式;

(3)出水水质的约束处理;

(4)利用自适应多目标粒子群优化目标函数获得pareto最优解,具体为:

①初始化粒子群的速度vi(0)、位置ai(0)、惯性权重ωi(0)、学习因子c1i(0)和c2i(0),并将各粒子的初始位置设为当前历史最优位置pi(0),同时设定种群规模s,最大进化代数m,变量维数d;

②根据目标函数计算每个粒子的适应度值,确定第t次迭代的个体最优解pi(t),其定义式为:

非支配解集a(t)通过a(t-1)更新,其公式为:

其中,a(t)=[a1(t),a2(t),…,aq(t)],q是知识库a(t)的最大容量,k是知识库中包含非支配解的数目,表示ai(t-1)和pi(t-1)互不支配;

③确定第t+1次迭代的全局最优解gbest(t+1)

其中,gbest(t+1)是第t+1次迭代的全局最优解,dgbest(t+1)是第t+1次迭代的多样性较好的全局最优解,其中,dgbest(t+1)的定义式为:

dgbest(t+1)=a(t),a(t)∈μtbest.(4)

cgbest(t+1)是第t+1次迭代的收敛性较好的全局最优解,其定义式为:

cgbest(t+1)=argmaxcdt(ai(t)),(5)

cdt(ai(t))是非支配解ai(t)的收敛度,i=1,2,3…k,e(t)是第t次迭代非支配解集的分布熵,其定义式为:

其中,un是知识库中非支配解的单元格,n=1,2,3…,n,n是单元格的总数目,pt(un)是第t次迭代单元格un的概率分布函数,其表达式为:

mt(un)是第t次迭代单元格un中非支配解的数目,mt(un)是第t次迭代单元格un中包含的所有非支配解的解集,其中,将最小的mt(un)记为mtbest,相应地,mtbest是第t次迭代最优的非支配解集,此外,cgbest(t+1)的选择通过能够反映支配关系的收敛度判断,收敛度cdt(ai(t))的定义为:

其中,是第j个能被非支配解ai(t)所支配的解,dst(ai(t))是第t次迭代非支配解ai(t)的支配强度,其定义式为:

dst(ai(t))是第t次迭代的非支配解ai(t)能够支配解的总数目,a(t)是第t次迭代的非支配解集,a(t-1)是第t-1次迭代的非支配解集;

④更新每个粒子的速度和位置;

⑤判断是否达到最大进化次数m,如达到,则结束,否则返回到②;

(5)从自适应多目标粒子群算法获得的一组pareto最优解中,找到当前状态下的一组满意优化解作为底层控制器的优化设定值;

(6)执行底层控制策略,溶解氧和硝态氮浓度分别通过曝气池氧气转换系数kla5和内回流量qa进行调节。

上述技术方案中,优选地,控制器为pid。

上述技术方案中,优选地,用于含油污水处理过程优化控制的目标函数:

其中,x(k)=[x1(k),x2(k)]t为k时刻的优化变量,x1(k)为k时刻的溶解氧do浓度设定值,x2(k)为k时刻的硝态氮浓度设定值,fae(x)和fpe(x)分别为曝气能耗、泵送能耗与优化变量的关系表达式,g1(x)和g2(x)分别为出水氨氮、总氮浓度与优化变量的关系表达式,c1和c2分别为出水氨氮和总氮的约束值,c1∈[0,4],c2∈[0,18];分别表示溶解氧优化设定值的下限和上限值,x1∈[0.4,3],分别表示硝态氮浓度优化设定值的下限和上限值,x2∈[0.5,2],优化周期为2h。

上述技术方案中,优选地,模糊神经网络的计算方式如下:

其中,x(k)=[x1(k),x2(k)]t为第k时刻模糊神经网络输入,cj=[c1j,c2j],σj=[σ1j,σ2j]分别是rbf层的第j个神经元的中心向量和宽度向量,为第k时刻第j个神经元的rbf层输出,p是rbf层和规则层的神经元数目;vl(k)为第k时刻对应的第l个规则层输出,v(k)=[v1(k),v2(k)…vp(k)]t为第k时刻规则层输出向量;w1=[w11,w21…wp1]和w2=[w12,w22…wp2]是输出神经元和规则层之间的权值向量,y(k)为该神经网络的输出,为污水处理系统实际物理量输出,基于bsm1模型数据获得;

设k时刻网络调整的目标函数为:

采用梯度下降算法,权重更新公式为:

式中,αq(k)=[θq(k)tcq(k)tσq(k)t]t为网络的学习参数向量,网络学习率η=0.01;

上述技术方案中,优选地,采用惩罚函数法对优化模型中的约束进行处理,其中,约束惩罚项定义为:

fpenalty(x)=max{g1(x)-4,0}+max{g2(x)-18,0},(16)

加入惩罚项的曝气能耗和泵送能耗目标函数为:

将建立的含油污水处理过程约束优化问题转化为无约束多目标优化问题,其中,ε为惩罚因子,设定为较大的正实数。

上述技术方案中,优选地,更新每个粒子的速度和位置:

vi(t+1)=ωi(t)vi(t)+c1r1(pi(t)-ai(t))+c2r2(gbestd(t)-ai(t));(18)

ai(t+1)=ai(t-1)+vi(t+1);

其中,r1和r2分别表示最好先前位置系数和全局最优位置系数,r1和r2取[0,1]的任意数。

上述技术方案中,优选地,步骤(6)中,执行底层控制策略,溶解氧和硝态氮浓度分别通过曝气池第5分区氧气转换系数kla5和内回流量qa进行调节。

上述技术方案中,优选地,步骤①中,初始化粒子群的速度vi(0)、位置ai(0)、惯性权重ωi(0)、学习因子c1i(0)和c2i(0),并将各粒子的初始位置设为当前历史最优位置pi(0),同时设定种群规模s=40,最大进化代数m=30,变量维数d=2。

本发明提供的一种基于自适应多目标粒子群的含油污水处理过程优化控制方法,包括以下步骤:

(1)设计用于含油污水处理过程优化控制的目标函数:

其中,x(k)=[x1(k),x2(k)]t为k时刻的优化变量,x1(k)为k时刻的溶解氧do浓度设定值,x2(k)为k时刻的硝态氮浓度设定值,fae(x)和fpe(x)分别为曝气能耗、泵送能耗与优化变量的关系表达式,g1(x)和g2(x)分别为出水氨氮、总氮浓度与优化变量的关系表达式,c1和c2分别为出水氨氮和总氮的约束值,c1∈[0,4],c2∈[0,18];分别表示溶解氧优化设定值的下限和上限值,x1∈[0.4,3],分别表示硝态氮浓度优化设定值的下限和上限值,x2∈[0.5,2],优化周期为2h。

(2)利用模糊神经网络建立溶解氧和硝态氮的优化设定值与曝气能耗和泵送能耗间的关系表达式,其中,模糊神经网络的计算方式如下:

其中,x(k)=[x1(k),x2(k)]t为第k时刻模糊神经网络输入,cj=[c1j,c2j],σj=[σ1j,σ2j]分别是rbf层的第j个神经元的中心向量和宽度向量,为第k时刻第j个神经元的rbf层输出,p是rbf层和规则层的神经元数目;vl(k)为第k时刻对应的第l个规则层输出,v(k)=[v1(k),v2(k)…vp(k)]t为第k时刻规则层输出向量;w1=[w11,w21…wp1]和w2=[w12,w22…wp2]是输出神经元和规则层之间的权值向量,y(k)为该神经网络的输出,为污水处理系统实际物理量输出,基于bsm1模型数据获得。

设k时刻网络调整的目标函数为:

采用梯度下降算法,权重更新公式为:

式中,αq(k)=[θq(k)tcq(k)tσq(k)t]t为网络的学习参数向量,网络学习率η=0.01。

(3)出水水质的约束处理

采用惩罚函数法对优化模型中的约束进行处理,其中,约束惩罚项定义为:

fpenalty(x)=max{g1(x)-4,0}+max{g2(x)-18,0},(7)

加入惩罚项的曝气能耗和泵送能耗目标函数为:

即将建立的含油污水处理过程约束优化问题转化为无约束多目标优化问题。其中,ε为惩罚因子,设定为较大的正实数。

(4)利用自适应多目标粒子群优化目标函数获得pareto最优解,具体为:

①初始化粒子群的速度vi(0)、位置ai(0)、惯性权重ωi(0)、学习因子c1i(0)和c2i(0),并将各粒子的初始位置设为当前历史最优位置pi(0)。同时,设定种群规模s=40,最大进化代数m=30,变量维数d=2;

②根据目标函数计算每个粒子的适应度值。确定第t次迭代的个体最优解pi(t),其定义式为:

非支配解集a(t)通过a(t-1)更新,其公式为:

其中,a(t)=[a1(t),a2(t),…,aq(t)],q是知识库a(t)的最大容量,k是知识库中包含非支配解的数目,表示ai(t-1)和pi(t-1)互不支配。

③确定第t+1次迭代的全局最优解gbest(t+1)。

其中,gbest(t+1)是第t+1次迭代的全局最优解,dgbest(t+1)是第t+1次迭代的多样性较好的全局最优解,其中,dgbest(t+1)的定义式为:

dgbest(t+1)=a(t),a(t)∈μtbest.(12)

cgbest(t+1)是第t+1次迭代的收敛性较好的全局最优解,其定义式为:

cgbest(t+1)=argmaxcdt(ai(t)),(13)

cdt(ai(t))是非支配解ai(t)的收敛度,i=1,2,3…k,e(t)是第t次迭代非支配解集的分布熵,其定义式为:

其中,un是知识库中非支配解的单元格,n=1,2,3…,n,n是单元格的总数目,pt(un)是第t次迭代单元格un的概率分布函数,其表达式为:

mt(un)是第t次迭代单元格un中非支配解的数目,mt(un)是第t次迭代单元格un中包含的所有非支配解的解集。其中,将最小的mt(un)记为mtbest,相应地,mtbest是第t次迭代最优的非支配解集。此外,cgbest(t+1)的选择通过能够反映支配关系的收敛度判断,收敛度cdt(ai(t))的定义为:

其中,是第j个能被非支配解ai(t)所支配的解,dst(ai(t))是第t次迭代非支配解ai(t)的支配强度,其定义式为:

dst(ai(t))是第t次迭代的非支配解ai(t)能够支配解的总数目,a(t)是第t次迭代的非支配解集,a(t-1)是第t-1次迭代的非支配解集;

④更新每个粒子的速度和位置:

vi(t+1)=ωi(t)vi(t)+c1r1(pi(t)-ai(t))+c2r2(gbestd(t)-ai(t));(18)

ai(t+1)=ai(t-1)+vi(t+1);(19)

其中,r1和r2分别表示最好先前位置系数和全局最优位置系数,r1和r2取[0,1]的任意数;

⑤判断是否达到最大进化次数m,如达到,则结束,否则返回到②。

(5)从自适应多目标粒子群算法获得的一组pareto最优解中,找到当前状态下的一组满意优化解作为底层pid控制器的优化设定值。

(6)执行底层pid控制策略,溶解氧和硝态氮浓度分别通过曝气池第5分区氧气转换系数kla5和内回流量qa进行调节。

本发明针对当前活性污泥法含油污水处理过程的复杂的、动态的、不稳定的生化反应过程,以及非线性、时变性、滞后性的特点;同时溶解氧do和硝态氮sno浓度之间存在强耦合关系,为了满足在出水水质达标的同时,降低运行能耗的需求,实现溶解氧do和硝态氮sno浓度的多目标控制,采用了基于多目标粒子群的含油污水处理模型预测控制方法实现溶解氧do和硝态氮sno浓度的控制,具有控制精度高、稳定性好等特点;本发明采用基于自适应多目标粒子群的含油污水处理模型对污水处理过程溶解氧do和硝态氮sno浓度进行优化控制,该优化控制方法解决了多个目标函数的优化求解问题,使控制器更好的满足当前环境的变化,实现了溶解氧do和硝态氮sno浓度实时闭环精确控制,避免了当前污水处理厂需要设计多个控制器进行控制的复杂过程,具有实时性强,结构简单等特点,取得了较好的技术效果。

附图说明

图1是本发明的自适应多目标粒子群优化控制系统整体结构图;

图2是本发明控制系统溶解氧do浓度结果图;

图3是本发明控制系统溶解氧do浓度误差图;

图4是本发明控制系统硝态氮sno浓度结果图;

图5是本发明控制系统硝态氮sno浓度结果误差图。

下面通过实施例对本发明作进一步的阐述,但不仅限于本实施例。

具体实施方式

【实施例1】

一种基于自适应多目标粒子群的含油污水处理过程优化控制方法,如图1所示,实现了含油污水处理过程中溶解氧do和硝态氮sno浓度的多目标优化控制;该控制方法通过模糊神经网络在线建模获得溶解氧do和硝态氮sno浓度的优化设定值与曝气能耗、泵送能耗、出水水质间的函数关系,将基于自适应多目标粒子群的优化控制方法应用在含油污水处理过程中,实现溶解氧do和硝态氮sno浓度的优化控制,在降低含油污水处理过程能耗的同时节约投资和运行成本,保证污水处理厂的平稳高效运行,既属于水处理领域,又属于智能控制领域。

本发明采用了如下的技术方案及实现步骤:

1.一种基于自适应多目标粒子群的含油污水处理过程优化控制方法设计包括以下步骤:针对序批式间歇活性污泥系统中溶解氧do浓度和硝态氮sno进行控制,以曝气能耗和泵送能耗为控制量,溶解氧do和硝态氮sno浓度为被控量,自适应多目标粒子群优化控制系统整体框架如图1;

(1)设计用于含油污水处理过程优化控制的目标函数:

其中,x(k)=[x1(k),x2(k)]t为k时刻的优化变量,x1(k)为k时刻的溶解氧浓度设定值,x2(k)为k时刻的硝态氮浓度设定值,fae(x)和fpe(x)分别为曝气能耗、泵送能耗与优化变量的关系表达式,g1(x)和g2(x)分别为出水氨氮、总氮浓度与优化变量的关系表达式,c1和c2分别为出水氨氮和总氮的约束值,c1∈[0,4],c2∈[0,18];分别表示溶解氧优化设定值的下限和上限值,x1∈[0.4,3],分别表示硝态氮浓度优化设定值的下限和上限值,x2∈[0.5,2],优化周期为2h。

(2)利用模糊神经网络建立溶解氧和硝态氮的优化设定值与曝气能耗和泵送能耗间的关系表达式,其中,模糊神经网络的计算方式如下:

其中,x(k)=[x1(k),x2(k)]t为第k时刻模糊神经网络输入,cj=[c1j,c2j],σj=[σ1j,σ2j]分别是rbf层的第j个神经元的中心向量和宽度向量,为第k时刻第j个神经元的rbf层输出,p是rbf层和规则层的神经元数目;vl(k)为第k时刻对应的第l个规则层输出,v(k)=[v1(k),v2(k)…vp(k)]t为第k时刻规则层输出向量;w1=[w11,w21…wp1]和w2=[w12,w22…wp2]是输出神经元和规则层之间的权值向量,y(k)为该神经网络的输出,为污水处理系统实际物理量输出,基于bsm1模型数据获得。

设k时刻网络调整的目标函数为:

采用梯度下降算法,权重更新公式为:

式中,αq(k)=[θq(k)tcq(k)tσq(k)t]t为网络的学习参数向量,网络学习率η=0.01。

(3)出水水质的约束处理

采用惩罚函数法对优化模型中的约束进行处理,其中,约束惩罚项定义为:

fpenalty(x)=max{g1(x)-4,0}+max{g2(x)-18,0},(7)

加入惩罚项的曝气能耗和泵送能耗目标函数为:

即将建立的含油污水处理过程约束优化问题转化为无约束多目标优化问题。其中,ε为惩罚因子,设定为较大的正实数。

(4)利用自适应多目标粒子群优化目标函数获得pareto最优解,具体为:

①初始化粒子群的速度vi(0)、位置ai(0)、惯性权重ωi(0)、学习因子c1i(0)和c2i(0),并将各粒子的初始位置设为当前历史最优位置pi(0)。同时,设定种群规模s=40,最大进化代数m=30,变量维数d=2;

②根据目标函数计算每个粒子的适应度值。确定第t次迭代的个体最优解pi(t),其定义式为:

非支配解集a(t)通过a(t-1)更新,其公式为:

其中,a(t)=[a1(t),a2(t),…,aq(t)],q是知识库a(t)的最大容量,k是知识库中包含非支配解的数目,表示ai(t-1)和pi(t-1)互不支配。

③确定第t+1次迭代的全局最优解gbest(t+1)。

其中,gbest(t+1)是第t+1次迭代的全局最优解,dgbest(t+1)是第t+1次迭代的多样性较好的全局最优解,其中,dgbest(t+1)的定义式为:

dgbest(t+1)=a(t),a(t)∈μtbest.(12)

cgbest(t+1)是第t+1次迭代的收敛性较好的全局最优解,其定义式为:

cgbest(t+1)=argmaxcdt(ai(t)),(13)

cdt(ai(t))是非支配解ai(t)的收敛度,i=1,2,3…k,e(t)是第t次迭代非支配解集的分布熵,其定义式为:

其中,un是知识库中非支配解的单元格,n=1,2,3…,n,n是单元格的总数目,pt(un)是第t次迭代单元格un的概率分布函数,其表达式为:

mt(un)是第t次迭代单元格un中非支配解的数目,mt(un)是第t次迭代单元格un中包含的所有非支配解的解集。其中,将最小的mt(un)记为mtbest,相应地,mtbest是第t次迭代最优的非支配解集。此外,cgbest(t+1)的选择通过能够反映支配关系的收敛度判断,收敛度cdt(ai(t))的定义为:

其中,是第j个能被非支配解ai(t)所支配的解,dst(ai(t))是第t次迭代非支配解ai(t)的支配强度,其定义式为:

dst(ai(t))是第t次迭代的非支配解ai(t)能够支配解的总数目,a(t)是第t次迭代的非支配解集,a(t-1)是第t-1次迭代的非支配解集;

④更新每个粒子的速度和位置:

vi(t+1)=ωi(t)vi(t)+c1r1(pi(t)-ai(t))+c2r2(gbestd(t)-ai(t));(18)

ai(t+1)=ai(t-1)+vi(t+1);(19)

其中,r1和r2分别表示最好先前位置系数和全局最优位置系数,r1和r2取[0,1]的任意数;⑤判断是否达到最大进化次数m,如达到,则结束,否则返回到②。

(5)从自适应多目标粒子群算法获得的一组pareto最优解中,找到当前状态下的一组满意优化解作为底层pid控制器的优化设定值。

(6)执行底层pid控制策略,pid控制器的参数设置为:kp,1=200,ki,1=15,kd,1=2和kp,2=20000,ki,2=5000,kd,2=400。溶解氧和硝态氮浓度分别通过曝气池第5分区氧气转换系数kla5和内回流量qa进行调节。

(7)利用pid控制器对溶解氧do和硝态氮sno浓度优化设定值进行追踪控制,整个控制系统的输出为溶解氧do和硝态氮sno浓度的值优化设定值和追踪控制值;图2显示系统的溶解氧do浓度优化设定值和追踪控制值,x轴:时间,单位是天数,y轴:溶解氧do的优化设定值和追踪控制值,单位是毫克/升,实线为溶解氧do浓度优化设定值,虚线是实际溶解氧do追踪控制值;优化设定溶解氧do浓度值与实际溶解氧do追踪控制浓度值的误差如图3,x轴:时间,单位是天数,y轴:溶解氧do浓度误差值,单位是毫克/升;图4显示系统的硝态氮sno优化设定值和追踪控制值,x轴:时间,单位是天数,y轴:硝态氮sno的优化设定值和追踪控制值,单位是毫克/升,实线为硝态氮sno浓度优化设定值,虚线是硝态氮sno追踪控制浓度值;优化设定硝态氮sno浓度值与实际硝态氮sno追踪控制浓度值的误差如图5,x轴:时间,单位是天数,y轴:硝态氮sno浓度误差值,单位是毫克/升,结果证明该方法的有效性。

特别要注意:本发明只是为了描述方便,采用的是对溶解氧do和硝态氮sno浓度的控制,同样该发明也可适用污水处理过程氨氮的控制等,只要采用了本发明的原理进行控制都应该属于本发明的范围。

当前第1页1 2 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1