Combined boundary conditions

hi all,

im creating a rod that has a combined boundary condition of neumann and direchlet.
usually when setting both of the boundary conditions, do i need to use a different solver that solve()? or i can use the same solve as i would use with only one kind of boundary condition ?

Do you mean that a part of the boundary has a Dirichlet condition, and another part of the boundary has Neumann conditions? That’s perfectly fine for solve() (even though I am not sure exactly to which function you are referring to), since Neumann boundary conditions result in an additional integral on the right-hand side term.

OR:

Do you mean that the same part of the boundary will have both a Dirichlet and a Neumann boundary condition? To my knowledge, that will typically result in an ill posed problem.

from dolfin import *
import numpy as np
import fenics as fe

#continue the code

class Omega_0(SubDomain):
    def inside(self, x, on_boundary):
        return x[0] <= length/2 + DOLFIN_EPS


class Omega_1(SubDomain):
    def inside(self, x, on_boundary):
        return x[0] > length/2 - DOLFIN_EPS


subdomain_0 = Omega_0()
subdomain_1 = Omega_1()
materials = MeshFunction('size_t', mesh, mesh.topology().dim())
materials.set_all(0)
subdomain_0.mark(materials, 1)
subdomain_1.mark(materials, 2)



# Boundary conditions
class Boundaryleft(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and near(x[0], 0, DOLFIN_EPS)


class RightBoundary(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and near(x[0], 1, DOLFIN_EPS)



V = VectorFunctionSpace(mesh, 'Lagrange', 2)
# ADDED

facets = MeshFunction("size_t", mesh, mesh.topology().dim() - 1)
facets.set_all(0)
Boundaryleft().mark(facets, 1)
RightBoundary().mark(facets, 2)

bc1 = DirichletBC(V.sub(0), Constant(0), facets, 1)


# Define measure for boundary integration
# applying neumann boundary condition on the right side
neumann = Expression("100000*x[0]", degree=1)
ds = Measure('ds', domain=mesh, subdomain_data=facets)
dx = Measure('dx', domain=mesh, subdomain_data=materials)
#applying direchla boundary contion on the left side
bc = bc1
 # rest of the code ....

u = TrialFunction(V)
d = u.geometric_dimension()  # space dimension
v = TestFunction(V)
f = Constant((0, 0, 0))
T = Constant((0, 0, 0))


a = inner(sigma(u, E_global, nu), epsilon(v)) * dx
n = FacetNormal(mesh)  # Normal vector to the boundary

L = dot(f, v) * dx + neumann*dot(n, v)*ds(2)


# Compute solution
u_solution = Function(V)
solve(a == L, u_solution, bc)

this is what i did, one facet is neumann and the other is direchlet

Looks reasonable to me, but it can’t be run because the mesh object is missing, and even if I could run it I wouldn’t know what to look for, considering that you haven’t stated what seems to be wrong with it.

this is my mesh

length, width, height = 1.0, 0.2, 0.2  # Dimensions of the rod
num_elements = 12  # Number of elements along each dimension
E1, E2 = 1E5, 0.5E5  # Young's moduli for the two halves of the rod

# Create a 3D mesh for the rod
mesh = BoxMesh(Point(0, 0, 0), Point(length, width, height), num_elements, num_elements, num_elements)

something else i though about, in the solution i dont take Neumann condition under consideration, is that an expected behavior ? or do i need to include it somehow ?

Please provide a minimal reproducible example.
This means that you cannot truncate the code with

and you have to express what is wrong.
Do you get an error message, if so, what is it?
If you do not get an error message, how are you interpreting that the results are wrong?