{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 3D Code for isothermal finite hyperelasticity\n", "\n", "Simple stretch\n", "\n", "Basic units:\n", "Length: mm\n", "Mass: kg\n", "Time: s\n", "Derived units:\n", "Force: milliNewtons\n", "Stress: kPa \n", "\n", "Eric Stewart and Lallit Anand \n", "ericstew@mit.edu and anand@mit.edu \n", "\n", "Converted to FEniCSx by Jorge Nin\n", "jorgenin@mit.edu\n", "September 2023\n", "\n" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import dolfinx\n", "\n", "from mpi4py import MPI\n", "from petsc4py import PETSc\n", "\n", "from dolfinx import fem, mesh, io, plot, log\n", "from dolfinx.fem import (Constant, dirichletbc, Function, functionspace, Expression )\n", "from dolfinx.fem.petsc import NonlinearProblem\n", "from dolfinx.nls.petsc import NewtonSolver\n", "from dolfinx.io import XDMFFile\n", "import ufl\n", "from ufl import (TestFunctions, TrialFunction, Identity, grad, det, div, dev, inv, tr, sqrt, conditional , gt, dx, inner, derivative, dot, ln, split)\n", "from datetime import datetime\n", "from basix.ufl import element,mixed_element\n", "from dolfinx.plot import vtk_mesh\n", "\n", "import pyvista\n", "pyvista.set_jupyter_backend('client')\n", "## Define temporal parameters" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Set level of detail for log messages (integer)\n", "Guide: \\\n", "CRITICAL = 50 errors that may lead to data corruption \\\n", "ERROR = 40 things that HAVE gone wrong \\\n", "WARNING = 30 things that MAY go wrong later \\\n", "INFO = 20 information of general interest (includes solver info) \\\n", "PROGRESS = 16 what's happening (broadly) \\\n", "TRACE = 13 what's happening (in detail) \\\n", "DBG = 10 sundry \n" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "log.set_log_level(log.LogLevel.WARNING)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# Define Geometry" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "length = 10.0 # mm\n", "domain = mesh.create_box(MPI.COMM_WORLD,[[0.0,0.0,0.0],[length,length,length]],[2,2,2],mesh.CellType.tetrahedron)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Visualize Gemometry" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "975f7f644a384df1b14a6ce11ee5c243", "version_major": 2, "version_minor": 0 }, "text/plain": [ "Widget(value='