版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领
文档简介
S系统辨识课程CHAPTER03SystemIdentification系统辨识第三章最小二乘参数辨识方法自动化系·系统辨识课程组2026年秋季学期S系统辨识CONTENTS本章内容概览Chapter3·最小二乘参数估计方法3.1最小二乘基本原理与几何意义LeastSquaresFundamentals3.2加权最小二乘法WeightedLeastSquares3.3递推最小二乘算法(RLS)RecursiveLeastSquares3.4增广最小二乘法(ELS)ExtendedLeastSquares3.5广义最小二乘法(GLS)GeneralizedLeastSquares3.6辅助变量法(IV)InstrumentalVariableMethod3.7Matlab仿真案例与本章小结SimulationCases&ChapterSummaryS系统辨识SECTIONDIVIDERSECTION3.1最小二乘基本原理从误差准则到正规方程LEASTSQUARESFORMULATIONLEASTSQUARES最小二乘问题的数学描述MathematicalFormulationofLeastSquaresProblem最小二乘估计的核心是最小化误差平方和准则函数J(θ)=||Y-Φθ||²,其解的存在唯一性依赖于数据矩阵Φ的列满秩条件,即ΦᵀΦ可逆。这一准则在统计学上对应于高斯噪声下的极大似然估计。01线性回归模型y(k)=φT(k)θ+e(k)其中φ(k)为k时刻的回归向量,θ为n维待估参数向量,e(k)为零均值随机噪声。02批量数据形式Y=Φθ+EY为N×1输出向量,Φ为N×n数据矩阵(每行为φT(k)),E为噪声向量,要求N≥n以保证信息充足。03最小二乘准则函数J(θ)=(Y−Φθ)T(Y−Φθ)04解的存在条件rank(Φ)=n→ΦTΦ可逆DERIVATION正规方程的推导过程最小二乘估计的解析解推导通过对二次型准则函数J(θ)求梯度∇J=‑2Φᵀ(Y‑Φθ)并令其为零,得到正规方程ΦᵀΦθ̂=ΦᵀY。该方程的物理意义是残差向量(Y‑Φθ̂)与数据矩阵Φ的列空间正交,即估计误差在数据张成的子空间上无分量。01展开准则函数J(θ)=YᵀY−2θᵀΦᵀY+θᵀΦᵀΦθ关于θ的二次凸函数,全局最小值点即为所求。02矩阵求导法则∂(aᵀθ)/∂θ=a∂(θᵀAθ)/∂θ=2Aθ据此得∇J(θ)=−2ΦᵀY+2ΦᵀΦθ。03令梯度为零−2ΦᵀY+2ΦᵀΦθ̂=0⇒ΦᵀΦθ̂=ΦᵀY此即正规方程(NormalEquation)。04最小二乘估计量θ̂LS=(ΦᵀΦ)⁻¹ΦᵀY前提是ΦᵀΦ非奇异;该式也称为伪逆解Φ⁺Y。S系统辨识LEASTSQUARES最小二乘估计的几何意义LS估计等价于将观测向量Y正交投影到数据矩阵Φ的列空间Col(Φ)上,预测值Ŷ=Φθ̂是Y在该子空间中的最佳逼近,残差ε=Y-Ŷ⊥Col(Φ)。列空间视角:Φ的n个列向量张成RN中的n维子空间S,Ŷ=Φθ̂∈S是Y在S上的正交投影点。正交条件:残差ε=Y-Ŷ与S中任意向量正交,即Φᵀε=0,这正是正规方程的几何表达。最优逼近性质:||Y-Φθ̂||≤||Y-Φθ||对任意θ成立,LS估计使预测误差的欧几里得范数最小化。投影矩阵:Ŷ=PY=Φ(ΦᵀΦ)⁻¹ΦᵀY,P=Φ(ΦᵀΦ)⁻¹Φᵀ称为投影算子,满足P²=P且Pᵀ=P。Y在Col(Φ)子空间上的正交投影示意图S系统辨识STATISTICALPROPERTIESProperties最小二乘估计的统计性质当噪声e(k)零均值且与回归向量φ(k)不相关时,LS估计θ̂是无偏的;若e(k)还为白噪声,则θ̂是有效估计(BLUE)。但若噪声有色或与输入相关,E[φ(k)e(k)]≠0导致偏差。四大核心统计性质01无偏性条件若E[e(k)]=0且E[φ(k)e(k)]=0,则E[θ̂LS]=θ₀,估计量无偏;否则存在渐近偏差,估计结果系统性偏离真值。02协方差矩阵Cov(θ̂)=σ²(ΦᵀΦ)⁻¹,其中σ²为噪声方差;数据越丰富(ΦᵀΦ越大),估计精度越高,参数不确定性越小。03Gauss-Markov定理在白噪声假设下,LS估计是所有线性无偏估计中方差最小的,即最佳线性无偏估计(BLUE),具有统计最优性。04一致性当N→∞时,若(1/N)ΦᵀΦ→R(正定)且(1/N)ΦᵀE→0,则θ̂依概率收敛于真值θ₀,样本量增大使估计渐近精确。CASESTUDY算例:直线拟合的最小二乘解LeastSquaresSolutionforLinearFitting以一元线性回归y=ax+b为例,构造Φ=[x1],Y=[y₁…y₅]ᵀ,直接套用θ̂=(ΦᵀΦ)⁻¹ΦᵀY即可手算出参数。该算例验证了LS公式的可操作性。01数据准备给定5组观测数据点:{(1,2.1),(2,3.9),(3,6.2),(4,7.8),(5,10.1)}设模型y=ax+b,待估参数θ=[a,b]ᵀ02构造矩阵信息矩阵与观测向量:Φ=[[1,1],[2,1],[3,1],[4,1],[5,1]]Y=[2.1,3.9,6.2,7.8,10.1]ᵀ计算ΦᵀΦ=[[55,15],[15,5]]03求解参数逆矩阵与中间量:(ΦᵀΦ)⁻¹=[[0.2,−0.6],[−0.6,2.2]]/10ΦᵀY=[82.3,30.1]ᵀθ̂=[1.99,0.08]ᵀ04结果验证拟合直线方程:y=1.99x+0.080.046残差平方和J<0.15最大偏差SECTION3.2加权最小二乘法处理非均匀噪声与数据可信度差异WEIGHTEDLEASTSQUARES加权最小二乘的原理与公式WLS通过引入正定权重矩阵W修正普通LS,准则函数JW=(Y-Φθ)TW(Y-Φθ)。当W=R-1(R为噪声协方差阵)时,WLS达到Cramér-Rao下界,成为最优线性无偏估计。WeightedLeastSquares—Theory&Derivation01加权准则函数JW(θ)=Σwk[y(k)−φT(k)θ]²=(Y−Φθ)TW(Y−Φθ),其中W=diag(w₁,...,wN)为正定对角权重矩阵,对各时刻残差赋予不同权重。02加权正规方程对准则函数求梯度∇JW=−2ΦTW(Y−Φθ)=0,得ΦTWΦθ̂WLS=ΦTWY,解为θ̂WLS=(ΦTWΦ)−1ΦTWY。03最优权重选择若噪声协方差Cov(E)=R已知,取W=R−1可使估计方差最小,此时WLS即为Markov估计,在线性无偏估计类中达到最优效率。04与普通LS的关系当W=I(单位阵)时WLS退化为普通最小二乘;当W≠I时,WLS对不同精度数据给予差异化信任度——方差越小的观测获得越大的权重。S系统辨识PROPERTIES加权最小二乘的性质与应用场景核心结论:WLS估计量仍为线性无偏,协方差Cov(θ̂WLS)=(ΦᵀWΦ)⁻¹ΦᵀWRWΦ(ΦᵀWΦ)⁻¹。当W=R⁻¹时简化为(ΦᵀR⁻¹Φ)⁻¹,达到理论下界。无偏性保持只要E[E]=0且E与Φ不相关,无论W如何选取,θ̂WLS始终是无偏估计。协方差表达式Cov(θ̂WLS)=(ΦᵀWΦ)⁻¹ΦᵀWRWΦ(ΦᵀWΦ)⁻¹,仅当W∝R⁻¹时取最小值。异方差处理当Var(e(k))=σk²不等时,取wk=1/σk²可有效抑制高噪声数据的影响。工程实例GPS定位按卫星仰角赋权;电力系统状态估计按量测精度加权;机器人多传感器融合。SECTION3.3递推最小二乘算法在线参数估计的核心引擎ALGORITHM为什么需要递推算法批量LS计算量O(N³)且需存储全部历史数据,不适用于在线场景。RLS将估计分解为θ̂(k)=θ̂(k-1)+K(k)[y(k)-φᵀ(k)θ̂(k-1)],每步仅需O(n²)运算和固定内存。01批量LS瓶颈每次新增数据需重新计算(ΦᵀΦ)⁻¹,计算量O(N³),存储O(Nn),无法满足实时性要求。02递推思想利用矩阵求逆引理,将k时刻的逆矩阵用k-1时刻的逆矩阵表示,避免重复求逆运算。03在线更新结构θ̂(k)=θ̂(k-1)+K(k)·ε(k),其中ε(k)=y(k)-φᵀ(k)θ̂(k-1)为新息,K(k)为增益向量。04适用场景自适应控制、实时信号处理、时变系统跟踪、嵌入式在线辨识等数据流式到达的场合。DERIVATIONSYSTEMIDENTIFICATIONRLS算法的数学推导利用矩阵求逆引理(A+BCD)⁻¹=A⁻¹−A⁻¹B(C⁻¹+DA⁻¹B)⁻¹DA⁻¹,导出P(k)=[I−K(k)φᵀ(k)]P(k−1)。配合K(k)=P(k−1)φ(k)/[1+φᵀ(k)P(k−1)φ(k)],构成完整RLS递推体系。01矩阵求逆引理(A+uvᵀ)⁻¹=A⁻¹−A⁻¹uvᵀA⁻¹/(1+vᵀA⁻¹u),这是RLS推导的数学基础,将矩阵逆运算转化为递推形式。02协方差递推P(k)=P(k−1)−P(k−1)φ(k)φᵀ(k)P(k−1)/[1+φᵀ(k)P(k−1)φ(k)],避免每步重新求逆,显著降低计算复杂度。03增益向量K(k)=P(k−1)φ(k)/[1+φᵀ(k)P(k−1)φ(k)],标量分母保证数值稳定性,增益随信息积累逐步衰减。04参数更新θ̂(k)=θ̂(k−1)+K(k)[y(k)−φᵀ(k)θ̂(k−1)],新息驱动修正,收敛速度由K(k)调节。H系统辨识ALGORITHMSUMMARYRLS算法总结与初始化RLS三方程构成闭环:K(k)=P(k-1)φ(k)/[λ+φᵀ(k)P(k-1)φ(k)],θ̂(k)=θ̂(k-1)+K(k)ε(k),P(k)=[I-K(k)φᵀ(k)]P(k-1)/λ。初值P(0)=αI(α≫1)加速收敛。算法三要素增益K(k)决定修正幅度,新息ε(k)提供修正方向,协方差P(k)反映估计不确定性。初值设定θ̂(0)=0或先验估计;P(0)=αI,α=10³~10⁶,大初值使早期数据权重高、收敛快。遗忘因子扩展引入λ∈(0,1],准则变为Σλ^{k-i}ε²(i),λ<1时旧数据指数衰减,可跟踪时变参数。数值稳定性直接RLS可能因舍入误差丧失P的正定性,工程上采用UD分解或平方根RLS保证鲁棒性。计算复杂度每步O(n²)浮点运算,n为参数维数;相比批量O(N³)显著降低,适合实时嵌入式应用。递推公式速查K(k)=P(k−1)φ(k)/[λ+φᵀ(k)P(k−1)φ(k)]θ̂(k)=θ̂(k−1)+K(k)·ε(k)P(k)=[I−K(k)φᵀ(k)]P(k−1)/λ初始化建议SECTION3.4增广最小二乘法处理有色噪声的偏差补偿策略SystemIdentificationTHEORY有色噪声导致的LS偏差问题在ARMAX模型中,等效噪声v(k)=C(z⁻¹)e(k)为有色噪声,且回归向量φ(k)含历史输出y(k-i),导致E[φ(k)v(k)]≠0。LS估计θ̂LS渐近有偏。ELS通过增广参数向量同时辨识系统与噪声模型来消除偏差。ARMAX模型:A(z⁻¹)y(k)=B(z⁻¹)u(k)+C(z⁻¹)e(k),e(k)~WN(0,σ²),C(z⁻¹)=1+c₁z⁻¹+…+cncz⁻ⁿᶜ偏差根源:将模型改写为y(k)=φᵀ(k)θ+v(k),其中v(k)=C(z⁻¹)e(k)为有色噪声,且φ(k)含y(k-i)与v(k)相关LS失效:E[φ(k)v(k)]≠0⇒plimθ̂LS≠θ₀,即使N→∞估计仍不收敛到真值ELS对策:将噪声参数ci纳入θ,构造增广回归向量φe(k)=[-y(k-1),…,u(k-nb),ê(k-1),…]ᵀ,使新噪声近似白噪声H系统辨识ALGORITHMRECURSIVEEXTENDEDLEASTSQUARES增广最小二乘算法实现RELS将噪声参数ci增广到θe,回归向量φe(k)包含残差估计ê(k−i)=y(k−i)−φ̂ᵀ(k−i)θ̂(k−i−1)。由于ê依赖θ̂,算法呈非线性耦合;采用一步滞后近似后仍可套用RLS框架递推更新。01增广参数向量θe=[a₁,…,anₐ,b₁,…,bnb,c₁,…,cnc]ᵀ,维度ne=na+nb+nc。将ARMA噪声模型参数与系统参数统一为一个向量进行辨识。02增广回归向量φe(k)=[−y(k−1),…,−y(k−na),u(k−1),…,u(k−nb),ê(k−1),…,ê(k−nc)]ᵀ。包含输入输出历史及残差历史三部分。03残差在线估计ê(k)=y(k)−φ̂eᵀ(k)·θ̂e(k−1),用当前预测误差近似真实噪声项,代入下一步回归向量,形成在线递推闭环。04递推更新完全套用RLS三方程,仅将φ替换为φe、θ替换为θe;初值Pe(0)=αI,α取大值加速收敛。SECTION3.5广义最小二乘法基于数据预滤波的两阶段估计ΣΦ⁻¹θ̂εALGORITHMSYSTEMIDENTIFICATION广义最小二乘的原理与步骤GLS采用两阶段解耦策略:Stage1用LS得θ̂¹并计算残差ê;Stage2对ê拟合AR模型得Ĉ(z⁻¹);Stage3用Ĉ⁻¹预滤波(y,u)后重做LS得θ̂²。迭代至收敛可获得一致估计。1Stage1·初步估计对原始数据(y,u)做普通LS,得θ̂¹,计算残差ê(k)=y(k)−φᵀ(k)θ̂¹2Stage2·噪声建模对残差序列ê(k)拟合AR模型Ĉ(z⁻¹)ê(k)≈e(k),用LS估计ĉᵢ参数3Stage3·预滤波+重估定义ỹ(k)=Ĉ⁻¹(z⁻¹)y(k)、ũ(k)=Ĉ⁻¹(z⁻¹)u(k),对(ỹ,ũ)做LS得θ̂²迭代收敛将θ̂²作为新初值重复Stage1-3,通常2-3次即可收敛;收敛解等价于极大似然估计GLS迭代过程中参数逐步收敛到真值H系统辨识COMPARISONAnalysisGLS与ELS的比较分析ELS将噪声参数增广到单一回归框架,每步O(ne²)在线更新,适合实时但存在非线性耦合风险;GLS通过预滤波解耦系统与噪声估计,理论性质更优且等价于MLE,但需多次批量LS,适合离线高精度场景。算法结构与计算特性ELS:单回路递推,每步更新增广参数向量,计算量O((na+nb+nc)²),内存固定,天然支持在线运行。GLS:多阶段迭代,每轮需两次批量LS(系统+噪声),计算量O(N·n²)每轮,通常需2-3轮收敛,适合离线批处理。理论性质与工程适用性ELS:在持续激励和适当条件下具有一致性,但因回归向量含估计残差,严格收敛证明较复杂;工程中广泛使用。GLS:收敛解等价于极大似然估计,具有渐近有效性和正态性;理论保证更强,常用于基准验证和高精度离线辨识。计算模式ELS在线递推vsGLS离线迭代理论保证ELS一致性vsGLS渐近有效性适用场景ELS实时控制vsGLS高精度辨识S系统辨识SECTIONDIVIDERSECTION3.6辅助变量法InstrumentalVariableEstimation22PRINCIPLESYSTEMIDENTIFICATION辅助变量法的基本原理IV估计量θ̂IV=(ZTΦ)−1ZTY通过引入辅助变量矩阵Z消除偏差。Z须满足两个条件:①与噪声渐近不相关plim(1/N)ZTE=0;②与回归矩阵渐近满秩plim(1/N)ZTΦ=RZΦ非奇异。01IV估计量定义θ̂IV=(ZTΦ)−1ZTY,其中Z为N×n辅助变量矩阵,结构与Φ相同但元素不含噪声污染。02条件一(正交性)plim(N→∞)(1/N)ZTE=0,即辅助变量与噪声渐近不相关,保证估计无偏。03条件二(相关性)plim(N→∞)(1/N)ZTΦ=RZΦ非奇异,即辅助变量与回归向量充分相关,保证可辨识。04一致性证明θ̂IV=θ₀+(ZTΦ)−1ZTE→θ₀(依概率),仅需上述两条件,不依赖噪声统计特性。S系统辨识METHODOLOGYInstrumentalVariables辅助变量的构造方法辅助变量构造是IV方法的关键。常用策略包括:①延迟输入z(k)=[u(k-d),...,u(k-d-nb)]ᵀ,简单可靠但要求输入充分激励;②模型预测输出ŷ(k)=φ̂ᵀ(k)θ̂(k-1),自适应性强但依赖初值;③外部参考信号r(k),适用于闭环系统。延迟输入法z(k)=[u(k-nk-1),...,u(k-nk-nb)]ᵀ,nk≥max(na,nb),利用输入与噪声独立性,最简单可靠。模型预测法z(k)=[-ŷ(k-1),...,-ŷ(k-na),u(k-1),...,u(k-nb)]ᵀ,ŷ由当前θ̂生成,自适应跟踪但需良好初值。外部参考信号闭环系统中用设定值r(k)或其延迟作为辅助变量,避免反馈导致的输入-噪声相关性。最优IV理论最优Z*=E[Φ|无噪声],不可得;实践中上述方法是其近似,性能取决于逼近程度。SECTION3.7仿真案例与本章小结从理论到实践的完整闭环SimulationSetupEXPERIMENT仿真实验设置与真实系统仿真采用二阶ARX模型A=[1,-1.5,0.7],B=[0,1,0.5],白噪声σ²=0.1,4阶M序列输入(周期15,幅值±1)确保持续激励。该设置覆盖了典型动态系统特征,适合全面评估各辨识算法性能。SystemModel真实系统A(z⁻¹)=1−1.5z⁻¹+0.7z⁻²B(z⁻¹)=z⁻¹+0.5z⁻²极点0.75±0.312j,阻尼比0.6,自然频率0.84rad/sNoiseConfig噪声设置e(k)~N(0,0.1),SNR≈15dB,模拟中等测量噪声环境另设有色噪声场景e_c(k)=0.8·e_c(k−1)+w(k)用于偏差测试InputSignal输入信号4阶M序列,长度L=15,幅值±1自相关近似δ函数,满足持续激励条件PE(n=4)DataScale数据规模采集N=500个样本,前50个初始化,后450个用于估计与验证重复50次MonteCarlo实验取平均ANALYSISLS与RLS参数收敛对比白噪声下LS与RLS均无偏一致:RLS在~100步内收敛至真值±2%范围,批量LS需完整数据集。两者稳态MSE差异<1%,验证理论等价性。RLS早期收敛速度受P(0)影响大。01收敛速度RLS在k≈80时参数误差降至5%以内,k≈150时稳定在1%以内;批量LS仅在N=500时给出单点估计。02稳态精度50次MC平均,LS的â₁=-1.498±0.012,RLS的â₁=-1.497±0.013,两者偏差<0.2%,方差相当。03初值敏感性P(0)=10⁶I时前20步出现明显超调;P(0)=10²I时收敛平滑但需200步;推荐10³~10⁴折中。04计算效率RLS每步0.02ms(n=4),批量LS每次12ms(N=500);在线场景RLS优势无可替代。COMPARISONSYSTEMIDENTIFICATION有色噪声下ELS/GLS/IV性能对比有色噪声(AR(1),a=0.8)下LS偏差达12%,三种改进方法均有效去偏:GLS精度最优(偏差<1%),ELS次之(偏差~2%,计算快5×),IV性能依赖辅助变量选择。工程建议:在线用ELS/IV,离线高精度用GLS。LS基线偏差有色噪声下â₁_LS=-1.32(真值-1.5),相对偏差12%;b̂₁_LS=0.88(真值1.0),偏差12%。普通LS在有噪声色时完全失效。ELS扩展最小二乘â₁=-1.47±0.02,偏差约2%,收敛需~200步。计算量仅为GLS的1/5,适合在线实时辨识场景。GLS广义最小二乘â₁=-1.495±0.008,偏差<1%,2轮迭代收敛。精度接近Cramér-Rao下界,但计算开销大,适合离线高精度分析。IV辅助变量法延迟输入IV偏差~3%;模型预测IV偏差~1.5%但依赖初值;最优IV(真实无噪声输出)偏差<0.5%。性能高度依赖辅助变量选取。在线场景:推荐ELS或IV,计算快、可递推离线场景:推荐GLS,精度逼近理论下界核心发现:有色噪声使LS完全失效,必须去偏VALIDATIONSYSTEMIDENTIFICATION模型验证与残差分析模型验证的核心是检验残差ε(k)=y(k)-ŷ(k)是否满足白噪声假设:①自相关R_ε(τ)≈0(τ≠0);②互相关R_uε(τ)≈0(∀τ)。若超出95%置信界±1.96/√N,表明模型未充分提取信息。残差白噪声检验计算Rε(τ)=Σε(k)ε(k−τ)/N,τ=1,…,20;95%置信界±1.96/√N,超限提示模型存在未捕获的动态结构缺陷。输入-残差不相关检验Ruε(τ)=Σu(k)ε(k−τ)/N应在全滞后范围内位于置信界内,否则表明存在未建模动态,模型对输入信息的提取不充分。拟合优度指标FIT%=100×(1−‖y−ŷ‖/‖y−mean(y)‖),>80%通常可接受;FPE/AIC准则用于模型阶次选择,
温馨提示
- 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
- 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
- 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
- 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
- 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
- 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
- 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。
最新文档
- 2026年中医耳鼻喉科肺经郁热鼻窒辨证技能测试卷及答案
- 林区火灾预防方法和扑救对策培训
- 起升机构操作要领及安全技术培训
- 汽车加油站火灾危险性及对策培训
- (2026年)医院放射防护管理规章制度
- (2026年)学校防欺凌信息员制度
- (2026年)物业管理公司走访制度
- 国内Geo优化服务商:2026年市场格局与三大代表平台深度解析
- 2025年河南省焦作市中站区数学四下期中检测试题含解析
- 2025年河南省南阳市方城县部分校数学四下期中复习检测模拟试题含解析
- 能源管理体系培训课件教学
- 元器件焊接技术
- 2025年黎明职业大学辅导员考试笔试题库附答案
- TJSTJXH5-2022高延性混凝土加固技术规程
- 医疗康复科操作礼仪要点
- DB31∕T 618-2022 电网电能计量装置配置技术规范
- 绿色企业能源公司企业管理制度
- GB/T 21387-2025供水系统用轴流式止回阀
- 设备除锈与刷漆标准规范手册
- 铁路工务安全教育课件
- T-ZZB 2977-2022 毛纺精梳机标准规范
评论
0/150
提交评论