Pulse boundary condition

variational formulation
1

Hi all,

I would like to impose a pulse function condition on a boundary.
For example, let consider that we have a rectangular box and we want to impose the following pressure condition on the left boundary:

P = 1 on left surface if t = 0
P = 0 on left surface if t > 0

Is there any easy way of doing this?

Thank you very much in advance!

2

Consider the following minimal example:

from dolfin import *

mesh = UnitSquareMesh(10, 10)

t = 0

expr = Expression("t<tol? 1: 0",degree=1, t=t, tol=1e-13)

V = FunctionSpace(mesh, "CG", 1)
u = Function(V)

out = File("pulse.pvd")
bc = DirichletBC(V, expr, "x[0]<DOLFIN_EPS and on_boundary")
bc.apply(u.vector())
out << u
expr.t = 5
bc.apply(u.vector())
out << u
1 Like
3

Thank you so much for the speedy response. But, is there an equivalent command for Expression in dolfinX?

Thanks again!

4

See: Test problem 2: Flow past a cylinder (DFG 2D-3 benchmark) — FEniCSx tutorial

1 Like
5

Thank you very much.
I have followed that tutorial and the minimal example that you provided, and I have managed to write a full working code for dolfinx.

I will leave it bellow in case that anyone else is interested in this problem.

import dolfinx
from dolfinx.io import XDMFFile
from dolfinx.fem import (Constant, DirichletBC, Function, FunctionSpace, apply_lifting, assemble_matrix, 
                         assemble_scalar, assemble_vector, create_vector, locate_dofs_topological, set_bc)
import gmsh
import ufl
from ufl import inner, div, dx
from mpi4py import MPI
import numpy as np
from petsc4py import PETSc

# Mesh
mesh = dolfinx.UnitSquareMesh(MPI.COMM_WORLD, 10, 10)
gdim = 2

# Velocity element space
v_cg2 = ufl.VectorElement('CG', mesh.ufl_cell(), 2)
V = dolfinx.FunctionSpace(mesh, v_cg2)
u = Function(V)

# Define pulse function for velocity in the left boundary
t = 0
class LeftVelocity():
    def __init__(self, t):
        self.t = t

    def __call__(self, x):
        values = np.zeros((gdim, x.shape[1]),dtype=PETSc.ScalarType)
        if t == 0:
            values[0] = 1
        else:
            values[0] = 0
        return values

# Get dofs of the left boundary
def left(x):
    return np.isclose(x[0], 0)
left_dofsp = dolfinx.fem.locate_dofs_geometrical(V, left)

# Define boundary condition
u_left = Function(V)
left_velocity = LeftVelocity(t)
u_left.interpolate(left_velocity)
bc = dolfinx.DirichletBC(u_left, left_dofsp)
dolfinx.fem.set_bc(u.vector, [bc])

# Output file
out = XDMFFile(MPI.COMM_WORLD, "Pulse.xdmf", "w")
out.write_mesh(mesh)
out.write_function(u, t)

dt = 0.1
for i in range(10):
    t += dt

    # Update left boundary condition
    left_velocity.t = t
    u_left.interpolate(left_velocity)
    bc = dolfinx.DirichletBC(u_left, left_dofsp)
    dolfinx.fem.set_bc(u.vector, [bc])

    # Output file
    out.write_function(u, t)

out.close()