{"id":253,"date":"2022-05-04T21:12:00","date_gmt":"2022-05-04T19:12:00","guid":{"rendered":"https:\/\/bowfinger.de\/blog\/?p=253"},"modified":"2022-05-04T18:13:15","modified_gmt":"2022-05-04T16:13:15","slug":"quick-transfer-functions-with-sympy","status":"publish","type":"post","link":"https:\/\/bowfinger.de\/blog\/2022\/05\/quick-transfer-functions-with-sympy\/","title":{"rendered":"Quick Transfer Functions with SymPy"},"content":{"rendered":"\n<p><a href=\"https:\/\/bowfinger.de\/blog\/category\/python\/sympy\/\" data-type=\"category\" data-id=\"39\">SymPy<\/a> allows quick and (relatively \ud83d\ude1c) easy manipulation of algebraic expressions. But it can do so much more! Using the <code>sympy.physics.control<\/code> toolbox, you can very inspect linear time-invariant systems in a very comfortable way.<\/p>\n\n\n\n<p>I was looking at a paper (<a href=\"http:\/\/dx.doi.org\/10.1109\/9.233177\" data-type=\"URL\" data-id=\"dx.doi.org\/10.1109\/9.233177\">Janiszowski, 1993<\/a>) which features three transfer functions of different properties and wanted to look at their pole-zero diagrams as well as their Bode plots in real frequency. The transfer functions where the following:<\/p>\n\n\n<p class=\"ql-center-displayed-equation\" style=\"line-height: 133px;\"><span class=\"ql-right-eqno\"> (1) <\/span><span class=\"ql-left-eqno\"> &nbsp; <\/span><img loading=\"lazy\" decoding=\"async\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/ql-cache\/quicklatex.com-ba6d7846c20a3f897dfcda1f39572929_l3.png\" height=\"133\" width=\"202\" class=\"ql-img-displayed-equation quicklatex-auto-format\" alt=\" &#92;&#98;&#101;&#103;&#105;&#110;&#123;&#97;&#108;&#105;&#103;&#110;&#42;&#125; &#71;&#95;&#49;&#40;&#115;&#41;&#32;&#38;&#61;&#32;&#92;&#102;&#114;&#97;&#99;&#123;&#32;&#50;&#32;&#43;&#32;&#52;&#50;&#115;&#32;&#125;&#123;&#32;&#40;&#49;&#32;&#43;&#32;&#50;&#115;&#41;&#40;&#49;&#32;&#43;&#32;&#52;&#48;&#115;&#41;&#32;&#125;&#32;&#92;&#92; &#71;&#95;&#50;&#40;&#115;&#41;&#32;&#38;&#61;&#32;&#92;&#102;&#114;&#97;&#99;&#123;&#32;&#53;&#32;&#45;&#32;&#54;&#48;&#115;&#32;&#125;&#123;&#32;&#40;&#49;&#32;&#43;&#32;&#52;&#115;&#41;&#40;&#49;&#32;&#43;&#32;&#52;&#48;&#115;&#41;&#32;&#125;&#32;&#92;&#92; &#71;&#95;&#51;&#40;&#115;&#41;&#32;&#38;&#61;&#32;&#92;&#102;&#114;&#97;&#99;&#123;&#32;&#52;&#40;&#49;&#32;&#43;&#32;&#115;&#41;&#32;&#125;&#123;&#32;&#49;&#32;&#43;&#32;&#52;&#115;&#32;&#43;&#32;&#56;&#115;&#94;&#50;&#32;&#125; &#92;&#101;&#110;&#100;&#123;&#97;&#108;&#105;&#103;&#110;&#42;&#125; \" title=\"Rendered by QuickLaTeX.com\"\/><\/p>\n\n\n\n<p>Drawing the pole-zero diagrams from this representation is easy:<\/p>\n\n\n\n<ul class=\"is-style-default wp-block-list\"><li>Zeros: Set numerator to zero<\/li><li>Poles: Set denominator (or its individual terms) to zero<\/li><\/ul>\n\n\n\n<p>Drawing the Bode diagram however would&#8217;ve involved some basic programming. But neither is necessary.<\/p>\n\n\n\n<pre class=\"wp-block-code\"><code lang=\"python\" class=\"language-python\">import sympy as sy\nfrom sympy.physics.control.lti import TransferFunction as TF\nfrom sympy.physics.control import bode_plot, pole_zero_plot\n\ns = sy.symbols(\"s\", complex=True)\n\nG1 = (2 + 42*s)\/((1 + 2*s)*(1 + 40*s))\ntf = TF(*G.as_numer_denom(), s)\nbode_plot(tf)\npole_zero_plot(tf)<\/code><\/pre>\n\n\n\n<div class=\"wp-block-image\"><figure class=\"aligncenter size-full\"><img loading=\"lazy\" decoding=\"async\" width=\"397\" height=\"301\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output.png\" alt=\"\" class=\"wp-image-256\" srcset=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output.png 397w, https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-300x227.png 300w\" sizes=\"auto, (max-width: 397px) 100vw, 397px\" \/><figcaption>Bode plot of the transfer function G<sub>1<\/sub>(s)<\/figcaption><\/figure><\/div>\n\n\n\n<div class=\"wp-block-image\"><figure class=\"aligncenter size-full\"><img loading=\"lazy\" decoding=\"async\" width=\"400\" height=\"297\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-1.png\" alt=\"\" class=\"wp-image-260\" srcset=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-1.png 400w, https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-1-300x223.png 300w\" sizes=\"auto, (max-width: 400px) 100vw, 400px\" \/><figcaption>Pole-Zero diagram of the transfer function G<sub>1<\/sub>(s)<\/figcaption><\/figure><\/div>\n\n\n\n<p>Well, that was easy. <code>TransferFunction<\/code> (abbreviated with <code>TF<\/code> in my examples) requires numerator and denominator to be passed as separate arguments. <code>sympy_expr.as_numer_denom()<\/code> is convenient, as it returns a <code>(numerator, denominator)<\/code> tuple. Using the asterisk expands this.<\/p>\n\n\n\n<p>Now, what else can we do with this? Look at RLC resonators, for example:<\/p>\n\n\n\n<div class=\"wp-block-image\"><figure class=\"aligncenter size-large is-resized\"><img loading=\"lazy\" decoding=\"async\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/RLC.svg\" alt=\"\" class=\"wp-image-257\" width=\"299\" height=\"248\"\/><figcaption>A parallel resonator of resistance R, inductance L and capacitance C<\/figcaption><\/figure><\/div>\n\n\n\n<p>The reactance of the inductor and capacitor are given by<\/p>\n\n\n<p class=\"ql-center-displayed-equation\" style=\"line-height: 36px;\"><span class=\"ql-right-eqno\"> (2) <\/span><span class=\"ql-left-eqno\"> &nbsp; <\/span><img loading=\"lazy\" decoding=\"async\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/ql-cache\/quicklatex.com-109804ec1c87eeac03e96d7be60ae404_l3.png\" height=\"36\" width=\"209\" class=\"ql-img-displayed-equation quicklatex-auto-format\" alt=\" &#92;&#98;&#101;&#103;&#105;&#110;&#123;&#101;&#113;&#117;&#97;&#116;&#105;&#111;&#110;&#42;&#125; &#88;&#95;&#67;&#32;&#61;&#32;&#92;&#102;&#114;&#97;&#99;&#123;&#49;&#125;&#123;&#115;&#67;&#125;&#32;&#92;&#113;&#117;&#97;&#100;&#32;&#92;&#116;&#101;&#120;&#116;&#123;&#97;&#110;&#100;&#125;&#32;&#92;&#113;&#117;&#97;&#100;&#32;&#88;&#95;&#76;&#32;&#61;&#32;&#115;&#76; &#92;&#101;&#110;&#100;&#123;&#101;&#113;&#117;&#97;&#116;&#105;&#111;&#110;&#42;&#125; \" title=\"Rendered by QuickLaTeX.com\"\/><\/p>\n\n\n\n<p>where <em>s<\/em> is basically the complex frequency (we&#8217;re looking at this with Laplace eyes). Now, off to <a href=\"https:\/\/bowfinger.de\/blog\/category\/python\/sympy\/\" data-type=\"category\" data-id=\"39\">SymPy<\/a> we go:<\/p>\n\n\n\n<pre class=\"wp-block-code\"><code lang=\"python\" class=\"language-python\">R, C, L = sy.symbols(\"R, C, L\", real=True, positive=True)\nXc, Xl = 1\/(s*C), s*L\n\ndef par(a, b):\n    # Shorthand for a parallel circuit\n    return a*b\/(a+b)\n\nG = par(Xc, par(R, Xl)).ratsimp()\nG = G.subs({R: 1e3, L: 1e-6, C: 1e-6})\ntf = TF(*G.as_numer_denom(), s)\nbode_plot(tf, 5, 7)\npole_zero_plot(tf)<\/code><\/pre>\n\n\n\n<p>For visualisation we choose R = 1 k\u03a9, L = 1 \u00b5H and C = 1 \u00b5F.<\/p>\n\n\n\n<div class=\"wp-block-image\"><figure class=\"aligncenter size-full\"><img loading=\"lazy\" decoding=\"async\" width=\"391\" height=\"301\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-2.png\" alt=\"\" class=\"wp-image-261\" srcset=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-2.png 391w, https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-2-300x231.png 300w\" sizes=\"auto, (max-width: 391px) 100vw, 391px\" \/><figcaption>Bode plot of the parallel RLC circuit with R = 1 k\u03a9, L = 1 \u00b5H and C = 1 \u00b5F, and hence a resonance frequency of 1 MHz<\/figcaption><\/figure><\/div>\n\n\n\n<div class=\"wp-block-image\"><figure class=\"aligncenter size-full\"><img loading=\"lazy\" decoding=\"async\" width=\"400\" height=\"297\" src=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-3.png\" alt=\"\" class=\"wp-image-262\" srcset=\"https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-3.png 400w, https:\/\/bowfinger.de\/blog\/wp-content\/uploads\/2022\/05\/output-3-300x223.png 300w\" sizes=\"auto, (max-width: 400px) 100vw, 400px\" \/><figcaption>And its pole-zero plot just for the lolz.<\/figcaption><\/figure><\/div>\n\n\n\n<p>Why is this useful? Well, because it is quick and robust and once you have your transfer function typed out, you get the impulse and step responses for free:<\/p>\n\n\n\n<pre class=\"wp-block-code\"><code lang=\"python\" class=\"language-python\">import sympy as sy\nfrom sympy.physics.control import impulse_response_plot\nfrom sympy.physics.control import step_response_plot\n\ns = sy.symbols(\"s\", complex=True)\nG1 = (2+42*s)\/((1+2*s)*(1+40*s))\nG2 = (5-60*s)\/((1+4*s)*(1+40*s))\nG3 = 4*(1+s)\/(1+4*s+8*s**2)\n\nfor G in [G1, G2, G3]:\n    tf = TF(*G.as_numer_denom(), s)\n    bode_plot(tf)\n    pole_zero_plot(tf)\n    step_response_plot(tf, upper_limit=60)\n    impulse_response_plot(tf, upper_limit=60)<\/code><\/pre>\n\n\n\n<p class=\"has-small-font-size\">Note: SymPy&#8217;s adaptive plotting is not particularly good at plotting oscillations, so the impulse response of the RLC circuit above will look ugly.<\/p>\n","protected":false},"excerpt":{"rendered":"<p>SymPy allows quick and (relatively \ud83d\ude1c) easy manipulation of algebraic expressions. But it can do so much more! Using the sympy.physics.control toolbox, you can very inspect linear time-invariant systems in a very comfortable way. I was looking at a paper (Janiszowski, 1993) which features three transfer functions of different properties and wanted to look at&hellip;<\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"closed","ping_status":"open","sticky":false,"template":"","format":"standard","meta":{"inline_featured_image":false,"footnotes":""},"categories":[39],"tags":[],"class_list":["post-253","post","type-post","status-publish","format-standard","hentry","category-sympy"],"_links":{"self":[{"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/posts\/253","targetHints":{"allow":["GET"]}}],"collection":[{"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/comments?post=253"}],"version-history":[{"count":12,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/posts\/253\/revisions"}],"predecessor-version":[{"id":270,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/posts\/253\/revisions\/270"}],"wp:attachment":[{"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/media?parent=253"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/categories?post=253"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/bowfinger.de\/blog\/wp-json\/wp\/v2\/tags?post=253"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}