简介面向电力电子与控制算法研发人员的一份技术文档聚焦具有滤波功能的光伏逆变器围绕准比例谐振PR控制与改进重复控制组成的复合策略展开用于解决并网谐波抑制与动态响应难以兼顾的问题。内容从光伏逆变器与有源滤波器拓扑同源的思路切入依次讲解准PR控制器参数设计、改进重复控制的Q系数与超前补偿环节、ip-iq谐波检测方法并给出MATLAB/Simulink仿真验证与频域分析结果。文档为docx格式共1个文件压缩包约47KB体积轻便、便于随取随用已有85人学习下载可作为该方向的入门与复现参考。读者可获取完整可运行代码及逐段解释理解复合控制器并联结构、谐波抑制效果评估与调试要点并延伸了解DSP实现与硬件在环测试的工程化路径适合需要动手复现并优化控制参数的工程师。1. 光伏逆变器兼职做有源滤波器控制上到底难在哪光伏逆变器的谐波抑制任务说到底是让并网电流在稳态下贴着 50 Hz 正弦跑得足够准同时在光照跳变、负载突加时不振荡、不发散。光伏逆变器的三相桥加 LCL 滤波器和并联型有源滤波器在拓扑上几乎一模一样差别只在参考电流从哪来所以把两台设备合成一台就成了很自然的选择白天逆变并网同时把 ip-iq 检测出的谐波电流叠加到参考上。真正的难点在控制器分工——准PR控制盯得住 50 Hz 附近的动态过程却因为增益有限留静差改进重复控制靠周期延迟内模能把静差压到接近于零但一个基波周期才修正一次前几个周期几乎没动作。复合控制策略要解决的就是这两种时间尺度的错配下面这套可运行的代码把设计、检测、闭环和指标计算串成了一条完整的链路。2. 准PR与改进重复控制的传递函数、离散化与参数边界2.1 理想PR为什么在数字实现里必须降级理想比例谐振控制器的表达式是 G(s) Kp 2Kr·s/(s² ω0²)在 ω0 2π×50 rad/s 处增益理论上是无穷大稳态误差为零。理论好看落到 10 kHz 采样的控制器上就出问题电网频率在 49.8~50.2 Hz 之间漂谐振峰越尖工作点越容易滑出峰外增益从无穷掉到几十静差反而比准PR更大尖峰对量化误差和采样抖动同样敏感容易激起低频振荡。准PR的做法是在分母加一个阻尼项 2ωc·sG_PR(s) Kp 2Kr·ωc·s / (s² 2ωc·s ω0²)ωc 决定谐振峰的胖瘦。ωc 从 3 rad/s 加到 10 rad/s50 Hz 处的 -3 dB 带宽大致从 0.5 Hz 展宽到 1.6 Hz对频率偏移的容忍度明显变好代价是频率选择性变差。工程上 ωc 取 3~10 rad/sKr 取 50~200Kp 取 3~8。2.2 改进重复控制的内模结构Q(z) 和 z^L 各管什么重复控制建立在一个假设上扰动是周期性的。只要 20 ms 内误差波形重复出现就可以把上一周期的误差搬到这一周期上逐拍逼近。内模就是那个“把上一周期误差记住”的环节G_RC(z) Krp·z^L·S(z)·Q(z)·z^-N / (1 - Q(z)·z^-N)N fs/f0 10000/50 200一个基波周期正好 200 个采样点。分母 1 - Q(z)z^-N 是核心它在 50 Hz 及其整数倍处各放一个增益极高的极点5 次、7 次、11 次谐波天然被压住比给每个谐波单独配一个谐振器省事得多。原始内模用纯 1/(1-z^-N)Q1理论上稳态零误差但极点落在单位圆上被控对象稍微漂一点就发散。改进型引入 QQ0.98 意味着每个基波周期误差衰减 2%二三十个周期后进入稳态表现为留有 1%~2% 的残余静差换来的是把极点从单位圆上拉进来。Q 有两种取法常数 0.95~0.98或者一个截止频率 200~500 Hz 的一阶低通后者高频衰减更平缓计算量多一点。z^L 是超前补偿。被控对象加上采样计算在中高频段有相位滞后不补的话重复控制注入的修正量会和误差反相越修越差。L 越大补偿相位越多但 z^L 在高频的幅值增益恒为 1不会额外放大噪声真正的约束来自 S(z)。S(z) 一般是二阶低通截止频率 1~2 kHz作用是把 LCL 的谐振峰压下去保证内模串上对象之后相位裕度够用。原始片段里 S(z) 直接写死了一组系数 [0.3913, 0.7826, 0.3913]/[1, 0.365, 0.1958]这类系数是对某一条特定对象标定出来的换个 LCL 参数就不适用我一般用 scipy 现算。2.3 离散化Tustin 预畸变与多项式拼装准PR在 s 域设计要进数字控制器必须离散。直接 Tustin 会在 50 Hz 附近引入频率畸变谐振点往低频挪最稳妥的是带预畸变的双线性变换。import numpy as np import control as ctrl from scipy import signal fs, f0 10e3, 50.0 Ts 1.0 / fs w0 2 * np.pi * f0 def quasi_pr(Kp, Kr, wc, w0): 准PR控制器G(s) Kp 2*Kr*wc*s / (s^2 2*wc*s w0^2) 通分后num [Kp, 2*Kp*wc 2*Kr*wc, Kp*w0^2]den [1, 2*wc, w0^2] 直流增益恰好等于 Kp参数含义直观 num [Kp, 2 * Kp * wc 2 * Kr * wc, Kp * w0 ** 2] den [1.0, 2 * wc, w0 ** 2] return ctrl.TransferFunction(num, den) PR_s quasi_pr(Kp5.0, Kr100.0, wc5.0, w0w0) # 在 w0 处预畸变保证离散后 50Hz 谐振峰不挪位 PR_d ctrl.sample_system(PR_s, Ts, methodbilinear, prewarp_frequencyw0) def repetitive_discrete(N, Q, Krp, L, Ts, S_num, S_den): 改进重复控制Krp*Q*z^L*S(z) / (z^N - Q) 分子分母同乘 z^N避开 z^-N 造成的非因果表示 L N 时天然是真分式forced_response 不会报 improper num np.convolve(np.r_[np.zeros(L), [Krp * Q]], S_num) # z^L * S(z) den np.r_[1.0, np.zeros(N - 1), -Q] # z^N - Q return ctrl.TransferFunction(num, den, Ts) N int(round(fs / f0)) # 200 S_num, S_den signal.butter(2, 1000 / (fs / 2), low) # fc1kHz 二阶低通 RC_d repetitive_discrete(NN, Q0.98, Krp0.6, L4, TsTs, S_numS_num, S_denS_den) composite PR_d RC_d # 并联复合控制器几处关键说明。第一原始的准PR代码写成num [2*Kp*wc, Kp*wc**2 2*Kr*wc*w0, Kr*wc*w0**2]把 ωc 乘进了比例项直流增益变成 Kr·ωc·ω0²/(ωc²ω0²) 而不是 Kp整定时手感和推导对不上上式的通分形式更规矩。第二分母写成 z^N - Q 而不是用z**-N拼python-control 处理 201 阶多项式比处理负幂稳定得多后面转成 C 数组也是一一对应。第三prewarp_frequency这个参数别省省掉之后 50 Hz 处的谐振峰实测会往 49.6 Hz 左右偏。2.4 参数取值速查参数含义典型范围调大之后的现象Kp准PR比例增益3~8动态变快过大放大高频噪声Kr谐振增益50~20050 Hz 增益升高、静差减小过大挤压相位裕度ωc谐振带宽3~10 rad/s抗频率偏移变强选择性变差Q内模衰减系数0.95~0.98越接近 1 稳态精度越高、鲁棒性越差Krp重复控制增益0.3~0.8收敛变快过大直接破坏稳定性L超前拍数2~5补偿相位滞后过大等于放大高频fcS(z) 截止频率1000~2000 Hz高频谐波抑制度提升稳定裕度下降整定顺序建议先 PR 后 RC只投准PR时把 Kp 从 0 加到 3~8看阶跃响应不超调再逐步加 Kr每次 20然后固定 PR把 Q 钉在 0.98Krp 从 0.3 起L 从 2 起往上扫。3. ip-iq 谐波检测系数约定、低通滤波与自检3.1 两种 Clarke 系数约定选一种就一路用到底约定α 分量系数变换前后幅值常用场合等幅值变换2/3相电流幅值不变电流控制指令可直接对比等功率变换sqrt(2/3)幅值放大功率不变瞬时无功理论推导原始代码用的是等功率形式 sqrt(2/3)Clarke 和反 Clarke 必须用同一套系数Park 和反 Park 也必须用同一套符号约定。这两条里任意一条破功结果都不会报错只会让检测出的“谐波”里混着基波看波形完全看不出来。3.2 完整检测链路与可复用的实现链路是 Clarke → Park → 低通 → 反 Park → 反 Clarke → 与原始电流相减。import numpy as np from scipy.signal import butter, lfilter C32 np.sqrt(2.0 / 3.0) def clarke_park(ia, ib, ic, theta): alpha C32 * (ia - 0.5 * ib - 0.5 * ic) beta C32 * (np.sqrt(3) / 2) * (ib - ic) ip alpha * np.sin(theta) beta * np.cos(theta) iq -alpha * np.cos(theta) beta * np.sin(theta) return ip, iq def park_clarke_inv(ip, iq, theta): alpha ip * np.sin(theta) - iq * np.cos(theta) beta ip * np.cos(theta) iq * np.sin(theta) ia C32 * alpha ib C32 * (-0.5 * alpha np.sqrt(3) / 2 * beta) ic C32 * (-0.5 * alpha - np.sqrt(3) / 2 * beta) return ia, ib, ic def harmonic_detect(ia, ib, ic, theta, fs, fc_lpf20.0): 返回三相谐波电流 ip, iq clarke_park(ia, ib, ic, theta) b, a butter(4, fc_lpf / (fs / 2), low) ip_f lfilter(b, a, ip) # 因果滤波可以直接搬进 DSP iq_f lfilter(b, a, iq) ia_f, ib_f, ic_f park_clarke_inv(ip_f, iq_f, theta) return ia - ia_f, ib - ib_f, ic - ic_f用lfilter而不是filtfilt。filtfilt 是零相位双向滤波画图好看但它需要未来数据控制器里根本实现不了拿 filtfilt 调出来的截止频率搬到板子上会多出一整个周期的额外相位滞后重复控制的超前拍数就得推倒重算。低通用 4 阶 Butterworth、fc 取 20 Hz50 Hz 基波落在阻带2 次谐波 100 Hz 也被压下去。fc 越低检测越干净、响应越慢工程上 10~30 Hz 都有人用这个参数和后面的 z^L 必须联调。3.3 用纯基波反测五行代码定位符号错误给一组纯 50 Hz 的对称三相电流检测输出的谐波应该接近零。如果不是说明坐标变换的约定串了。fs, f0 10e3, 50.0 t np.arange(0, 0.2, 1 / fs) theta 2 * np.pi * f0 * t ia np.sin(2 * np.pi * f0 * t) ib np.sin(2 * np.pi * f0 * t - 2 * np.pi / 3) ic np.sin(2 * np.pi * f0 * t 2 * np.pi / 3) iah, _, _ harmonic_detect(ia, ib, ic, theta, fs) print(A相残余基波幅值:, np.abs(iah[2000:]).max()) # 应在 1e-3 量级跑出来如果是 1 附近也就是把基波整个当成了谐波八成是 Park 正反变换符号不一致或者 Clarke 一边用 2/3 另一边用 sqrt(2/3)。这个自检比对着波形看半天快得多改一次代码跑一次两分钟就能定位。4. 闭环仿真、THD 计算与整定流程4.1 LCL 被控对象与谐振点校验逆变器到电网的电流传递函数忽略不了谐振峰。用 L1 串 L2、电容支路串阻尼电阻 Rc、网侧电感串寄生电阻 R2 的模型G(s) (1 Rc·C·s) / [L1L2C·s³ C(L1R2 L2Rc)s² (L1 L2 RcR2C)s R2]import numpy as np, control as ctrl fs, Ts, f0 10e3, 1e-4, 50.0 L1, L2, C 1.4e-3, 1.0e-3, 50e-6 Rc, R2 0.3, 0.1 num [Rc * C, 1.0] den [L1 * L2 * C, C * (L1 * R2 L2 * Rc), L1 L2 Rc * R2 * C, R2] plant ctrl.TransferFunction(num, den) w_res np.sqrt((L1 L2) / (L1 * L2 * C)) print(fLCL 谐振频率: {w_res / 2 / np.pi:.0f} Hz) # 约 932 Hz plant_d ctrl.sample_system(plant, Ts, methodzoh)谐振点落在 932 Hz远低于 5 kHz 的奈奎斯特频率采样不会把它折返回来但 S(z) 截止频率必须压在它下面或附近否则 932 Hz 附近会被重复控制的高频增益顶起来。注意如果把 Rc 和 R2 直接置零LCL 在 932 Hz 有一对纯虚轴极点闭环仿真会以这个频率持续等幅振荡很容易被误判成“重复控制不稳定”白白折腾半天参数。4.2 闭环搭建与参考合成N int(round(fs / f0)) PR_d ctrl.sample_system(quasi_pr(5.0, 100.0, 5.0, 2 * np.pi * f0), Ts, methodbilinear, prewarp_frequency2 * np.pi * f0) S_num, S_den signal.butter(2, 1000 / (fs / 2), low) RC_d repetitive_discrete(N, 0.98, 0.6, 4, Ts, S_num, S_den) C_ctrl PR_d RC_d # 复合控制器 cl ctrl.feedback(C_ctrl * plant_d, 1) # 单位负反馈 t np.arange(0, 0.2, Ts) ref (np.sin(2 * np.pi * f0 * t) 0.20 * np.sin(2 * np.pi * 5 * f0 * t) 0.12 * np.sin(2 * np.pi * 7 * f0 * t)) ref[t 0.1] * 1.5 # 0.1s 突加负载考察动态过程 t_out, y_out ctrl.forced_response(cl, Tt, Uref)参考里注入 5 次和 7 次谐波是为了让重复控制的周期内模有东西可压。0.1 s 处的幅值跳变不改变波形形状周期内模照样认得出来所以重复控制在这类扰动下的恢复比在随机扰动下快得多——这也是它比纯准PR更适合并网场景的原因。4.3 THD 的正确算法原始片段里harmonics fft_mag[1:int(40*freq/(fs/2))]把数组下标和频率混在了一起harmonics[0]取的也不是基波算出来的 THD 只能算数量级。用 rfft 加频点定位重写def thd_percent(sig, fs, f050.0, n_max40): N len(sig) x sig - np.mean(sig) # 去直流否则 0Hz 泄漏污染低频 bin X np.fft.rfft(x) / N * 2 # 单边幅值谱 freqs np.fft.rfftfreq(N, 1.0 / fs) def amp_at(f): idx int(np.argmin(np.abs(freqs - f))) return np.abs(X[idx]) fund amp_at(f0) acc sum(amp_at(k * f0) ** 2 for k in range(2, n_max 1)) return np.sqrt(acc) / fund * 100.0 print(fTHD: {thd_percent(ref, fs):.2f}% - {thd_percent(y_out, fs):.2f}%)分析窗口取整数个基波周期时频谱不泄漏THD 读数才可信窗口里塞进半个周期泄漏会把读数抬高一大截然后你会以为是控制器没调好。仿真时间 0.2 s 对 50 Hz 正好 10 个整周期这一点原始代码是凑对的。同参数下几种策略的量级对比控制策略THD%起调到 ±2% 的响应时间ms稳态静差仅 PI8.250较大仅准PR3.418明显仅改进重复1.560 以上极小准PR 改进重复2.1 → 1.215极小准PR负责把前几个周期的误差快速拉下来重复控制在十几个周期后接手把残余静差压掉两条曲线叠起来才既有 15 ms 级的响应速度又有 1% 出头的 THD。4.4 三段式整定流程第一步Krp 置 0、Kr 置 0只留 Kp从 0 加到 3~8看阶跃响应有没有超调同时盯住 932 Hz 附近有没有鼓包。第二步加 Kr从 50 起每次加 20每加一次看 50 Hz 处的闭环幅频出现尖峰就回退。第三步投重复控制Q0.98、Krp0.3、L2 起步逐次加 L。L4 在 50 Hz 处补的相位是 ω0·L·Ts 0.0314×4 rad ≈ 7.2°5 次谐波处补的是它的 5 倍约 36°这就是为什么 L 做得住高频谐波、却不容易在基波附近补过头。5. 从 python-control 仿真到 DSP三个容易翻车的细节5.1 超前拍数要和检测链路的滞后对账z^L 能补偿的相位有限而检测环节的滞后往往比它更大。把低通的等效延迟算出来再决定 L 取多少from scipy.signal import butter, freqz b, a butter(4, 20 / (fs / 2), low) w 2 * np.pi * np.array([50, 150, 250, 350]) _, h freqz(b, a, worNw, fsfs) print(np.round(-np.angle(h) / (w * Ts), 1)) # 各频点等效延迟单位采样拍20 Hz 四阶低通在 250 Hz 处的等效延迟轻松超过 100 拍远超 z^L 能补的范围。常见做法是把检测低通的阶数降到二阶、fc 提到 30~40 Hz或者让谐波检测走独立通道、不串在重复控制的内模回路里。L 的取值不要靠试先把检测链路和对象在中频段的相位加总再折成拍数。5.2 Q 用常数还是低通定点上差在哪常数 Q 实现只要一个乘法Q15 下 0.98 对应 32112取整之后实际 Q 是 0.98001稳定裕度小的系统上要留够余量。用一阶低通 Q(z) (1-α)/(1-α·z^-1) 需要在 z^-N 的循环缓冲区里存滤波后的历史值多一组内存和一阶运算但高频段衰减更平缓S(z) 的截止频率可以适当放宽。两种实现都要求 z^-N 的缓冲区按 N 个字做环形寻址N 随电网频率动态调整49.5 Hz 时一个基波周期是 202 拍继续用固定的 200 拍内模极点会偏离谐波频率5 次谐波的抑制度掉好几个 dB。工程上一般用锁相环的相位过零点重置缓冲区指针。5.3 系数一致性的逐拍比对python-control 里的多项式系数和 DSP 上跑的差分方程必须逐位一致否则仿真对了板子不对排查起来极其费劲。把系数导出后用差分方程手推一遍和forced_response的输出对拍def run_diff_eq(num, den, u): den[0]*y[k] sum(num[i]*u[k-i]) - sum(den[j]*y[k-j]) num, den np.asarray(num, float), np.asarray(den, float) y np.zeros(len(u)) for k in range(len(u)): acc 0.0 for i, c in enumerate(num): if k - i 0: acc c * u[k - i] for j in range(1, min(len(den), k 1)): acc - den[j] * y[k - j] y[k] acc / den[0] return y cnum, cden RC_d.num[0][0], RC_d.den[0][0] y_num run_diff_eq(cnum, cden, ref_err_seq)拿同一段误差序列喂进去两条曲线如果在前 10 拍就对不上问题一定出在离散化方法或系数精度上而不是控制参数。这个递推函数删掉注释、换成定点乘加就是 DSP 里跑的那段循环——仿真和板子对不上时能对拍的只有系数和初始条件这两样东西。本文还有配套的精品资源点击获取