# Point load and body force

**URL:** <https://fenicsproject.discourse.group/t/point-load-and-body-force/2556>\
**Category:** variational formulation\
**Created:** [March 3, 2020, 3:08pm UTC](https://fenicsproject.discourse.group/t/point-load-and-body-force/2556 "2020-03-03T15:08:40Z")\
**Posts on this page:** 1\
**Showing post:** 7

<div class="post-metadata">

**Author:** ![amy](https://avatars.discourse-cdn.com/v4/letter/a/4da419/32.png) [@amy](https://fenicsproject.discourse.group/u/amy)\
**Post date:** [March 3, 2020, 11:28pm UTC](https://fenicsproject.discourse.group/t/point-load-and-body-force/2556/7 "2020-03-03T23:28:30Z")

</div>

I am running this code now, and getting the error above

from dolfin import \*  
import mshr  
import numpy as np  
from **future** import print\_function  
from fenics import \*  
from **future** import print\_function  
from dolfin import \*  
from **future** import division  
import matplotlib.pyplot as plt  
from ufl import nabla\_div

from mshr import \*  
from dolfin import mesh as module

# Load mesh and define function space

L = 8.  
W = 1.

n0 = int(L_6)  
n1 = int(W_6)  
n2 = int(W\*6)

# Mesh for cantilever beam

Omega = BoxMesh(Point(0, 0, 0), Point(L, W, W), 10, 3, 3)

# = Box(0., 0., 0., d0, d1, d2, n0, n1, n2)

Omega.init()

V = VectorFunctionSpace(Omega, “CG”, 1)

tol = 1.e-8 # tolerance

# Left boundary

class Left(SubDomain):  
def inside(self, x, on\_boundary):  
return near(x[0], 0.0)

# Right boundary

class Right(SubDomain):  
def inside(self, x, on\_boundary):  
return near(x[0], 8.0)

# Point (8,0.5,0.5) at which load is added

class point(SubDomain):  
def inside(self, x, on\_boundary):  
return near(x[0],8.) and near(x[1],0.5,1e-2) and near(x[2],0.5,1e-2)

#Boundary segments  
left = Left()  
pt = point()  
right = Right()

# Initialize mesh function for boundary domains

#boundaries = MeshFunction(“size\_t”, Omega, file)#“size\_t”, mesh, 2, mesh.domains()  
#boundaries.set\_all(0)  
#left.mark(boundaries, 1)  
#pt.mark(boundaries, 2)  
boundaries = MeshFunction(“size\_t”, Omega,0) #mesh omega yerinde kalmisti  
#boundaries = CellFunction(‘size\_t’, mesh)  
boundaries.set\_all(0)  
left.mark(boundaries, 1)  
pt.mark(boundaries, 2)  
#all\_boundary.mark(boundaries, 1)

# Define variational problem

f = Constant((0.0, 0.0, 00.0))  
E = 210000.0  
nu = 0.3

mu, lmbda = Constant(E/(2\*(1 + nu))), Constant(E_nu/((1 + nu)_(1 - 2\*nu)))

# Strain

def eps(v):  
return sym(grad(v))  
model = “plane\_stress”

mu = E/2/(1+nu)  
lmbda = E_nu/(1+nu)/(1-2_nu)  
if model == “plane\_stress”:  
lmbda = 2_mu_lmbda/(lmbda+2\*mu)

def sigma(v):  
return lmbda_tr(eps(v))Identity(3) + 2.0mu_eps(v)

ds = ds(subdomain\_data=boundaries)

# Define variational problem

ds= Measure(“ds”)(subdomain\_data=boundaries)  
g\_T= Constant((0.0, 0.0, -0.3))  
du = TrialFunction(V)  
u\_ = TestFunction(V)  
a = inner(sigma(du), eps(u\_))\*dx  
left\_end\_displacement = Constant((0.0, 0.0, 0.0))

#bc = [DirichletBC(V, left\_end\_displacement, left)]  
bc = [DirichletBC(V, left\_end\_displacement, left),DirichletBC(V, g\_T, pt, method=“pointwise”)]  
l= dot(f, v)\*dx +dot(v,g\_T)\*ds(2)  
u = Function(V)  
solve(a ==l, u, bc);

---

_[View the full topic](https://fenicsproject.discourse.group/t/point-load-and-body-force/2556)._