{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Elasticity solver examples\n", "\n", "Author(s): Jukka Aho " ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Docile: upgrading cache to 0.0.2.\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\n", "WARNING: deprecated syntax \"{a=>b, ...}\" at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:11.\n", "Use \"Dict{Any,Any}(a=>b, ...)\" instead.\n", "\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", "Warning: requiring \"JuliaFEM\" did not define a corresponding module.\n", "Warning: requiring \"JuliaFEM\" did not define a corresponding module.\n", "Warning: requiring \"JuliaFEM\" did not define a corresponding module.\n" ] }, { "data": { "text/plain": [ "Logger(root,DEBUG,Pipe(open, 0 bytes waiting),root)" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# These are internal module functions and not intended to use like this.\n", "using JuliaFEM.elasticity_solver\n", "using JuliaFEM.xdmf\n", "using JuliaFEM.abaqus_reader\n", "using Logging\n", "Logging.configure(level=DEBUG)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [], "source": [ "using LightXML" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 2d beam with linear elements" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "Array{Float64,2}" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "typeof(Float64[1 2; 2 3])" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Converged\n" ] }, { "data": { "text/plain": [ "2x4 Array{Float64,2}:\n", " 0.0 -0.399145 -0.0722858 0.0\n", " 0.0 -2.17799 -2.22224 0.0" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "\"\"\"\n", "#### Parameters\n", "X : Array{Float64, 2}\n", " \n", "\"\"\"\n", "function solve_elasticity_increment!(X, u, du, elmap, nodalloads,\n", " dirichletbc, lambda, mu, N, dNdchi, ipoints,\n", " iweights)\n", " if length(size(elmap)) == 1\n", " # quick hack for just one element\n", " elmap = elmap''\n", " end\n", " nelnodes, nelements = size(elmap)\n", " dim, nnodes = size(u)\n", " dofs = dim*nelnodes\n", "\n", " Imat = Int64[]\n", " Jmat = Int64[]\n", " Vmat = Float64[]\n", " Ivec = Int64[]\n", " Vvec = Float64[]\n", "\n", " # FIXME: different number of nodes/element\n", " R = zeros(dim, nelnodes)\n", " Kt = zeros(dofs, dofs)\n", "\n", " # this can be parallelized\n", " for i in 1:nelements\n", " eldofs = elmap[:,i]\n", " calc_local_matrices!(X[:, eldofs], u[:, eldofs], R, Kt, N, dNdchi,\n", " lambda[eldofs], mu[eldofs], ipoints, iweights)\n", " assemble!(Kt, eldofs, Imat, Jmat, Vmat)\n", " assemble!(R, eldofs, Ivec, Vvec)\n", " end\n", "\n", " # add additional neumann boundary conditions to force vector\n", " for (i, nodal_load) in enumerate(nodalloads)\n", " if nodal_load == 0\n", " continue\n", " end\n", " push!(Ivec, i)\n", " push!(Vvec, -nodal_load)\n", " end\n", "\n", " # Create sparse matrix and vector\n", " A = sparse(Imat, Jmat, Vmat)\n", " b = sparsevec(Ivec, Vvec)\n", "\n", " # Remove dirichlet boundary conditions\n", " free_dofs = find(isnan(dirichletbc))\n", " #Imat, Jmat, Vmat = eliminate_boundary_conditions(dirichletbc, Imat, Jmat, Vmat)\n", " b = b[free_dofs]\n", " A = A[free_dofs, free_dofs]\n", "\n", " # solution\n", " du[free_dofs] = lufact(A) \\ -full(b)\n", "end\n", "\n", "function one_elem_fixture()\n", " X = [0.0 0.0; 10.0 0.0; 10.0 1.0; 0.0 1.0]'\n", " elmap = [1; 2; 3; 4]\n", " nodalloads = [0 0; 0 0; 0 -2; 0 0]'\n", " dirichletbc = [0 0; NaN NaN; NaN NaN; 0 0]'\n", "\n", " E = 90\n", " nu = 0.25\n", " mu = E/(2*(1+nu))\n", " la = E*nu/((1+nu)*(1-2*nu))\n", " la = 2*la*mu/(la + 2*mu)\n", "\n", " la = la*ones(1, 4)\n", " mu = mu*ones(1, 4)\n", " u = zeros(2, 4)\n", " du = zeros(2, 4)\n", "\n", " N(xi) = [\n", " (1-xi[1])*(1-xi[2])/4\n", " (1+xi[1])*(1-xi[2])/4\n", " (1+xi[1])*(1+xi[2])/4\n", " (1-xi[1])*(1+xi[2])/4\n", " ]\n", "\n", " dNdξ(ξ) = [-(1-ξ[2])/4.0 -(1-ξ[1])/4.0\n", " (1-ξ[2])/4.0 -(1+ξ[1])/4.0\n", " (1+ξ[2])/4.0 (1+ξ[1])/4.0\n", " -(1+ξ[2])/4.0 (1-ξ[1])/4.0]\n", "\n", " ipoints = 1/sqrt(3)*[-1 -1; 1 -1; 1 1; -1 1]\n", " iweights = [1 1 1 1]\n", "\n", " return (X, u, du, elmap, nodalloads, dirichletbc,\n", " la, mu, N, dNdξ, ipoints, iweights)\n", "end\n", "\n", "#(X, u, du, elmap, nodalloads, dirichletbc,\n", "# la, mu, N, dNdξ, ipoints, iweights) = one_elem_fixture()\n", "function solve_one_element()\n", "\n", " m = new_model()\n", "\n", " # create nodes separately ...\n", " n1 = new_node()\n", " set_node_id(n1, 1)\n", " set_node_coords(n1, [0.0, 0.0])\n", "\n", " n2 = new_node()\n", " set_node_id(n2, 2)\n", " set_node_coords(n2, [10.0, 0.0])\n", " \n", " n3 = new_node()\n", " set_node_id(n3, 3)\n", " set_node_coords(n3, [10.0, 1.0])\n", "\n", " n4 = new_node()\n", " set_node_id(n4, 4)\n", " set_node_coords(n4, [0.0, 1.0])\n", "\n", " nodes = [n1, n2, n3, n4]\n", " add_nodes(m, \"NALL\", nodes)\n", "\n", " # .. or create somewhat simpler syntax\n", " # nodes = Dict(1 => [0.0, 0.0], 2 => [10.0, 0.0], 3 => [10.0, 1.0], 4 => [0.0, 1.0])\n", " # add_nodes(m, \"NALL\", nodes)\n", "\n", " # create element (hard way)\n", " e = new_element()\n", " set_element_id(1)\n", " set_node_ids(e, [1, 2, 3, 4])\n", "\n", " # we can set function spaces and integration schema for element-wise ...\n", " set_function_space(e, \"Lagrange(1)\") # use linear Lagrange function space to approximate unknown field\n", " set_integration_schema(e, \"FPG4\") # use four integration points\n", " # or for model as a \"default value\"\n", " # set_function_space(m, \"Lagrange\", 1)\n", " # set_integration_schema(m, \"FPG4\")\n", "\n", " # set necessary material parameters\n", " E = 90\n", " nu = 0.25\n", " mu = E/(2*(1+nu))\n", " la = E*nu/((1+nu)*(1-2*nu))\n", " la = 2*la*mu/(la + 2*mu)\n", "\n", " # we can add attribute for each element node ...\n", " #add_attribute(e, \"lambda\", [la, la, la, la])\n", " #add_attribute(e, \"mu\", [mu, mu, mu, mu])\n", " # or just one for each element ...\n", " add_attribute(e, \"lambda\", lambda)\n", " add_attribute(e, \"mu\", mu)\n", " # .. or just for model as a default value\n", " # add_attribute(m, \"lambda\", lambda)\n", " # add_attribute(m, \"mu\", mu)\n", " # note that we don't assign attributes to nodes here so we can describe discontinous fields\n", " # in attributes\n", "\n", " # or we can use convenient syntax\n", " # e = new_element(element_id=1, node_ids=[1, 2, 3, 4], function_space=\"Lagrange(1)\",\n", " # integration_schema=\"FPG4\", attributes=Dict(\"lambda\" => lambda, \"mu\" => mu))\n", "\n", " elements = [e]\n", " add_elements(m, \"EALL\", elements)\n", "\n", " # boundary conditions are always assignet to sets\n", "\n", " # element boundary conditions\n", " \n", " # add load for element surface S1 in normal-tangential coordinate system\n", " # add_attribute(e, \"displacement S1 normal load\", 1)\n", " # add load for element surface S1 from -1 to -2\n", " # add_attribute(e, \"displacement S1 load\", [-1, -2])\n", "\n", " # nodal boundary conditions\n", " # add neumann boundary condition to node 3, in y direction\n", " add_attribute(n3, \"displacement 2 load\", -2)\n", " # add dirichlet boundary condition to nodes 1 and 2 (encastre)\n", " add_attribute(n1, \"displacement\", [0.0, 0.0])\n", " # or\n", " # add_attribute(n1, \"displacement 1\", 0.0)\n", " # add_attribute(n1, \"displacement 2\", 0.0)\n", " add_attribute(n2, \"displacement\", [0.0, 0.0])\n", " # or, for example, add dirichlet boundary conditions in directions 1 and 3 for node\n", " #add_attribute(n2, \"displacement 1,3\", [0.0, 1.0])\n", " # or fix all dofs of a nodes in nodeset \"SUPPORT\"\n", " # support_bc = add_nodeset(m, \"SUPPORT\", [n1, n2])\n", " # add_attribute(support_bc, \"displacement\", 0)\n", "\n", "\n", "\n", " # Everything is very general so far. Datamodel is well defined and in this point\n", " # we can save or load it to disk\n", "\n", " # save_model(m, \"mymodel\") # save model to xml/h5 (Xdmf)\n", " # m = load_model(\"mymodel\") # load model from file\n", "\n", " \n", " # END OF MODEL DEFINITON\n", " \n", " # Now we kick in elasticity iterations and solve displacement\n", " # field but of course it could be something else too\n", "\n", "\n", " # m1, m2 = make_domain_decomposition(parts=2, keep_in_one_domain=[\"CONTACT_BOUNDARY\"])\n", "\n", "\n", " \n", " # In newton iteration, there might be situations where we just want to update RHS\n", " # like in radiation problems, it makes no sense to update stiffness matrix in that\n", " # case. For this reason local element matrix and force vector can be updated both\n", " # or just one of them. This time we don't have anything nonlinear in RHS so we can\n", " # calculate it outside of iteration loop\n", " p = new_problem(\"displacement\")\n", " \n", " RHS = zeros(length(nodes)*dofs)\n", " \n", " for i=1:10\n", " \n", " solve_elasticity_increment!(X, u, du, elmap, nodalloads, dirichletbc,\n", " la, mu, N, dNdξ, ipoints, iweights)\n", " u += du\n", " if norm(du) < 1.0e-9\n", " println(\"Converged\")\n", " break\n", " end\n", " end\n", " return u\n", "end\n", "\n", "u" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [], "source": [ "u3d = [u; 0 0 0 0] # extend to 3d vector field\n", "X3d = [X; 0 0 0 0]\n", "elmap2 = [0x5; elmap]'';" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", "\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "WARNING: int(x) is deprecated, use Int(x) instead.\n" ] }, { "data": { "text/plain": [ "1015" ] }, "execution_count": 4, "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 int at deprecated.jl:49\n", " in save_file at /Users/jukka/.julia/v0.4/LightXML/src/document.jl:108\n", " in xdmf_save_model at /Users/jukka/.julia/v0.4/JuliaFEM/src/xdmf.jl:106\n", " in include_string at loading.jl:99\n", " in execute_request_0x535c5df2 at /Users/jukka/.julia/v0.4/IJulia/src/execute_request.jl:157\n", " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[4], in expression starting on line 7\n" ] } ], "source": [ "xdoc, model = JuliaFEM.xdmf.xdmf_new_model()\n", "temporal_collection = JuliaFEM.xdmf.xdmf_new_temporal_collection(model)\n", "grid = JuliaFEM.xdmf.xdmf_new_grid(temporal_collection; time=0)\n", "JuliaFEM.xdmf.xdmf_new_mesh(grid, X3d, elmap2)\n", "JuliaFEM.xdmf.xdmf_new_field(grid, \"Displacement\", \"nodes\", u3d)\n", "print(xdoc)\n", "JuliaFEM.xdmf.xdmf_save_model(xdoc, \"/tmp/foo.xmf\")" ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": [ "iVBORw0KGgoAAAANSUhEUgAAA5oAAAHhCAIAAACA7Vb3AAAgAElEQVR4Xu3dfaxteVkf8GefamxtjQHFEiiCgFNAcBheh5m5d15AUBRfqCIaVDIIWKqkI4oYLA4KAQShgKADDDPMDPNyQRBEyKi8pS1pmqaNbZqmjWnSpklbkzY1sa1G5+z+sfd62WvttX6/tc7Z5+zfOZ/PH961137Wb6+958r53uc8e63FcrkMAAAo00GqAAAA9pc4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLcDyWsUiVAHD8xFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCfVWqgC0Wi5tXG8vlzWN1AADsmDg707df/fXRyrUh2gIAnAZx9kguv+ZBq43FwaIdbUO6BQA4EeLssXnyxQdHRMRi9VDjFgDgBIizu3LFdQ+JiMUiYjPahnQLAHB8xNmT8JTrH3JwUF1EYrHQuAUAOC7i7Cl42nMeuliPJGjcAgAciTh7yp7xvIetNhaLRZi4BQCYSJzdL8/8zoevc61LJQAAZBBn99ezvvsRUXVtQ+MWAGAbcbYYV7/gkRERB92ZhJBuAYBzTJwt0tXf96iDumt74FIJAMD5Jc6eBRdf+OiIWITGLQBw7oizZ821P/iY1cbqSrcatwDA2SbOnnHXv+hbV13bONC4BQDOIHH2HLnhxZeFSyUAAGeLOHt+fcePPm614Rq3AEC5xFkiIp77ksdH07V1qQQAoBjiLFs878e/LSIOete4FW0BgH0jzpLw/Jc+cb1lJgEA2D/iLBM8/8YnRcSBL5MBAHtDnGW+F/zk5RERBxFu3wAAnBJxluPxgldcfuDLZADAiRNn2Ynv/6krtl7gNqRbAOBYibPs3Atf9ZTVxsKlEgCA4ybOctL+3k8/ddW37cwkhHQLAEwnznKafujVTwv33QUAjkCcZY/88E1Pj4hFdGcSQroFAAaIs+ypF9/0jMXBetulEgCAIeIsZfjRn3tmVHdw0LgFAGriLOV5yWuvXG1UXym7uX5KtAWA80acpXg/9gtXbr3vbki3AHAOiLOcKT/xi8+K6gK3oXELAOeAOMtZduPrr1ptLBaiLQCcTeIs58XLfumq1qUSbm4/Jd0CQLnEWc6jl7/h6tXGKuAeQ+O2Gt7lXFv6m0BERCyXqQrgOImzEK+4+erVnckWsxu3fnoREbHwNwHg5ImzsOGn3nhNRBxsG0vIjbYAwAkSZ2HMq950IbZd4DakWwDYD+Is5PoHb75QNW1jcaBxCwB7QZyFmV79lgsRERq3AHCqxFk4Bv/wrRdWG+67CwAnTJyF43fT2y5UlwBbaNwCwE6Js7Bbr3nHxViPJERo3ALAcRNn4UT9/K9fXG107rsb0i0AzCLOwql57TsvHjRtW41bAJhDnIV98bp/fDFc4xYAJhJnYR+9/j3VTMLq/2rcAsAAcRYK8EvvuVh/m0y0BYA2cRYK84bfMJMAAA1xFgr2y++7uFgsV9sLjVsAziVxFs6ON77/wmLbTEJItwCcXeIsnE2/+lvuuwvAuSDOwrnwplsurPq2B27fAMDZIs7CufOWD1yICJdKAOBsEGfhvHvrB6+JiAOXSgCgTOIs0Pi1D13TdG1jqXELwP4TZ4FBb7/16nCNWwD2mzgLZHnnbVe1H2rcArAnxFlgjnfedpVr3AKwD8RZ4KjefftVEbGI6v5kGrcAnCBxFjhm7/nIs1Ybi8VS4xaAXRNngR167x1XHtQPFhq3ABw/cRY4Oe+788qoxhJEWwCOhTgLnI733/nMVa51310AjkKcBU7fLXc9IyLCl8kAmE6cBfbOBz769IhYLLozCSHdAtAjzgJ77YN3P626vq377gKwhTgLlOTWu58aGrcAtIizQKluu+cpqw2XSgA4z8TZmf7NP/vTiLj8mgelCoETcvs9V8Q617p9A8A5Is7OUf9obP/IfPLFB28tBk7eR+69PCLqOzho3AKcYeLskbR/LrZ/Xl5x3UP6xcBpueveb19vadwCnDni7LEZirYR8ZTrq3S7qL+iDZyOu+570qK+wK3GLUD5xNmd6PxQbP+8fNpzHhrDFvIunKy773tibPsyWUi3AIUQZ0/CSOP2Gc97WEixsB/uufSE1YZLJQAURJw9aSON22d+58M7xXXMlXfh5N176fHNHRw0bgH2lTh7ykYat8/67kdEyuKgyrsh78IOXbr0uKi6tqFxC7BPxNk9MtK4vfoFj2w/dZDdrBVzYRc+dumy1Ub7vruHy8F6AHZHnN1fQ43bq7/vUb3aSnWZzYOD+oKbCXV/F5jn45e+dWvXNjRuAU6EOFuGkZmEiy98dH4Ltp7Bze/vAvk+cekxq43F0pfJAE6IOFuekZmEa3/wMbFVbq+2lXd1beHIPnnfo9f33V26fQPAroizxRtp3N7wostiVDNpkN2sPXCxBZjld+59VPgyGcAOiLNnykjj9oYXN9E2f9Kgqctu1oq5kOPT93xzRITbNwAcmTh7lg01br/jRx/XL26bMGlQf/ksO8VmF8J58el7HrGIw9X2YqlxCzCNOHtejMwkPPfHqjshZSfNKZXVRioiL7IHfOFs+8zdD49qLEHjFiBJnD2PRmYSnvfj3xZ9VRJNtmBn3MYsvzL56nD2fPajD11vLSNM3AJsI84y1rh9/kufGBnyg2bTrJ1wSG5l9vXKoFSfveuhzViCxi1ARIizdIw0bp9/45OqnVVszB4PSE4a1Opr6ObPHuSHaThLPnfnQyJcKgFAnGXUUOP2BT95eb+41ho5GKlaFUwOr3ULNrl4vWSyEs6A++/4hmrTNW6B80WcJdfITML3vfLJMWUqYEZlfgt2Qie4Crz55wNFuP+OB8ey6touDzVugbNtsaz+Jw9ma/+w/P6fuqLTna3zZX9wtt+dXQ0b9INmE2c3u7Nb9tdfXIuNyqiW7c8zNKdRHXJwsHFkf7yiaScvNvd3Xqu3cutDiGhfz3dgfxO063NYNGURsVhUqSU6+3sbrcObx7H+bfXQ/mhOqfu/FZ2Ppa5vvZdltD6cLSt3XnrzFTtvbePAxXLr/uoyruuHrf9K9bkt2w+7KzSvuL2g80L9/c3GcqOyftgrqPc3F+pqVy6Wh61lNipXgXXoFdtxtjowIuI7Xvqn0SLdAmeAOMsx6zRuX/iqpyTjbCfLbuzpxNleREvG2U6+jOFXqYPLQS+cDsbZ/v7NlxNnI8TZZfvhqcfZiIjlYV3/nBv/rNkt2gJlEmfZrXa6/cGfeWpEN8tGpIPmUGu2/1SyNRu9VxlszbYOXue5evdmazaGX06cjRBnl+2H+xZn2ys8+8b/Ey3SLVAKcZaT0462P/Tqpx190qB5qt6firMTJg2qgzuhLWIwzm6JzuKsOFtOnI3NU3r2y/5fc4RoC+wxcZbT0ZlJ+OGbnp6eNIj1D+Rulq32RwzG2WQDOKo1jnHSIHqL9yOmOBvi7Hr/3sXZWDandMPL/6I5WrQF9ow4y15op9sfec0zqp2LyGnNtp86WMR4vqxKp04aRB3aBlqz7crMSYPohcj+/iaCbD4lzm7dL842pccaZ9t7bnjFX0aLdAucOnGWvdNp3L7k56+snqgL1htTJw0iYmieof5JP3vSoC4++qRB+6kmgmweIs5u3S/ONqU7i7OdF7r+lQ9E/Vi0BU6DOMu+a6fbl7z2yog6XXUT39EnDSKO7RJdsW3xKqWtHybjbJM/eod04mznJbZvRIQ4K87GMcfZ1v647u/XT0VIt8BJEWcpSadx+xOve1ZE81N66qRBREy9RFeyNRvDi/fz5exJg7p4S+Drvcq6sioQZ9sPuyuIs/XhVVm0TikZZzsHXveq1muKtsDOiLMUrJ1uX/r6q1YbyTibbs22Dl4nuXp3Ks4mW7PtjXWMG9rfzh+Lzf2bebFd2U2xvaXE2fbD7gribH14VRatU5oaZ+tziOXyup/+a9Ei3QLHSJzljGhH2xtff1VEN87uYtKgfmqoNdsq7ObL5KRB80x//2ZerEuGWrMRrfe42HhcRxlxtn52pECcjf4LNfs3yiI24mzn4XWv/uqqWrQFjkqc5QzqzCT85BuuiohknE1OGkR9bH9/HdTmThrUT+1i0qCuGcms3cmHzWQZ4qw42zvwKHG2Kj2MiOte/TXRIt0CU4mznH3tdPvyN1w99RJdydZs64jpl+ga2t8OH5tPDbVmN5/afFgVpONsb/9QnN2S9cXZ9oY4G7lxtlnwcBkR1970N6pnRVsgizjL+dJp3L7yjdc0UW8znA62ZqPJBVMnDeqNXUwa1E8NtWZjuAWbnDSI6XG2nym3LC7OirPtBQ+7K19709dGi3QLbCXOcq610+2rfvWa1a7qqfX+oTibnDSoCrfky2ScbZJH75BOnG1W2lJQbURETmu299RQlt14qlW/2qoKqsWbY+o/N3KYOBsR4mzEljjbOfDa1/ytqkC0BRriLKxtRNs3XejGtV1MGlRFg63Z1tbUSYN6oy5Ix9ne/mScTU4aNDX9xcXZjXXE2XScba987c9/fbRIt3CeibOwRWcm4WfefCEZZ3cyaVBt7WLSIKocs4tJg6jOMDlpEM0JdF9UnI0QZzsHbq7cOvDa1z2oOla0hXNHnIW0drp99VsvxHBrNmIwzg61ZvtPNbGjd8hQ1BtqzUZrtaEWbDLOJicNmkN6pzd70iD673Ezy0azZDdlirP1s+0VWp/0GYyz7QOvfd2Do0W6hTNPnIVpOo3bm952Yb1/M852smx7Y/akQV28k0mD6qlOlo1jmTSI3uKpOJtszW57qlsgzrZXOD9xtnlqdamE139j9axoC2eTOAtH0k63P/v2C8lJg+h3SVNxNtma3Xxq82FVsPeTBrGOfSPvcTOt9ven42wqrSYLxNnov1Czf6MsYl/i7Kpi9ce1r/+mao9oC2eHOAvHph1tf+4dF6udzf+NjEmD5pn+/oGoN2HSoNq1JU0OxNkJkwZRLT7Qmm1XdsP05otG/z32omoyzg61Ztt7OouLs83hVVm0TukMxNmNru0bHhot0i2US5yFnejMJLz2net0OxRn8ycNInrtyVSc3cWkQVR5sZNlozq9iG6c7We+ZJztZ9ajTxrUe7YE4lTeFWej/0LN/o2yiH2Ps83DiFgur735YVWxaAuFEWfhJLTT7evedbGV/NYbyTi7i0mD5qne/qNPGjQ1/cUH4mz+pEH91JZ8ORRnU53Xfs2WVxdn+y/U7N8oi+jl1/2Os9WfhxFx3Rv/TrRIt7DnxFk4aZ3G7S+++2IcYdKgfmoXkwYxPc7mTxpEcwLdFx2Ks/mTBrGledx99e4iyYJmo1/Z2xBnoxcWy4mzrb9WhxFx3a9+c9Q7RFvYP+IsnLJ2uv2l91ys9q7/PPqkQQy3YIfibP6kQX1Is0Iqzm559YE4uyVNpuLsLiYNWod0M9NQnO0HPnF2y8OIUuJs+5Dr3vzIaJFuYR+Is7BHNqLte5uubbI129/oZtnWrm6g7O1PxtldTBpEfebDmXX+pEGs09ZQVK0L6j1bXj0VZweTZWvPYN4VZ6OMOLt+7nC9/7q3PLraL9rCqRFnYU91ZhJ++X0XIwbjbLo129q1ii+7mDSI6gx3MWlQP7UliW4WxIxJg6pmpKDzcv39yTg7VrAZZ1uV4mz7wL2Ls60XPLz+bY+NFukWTow4C2Vop9ub339h3Z7sp9ihONuLlck4u4tJg3ojOWlQH9rJstEs2cuXqdZse09n8WTe3RKdRzY2c+HwpMG6Zqg1GzEcZ1NrirOb+08izm4ceHh4/Tsuqx6ItrBb4iyUp9O4/ZXfvDAUZ9OTBtVTnSwbxzJpEL3F+6+++brJ1uy2p7oFR580qGv6Zx69PZ28m5w0aO2p96fibLO/mxTF2bEDTzXOVvuXEXHDr//daJFu4XiJs1C8drp90y0XohVc9n7SIFY//3cxadDUpDqvdc1IwVCc7RccfdKgLk5OGrQqu2uKs5v7TznOrrbqh89+1+OrOtEWjoE4C2dKO9q++ZYLsy/RNWHSIKrFdzBpUB+SnDSoa/InDZo9yYJmo1/Z29gMhclJg4h1SD3OSYN2zSprVk+3Qmp3EXE24oTibESz4LPf/YRokW5hBnEWzqzOTMJbPnghGWeTkwax5VtWdUG1eHNM/eeyfUgyzg61Zrc91S1Ixtmh1mx7z6L7cDPKbHtqvXg/qg5Gz3p/NzUOxdljnDRoH7LYPKQV0MTZOIE4W72Xdf2z3/ukagnRFnKJs3BetNPtr33omogm47RS43rP7EmDpmagNRtHmDSon9qSL488aVDvGW7NRkSnsrd4Ks7uYtIgmmWr0lScHW7NRmuRzRPr5eBOpTi7ceDcOFu/YCwPn/P+y6NFuoUh4iycR53G7dtvvbrav94zNc4mW7PtjamX6MqfNIgtzePuq3cXSRb0Xm5Liu2F19gMr2N599gv0VUVRLVIujXbOmS9yEBrNvovN1wpzkYcKc62V1geLp97yxXVc6ItbBBngY10+44PX52cNIgmC9YFy05BJ872smzUP/nPzCW6kq3Zfs1OJg2qml1MGrSOTVcOplhxNubE2Xb9cz/wlKpOtAVxFtjUady+87arImLyJbrSrdmof/JPvUTXlnyZirNDUbUuqPckW7P9p0Zi5exJg4gqjDb7k3G2Kl1urtB6U51Iuugd0lpk88SGWrP9ynrPSIoVZ+OocbbeH4fL5936tGiRbjmHxFlgTDvdvvv2q9Y7WyGrKlvGbiYN6qd2MmlQ1YwUdF6uvz8ZZ8cK5k4a1MXJSYOo39RQa7Z1SGuRVdCM9sPY/nLbKwdT7GaWjfaa4mzrBfPjbKf+eR9+RrWEaMt5Ic4CudrR9j23Pyt6STQZZ3cxaVDX7GLSoLXRr+xtLDcre6lx9qRBU5xas65Jtmaj/TlsnvnI/MDkSYN6o59uB+Nsr1KcjQlxtvV35/C7br8yWqRbzipxFpijM5Pw3juujF6c7Qe+ZJxNThq0nuoWHH3SoK7pn3n09qybo8OxcjjO1vtTcbbZ382Ig6+73Fyh9a47kbSXZaMVzbaH1F1MGkS97FAk3dgYWFyc3bbCcvWuD5v659+5/gVLhGjLmSLOAsegnW7fd+e6IZSMs7uYNGhqMjqvybw7FGf7BUefNKiLdzJpENGJpMk4mz9p0BQPtWa3bPTWTMbZkcU7x27ZL852X/r5H11fz6R69uaAYomzwDHrNG7ff+czV8kp2ZqNXqDcxaRBsydZ0Gz0K3sbm1EvOWkQEZNvBhbrKDOyZifO7mbSoCoeruwl0e7ppeNs/5ChxUdysDjbKui89HJz//fceyEqoi3FEWeB3Wqn21vuWn1JZVk9VdcsY2uaTMXZ5KRBvSd/0qB1yMZrbVZuLt6PqoNxtt7fzYJDcXYXkwZRL9I786GQmm7NtvYMptihOJtszcbw4uk428vB4mxrwXp/vecFly5Gi3TL/hNngZPTjrYf+OjTq53r/xVKxtmh1mwcYdKg3rOLSYPmkH7BkScNolm2Kk3F2fxJg6hPbKg1O1w5lmI34+xOJg3qp4Zas+2n1osvq0fdD78XZ+sVznKcXb2N+uH3fvy6ar9oy576qlQBwLFp/yzszCR88O6Na2dGbIm5Sa0AmtLLu0kTKvupMWVGZStfpjQRKaXfcE3Jr+xk2TFb4m9S9iH91mzSYXblWfTpF35xtbFcLheb/98o3bInxFngdHR+ELbT7a13PzUybLlAWMqM1BjZh+RX9luzQ/qt2aQtrdkh/YZrWp0FdxlJk7a0ZlPyT6P5VFLvsbLMrizdp37gC+22bjvdiracInEW2Asjjdvb7mnu57mqjTz9SYMh22ZzE/K7pFt+TZ+UXTihK1nJf4/bJg1Ssiu3TRqkZC++ZdIgZUIk7U0apOWfeVE++b1/GLH+tDVuOUXiLLB3Rhq3t99zRQzoD8Um5Qe7bQO1CTMqF9mJakq+nBxJsz+/VmX2mefny1Z/d/ohSfnN4/5sblJ/cHZIsqAcn/yeP2iP+WrccpLEWWDfDTVuP3Lv5TFv0iA/r+UnmOxmbW3KWcwI0/nyU2N+ZSU/NU5vM08ImnM6wbnvccsXyJKyTyMdiPfVb3/X/fXJa9yya+IsUJKRmYS77v322CY/2M1pwU5IG/XiqZxUp+5kZa3ODakQtu2aBgkT3uOE8FrJDnZzFs/+AOdE0gmVy+rP7GNS/x3L8rHnfi6i9bdU45bjJs4CpRqZSbjrviclI2n+ZG0tP+/mT9bWkkm0NmMYd8ppVFvZ5zOhsn/3hKT8xee0YLMr64/6MHU++ZMGtf5VuobMSMb7574bfm+14ctkHBdxFjgjRhq3d9/3xObBlN/xr/5Ihtdafkrr3z0haUplJZ2PKtMrp3yQ2aX5Ldj8SYNa/iW66k9jQjJOFdTyT6OyzD6N7kVn9969138mqjdoJoHZxFngDBpp3N5z6QkxakJqnNElzQ4ZR7lEV1r2mc94j7uNpMk1G5M/wPyzaI1s5H7mE/qpdbN2+iFp2cn45N1z8dOrDV8mYypxFjj7hhq39156fP4luvInDWr5WbB/M7AMualxxm0dJsSj5eTUmN/DnrH4hLjW/OI++80mJw0qnZuB5cg/iy23FhvSvyVYWn7lbt19zad8mYxM4ixwvozMJFy69LjoyY+AM1LjjMo5WTBp+vxAfmUrz+Wfz4TVqz+zD8mvzG8eN29x+ntM/Wfq3942KX/SYNsdbvfRXVd+Mlpnq3FLhzgLnF8jMwkfu3RZ5MrNAfk3A2tkJ5gJNwNr5C7eSo2pnDS9KzkhkuZ3JWv5hzTZNfUeK/lt3Qn91MqEfmr+N8kqE8LrnP7uzt3x9E+sNtx3lxVxFmBtpHH78Uvfut4/o0uanxqrwowR2+wIWJlwia6mWZuqrGVXuhnYVs1ZTDif7MWnf+YF5d2PPO23lw8056Bxez6JswBbjDRuP3HpMTEsP++6Gdh22b+Fb2S/xyM1j5OOcImutHoqILV4M2mQHUlb4TW1+PTwOiEZH4fbnvyxiFgFXI3b80OcBUgbadx+8r5vid2kxjrBTBmxzTcjNeZXTk+N+ZXZ2atZMz9RZV+ia86dF7KnYGeMtE74jzMjkuZXzlg8+z1O9eEnXXLf3XNCnAWYZqRx+zv3Piq66kiaGzfyK5u8OyHL5EaH/MsyzAmvMyJp/uLZH+CcSDqhcnKwm5BJm2Q8/ZCUJl/mH5L9Hg8fmPyx1LMER/ShJ9xbf1oat2eMOAtwJEON20/f88394r78S3RNyJeV/IstzJkfyK88DzcDq+SntDlf4cpfPLtLmj9pUMsfxp1zWYb8M08VjLjlsrsjIqqgrHFbOnEW4NiMzCT87j0Pj1ZqTObL2pTKyoR8NLlyxrBE2pwWbHbl9ESVTGmNujCZd6efRn4WbP2bKHUalfwsOGfyIX/x7PdYm/B3NnvN3/yWj642loculVAkcRZgJ0ZmEj5z98NjwJ7cDKwlu3KnkTS7HRjTP8DmLJKfZPMWU5WV/ETVNGunH5JWtzzTH2Alu3K3effEhxPe/8i7onovywdcKqEM4izASRhp3H72ow/Nz3XtSDVaNmHSoJF/PYTl5NQ44RJd0xdPr1mrs9eELJhb2USu7DPPP4v8X/HX8iNg/fcqP5Lmn0b9aRymgmb+pEEtP+82X/PL/g8aEb/xsDsiYvlAhInbPSbOApy0kcbt5+76po2nmot55f50n5KPJs8P5Fe28lz++UxYvfoz+5D8yumzpPl5t9XDnn5ISv4luub0U7OD5ozFJ7zH5l86+YdMrhw65L1/+45Vel4+cBgmbveJOAtwyoYat5+78yH94o59uRlYLTtkzImk2bmkWTN5SJNdc99j/vUQjtJPTYewGd8kS65Zyw6vtfzFJ7zHyoTTyG7W1ib8Y2Szu/zuB912+FfrbdH2dImzAHtkZCbh/ju+oV0XeSbcDKwypROcu/iESYNa9uKtNXODyZxJg9RZNJInXGnOYsL5ZC8+/TPPz5cz8m5y0qDW5N3phyQ1kwapxZtk/MB4YUTEO//6hyNi+ZfLqKKtUHuSxFmAPTUyk3D/Rx4UW+UPv1byL9GVX9mSXbnL38JPSMaN7Mr9uBlYIzuS5l+iK+2tavUAAA81SURBVH/SoDbhkgVHaNYm825rfmC8sJF/GnVrdpVi22TZEybOApRhpHH7B7d/feSakRrzK7MTTC2/Mjt7NbJzSfMeJxySKqhlB7s5v4VveseJ/0xzIml+5YzF899j87GMFzbyT6NpvuefzwPZ58FJEWcByjPSuP2D274u+rJ/us+4WUP+4nMiaf7i2Zfoyh9+bTSJNBVlZvRTk2vWmmQ8/ZCUJs/lH5L9Hnd6M7D8Zm0ju7KeNJiwOCdOnAUo3lDj9vMf/pv94rY58wPL3NRY58v8yYcJi9cF2dkrvWajXjz3Y8lPafV7nPCR5y8+vUuaH+pnDOPmf+b5/d0Z4bVpvqf+gzankb94NWmw2njNAy8zaXDyxFmAM2VkJuHzt35tVZSfpOpIOl7Xkoojjboy/5Ds7JX/K/7GhJRWbaTiUSP7NCZMnVaSkwa1/HmG/Mpa/iFz3mPuW5xwGrUJldWkQT04yz4QZwHOrJGZhM9/6GtiTO5P991G0l22Aydcoqt5i6nKyoRGZt2snXDm2f3d5otQ2YtnV85Jjdnvccb1tvL7qfWkwYR/jOQvzmkQZwHOi5HG7Rc+8NUxbSogt3LCJbryJw1qzZqpYFJnr3QGrGSHnSZyZZ95/llMSMaV/AhYR/X8SJp/GvWnMeMSXUkz8m5y0qCWf4kukwZ7QpwFOI9GGrdfvKX7o2HCiG2TLVKVtew4Eq026WhZS37ljK+pZcejVg97+iEp+TcDq+WH1/yR1lr+4s2/MVKHtP6lk7t4/iW6WpW5i5s02DfiLACDjdsv/lbe7cayQ8acSJodMpo1JxySCjuVOddDyD6NCb+4z580qEyIpPmVlfzFJ7zHSnLN2qxmbaqi0mrWZh/DyRJnAdgwMpPwpfevK+rSyJSdYOYsnt0JnhBJj3AzsGSiapJofibN/gDnDCdkV87IuzMmDfIHVSecRvbNwGrJSYOaSYNTJ84CMGhgJmEREV9630gUWz+VHsY9ym/hkyFsTjLOrjwHNwOr5Z/GrEsWZFdm5938SYNa/iW6Rm4GxmkRZwHINda4fe9f5SeYCZXTs1d+sGvCzoRDUgW1nQa7pnecm9cmRNLsqYD8SYNaft5tncZ4YSP/NNwM7IwRZwGYY+TLZF96z19slM6IpNm5pJk0yD8ku/BIl+hKyr9EVxMBc08j/xJdTUGqspY+4cr5uRmYSYPTJc4CcAyGGrdffvef94u78i/RlT9pUGt+cZ/MgvXiycq1JqWlz3z6V7iSa9amd0nzQ/2MYdz0p1GZ09/Ny5fR+qhnXKIrqXOJLk6XOAvAMRuZSfjyu/5v5EfARnZl9u+yGxNSWrWRHY/yTyP/t/C1GZMGyWZtbULl9MUnvMfmXzqJQ2acRv4lutwMbM+JswDs0MhMwpd//c9iq122AyfcDKyS0dZdm9HITKa0Rn5/NzulNbK7pHNSY/YhM663ld9P3d3NwEwanDpxFoCTM9a4ffufRqYmFeUGk/xImn+JrgmTBpX8S3TtNhm7Gdg2+TcDq5k02BPiLACnY6xx+2v/u11Y/ZkdHfKD5oyvqWXHo9ZpTD8k5TzcDKw24zSSH3n+pEHNJbr2ljgLwF4Ya9y+9X/FiOmNzHTYqUy480It+zTyfwtfyz7xWZE0vzJ78TnvsTlkvHBeszZVUcm8GZhJg30gzgKwd8Yat2/5n9VmbjCZE0knVE5OVOmYVqsjYOqQOcMJh9mHZIfX2oxJg8xB1Zh0Gru/GRj7QJwFYN8NXgXszX/SL54wadCog10qaB7hZmCpRDorGSeTaKXVHM0+JPs05lyWIf/MZ+Td5EddadrM+YtLsftHnAWgJGMzCb/y3yOHm4Ftc5SpgKRZeTdVUck/86PcDGzrJbpMGuwJcRaAUo3NJLzxv3WKJ0TSJnrlRqr8YNc0a6cfkpb99aYJkwaVCalx+khr/hTsjGbtjEmDCYuzB8RZAM6Ikcbtl375v0ZSctKg0kSuVPaaMGlQmTAFO32kNT/UTziNWvZp5H+TrJafL0/sZmBas/tDnAXgDBpp3H7pH/2XiCmX6Kqfz45HM4LdhEOmTxok+6m1CZXTF5/wHuupgAmHTK7MP8TNwPacOAvA2TfUuP3S6/9zv3hlzqRBdjzKT2lz+rvZpzGjv5u/eH7ezZ80qE045Ag3Axu/RBf7Q5wF4HwZm0n4xf+U7NU2kkGqkn8zsEYqAtamhOncfNnIz7tVwYxLdCVNCK+VGZMGky7RZdJgr4izAJxfIzMJX/yFP94ozb9EVzNZm5uo6kOSi7sZ2Fat+YHxwkb+abgZ2P4TZwFgbaRx+8Wf+48xLjse5f8WvjYjpU3IgvmVMxbPf4/NxzJe2Mg/jfxLdJk0KJE4CwBbjDRuv/Ca/5DfT23kx7Tm+2HTD0k5yiW6kvIv0TXjelv5h+z6ZmAmDfaNOAsAaUON28/f9O/7xRFNEk0n0hnJOLtL2mpG5i6eP4w757IMqTVr+eG1lv6oK02bOX9xkwZ7TJwFgGlGZhI+/+p/F3kmBLv6F+XTL9GVNGfyIX/xOXk3VVGZc+b5laM3A2PfiLMAMN/ITMLnf+bftuqyv0lWy092h5ODXbLzWpuTGrPfY/4lC2Y0a/Mv0TXpZmAmDfaQOAsAx2akcfuHr/qjGJb/K/5aMgK25EbSGadRd16Tl+ia802y1Jq1E7sZ2Hglp0KcBYCdGGnc/v4r/3VdFHnyL9E1p5+aXzl98QnvsR6rmHDI5Mr8Q0walEKcBYCTMNS4/f1X/Kt+8YRvklXyU1r+N8lq+YvPyLsTTiN7OKGW3avNukSXSYP9JM4CwEkbmUm4/2X/MjLMmArIz5cz8m5y0qCWPwU7K7zmXqJr3s3AUiWcDnEWAE7TyEzC/Tf+i3VN9qRBLf8SXXNGWvMvWbDTZm0zPzBe2Mg/DTcDK4g4CwB7ZKRx+7mX/vMYMCeS5lfOWDw7NbaGE8YLG/mnkX8zsJpJgxKJswCwp0Yat5/9sa9EzGrW5h+SnxqnTwUkJw1q+cMJ+ZMGtfxLdOnR7jNxFgDKMNK4/b0f+aexXW4Iyx/GnTBpUMnv7+aH19qEtu7h9MWl2BIskn+xAIA91063n3nxP+kPzna6s/3rD/TjbKdmyyULqj2d7my/oP8qy+4h0X7YL6j3bDnz/hVnN2+g0L9EV/97YOvFezcDq+Psz/75jSYN9pbuLAAUb6hx+7sv+vJGXarzWtu3623l91PzbwZWm7A4e0mcBYAzZWQm4dMv/GKMyw6v/dZs0ozwmjzkKDcDS16iy83ASiHOAsCZ1fsyWbP9qR/4wpyR1ul5N3nI/t8MzDUN9pw4CwDnxWbjduOpT37vH0ZPfgSckXcnhNfsZm0tu1ebdTMw9pw4CwDn0Ujj9hPf/fuRZ8akQf6g6pTwOv8SXUkmDfafOAsADDZuf/u77m9q8i/RdYRmbTLvHuVmYMnF+9c0MGmw/8RZAGDDyEzCx577uRg2J7zmR9Jkhq4c783A2H/iLAAwaGQm4b4bfq+qyU+N2ZXZzdpGdmX+zcAogjgLAOQaadzee/1nusXZV06o5efL/Et05U8a1Nw9oSziLAAwx0jj9p6Ln44M+ZfomnPlhPzK3s3AKIs4CwAcg6HG7d3XfGpdMP16W/n9VDcDO8/EWQDgmI3MJNx15SdjwIy8m5w0qOXfDKxm0qAU4iwAsEMjMwl3PP0TMSW85l+i6yg3A3Oh2eKIswDAyRlp3N5+xcc3Kuc0a1MVFTcDO0vEWQDgdIw0bm978sciZac3AzNpUBBxFgDYC0ON21ufeN9GWXY7Nf8SXSYNiibOAgB7Z2Qm4YOPuze2cTOwc0ucBQD22shMwi2X3T1j0iDZrDVpUBZxFgAoyUjj9je/5aPRkwyvNZMGhRJnAYBSjTRu3/eIO2PU0M3AtGaLI84CAGfEUOP2Nx52R1OT3aylFOIsAHAGjcwkvOcbPxKcIeIsAHDGjcwkvOvrblvX/OUyTBqUSZwFAM6XkcYtJRJnAYDzq9e43XhIERb590EGAIB9c5AqAACA/SXOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABSs1Dh7eHjYedjZAwDAeVBqnD04OOjk14ODUt8LAACzlRcBv/KVxVe+slhtrxLt4eGhLAsAcD4VnAJFWAAACkuEl1++eOxj47GPjcsvX0Q1ciDXAgCcW4UFwRe9qLsBAMB5Vlicvf767gYAAOfZYrlcpmr2zWKxiPq0DRsAAJxnRcbZiOLOGQCAndDXBACgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZ9cODw8PDw/bD0eKAQDYE+Ls2sFB81EcHh62HwIAsLe+KlVwjhwcHGjKAgCURQ9yC61ZAIBSFNadXfz4YvVnRCzvWI4XAwBw5pXUhqyy7PaHR/GVxeLNi8Xli8XBwYGRAwCAghTWnd2FxWLxPyIeG/HyiK8sFhFx1XLp22AAAEU4W3F2Madfu4z4k95OWRYAoAhnK84uZ07T/vFi8cWISxF/NHcFAABOxWJZVIBrz8v6KhgAAIXFWQAAaDMhCgBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACiYOAsAQMHEWQAACibOAgBQMHEWAICCibMAABRMnAUAoGDiLAAABRNnAQAomDgLAEDBxFkAAAomzgIAUDBxFgCAgomzAAAUTJwFAKBg4iwAAAUTZwEAKJg4CwBAwcRZAAAKJs4CAFAwcRYAgIKJswAAFEycBQCgYOIsAAAFE2cBACjY/weRifNSE6DioAAAAABJRU5ErkJggg==" ], "text/plain": [ "PyObject " ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using PyCall\n", "@pyimport IPython.display as d\n", "d.Image(\"/tmp/displacement.png\")" ] }, { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "## 3d beam with quadratic elements" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "27-Jun 23:44:05: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", " in parse_abaqus at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:113\n", " in include_string at loading.jl:99\n", " in execute_request_0x535c5df2 at /Users/jukka/.julia/v0.4/IJulia/src/execute_request.jl:157\n", " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[4], in expression starting on line 2\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", " in parse_abaqus at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:113\n", " in include_string at loading.jl:99\n", " in execute_request_0x535c5df2 at /Users/jukka/.julia/v0.4/IJulia/src/execute_request.jl:157\n", " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[4], in expression starting on line 2\n", "27-Jun 23:44:06:DEBUG:root:Found NODE section\n", "27-Jun 23:44:07: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", " in map at abstractarray.jl:1251\n", " in parse_element_section at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:59\n", " in process_section at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:108\n", " in parse_abaqus at /Users/jukka/.julia/v0.4/JuliaFEM/src/abaqus_reader.jl:117\n", " in include_string at loading.jl:99\n", " in execute_request_0x535c5df2 at /Users/jukka/.julia/v0.4/IJulia/src/execute_request.jl:157\n", " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[4], in expression starting on line 2\n", "27-Jun 23:44:08:DEBUG:root:120 elements found\n", "27-Jun 23:44:08:INFO:root:Creating ELSET Body1\n", "27-Jun 23:44:08:DEBUG:root:Found NSET section\n", "27-Jun 23:44:08:DEBUG:root:Creating node set SUPPORT\n", "27-Jun 23:44:08:DEBUG:root:Found NSET section\n", "27-Jun 23:44:08:DEBUG:root:Creating node set LOAD\n", "27-Jun 23:44:08:DEBUG:root:Found NSET section\n", "27-Jun 23:44:08:DEBUG:root:Creating node set TOP\n" ] }, { "data": { "text/plain": [ "Dict{Any,Any} with 4 entries:\n", " \"nodes\" => Dict{Any,Any}(288=>[97.5,7.5,10.0],11=>[92.5,2.5,5.0],134=>[45.…\n", " \"elements\" => Dict{Any,Any}(68=>[71,144,149,198,51,150,57,43,50,214],2=>[204,…\n", " \"elsets\" => Dict{Any,Any}(\"Body1\"=>[1,2,3,4,5,6,7,8,9,10 … 111,112,113,11…\n", " \"nsets\" => Dict{Any,Any}(\"LOAD\"=>[82,84,87,179,197,246,249,256,257],\"SUPPO…" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "fid = open(\"../geometry/3d_beam/palkki.inp\")\n", "model = JuliaFEM.abaqus_reader.parse_abaqus(fid)\n", "close(fid)\n", "model" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "10x120 Array{Int64,2}:\n", " 243 204 259 145 96 96 236 285 … 217 69 154 179 203 96 259\n", " 240 199 70 175 88 101 88 179 216 144 114 91 204 267 199\n", " 191 175 69 199 236 164 285 178 278 78 278 178 259 95 204\n", " 117 130 130 130 178 97 178 83 155 71 218 83 199 97 130\n", " 245 207 265 177 141 102 290 12 219 146 20 181 206 268 39\n", " 242 208 72 208 290 171 289 182 … 282 152 280 180 263 272 207\n", " 244 209 5 202 291 9 287 11 33 79 32 182 262 98 263\n", " 1 3 6 174 7 99 237 13 224 74 223 14 205 99 6\n", " 2 4 132 176 8 103 8 14 225 51 284 93 207 24 4\n", " 196 176 134 4 237 10 11 15 17 80 283 15 39 100 3" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nnodes = length(model[\"nodes\"])\n", "nelements = length(model[\"elements\"])\n", "dim = 3\n", "E = 90\n", "nu = 0.25\n", "mu = E/(2*(1+nu))\n", "la = E*nu/((1+nu)*(1-2*nu))\n", "\n", "X = zeros(dim, nnodes)\n", "u = zeros(dim, nnodes)\n", "du = zeros(dim, nnodes)\n", "elmap = zeros(Int, 10, nelements)\n", "nodalloads = zeros(3, nnodes)\n", "dirichletbc = NaN*ones(3, nnodes)\n", "la = la*ones(1, nnodes)\n", "mu = mu*ones(1, nnodes)\n", "\n", "# calculate permutation which maps node ids to matrix indices\n", "perm = Dict()\n", "for (j, k) in enumerate(keys(model[\"nodes\"]))\n", " perm[k] = j\n", "end\n", "\n", "for j=1:nnodes\n", " #X[:,j] = model[\"nodes\"][perm[j]]\n", " X[:,j] = model[\"nodes\"][j]\n", "end\n", "\n", "for (j, k) in enumerate(keys(model[\"elements\"]))\n", " #node_ids = model[\"elements\"][j]\n", " elmap[:,j] = model[\"elements\"][j]\n", " #for l=1:10\n", " # elmap[l, j] = perm[node_ids[l]]\n", " #end\n", "end\n", "\n", "elmap" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "3x298 Array{Float64,2}:\n", " NaN NaN NaN NaN NaN NaN NaN NaN … NaN NaN NaN NaN NaN NaN NaN\n", " NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN\n", " NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN NaN" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Handle dirichlet boundaries on SUPPORT\n", "for j in model[\"nsets\"][\"SUPPORT\"]\n", " #dirichletbc[perm[j]] = 0.0\n", " dirichletbc[j] = 0.0\n", "end\n", "dirichletbc" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Shape functions and integration points" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "dNtet (generic function with 1 method)" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Ntet(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", "dNtet(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]" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "1" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sum(Ntet([0, 1, 1]))" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "1x4 Array{Float64,2}:\n", " 0.0416667 0.0416667 0.0416667 0.0416667" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# from code aster documentation\n", "a = 1/20*(5-sqrt(5))\n", "b = 1/20*(5+3*sqrt(5))\n", "ipoints = [a a a; a a b; a b a; b a a]\n", "iweights = 1/24*[1 1 1 1]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Add point force to LOAD nodeset" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "9-element Array{Int64,1}:\n", " 82\n", " 84\n", " 87\n", " 179\n", " 197\n", " 246\n", " 249\n", " 256\n", " 257" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "model[\"nsets\"][\"LOAD\"]" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[100.0,10.0,0.0]\n" ] } ], "source": [ "#nodalloads[3, perm[82]] = -0.06\n", "nodalloads[3, 82] = -50.0\n", "println(model[\"nodes\"][82])" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "iteration 1, norm = 93.13879357262775\n", "iteration 2, norm = 14.241640589907869\n", "iteration 3, norm = 2.3999534862146534\n", "iteration 4, norm = 0.8191634532858331\n", "iteration 5, norm = 0.08574651168344295\n", "iteration 6, norm = 0.00032306479346201524\n", "iteration 7, norm = 8.052462405199464e-9\n", "iteration 8, norm = 3.0878906643174355e-14\n" ] }, { "data": { "text/plain": [ "3x298 Array{Float64,2}:\n", " -0.0479657 -0.0483971 -0.0497284 … -0.0155295 -0.0150513 -0.0462958\n", " -0.474054 -0.418817 -0.254176 -0.229641 -0.283098 -0.694379 \n", " 0.232677 0.198051 0.0948355 0.128246 0.162511 0.371043 " ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" }, { "name": "stdout", "output_type": "stream", "text": [ "Converged\n" ] } ], "source": [ "u = zeros(dim, nnodes)\n", "du = zeros(dim, nnodes)\n", "\n", "for i=1:10\n", " JuliaFEM.elasticity_solver.solve_elasticity_increment!(X, u, du, elmap, nodalloads, dirichletbc,\n", " la, mu, Ntet, dNtet, ipoints, iweights)\n", " u += du\n", " println(\"iteration $i, norm = $(norm(du))\")\n", " if norm(du) < 1.0e-9\n", " println(\"Converged\")\n", " break\n", " end\n", "end\n", "u" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "15.449170689704438" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "maximum(abs(u))" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(10,120)" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "size(elmap)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "120" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "nelements" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "11x120 Array{Int64,2}:\n", " 38 38 38 38 38 38 38 38 … 38 38 38 38 38 38 38\n", " 243 204 259 145 96 96 236 285 217 69 154 179 203 96 259\n", " 240 199 70 175 88 101 88 179 216 144 114 91 204 267 199\n", " 191 175 69 199 236 164 285 178 278 78 278 178 259 95 204\n", " 117 130 130 130 178 97 178 83 155 71 218 83 199 97 130\n", " 245 207 265 177 141 102 290 12 … 219 146 20 181 206 268 39\n", " 242 208 72 208 290 171 289 182 282 152 280 180 263 272 207\n", " 244 209 5 202 291 9 287 11 33 79 32 182 262 98 263\n", " 1 3 6 174 7 99 237 13 224 74 223 14 205 99 6\n", " 2 4 132 176 8 103 8 14 225 51 284 93 207 24 4\n", " 196 176 134 4 237 10 11 15 … 17 80 283 15 39 100 3" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "elcodes = 0x0026*ones(Int, nelements)\n", "elmap2 = [elcodes'\n", " elmap]" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "WARNING: int(x) is deprecated, use Int(x) instead.\n" ] }, { "data": { "text/plain": [ "9568" ] }, "execution_count": 17, "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 int at deprecated.jl:49\n", " in save_file at /Users/jukka/.julia/v0.4/LightXML/src/document.jl:108\n", " in xdmf_save_model at /Users/jukka/.julia/v0.4/JuliaFEM/src/xdmf.jl:140\n", " in include_string at loading.jl:99\n", " in execute_request_0x535c5df2 at /Users/jukka/.julia/v0.4/IJulia/src/execute_request.jl:157\n", " in eventloop at /Users/jukka/.julia/v0.4/IJulia/src/IJulia.jl:123\n", " in anonymous at task.jl:365\n", "while loading In[17], in expression starting on line 7\n" ] } ], "source": [ "xdoc, xmodel = JuliaFEM.xdmf.xdmf_new_model()\n", "temporal_collection = JuliaFEM.xdmf.xdmf_new_temporal_collection(xmodel)\n", "grid = JuliaFEM.xdmf.xdmf_new_grid(temporal_collection; time=0)\n", "JuliaFEM.xdmf.xdmf_new_mesh(grid, X, elmap2)\n", "#JuliaFEM.xdmf.xdmf_new_field(grid, \"Displacement\", \"nodes\", u)\n", "#print(xdoc)\n", "JuliaFEM.xdmf.xdmf_save_model(xdoc, \"/tmp/foo3d2.xmf\")" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "10x2 Array{Int64,2}:\n", " 243 145\n", " 240 199\n", " 191 69\n", " 117 130\n", " 245 202\n", " 242 47\n", " 244 148\n", " 1 174\n", " 2 4\n", " 196 134" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "elmap[:,[1, 101]]" ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "3x10 Array{Float64,2}:\n", " 20.0 30.0 20.0 20.0 25.0 25.0 20.0 20.0 25.0 20.0\n", " 0.0 0.0 0.0 10.0 0.0 0.0 0.0 5.0 5.0 5.0\n", " 10.0 10.0 0.0 0.0 10.0 5.0 5.0 5.0 5.0 0.0" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "tmp = elmap[:,[1]]\n", "tmp = reshape(tmp, length(tmp))\n", "X[:, tmp]" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "11x1 Array{Int64,2}:\n", " 38\n", " 1\n", " 2\n", " 3\n", " 4\n", " 5\n", " 6\n", " 7\n", " 8\n", " 9\n", " 10" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#tmpelmap = [38 1 2 3 4 5 6 7 8 9 10; 38 11 12 13 14 15 16 17 18 19 20]'\n", "tmpelmap = [38 1 2 3 4 5 6 7 8 9 10]'" ] }, { "cell_type": "code", "execution_count": 48, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "\n", "\n" ] }, { "data": { "text/plain": [ "28476" ] }, "execution_count": 48, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using LightXML\n", "\n", "function xdmf_new_mesh(grid, X, elmap)\n", " dim, nnodes = size(X)\n", " geometry = new_child(grid, \"Geometry\")\n", " set_attribute(geometry, \"Type\", \"XYZ\")\n", " dataitem = new_child(geometry, \"DataItem\")\n", " set_attribute(dataitem, \"DataType\", \"Float\")\n", " set_attribute(dataitem, \"Dimensions\", \"$nnodes $dim\")\n", " set_attribute(dataitem, \"Format\", \"XML\")\n", " set_attribute(dataitem, \"Precision\", 8)\n", " #add_text(dataitem, join(X, \" \"))\n", " s = \"\\n\"\n", " \n", " for i=1:nnodes\n", " s *= \"\\t\\t\" * join(X[:,i], \" \") * \"\\n\"\n", " end\n", " s *= \" \"\n", " add_text(dataitem, s)\n", "\n", " elmap2 = copy(elmap)\n", " elmap2[2:end,:] -= 1\n", " dim, nelements = size(elmap2)\n", "\n", " topology = new_child(grid, \"Topology\")\n", " #set_attribute(topology, \"Dimensions\", \"1\")\n", " set_attribute(topology, \"TopologyType\", \"Mixed\")\n", " set_attribute(topology, \"NumberOfElements\", nelements)\n", " dataitem = new_child(topology, \"DataItem\")\n", " set_attribute(dataitem, \"DataType\", \"Int\")\n", " set_attribute(dataitem, \"Dimensions\", \"$nelements $dim\")\n", " set_attribute(dataitem, \"Format\", \"XML\")\n", " set_attribute(dataitem, \"Precision\", 8)\n", " s = \"\\n\"\n", " for i=1:nelements\n", " s *= \"\\t\\t\" * join(elmap2[:,i], \" \") * \"\\n\"\n", " end\n", " add_text(dataitem, s)\n", " #add_text(dataitem, join(elmap2, \" \"))\n", " \n", "end\n", "\n", "xdoc, xmodel = JuliaFEM.xdmf.xdmf_new_model()\n", "temporal_collection = JuliaFEM.xdmf.xdmf_new_temporal_collection(xmodel)\n", "grid = JuliaFEM.xdmf.xdmf_new_grid(temporal_collection; time=0)\n", "xdmf_new_mesh(grid, X, elmap2)\n", "JuliaFEM.xdmf.xdmf_new_field(grid, \"Displacement\", \"nodes\", u)\n", "print(xdoc)\n", "JuliaFEM.xdmf.xdmf_save_model(xdoc, \"/tmp/foo3d.xmf\")" ] } ], "metadata": { "kernelspec": { "display_name": "Julia 0.4.0-dev", "language": "julia", "name": "julia-0.4" }, "language_info": { "name": "julia", "version": "0.4.0" } }, "nbformat": 4, "nbformat_minor": 0 }