Point load and body force

I have read it but I could not be sure because when I try this in the code above, the deformation is comes as 2.439431848880892e-07.
With linear elements I got this value 0.00012453125397051366. I did not understand this, why are they that different? I have read in the document that if the deformation is large, quad elements give more accurate answer.What is the reason, could you explain?

The reason is because you are trying to apply a point load at a place where there are no degrees of freedom. When you apply a pointwise boundary there has to be a dof there. Using:

class point(SubDomain):
    def inside(self, x, on_boundary):
        return near(x[0],8.) and near(x[1],1,1e-5) and near(x[2],1,1e-5)

Will give you similar results for both P1 and P2, as there is always a dof at the corner of the geometry.

Ok thank you, I am trying to do large deformation version of this problem using this demo
(https://fenicsproject.org/docs/dolfin/1.4.0/python/_downloads/demo_hyperelasticity.py)
but I am getting RuntimeError: Unable to compile C++ code with dijitso
Why would it be?

Again, you Need to make a minimal example and supply the full error message.

dolfin 1.4.0 is a very old version. Use the demos for the version of dolfin that you’re running.

1 Like

I have tried the last version but still I am getting an error https://fenicsproject.org/docs/dolfin/latest/python/demos/hyperelasticity/demo_hyperelasticity.py.html

    # .. _demo_hyperelasticity:
    # 
    # Hyperelasticity
    # ===============
    # 
    # This demo is implemented in a single Python file,
    # :download:`demo_hyperelasticity.py`, which contains both the
    # variational forms and the solver.
    # 
    # 
    # Background
    # ----------
    # 
    # This example demonstrates the solution of a three-dimensional
    # elasticity problem. In addition to illustrating how to use
    # FunctionSpaces, Expressions and how to apply Dirichlet boundary
    # conditions, it focuses on how to:
    # 
    # * Minimise a non-quadratic functional
    # * Use automatic computation of the directional derivative
    # * Solve a nonlinear variational problem
    # * Define compiled sub-domains
    # * Use specific form compiler optimization options
    # 
    # Equation and problem definition
    # -------------------------------
    # 
    # By definition, boundary value problems for hyperelastic media can be
    # expressed as minimisation problems, and the minimization approach is
    # adopted in this example. For a domain :math:`\Omega \subset
    # \mathbb{R}^{d}`, where :math:`d` denotes the spatial dimension, the
    # task is to find the displacement field :math:`u: \Omega \rightarrow
    # \mathbb{R}^{d}` that minimises the total potential energy :math:`\Pi`:
    # 
    # .. math::
    #    \min_{u \in V} \Pi,
    # 
    # where :math:`V` is a suitable function space that satisfies boundary
    # conditions on :math:`u`.  The total potential energy is given by
    # 
    # .. math::
    #    \Pi = \int_{\Omega} \psi(u) \, {\rm d} x
    #    - \int_{\Omega} B \cdot u \, {\rm d} x
    #    - \int_{\partial\Omega} T \cdot u \, {\rm d} s,
    # 
    # where :math:`\psi` is the elastic stored energy density, :math:`B` is a
    # body force (per unit reference volume) and :math:`T` is a traction force
    # (per unit reference area).
    # 
    # At minimum points of :math:`\Pi`, the directional derivative of :math:`\Pi`
    # with respect to change in :math:`u`
    # 
    # .. math::
    #    :label: first_variation
    # 
    #    L(u; v) = D_{v} \Pi = \left. \frac{d \Pi(u + \epsilon v)}{d\epsilon} \right|_{\epsilon = 0}
    # 
    # is equal to zero for all :math:`v \in V`:
    # 
    # .. math::
    #    L(u; v) = 0 \quad \forall \ v \in V.
    # 
    # To minimise the potential energy, a solution to the variational
    # equation above is sought. Depending on the potential energy
    # :math:`\psi`, :math:`L(u; v)` can be nonlinear in :math:`u`. In such a
    # case, the Jacobian of :math:`L` is required in order to solve this
    # problem using Newton's method. The Jacobian of :math:`L` is defined as
    # 
    # .. math::
    #    :label: second_variation
    # 
    #    a(u; du, v) = D_{du} L = \left. \frac{d L(u + \epsilon du; v)}{d\epsilon} \right|_{\epsilon = 0} .
    # 
    # 
    # 
    # Elastic stored energy density
    # ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
    # 
    # To define the elastic stored energy density, consider the deformation
    # gradient :math:`F`
    # 
    # .. math::
    # 
    #    F = I + \nabla u,
    # 
    # the right Cauchy-Green tensor :math:`C`
    # 
    # .. math::
    # 
    #    C = F^{T} F,
    # 
    # and the scalars  :math:`J` and :math:`I_{c}`
    # 
    # .. math::
    #    J     &= \det(F), \\
    #    I_{c} &= {\rm trace}(C).
    # 
    # This demo considers a common neo-Hookean stored energy model of the form
    # 
    # .. math::
    #    \psi =  \frac{\mu}{2} (I_{c} - 3) - \mu \ln(J) + \frac{\lambda}{2}\ln(J)^{2},
    # 
    # where :math:`\mu` and :math:`\lambda` are the Lame parameters. These
    # can be expressed in terms of the more common Young's modulus :math:`E`
    # and Poisson ratio :math:`\nu` by:
    # 
    # .. math::
    #     \lambda = \frac{E \nu}{(1 + \nu)(1 - 2\nu)}, \quad  \quad
    #     \mu     =  \frac{E}{2(1 + \nu)} .
    # 
    # 
# Demo parameters
# ^^^^^^^^^^^^^^^
# 
# We consider a unit cube domain:
# 
# * :math:`\Omega = (0, 1) \times (0, 1) \times (0, 1)` (unit cube)
# 
# We use the following definitions of the boundary and boundary conditions:
# 
# * :math:`\Gamma_{D_{0}} = 0 \times (0, 1) \times (0, 1)` (Dirichlet boundary)
# 
# * :math:`\Gamma_{D_{1}} = 1 \times (0, 1) \times (0, 1)` (Dirichlet boundary)
# 
# * :math:`\Gamma_{N} = \partial \Omega \backslash \Gamma_{D}` (Neumann boundary)
# 
# * On :math:`\Gamma_{D_{0}}`:  :math:`u = (0, 0, 0)`
# 
# * On  :math:`\Gamma_{D_{1}}`
#     .. math::
#        u = (&0, \\
#        &(0.5 + (y - 0.5)\cos(\pi/3) - (z - 0.5)\sin(\pi/3) - y)/2, \\
#        &(0.5 + (y - 0.5)\sin(\pi/3) + (z - 0.5)\cos(\pi/3) - z))/2)
# 
# * On :math:`\Gamma_{N}`: :math:`T = (0.1, 0, 0)`
# 
# These are the body forces and material parameters used:
# 
# * :math:`B = (0, -0.5, 0)`
# 
# * :math:`E    = 10.0`
# 
# * :math:`\nu  = 0.3`
# 
# With the above input the solution for :math:`u` will look as follows:
# 
# .. image:: hyperelasticity_u0.png
#     :scale: 75
#     :align: center
# 
# .. image:: hyperelasticity_u1.png
#     :scale: 75
#     :align: center
# 
# 
# Implementation
# --------------
# 
# This demo is implemented in the :download:`demo_hyperelasticity.py`
# file.
# 
# First, the required modules are imported::

import matplotlib.pyplot as plt
from dolfin import *

# The behavior of the form compiler FFC can be adjusted by prescribing
# various parameters. Here, we want to use the UFLACS backend of FFC::

# Optimization options for the form compiler
parameters["form_compiler"]["cpp_optimize"] = True
parameters["form_compiler"]["representation"] = "uflacs"

# The first line tells the form compiler to use C++ compiler optimizations when
# compiling the generated code. The remainder is a dictionary of options which
# will be passed to the form compiler. It lists the optimizations strategies
# that we wish the form compiler to use when generating code.
# 
# .. index:: VectorFunctionSpace
# 
# First, we need a tetrahedral mesh of the domain and a function space
# on this mesh. Here, we choose to create a unit cube mesh with 25 ( =
# 24 + 1) vertices in one direction and 17 ( = 16 + 1) vertices in the
# other two direction. On this mesh, we define a function space of
# continuous piecewise linear vector polynomials (a Lagrange vector
# element space)::

# Create mesh and define function space
mesh = UnitCubeMesh(24, 16, 16)
V = VectorFunctionSpace(mesh, "Lagrange", 1)

# Note that :py:class:`VectorFunctionSpace
# <dolfin.functions.functionspace.VectorFunctionSpace>` creates a
# function space of vector fields. The dimension of the vector field
# (the number of components) is assumed to be the same as the spatial
# dimension, unless otherwise specified.
# 
# .. index:: compiled subdomain
# 
# The portions of the boundary on which Dirichlet boundary conditions
# will be applied are now defined::

# Mark boundary subdomians
left =  CompiledSubDomain("near(x[0], side) && on_boundary", side = 0.0)
right = CompiledSubDomain("near(x[0], side) && on_boundary", side = 1.0)

# The boundary subdomain ``left`` corresponds to the part of the
# boundary on which :math:`x=0` and the boundary subdomain ``right``
# corresponds to the part of the boundary on which :math:`x=1`. Note
# that C++ syntax is used in the :py:func:`CompiledSubDomain`
# <dolfin.compilemodules.subdomains.CompiledSubDomain>` function since
# the function will be automatically compiled into C++ code for
# efficiency. The (built-in) variable ``on_boundary`` is true for points
# on the boundary of a domain, and false otherwise.
# 
# .. index:: compiled expression
# 
# The Dirichlet boundary values are defined using compiled expressions::

# Define Dirichlet boundary (x = 0 or x = 1)
c = Constant((0.0, 0.0, 0.0))
r = Expression(("scale*0.0",
                "scale*(y0 + (x[1] - y0)*cos(theta) - (x[2] - z0)*sin(theta) - x[1])",
                "scale*(z0 + (x[1] - y0)*sin(theta) + (x[2] - z0)*cos(theta) - x[2])"),
                scale = 0.5, y0 = 0.5, z0 = 0.5, theta = pi/3, degree=2)

# Note the use of setting named parameters in the :py:class:`Expression
# <dolfin.functions.expression.Expression>` for ``r``.
# 
# The boundary subdomains and the boundary condition expressions are
# collected together in two :py:class:`DirichletBC
# <dolfin.fem.bcs.DirichletBC>` objects, one for each part of the
# Dirichlet boundary::

bcl = DirichletBC(V, c, left)
bcr = DirichletBC(V, r, right)
bcs = [bcl, bcr]

# The Dirichlet (essential) boundary conditions are constraints on the
# function space :math:`V`. The function space is therefore required as
# an argument to :py:class:`DirichletBC <dolfin.fem.bcs.DirichletBC>`.
# 
# .. index:: TestFunction, TrialFunction, Constant
# 
# Trial and test functions, and the most recent approximate displacement
# ``u`` are defined on the finite element space ``V``, and two objects
# of type :py:class:`Constant <dolfin.functions.constant.Constant>` are
# declared for the body force (``B``) and traction (``T``) terms::

# Define functions
du = TrialFunction(V)            # Incremental displacement
v  = TestFunction(V)             # Test function
u  = Function(V)                 # Displacement from previous iteration
B  = Constant((0.0, -0.5, 0.0))  # Body force per unit volume
T  = Constant((0.1,  0.0, 0.0))  # Traction force on the boundary

# In place of :py:class:`Constant <dolfin.functions.constant.Constant>`,
# it is also possible to use ``as_vector``, e.g.  ``B = as_vector( [0.0,
# -0.5, 0.0] )``. The advantage of Constant is that its values can be
# changed without requiring re-generation and re-compilation of C++
# code. On the other hand, using ``as_vector`` can eliminate some
# function calls during assembly.
# 
# With the functions defined, the kinematic quantities involved in the model
# are defined using UFL syntax::

# Kinematics
d = len(u)
I = Identity(d)             # Identity tensor
F = I + grad(u)             # Deformation gradient
C = F.T*F                   # Right Cauchy-Green tensor

# Invariants of deformation tensors
Ic = tr(C)
J  = det(F)

# Next, the material parameters are set and the strain energy density
# and the total potential energy are defined, again using UFL syntax::

# Elasticity parameters
E, nu = 10.0, 0.3
mu, lmbda = Constant(E/(2*(1 + nu))), Constant(E*nu/((1 + nu)*(1 - 2*nu)))

# Stored strain energy density (compressible neo-Hookean model)
psi = (mu/2)*(Ic - 3) - mu*ln(J) + (lmbda/2)*(ln(J))**2

# Total potential energy
Pi = psi*dx - dot(B, u)*dx - dot(T, u)*ds

# Just as for the body force and traction vectors, :py:class:`Constant
# <dolfin.functions.constant.Constant>` has been used for the model
# parameters ``mu`` and ``lmbda`` to avoid re-generation of C++ code
# when changing model parameters. Note that ``lambda`` is a reserved
# keyword in Python, hence the misspelling ``lmbda``.
# 
# .. index:: directional derivative; derivative, taking variations; derivative, automatic differentiation; derivative
# 
# Directional derivatives are now computed of :math:`\Pi` and :math:`L`
# (see :eq:`first_variation` and :eq:`second_variation`)::

# Compute first variation of Pi (directional derivative about u in the direction of v)
F = derivative(Pi, u, v)

# Compute Jacobian of F
J = derivative(F, u, du)

# The complete variational problem can now be solved by a single call to
# :py:func:`solve <dolfin.fem.solving.solve>`::

# Solve variational problem
solve(F == 0, u, bcs, J=J)

# Finally, the solution ``u`` is saved to a file named
# ``displacement.pvd`` in VTK format, and the deformed mesh is plotted
# to the screen::

# Save solution in VTK format
file = File("displacement.pvd");
file << u;

# Plot solution
plot(u)
plt.show()

error-------------------------
Moving new file over differing existing file:
src: /home/fenics/jitfailure-dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2/error.log.74ff57ab73434456b2197958b6e8d895
dst: /home/fenics/jitfailure-dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2/error.log
backup: /home/fenics/jitfailure-dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2/error.log.old
------------------- Start compiler output ------------------------
/tmp/tmpobjiduqo/dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2.cpp:13:10: fatal error: dolfin/common/Array.h: No such file or directory
#include <dolfin/common/Array.h>
^~~~~~~~~~~~~~~~~~~~~~~
compilation terminated.

------------------- End compiler output ------------------------
Compilation failed! Sources, command, and errors have been written to: /home/fenics/jitfailure-dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2

DijitsoError Traceback (most recent call last)
/usr/lib/python3/dist-packages/dolfin/jit/jit.py in compile_class(cpp_data, mpi_comm)
166 generate=cpp_data[‘jit_generate’],
–> 167 mpi_comm=mpi_comm)
168 submodule = dijitso.extract_factory_function(module, “create_” + module_name)()

/usr/lib/python3/dist-packages/dolfin/jit/jit.py in mpi_jit(*args, **kwargs)
46 if MPI.size(mpi_comm) == 1:
—> 47 return local_jit(*args, **kwargs)
48

/usr/lib/python3/dist-packages/dolfin/jit/jit.py in dijitso_jit(*args, **kwargs)
102 def dijitso_jit(*args, **kwargs):
–> 103 return dijitso.jit(*args, **kwargs)
104

/home/fenics/.local/lib/python3.6/site-packages/dijitso/jit.py in jit(jitable, name, params, generate, send, receive, wait)
216 raise DijitsoError(“Dijitso JIT compilation failed, see ‘%s’ for details”
–> 217 % err_info[‘fail_dir’], err_info)
218

DijitsoError: Dijitso JIT compilation failed, see ‘/home/fenics/jitfailure-dolfin_subdomain_daf8754d301c8cd9a0f303f23d7155a2’ for details

During handling of the above exception, another exception occurred:

RuntimeError Traceback (most recent call last)
in ()
203
204 # Mark boundary subdomians
–> 205 left = CompiledSubDomain(“near(x[0], side) && on_boundary”, side = 0.0)
206 right = CompiledSubDomain(“near(x[0], side) && on_boundary”, side = 1.0)
207

/usr/lib/python3/dist-packages/dolfin/mesh/subdomain.py in new(cls, inside_code, **kwargs)
123 def new(cls, inside_code, **kwargs):
124 properties = kwargs
–> 125 return compile_subdomain(inside_code, properties)

/usr/lib/python3/dist-packages/dolfin/mesh/subdomain.py in compile_subdomain(inside_code, properties)
116 ‘name’: ‘subdomain’, ‘jit_generate’: jit_generate}
117
–> 118 subdomain = compile_class(cpp_data, mpi_comm=mpi_comm)
119 return subdomain
120

/usr/lib/python3/dist-packages/dolfin/jit/jit.py in compile_class(cpp_data, mpi_comm)
168 submodule = dijitso.extract_factory_function(module, “create_” + module_name)()
169 except Exception:
–> 170 raise RuntimeError(“Unable to compile C++ code with dijitso”)
171
172 if name == ‘expression’:

RuntimeError: Unable to compile C++ code with dijitso

I cannot reproduce this error. Running it with docker does not cause this error:
docker run -ti -v $(pwd):/home/fenics/shared -w /home/fenics/shared/ --rm quay.io/fenicsproject/dev:latest

ok I could not install docker. I am using windows 10 home and the one I downloaded requires windows 10 pro or student? Is there a way to install it to my PC?

See for instance: https://forums.docker.com/t/installing-docker-on-windows-10-home/11722/67
Could you Also tell me Which version of dolfin you are using. Can be obtained with

import dolfin
print(dolfin.__version__)

ok I will try,this is my version 2019.1.0

I still can’t reproduce your issue using docker:

docker run -ti -v $(pwd):/home/fenics/shared -w /home/fenics/shared/ --rm quay.io/fenicsproject/stable:latest

Could you try to run dijitso clean to clear your cache?

How did you install dolfin?

for the installation I followed the instructions on fenics.org. I don’t remember the details much
I run this command on the terminal it says
docker:invalid reference format

is it certain that it is elated with docker?

If you dont Remember how you installed it I cant really help you any further. You keep on posting half of an error message and not all of the commands you run, so there is no way to reproduce or debug the behavior.

Docker or anaconda are the simplest ways of running FEniCS, and i would strongly suggest using either.

Did you mean how did I install docker or dolfin earlier? This is how I installed fenics and Ubuntu, Anconda but as I said docker asks for windows 10 pro or student .
https://fenicsproject.org/download/

In general your messages does not contain enough information for me to help you.

For your installation of FEniCS, there Seems to be some files missing. Since you are using Windows, i would install ubuntu through Windows subsystem, and run docker on this subsystem.

I have already Ubuntu in my PC, do you think that my problem is gonna be solved if I have docker in my PC. I followed this https://fenicsproject.org/download/ for running fenics . Docker does not work in my PC since I have windows 10 home. If I upgrade to windows 10 pro docker would work. What do you think? I am not sure also if upgrading will affect the existing files?
Thanks

If you have ubuntu, you should be able to install docker there and use the images described further up

docker does not work on windows10 home as far as I read,

But you are running ubuntu aren’t you? In ubuntu docker works…