首页> 中国专利> 散堆填料中液泛气速的建模和预报方法

散堆填料中液泛气速的建模和预报方法

摘要

用于散堆填料塔中液泛气速的非线性建模和预报方法,首先采用高斯混合模型方法对散堆填料塔的液泛数据进行有效的聚类;根据聚类后得到的子集分别建立相应的支持向量机子模型;再利用加权最小二乘方法,对这些液泛气速的子模型进行融合,最终,得到液泛气速的预报结果。使用该非线性集成模型,可以为各种类型的散堆填料塔中的液泛气速提供比较准确的模型,用于液泛气速的预报,有效提高相应塔设备的操作效率和安全性能。

著录项

  • 公开/公告号CN103279635A

    专利类型发明专利

  • 公开/公告日2013-09-04

    原文格式PDF

  • 申请/专利权人 浙江工业大学;

    申请/专利号CN201310113831.5

  • 发明设计人 刘毅;高增梁;李常跞;

    申请日2013-04-02

  • 分类号G06F19/00;

  • 代理机构杭州天正专利事务所有限公司;

  • 代理人王兵

  • 地址 310014 浙江省杭州市下城区潮王路18号

  • 入库时间 2024-02-19 20:08:03

法律信息

  • 法律状态公告日

    法律状态信息

    法律状态

  • 2016-05-18

    授权

    授权

  • 2013-10-09

    实质审查的生效 IPC(主分类):G06F19/00 申请日:20130402

    实质审查的生效

  • 2013-09-04

    公开

    公开

说明书

技术领域

本发明涉及填料塔设计阶段重要参数建模和预报方法的技术领域,特别涉及一种通用的适合散堆(乱堆)填料塔中液泛气速的非线性建模和预报方法。 

背景技术

塔设备常用于工业生产中的精馏、吸收、萃取以及冷却过程中的传热与传质。填料塔具有生产能力大、分离效率高、压力降小、操作弹性大和持液量小等优点,在精细化工、石油化工、医药、环保等行业都得到广泛应用。尤其在石油化工中,与精馏和吸收塔相关的能耗占了很大的比例。因此,填料塔的操作越接近其最大的有效生产能力,装置的生产效率越高,相应的能耗越少。 

塔内液泛的发生是现有工业填料塔在操作时主要存在的问题,它破坏了塔的正常操作,也是限制其运行效率不高、能耗大的重要原因。塔中发生液泛时的速度称为液泛气速,也称泛点气速。一般来说,液泛点是填料塔的操作极限。填料塔只有在泛点气速下才能稳定操作,但如果操作气速太低又会造成气液的分布不均匀及投资浪费。同时,液泛气速也是计算在给定液相负荷下塔的最大承载能力的必要参数,是填料塔设计的重要依据。因此,确定各种类型填料的液泛气速也变得尤为重要,准确的预报液泛气速对填料塔的设计与操作都具有重要的意义。 

迄今为止,利用各种机理和操作条件建立了不同形式的液泛气速预报模型,在工业应用上主要依靠传统的经验公式(或相应的关联图)。但是这些经验模型都需特定的填料常数,其通用性较差。另外,随着各种新型填料的出现,许多条件与经验公式并不相符,使模型的预报效果降低。因此,由于各种限制条件,至今还没有较准确且通用的模型适合某类工业填料塔的应用。在实际生产操作中,考虑到传统液泛预报模型的准确性低和适用范围等因素,大多数塔的操作指标都远低于最大有效能力。为防止液泛的发生,实际操作速度一般在液泛气速50%以下,即实际工作点远离最佳工作点。这种保守操作导致了塔生产能力和分离效率不高,经济效益较差,能耗显著增加。因此,这些传统模型在工业应用中都还有较大的局限,不能满足工程应用以及激烈的市场竞争和需求,有必要建立一种准确度较高,通用性较强的液泛气速模型以适应工程应用。 

随着过程数据能够及时获得,各种数据驱动的建模方法得到广泛研究和应用,但用于液泛气速预报的却很少。经文献检索发现,一种反向传播神经网络(BP-NN)曾用于液泛气速的建模和预报,虽然BP-NN方法的非线性建模能力较强,但存在训练费时、拓扑结构难确定和推广能力差等问题,且当样本较少时,易出现过拟合问题。 

最近,最小二乘支持向量机(LSSVM)在不少化工过程建模中获得了较好的效果。LSSVM采用结构风险最小化原则和核技巧,较NN等方法能提高模型在限样本情况下的建模能力和泛化能力。同时,LSSVM采用等式约束,避免了求解耗时的二次规划问题,训练速度 很快。这为液泛重要参数,如液泛气速的建模和预报提供了新方法。 

同时,考虑到实际工业生产中的填料种类繁多,加上相关的液泛数据并不好获得,使得采集的液泛数据出现多样性和不平衡性的特点。此时,单一模型并不能很好地提取数据中的相关特征信息。因此,本发明提供了一种集成方法,能够针对不同类型填料的气速进行较准确的预报。 

发明内容

本发明提供了一种散堆(乱堆)填料塔液泛气速的非线性建模和预报方法,针对以上所述的现有填料塔中液泛气速预报技术中存在的不足和缺陷,提出一种结合高斯混合模型(GMM)聚类和最小二乘支持向量机(LSSVM)方法的集成预报方法,可以较有效的针对液泛数据的特点,提取相关特征信息,提高液泛气速模型的预报准确度和可靠程度。 

一种散堆填料塔液泛气速的预报方法,包括: 

(1)将从以往的填料塔液泛实验中采集的相关液泛数据集作为一个样本,分析填料塔液泛的特征及影响因素,确定液泛气速预报模型的输入和输出。 

首先,选择影响液泛的流动参数作为模型的输入变量,主要有液相雷诺数ReL、斯托克斯数StL、伽利略数GaL及填料层厚度校正系数SB和填料球形度。输出变量为Lockhart-Martinelli参数(χ)。这些流动参数组成一个样本,用[xi,yi]表示,其中xi表示第i个样本的输入变量(每个样本包含5个输入变量),yi表示该样本的输出变量。模型的输出和输入变量的函数关系式可表示成下式: 

液泛数据来自国内1976~2007年的液泛相关文献,这些液泛数据源自不同塔径的填料塔,基本涵盖了目前工业现场常用的散堆填料。 

(2)采用高斯混合模型(GMM)聚类算法将液泛样本数据集聚类成若干个样本子集。 

假设n组液泛数据服从混合高斯分布,即样本数据可看成是从数个高斯分布中生成出来的,则GMM也可看成是由k个单高斯模型(SGM)组成的。则液泛数据样本xs的概率分布密度函数p(xs)为: 

p(xs)=Σi=1kp(i)p(xs|i)=Σi=1kαiNi(xs|μi,Σi)---(2)

其中,αi为混合权值,服从μi为均值,Σi为方差矩阵; 

而Ni(xsii)则表示第i个SGM的概率密度函数: 

Ni(xs|μi,Σi)=1(2π)m|Σi|exp[-12(xs-μi)TΣi-1(xs-μi)]---(3)

GMM聚类的核心则是找到一组这样的参数使它们所确定的概率分布生成这些给定的数据点的概率最大,这个概率为也称为似然函数。为计算方便,取其对数,得到对数似然函数: 

L(X|Θ)=logΠs=1nΣi=1kαiNi(xs|μi,Σi)

=Σs=1nlogΣi=1kαiNi(xs|μi,Σi)---(4)

其中,Θ=(θ1,...θk),θi=(αiii)。由于对数似然函数里有加和,不能直接用求导方法直接求得最大值。因此本发明中采用期望最大(Expectation Maximum,EM)算法根据式(4)对GMM的参数αi, μi和Σi进行估计。 

EM算法分为E步和M步。E步是计算对数似然函数的期望,M步选择使期望最大的参数,并将所选择的参数代入E步,进行反复计算。 

E步:估计数据由第k个成分生成的概率为,即样本属于第k个成分的后验概率: 

ηsi=αiNi(xs|μi,Σi)Σi=1kαiNi(xs|μi,Σi)---(5)

M步:估计参数αi,μi和Σi的计算公式如下: 

αi=Σs=1nηsin---(6)

μi=Σs=1nηsixjΣi=1nηsi---(7)

Σi=Σs=1nηsi(xs-μi)(xs-μi)TΣs=1nηsi---(8)

GMM是一个迭代过程,通过不断更新式(6)、式(7)与式(8),使式(4)计算的对数似然函数达到目标要求,然后停止迭代,输出此次迭代的概率分布密度函数p(x),即为液泛数据训练集最终的概率分布密度函数。基于概率分布密度函数使液泛数据训练集聚类成k个子集。 

(3)对聚类后的每个样本子集单独进行学习训练,建立基于最小二乘支持向量机(LSSVM)的预报液泛气速的子模型。 

聚类后的k个液泛数据样本子集单独进行训练,分别建立液泛气速模型。本发明采用核学习建模技术,先采用一个非线性映射φ:x→H将数据集的输入映射到高维特征空间中,后在此高维特征空间中估计回归函数。对于液泛数据集{xi,yi}(i=1,2,...,n),其模型可表示成: 

f(xi)=wTφ(xi)+b         (9) 

式中xi∈Rd,yi∈R,为模型的权值向量,b为模型的偏置。根据结构风险最小化原则,使LSSVM在优化目标中的损失函数为误差ξi的二次项,构造其最优化问题: 

minJ(w,b,ξi)=12||w||2+γ2Σi=1nξi2---(10)

s.t.yi=wTφ(xi)+b+ξi,i=1,2,...,n              (11) 

其目标函数由拟合误差平方和项以及正则化项组成。其中,通过调节γ,可控制模型的拟合能力。 

用拉格朗日法构造优化问题,并根据优化条件求解优化问题,可得如下方程组: 

0vTvΩ+1γIbβ=0y---(12)

式中:y=[y1,…,yn]T,v=[1,…,1]T,β=[β1,…,βn]T,Ωik=φ(xi)Tφ(xk)(i,k=1,...,n)是满足Mercer条件的核函数。由式(15)求解可得到β与b,最终可得函数估计的液泛气速LSSVM模型: 

f(x)=Σi=1nβiK(x,xi)+b---(13)

本发明中核函数采用径向基(RBF)核函数,其数学表达式为 通过引入核函数,将高维特征空间中的点积运算转化为低维空间的核函数运算,避免了在高维空间中直接计算非线性映射φ的困难。 

(4)分别使用这些子模型单独进行液泛气速的预报,得到若干组预报结果,这些结果由一种加权最小二乘法(WLS)进行连接,得到的结果为填料液泛气速的最终预报值。 

由步骤(3)获得的k个液泛气速的子模型,分别对需要预报的填料进行预报,k组预报结果由一种连接方法进行连接融合,得到最后的预报结果。本发明采用WLS进行子模型结果的连接。k个子模型的预报结果的加权融合的实现过程为:根据最优原理推导出液泛气速子模型的输出残差与权系数的关系,建立权系数的优化模型: 

y^=hy+e---(14)

式中:为k个子模型输出分量构成的k维矢量;y为被估计的实际值;e=[e1,…,ek]T为各子模型输出分量的预报误差所构成的误差矢量;h=[1,…,1]T为k维常向量。 

WLS的估计准则是使加权误差平方和最小: 

Jr(y~)=(y^-hy~)TR(y^-hy~)---(15)

式中:R=diag(r1,r2,…,rc)为正定对角加权矩阵。即对求极小值,有: 

Jr(y~)y~=-hT(R+RT)(y^-hy~)---(16)

由此可得到y的最小二乘估计,即子模型的预报结果的融合结果: 

y~=(hTRh)-1hTRy^=Σi=1kriy^i/Σi=1kri---(17)

本发明中采用预报误差计算各组的权系数,即式中权系数δ∈(0,∞)为误差指数,控制子模型对集成模型的影响。从权系数计算式可知,子模型的预报误差越大,对应的权系数越小,该输出量对估计解的影响越低,可降低噪声数据的影响程度。 

本发明提出LSSVM的集成建模方法,对散堆填料塔中的液泛气速进行建模和预报,克服了传统液泛气速模型预报通用性差,准确度低的缺点。通过对液泛数据样本进行GMM聚类后再建模,能较好的提取液泛数据的特征信息,提高模型精度。采用WLS进行子模型的连接,能提高模型的通用性,有利于新型填料的液泛气速的预报。 

附图说明

图1是本发明对散堆填料塔中液泛气速的建模及预报的流程图; 

图2是本发明中GMM聚类算法的步骤流程图; 

具体实施方式

参照附图: 

一种散堆填料塔液泛气速的预报方法,包括: 

(1)将从以往的填料塔液泛实验中采集的相关液泛数据集作为一个样本,分析填料塔液泛的特征及影响因素,确定液泛气速预报模型的输入和输出。 

首先,选择影响液泛的流动参数作为模型的输入变量,主要有液相雷诺数ReL、斯托克斯数StL、伽利略数GaL及填料层厚度校正系数SB和填料球形度。输出变量为Lockhart-Martinelli参数(χ)。这些流动参数组成一个样本,用[xi,yi]表示,其中xi表示第i个样本的输入变量(每个样本包含5个输入变量),yi表示该样本的输出变 量。模型的输出和输入变量的函数关系式可表示成下式: 

液泛数据来自国内1976~2007年的液泛相关文献,这些液泛数据源自不同塔径的填料塔,基本涵盖了目前工业现场常用的散堆填料。 

(2)采用高斯混合模型(GMM)聚类算法将液泛样本数据集聚类成若干个样本子集。 

假设n组液泛数据服从混合高斯分布,即样本数据可看成是从数个高斯分布中生成出来的,则GMM也可看成是由k个单高斯模型(SGM)组成的。则液泛数据样本的概率分布密度函数为: 

p(x)=Σi=1kp(i)p(x|i)=Σi=1kαiNi(x|μi,Σi)---(2)

其中,αi为混合权值,服从μi为均值,Σi为方差矩阵;而Ni(xsii)则表示第i个SGM的概率密度函数: 

Ni(x|μi,Σi)=1(2π)m|Σi|exp[-12(x-μi)TΣi-1(x-μi)]---(3)

GMM聚类的核心则是找到一组这样的参数使它们所确定的概率分布生成这些给定的数据点的概率最大,这个概率为也称为似然函数。为计算方便,取其对数,得到对数似然函数: 

L(X|Θ)=logΠs=1nΣi=1kαiNi(xs|μi,Σi)

=Σs=1nlogΣi=1kαiNi(xs|μi,Σi)---(4)

其中,Θ=(θ1,...θk),θi=(αiii)。由于对数似然函数里有加和,不能直接用求导方法直接求得最大值。因此本发明中采用期望最大 (Expectation Maximum,EM)算法根据式(4)对GMM的参数αi,μi和Σi进行估计。 

EM算法分为E步和M步。E步是计算对数似然函数的期望,M步选择使期望最大的参数,并将所选择的参数代入E步,进行反复计算。 

E步:估计数据由第k个成分生成的概率为,即样本属于第k个成分的后验概率: 

ηsi=αiNi(xs|μi,Σi)Σi=1kαiNi(xs|μi,Σi)---(5)

M步:估计参数αi,μi和Σi的计算公式如下: 

αi=Σs=1nηsin---(6)

μi=Σs=1nηsixjΣi=1nηsi---(7)

Σi=Σs=1nηsi(xs-μi)(xs-μi)TΣs=1nηsi---(8)

GMM是一个迭代过程,通过不断更新式(6)、式(7)与式(8),使式(4)计算的对数似然函数达到目标要求,然后停止迭代,输出此次迭代的概率分布密度函数p(x),即为液泛数据训练集最终的概率分布密度函数。基于概率分布密度函数使液泛数据训练集聚类成k个子集。 

(3)对聚类后的每个样本子集单独进行学习训练,建立基于最 小二乘支持向量机(LSSVM)的预报液泛气速的子模型。 

聚类后的k个液泛数据样本子集单独进行训练,分别建立液泛气速模型。本发明采用核学习建模技术,先采用一个非线性映射φ:x→H将数据集的输入映射到高维特征空间中,后在此高维特征空间中估计回归函数。对于液泛数据集{xi,yi}(i=1,2,...,n),其模型可表示成: 

f(xi)=wTφ(xi)+b        (9) 

式中xi∈Rd,yi∈R,为模型的权值向量,b为模型的偏置。根据结构风险最小化原则,使LSSVM在优化目标中的损失函数为误差ξi的二次项,构造其最优化问题: 

minJ(w,b,ξi)=12||w||2+γ2Σi=1nξi2---(10)

s.t.yi=wTφ(xi)+b+ξi,i=1,2,...,n     (11) 

其目标函数由拟合误差平方和项以及正则化项组成。其中,通过调节γ,可控制模型的拟合能力。 

用拉格朗日法构造优化问题,并根据优化条件求解优化问题,可得如下方程组: 

0vTvΩ+1γIbβ=0y---(12)

式中:y=[y1,…,yn]T,v=[1,…,1]T,β=[β1,…,βn]T,Ωik=φ(xi)Tφ(xk)(i,k=1,...,n)是满足Mercer条件的核函数。由式(15)求解可得到β与b,最终可得函数估计的液泛气速LSSVM模型: 

f(x)=Σi=1nβiK(x,xi)+b---(13)

本发明中核函数采用径向基(RBF)核函数,其数学表达式为 通过引入核函数,将高维特征空间中的点积运算转化为低维空间的核函数运算,避免了在高维空间中直接计算非线性映射φ的困难。 

(4)分别使用这些子模型单独进行液泛气速的预报,得到若干组预报结果,这些结果由一种加权最小二乘法(WLS)进行连接,得到的结果为填料液泛气速的最终预报值。 

由步骤(3)获得的k个液泛气速的子模型,分别对需要预报的填料进行预报,k组预报结果由一种连接方法进行连接融合,得到最后的预报结果。本发明采用WLS进行子模型结果的连接。k个子模型的预报结果的加权融合的实现过程为:根据最优原理推导出液泛气速子模型的输出残差与权系数的关系,建立权系数的优化模型: 

y^=hy+e---(14)

式中:为k个子模型输出分量构成的k维矢量;y为被估计的实际值;e=[e1,…,ek]T为各子模型输出分量的预报误差所构成的误差矢量;h=[1,…,1]T为k维常向量。 

WLS的估计准则是使加权误差平方和最小: 

Jr(y~)=(y^-hy~)TR(y^-hy~)---(15)

式中:R=diag(r1,r2,…,rc)为正定对角加权矩阵。即对求极小值,有: 

Jr(y~)y~=-hT(R+RT)(y^-hy~)---(16)

由此可得到y的最小二乘估计,即子模型的预报结果的融合结果: 

y~=(hTRh)-1hTRy^=Σi=1kriy^i/Σi=1kri---(17)

本发明中采用预报误差计算各组的权系数,即式中权系数δ∈(0,∞)为误差指数,控制子模型对集成模型的影响。从权系数计算式可知,子模型的预报误差越大,对应的权系数越小,该输出量对估计解的影响越低,可降低噪声数据的影响程度。 

如图1所示,一种散堆填料塔中液泛气速的预报方法包括以下步骤: 

(1)选择合适的模型输入与输出变量,建立数据样本X; 

由于液泛气速受多种因素的影响,如操作介质的物理性质,填料塔的操作条件以及填料的几何特征和物理性质等。因此,需要考虑这些影响因素,选出能充分描述液泛特征的各种参数作为模型的输入与输出变量。除此之外,为取得好的预报效果,需要从模型方法考虑以下几方面:(1)所选变量能产生良好的预报效果,尽量使预报误差较小;(2)所选变量应与液泛现象有重要关联;(3)所选变量在达到以上两点要求的前提下,应选择最少数目的变量。综合考虑这些影响因素,最终选择以下特征变量为模型的输入与输出: 

输入变量:液相雷诺数ReL、斯托克斯数StL、伽利略数GaL及填料层厚度校正系数SB和填料球形度; 

输出变量:Lockhart-Martinelli参数(χ)。 

本发明的液泛数据从国内1976~2007年的液泛相关文献中采集而来,这些液泛数据源自不同塔径,不同填料的液泛实验,基本涵盖了目前工业现场常用的散堆填料,包含了建模和预报所用的共440组液泛数据。为了使所建模型能克服传统检验公式的缺陷,提高其通用性,本发明采用传统的填料数据为液泛气速进行建模,并预报某些新 型填料。模型的训练数据取数据集中的333组,剩余的107组液泛数据则为测试模型所用。 

(2)对液泛数据的训练集进行GMM聚类,将训练集分为若干子集; 

由于采集的液泛数据的多样性与不平衡性,本发明用GMM对液泛数据进行聚类。采用GMM聚类方法对训练数据进行聚类时,经过对液泛数据特点的分析,鉴于液泛数据的训练集主要包含了四种典型填料,故取聚类个数k=4,即最后训练集将被聚类成四个子集。GMM聚类的步骤为: 

步骤1:初始化参数αi,μi和Σi; 

步骤2:根据式(5)计算后验概率ηsi; 

步骤3:根据式(6)和式(7)分别更新αi和μi,后根据式(8)更新方差矩阵Σi; 

步骤4:根据式(2)计算概率分布密度函数p(x),再计算对数似然函数,当此次的函数值与前次的差值小于某个阈值ε时,停止迭代,输出p(x),同时输出各个SGM的概率密度函数N(x|μ,Σ);否则返回步骤2继续进行迭代计算。 

其流程如附图2。液泛样本数据经GMM聚类后,以不同的概率归为四类,使每一类集中的数据具有较为相似的信息特征。即GMM聚类较好的提取了液泛数据样本中的特征信息,能更好的建立液泛气速的预报模型。 

(3)聚类后的每个训练子集分别建立LSSVM模型; 

聚类后对4个子集建立模型,其模型如式(13)选用RBF核函 数K(x,xi)=exp(-||x-xi||22σ2).

为了使模型具有更高的通用性,本发明中正则化参数γ和宽度系数σ的确定采用交叉验证进行选择。首先确定正则化参数γ和宽度系数σ的选择范围分别为[10,1000]和[0.1,10],然后确定这两个参数的交叉计算的间距,在选择范围内的进行两参数的交叉计算。计算每次参数变化后模型的训练误差,具有最小训练误差的模型参数(即正则化参数γ和宽度系数σ)为最佳模型参数,以此建立LSSVM液泛气速模型。 

(4)采用WLS对子模型进行连接融合 

经过聚类后建立的子模型分别对测试数据进行预报,得到4组液泛气速预报结果。进一步,针对每个子模型的预报结果,本发明引进WLS进行连接融合,以提高单一LSSVM模型的预报精度和通用性。 

采用WLS进行连接时,需对权系数计算式中的误差指数δ进行正确的选择。若δ取太小,会放大预报性能差的子模型对集成模型的贡献,使集成模型的预报精度降低;δ太大,则使预报性能较好的子模型对集成模型贡献减少,同样会降低模型的预报效果。因此,结合步骤(3)中的正则化参数γ和宽度系数σ的选择,采用交叉验证选择LSSVM集成模型参数γ、σ以及δ,确定其范围分别为[1,1000]和[0.1,100]以及[1,5]。 

确定好参数后,根据式(17)对子模型的预报结果进行连接,得到最终的预报结果。 

(5)由集成模型得到的预报结果,通过液泛气速与Lockhart-Martinelli参数的关系式计算得到液泛气速的预报值。 

将本发明方法与文献中的两种典型的数据驱动方法比较: 

采用LSSVM的单模型方法,即用所有训练样本数据建立单一的LSSVM模型,参考文献为(Li,C.L.;Liu,Y.;Yang,J.;Gao,Z.L.Prediction of flooding velocity in packed towers using least squares support vector machine.In Proceedings of the10th World Congress on Intelligent Control and Automation,Beijing,China,July2012,pp:3226-3231.); 

采用BP-NN的单模型方法,即用所有训练样本数据建立单一的RBF-NN模型,参考文献为(Piche,S.;Larachi,F.;Grandjean,B.P.A.Flooding capacity in packed towers:database,correlations,and analysis.Ind.Eng.Chem.Res.2001,40,476–487.); 

比较结果如下表: 

表中,RMSE表示预报均方误差,AARE表示平均绝对相对误差,两个技术指标均是越小越好。从结果可知,本发明方法建立的液泛气速模型具有较高的预报精度,其预报效果优于目前已有的经验关联式等其它方法。 

上述要求的条件在填料塔的工业生产中均能够满足,而且实施简单,成本低。因此,本发明具有普遍性和通用性,能够为散堆(乱堆)填料塔的液泛气速提供较准确的模型和预报,供其设计和操作阶段使用,在保障塔设备的安全操作下,对提高过程效率具有重要的应用价值。 

去获取专利,查看全文>

相似文献

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

客服邮箱:kefu@zhangqiaokeyan.com

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

  • 服务号