diff --git a/notebooks/2015-06-25-elasticity-solver-example.ipynb b/notebooks/2015-06-25-elasticity-solver-example.ipynb index c4b82ae..361a2d5 100644 --- a/notebooks/2015-06-25-elasticity-solver-example.ipynb +++ b/notebooks/2015-06-25-elasticity-solver-example.ipynb @@ -50,6 +50,27 @@ "execution_count": 1, "metadata": {}, "output_type": "execute_result" + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ + " in depwarn at /Applications/Julia-0.4.0-dev-539c818c4e.app/Contents/Resources/julia/lib/julia/sys.dylib\n", + " in int32 at deprecated.jl:49\n", + " in recv at /Users/jukka/.julia/v0.4/ZMQ/src/ZMQ.jl:617\n", + " in recv_ipython at /Users/jukka/.julia/v0.4/IJulia/src/msg.jl:63\n", + " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:120\n", + " in anonymous at task.jl:365\n", + "while loading /Users/jukka/.julia/v0.4/IJulia/src/kernel.jl, in expression starting on line 35\n", + "WARNING: int32(x) is deprecated, use Int32(x) instead.\n", + " in depwarn at /Applications/Julia-0.4.0-dev-539c818c4e.app/Contents/Resources/julia/lib/julia/sys.dylib\n", + " in int32 at deprecated.jl:49\n", + " in recv at /Users/jukka/.julia/v0.4/ZMQ/src/ZMQ.jl:617\n", + " in recv_ipython at /Users/jukka/.julia/v0.4/IJulia/src/msg.jl:63\n", + " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:120\n", + " in anonymous at task.jl:365\n", + "while loading /Users/jukka/.julia/v0.4/IJulia/src/kernel.jl, in expression starting on line 35\n" + ] } ], "source": [ @@ -73,7 +94,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": 2, "metadata": { "collapsed": false }, @@ -84,7 +105,7 @@ "calc_local_matrices! (generic function with 1 method)" ] }, - "execution_count": 18, + "execution_count": 2, "metadata": {}, "output_type": "execute_result" } @@ -218,7 +239,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": 3, "metadata": { "collapsed": true }, @@ -229,7 +250,7 @@ }, { "cell_type": "code", - "execution_count": 19, + "execution_count": 4, "metadata": { "collapsed": false, "scrolled": false @@ -246,17 +267,17 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 21:42:03:DEBUG:root:Converged in 6 iterations.\n", - "11-Aug 21:42:03:DEBUG:root:solution vector: \n", + "13-Aug 00:37:26:DEBUG:root:Converged in 6 iterations.\n", + "13-Aug 00:37:26:DEBUG:root:solution vector: \n", " [0.0 -0.3991450609547433 -0.07228582695592461 0.0\n", " 0.0 -2.1779892317073504 -2.222244754401764 0.0]\n", - "11-Aug 21:42:03:DEBUG:root:norm of u: 3.1292483947150047\n", - "11-Aug 21:42:03:DEBUG:root:Converged in 6 iterations.\n", - "11-Aug 21:42:03:DEBUG:root:solution vector: \n", + "13-Aug 00:37:26:DEBUG:root:norm of u: 3.1292483947150047\n", + "13-Aug 00:37:27:DEBUG:root:Converged in 6 iterations.\n", + "13-Aug 00:37:27:DEBUG:root:solution vector: \n", " [0.0 0.7433248532717796 1.048521014723486 0.0\n", " 0.0 -2.085766534304891 -1.9606633242166027 0.0]\n", - "11-Aug 21:42:03:DEBUG:root:norm of u: 3.1292483947150056\n", - "11-Aug 21:42:03:DEBUG:root:Iteration 1\n" + "13-Aug 00:37:27:DEBUG:root:norm of u: 3.1292483947150056\n", + "13-Aug 00:37:27:DEBUG:root:Iteration 1\n" ] }, { @@ -270,15 +291,38 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 21:42:04:DEBUG:root:Norm of du: 0.5992228342549063\n", - "11-Aug 21:42:04:DEBUG:root:Element displacement: [0.0 -0.02264423092574128 0.022536491965822688 0.0\n", + "13-Aug 00:37:27:DEBUG:root:Norm of du: 0.5992228342549063\n", + "13-Aug 00:37:27:DEBUG:root:Element displacement: [0.0 -0.02264423092574128 0.022536491965822688 0.0\n", " 0.0 -0.12688379170176511 -0.12679760053383046 0.0]\n", - "11-Aug 21:42:04:DEBUG:root:Element displacement: [-0.02264423092574128 -0.029998807330416523 0.030242156525002267 0.022536491965822688\n", + "13-Aug 00:37:27:DEBUG:root:Element displacement: [-0.02264423092574128 -0.029998807330416523 0.030242156525002267 0.022536491965822688\n", " -0.12688379170176511 -0.4041002710283739 -0.40446733271991214 -0.12679760053383046]\n", - "11-Aug 21:42:04:DEBUG:root:solution vector: \n", + "13-Aug 00:37:27:DEBUG:root:solution vector: \n", " [-0.02264423092574128 -0.029998807330416523 0.022536491965822688 0.030242156525002267 0.0 0.0\n", " -0.12688379170176511 -0.4041002710283739 -0.12679760053383046 -0.40446733271991214 0.0 0.0]\n", - "11-Aug 21:42:04:DEBUG:root:norm of u: 0.5992228342549063\n" + "13-Aug 00:37:27:DEBUG:root:norm of u: 0.5992228342549063\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + ":\n", + " 132.8 0.0 23.6 -3.0 -113.6 … 23.6 3.0 -33.2 15.0\n", + " 0.0 324.8 3.0 77.6 0.0 -3.0 77.6 15.0 -81.2\n", + " 23.6 3.0 66.4 -15.0 -33.2 0.0 0.0 0.0 0.0\n", + " -3.0 77.6 -15.0 162.4 15.0 0.0 0.0 0.0 0.0\n", + " -113.6 0.0 -33.2 15.0 132.8 -33.2 -15.0 23.6 -3.0\n", + " 0.0 -317.6 15.0 -81.2 0.0 … -15.0 -81.2 3.0 77.6\n", + " -33.2 -15.0 -56.8 3.0 23.6 0.0 0.0 0.0 0.0\n", + " -15.0 -81.2 -3.0 -158.8 3.0 0.0 0.0 0.0 0.0\n", + " 23.6 -3.0 0.0 0.0 -33.2 66.4 15.0 -56.8 3.0\n", + " 3.0 77.6 0.0 0.0 -15.0 15.0 162.4 -3.0 -158.8\n", + " -33.2 15.0 0.0 0.0 23.6 … -56.8 -3.0 66.4 -15.0\n", + " 15.0 -81.2 0.0 0.0 -3.0 3.0 -158.8 -15.0 162.4\n", + "K norm = 708.0378644377365\n", + "du = [-0.02264423092574128 -0.029998807330416523 0.022536491965822688 0.030242156525002267 0.0 0.0\n", + " -0.12688379170176511 -0.4041002710283739 -0.12679760053383046 -0.40446733271991214 0.0 0.0]\n", + "Out of 3 total facts:" ] }, { @@ -287,34 +331,9 @@ "delayed_handler (generic function with 4 methods)" ] }, - "execution_count": 19, + "execution_count": 4, "metadata": {}, "output_type": "execute_result" - }, - { - "name": "stdout", - "output_type": "stream", - "text": [ - ":\n", - " 132.8 0.0 23.6 -3.0 -113.6 0.0 -33.2 -15.0 23.6 3.0 -33.2 15.0\n", - " 0.0 324.8 3.0 77.6 0.0 -317.6 -15.0 -81.2 -3.0 77.6 15.0 -81.2\n", - " 23.6 3.0 66.4 -15.0 -33.2 15.0 -56.8 -3.0 0.0 0.0 0.0 0.0\n", - " -3.0 77.6 -15.0 162.4 15.0 -81.2 3.0 -158.8 0.0 0.0 0.0 0.0\n", - " -113.6 0.0 -33.2 15.0 132.8 0.0 23.6 3.0 -33.2 -15.0 23.6 -3.0\n", - " 0.0 -317.6 15.0 -81.2 0.0 324.8 -3.0 77.6 -15.0 -81.2 3.0 77.6\n", - " -33.2 -15.0 -56.8 3.0 23.6 -3.0 66.4 15.0 0.0 0.0 0.0 0.0\n", - " -15.0 -81.2 -3.0 -158.8 3.0 77.6 15.0 162.4 0.0 0.0 0.0 0.0\n", - " 23.6 -3.0 0.0 0.0 -33.2 -15.0 0.0 0.0 66.4 15.0 -56.8 3.0\n", - " 3.0 77.6 0.0 0.0 -15.0 -81.2 0.0 0.0 15.0 162.4 -3.0 -158.8\n", - " -33.2 15.0 0.0 0.0 23.6 3.0 0.0 0.0 -56.8 -3.0 66.4 -15.0\n", - " 15.0 -81.2 0.0 0.0 -3.0 77.6 0.0 0.0 3.0 -158.8 -15.0 162.4\n", - "K norm = 708.0378644377365\n", - "du = [-0.02264423092574128 -0.029998807330416523 0.022536491965822688 0.030242156525002267 0.0 0.0\n", - " -0.12688379170176511 -0.4041002710283739 -0.12679760053383046 -0.40446733271991214 0.0 0.0]\n", - "Out of 3 total facts:\n", - " Verified: 2\n", - " Pending: 1\n" - ] } ], "source": [ @@ -473,7 +492,17 @@ "metadata": { "collapsed": false }, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n", + " Verified: 2\n", + " Pending: 1\n" + ] + } + ], "source": [ "type Assembly\n", " # LHS\n", @@ -533,28 +562,30 @@ " -(1+xi[2])/4.0 (1-xi[1])/4.0]\n", " return basis, dbasis\n", " elseif (nnodes == 10) & (ndim == 3)\n", - " basis(xi) = [(xi[1] + xi[2] + xi[3] - 1)*(2*xi[1] + 2*xi[2] + 2*xi[3] - 1)\n", - " -xi[1]*(-2*xi[1] + 1)\n", - " -xi[2]*(-2*xi[2] + 1)\n", - " -xi[3]*(-2*xi[3] + 1)\n", - " 4*xi[1]*(-xi[1] - xi[2] - xi[3] + 1)\n", - " 4*xi[1]*xi[2]\n", - " 4*xi[2]*(-xi[1] - xi[2] - xi[3] + 1)\n", - " 4*xi[1]*xi[3]\n", - " 4*xi[2]*xi[3]\n", - " 4*xi[3]*(-xi[1] - xi[2] - xi[3] + 1)]\n", - "\n", + " basis(xi) = [\n", + " (xi[1] + xi[2] + xi[3] - 1)*(2*xi[1] + 2*xi[2] + 2*xi[3] - 1)\n", + " xi[1]*(2*xi[1] - 1)\n", + " xi[2]*(2*xi[2] - 1)\n", + " xi[3]*(2*xi[3] - 1)\n", + " -4*xi[1]*(xi[1] + xi[2] + xi[3] - 1)\n", + " 4*xi[1]*xi[2]\n", + " -4*xi[2]*(xi[1] + xi[2] + xi[3] - 1)\n", + " -4*xi[3]*(xi[1] + xi[2] + xi[3] - 1)\n", + " 4*xi[1]*xi[3]\n", + " 4*xi[2]*xi[3]\n", + " ]\n", " dbasis(xi) = [\n", - " 4*xi[1] + 4*xi[2] + 4*xi[3] - 3 4*xi[1] + 4*xi[2] + 4*xi[3] - 3 4*xi[1] + 4*xi[2] + 4*xi[3] - 3\n", - " 4*xi[1] - 1 0 0\n", - " 0 4*xi[2] - 1 0\n", - " 0 0 4*xi[3] - 1\n", - " -8*xi[1] - 4*xi[2] - 4*xi[3] + 4 -4*xi[1] -4*xi[1]\n", - " 4*xi[2] 4*xi[1] 0\n", - " -4*xi[2] -4*xi[1] - 8*xi[2] - 4*xi[3] + 4 -4*xi[2]\n", - " 4*xi[3] 0 4*xi[1]\n", - " 0 4*xi[3] 4*xi[2]\n", - " -4*xi[3] -4*xi[3] -4*xi[1] - 4*xi[2] - 8*xi[3] + 4]\n", + " 4*xi[1] + 4*xi[2] + 4*xi[3] - 3 4*xi[1] + 4*xi[2] + 4*xi[3] - 3 4*xi[1] + 4*xi[2] + 4*xi[3] - 3\n", + " 4*xi[1] - 1 0 0\n", + " 0 4*xi[2] - 1 0\n", + " 0 0 4*xi[3] - 1\n", + " -4*(2*xi[1] + xi[2] + xi[3] - 1) -4*xi[1] -4*xi[1]\n", + " 4*xi[2] 4*xi[1] 0\n", + " -4*xi[2] -4*(xi[1] + 2*xi[2] + xi[3] - 1) -4*xi[2]\n", + " -4*xi[3] -4*xi[3] -4*(xi[1] + xi[2] + 2*xi[3] - 1)\n", + " 4*xi[3] 0 4*xi[1]\n", + " 0 4*xi[3] 4*xi[2]\n", + " ]\n", " return basis, dbasis\n", " end\n", " throw(\"Unknown function space, ndim=$ndim, nnodes=$nnodes\")\n", @@ -668,28 +699,28 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 21:39:42:DEBUG:root:Adding nodes to array\n", - "11-Aug 21:39:42:DEBUG:root:Creating elements\n", - "11-Aug 21:39:42:DEBUG:root:Starting iteration 1\n", - "11-Aug 21:39:42:DEBUG:root:Assembling\n", - "11-Aug 21:39:42:DEBUG:root:Solution norm = 3.090022136728999\n", - "11-Aug 21:39:42:DEBUG:root:Starting iteration 2\n", - "11-Aug 21:39:42:DEBUG:root:Assembling\n", - "11-Aug 21:39:43:DEBUG:root:Solution norm = 0.32121316021535135\n", - "11-Aug 21:39:43:DEBUG:root:Starting iteration 3\n", - "11-Aug 21:39:43:DEBUG:root:Assembling\n", - "11-Aug 21:39:43:DEBUG:root:Solution norm = 0.040431781939994194\n", - "11-Aug 21:39:43:DEBUG:root:Starting iteration 4\n", - "11-Aug 21:39:43:DEBUG:root:Assembling\n", - "11-Aug 21:39:43:DEBUG:root:Solution norm = 0.0009291101052124042\n", - "11-Aug 21:39:43:DEBUG:root:Starting iteration 5\n", - "11-Aug 21:39:43:DEBUG:root:Assembling\n", - "11-Aug 21:39:43:DEBUG:root:Solution norm = 1.5638899022213175e-7\n", - "11-Aug 21:39:43:DEBUG:root:Starting iteration 6\n", - "11-Aug 21:39:43:DEBUG:root:Assembling\n", - "11-Aug 21:39:43:DEBUG:root:Solution norm = 1.0118539067290854e-14\n", - "11-Aug 21:39:43:DEBUG:root:Converged in 6 iterations.\n", - "11-Aug 21:39:43:DEBUG:root:Displacement of element = \n", + "13-Aug 00:37:30:DEBUG:root:Adding nodes to array\n", + "13-Aug 00:37:30:DEBUG:root:Creating elements\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 1\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 3.090022136728999\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 2\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 0.32121316021535135\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 3\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 0.040431781939994194\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 4\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 0.0009291101052124042\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 5\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 1.5638899022213175e-7\n", + "13-Aug 00:37:30:DEBUG:root:Starting iteration 6\n", + "13-Aug 00:37:30:DEBUG:root:Assembling\n", + "13-Aug 00:37:30:DEBUG:root:Solution norm = 1.0118539067290854e-14\n", + "13-Aug 00:37:30:DEBUG:root:Converged in 6 iterations.\n", + "13-Aug 00:37:30:DEBUG:root:Displacement of element = \n", "[-0.39914506095474334 -0.0722858269559246 0.0 0.0\n", " -2.1779892317073504 -2.2222447544017645 0.0 0.0]\n" ] @@ -813,48 +844,48 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 21:39:49:DEBUG:root:Creating elements\n", - "11-Aug 21:39:49:INFO:root:create_ldof2gdofmap: dofs per node: 2\n", - "11-Aug 21:39:49:DEBUG:root:Dict(4=>[7,8],2=>[3,4],3=>[5,6],5=>[9,10],6=>[11,12],1=>[1,2])\n", - "11-Aug 21:39:50:INFO:root:solve!: dofs per node: 2\n", - "11-Aug 21:39:50:DEBUG:root:Problem size = 12\n", - "11-Aug 21:39:50:DEBUG:root:Starting iteration 1\n", - "11-Aug 21:39:50:DEBUG:root:Assembling\n", - "11-Aug 21:39:50:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 21:39:50:DEBUG:root:Added 5 Lagrange multipliers to model\n", - "11-Aug 21:39:50:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 21:39:50:DEBUG:root:Solving system of equations. Total size = 16\n", - "11-Aug 21:39:51:DEBUG:root:Solution norm du = 0.6015838690633517\n", - "11-Aug 21:39:51:DEBUG:root:Starting iteration 2\n", - "11-Aug 21:39:51:DEBUG:root:Assembling\n", - "11-Aug 21:39:51:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 21:39:51:DEBUG:root:Added 5 Lagrange multipliers to model\n", - "11-Aug 21:39:51:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 21:39:51:DEBUG:root:Solving system of equations. Total size = 16\n", - "11-Aug 21:39:51:DEBUG:root:Solution norm du = 0.013420417380980414\n", - "11-Aug 21:39:51:DEBUG:root:Starting iteration 3\n", - "11-Aug 21:39:51:DEBUG:root:Assembling\n", - "11-Aug 21:39:51:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 21:39:51:DEBUG:root:Added 5 Lagrange multipliers to model\n", - "11-Aug 21:39:51:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 21:39:51:DEBUG:root:Solving system of equations. Total size = 16\n", - "11-Aug 21:39:51:DEBUG:root:Solution norm du = 0.00032202957936854873\n", - "11-Aug 21:39:51:DEBUG:root:Starting iteration 4\n", - "11-Aug 21:39:51:DEBUG:root:Assembling\n", - "11-Aug 21:39:51:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 21:39:51:DEBUG:root:Added 5 Lagrange multipliers to model\n", - "11-Aug 21:39:51:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 21:39:51:DEBUG:root:Solving system of equations. Total size = 16\n", - "11-Aug 21:39:51:DEBUG:root:Solution norm du = 9.906677094801476e-8\n", - "11-Aug 21:39:51:DEBUG:root:Starting iteration 5\n", - "11-Aug 21:39:51:DEBUG:root:Assembling\n", - "11-Aug 21:39:51:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 21:39:51:DEBUG:root:Added 5 Lagrange multipliers to model\n", - "11-Aug 21:39:51:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 21:39:51:DEBUG:root:Solving system of equations. Total size = 16\n", - "11-Aug 21:39:51:DEBUG:root:Solution norm du = 5.027624635820698e-15\n", - "11-Aug 21:39:51:DEBUG:root:Converged in 5 iterations.\n", - "11-Aug 21:39:51:DEBUG:root:Displacement of element = \n", + "13-Aug 00:37:31:DEBUG:root:Creating elements\n", + "13-Aug 00:37:31:INFO:root:create_ldof2gdofmap: dofs per node: 2\n", + "13-Aug 00:37:31:DEBUG:root:Dict(4=>[7,8],2=>[3,4],3=>[5,6],5=>[9,10],6=>[11,12],1=>[1,2])\n", + "13-Aug 00:37:32:INFO:root:solve!: dofs per node: 2\n", + "13-Aug 00:37:32:DEBUG:root:Problem size = 12\n", + "13-Aug 00:37:32:DEBUG:root:Starting iteration 1\n", + "13-Aug 00:37:32:DEBUG:root:Assembling\n", + "13-Aug 00:37:32:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:32:DEBUG:root:Added 5 Lagrange multipliers to model\n", + "13-Aug 00:37:32:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:32:DEBUG:root:Solving system of equations. Total size = 16\n", + "13-Aug 00:37:32:DEBUG:root:Solution norm du = 0.6015838690633517\n", + "13-Aug 00:37:32:DEBUG:root:Starting iteration 2\n", + "13-Aug 00:37:32:DEBUG:root:Assembling\n", + "13-Aug 00:37:32:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:32:DEBUG:root:Added 5 Lagrange multipliers to model\n", + "13-Aug 00:37:32:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:32:DEBUG:root:Solving system of equations. Total size = 16\n", + "13-Aug 00:37:32:DEBUG:root:Solution norm du = 0.013420417380980414\n", + "13-Aug 00:37:33:DEBUG:root:Starting iteration 3\n", + "13-Aug 00:37:33:DEBUG:root:Assembling\n", + "13-Aug 00:37:33:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:33:DEBUG:root:Added 5 Lagrange multipliers to model\n", + "13-Aug 00:37:33:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:33:DEBUG:root:Solving system of equations. Total size = 16\n", + "13-Aug 00:37:33:DEBUG:root:Solution norm du = 0.00032202957936854873\n", + "13-Aug 00:37:33:DEBUG:root:Starting iteration 4\n", + "13-Aug 00:37:33:DEBUG:root:Assembling\n", + "13-Aug 00:37:33:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:33:DEBUG:root:Added 5 Lagrange multipliers to model\n", + "13-Aug 00:37:33:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:33:DEBUG:root:Solving system of equations. Total size = 16\n", + "13-Aug 00:37:33:DEBUG:root:Solution norm du = 9.906677094801476e-8\n", + "13-Aug 00:37:33:DEBUG:root:Starting iteration 5\n", + "13-Aug 00:37:33:DEBUG:root:Assembling\n", + "13-Aug 00:37:33:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:33:DEBUG:root:Added 5 Lagrange multipliers to model\n", + "13-Aug 00:37:33:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:33:DEBUG:root:Solving system of equations. Total size = 16\n", + "13-Aug 00:37:33:DEBUG:root:Solution norm du = 5.027624635820698e-15\n", + "13-Aug 00:37:33:DEBUG:root:Converged in 5 iterations.\n", + "13-Aug 00:37:33:DEBUG:root:Displacement of element = \n", "[-0.02442313597467864 -0.039356000063335075 0.021097993207232584 0.020877031423993653\n", " -0.12673626841485705 -0.40433021969759375 -0.40656320177872923 -0.1275776940913048]\n" ] @@ -1105,7 +1136,7 @@ "\n", "WARNING: deprecated syntax \"{a=>b, ...}\" at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:32.\n", "Use \"Dict{Any,Any}(a=>b, ...)\" instead.\n", - "11-Aug 21:40:01:INFO:root:Registered handlers: Any[\"ELEMENT\",\"NODE\",\"NSET\"]\n", + "13-Aug 00:37:37:INFO:root:Registered handlers: Any[\"ELEMENT\",\"NODE\",\"NSET\"]\n", "WARNING: beginswith is deprecated, use startswith instead.\n", " in depwarn at /Applications/Julia-0.4.0-dev-539c818c4e.app/Contents/Resources/julia/lib/julia/sys.dylib\n", " in beginswith at deprecated.jl:30\n", @@ -1124,8 +1155,8 @@ " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[13], in expression starting on line 3\n", - "11-Aug 21:40:01:DEBUG:root:Found NODE section\n", - "11-Aug 21:40:01:DEBUG:root:Found ELEMENT section\n", + "13-Aug 00:37:38:DEBUG:root:Found NODE section\n", + "13-Aug 00:37:38:DEBUG:root:Found ELEMENT section\n", "WARNING: integer(s::AbstractString) is deprecated, use parse(Int,s) instead.\n", " in depwarn at /Applications/Julia-0.4.0-dev-539c818c4e.app/Contents/Resources/julia/lib/julia/sys.dylib\n", " in integer at deprecated.jl:49\n", @@ -1138,14 +1169,14 @@ " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[13], in expression starting on line 3\n", - "11-Aug 21:40:03:DEBUG:root:120 elements found\n", - "11-Aug 21:40:03:INFO:root:Creating ELSET Body1\n", - "11-Aug 21:40:03:DEBUG:root:Found NSET section\n", - "11-Aug 21:40:03:DEBUG:root:Creating node set SUPPORT\n", - "11-Aug 21:40:03:DEBUG:root:Found NSET section\n", - "11-Aug 21:40:03:DEBUG:root:Creating node set LOAD\n", - "11-Aug 21:40:03:DEBUG:root:Found NSET section\n", - "11-Aug 21:40:03:DEBUG:root:Creating node set TOP\n" + "13-Aug 00:37:39:DEBUG:root:120 elements found\n", + "13-Aug 00:37:40:INFO:root:Creating ELSET Body1\n", + "13-Aug 00:37:40:DEBUG:root:Found NSET section\n", + "13-Aug 00:37:40:DEBUG:root:Creating node set SUPPORT\n", + "13-Aug 00:37:40:DEBUG:root:Found NSET section\n", + "13-Aug 00:37:40:DEBUG:root:Creating node set LOAD\n", + "13-Aug 00:37:40:DEBUG:root:Found NSET section\n", + "13-Aug 00:37:40:DEBUG:root:Creating node set TOP\n" ] }, { @@ -1173,7 +1204,7 @@ }, { "cell_type": "code", - "execution_count": 34, + "execution_count": 14, "metadata": { "collapsed": false, "scrolled": false @@ -1183,47 +1214,47 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 22:15:29:DEBUG:root:Creating elements\n", - "11-Aug 22:15:29:INFO:root:create_ldof2gdofmap: dofs per node: 3\n", - "11-Aug 22:15:29:INFO:root:solve!: dofs per node: 3\n", - "11-Aug 22:15:29:DEBUG:root:Problem size = 894\n", - "11-Aug 22:15:29:DEBUG:root:Starting iteration 1\n", - "11-Aug 22:15:29:DEBUG:root:Assembling\n", - "11-Aug 22:15:30:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 22:15:30:DEBUG:root:Added 28 Lagrange multipliers to model\n", - "11-Aug 22:15:30:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 22:15:30:DEBUG:root:Solving system of equations. Total size = 921\n", - "11-Aug 22:15:30:DEBUG:root:Solution norm du = 72.87727091053921\n", - "11-Aug 22:15:30:DEBUG:root:Starting iteration 2\n", - "11-Aug 22:15:31:DEBUG:root:Assembling\n", - "11-Aug 22:15:32:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 22:15:32:DEBUG:root:Added 28 Lagrange multipliers to model\n", - "11-Aug 22:15:32:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 22:15:32:DEBUG:root:Solving system of equations. Total size = 921\n", - "11-Aug 22:15:32:DEBUG:root:Solution norm du = 2.737285770100854\n", - "11-Aug 22:15:32:DEBUG:root:Starting iteration 3\n", - "11-Aug 22:15:32:DEBUG:root:Assembling\n", - "11-Aug 22:15:33:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 22:15:34:DEBUG:root:Added 28 Lagrange multipliers to model\n", - "11-Aug 22:15:34:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 22:15:34:DEBUG:root:Solving system of equations. Total size = 921\n", - "11-Aug 22:15:34:DEBUG:root:Solution norm du = 0.07997112801214978\n", - "11-Aug 22:15:34:DEBUG:root:Starting iteration 4\n", - "11-Aug 22:15:34:DEBUG:root:Assembling\n", - "11-Aug 22:15:35:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 22:15:35:DEBUG:root:Added 28 Lagrange multipliers to model\n", - "11-Aug 22:15:35:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 22:15:35:DEBUG:root:Solving system of equations. Total size = 921\n", - "11-Aug 22:15:35:DEBUG:root:Solution norm du = 6.406557430235748e-5\n", - "11-Aug 22:15:35:DEBUG:root:Starting iteration 5\n", - "11-Aug 22:15:35:DEBUG:root:Assembling\n", - "11-Aug 22:15:37:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", - "11-Aug 22:15:37:DEBUG:root:Added 28 Lagrange multipliers to model\n", - "11-Aug 22:15:37:DEBUG:root:Adding Neumann boundary conditions\n", - "11-Aug 22:15:37:DEBUG:root:Solving system of equations. Total size = 921\n", - "11-Aug 22:15:37:DEBUG:root:Solution norm du = 6.870072007793185e-11\n", - "11-Aug 22:15:37:DEBUG:root:Converged in 5 iterations.\n", - "11-Aug 22:15:37:INFO:root:Maximum absolute displacement in y direction: 6.9925206227884695\n" + "13-Aug 00:37:41:DEBUG:root:Creating elements\n", + "13-Aug 00:37:41:INFO:root:create_ldof2gdofmap: dofs per node: 3\n", + "13-Aug 00:37:41:INFO:root:solve!: dofs per node: 3\n", + "13-Aug 00:37:41:DEBUG:root:Problem size = 894\n", + "13-Aug 00:37:41:DEBUG:root:Starting iteration 1\n", + "13-Aug 00:37:41:DEBUG:root:Assembling\n", + "13-Aug 00:37:43:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:43:DEBUG:root:Added 28 Lagrange multipliers to model\n", + "13-Aug 00:37:43:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:43:DEBUG:root:Solving system of equations. Total size = 921\n", + "13-Aug 00:37:43:DEBUG:root:Solution norm du = 54.19642700242575\n", + "13-Aug 00:37:43:DEBUG:root:Starting iteration 2\n", + "13-Aug 00:37:43:DEBUG:root:Assembling\n", + "13-Aug 00:37:45:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:45:DEBUG:root:Added 28 Lagrange multipliers to model\n", + "13-Aug 00:37:45:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:45:DEBUG:root:Solving system of equations. Total size = 921\n", + "13-Aug 00:37:45:DEBUG:root:Solution norm du = 1.68729144400063\n", + "13-Aug 00:37:45:DEBUG:root:Starting iteration 3\n", + "13-Aug 00:37:45:DEBUG:root:Assembling\n", + "13-Aug 00:37:47:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:47:DEBUG:root:Added 28 Lagrange multipliers to model\n", + "13-Aug 00:37:47:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:47:DEBUG:root:Solving system of equations. Total size = 921\n", + "13-Aug 00:37:47:DEBUG:root:Solution norm du = 0.04175080278208098\n", + "13-Aug 00:37:47:DEBUG:root:Starting iteration 4\n", + "13-Aug 00:37:47:DEBUG:root:Assembling\n", + "13-Aug 00:37:48:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:48:DEBUG:root:Added 28 Lagrange multipliers to model\n", + "13-Aug 00:37:48:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:48:DEBUG:root:Solving system of equations. Total size = 921\n", + "13-Aug 00:37:49:DEBUG:root:Solution norm du = 2.4176157836844963e-5\n", + "13-Aug 00:37:49:DEBUG:root:Starting iteration 5\n", + "13-Aug 00:37:49:DEBUG:root:Assembling\n", + "13-Aug 00:37:50:DEBUG:root:Adding Dirichlet boundary conditions using Lagrange multipliers\n", + "13-Aug 00:37:50:DEBUG:root:Added 28 Lagrange multipliers to model\n", + "13-Aug 00:37:50:DEBUG:root:Adding Neumann boundary conditions\n", + "13-Aug 00:37:50:DEBUG:root:Solving system of equations. Total size = 921\n", + "13-Aug 00:37:50:DEBUG:root:Solution norm du = 1.4448231751500831e-11\n", + "13-Aug 00:37:50:DEBUG:root:Converged in 5 iterations.\n", + "13-Aug 00:37:50:INFO:root:Maximum absolute displacement in y direction: 5.245400568184194\n" ] } ], @@ -1294,7 +1325,7 @@ }, { "cell_type": "code", - "execution_count": 49, + "execution_count": 15, "metadata": { "collapsed": false }, @@ -1307,7 +1338,7 @@ "\n" ] }, - "execution_count": 49, + "execution_count": 15, "metadata": {}, "output_type": "execute_result" } @@ -1327,7 +1358,7 @@ }, { "cell_type": "code", - "execution_count": 50, + "execution_count": 16, "metadata": { "collapsed": false }, @@ -1336,8 +1367,8 @@ "name": "stderr", "output_type": "stream", "text": [ - "11-Aug 22:38:28:INFO:root:Number of nodes in model: 298\n", - "11-Aug 22:38:28:INFO:root:Number of elements in model: 120\n" + "13-Aug 00:37:51:INFO:root:Number of nodes in model: 298\n", + "13-Aug 00:37:51:INFO:root:Number of elements in model: 120\n" ] } ], @@ -1365,7 +1396,7 @@ }, { "cell_type": "code", - "execution_count": 51, + "execution_count": 17, "metadata": { "collapsed": false }, @@ -1376,7 +1407,7 @@ "true" ] }, - "execution_count": 51, + "execution_count": 17, "metadata": {}, "output_type": "execute_result" } @@ -1429,7 +1460,7 @@ }, { "cell_type": "code", - "execution_count": 52, + "execution_count": 18, "metadata": { "collapsed": false }, @@ -1440,7 +1471,7 @@ "10435" ] }, - "execution_count": 52, + "execution_count": 18, "metadata": {}, "output_type": "execute_result" } @@ -1458,7 +1489,7 @@ }, { "cell_type": "code", - "execution_count": 53, + "execution_count": 19, "metadata": { "collapsed": false }, @@ -1467,35 +1498,35 @@ "data": { "text/plain": [ "Dict{Any,Any} with 298 entries:\n", - " 288 => [-0.21784897897079505,-6.64256844921929,0.18339325827190953]\n", - " 11 => [-0.17342410820059667,-6.343125806991628,0.17544367431766225]\n", - " 158 => [-0.0007723334905839533,-0.2321904799732306,-0.03883691791045858]\n", - " 215 => [-0.03299227384130259,-4.2236276356937505,-0.0017757488281830807]\n", - " 134 => [-0.06336705102703495,-3.3332586850719337,0.06555932581531249]\n", - " 160 => [0.03205303832663507,-0.536361456950779,-0.07291283371223893]\n", - " 29 => [-0.01620524134284992,-0.19455583841608232,-0.03838421619811468]\n", - " 131 => [0.04015396763870081,-3.5101776818531523,0.10045744211179368]\n", - " 249 => [-0.2491004521409732,-6.570813533283954,0.06837890174506163]\n", - " 207 => [0.03974407267311282,-3.48588780058722,0.007971846403163977]\n", - " 173 => [-0.04472757590353555,-4.582185697159958,-0.04434348606400347]\n", - " 289 => [-0.16494613158592017,-6.360857515886982,0.33435404817411213]\n", - " 74 => [-0.03449870789750817,-3.947284664599636,0.0692691027146262]\n", - " 201 => [-0.00128347996939122,-3.8527236749372054,0.054216444201382656]\n", - " 176 => [-0.004729440266412542,-3.132546258410676,0.0036604388995965386]\n", - " 57 => [-0.021453552968802154,-4.134921994863808,0.04888988010136468]\n", - " 31 => [-0.005045439296021616,-0.3424143474221696,-0.0365148544462663]\n", - " 285 => [-0.19371941654927557,-6.446502168659889,0.18022053743770028]\n", - " 70 => [0.006277500845425116,-3.802694347992476,0.06251589929767469]\n", - " 33 => [-0.02778798411401373,-0.5091869603980016,0.03538255963543149]\n", - " 252 => [-0.11915625386204923,-1.2182594899567352,-0.025975085489623087]\n", - " 114 => [-0.14643901527350212,-0.4628626545501056,-0.0057156433415159625]\n", - " 165 => [-0.07462160313837857,-4.426674712264929,-0.04391022455570823]\n", - " 96 => [0.0541064545893336,-5.057801826440659,-0.01669155772495867]\n", - " 133 => [0.01615115174975327,-3.140610327091724,0.028951801068022566]\n", + " 288 => [-0.11672879201601734,-4.945796660856531,0.34641446977739976]\n", + " 11 => [-0.10351195894967631,-4.756564623966246,0.3419475317395125]\n", + " 158 => [-0.004781967238459207,-0.15800625216143013,-0.03641750463465448]\n", + " 215 => [-0.009679264390442352,-3.1283115653077243,0.1060066992444004]\n", + " 134 => [-0.012153531002752209,-2.568061461984137,0.1358780416251414]\n", + " 160 => [0.045841800099574996,-0.3952996715282528,-0.04337647899598485]\n", + " 29 => [-0.019610173840903154,-0.13106872577808154,-0.04187494739380516]\n", + " 131 => [0.04911611452200464,-2.5861023637514955,0.15290704689462717]\n", + " 249 => [-0.18589774363886435,-4.932931637644884,0.27789363131184]\n", + " 207 => [0.05588675730587902,-2.563471323701039,0.0669746466494972]\n", + " 173 => [-0.017826777186855432,-3.4323824666471117,0.08923268461041232]\n", + " 289 => [-0.08168560138254058,-4.782381469728047,0.4566974667444831]\n", + " 74 => [-0.013759506382279562,-2.9054222587693794,0.1389413161160111]\n", + " 201 => [0.006183182231407451,-2.8511450493652664,0.14659955727801022]\n", + " 176 => [0.04229935999959756,-2.349387641785469,0.09405299762367403]\n", + " 57 => [-0.0022029054234204426,-3.062804411311773,0.14846103116563106]\n", + " 31 => [0.005313697788644969,-0.25237461591643895,-0.030307540603490793]\n", + " 285 => [-0.10859185872935265,-4.832314881155342,0.35832378735042353]\n", + " 70 => [0.042811232565748765,-2.7950866978216338,0.1421452814276251]\n", + " 33 => [-0.016450709743293015,-0.379272230318473,0.027549328283933208]\n", + " 252 => [-0.14470395638682446,-0.9284887677850233,0.043594035381458736]\n", + " 114 => [-0.14842327883584733,-0.3532980011491899,-0.018641651848308006]\n", + " 165 => [-0.04346921888336824,-3.2745129973662452,0.07316785266950793]\n", + " 96 => [0.1501458469524447,-3.8471140528652787,0.1489845291162686]\n", + " 133 => [0.09045793517058681,-2.3950191788705957,0.12799724500423504]\n", " ⋮ => ⋮" ] }, - "execution_count": 53, + "execution_count": 19, "metadata": {}, "output_type": "execute_result" } @@ -1512,7 +1543,7 @@ }, { "cell_type": "code", - "execution_count": 54, + "execution_count": 20, "metadata": { "collapsed": false }, @@ -1521,12 +1552,12 @@ "data": { "text/plain": [ "3x298 Array{Float64,2}:\n", - " -0.202475 -0.17437 0.0131267 -0.0340868 -0.0252162 0.0521758 … -0.054781 -0.169903 0.186936 0.102666 -0.0970747 -0.164397 0.0\n", - " -1.34204 -1.66404 -3.20587 -3.27673 -3.68732 -3.36789 -4.33822 -1.03943 -2.27466 -1.67467 -2.91213 -2.52353 0.0\n", - " 0.0211856 0.0324162 0.00213878 0.0377156 0.0683502 0.030081 -0.00684333 -0.035509 0.304867 0.26033 0.113802 0.152433 0.0" + " -0.174259 -0.175081 0.0359566 0.0279021 0.00925611 0.0379208 … -0.0248959 -0.116472 0.191156 0.193018 -0.0153311 -0.123282 0.0\n", + " -0.978534 -1.32327 -2.42646 -2.50175 -2.72095 -2.59357 -3.20696 -0.707364 -1.76025 -1.33439 -2.22191 -1.93398 0.0\n", + " 0.0535874 0.111758 0.125464 0.115461 0.145374 0.1251 0.0902774 0.00418633 0.368958 0.387435 0.185441 0.230546 0.0" ] }, - "execution_count": 54, + "execution_count": 20, "metadata": {}, "output_type": "execute_result" } @@ -1541,7 +1572,7 @@ }, { "cell_type": "code", - "execution_count": 55, + "execution_count": 21, "metadata": { "collapsed": false }, @@ -1552,7 +1583,7 @@ "true" ] }, - "execution_count": 55, + "execution_count": 21, "metadata": {}, "output_type": "execute_result" } @@ -1563,7 +1594,7 @@ }, { "cell_type": "code", - "execution_count": 56, + "execution_count": 22, "metadata": { "collapsed": false }, @@ -1571,10 +1602,10 @@ { "data": { "text/plain": [ - "28330" + "28215" ] }, - "execution_count": 56, + "execution_count": 22, "metadata": {}, "output_type": "execute_result" } diff --git a/notebooks/2015-06-25-shape-functions.ipynb b/notebooks/2015-06-25-shape-functions.ipynb index 6b87a27..3113e17 100644 --- a/notebooks/2015-06-25-shape-functions.ipynb +++ b/notebooks/2015-06-25-shape-functions.ipynb @@ -1,5 +1,16 @@ { "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Shape function and integration points\n", + "\n", + "Author(s): Jukka Aho\n", + "\n", + "**Abstract**: Shape functions and element descriptions used in JuliaFEM." + ] + }, { "cell_type": "code", "execution_count": 1, @@ -8,18 +19,325 @@ }, "outputs": [], "source": [ - "from sympy import *" + "from sympy import *\n", + "#init_printing()" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { - "collapsed": true + "collapsed": false }, "outputs": [], "source": [ - "xi = DeferredVector(\"xi\")" + "xi = DeferredVector(r\"xi\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 1D shape function" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Linear 2-node segment (Lagrange family)\n", + "\n", + "| | $\\xi_1$ |\n", + "| ----- | -------:|\n", + "| $N_1$ | -1 |\n", + "| $N_2$ | 1 |\n", + "\n", + "\\begin{equation}\n", + " \\left(\\mathbf{P}\\boldsymbol{\\alpha}\\right)\\left(\\xi_1\\right) = \\alpha_{1}+\\alpha_{2}\\xi_{1}\n", + "\\end{equation}" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "(Matrix([\n", + " [-xi[1]/2 + 1/2],\n", + " [ xi[1]/2 + 1/2]]), Matrix([\n", + " [-1/2],\n", + " [ 1/2]]))" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "A = Matrix([[1, -1], [1, 1]])\n", + "P = Matrix([1, xi[1]]).T\n", + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T]).T\n", + "N, dN" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Quadratic 3-node segment (Lagrange family)\n", + "\n", + "| | $\\xi_1$ |\n", + "| ----- | -------:|\n", + "| $N_1$ | -1 |\n", + "| $N_2$ | 1 |\n", + "| $N_3$ | 0 |\n", + "\n", + "\\begin{equation}\n", + " \\left(\\mathbf{P}\\boldsymbol{\\alpha}\\right)\\left(\\xi_1\\right) = \\alpha_1 + \\alpha_2\\xi_1 + \\alpha_3\\xi_1^2\n", + "\\end{equation}" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "(Matrix([\n", + " [xi[1]**2/2 - xi[1]/2],\n", + " [xi[1]**2/2 + xi[1]/2],\n", + " [ -xi[1]**2 + 1]]), Matrix([\n", + " [xi[1] - 1/2],\n", + " [xi[1] + 1/2],\n", + " [ -2*xi[1]]]))" + ] + }, + "execution_count": 4, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "A = Matrix([[1, -1, (-1)**2],\n", + " [1, 1, 1**2],\n", + " [1, 0, 0**2]])\n", + "P = Matrix([1, xi[1], xi[1]**2]).T\n", + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T]).T\n", + "N, dN" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### P-elements" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 2D shape functions" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Linear triangle\n", + "\n", + "| | $\\xi_1$ | $\\xi_2$ |\n", + "| ----- | -------:| -------:|\n", + "| $N_1$ | 0 | 0 |\n", + "| $N_2$ | 1 | 0 |\n", + "| $N_3$ | 0 | 1 |" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "(Matrix([\n", + " [-xi[1] - xi[2] + 1],\n", + " [ xi[1]],\n", + " [ xi[2]]]), Matrix([\n", + " [-1, -1],\n", + " [ 1, 0],\n", + " [ 0, 1]]))" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "A = Matrix([[1, 0, 0], [1, 1, 0], [1, 0, 1]])\n", + "P = Matrix([1, xi[1], xi[2]]).T\n", + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T, N.diff(xi[2]).T]).T\n", + "N, dN" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Quadratic triangle\n", + "\n", + "| | $\\xi_1$ | $\\xi_2$ |\n", + "| ----- | -------:| -------:|\n", + "| $N_1$ | 0 | 0 |\n", + "| $N_2$ | 1 | 0 |\n", + "| $N_3$ | 0 | 1 |\n", + "| $N_4$ | 1/2 | 0 |\n", + "| $N_5$ | 1/2 | 1/2 |\n", + "| $N_6$ | 0 | 1/2 |" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "Matrix([\n", + "[1, 0, 0, 0, 0, 0],\n", + "[1, 1, 0, 1, 0, 0],\n", + "[1, 0, 1, 0, 1, 0],\n", + "[1, 1/2, 0, 1/4, 0, 0],\n", + "[1, 1/2, 1/2, 1/4, 1/4, 1/4],\n", + "[1, 0, 1/2, 0, 1/4, 0]])" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "P = Matrix([1, xi[1], xi[2], xi[1]**2, xi[2]**2, xi[1]*xi[2]]).T\n", + "A = Matrix([\n", + " P.subs({xi[1]: 0, xi[2]: 0}),\n", + " P.subs({xi[1]: 1, xi[2]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: 1}),\n", + " P.subs({xi[1]: Rational(1,2), xi[2]: 0}),\n", + " P.subs({xi[1]: Rational(1,2), xi[2]: Rational(1,2)}),\n", + " P.subs({xi[1]: 0, xi[2]: Rational(1,2)}),\n", + " ])\n", + "A" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "(Matrix([\n", + " [2*xi[1]**2 + 4*xi[1]*xi[2] - 3*xi[1] + 2*xi[2]**2 - 3*xi[2] + 1],\n", + " [ 2*xi[1]**2 - xi[1]],\n", + " [ 2*xi[2]**2 - xi[2]],\n", + " [ -4*xi[1]**2 - 4*xi[1]*xi[2] + 4*xi[1]],\n", + " [ 4*xi[1]*xi[2]],\n", + " [ -4*xi[1]*xi[2] - 4*xi[2]**2 + 4*xi[2]]]), Matrix([\n", + " [ 4*xi[1] + 4*xi[2] - 3, 4*xi[1] + 4*xi[2] - 3],\n", + " [ 4*xi[1] - 1, 0],\n", + " [ 0, 4*xi[2] - 1],\n", + " [-8*xi[1] - 4*xi[2] + 4, -4*xi[1]],\n", + " [ 4*xi[2], 4*xi[1]],\n", + " [ -4*xi[2], -4*xi[1] - 8*xi[2] + 4]]))" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T, N.diff(xi[2]).T]).T\n", + "N, dN" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 3D shape functions" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Linear tetrahedra, **tet4**\n", + "\n", + "| | $\\xi_1$ | $\\xi_2$ | $\\xi_2$ |\n", + "| ----- | -------:| -------:| -------:|\n", + "| $N_1$ | 0 | 0 | 0 |\n", + "| $N_2$ | 1 | 0 | 0 |\n", + "| $N_3$ | 0 | 1 | 0 |\n", + "| $N_4$ | 0 | 0 | 1 |" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "Matrix([\n", + "[1, 0, 0, 0],\n", + "[1, 1, 0, 0],\n", + "[1, 0, 1, 0],\n", + "[1, 0, 0, 1]])" + ] + }, + "execution_count": 8, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "P = Matrix([1, xi[1], xi[2], xi[3]]).T\n", + "A = Matrix([\n", + " P.subs({xi[1]: 0, xi[2]: 0, xi[3]: 0}),\n", + " P.subs({xi[1]: 1, xi[2]: 0, xi[3]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: 1, xi[3]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: 0, xi[3]: 1}),\n", + " ])\n", + "A" ] }, { @@ -32,17 +350,15 @@ { "data": { "text/plain": [ - "Matrix([\n", - "[ -xi[1] - xi[2] - xi[3] + 1],\n", - "[ xi[1]],\n", - "[ xi[2]],\n", - "[ xi[3]],\n", - "[4*xi[1]*(-xi[1] - xi[2] - xi[3] + 1)],\n", - "[ 4*xi[1]*xi[2]],\n", - "[4*xi[2]*(-xi[1] - xi[2] - xi[3] + 1)],\n", - "[4*xi[3]*(-xi[1] - xi[2] - xi[3] + 1)],\n", - "[ 4*xi[1]*xi[3]],\n", - "[ 4*xi[2]*xi[3]]])" + "(Matrix([\n", + " [-xi[1] - xi[2] - xi[3] + 1],\n", + " [ xi[1]],\n", + " [ xi[2]],\n", + " [ xi[3]]]), Matrix([\n", + " [-1, -1, -1],\n", + " [ 1, 0, 0],\n", + " [ 0, 1, 0],\n", + " [ 0, 0, 1]]))" ] }, "execution_count": 9, @@ -51,23 +367,29 @@ } ], "source": [ - "def c3d10():\n", - " N1 = 1 - xi[1] - xi[2] - xi[3]\n", - " N2 = xi[1]\n", - " N3 = xi[2]\n", - " N4 = xi[3]\n", - " N5 = 4*xi[1]*(1-xi[1]-xi[2]-xi[3])\n", - " N6 = 4*xi[1]*xi[2]\n", - " N7 = 4*xi[2]*(1-xi[1]-xi[2]-xi[3])\n", - " N8 = 4*xi[3]*(1-xi[1]-xi[2]-xi[3])\n", - " N9 = 4*xi[1]*xi[3]\n", - " N10 = 4*xi[2]*xi[3]\n", - " N = Matrix([N1, N2, N3, N4, N5, N6, N7, N8, N9, N10])\n", - " dN = Matrix([N.diff(xi[1]).T, N.diff(xi[2]).T, N.diff(xi[3]).T]).T\n", - " return N, dN\n", + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T, N.diff(xi[2]).T, N.diff(xi[3]).T]).T\n", + "N, dN" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Quadratic Lagrange tetrahedral element, 10 nodes, **tet10**\n", "\n", - "N, dN = c3d10()\n", - "N" + "| | $\\xi_1$ | $\\xi_2$ | $\\xi_2$ |\n", + "| ----- | -------:| -------:| -------:|\n", + "| $N_1$ | 0 | 0 | 0 |\n", + "| $N_2$ | 1 | 0 | 0 |\n", + "| $N_3$ | 0 | 1 | 0 |\n", + "| $N_4$ | 0 | 0 | 1 |\n", + "| $N_5$ | 1/2 | 0 | 0 |\n", + "| $N_6$ | 1/2 | 1/2 | 0 |\n", + "| $N_7$ | 0 | 1/2 | 0 |\n", + "| $N_8$ | 0 | 0 | 1/2 |\n", + "| $N_9$ | 1/2 | 0 | 1/2 |\n", + "| $N_{10}$ | 0 | 1/2 | 1/2 |" ] }, { @@ -81,16 +403,16 @@ "data": { "text/plain": [ "Matrix([\n", - "[ -1, -1, -1],\n", - "[ 1, 0, 0],\n", - "[ 0, 1, 0],\n", - "[ 0, 0, 1],\n", - "[-8*xi[1] - 4*xi[2] - 4*xi[3] + 4, -4*xi[1], -4*xi[1]],\n", - "[ 4*xi[2], 4*xi[1], 0],\n", - "[ -4*xi[2], -4*xi[1] - 8*xi[2] - 4*xi[3] + 4, -4*xi[2]],\n", - "[ -4*xi[3], -4*xi[3], -4*xi[1] - 4*xi[2] - 8*xi[3] + 4],\n", - "[ 4*xi[3], 0, 4*xi[1]],\n", - "[ 0, 4*xi[3], 4*xi[2]]])" + "[1, 0, 0, 0, 0, 0, 0, 0, 0, 0],\n", + "[1, 1, 0, 0, 1, 0, 0, 0, 0, 0],\n", + "[1, 0, 1, 0, 0, 0, 1, 0, 0, 0],\n", + "[1, 0, 0, 1, 0, 0, 0, 0, 1, 0],\n", + "[1, 1/2, 0, 0, 1/4, 0, 0, 0, 0, 0],\n", + "[1, 1/2, 1/2, 0, 1/4, 1/4, 1/4, 0, 0, 0],\n", + "[1, 0, 1/2, 0, 0, 0, 1/4, 0, 0, 0],\n", + "[1, 0, 0, 1/2, 0, 0, 0, 0, 1/4, 0],\n", + "[1, 1/2, 0, 1/2, 1/4, 0, 0, 0, 1/4, 1/4],\n", + "[1, 0, 1/2, 1/2, 0, 0, 1/4, 1/4, 1/4, 0]])" ] }, "execution_count": 10, @@ -99,17 +421,89 @@ } ], "source": [ - "dN" + "P = Matrix([1, xi[1], xi[2], xi[3], xi[1]**2, xi[1]*xi[2], xi[2]**2, xi[2]*xi[3], xi[3]**2, xi[1]*xi[3]]).T\n", + "A = Matrix([\n", + " P.subs({xi[1]: 0, xi[2]: 0, xi[3]: 0}),\n", + " P.subs({xi[1]: 1, xi[2]: 0, xi[3]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: 1, xi[3]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: 0, xi[3]: 1}),\n", + "\n", + " P.subs({xi[1]: Rational(1,2), xi[2]: 0, xi[3]: 0}),\n", + " P.subs({xi[1]: Rational(1,2), xi[2]: Rational(1,2), xi[3]: 0}),\n", + " P.subs({xi[1]: 0, xi[2]: Rational(1,2), xi[3]: 0}),\n", + "\n", + " P.subs({xi[1]: 0, xi[2]: 0, xi[3]: Rational(1,2)}),\n", + " P.subs({xi[1]: Rational(1,2), xi[2]: 0, xi[3]: Rational(1,2)}),\n", + " P.subs({xi[1]: 0, xi[2]: Rational(1,2), xi[3]: Rational(1,2)}),\n", + " ])\n", + "A" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 11, "metadata": { - "collapsed": true + "collapsed": false }, - "outputs": [], - "source": [] + "outputs": [ + { + "data": { + "text/plain": [ + "Matrix([\n", + "[(xi[1] + xi[2] + xi[3] - 1)*(2*xi[1] + 2*xi[2] + 2*xi[3] - 1)],\n", + "[ xi[1]*(2*xi[1] - 1)],\n", + "[ xi[2]*(2*xi[2] - 1)],\n", + "[ xi[3]*(2*xi[3] - 1)],\n", + "[ -4*xi[1]*(xi[1] + xi[2] + xi[3] - 1)],\n", + "[ 4*xi[1]*xi[2]],\n", + "[ -4*xi[2]*(xi[1] + xi[2] + xi[3] - 1)],\n", + "[ -4*xi[3]*(xi[1] + xi[2] + xi[3] - 1)],\n", + "[ 4*xi[1]*xi[3]],\n", + "[ 4*xi[2]*xi[3]]])" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "N = (P*A.inv()).T\n", + "dN = Matrix([N.diff(xi[1]).T, N.diff(xi[2]).T, N.diff(xi[3]).T]).T\n", + "factor(N)" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": { + "collapsed": false + }, + "outputs": [ + { + "data": { + "text/plain": [ + "Matrix([\n", + "[ 4*xi[1] + 4*xi[2] + 4*xi[3] - 3, 4*xi[1] + 4*xi[2] + 4*xi[3] - 3, 4*xi[1] + 4*xi[2] + 4*xi[3] - 3],\n", + "[ 4*xi[1] - 1, 0, 0],\n", + "[ 0, 4*xi[2] - 1, 0],\n", + "[ 0, 0, 4*xi[3] - 1],\n", + "[-4*(2*xi[1] + xi[2] + xi[3] - 1), -4*xi[1], -4*xi[1]],\n", + "[ 4*xi[2], 4*xi[1], 0],\n", + "[ -4*xi[2], -4*(xi[1] + 2*xi[2] + xi[3] - 1), -4*xi[2]],\n", + "[ -4*xi[3], -4*xi[3], -4*(xi[1] + xi[2] + 2*xi[3] - 1)],\n", + "[ 4*xi[3], 0, 4*xi[1]],\n", + "[ 0, 4*xi[3], 4*xi[2]]])" + ] + }, + "execution_count": 12, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "factor(dN)" + ] } ], "metadata": { @@ -128,7 +522,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", - "version": "2.7.9" + "version": "2.7.10" } }, "nbformat": 4,