{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Finite differences for linear elasticity" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Definitions" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "from sympy import *\n", "from continuum_mechanics.vector import lap, sym_grad\n", "from continuum_mechanics.solids import navier_cauchy, strain_stress" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "init_printing()" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "x, y = symbols(\"x y\")\n", "lamda, mu, h = symbols(\"lamda mu h\")" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "def construct_poly(pts, terms, var=\"u\"):\n", " npts = len(pts)\n", " u = symbols(\"{}:{}\".format(var, npts))\n", " vander = Matrix(npts, npts,\n", " lambda i, j: (terms[j]).subs({x: pts[i][0],\n", " y: pts[i][1]}))\n", " inv_vander = simplify(vander.inv())\n", " shape_funs = simplify(inv_vander.T * Matrix(terms))\n", " poly = sum(Matrix(u).T * shape_funs)\n", " return poly\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Nine-point stencil" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can compute the finite difference for elasticity applying\n", "the Navier-Cauchy operator to a polynomial interpolator for\n", "our stencil.\n", "\n", "