reissner_mindlin_quads.py.html 26.3 KB
Newer Older
Jeremy BLEYER's avatar
Jeremy BLEYER committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27

<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
  "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">

<html xmlns="http://www.w3.org/1999/xhtml">
  <head>
    <meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
    <title>Reissner-Mindlin plate with Quadrilaterals &#8212; Numerical tours of continuum mechanics using FEniCS master documentation</title>
    <link rel="stylesheet" href="../../_static/agogo.css" type="text/css" />
    <link rel="stylesheet" href="../../_static/pygments.css" type="text/css" />
    <script type="text/javascript">
      var DOCUMENTATION_OPTIONS = {
        URL_ROOT:    '../../',
        VERSION:     'master',
        COLLAPSE_INDEX: false,
        FILE_SUFFIX: '.html',
        HAS_SOURCE:  true,
        SOURCELINK_SUFFIX: '.txt'
      };
    </script>
    <script type="text/javascript" src="../../_static/jquery.js"></script>
    <script type="text/javascript" src="../../_static/underscore.js"></script>
    <script type="text/javascript" src="../../_static/doctools.js"></script>
    <script type="text/javascript" src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.1/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>
    <link rel="index" title="Index" href="../../genindex.html" />
    <link rel="search" title="Search" href="../../search.html" />
    <link rel="next" title="Reissner-Mindlin plate with a Discontinuous-Galerkin approach" href="reissner_mindlin_dg.py.html" />
Jeremy BLEYER's avatar
Jeremy BLEYER committed
28
    <link rel="prev" title="Eulerian buckling of a beam" href="../timoshenko/beam_buckling.html" /> 
Jeremy BLEYER's avatar
Jeremy BLEYER committed
29 30 31 32 33 34 35
  </head>
  <body>
    <div class="header-wrapper" role="banner">
      <div class="header">
        <div class="headertitle"><a
          href="../../index.html">Numerical tours of continuum mechanics using FEniCS</a></div>
        <div class="rel" role="navigation" aria-label="related navigation">
Jeremy BLEYER's avatar
Jeremy BLEYER committed
36
          <a href="../timoshenko/beam_buckling.html" title="Eulerian buckling of a beam"
Jeremy BLEYER's avatar
Jeremy BLEYER committed
37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201
             accesskey="P">previous</a> |
          <a href="reissner_mindlin_dg.py.html" title="Reissner-Mindlin plate with a Discontinuous-Galerkin approach"
             accesskey="N">next</a> |
          <a href="../../genindex.html" title="General Index"
             accesskey="I">index</a>
        </div>
       </div>
    </div>

    <div class="content-wrapper">
      <div class="content">
        <div class="document">
            
      <div class="documentwrapper">
        <div class="bodywrapper">
          <div class="body" role="main">
            
  <div class="section" id="reissner-mindlin-plate-with-quadrilaterals">
<span id="reissnermindlinquads"></span><h1>Reissner-Mindlin plate with Quadrilaterals<a class="headerlink" href="#reissner-mindlin-plate-with-quadrilaterals" title="Permalink to this headline"></a></h1>
<div class="section" id="introduction">
<h2>Introduction<a class="headerlink" href="#introduction" title="Permalink to this headline"></a></h2>
<p>This program solves the Reissner-Mindlin plate equations on the unit
square with uniform transverse loading and fully clamped boundary conditions.
The corresponding file can be obtained from <a class="reference download internal" href="../../_downloads/reissner_mindlin_quads.py" download=""><code class="xref download docutils literal"><span class="pre">reissner_mindlin_quads.py</span></code></a>.</p>
<p>It uses quadrilateral cells and selective reduced integration (SRI) to
remove shear-locking issues in the thin plate limit. Both linear and
quadratic interpolation are considered for the transverse deflection
<span class="math">\(w\)</span> and rotation <span class="math">\(\underline{\theta}\)</span>.</p>
<div class="admonition note">
<p class="first admonition-title">Note</p>
<p class="last">Note that for a structured square grid such as this example, quadratic
quadrangles will not exhibit shear locking because of the strong symmetry (similar
to the criss-crossed configuration which does not lock). However, perturbating
the mesh coordinates to generate skewed elements suffice to exhibit shear locking.</p>
</div>
<p>The solution for <span class="math">\(w\)</span> in this demo will look as follows:</p>
<a class="reference internal image-reference" href="../../_images/clamped_40x40.png"><img alt="../../_images/clamped_40x40.png" src="../../_images/clamped_40x40.png" style="width: 635.2px; height: 374.4px;" /></a>
</div>
<div class="section" id="implementation">
<h2>Implementation<a class="headerlink" href="#implementation" title="Permalink to this headline"></a></h2>
<p>Material parameters for isotropic linear elastic behavior are first defined:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="kn">from</span> <span class="nn">__future__</span> <span class="k">import</span> <span class="n">print_function</span>
<span class="kn">from</span> <span class="nn">fenics</span> <span class="k">import</span> <span class="o">*</span>

<span class="n">E</span> <span class="o">=</span> <span class="n">Constant</span><span class="p">(</span><span class="mf">1e3</span><span class="p">)</span>
<span class="n">nu</span> <span class="o">=</span> <span class="n">Constant</span><span class="p">(</span><span class="mf">0.3</span><span class="p">)</span>
</pre></div>
</div>
<p>Plate bending stiffness <span class="math">\(D=\dfrac{Eh^3}{12(1-\nu^2)}\)</span> and shear stiffness <span class="math">\(F = \kappa Gh\)</span>
with a shear correction factor <span class="math">\(\kappa = 5/6\)</span> for a homogeneous plate
of thickness <span class="math">\(h\)</span>:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">thick</span> <span class="o">=</span> <span class="n">Constant</span><span class="p">(</span><span class="mf">1e-3</span><span class="p">)</span>
<span class="n">D</span> <span class="o">=</span> <span class="n">E</span><span class="o">*</span><span class="n">thick</span><span class="o">**</span><span class="mi">3</span><span class="o">/</span><span class="p">(</span><span class="mi">1</span><span class="o">-</span><span class="n">nu</span><span class="o">**</span><span class="mi">2</span><span class="p">)</span><span class="o">/</span><span class="mf">12.</span>
<span class="n">F</span> <span class="o">=</span> <span class="n">E</span><span class="o">/</span><span class="mi">2</span><span class="o">/</span><span class="p">(</span><span class="mi">1</span><span class="o">+</span><span class="n">nu</span><span class="p">)</span><span class="o">*</span><span class="n">thick</span><span class="o">*</span><span class="mf">5.</span><span class="o">/</span><span class="mf">6.</span>
</pre></div>
</div>
<p>The uniform loading <span class="math">\(f\)</span> is scaled by the plate thickness so that the deflection converges to a
constant value in the thin plate Love-Kirchhoff limit:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">f</span> <span class="o">=</span> <span class="n">Constant</span><span class="p">(</span><span class="o">-</span><span class="n">thick</span><span class="o">**</span><span class="mi">3</span><span class="p">)</span>
</pre></div>
</div>
<p>The unit square mesh is divided in <span class="math">\(N\times N\)</span> quadrilaterals:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">N</span> <span class="o">=</span> <span class="mi">50</span>
<span class="n">mesh</span> <span class="o">=</span> <span class="n">UnitSquareMesh</span><span class="o">.</span><span class="n">create</span><span class="p">(</span><span class="n">N</span><span class="p">,</span> <span class="n">N</span><span class="p">,</span> <span class="n">CellType</span><span class="o">.</span><span class="n">Type_quadrilateral</span><span class="p">)</span>
</pre></div>
</div>
<p>Continuous interpolation using of degree <span class="math">\(d=\texttt{deg}\)</span> is chosen for both deflection and rotation:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">deg</span> <span class="o">=</span> <span class="mi">1</span>
<span class="n">We</span> <span class="o">=</span> <span class="n">FiniteElement</span><span class="p">(</span><span class="s2">&quot;Lagrange&quot;</span><span class="p">,</span> <span class="n">mesh</span><span class="o">.</span><span class="n">ufl_cell</span><span class="p">(),</span> <span class="n">deg</span><span class="p">)</span>
<span class="n">Te</span> <span class="o">=</span> <span class="n">VectorElement</span><span class="p">(</span><span class="s2">&quot;Lagrange&quot;</span><span class="p">,</span> <span class="n">mesh</span><span class="o">.</span><span class="n">ufl_cell</span><span class="p">(),</span> <span class="n">deg</span><span class="p">)</span>
<span class="n">V</span> <span class="o">=</span> <span class="n">FunctionSpace</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">MixedElement</span><span class="p">([</span><span class="n">We</span><span class="p">,</span> <span class="n">Te</span><span class="p">]))</span>
</pre></div>
</div>
<p>Clamped boundary conditions on the lateral boundary are defined as:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="k">def</span> <span class="nf">border</span><span class="p">(</span><span class="n">x</span><span class="p">,</span> <span class="n">on_boundary</span><span class="p">):</span>
    <span class="k">return</span> <span class="n">on_boundary</span>

<span class="n">bc</span> <span class="o">=</span>  <span class="p">[</span><span class="n">DirichletBC</span><span class="p">(</span><span class="n">V</span><span class="p">,</span> <span class="n">Constant</span><span class="p">((</span><span class="mf">0.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">,</span> <span class="mf">0.</span><span class="p">)),</span> <span class="n">border</span><span class="p">)]</span>
</pre></div>
</div>
<p>Some useful functions for implementing generalized constitutive relations are now
defined:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="k">def</span> <span class="nf">strain2voigt</span><span class="p">(</span><span class="n">eps</span><span class="p">):</span>
    <span class="k">return</span> <span class="n">as_vector</span><span class="p">([</span><span class="n">eps</span><span class="p">[</span><span class="mi">0</span><span class="p">,</span> <span class="mi">0</span><span class="p">],</span> <span class="n">eps</span><span class="p">[</span><span class="mi">1</span><span class="p">,</span> <span class="mi">1</span><span class="p">],</span> <span class="mi">2</span><span class="o">*</span><span class="n">eps</span><span class="p">[</span><span class="mi">0</span><span class="p">,</span> <span class="mi">1</span><span class="p">]])</span>
<span class="k">def</span> <span class="nf">voigt2stress</span><span class="p">(</span><span class="n">S</span><span class="p">):</span>
    <span class="k">return</span> <span class="n">as_tensor</span><span class="p">([[</span><span class="n">S</span><span class="p">[</span><span class="mi">0</span><span class="p">],</span> <span class="n">S</span><span class="p">[</span><span class="mi">2</span><span class="p">]],</span> <span class="p">[</span><span class="n">S</span><span class="p">[</span><span class="mi">2</span><span class="p">],</span> <span class="n">S</span><span class="p">[</span><span class="mi">1</span><span class="p">]]])</span>
<span class="k">def</span> <span class="nf">curv</span><span class="p">(</span><span class="n">u</span><span class="p">):</span>
    <span class="p">(</span><span class="n">w</span><span class="p">,</span> <span class="n">theta</span><span class="p">)</span> <span class="o">=</span> <span class="n">split</span><span class="p">(</span><span class="n">u</span><span class="p">)</span>
    <span class="k">return</span> <span class="n">sym</span><span class="p">(</span><span class="n">grad</span><span class="p">(</span><span class="n">theta</span><span class="p">))</span>
<span class="k">def</span> <span class="nf">shear_strain</span><span class="p">(</span><span class="n">u</span><span class="p">):</span>
    <span class="p">(</span><span class="n">w</span><span class="p">,</span> <span class="n">theta</span><span class="p">)</span> <span class="o">=</span> <span class="n">split</span><span class="p">(</span><span class="n">u</span><span class="p">)</span>
    <span class="k">return</span> <span class="n">theta</span><span class="o">-</span><span class="n">grad</span><span class="p">(</span><span class="n">w</span><span class="p">)</span>
<span class="k">def</span> <span class="nf">bending_moment</span><span class="p">(</span><span class="n">u</span><span class="p">):</span>
    <span class="n">DD</span> <span class="o">=</span> <span class="n">as_tensor</span><span class="p">([[</span><span class="n">D</span><span class="p">,</span> <span class="n">nu</span><span class="o">*</span><span class="n">D</span><span class="p">,</span> <span class="mi">0</span><span class="p">],</span> <span class="p">[</span><span class="n">nu</span><span class="o">*</span><span class="n">D</span><span class="p">,</span> <span class="n">D</span><span class="p">,</span> <span class="mi">0</span><span class="p">],[</span><span class="mi">0</span><span class="p">,</span> <span class="mi">0</span><span class="p">,</span> <span class="n">D</span><span class="o">*</span><span class="p">(</span><span class="mi">1</span><span class="o">-</span><span class="n">nu</span><span class="p">)</span><span class="o">/</span><span class="mf">2.</span><span class="p">]])</span>
    <span class="k">return</span> <span class="n">voigt2stress</span><span class="p">(</span><span class="n">dot</span><span class="p">(</span><span class="n">DD</span><span class="p">,</span><span class="n">strain2voigt</span><span class="p">(</span><span class="n">curv</span><span class="p">(</span><span class="n">u</span><span class="p">))))</span>
<span class="k">def</span> <span class="nf">shear_force</span><span class="p">(</span><span class="n">u</span><span class="p">):</span>
    <span class="k">return</span> <span class="n">F</span><span class="o">*</span><span class="n">shear_strain</span><span class="p">(</span><span class="n">u</span><span class="p">)</span>
</pre></div>
</div>
<p>The contribution of shear forces to the total energy is under-integrated using
a custom quadrature rule of degree <span class="math">\(2d-2\)</span> i.e. for linear (<span class="math">\(d=1\)</span>)
quadrilaterals, the shear energy is integrated as if it were constant (1 Gauss point instead of 2x2)
and for quadratic (<span class="math">\(d=2\)</span>) quadrilaterals, as if it were quadratic (2x2 Gauss points instead of 3x3):</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">u</span> <span class="o">=</span> <span class="n">Function</span><span class="p">(</span><span class="n">V</span><span class="p">)</span>
<span class="n">u_</span> <span class="o">=</span> <span class="n">TestFunction</span><span class="p">(</span><span class="n">V</span><span class="p">)</span>
<span class="n">du</span> <span class="o">=</span> <span class="n">TrialFunction</span><span class="p">(</span><span class="n">V</span><span class="p">)</span>

<span class="n">dx_shear</span> <span class="o">=</span> <span class="n">dx</span><span class="p">(</span><span class="n">metadata</span><span class="o">=</span><span class="p">{</span><span class="s2">&quot;quadrature_degree&quot;</span><span class="p">:</span> <span class="mi">2</span><span class="o">*</span><span class="n">deg</span><span class="o">-</span><span class="mi">2</span><span class="p">})</span>

<span class="n">L</span> <span class="o">=</span> <span class="n">f</span><span class="o">*</span><span class="n">u_</span><span class="p">[</span><span class="mi">0</span><span class="p">]</span><span class="o">*</span><span class="n">dx</span>
<span class="n">a</span> <span class="o">=</span> <span class="n">inner</span><span class="p">(</span><span class="n">bending_moment</span><span class="p">(</span><span class="n">u_</span><span class="p">),</span> <span class="n">curv</span><span class="p">(</span><span class="n">du</span><span class="p">))</span><span class="o">*</span><span class="n">dx</span> <span class="o">+</span> <span class="n">dot</span><span class="p">(</span><span class="n">shear_force</span><span class="p">(</span><span class="n">u_</span><span class="p">),</span> <span class="n">shear_strain</span><span class="p">(</span><span class="n">du</span><span class="p">))</span><span class="o">*</span><span class="n">dx_shear</span>
</pre></div>
</div>
<p>We then solve for the solution and export the relevant fields to XDMF files</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="n">solve</span><span class="p">(</span><span class="n">a</span> <span class="o">==</span> <span class="n">L</span><span class="p">,</span> <span class="n">u</span><span class="p">,</span> <span class="n">bc</span><span class="p">)</span>

<span class="p">(</span><span class="n">w</span><span class="p">,</span> <span class="n">theta</span><span class="p">)</span> <span class="o">=</span> <span class="n">split</span><span class="p">(</span><span class="n">u</span><span class="p">)</span>

<span class="n">Vw</span> <span class="o">=</span> <span class="n">FunctionSpace</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">We</span><span class="p">)</span>
<span class="n">Vt</span> <span class="o">=</span> <span class="n">FunctionSpace</span><span class="p">(</span><span class="n">mesh</span><span class="p">,</span> <span class="n">Te</span><span class="p">)</span>
<span class="n">ww</span> <span class="o">=</span> <span class="n">Function</span><span class="p">(</span><span class="n">Vw</span><span class="p">,</span> <span class="n">name</span><span class="o">=</span><span class="s2">&quot;Deflection&quot;</span><span class="p">)</span>
<span class="n">tt</span> <span class="o">=</span> <span class="n">Function</span><span class="p">(</span><span class="n">Vt</span><span class="p">,</span> <span class="n">name</span><span class="o">=</span><span class="s2">&quot;Rotation&quot;</span><span class="p">)</span>
<span class="n">ww</span><span class="o">.</span><span class="n">assign</span><span class="p">(</span><span class="n">project</span><span class="p">(</span><span class="n">w</span><span class="p">,</span> <span class="n">Vw</span><span class="p">))</span>
<span class="n">tt</span><span class="o">.</span><span class="n">assign</span><span class="p">(</span><span class="n">project</span><span class="p">(</span><span class="n">theta</span><span class="p">,</span> <span class="n">Vt</span><span class="p">))</span>

<span class="n">file_results</span> <span class="o">=</span> <span class="n">XDMFFile</span><span class="p">(</span><span class="s2">&quot;RM_results.xdmf&quot;</span><span class="p">)</span>
<span class="n">file_results</span><span class="o">.</span><span class="n">parameters</span><span class="p">[</span><span class="s2">&quot;flush_output&quot;</span><span class="p">]</span> <span class="o">=</span> <span class="kc">True</span>
<span class="n">file_results</span><span class="o">.</span><span class="n">parameters</span><span class="p">[</span><span class="s2">&quot;functions_share_mesh&quot;</span><span class="p">]</span> <span class="o">=</span> <span class="kc">True</span>
<span class="n">file_results</span><span class="o">.</span><span class="n">write</span><span class="p">(</span><span class="n">ww</span><span class="p">,</span> <span class="mf">0.</span><span class="p">)</span>
<span class="n">file_results</span><span class="o">.</span><span class="n">write</span><span class="p">(</span><span class="n">tt</span><span class="p">,</span> <span class="mf">0.</span><span class="p">)</span>
</pre></div>
</div>
<p>The solution is compared to the Kirchhoff analytical solution:</p>
<div class="highlight-default"><div class="highlight"><pre><span></span><span class="nb">print</span><span class="p">(</span><span class="s2">&quot;Kirchhoff deflection:&quot;</span><span class="p">,</span> <span class="o">-</span><span class="mf">1.265319087e-3</span><span class="o">*</span><span class="nb">float</span><span class="p">(</span><span class="n">f</span><span class="o">/</span><span class="n">D</span><span class="p">))</span>
<span class="nb">print</span><span class="p">(</span><span class="s2">&quot;Reissner-Mindlin FE deflection:&quot;</span><span class="p">,</span> <span class="o">-</span><span class="nb">min</span><span class="p">(</span><span class="n">ww</span><span class="o">.</span><span class="n">vector</span><span class="p">()</span><span class="o">.</span><span class="n">get_local</span><span class="p">()))</span> <span class="c1"># point evaluation for quads</span>
                                                                        <span class="c1"># is not implemented in fenics 2017.2</span>
</pre></div>
</div>
<p>For <span class="math">\(h=0.001\)</span> and 50 quads per side, one finds <span class="math">\(w_{FE} = 1.38182\text{e-5}\)</span> for linear quads
and <span class="math">\(w_{FE} = 1.38176\text{e-5}\)</span> for quadratic quads against <span class="math">\(w_{\text{Kirchhoff}} = 1.38173\text{e-5}\)</span> for
the thin plate solution.</p>
</div>
</div>


          </div>
        </div>
      </div>
        </div>
        <div class="sidebar">
          <h3>Table Of Contents</h3>
          <ul class="current">
<li class="toctree-l1"><a class="reference internal" href="../../intro.html">Introduction</a><ul>
<li class="toctree-l2"><a class="reference internal" href="../../intro.html#what-is-it-about">What is it about ?</a></li>
<li class="toctree-l2"><a class="reference internal" href="../../intro.html#how-do-i-get-started">How do I get started ?</a></li>
<li class="toctree-l2"><a class="reference internal" href="../../intro.html#about-the-author">About the author</a></li>
</ul>
</li>
<li class="toctree-l1"><a class="reference internal" href="../../linear_problems.html">Linear problems in solid mechanics</a><ul>
<li class="toctree-l2"><a class="reference internal" href="../elasticity/2D_elasticity.py.html">2D linear elasticity</a></li>
<li class="toctree-l2"><a class="reference internal" href="../elasticity/orthotropic_elasticity.py.html">Orthotropic linear elasticity</a></li>
<li class="toctree-l2"><a class="reference internal" href="../elasticity/axisymmetric_elasticity.html">Axisymmetric formulation for elastic structures of revolution</a></li>
<li class="toctree-l2"><a class="reference internal" href="../thermoelasticity/thermoelasticity.html">Linear thermoelasticity (weak coupling)</a></li>
<li class="toctree-l2"><a class="reference internal" href="../thermoelasticity/thermoelasticity_transient.html">Thermo-elastic evolution problem (full coupling)</a></li>
<li class="toctree-l2"><a class="reference internal" href="../modal_analysis_dynamics/cantilever_modal.py.html">Modal analysis of an elastic structure</a></li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
202
<li class="toctree-l2"><a class="reference internal" href="../elastodynamics/demo_elastodynamics.py.html">Time-integration of elastodynamics equation</a></li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
203 204 205 206 207 208 209
</ul>
</li>
<li class="toctree-l1"><a class="reference internal" href="../../homogenization.html">Homogenization of heterogeneous materials</a><ul>
<li class="toctree-l2"><a class="reference internal" href="../periodic_homog_elas/periodic_homog_elas.html">Periodic homogenization of linear elastic materials</a></li>
</ul>
</li>
<li class="toctree-l1"><a class="reference internal" href="../../nonlinear_problems.html">Nonlinear problems in solid mechanics</a><ul>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
210
<li class="toctree-l2"><a class="reference internal" href="../viscoelasticity/linear_viscoelasticity.html">Linear viscoelasticity</a></li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
211
<li class="toctree-l2"><a class="reference internal" href="../2D_plasticity/vonMises_plasticity.py.html">Elasto-plastic analysis of a 2D von Mises material</a></li>
212
<li class="toctree-l2"><a class="reference internal" href="../contact/penalty.html">Hertzian contact with a rigid indenter using a penalty approach</a></li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
213 214
</ul>
</li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
215 216
<li class="toctree-l1 current"><a class="reference internal" href="../../beams_and_plates.html">Beams and plates</a><ul class="current">
<li class="toctree-l2"><a class="reference internal" href="../timoshenko/beam_buckling.html">Eulerian buckling of a beam</a></li>
Jeremy BLEYER's avatar
Jeremy BLEYER committed
217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244
<li class="toctree-l2 current"><a class="current reference internal" href="#">Reissner-Mindlin plate with Quadrilaterals</a></li>
<li class="toctree-l2"><a class="reference internal" href="reissner_mindlin_dg.py.html">Reissner-Mindlin plate with a Discontinuous-Galerkin approach</a></li>
</ul>
</li>
<li class="toctree-l1"><a class="reference internal" href="../../tips_and_tricks.html">Tips and Tricks</a><ul>
<li class="toctree-l2"><a class="reference internal" href="../../tips_and_tricks.html#efficient-projection-on-dg-or-quadrature-spaces">Efficient projection on DG or Quadrature spaces</a></li>
</ul>
</li>
</ul>

          <div role="search">
            <h3 style="margin-top: 1.5em;">Search</h3>
            <form class="search" action="../../search.html" method="get">
                <input type="text" name="q" />
                <input type="submit" value="Go" />
                <input type="hidden" name="check_keywords" value="yes" />
                <input type="hidden" name="area" value="default" />
            </form>
          </div>
        </div>
        <div class="clearer"></div>
      </div>
    </div>

    <div class="footer-wrapper">
      <div class="footer">
        <div class="left">
          <div role="navigation" aria-label="related navigaton">
Jeremy BLEYER's avatar
Jeremy BLEYER committed
245
            <a href="../timoshenko/beam_buckling.html" title="Eulerian buckling of a beam"
Jeremy BLEYER's avatar
Jeremy BLEYER committed
246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271
              >previous</a> |
            <a href="reissner_mindlin_dg.py.html" title="Reissner-Mindlin plate with a Discontinuous-Galerkin approach"
              >next</a> |
            <a href="../../genindex.html" title="General Index"
              >index</a>
          </div>
          <div role="note" aria-label="source link">
              <br/>
              <a href="../../_sources/demo/reissner_mindlin/reissner_mindlin_quads.py.rst.txt"
                rel="nofollow">Show Source</a>
          </div>
        </div>

        <div class="right">
          
    <div class="footer" role="contentinfo">
        &#169; Copyright 2016, Jeremy Bleyer.
      Created using <a href="http://sphinx-doc.org/">Sphinx</a> 1.6.5.
    </div>
        </div>
        <div class="clearer"></div>
      </div>
    </div>

  </body>
</html>