index e937db3..286120c 100644 (file)
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2010 - INRIA - Serge Steer
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+
+// <-- ENGLISH IMPOSED -->
+
//gain margin
//-----------
+function ok=check_gmargin(h,fref,gref)
+  eps=1e-10
+  [g,f]=g_margin(h)
+  if fref==0 then
+     ok=abs(f - fref)<eps
+  else
+    ok=abs(f - fref)/abs(fref)<eps
+  end
+  if gref==0 then
+     ok=ok&abs(g - gref)<eps
+  else
+    ok=ok&abs(g - gref)/abs(gref)<eps
+  end
+endfunction
+
//discrete time case
-h = syslin(0.1,0.04798*%z+0.0464,%z^2-1.81*%z+0.9048);//pk
-[g ,f]=g_margin(h);
-if norm(g - 6.2423793565916954)>1d-5  then pause,end
-if norm(f - 0.86539118158628991)>1d-5  then pause,end
+z=poly(0,'z');
+//the references solutions are computed using the following Maple instructions
+//Digits:=50:
+//z:=exp(I*w):
+//assume(w, 'real', w > 0)
+//h:= rational fraction in z
+//M := solve(Im(h) = 0, w)
+//evalf(M/(2*Pi))
+//evalf(eval(20*log10(1/abs(h)), w = M[k]))
+
+h = syslin(0.1,((29/625)+(2399/50000)*z)/((1131/1250)+(-181/100)*z+z^2));
+f_ref=atan(sqrt(9003479136639)/4963519)/(0.1*2*%pi);
+g_ref=6.242379356591869534116;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end
+

-h=syslin(0.1,0.086-0.161*%z+0.078*%z^2,-0.67+2.286*%z-2.61*%z^2+%z^3);//ok
-[g ,f]=g_margin(h);
-if norm(g- 12.3470513198133762)>1d-5  then pause,end
-if norm(f-0.26109869168929700)>1d-5  then pause,end
+h=syslin(0.1,((43/500)+(-161/1000)*z+(39/500)*z^2)/((-67/100)+(1143/500)*z+(-261/100)*z^2+z^3));
+f_ref=atan((1/215)*sqrt(1136805-8330*sqrt(5970))/(833/43+(1/215)*sqrt(5970)))/(0.1*2*%pi);
+g_ref=12.3470513198103944037870;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end

-h=syslin(0.1,3*(0.086-0.161*%z+0.078*%z^2),-0.67+2.286*%z-2.61*%z^2+%z^3);//ok
-[g ,f]=g_margin(h);
-if norm(g-2.80462622542011131)>1d-5  then pause,end
-if norm(f-0.26109869168929828)>1d-5  then pause,end
+h=syslin(0.1,((129/500)+(-483/1000)*z+(117/500)*z^2)/((-67/100)+(1143/500)*z+(-261/100)*z^2+z^3));
+f_ref=atan((1/215)*sqrt(1136805-8330*sqrt(5970))/(833/43+(1/215)*sqrt(5970)))/(0.1*2*%pi);
+g_ref=2.8046262254171456578864840;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end

+h=syslin(1,((21/500)+(-3933/100000)*z+(-15407/100000)*z^2+(9259/50000)*z^3+(6939/(10^12))*z^4)/((-21/500)+(11/125)*z+(57/500)*z^2+(-9/25)*z^3+(1/5)*z^4));
+f_ref=1/2;
+g_ref=6.48227782514616029706190;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end
+
//continuuous time case
-h=syslin('c',-1+%s,3+2*%s+%s^2);//ok
-[g,f]=g_margin(h);
-if norm(g-20*log(3)/log(10))>1d-5  then pause,end
-if norm(f)>1d-5  then pause,end
+s=poly(0,'s');
+//the reference solutions are computed using  the following Maple instructions
+//s:=I*w:
+//assume(w, 'real', w > 0)
+//h:= rational fraction in z
+//M := solve(Im(h) = 0, w)
+//evalf(M/(2*Pi))
+//evalf(eval(20*log10(1/abs(h)), w = M[k]))
+
+h=syslin('c',(-1+s)/(3+2*s+s^2));
+f_ref=0;
+g_ref=9.5424250943932487459005580;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end
+
+
+h = syslin('c',((1/5)+(4/5)*s+(0/1)*s^2+(3/10)*s^3)/((409/10000)+(1827/10000)*s+(5129/4000)*s^2+(31909/10000)*s^3+(64/25)*s^4 +s^5));
+f_ref=0.0714555258202006740373134;
+g_ref=-4.91687406933815400242335;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end
+

+h=syslin('c',485000/(10000*s+200*s^2+s^3));
+f_ref=100/(2*%pi);
+g_ref=12.305765141234350772862319;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end

- h = syslin('c',0.2+0.8*%s+0.3*%s^3,0.0409+0.1827*%s+1.28225*%s^2+3.1909*%s^3+2.56*%s^4+%s^5);
-[g,f]=g_margin(h);
-if norm(g + 4.91687406933815385)>1d-5  then pause,end
-if norm(f-0.07145552582020068)>1d-5  then pause,end
+
+h = syslin('c',1/(s+2*s^2+s^3));
+f_ref=1/(2*%pi);
+g_ref=6.0205999132796239042747779;
+if ~check_gmargin(h,f_ref,g_ref) then pause,end

//phase margin
//-----------
+function ok=check_pmargin(h,f_ref,p_ref)
+  eps=1e-9
+  [p,f]=p_margin(h)
+  if f_ref==0 then
+     ok=abs(f - f_ref)<eps
+  else
+    ok=abs(f - f_ref)/abs(f_ref)<eps
+  end
+  if p_ref==0 then
+     ok=ok&abs(p - p_ref)<eps
+  else
+    ok=ok&abs(p - p_ref)/abs(p_ref)<eps
+  end
+endfunction
+
//discrete time case
-h = syslin(0.1,0.04798*%z+0.0464,%z^2-1.81*%z+0.9048);//ok
-[p ,f]=p_margin(h);
-if norm(p-13.5711556361199754)>1d-5  then pause,end
-if norm(f-0.69301660031536083)>1d-5  then pause,end
+//the reference solutions are computed using  the following Maple instructions
+//z:=exp(I*w):
+//assume(w, 'real', w > 0)
+//P := solve(abs(h) = 1, w)
+//evalf(-(eval(180-180*argument(h)/Pi, w = Re(P[k]))))+360

+h = syslin(0.1,((29/625)+(2399/50000)*z)/((1131/1250)+(-181/100)*z+z^2));
+f_ref=0.693016600315284442350578876;
+p_ref=13.57115563612946355428439468;
+if  ~check_pmargin(h,f_ref,p_ref) then pause,end

-h=syslin(0.1,0.086-0.161*%z+0.078*%z^2,-0.67+2.286*%z-2.61*%z^2+%z^3);//ok
+
+h=syslin(0.1,((43/500)+(-161/1000)*z+(39/500)*z^2)/((-67/100)+(1143/500)*z+(-261/100)*z^2+z^3));
[p ,f]=p_margin(h);
if p<>%inf  then pause,end
if f<>[]  then pause,end

-h=syslin(0.1,3*(0.086-0.161*%z+0.078*%z^2),-0.67+2.286*%z-2.61*%z^2+%z^3);//ok
-[p ,f]=p_margin(h);
-if norm(p-52.9496741629204308 )>1d-5  then pause,end
-if norm(f-0.21233648894721618)>1d-5  then pause,end
-
+h=syslin(0.1,3*(0.086-0.161*%z+0.078*%z^2),-0.67+2.286*%z-2.61*%z^2+%z^3);
+h=syslin(0.1,((129/500)+(-483/1000)*z+(117/500)*z^2)/((-67/100)+(1143/500)*z+(-261/100)*z^2+z^3));
+f_ref=0.212336488950669705771059018;
+p_ref=52.94967415965772478856630911;
+if  ~check_pmargin(h,f_ref,p_ref) then pause,end

//continuous case
-h=syslin('c',1.1+2.4*%s+0.7*%s^2,3+2*%s+%s^2);
-[p,f]=p_margin(h);
-if norm(p-(-148.547076202317413))>1d-5  then pause,end
-if norm(f-0.18945852512301298)>1d-5  then pause,end
-
- h = syslin('c',0.2+0.8*%s+0.3*%s^3,0.0409+0.1827*%s+1.28225*%s^2+3.1909*%s^3+2.56*%s^4+%s^5);
-[p,f]=p_margin(h);
-if norm(p + 13.1128497150069734 )>1d-5  then pause,end
-if norm(f-0.09144216563554157  )>1d-5  then pause,end
+//the reference solutions are computed using  the following Maple instructions
+//z:=I*w:
+//assume(w, 'real', w > 0)
+//P := solve(abs(h) = 1, w)
+//evalf(-(eval(180-180*argument(h)/Pi, w = Re(P[k]))))+360
+
+h=syslin('c',((11/10)+(12/5)*s+(7/10)*s^2)/(3+2*s+s^2));
+f_ref=(1/51)*sqrt(15861-204*sqrt(3562))/(2*%pi);
+p_ref=-148.547076202317410601324666;
+if  ~check_pmargin(h,f_ref,p_ref) then pause,end
+
+
+h = syslin('c',((1/5)+(4/5)*s+(3/10)*s^3)/((409/10000)+(1827/10000)*s+(5129/4000)*s^2+(31909/10000)*s^3+(64/25)*s^4+s^5));
+f_ref=0.09144216563554157543991;
+p_ref=-13.1128497150069802772313;
+if  ~check_pmargin(h,f_ref,p_ref) then pause,end