{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Strain and stress tensors in Cartesian coordinates\n", "\n", "This worksheet demonstrates a few capabilities of [SageManifolds](http://sagemanifolds.obspm.fr/) (version 1.0, as included in SageMath 7.5) in computations regarding elasticity theory in Cartesian coordinates.\n", "\n", "Click [here](https://raw.githubusercontent.com/sagemanifolds/SageManifolds/master/Worksheets/v1.0/SM_elasticity_Cartesian.ipynb) to download the worksheet file (ipynb format). To run it, you must start SageMath with the Jupyter notebook, via the command `sage -n jupyter`" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "*NB:* a version of SageMath at least equal to 7.5 is required to run this worksheet:" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "'SageMath version 7.5.1, Release Date: 2017-01-15'" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "version()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First we set up the notebook to display mathematical objects using LaTeX rendering:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "%display latex" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Euclidean 3-space and Cartesian coordinates" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We introduce the Euclidean space as a 3-dimensional differentiable manifold:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3-dimensional differentiable manifold M\n" ] } ], "source": [ "M = Manifold(3, 'M', start_index=1)\n", "print(M)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We then introduce the Cartesian coordinates $(x,y,z)$ as a chart $X$ on $M$:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Chart (M, (x, y, z))\n" ] }, { "data": { "text/html": [ "" ], "text/plain": [ "Chart (M, (x, y, z))" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "X. = M.chart()\n", "print(X)\n", "X" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The associated vector frame is" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "Coordinate frame (M, (d/dx,d/dy,d/dz))" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "X.frame()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We shall expand vector and tensor fields not on this frame, which is the default one on $M$:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "Coordinate frame (M, (d/dx,d/dy,d/dz))" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "M.default_frame()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Displacement vector and strain tensor\n", "\n", "Let us define the **displacement vector** $U$ in terms of its components w.r.t. the orthonormal Cartesian frame:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "U = U_x(x, y, z) d/dx + U_y(x, y, z) d/dy + U_z(x, y, z) d/dz" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "U = M.vector_field(name='U')\n", "U[:] = [function('U_x')(x,y,z), function('U_y')(x,y,z), \n", " function('U_z')(x,y,z)]\n", "U.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The following computations will involve the metric $g$ of the Euclidean space. At the current stage of SageManifolds, we need to introduce it explicitly, as a Riemannian metric on the manifold $M$ (in a future version of SageManifolds, one shall to declare $M$ as an Euclidean space, and not merely as a manifold, so that it will come equipped with $g$):" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Riemannian metric g on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "g = M.riemannian_metric('g')\n", "print(g)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We initialize $g$ by declaring that its components with respect to the frame of Cartesian coordinates are\n", "$\\mathrm{diag}(1,1,1)$:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "g = dx*dx + dy*dy + dz*dz" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g[1,1], g[2,2], g[3,3] = 1, 1, 1\n", "g.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The covariant derivative operator $\\nabla$ is introduced as the (Levi-Civita) connection associated with $g$: " ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Levi-Civita connection nabla_g associated with the Riemannian metric g on the 3-dimensional differentiable manifold M\n" ] }, { "data": { "text/html": [ "" ], "text/plain": [ "Levi-Civita connection nabla_g associated with the Riemannian metric g on the 3-dimensional differentiable manifold M" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nabla = g.connection()\n", "print(nabla)\n", "nabla" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The covariant derivative of the displacement vector $U$ is" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Tensor field nabla_g(U) of type (1,1) on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "nabU = nabla(U)\n", "print(nabU)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "nabla_g(U) = d(U_x)/dx d/dx*dx + d(U_x)/dy d/dx*dy + d(U_x)/dz d/dx*dz + d(U_y)/dx d/dy*dx + d(U_y)/dy d/dy*dy + d(U_y)/dz d/dy*dz + d(U_z)/dx d/dz*dx + d(U_z)/dy d/dz*dy + d(U_z)/dz d/dz*dz" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nabU.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We convert it to a tensor field of type (0,2) (i.e. a bilinear form) by lowering the upper index with $g$:" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Tensor field of type (0,2) on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "nabU_form = nabU.down(g)\n", "print(nabU_form)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "d(U_x)/dx dx*dx + d(U_x)/dy dx*dy + d(U_x)/dz dx*dz + d(U_y)/dx dy*dx + d(U_y)/dy dy*dy + d(U_y)/dz dy*dz + d(U_z)/dx dz*dx + d(U_z)/dy dz*dy + d(U_z)/dz dz*dz" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nabU_form.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The **strain tensor** $\\varepsilon$ is defined as the symmetrized part of this tensor:" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Field of symmetric bilinear forms on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "E = nabU_form.symmetrize()\n", "print(E)" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "E = d(U_x)/dx dx*dx + (1/2*d(U_x)/dy + 1/2*d(U_y)/dx) dx*dy + (1/2*d(U_x)/dz + 1/2*d(U_z)/dx) dx*dz + (1/2*d(U_x)/dy + 1/2*d(U_y)/dx) dy*dx + d(U_y)/dy dy*dy + (1/2*d(U_y)/dz + 1/2*d(U_z)/dy) dy*dz + (1/2*d(U_x)/dz + 1/2*d(U_z)/dx) dz*dx + (1/2*d(U_y)/dz + 1/2*d(U_z)/dy) dz*dy + d(U_z)/dz dz*dz" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E.set_name('E', latex_name=r'\\varepsilon')\n", "E.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us display the components of $\\varepsilon$, skipping those that can be deduced by symmetry:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "E_xx = d(U_x)/dx \n", "E_xy = 1/2*d(U_x)/dy + 1/2*d(U_y)/dx \n", "E_xz = 1/2*d(U_x)/dz + 1/2*d(U_z)/dx \n", "E_yy = d(U_y)/dy \n", "E_yz = 1/2*d(U_y)/dz + 1/2*d(U_z)/dy \n", "E_zz = d(U_z)/dz " ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "E.display_comp(only_nonredundant=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Stress tensor and Hooke's law\n", "\n", "To form the stress tensor according to Hooke's law, we introduce first the Lamé constants:" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "ll" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "var('ll', latex_name=r'\\lambda')" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "mu" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "var('mu', latex_name=r'\\mu')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The trace (with respect to $g$) of the bilinear form $\\varepsilon$ is obtained by (i) raising the first index (`pos=0`) by means of $g$ and (ii) by taking the trace of the resulting endomorphism:" ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Scalar field on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "trE = E.up(g, pos=0).trace()\n", "print(trE)" ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "M --> R\n", "(x, y, z) |--> d(U_x)/dx + d(U_y)/dy + d(U_z)/dz" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "trE.display()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The **stress tensor** $S$ is obtained via Hooke's law for isotropic material:\n", "$$ S = \\lambda \\, \\mathrm{tr}\\varepsilon \\; g + 2\\mu \\, \\varepsilon$$" ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Field of symmetric bilinear forms on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "S = ll*trE*g + 2*mu*E\n", "print(S)" ] }, { "cell_type": "code", "execution_count": 23, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "S = ((ll + 2*mu)*d(U_x)/dx + ll*d(U_y)/dy + ll*d(U_z)/dz) dx*dx + (mu*d(U_x)/dy + mu*d(U_y)/dx) dx*dy + (mu*d(U_x)/dz + mu*d(U_z)/dx) dx*dz + (mu*d(U_x)/dy + mu*d(U_y)/dx) dy*dx + (ll*d(U_x)/dx + (ll + 2*mu)*d(U_y)/dy + ll*d(U_z)/dz) dy*dy + (mu*d(U_y)/dz + mu*d(U_z)/dy) dy*dz + (mu*d(U_x)/dz + mu*d(U_z)/dx) dz*dx + (mu*d(U_y)/dz + mu*d(U_z)/dy) dz*dy + (ll*d(U_x)/dx + ll*d(U_y)/dy + (ll + 2*mu)*d(U_z)/dz) dz*dz" ] }, "execution_count": 23, "metadata": {}, "output_type": "execute_result" } ], "source": [ "S.set_name('S')\n", "S.display()" ] }, { "cell_type": "code", "execution_count": 24, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "S_xx = (ll + 2*mu)*d(U_x)/dx + ll*d(U_y)/dy + ll*d(U_z)/dz \n", "S_xy = mu*d(U_x)/dy + mu*d(U_y)/dx \n", "S_xz = mu*d(U_x)/dz + mu*d(U_z)/dx \n", "S_yy = ll*d(U_x)/dx + (ll + 2*mu)*d(U_y)/dy + ll*d(U_z)/dz \n", "S_yz = mu*d(U_y)/dz + mu*d(U_z)/dy \n", "S_zz = ll*d(U_x)/dx + ll*d(U_y)/dy + (ll + 2*mu)*d(U_z)/dz " ] }, "execution_count": 24, "metadata": {}, "output_type": "execute_result" } ], "source": [ "S.display_comp(only_nonredundant=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Each component can be accessed individually:" ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "mu*d(U_x)/dy + mu*d(U_y)/dx" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "S[1,2]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Divergence of the stress tensor\n", "\n", "The divergence of the stress tensor is the 1-form:\n", "$$ f_i = \\nabla_j S^j_{\\ \\, i} $$\n", "In a next version of SageManifolds, there will be a function `divergence()`. For the moment, to evaluate $f$, \n", "we first form the tensor $S^j_{\\ \\, i}$ by raising the first index (`pos=0`) of $S$ with $g$:" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Tensor field of type (1,1) on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "SU = S.up(g, pos=0)\n", "print(SU)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The divergence is obtained by taking the trace on the first index (`0`) and the third one (`2`) of the tensor\n", "$(\\nabla S)^j_{\\ \\, ik} = \\nabla_k S^j_{\\ \\, i}$:" ] }, { "cell_type": "code", "execution_count": 27, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1-form on the 3-dimensional differentiable manifold M\n" ] } ], "source": [ "divS = nabla(SU).trace(0,2)\n", "print(divS)" ] }, { "cell_type": "code", "execution_count": 28, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "f = ((ll + 2*mu)*d^2(U_x)/dx^2 + mu*d^2(U_x)/dy^2 + mu*d^2(U_x)/dz^2 + (ll + mu)*d^2(U_y)/dxdy + (ll + mu)*d^2(U_z)/dxdz) dx + ((ll + mu)*d^2(U_x)/dxdy + mu*d^2(U_y)/dx^2 + (ll + 2*mu)*d^2(U_y)/dy^2 + mu*d^2(U_y)/dz^2 + (ll + mu)*d^2(U_z)/dydz) dy + ((ll + mu)*d^2(U_x)/dxdz + (ll + mu)*d^2(U_y)/dydz + mu*d^2(U_z)/dx^2 + mu*d^2(U_z)/dy^2 + (ll + 2*mu)*d^2(U_z)/dz^2) dz" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "divS.set_name('f')\n", "divS.display()" ] }, { "cell_type": "code", "execution_count": 29, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "f_x = (ll + 2*mu)*d^2(U_x)/dx^2 + mu*d^2(U_x)/dy^2 + mu*d^2(U_x)/dz^2 + (ll + mu)*d^2(U_y)/dxdy + (ll + mu)*d^2(U_z)/dxdz \n", "f_y = (ll + mu)*d^2(U_x)/dxdy + mu*d^2(U_y)/dx^2 + (ll + 2*mu)*d^2(U_y)/dy^2 + mu*d^2(U_y)/dz^2 + (ll + mu)*d^2(U_z)/dydz \n", "f_z = (ll + mu)*d^2(U_x)/dxdz + (ll + mu)*d^2(U_y)/dydz + mu*d^2(U_z)/dx^2 + mu*d^2(U_z)/dy^2 + (ll + 2*mu)*d^2(U_z)/dz^2 " ] }, "execution_count": 29, "metadata": {}, "output_type": "execute_result" } ], "source": [ "divS.display_comp()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Displaying the components one by one:" ] }, { "cell_type": "code", "execution_count": 30, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "(ll + 2*mu)*d^2(U_x)/dx^2 + mu*d^2(U_x)/dy^2 + mu*d^2(U_x)/dz^2 + (ll + mu)*d^2(U_y)/dxdy + (ll + mu)*d^2(U_z)/dxdz" ] }, "execution_count": 30, "metadata": {}, "output_type": "execute_result" } ], "source": [ "divS[1]" ] }, { "cell_type": "code", "execution_count": 31, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "(ll + mu)*d^2(U_x)/dxdy + mu*d^2(U_y)/dx^2 + (ll + 2*mu)*d^2(U_y)/dy^2 + mu*d^2(U_y)/dz^2 + (ll + mu)*d^2(U_z)/dydz" ] }, "execution_count": 31, "metadata": {}, "output_type": "execute_result" } ], "source": [ "divS[2]" ] }, { "cell_type": "code", "execution_count": 32, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "(ll + mu)*d^2(U_x)/dxdz + (ll + mu)*d^2(U_y)/dydz + mu*d^2(U_z)/dx^2 + mu*d^2(U_z)/dy^2 + (ll + 2*mu)*d^2(U_z)/dz^2" ] }, "execution_count": 32, "metadata": {}, "output_type": "execute_result" } ], "source": [ "divS[3]" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 7.5.1", "language": "", "name": "sagemath" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.13" } }, "nbformat": 4, "nbformat_minor": 0 }