 # Computing the self-weight of a truss

Hello community,
I have a truss like structure. I am trying to apply the self-weight (-g * rho) to the back cells only. My current implementation only applies force on the hollow areas (white cells) of the section (shown in Figure 1).
Could someone help me with properly applying the self-weight given the density? Fig. 1 The black areas is the actual truss (where the material is and forces should be applied here) and the yellow areas are the hollow part of the truss where the forces are being applied (where there is actually no material)

Here is a working example:

``````from dolfin import *
import torch_fenics
from ufl import nabla_div

def left(x, on_boundary):
return near(x, -2) and on_boundary

def eps(v):

def sigma(v):
return coeff*(lamda*div(v)*Identity(2)+2*mu*eps(v))

# Problem parameters
p = Constant(5)
E = Constant(203 * 1e8) # kg/m2
nu = Constant(0.265)

g = Constant(9.8) # m/s2
rho = Constant(7985) # kg/m3

lamda = E*nu/(1+nu)/(1-2*nu)
mu = E/(2*(1+nu))

# Mesh
mesh = RectangleMesh(Point(-2, 0), Point(2, 1), 20, 10, "crossed")
facets = MeshFunction("size_t", mesh, 1)

V0 = FunctionSpace(mesh, "DG", 0)        # Function space for density field
V2 = VectorFunctionSpace(mesh, "CG", 2)  # Function space for displacement

# Assigning the density
# Peta is a numpy array containing the values of the truss; here i used a random value to have a working example
peta = np.ones((800)) * 0.5
theta = Function(V0, name = 'Density')
theta.vector().set_local(peta)

# Inhomogeneous elastic variational problem
u_ = TestFunction(V2)
du = TrialFunction(V2)

u = Function(V2, name="Displacement")
coeff = theta ** p

f = theta * Constant((0, -rho*g))

a = inner(sigma(u_), eps(du))*dx
L = inner(f, u_)* dx

# Applying clamped end
bc = DirichletBC(V2, Constant((0, 0)), left)

# solve
solve(a == L, u, bc)

# Plotting the displacement
plot(u/1e9, mode="displacement")
``````

Feel free to message me so I can also send you the mesh file if needed.