Dirichletbc on only one coordinate

variational formulation
1

Hello,
I want to apply the following boundary condition on a 2D face : x[0]=0 and no constrain on x[1] and x[2].
I wanna do it on this hyperelasticity example, on the right face.

If i do :

r1=Expression(("0","x[1]","x[2]"),degree=2)

it doesn’t work.
Here is my code :

from dolfin import *
import matplotlib.pyplot as plt
import numpy as np

E=1e6
nu=0.3

#optimisation de la compilation
parameters["form_compiler"]["cpp_optimize"] = True
parameters["form_compiler"]["representation"] = "uflacs"


# import mesh
mesh = BoxMesh(Point(0,0,0),Point(1,1,0.1),20,20,2)

# 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)

class PeriodicBoundary(SubDomain):
    def __init__(self, tolerance=DOLFIN_EPS):
        """ vertices stores the coordinates of the 4 unit cell corners"""
        SubDomain.__init__(self, tolerance)
        self.tol = tolerance
        
    def inside(self, x, on_boundary):
        # return True if on bottom boundary
        return bool((near(x[1], 0, self.tol)) and on_boundary)

    def map(self, x, y):
        y[0] = x[0]
        y[1] = x[1] - 1
        y[2] = x[2]


PBC = PeriodicBoundary(tolerance=1e-8)
V = VectorFunctionSpace(mesh, "Lagrange", 1,constrained_domain=PBC)

v = TestFunction(V)
du = TrialFunction(V)

u = Function(V)                 # Displacement from previous iteration
B  = Constant((0.0, 0.0, -100000.0))  # Body force per unit volume
T  = Constant((0.0,  0.0, 0.0))  # Traction force on the boundary


r0= Constant((0.,0.,0.))
r1= Constant((0.,0.,0.))


bc_L = DirichletBC(V,r0,left)
bc_R = DirichletBC(V,r1,right)

bcs =[bc_L,bc_R]


# 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)

# Lame parameters
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

# 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)

# Solve variational problem

problem = NonlinearVariationalProblem(F, u, bcs, J)
solver  = NonlinearVariationalSolver(problem)
prm = solver.parameters
prm['newton_solver']['absolute_tolerance'] = 1E-8
prm['newton_solver']['relative_tolerance'] = 1E-7
prm['newton_solver']['maximum_iterations'] = 20
prm['newton_solver']['relaxation_parameter'] = 1.0
solver.solve()

# Save solution in VTK format

file = File("displacement.pvd");
file << u;

Thanks

2

Use V.sub(0) in your dirichlet condition, as for instance done in Boundary conditions on edges or nodes - #6 by dokken

1 Like