带通滤波器电路分析.doc

带通滤波器电路分析与优化
一、实验目的
1．熟悉利用修正节点法求解电路传递函数的方法；
2．理解各种电路优化的方法，并能将其灵活应用于具体电路的优化；
3．掌握利用MATLAB工具来优化电路的方法，并能够用程序来实现优化方法，以利于对电路的研究。
二、实验要求
分析下面的带通滤波器电路

图1. 高通滤波器电路图
考虑电路的频率响应：
要求优化上面各个元件参数，使得有以下的理想响应：

图2. 理想带通滤波器的幅频响应
设计过程可以看成对下面的代价函数进行优化，求最小值：

（其中是上图给出的理想的频率响应）
可选的优化内容包括：
1 在元件值没有约束的情况下寻找最优的元件参数，使得归一化后的频率响应满足上图。
2 假设元件的值限制在下面的范围内：
,,
如何选择元件，满足设计要求
3 假设元件的取值是离散的，即：
1、1.1、1.2、1.3、1.5、1.6、1.8、2.0、2.2、2.4、2.7、3.0、3.3、3.6、3.9、4.3、4.7、5.1、5.6、6.2、6.8、7.5、8.2、9.1此时如何选择元件，满足设计要求
4 假设所有元件值有的误差，（简单起见，假设元件值在这个范围内均匀分布）分析此时代价函数的变化（代价函数成为一个概率分布）
5 设定一个代价函数变化极限，然后尝试反过来确定允许的元件误差范围
三、实验内容
1．问题分析
首先，利用修正节点法求解该电路的传递函数，其传递函数为，当输入为1时，所得到的输出即为传递函数的值。为了求得与理想带通滤波器的幅频响应较一致的滤波器，我们设定了如下代价函数：，由于带通滤波器在滤波器边缘变化较快，而其他地方变化较慢，我们应当对边缘处进行着重考虑。因此，我们最终选定如下代价函数：，其中为频点处的权重，经过基本的筛选和分析，最后采取在之间的权重为15，之间权重为50，之间权重为100，其余频率段的权重为1。
2．代价函数的求取
首先，要求得该电路传递函数的值，可以有两种方法供选择，第一是通过修改第一次大作业的读网表文件来求取该值，此方法通用性较好，可以在分析不同的电路时较为方便的得到结果。第二是直接列出修正节点法方程进行求取，此方法通用性不好，只能针对此电路，但是对于本次所研究的对象，表达较为简单，实现起来较为容易。综合考虑，我们选择用第二种方法来求取。
在实现时，首先进行数据转换，即电阻转为电导，电容和电感也进行相应的转换。将转换后的数据保存在xteg中，而修正节点法电路表达式中T矩阵的元素保存在数组T中，通过p=inv(T)*b;求得电路参数的值。为了得到多个频点的值，可用如下方法得到不同频率的传递函数的值：

10

f_start=10000;
f_step=10000;
f_end=8000000;
for f=fstart:f_step:f_end
s=2*pi*f*j;
i=f/10000;
……
H(i)=abs(p(4))/abs(p(1));

其中，H(i)为传递函数。
代价函数的求取要用到所求得的传递函数和理想幅频响应，根据以上结论，可通过如下方法实现：

for i=1:800
w(i)=1;
end
for i=1:800
……%此处为权重的设置
end
end
for i=1:800
a(i)=w(i)*(H(i)-Hteg(i))^2;
end
for i=1:800
G=G+a(i);
end

该程序主要实现了权重的设置，和通过a(i)=w(i)*(H(i)-Hteg(i))^2;和G=G+a(i);语句进行代价函数的求取。完整的代码见附录fun.m。
3．在元件值没有约束的情况下寻找最优的元件参数，使得归一化后的频率响应满足要求。
1） 方法：此部分的关键在于最优元件参数的求取，使用MATLAB自带的优化工具箱的fminsearch函数进行搜索，使用途中所给的参数做初值，通过如下语句进行搜索：
options=optimset('MaxFunEvals',100000,'MaxIter',100000);%函数参数设定
x=fminsearch(@(x)fun(x),x0,options);
具体的程序见附录中main.m。
2） 结果：得到的x可通过调用G=fun(x);求取代价函数的值，得到的最优元件值和幅频响应如下：
表1. 无约束优化下最优元件参数
元件
R1
R2
C1
C2
数值
20．6918469451
21．9290683997
30304.4590029536
225.4344914925
元件
C3
L1
L2
L3
数值
21351.3454683340
0.0592719030
8.0145017495
0.0844666000
其中，电阻的单位为Ω，电容的单位为pF，电感为μH.

图3. 无约束优化下的幅频响应
3） 结论：通过电路的优化，得到了较为理想的幅频响应，其边缘的陡峭性也较为理想，但是中间出现了一定的纹波，为了使纹波得到改善，我们尝试了在代价函数中将抖动大的频率段加更大的的权重或者取更多的频点，这些方法都能使纹波得到改善。但抖动和上升沿缓慢是不能同时避免的，最后综合考虑了各种因素，选定了该权重。
4．假设元件的值限制在下面的范围内：，，。如何选择元件，满足设计要求
1） 方法：思路类似无约束条件，但是需要对元件的取值做一定的约束。设定元件取值的上下限：
lb=[0.1,0.1,0.01,0.01,0.01,0.01,0.01,0.01];%下限设定
ub=[10^6,10^6,10^7,10^7,10^7,100,100,100];%上限设定
选用优化工具箱的有约束条件优化函数fmincon来完成该优化：
options=optimset('MaxFunEvals',100000,'MaxIter',100000,'TolCon',1e-9,'TolFun',1e-9);
x=fmincon(@(x)fun(x),x0,A,b,Aeq,Beq,lb,ub,nonlcon,options);
2） 结果：得到的x可通过调用G=fun(x);求取代价函数的值，得到的最优元件值和幅频响应如下：
表2. 元件参数有约束时的最优元件参数
元件
R1
R2
C1
C2
数值
184.28110122931
193.88084852187
3006.24782995410
25.76876139846
元件
C3
L1
L2
L3
数值
3006.77276419080
0.59849245686
70.1456325938
0.59823789471

图4. 有取值范围约束下的最优幅频响应
3） 结论：通过对比图3和图4的结果可以看出来，有约束时的幅频响应不如无约束时的结果，但是两者相差不大，可见此时对传递函数的幅频响应没有太大的影响。但是通过表1和表2的对比可以看出来，此时元件的取值有较大的差别。
5．假设元件的取值是离散的，即：1、1.1、1.2、1.3、1.5、1.6、1.8、2.0、2.2、2.4、2.7、3.0、3.3、3.6、3.9、4.3、4.7、5.1、5.6、6.2、6.8、7.5、8.2、9.1此时如何选择元件，满足设计要求
1） 方法：离散优化的主要想法是模仿模拟退火法，但由于模拟退火法中关于温度的设定更多时候是一个经验值，在此，我们就略作简化。整个离散优值的选取是在连续优值确定的情况之下，依据可选优值附近的离散值作为初值，对每一元素X（电阻电容电感）实施线性搜索迭代，不断在更优值的基础上获取，最终获取在这一搜索范围内的最优离散值。核心实现方式主体如下：

for n=1:3 %综合考虑时间与效率设定循环搜索三周
x=y;
for k=1:8
for i=1:8
for j=1:24
x(k)=model(j)*gain(i);
cur_val=fun1(x);
if(cur_val end
end
end
x(k)=y(k);%很重要
end

在此值得注意的是，离散最优值的搜索往往会陷入局部最优的僵局而造成无法获取整体结果最优的优化值，对此主要有两种解决办法：1）可选取更加良好的搜索方式2）对于搜索初值首先进行一定的优化，保证足够好的搜索方向和起点。
另外，对于优化函数频带特性和纹波要求显然是较为矛盾的要求，因而，我们选择了对纹波要求较高的初值再配合对频带要求较高的目标函数获取兼顾两者的优化值。
2） 结果：得到的x可通过调用G=fun(x);求取代价函数的值，得到的最优元件值和幅频响应如下：
表3. 元件参数离散时的最优元件参数
元件
R1
R2
C1
C2
数值
13．0000000000
16.0000000000
56000.0000000000
270.0000000000
元件
C3
L1
L2
L3
数值
11000.0000000000
0.0330000000
6.8000000000
0.2000000000

图5. 元件取值离散时的最优幅频响应
3） 结论：通过与图3和图4的结果对比可以看出来，元件取值离散时的幅频响应不如前两种情况时的结果，主要体现在幅频响应曲线的陡峭性不够好，但是纹波稍有改善。通过与表1和表2的对比可以看出来，此时元件的取值与前两种情况也有较大的差别。
6．假设所有元件值有的误差，（简单起见，假设元件值在这个范围内均匀分布）分析此时代价函数的变化（代价函数成为一个概率分布）
1） 方法：容差设计的任务是设计电路标称值与分配电路中元器件的参数容差，使电路性能的偏差最小，或者在保证电路性能满足指标要求的条件下允许元器件参数的容差范围最大。
统计过程中，我们采用蒙特卡洛分析法，先是使用随机数发生器按元件值的概率分布来选择元件，然后对整个输出量进行统计，得到最后的模拟分析结果。为了方便评估，我们选定无约束优化所得到的结果作为偏差的基准：
x1=1.0e+004 *[0.00206918469451 0.00219290683997 3.03044590029536 0.02254344914925 2.13513454683340 0.00000592719030 0.00080145017495 0.00000844666000];%无约束优化最佳值
主要通过如下语句来实现偏差的随机性：
x(j)=x1(j)*0.95+x1(j)*(1.05-0.95)*rand(1);
完整的代码见附件main4.m
2） 结果：求得代价函数并通过i=0:10:100;hist(G,i);语句做得如下分布图：

图6. 代价函数概率分布图
3 结论：由代价函数的概率分布可以看到，若元件值在最优解的基础上有一定的随机误差时，代价函数值将会增大。而且可以看到，代价函数的概率分布在大于10之后呈指数下降，也就是说当元件已取得最优解后，其本身的误差对电路造成极大影响的概率很小，并且系统误差越大，发生的概率越小，这也是我们进行电路优化的意义所在。
6．设定一个代价函数变化极限，然后尝试反过来确定允许的元件误差范围
1） 方法：由容差分析得出如下结论：当各独立元件存在一定偏差的情况下，电路的整体性能也将存在一个偏差。这在实际工业生产中十分有研究意义。我们为了保证某个大规模生产的产品存在一定的正品率，即代价函数应存在某个置信区间内。为了达到这样的效果，必须给所采用的独立元件设置一定的误差范围。
相对应于本次任务分析，我们分别对无约束条件和离散参数条件时所得到的最优解进行讨论。对于无约束条件，假设次品率为5%，即正品率为95%。定义正品率的范围为代价函数应在[0，15]范围内；对于离散参数，假设次品率为10%，即正品率为90%。定义正品率的范围为代价函数应在[0，25]范围内。为了保证数据的有效可靠，我进行了1000次统计，那么在允许的元件误差范围内就应该分别至少有950和900次满足上述条件。
这两种情况下分别使rate从3.5%和2%开始进行评估，如果满足条件，则减0.1%之后再次进行评估，直至条件不能够满足，所得到的rate即是所求的允许的原件误差范围。其中在误差范围内的元件值由如下代码实现：
for j=1:8
y(j)=x1(j)*(1-rate)+x1(j)*2*rate*rand(1);
end
而1000次测试进行评估由
if k<25
m=m+1;
end
进行统计，如果m大于950或者900则满足条件，进行下一次测试，否则得到结果。具体代码见附录main5.m。
2） 结果：经运算，无约束条件时允许的元件误差范围为±3.1%，离散取值条件下允许的元件误差范围为±1.6%。画出容差范围在±3.1%和±1.6%时的代价函数概率分布图分别如图7、图8所示：

图7. 容差范围±3.1%代价函数概率分布图 图8. 容差范围在±1.6%代价函数概率分布图
四、实验感想
经过本次试验，我们更加深入地了解了电路优化的各种方法和技巧，对各种方案的具体实现有了更深的理解。同时我们也掌握了利用MATLAB来实现电路优化的方法，并从多个角度来看电路优化的方法和成效，懂得如何评估电路的偏差。这次试验不仅让我们学到了很多知识，更重要的是让我们学到了考虑问题的角度和方法。最后感谢应忍冬老师课上认真的讲解和指导。
附录（源代码）
fun.m

function G=fun(x) %x=[R1,R2,C1,C2,C3,L1,L2,L3];
%H=Fr(x);%实际频响曲线 实际运行为加快速度合并调用函数
R1=x(1);R2=x(2);C1=x(3);C2=x(4);C3=x(5);L1=x(6);L2=x(7);L3=x(8);
xteg=[1/R1,1/R2,C1*10^(-12),C2*10^(-12),C3*10^(-12),L1*10^(-6),L2*10^(-6),L3*10^(-6)];% 数据转化 注意关于short,long型的设置
%初始化
f_start=10000;
f_step=10000;
f_end=8000000;
T=zeros(8,8);
b=zeros(8,1);b(8)=1;
H=zeros(1,800);
for f=f_start:f_step:f_end
s=2*pi*f*j;
i=f/10000;
T(1,1)=xteg(1);T(1,2)=-xteg(1);T(1,8)=1;
T(2,1)=-xteg(1);T(2,2)=xteg(1)+s*xteg(3);T(2,5)=1;T(2,6)=1;
T(3,3)=s*xteg(4); T(3,4)=-s*xteg(4);T(3,6)=-1;
T(4,3)=-s*xteg(4);T(4,4)=s*(xteg(4)+xteg(5))+xteg(2);T(4,7)=1;
T(5,2)=1;T(5,5)=-s*xteg(6);
T(6,2)=1;T(6,3)=-1;T(6,6)=-s*xteg(7);
T(7,4)=1;T(7,7)=-s*xteg(8);
T(8,1)=1;
p=inv(T)*b;
H(i)=abs(p(4))/abs(p(1));
end
Hteg=zeros(1,800);%理想频响曲线
w=zeros(1,800);%权重系数 权重的设置很大程度会影响程序速度
a=zeros(1,800);%暂时寄存
G=0;
for i=1:800
w(i)=1;
end
for i=1:800
if(i>=350)&&(i<=400)
Hteg(i)=0.5;
end
if(i>=350)&&(i<=370)
w(i)=15;
end
if(i>=380)&&(i<=400)
w(i)=50;
end
if(i>=370)&&(i<=380)
w(i)=100;
end
end
for i=1:800
a(i)=w(i)*(H(i)-Hteg(i))^2;
end
for i=1:800
G=G+a(i);
end
return

main1.m

tic;
%无约束优化
x0=[50,50,3000,150,3000,0.5,10,0.5];%初始值设定
options=optimset('MaxFunEvals',100000,'MaxIter',100000);%函数参数设定
x=fminsearch(@(x)fun(x),x0,options);
G=fun(x);
plot(Fr(x));axis([0 8000 0 0.55]);
disp(' ');disp(' ');
toc;

main2.m

tic;
%约束优化
lb=[0.1,0.1,0.01,0.01,0.01,0.01,0.01,0.01];%下限设定
ub=[10^6,10^6,10^7,10^7,10^7,100,100,100];%上限设定
A=[];b=[];Aeq=[];Beq=[];nonlcon = [];%函数参数设定
options=optimset('MaxFunEvals',100000,'MaxIter',100000,'TolCon',1e-9,'TolFun',1e-9);
x=fmincon(@(x)fun(x),x0,A,b,Aeq,Beq,lb,ub,nonlconoptions);
G=fun(x);
plot(Fr(x));axis([0 8000 0 0.55]);
toc;

main3.m

%离散取值优化 原数据有误
model=[1,1.1,1.2,1.3,1.5,1.6,1.8,2.0,2.2,2.4,2.7,3.0,3.3,3.6,3.9,4.3,4.7,5.1,5.6,6.2,6.8,7.5,8.2,9.1];%离散取值标称值
gain=[10^-1,10^0,10^1,10^2,10^3,10^4,10^5,10^6];%标称值放大倍数
x1=1.0e+004 *[0.00092482722817 0.00162222092086 5.96077986659246 0.02583283527855 1.02955757469816 0.00000320091835 0.00068597370239 0.00002039678512];%无约束优化最佳值
y1=1.0e+004 *[0.00091 0.0016 5.6 0.024 1 0.0000033 0.00068 0.00002];%最佳值估计
pre_val=1000;%初始容忍度设定
cur_val=1000;
x=zeros(1,8);%参数值
y=zeros(1,8);%优化值
y=y1;
for n=1:3 %综合考虑时间与效率设定循环搜索三周
x=y;
for k=1:8
for i=1:8
for j=1:24
x(k)=model(j)*gain(i);
cur_val=fun1(x);
if(cur_val end
end
x(k)=y(k);%很重要
end
G=fun1(y);
end
plot(Fr(y));axis([0 8000 0 0.55]);
toc;

main4.m

tic;
%代价函数概率分布
G=zeros(1,1000);
x1=1.0e+004 *[0.00206918469451 0.00219290683997 3.03044590029536 0.02254344914925 2.13513454683340 0.00000592719030 0.00080145017495 0.00000844666000];%无约束优化最佳值
for i=1:1000
for j=1:8
x(j)=x1(j)*0.95+x1(j)*(1.05-0.95)*rand(1);
end
G(i)=fun(x);
end
i=0:10:100;
hist(G,i);
toc;

main5.m
tic;
% x1=1.0e+004 *[0.00206918469451 0.00219290683997 3.03044590029536 0.02254344914925 2.13513454683340 0.00000592719030 0.00080145017495 0.00000844666000];%无约束优化最佳值
x1=1.0e+004 *[0.00130000000000 0.00160000000000 5.60000000000000 0.02700000000000 1.10000000000000 0.00000330000000 0.00068000000000 0.00002000000000];%离散取值时的最佳值
y=zeros(1,8);
y=x1;
rate=0.033;
m=0;
while (m<900)
rate=rate-0.001;
m=0;
for i=1:1000
for j=1:8
y(j)=x1(j)*(1-rate)+x1(j)*2*rate*rand(1);
end
k=fun(y);
if k<25
m=m+1;
end
end
end
rate
for i=1:1000
for j=1:8
y(j)=x1(j)*(1-rate)+x1(j)*2*rate*rand(1);
end
G(i)=fun(y);
end
i=0:10:100;
hist(G,i);
toc;