From ab5b5fb8600736dc0b32ccd62b44a95b296164bf Mon Sep 17 00:00:00 2001 From: Christian Kral Date: Sat, 23 Jan 2021 22:39:12 +0100 Subject: [PATCH] Fix bug in calculation of complex transfer function --- .../ComplexMath/TransferFunction.mo | 14 ++++++++----- .../Examples/ShowTransferFunction.mo | 20 ++++++++++--------- 2 files changed, 20 insertions(+), 14 deletions(-) diff --git a/Modelica/ComplexBlocks/ComplexMath/TransferFunction.mo b/Modelica/ComplexBlocks/ComplexMath/TransferFunction.mo index cc9d344007..435b506f2a 100644 --- a/Modelica/ComplexBlocks/ComplexMath/TransferFunction.mo +++ b/Modelica/ComplexBlocks/ComplexMath/TransferFunction.mo @@ -3,9 +3,9 @@ block TransferFunction "Complex Transfer Function" extends Modelica.ComplexBlocks.Interfaces.ComplexSISO; import Modelica.ComplexMath.j; parameter Real b[:]={1} - "Numerator coefficients of transfer function (e.g., 2*s+3 is specified as {2,3})"; + "Numerator coefficients of transfer function (e.g., 2*(jw)+3 is specified as {2,3})"; parameter Real a[:]={1} - "Denominator coefficients of transfer function (e.g., 5*s+6 is specified as {5,6})"; + "Denominator coefficients of transfer function (e.g., 5*(jw)+6 is specified as {5,6})"; Modelica.Blocks.Interfaces.RealInput w(unit="rad/s") "Angular frequency input" annotation (Placement(transformation( extent={{-20,-20},{20,20}}, rotation=90, @@ -19,15 +19,15 @@ protected input Integer k; output Complex x; protected - Integer m=mod(k,4); + Integer m=mod(k,4); algorithm x:=if m==0 then Complex(1) elseif m==1 then j elseif m==2 then Complex(-1) else -j; annotation(Inline=true); end powerOfJ; equation // Avoid computing power of Complex numbers - since it fails for w==0 - bw = {b[i]*(w)^(i-1)*powerOfJ(i-1) for i in 1:size(b,1)}; - aw = {a[i]*(w)^(i-1)*powerOfJ(i-1) for i in 1:size(a,1)}; + bw = {b[i]*(w)^(size(b,1)-i)*powerOfJ(size(b,1)-i) for i in 1:size(b,1)}; + aw = {a[i]*(w)^(size(a,1)-i)*powerOfJ(size(a,1)-i) for i in 1:size(a,1)}; bSum = Complex(sum(bw.re), sum(bw.im)); aSum = Complex(sum(aw.re), sum(aw.im)); y = u*bSum/aSum; @@ -50,5 +50,9 @@ The complex input u is multiplied by the complex transfer function (depending on y(jw) = ------------------------------------------------- * u(jw) a[1]*(jw)^[na-1] + a[2]*(jw)^[na-2] + ... + a[na] +", revisions=" + ")); end TransferFunction; diff --git a/Modelica/ComplexBlocks/Examples/ShowTransferFunction.mo b/Modelica/ComplexBlocks/Examples/ShowTransferFunction.mo index a0425b6642..45a5168eda 100644 --- a/Modelica/ComplexBlocks/Examples/ShowTransferFunction.mo +++ b/Modelica/ComplexBlocks/Examples/ShowTransferFunction.mo @@ -1,9 +1,11 @@ within Modelica.ComplexBlocks.Examples; -model ShowTransferFunction "Test Complex Transfer Function Block" +model ShowTransferFunction "Test complex transfer function block" extends Modelica.Icons.Example; - parameter Real d=1/sqrt(2) "Damping coefficient"; - parameter Real b[:]={1} "Numerator polynomial coefficients of the transfer function"; - parameter Real a[:]={1,2*d,1} "Denominator polynomial coefficients of the transfer function"; + parameter Real d=0.01 "Damping coefficient in kg/s (not the damping ratio)"; + parameter Real m=0.2 "Mass in kg"; + parameter Real c=0.1 "Stiffness in N/m"; + parameter Real b[:]={-m} "Numerator polynomial coefficients {-m} of the transfer function"; + parameter Real a[:]={m,d,c} "Denominator polynomial coefficients {m,d,c} of the transfer function"; parameter Real wMin=0.01 "Lower bound for frequency sweep"; parameter Real wMax=100 "Upper bound for frequency sweep"; Real lg_w=log10(logFrequencySweep.y) "Logarithm of frequency"; @@ -28,13 +30,13 @@ equation connect(transferFunction.y, complexToPolar.u) annotation (Line(points={{-19,0},{-2,0}}, color={85,170,255})); annotation ( - experiment(StopTime=1, Interval=0.001), Documentation(info= - " -

This example shows the response of a PT2 defined by its transfer function

+ experiment(StopTime=1, Interval=0.001), Documentation(info=" +

This example shows the response of a PT2 (mechanical spring-mass-damper system with +an acceleration acting on the mass) defined by its transfer function

-            1
+              -m
 H(jw)=-------------------
-      1 + 2 d jw + (jw)^2
+      m*(jw)^2 + d*jw + c

Frequency performs a logarithmic ramp from 0.01 to 100 s^-1.