{"id":684,"date":"2013-04-26T06:14:21","date_gmt":"2013-04-26T11:14:21","guid":{"rendered":"https:\/\/blogs.mathworks.com\/loren\/?p=684"},"modified":"2017-01-06T11:12:30","modified_gmt":"2017-01-06T16:12:30","slug":"using-symbolic-math-toolbox-to-compute-area-moments-of-inertia","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/loren\/2013\/04\/26\/using-symbolic-math-toolbox-to-compute-area-moments-of-inertia\/","title":{"rendered":"Using Symbolic Math Toolbox to Compute Area Moments of Inertia"},"content":{"rendered":"<div class=\"content\">\n<p><!--introduction-->Once more I am pleased to introduce guest blogger Kai Gehrs. Kai has been a Software Engineer at MathWorks for the past five years working on the <a href=\"https:\/\/www.mathworks.com\/help\/toolbox\/symbolic\">Symbolic Math Toolbox<\/a>. He has a background in mathematics and computer science. He already contributed to my BLOG in the past writing about <a href=\"https:\/\/blogs.mathworks.com\/loren\/2012\/07\/27\/using-symbolic-equations-and-symbolic-functions-in-matlab\">Using Symbolic Equations And Symbolic Functions In MATLAB<\/a> as well as on approaches for <a href=\"https:\/\/blogs.mathworks.com\/loren\/2011\/10\/25\/simplifying-symbolic-results\">Simplifying Symbolic Results<\/a>.<\/p>\n<p><!--\/introduction--><\/p>\n<h3>Contents<\/h3>\n<div>\n<ul>\n<li><a href=\"#05e69b68-4591-4e86-a7a3-4366fbda2fad\">In a Nutshell: What Is This Article About?<\/a><\/li>\n<li><a href=\"#510e5d8c-2254-4f99-9565-cc9fb1f1338c\">Basic Example: Cross Section of an Elliptical Tube<\/a><\/li>\n<li><a href=\"#319f5ba4-64a6-4e2e-b648-831a0818d6ff\">Area Moment Of Inertia<\/a><\/li>\n<li><a href=\"#3b715d35-70da-4013-bb0b-470df2df3da5\">Math Behind This Example<\/a><\/li>\n<li><a href=\"#19841dd9-5c50-4448-ba67-76dd636545d3\">Symbolic Integration<\/a><\/li>\n<li><a href=\"#6c55693e-8593-48c1-90be-48129f5b3b0e\">Advanced Example: Cross Section of an Elliptical Tube Defined by Symbolic Parameters<\/a><\/li>\n<li><a href=\"#2415b2ac-5410-444a-8c2b-f15693021763\">How to Create MuPAD Graphics Shown in This Article<\/a><\/li>\n<li><a href=\"#ebec3076-fed9-49e8-967d-2c5a6a31fea0\">Have You Tried Symbolic Math Toolbox?<\/a><\/li>\n<\/ul>\n<\/div>\n<h4>In a Nutshell: What Is This Article About?<a name=\"05e69b68-4591-4e86-a7a3-4366fbda2fad\"><\/a><\/h4>\n<p>If you are interested in using MATLAB and the Symbolic Math Toolbox in teaching some basics in mechanical engineering, this might be of interest to you.<\/p>\n<p>Computing <i>area moments of inertia<\/i> is an important task in mechanics. For example, area moments of inertia play a critical role in stress, dynamic, and stability analysis of structures. In this article, we use capabilities of the <a href=\"www.mathworks.com\/products\/symbolic\">Symbolic Math Toolbox<\/a> to compute area moments for cross sections of elliptical tubes.<\/p>\n<p>We start with a basic case involving only numeric parameters, and then make the computations more general by introducing symbolic parameters.<\/p>\n<p>All plots used in this article have been created in the <a href=\"https:\/\/www.mathworks.com\/discovery\/mupad.html\">MuPAD Notebook app<\/a>. You can find the source code for these plots at the end of the article.<\/p>\n<h4>Basic Example: Cross Section of an Elliptical Tube<a name=\"510e5d8c-2254-4f99-9565-cc9fb1f1338c\"><\/a><\/h4>\n<p>The following picture shows the cross section of an elliptical tube.<\/p>\n<p><img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2013\/lateralCut.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/><\/p>\n<p>The following ellipses define outer and inner contours of the section:<\/p>\n<div>\n<ul>\n<li>$y = \\pm 2 \\, \\sqrt{1 - \\frac{x^2}{9}}$ describe the <i>outer<\/i> contour line for $x \\in [-3,3]$.<\/li>\n<li>$y = \\pm \\sqrt{1 - \\frac{x^2}{4}}$ describe the <i>inner<\/i> contour line for $x \\in [-2,2]$.<\/li>\n<\/ul>\n<\/div>\n<p>We will compute the area moment of inertia of this section with respect to the $x$-axis.<\/p>\n<h4>Area Moment Of Inertia<a name=\"319f5ba4-64a6-4e2e-b648-831a0818d6ff\"><\/a><\/h4>\n<p>The moment of inertia of an area $A$ with respect to the $x$-axis is defined in terms of the double integral<\/p>\n<p>$$I_x = {\\int \\int}_A y^2\\, \\mathrm{d}A.$$<\/p>\n<p>You can find this definition on <a href=\"http:\/\/en.wikipedia.org\/wiki\/Second_moment_of_area\">Wikipedia<\/a>.<\/p>\n<h4>Math Behind This Example<a name=\"3b715d35-70da-4013-bb0b-470df2df3da5\"><\/a><\/h4>\n<p>In our example, the area $A$ is the hatched area shown in the previous plot. Taking advantage of the symmetry of the section with respect to both the $x$- and $y$-axes, we can restrict our computations to the first quadrant of the $x-y$-plane.<\/p>\n<p><img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2013\/UpperQuarterLateralCut.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/><\/p>\n<p>To compute the necessary double integral over the hatched area, we divide this area into two separate areas $A1$ and $A2$.<\/p>\n<p><img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2013\/IntegrationAreas.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/><\/p>\n<p>To compute the final area moment of inertia about the x-axis, we multiply the sum of the double integrals over $A_1$ and $A_2$ by four.<\/p>\n<p>$$I_x = {\\int\\int}_A y^2 \\, \\mathrm dA = 4 \\cdot \\Big({\\int\\int}_{A_1}<br \/>\ny^2 \\, \\mathrm dA_1 + {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2\\Big).$$<\/p>\n<p>Now, we only need to compute these two double integrals. The mathematical theory behind multi-dimensional integration tells us that we can rewrite each of these double integrals in terms of two integrals:<\/p>\n<p>$$I_1 = {\\int\\int}_{A_1} y^2 \\, \\mathrm dA_1 = \\int_0^2 \\int_{\\sqrt{1 -<br \/>\n\\frac{x^2}{4}}}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}} y^2 \\, \\mathrm dy \\,<br \/>\n\\mathrm dx,$$<\/p>\n<p>$$I_2 = {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2 = \\int_2^3 \\int_{0}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}}<br \/>\ny^2\\, \\mathrm dy\\, \\mathrm dx.$$<\/p>\n<p>The <a title=\"https:\/\/www.mathworks.com\/help\/symbolic\/int.html (link no longer works)\"><tt>int<\/tt><\/a> function from the Symbolic Math Toolbox can compute these integrals.<\/p>\n<h4>Symbolic Integration<a name=\"19841dd9-5c50-4448-ba67-76dd636545d3\"><\/a><\/h4>\n<p>We define $x$ and $y$ as symbolic variables:<\/p>\n<pre class=\"codeinput\">syms <span class=\"string\">x<\/span> <span class=\"string\">y<\/span>;\r\n<\/pre>\n<p>We start computing the inner integral $\\int_{\\sqrt{1 - \\frac{x^2}{4}}}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}} y^2 \\, \\mathrm dy$ of $I_1$:<\/p>\n<pre class=\"codeinput\">innerI1 = int(y^2, y, sqrt(1-x^2\/4), 2*sqrt(1-x^2\/9))\r\n<\/pre>\n<pre class=\"codeoutput\">innerI1 =\r\n(8*(9 - x^2)^(3\/2))\/81 - (4 - x^2)^(3\/2)\/24\r\n<\/pre>\n<p>Next, we integrate <tt>innerI1<\/tt> with respect to $x$ from $0$ to $2$. This provides $I_1$:<\/p>\n<pre class=\"codeinput\">I1 = int(innerI1, x, 0, 2)\r\n<\/pre>\n<pre class=\"codeoutput\">I1 =\r\n3*asin(2\/3) - pi\/8 + (74*5^(1\/2))\/81\r\n<\/pre>\n<p>We use the same strategy to compute $I_2$. First, we compute the inner integral $\\int_{0}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}} y^2\\, \\mathrm dy$. Then, we integrate the resulting expression with respect to $x$ from $2$ to $3$:<\/p>\n<pre class=\"codeinput\">I2 = int(int(y^2, y, 0, 2*sqrt(1-x^2\/9)), x, 2, 3);\r\n<\/pre>\n<p>Hence, the area moment of inertia of the elliptical tube with respect to the $x$-axis is<\/p>\n<pre class=\"codeinput\">Ix = 4 * (I1 + I2)\r\n<\/pre>\n<pre class=\"codeoutput\">Ix =\r\n(11*pi)\/2\r\n<\/pre>\n<p>We can approximate the result numerically by using the <a href=\"https:\/\/www.mathworks.com\/help\/toolbox\/symbolic\/vpa.html\"><tt>vpa<\/tt><\/a> function. For example, we approximate the result using 5 significant (nonzero) digits:<\/p>\n<pre class=\"codeinput\">vpa(Ix, 5)\r\n<\/pre>\n<pre class=\"codeoutput\">ans =\r\n17.279\r\n<\/pre>\n<h4>Advanced Example: Cross Section of an Elliptical Tube Defined by Symbolic Parameters<a name=\"6c55693e-8593-48c1-90be-48129f5b3b0e\"><\/a><\/h4>\n<p>We can use a numerical integrator, such as MATLAB's <tt>integral2<\/tt>, to compute the area moment of inertia in the previous example. A numerical integrator might return slightly less accurate results, but other than that there is not much benefit from using symbolic integration there.<\/p>\n<p>But what if we consider the cross section of a more general elliptical tube whose shape is defined by arbitrary symbolic parameters?<\/p>\n<p>For example, we want to compute the area moment of inertia of this elliptical tube.<\/p>\n<p><img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2013\/LateralCutMoreGeneral.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/><\/p>\n<p>Its contour lines are still ellipses, but they are defined by the symbolic parameters $r_1$, $r_2$, $R_1$ and $R_2$:<\/p>\n<div>\n<ul>\n<li>$y = \\pm R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}$ describe the outer contour line for $x \\in [-R_1,R_1]$.<\/li>\n<li>$y = \\pm r_2 \\, \\sqrt{1 - \\frac{x^2}{r_1^2}}$ describe the inner contour line for $x \\in [-r_1,r_1]$.<\/li>\n<\/ul>\n<\/div>\n<p>Now we need to compute these integrals:<\/p>\n<p>$$I_1 = {\\int\\int}_{A_1} y^2 \\, \\mathrm dA_1 = \\int_0^{r_1} \\int_{r_2 \\, \\sqrt{1 -<br \/>\n\\frac{x^2}{r_1^2}}}^{R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}} y^2 \\, \\mathrm dy \\,<br \/>\n\\mathrm dx,$$<\/p>\n<p>$$I_2 = {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2 = \\int_{r_1}^{R_1}<br \/>\n\\int_{0}^{R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}} y^2\\, \\mathrm dy\\, \\mathrm dx.$$<\/p>\n<p>Although the situation is more complicated now, we can still use the same strategy as above, using symbolic variables instead of numbers. Nevertheless, we need to add one more step: specify relationships between symbolic parameters. This is because the variables $r_1$, $r_2$, $R_1$, and $R_2$ can be any complex number, unless we explicitly restrict their values. For example, the system does not know if these variables are positive or negative, if $r_1 &lt; R_1$ or otherwise. In this example, we want to specify that $r_1&gt;0$, $r_2&gt;0$, $R_1&gt;r_1$, and $R_2&gt;r_2$. In the Symbolic Math Toolbox, such restrictions are called <i>assumptions on variables<\/i>. They can be set by using the <tt>assume<\/tt> and <tt>assumeAlso<\/tt> functions:<\/p>\n<pre class=\"codeinput\">syms <span class=\"string\">x<\/span> <span class=\"string\">y<\/span> <span class=\"string\">r1<\/span> <span class=\"string\">r2<\/span> <span class=\"string\">R1<\/span> <span class=\"string\">R2<\/span>;\r\nassume(r1 &gt; 0);\r\nassume(r2 &gt; 0);\r\nassumeAlso(r1 &lt; R1);\r\nassumeAlso(r2 &lt; R2);\r\nI1 = int(int(y^2, y, r2*sqrt(1-x^2\/r1^2), R2*sqrt(1-x^2\/R1^2)), x, 0, r1);\r\nI2 = int(int(y^2, y, 0, R2*sqrt(1-x^2\/R1^2)), x, r1, R1);\r\nIxgeneral = 4 * (I1 + I2);\r\npretty(Ixgeneral)\r\n<\/pre>\n<pre class=\"codeoutput\"> \r\n        \/          4     \/ r1 \\ \\         \/                     4     \/ r1 \\ \\ \r\n        |      3 R1  asin| -- | |         |             4   3 R1  asin| -- | | \r\n      3 |                \\ R1 \/ |       3 |      3 pi R1              \\ R1 \/ | \r\n  4 R2  | #1 + ---------------- |   4 R2  | #1 - -------- + ---------------- | \r\n        \\             8         \/         \\         16             8         \/ \r\n  ------------------------------- - ------------------------------------------ \r\n                   3                                      3 \r\n               3 R1                                   3 R1 \r\n   \r\n                3 \r\n        pi r1 r2 \r\n      - --------- \r\n            4 \r\n   \r\n  where \r\n   \r\n           \/     2        3 \\ \r\n           | 5 R1  r1   r1  |    2     2 1\/2 \r\n     #1 == | -------- - --- | (R1  - r1 ) \r\n           \\    8        4  \/\r\n<\/pre>\n<p>Now, let's substitute our symbolic parameters with the same values that we had in the first example. Do we get the same result? Let's double-check: when substituting $r_1=2,r_2,=1,R_1=3,R_2=2$ in <tt>Ixgeneral<\/tt>, do we get the same result for <tt>Ix<\/tt> as obtained for the section initially considered?<\/p>\n<p>As expected, the general solution <tt>Ixgeneral<\/tt> reduces to the specific solution <tt>Ix<\/tt> when we substitute our original dimensions:<\/p>\n<pre class=\"codeinput\">vpa(subs(Ixgeneral, [r1, R1, r2, R2],[2, 3, 1, 2]), 5)\r\n<\/pre>\n<pre class=\"codeoutput\">ans =\r\n17.279\r\n<\/pre>\n<p>Note that it is essential to set the assumptions on $r_1$, $r_2$, $R_1$ and $R_2$. By default, the Symbolic Math Toolbox assumes that all variables including symbolic parameters are arbitrary complex numbers. This makes computing integrals much more complicated and time-consuming. For example, try computing <tt>Ixgeneral<\/tt> without setting the assumptions on $r_1$, $r_2$, $R_1$ and $R_2$. In general, it can be impossible to even get the result without any additional assumptions.<\/p>\n<h4>How to Create MuPAD Graphics Shown in This Article<a name=\"2415b2ac-5410-444a-8c2b-f15693021763\"><\/a><\/h4>\n<p>To generate the graphics shown in this article, use <a href=\"https:\/\/www.mathworks.com\/discovery\/mupad.html\">MuPAD Notebook app<\/a>. Note that the code in this section does not run in the MATLAB Command Window.<\/p>\n<p>To open the MuPAD Notebook app:<\/p>\n<div>\n<ol>\n<li>On the MATLAB Toolstrip, click the Apps tab.<\/li>\n<li>On the Apps tab, click the down arrow at the end of the Apps section.<\/li>\n<li>Under Math, Statistics and Optimization, click the MuPAD Notebook button.<\/li>\n<\/ol>\n<\/div>\n<p>Alternatively, type <tt>mupad<\/tt> in the MATLAB Command Window.<\/p>\n<p>Paste the following commands to a MuPAD notebook to obtain the four graphics used above:<\/p>\n<pre>Ellipse:= (x,a,b) -&gt; sqrt(b^2*(1-x^2\/a^2)):<\/pre>\n<pre>f1 := plot::Function2d(Ellipse(x,2,1), x = -3..3):\r\nf2 := plot::Function2d(Ellipse(x,3,2), x = -3..3):\r\nf3 := plot::Function2d(Ellipse(x,3,2), x = -3..-2):\r\nf4 := plot::Function2d(Ellipse(x,3,2), x = 2..3):<\/pre>\n<pre>f5 := plot::Function2d(-Ellipse(x,2,1), x = -3..3):\r\nf6 := plot::Function2d(-Ellipse(x,3,2), x = -3..3):\r\nf7 := plot::Function2d(-Ellipse(x,3,2), x = -3..-2):\r\nf8 := plot::Function2d(-Ellipse(x,3,2), x = 2..3):<\/pre>\n<pre>f9 := plot::Function2d(Ellipse(x,2,1), x = 0..2, VerticalAsymptotesVisible = FALSE):\r\nf10 := plot::Function2d(Ellipse(x,3,2), x = 0..3, VerticalAsymptotesVisible = FALSE):\r\nf11 := plot::Function2d(Ellipse(x,3,2), x = 2..3, VerticalAsymptotesVisible = FALSE):<\/pre>\n<pre>plot(plot::Hatch(f1,f2),plot::Hatch(f3),plot::Hatch(f4),f1,f2,\r\n     plot::Hatch(f5,f6),plot::Hatch(f7),plot::Hatch(f8),f5,f6,\r\n     Scaling = Constrained,\r\n     GridVisible = TRUE,\r\n     VerticalAsymptotesVisible = FALSE,\r\n     Height = 120, Width = 200,\r\n     Header = \"Cross section of elliptical tube (x-y-plane)\"):<\/pre>\n<pre>plot(plot::Hatch(f9,f10),plot::Hatch(f11),f9,f10,\r\n     Scaling = Constrained,\r\n     GridVisible = TRUE,\r\n     Height = 120, Width = 200,\r\n     Header = \"Restriction to first quadrant (x-y-plane)\"):<\/pre>\n<pre>plot(plot::Hatch(f9,f10,Color=RGB::Magenta),plot::Hatch(f11,Color = RGB::Green),f9,f10,\r\n     Scaling=Constrained,\r\n     GridVisible = TRUE,\r\n     plot::Text2d(\"A1\",[1.0,1.25]),\r\n     plot::Text2d(\"A2\",[2.5,0.25]),\r\n     Height = 120, Width = 200,\r\n     Header = \"Integration areas\"):<\/pre>\n<pre>plot(plot::Hatch(f1,f2),plot::Hatch(f3),plot::Hatch(f4),f1,f2,\r\n     plot::Hatch(f5,f6),plot::Hatch(f7),plot::Hatch(f8),f5,f6,\r\n     plot::Text2d(\"r1\",[1.8,0.05]),\r\n     plot::Text2d(\"r2\",[0.05,0.85]),\r\n     plot::Text2d(\"R1\",[2.74,0.05]),\r\n     plot::Text2d(\"R2\",[0.05,1.84]),\r\n     Scaling = Constrained,\r\n     GridVisible = TRUE,\r\n     VerticalAsymptotesVisible = FALSE,\r\n     XTicksLabelsVisible = FALSE,\r\n     YTicksLabelsVisible = FALSE,\r\n     ViewingBox = [-3.5..3.5,-2.5..2.5],\r\n     Height = 120, Width = 200,\r\n     Header = \"Cross section of elliptical tube (more general situation)\"):<\/pre>\n<h4>Have You Tried Symbolic Math Toolbox?<a name=\"ebec3076-fed9-49e8-967d-2c5a6a31fea0\"><\/a><\/h4>\n<p>Have you used the Symbolic Math Toolbox in computations related to mechanics? Let me know <a href=\"https:\/\/blogs.mathworks.com\/loren\/?p=684#respond\">here<\/a>.<\/p>\n<p><script>\/\/ <![CDATA[\nfunction grabCode_75b0bee11e3e4c56aafb679a67295a76() {\n        \/\/ Remember the title so we can use it in the new page\n        title = document.title;\n\n        \/\/ Break up these strings so that their presence\n        \/\/ in the Javascript doesn't mess up the search for\n        \/\/ the MATLAB code.\n        t1='75b0bee11e3e4c56aafb679a67295a76 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 75b0bee11e3e4c56aafb679a67295a76';\n    \n        b=document.getElementsByTagName('body')[0];\n        i1=b.innerHTML.indexOf(t1)+t1.length;\n        i2=b.innerHTML.indexOf(t2);\n \n        code_string = b.innerHTML.substring(i1, i2);\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\n\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \n        \/\/ in the XML parser.\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\n        \/\/ doesn't go ahead and substitute the less-than character. \n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\n\n        copyright = 'Copyright 2013 The MathWorks, Inc.';\n\n        w = window.open();\n        d = w.document;\n        d.write('\n\n\n\n\n\n\n\n<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\n\n\n\n\n\n\n\n\n\\n');\n\n        d.title = title + ' (MATLAB code)';\n        d.close();\n    }\n\/\/ ]]><\/script><\/p>\n<p style=\"text-align: right; font-size: xx-small; font-weight: lighter; font-style: italic; color: gray;\"><a><span style=\"font-size: x-small; font-style: italic;\">Get<br \/>\nthe MATLAB code<noscript>(requires JavaScript)<\/noscript><\/span><\/a><\/p>\n<p>Published with MATLAB\u00ae R2013a<\/p>\n<p class=\"footer\">Published with MATLAB\u00ae R2013a<\/p>\n<\/div>\n<p><!--\n75b0bee11e3e4c56aafb679a67295a76 ##### SOURCE BEGIN #####\n%% Using Symbolic Math Toolbox to Compute Area Moments of Inertia\n% Once more I am pleased to introduce guest blogger Kai Gehrs. Kai has been\n% a Software Engineer at MathWorks for the past five years working on the\n% <https:\/\/www.mathworks.com\/help\/toolbox\/symbolic Symbolic Math Toolbox>.\n% He has a background in mathematics and computer science. He already\n% contributed to my BLOG in the past writing about\n% <https:\/\/blogs.mathworks.com\/loren\/2012\/07\/27\/using-symbolic-equations-and-symbolic-functions-in-matlab Using Symbolic Equations And Symbolic Functions In MATLAB>\n% as well as on approaches for <https:\/\/blogs.mathworks.com\/loren\/2011\/10\/25\/simplifying-symbolic-results Simplifying Symbolic Results>.\n%\n%% In a Nutshell: What Is This Article About?\n% If you are interested in using MATLAB and the Symbolic Math Toolbox\n% in teaching some basics in mechanical engineering, this might be of\n% interest to you.\n%\n% Computing _area moments of inertia_ is an important task in\n% mechanics. For example, area moments of inertia play a critical role in stress, dynamic, and stability\n% analysis of structures. In this article, we use capabilities of the\n% <www.mathworks.com\/products\/symbolic Symbolic Math Toolbox> to compute\n% area moments for cross sections of elliptical tubes.\n%\n% We start with a basic case involving only numeric parameters, and\n% then make the computations more general by introducing symbolic\n% parameters.\n%\n% All plots used in this article have been created in the\n% <https:\/\/www.mathworks.com\/discovery\/mupad.html MuPAD Notebook app>.\n% You can find the source code for these plots at the end of the article.\n%\n%% Basic Example: Cross Section of an Elliptical Tube\n% The following picture shows the cross section of an elliptical tube.\n%\n% <<lateralCut.png>>\n%\n% The following ellipses define outer and inner contours of the section:\n%\n% * $y = \\pm 2 \\, \\sqrt{1 - \\frac{x^2}{9}}$ describe the _outer_\n% contour line for $x \\in [-3,3]$.\n% * $y = \\pm \\sqrt{1 - \\frac{x^2}{4}}$ describe the _inner_\n% contour line for $x \\in [-2,2]$.\n%\n% We will compute the area moment of inertia of this section with respect\n% to the $x$-axis.\n%\n%% Area Moment Of Inertia\n% The moment of inertia of an area $A$ with respect to the $x$-axis\n% is defined in terms of the double integral\n%\n% $$I_x = {\\int \\int}_A y^2\\, \\mathrm{d}A.$$\n%\n% You can find this definition on\n% <http:\/\/en.wikipedia.org\/wiki\/Second_moment_of_area Wikipedia>.\n%\n%% Math Behind This Example\n% In our example, the area $A$ is the hatched area shown in the\n% previous plot. Taking advantage of the symmetry of the section with respect\n% to both the $x$- and $y$-axes, we can restrict our computations to\n% the first quadrant of the $x-y$-plane.\n%\n% <<UpperQuarterLateralCut.png>>\n%\n% To compute the necessary double integral over the hatched area,\n% we divide this area into two separate areas $A1$ and $A2$.\n%\n% <<IntegrationAreas.png>>\n%\n% To compute the final area moment of inertia about the x-axis, we multiply\n% the sum of the double integrals over $A_1$ and $A_2$ by four.\n%\n% $$I_x = {\\int\\int}_A y^2 \\, \\mathrm dA = 4 \\cdot \\Big({\\int\\int}_{A_1}\n% y^2 \\, \\mathrm dA_1 + {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2\\Big).$$\n%\n% Now, we only need to compute these two double integrals. The mathematical\n% theory behind multi-dimensional integration tells us that we can rewrite\n% each of these double integrals in terms of two integrals:\n%\n% $$I_1 = {\\int\\int}_{A_1} y^2 \\, \\mathrm dA_1 = \\int_0^2 \\int_{\\sqrt{1 -\n% \\frac{x^2}{4}}}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}} y^2 \\, \\mathrm dy \\,\n% \\mathrm dx,$$\n%\n% $$I_2 = {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2 = \\int_2^3 \\int_{0}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}}\n% y^2\\, \\mathrm dy\\, \\mathrm dx.$$\n%\n% The <https:\/\/www.mathworks.com\/help\/toolbox\/symbolic\/int.html |int|>\n% function from the Symbolic Math Toolbox can compute these integrals.\n%% Symbolic Integration\n%\n% We define $x$ and $y$ as symbolic variables:\nsyms x y;\n%%\n% We start computing the inner integral $\\int_{\\sqrt{1 - \\frac{x^2}{4}}}^{2 \\, \\sqrt{1 -\n% \\frac{x^2}{9}}} y^2 \\, \\mathrm dy$ of $I_1$:\ninnerI1 = int(y^2, y, sqrt(1-x^2\/4), 2*sqrt(1-x^2\/9))\n%%\n% Next, we integrate |innerI1| with respect to $x$ from $0$ to $2$.\n% This provides $I_1$:\nI1 = int(innerI1, x, 0, 2)\n%%\n% We use the same strategy to compute $I_2$. First, we compute the inner\n% integral $\\int_{0}^{2 \\, \\sqrt{1 - \\frac{x^2}{9}}} y^2\\, \\mathrm dy$.\n% Then, we integrate the resulting expression with respect to $x$ from $2$\n% to $3$:\nI2 = int(int(y^2, y, 0, 2*sqrt(1-x^2\/9)), x, 2, 3);\n%%\n% Hence, the area moment of inertia of the elliptical tube with respect to\n% the $x$-axis is\nIx = 4 * (I1 + I2)\n%%\n% We can approximate the result numerically by using the\n% <https:\/\/www.mathworks.com\/help\/toolbox\/symbolic\/vpa.html |vpa|>\n% function. For example, we approximate the result using 5 significant\n% (nonzero) digits:\nvpa(Ix, 5)\n%% Advanced Example: Cross Section of an Elliptical Tube Defined by Symbolic Parameters\n%\n% We can use a numerical integrator, such as MATLAB's\n% <http:\/\/www.mathworks.de\/de\/help\/matlab\/ref\/integral2.html |integral2|>,\n% to compute the area moment of inertia in the previous example. A\n% numerical integrator might return slightly less accurate results, but\n% other than that there is not much benefit from using symbolic\n% integration there.\n%\n% But what if we consider the cross section of a more general elliptical\n% tube whose shape is defined by arbitrary symbolic parameters?\n%\n% For example, we want to compute the area moment of inertia of this\n% elliptical tube.\n%\n% <<LateralCutMoreGeneral.png>>\n%\n% Its contour lines are still ellipses, but they are defined by the\n% symbolic parameters $r_1$, $r_2$, $R_1$ and $R_2$:\n%\n% * $y = \\pm R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}$ describe the outer\n% contour line for $x \\in [-R_1,R_1]$.\n% * $y = \\pm r_2 \\, \\sqrt{1 - \\frac{x^2}{r_1^2}}$ describe the inner\n% contour line for $x \\in [-r_1,r_1]$.\n%\n% Now we need to compute these integrals:\n%\n% $$I_1 = {\\int\\int}_{A_1} y^2 \\, \\mathrm dA_1 = \\int_0^{r_1} \\int_{r_2 \\, \\sqrt{1 -\n% \\frac{x^2}{r_1^2}}}^{R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}} y^2 \\, \\mathrm dy \\,\n% \\mathrm dx,$$\n%\n% $$I_2 = {\\int\\int}_{A_2} y^2 \\, \\mathrm dA_2 = \\int_{r_1}^{R_1}\n% \\int_{0}^{R_2 \\, \\sqrt{1 - \\frac{x^2}{R_1^2}}} y^2\\, \\mathrm dy\\, \\mathrm dx.$$\n%\n% Although the situation is more complicated now, we can still use the same\n% strategy as above, using symbolic variables instead of numbers.\n% Nevertheless, we need to add one more step: specify relationships between\n% symbolic parameters. This is because the variables $r_1$, $r_2$, $R_1$,\n% and $R_2$ can be any complex number, unless we explicitly restrict their\n% values. For example, the system does not know if these variables are\n% positive or negative, if $r_1 < R_1$ or otherwise. In this example, we % want to specify that $r_1>0$, $r_2>0$, $R_1>r_1$, and $R_2>r_2$. In the\n% Symbolic Math Toolbox, such restrictions are called _assumptions on\n% variables_. They can be set by using the\n% <http:\/\/www.mathworks.de\/de\/help\/matlab\/ref\/assume.html |assume|>\n% and\n% <http:\/\/www.mathworks.de\/de\/help\/matlab\/ref\/assumeAlso.html |assumeAlso|>\n% functions:\nsyms x y r1 r2 R1 R2;\nassume(r1 > 0);\nassume(r2 > 0);\nassumeAlso(r1 < R1);\nassumeAlso(r2 < R2);\nI1 = int(int(y^2, y, r2*sqrt(1-x^2\/r1^2), R2*sqrt(1-x^2\/R1^2)), x, 0, r1);\nI2 = int(int(y^2, y, 0, R2*sqrt(1-x^2\/R1^2)), x, r1, R1);\nIxgeneral = 4 * (I1 + I2);\npretty(Ixgeneral)\n%%\n% Now, let's substitute our symbolic parameters  with the same values that\n% we had in the first example. Do we get the same result? Let's double-check:\n% when substituting $r_1=2,r_2,=1,R_1=3,R_2=2$ in |Ixgeneral|, do we get the\n% same result for |Ix| as obtained for the section initially considered?\n%\n% As expected, the general solution |Ixgeneral| reduces to the specific\n% solution |Ix| when we substitute our original dimensions:\nvpa(subs(Ixgeneral, [r1, R1, r2, R2],[2, 3, 1, 2]), 5)\n%%\n% Note that it is essential to set the assumptions on $r_1$, $r_2$, $R_1$\n% and $R_2$. By default, the Symbolic Math Toolbox assumes that all variables\n% including symbolic parameters are arbitrary complex numbers. This makes\n% computing integrals much more complicated and time-consuming. For example,\n% try computing |Ixgeneral| without setting the assumptions on $r_1$, $r_2$,\n% $R_1$ and $R_2$. In general, it can be impossible to even get the\n% result without any additional assumptions.\n%\n%% How to Create MuPAD Graphics Shown in This Article\n% To generate the graphics shown in this article, use\n% <https:\/\/www.mathworks.com\/discovery\/mupad.html MuPAD Notebook app>.\n% Note that the code in this section does not run in the MATLAB Command\n% Window.\n%\n% To open the MuPAD Notebook app:\n%\n% # On the MATLAB Toolstrip, click the Apps tab.\n% # On the Apps tab, click the down arrow at the end of the Apps section.\n% # Under Math, Statistics and Optimization, click the MuPAD Notebook button.\n%\n% Alternatively, type |mupad| in the MATLAB Command Window.\n%\n% Paste the following commands to a MuPAD notebook to obtain the four\n% graphics used above:\n%\n%  Ellipse:= (x,a,b) -> sqrt(b^2*(1-x^2\/a^2)):\n%\n%  f1 := plot::Function2d(Ellipse(x,2,1), x = -3..3):\n%  f2 := plot::Function2d(Ellipse(x,3,2), x = -3..3):\n%  f3 := plot::Function2d(Ellipse(x,3,2), x = -3..-2):\n%  f4 := plot::Function2d(Ellipse(x,3,2), x = 2..3):\n%\n%  f5 := plot::Function2d(-Ellipse(x,2,1), x = -3..3):\n%  f6 := plot::Function2d(-Ellipse(x,3,2), x = -3..3):\n%  f7 := plot::Function2d(-Ellipse(x,3,2), x = -3..-2):\n%  f8 := plot::Function2d(-Ellipse(x,3,2), x = 2..3):\n%\n%  f9 := plot::Function2d(Ellipse(x,2,1), x = 0..2, VerticalAsymptotesVisible = FALSE):\n%  f10 := plot::Function2d(Ellipse(x,3,2), x = 0..3, VerticalAsymptotesVisible = FALSE):\n%  f11 := plot::Function2d(Ellipse(x,3,2), x = 2..3, VerticalAsymptotesVisible = FALSE):\n%\n%  plot(plot::Hatch(f1,f2),plot::Hatch(f3),plot::Hatch(f4),f1,f2,\n%       plot::Hatch(f5,f6),plot::Hatch(f7),plot::Hatch(f8),f5,f6,\n%       Scaling = Constrained,\n%       GridVisible = TRUE,\n%       VerticalAsymptotesVisible = FALSE,\n%       Height = 120, Width = 200,\n%       Header = \"Cross section of elliptical tube (x-y-plane)\"):\n%\n%  plot(plot::Hatch(f9,f10),plot::Hatch(f11),f9,f10,\n%       Scaling = Constrained,\n%       GridVisible = TRUE,\n%       Height = 120, Width = 200,\n%       Header = \"Restriction to first quadrant (x-y-plane)\"):\n%\n%  plot(plot::Hatch(f9,f10,Color=RGB::Magenta),plot::Hatch(f11,Color = RGB::Green),f9,f10,\n%       Scaling=Constrained,\n%       GridVisible = TRUE,\n%       plot::Text2d(\"A1\",[1.0,1.25]),\n%       plot::Text2d(\"A2\",[2.5,0.25]),\n%       Height = 120, Width = 200,\n%       Header = \"Integration areas\"):\n%\n%  plot(plot::Hatch(f1,f2),plot::Hatch(f3),plot::Hatch(f4),f1,f2,\n%       plot::Hatch(f5,f6),plot::Hatch(f7),plot::Hatch(f8),f5,f6,\n%       plot::Text2d(\"r1\",[1.8,0.05]),\n%       plot::Text2d(\"r2\",[0.05,0.85]),\n%       plot::Text2d(\"R1\",[2.74,0.05]),\n%       plot::Text2d(\"R2\",[0.05,1.84]),\n%       Scaling = Constrained,\n%       GridVisible = TRUE,\n%       VerticalAsymptotesVisible = FALSE,\n%       XTicksLabelsVisible = FALSE,\n%       YTicksLabelsVisible = FALSE,\n%       ViewingBox = [-3.5..3.5,-2.5..2.5],\n%       Height = 120, Width = 200,\n%       Header = \"Cross section of elliptical tube (more general situation)\"):\n%% Have You Tried Symbolic Math Toolbox?\n% Have you used the Symbolic Math Toolbox in computations related to\n% mechanics? Let me know <https:\/\/blogs.mathworks.com\/loren\/?p=684#respond here>.\n##### SOURCE END ##### 75b0bee11e3e4c56aafb679a67295a76\n--><\/p>\n","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/loren\/2013\/LateralCutMoreGeneral.png\" onError=\"this.style.display ='none';\" \/><\/div>\n<p><!--introduction-->Once more I am pleased to introduce guest blogger Kai Gehrs. Kai has been a Software Engineer at MathWorks for the past five years working on the <a href=\"https:\/\/www.mathworks.com\/help\/toolbox\/symbolic\">Symbolic Math Toolbox<\/a>. He has a background in mathematics and computer science. He already contributed to my BLOG in the past writing about <a href=\"https:\/\/blogs.mathworks.com\/loren\/2012\/07\/27\/using-symbolic-equations-and-symbolic-functions-in-matlab\">Using Symbolic Equations And Symbolic Functions In MATLAB<\/a> as well as on approaches for <a href=\"https:\/\/blogs.mathworks.com\/loren\/2011\/10\/25\/simplifying-symbolic-results\">Simplifying Symbolic Results<\/a>.<\/p>\n<p><!--\/introduction-->... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/loren\/2013\/04\/26\/using-symbolic-math-toolbox-to-compute-area-moments-of-inertia\/\">read more >><\/a><\/p>\n","protected":false},"author":39,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[38],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/684"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/users\/39"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/comments?post=684"}],"version-history":[{"count":7,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/684\/revisions"}],"predecessor-version":[{"id":2185,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/posts\/684\/revisions\/2185"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/media?parent=684"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/categories?post=684"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/loren\/wp-json\/wp\/v2\/tags?post=684"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}