考尔型电路设计中求连分式商的Matlab实现

考尔型电路设计中求连分式商的Matlab实现

背景

搞电气工程的，研究生阶段有一门必修课《电网络理论》。当给出一个网络函数时F(s)时，需要用实际的电路去实现那个网络函数。一般用到的是考尔型电路去实现。下图所示就是考尔型电路。

图中电路参数与数学表达式的关系为：

而计算这些参数要用到“辗转相除”的方法，从而才能的到上面的连分式形式。
比如一个网络函数:
$$F(s)=\frac{s^4+10s^2+9}{s^3+2s}$$
下图是课本用辗转相除法的解体过程：

目的就是求出图中右边的那几个数，$\frac{9}{2s}$，$\frac{4}{11s}$，$\frac{121}{14s}$，$\frac{7}{11s}$。这几个参数分别对应上图电路中原件的参数。也就是通过数学方法把传递函数形式转化成连分式的形式，由连分式的形式就可以得到实际的电路了。

解决方案

此计算还是比较繁琐。于是我就想用Matlab编程实现。
程序如下：

function g = lianfenshi(N,D,H)
    clc;
    x=size(N);y=size(D);
    L=size(N)-size(D);
    if (L(1,2)>0)
         D((y(1,2)+1):x(1,2))=0;
    end
    if (L(1,2)<0)
        N((x(1,2)+1):y(1,2))=0;
    end
    clear x y L;

    syms s;
    L=size(D);
    P=L(1,2);
    j=1;q=1;

    if H==1%升幂
    for i=1:P
        if N(i)~=0
            A(j)=N(i)*s^(i-1);
            j=j+1;
        end
        if D(i)~=0
            B(q)=D(i)*s^(i-1);
            q=q+1;
        end

    end
    end

    if H==0%降幂

    for i=P:-1:1
        if N(i)~=0
            A(j)=N(i)*s^(i-1);
            j=j+1;
        end
        if D(i)~=0
            B(q)=D(i)*s^(i-1);
            q=q+1;
        end

    end
    end

    M(1,:)=A(1,:);   %第1次计算时的 被除数的多项式的因子构成的矩阵
    U(1,:)=B(1,:);%第1次计算时的 除数的多项式的因子构成的矩阵
    w(1)=M(1,1)/U(1,1);
    W(1)=w(1,1);%第1次计算时的 商
    W(1)

    x=size(M);y=size(U);
    L=size(M)-size(U);
    if (L(1,2)>0)
         U((y(1,2)+1):x(1,2))=0;
    end
    if (L(1,2)<0)
        M((x(1,2)+1):y(1,2))=0;
    end
    clear x y L;

    Z(1,:)=M(1,:)-U(1,:)*W(1);%第1次计算时的 余数多项式的因子构成的矩阵
    %z(1)=sum(M(1,:)')'-sum(U(1,:)')'*W(1);
    k=2;

    while sum(Z(k-1,:)')~=0
        M(k,:)=U(k-1,:);
        L=size(Z);
        clear U;
        U(k,:)=Z(k-1,2:L(1,2));
        clear L;
        w(k)=M(k,1)/U(k,1);
        W(k)=w(k);

        x=size(M);y=size(U);
        L=size(M)-size(U);
        if (L(1,2)>0)
         U(:,(y(1,2)+1):x(1,2))=0;
        end
        if (L(1,2)<0)
        M(:,(x(1,2)+1):y(1,2))=0;
        end
        clear x y L;

        Z(k,:)=M(k,:)-U(k,:)*W(k);
        W(k)
        k=k+1;
    end
end

验证

当在命令窗口输入如下代码时：

>> D=[0 2 0 1];
>> N=[9 0 10 0 1];
>> lianfenshi(N,D,1);

计算截图如下：