diff --git a/scripts/run_unit_tests.sh b/scripts/run_unit_tests.sh
deleted file mode 100755
index 2727ccb..0000000
--- a/scripts/run_unit_tests.sh
+++ /dev/null
@@ -1 +0,0 @@
-HIP_VISIBLE_DEVICES=0 julia +1.12 --project=@. -e 'using Pkg; Pkg.test()'
diff --git a/src/FiniteElementContainers.jl b/src/FiniteElementContainers.jl
index 87bb075..bfad648 100644
--- a/src/FiniteElementContainers.jl
+++ b/src/FiniteElementContainers.jl
@@ -72,11 +72,11 @@ export update_field_ics!
 export update_ic_values!
 
 # Integrals
-export MatrixIntegral
-export ScalarIntegral
-export VectorIntegral
-export integrate
-export remove_fixed_dofs!
+# export MatrixIntegral
+# export ScalarIntegral
+# export VectorIntegral
+# export integrate
+# export remove_fixed_dofs!
 
 # Integrators
 export AbstractIntegrator
diff --git a/src/assemblers/Constraints.jl b/src/assemblers/Constraints.jl
index 664d8cc..a52a41b 100644
--- a/src/assemblers/Constraints.jl
+++ b/src/assemblers/Constraints.jl
@@ -135,7 +135,7 @@ function _adjust_vector_entries_for_constraints!(b, constraint_storage, ::KA.CPU
     # modify b => (I - G) * b + (Gu - g)
     # but Gu = g, so we don't need that here
     # unless we want to modify this to support weakly
-    # enforced BCs later
+    # enforced BCs later TODO
     for i in 1:length(constraint_storage)
         @inbounds b[i] = (1. - constraint_storage[i]) * b[i]
     end
diff --git a/src/bcs/DirichletBCs.jl b/src/bcs/DirichletBCs.jl
index f9fbcb9..77556d2 100644
--- a/src/bcs/DirichletBCs.jl
+++ b/src/bcs/DirichletBCs.jl
@@ -1,5 +1,7 @@
 abstract type AbstractDirichletBC{F} <: AbstractBC{F} end
 
+const EntityName = Union{Nothing, Symbol}
+
 """
 $(TYPEDEF)
 $(TYPEDSIGNATURES)
@@ -8,7 +10,9 @@ User facing API to define a ```DirichletBC````.
 """
 struct DirichletBC{F} <: AbstractDirichletBC{F}
   func::F
-  sset_name::Symbol
+  block_name::EntityName
+  nset_name::EntityName
+  sset_name::EntityName
   var_name::Symbol
 end
 
@@ -17,8 +21,13 @@ $(TYPEDEF)
 $(TYPEDSIGNATURES)
 $(TYPEDFIELDS)
 """
-function DirichletBC(var_name::Symbol, sset_name::Symbol, func::Function)
-  return DirichletBC(func, sset_name, var_name)
+function DirichletBC(
+  var_name::Symbol, func::Function;
+  block_name::EntityName = nothing,
+  nodeset_name::EntityName = nothing,
+  sideset_name::EntityName = nothing
+)
+  return DirichletBC(func, block_name, nodeset_name, sideset_name, var_name)
 end
 
 """
@@ -26,8 +35,25 @@ $(TYPEDEF)
 $(TYPEDSIGNATURES)
 $(TYPEDFIELDS)
 """
-function DirichletBC(var_name::String, sset_name::String, func::Function)
-  return DirichletBC(Symbol(var_name), Symbol(sset_name), func)
+function DirichletBC(
+  var_name::String, func::Function;
+  block_name::Union{Nothing, String} = nothing,
+  nodeset_name::Union{Nothing, String} = nothing,
+  sideset_name::Union{Nothing, String} = nothing
+)
+  # return DirichletBC(Symbol(var_name), Symbol(sset_name), func)
+  if block_name !== nothing
+    block_name = Symbol(block_name)
+  end
+  
+  if nodeset_name !== nothing
+    nodeset_name = Symbol(nodeset_name)
+  end
+
+  if sideset_name !== nothing
+    sideset_name = Symbol(sideset_name)
+  end
+  return DirichletBC(func, block_name, nodeset_name, sideset_name, Symbol(var_name))
 end
 
 """
@@ -53,7 +79,17 @@ $(TYPEDSIGNATURES)
 $(TYPEDFIELDS)
 """
 function DirichletBCContainer(mesh, dof::DofManager, dbc::DirichletBC)
-  bk = BCBookKeeping(mesh, dof, dbc.var_name, sset_name=dbc.sset_name)
+  if dbc.block_name !== nothing
+    bk = BCBookKeeping(mesh, dof, dbc.var_name, block_name=dbc.block_name)
+  elseif dbc.nset_name !== nothing
+    bk = BCBookKeeping(mesh, dof, dbc.var_name, nset_name=dbc.nset_name)
+  elseif dbc.sset_name !== nothing
+    bk = BCBookKeeping(mesh, dof, dbc.var_name, sset_name=dbc.sset_name)
+  else
+    @assert false
+  end
+
+  # bk = BCBookKeeping(mesh, dof, dbc.var_name, sset_name=dbc.sset_name)
 
   # sort nodes and dofs for dirichlet bc
   dof_perm = _unique_sort_perm(bk.dofs)
diff --git a/src/integrals/Integrals.jl b/src/integrals/Integrals.jl
index bc7eb69..140228e 100644
--- a/src/integrals/Integrals.jl
+++ b/src/integrals/Integrals.jl
@@ -12,53 +12,86 @@ function Adapt.adapt_storage(to, int::AbstractIntegral)
     )
 end
 
+function Base.show(io::IO, int::AbstractIntegral)
+    println(io, typeof(int).name.name, ":")
+    println(io, "  Assembler type     = $(typeof(int.assembler))")
+    println(io, "  Cache type         = $(typeof(int.cache))")
+    println(io, "  Integrand function = $(int.integrand)")
+end
+
 KA.get_backend(int::AbstractIntegral) = KA.get_backend(int.cache)
 
+function create_field(int::AbstractIntegral)
+    return create_field(int.assembler)
+end
+
+function create_unknowns(int::AbstractIntegral)
+    return create_unknowns(int.assembler)
+end
+
 function integrate end
 
 function (integral::AbstractIntegral)(U, p)
     return integrate(integral, U, p)
 end
 
-struct ScalarIntegral{A, C <: NamedTuple, I} <: AbstractIntegral{A, C, I}
+# abstract type AbstractIntegral{A, C, I} <: AbstractIntegral{A, C, I} end
+# abstract type AbstractBoundaryIntegral{A, C, I} <: AbstractIntegral{A, C, I} end
+
+abstract type AbstractMatrixIntegral{A, C, I} <: AbstractIntegral{A, C, I} end
+abstract type AbstractScalarIntegral{A, C, I} <: AbstractIntegral{A, C, I} end
+abstract type AbstractVectorIntegral{A, C, I} <: AbstractIntegral{A, C, I} end
+
+struct ScalarCellIntegral{A, C <: NamedTuple, I} <: AbstractIntegral{A, C, I}
     assembler::A
     cache::C
     integrand::I
 end
 
-function ScalarIntegral(asm, integrand)
+function ScalarCellIntegral(asm, integrand)
     cache = create_assembler_cache(asm, AssembledScalar())
-    return ScalarIntegral(asm, cache, integrand)
+    return ScalarCellIntegral(asm, cache, integrand)
 end
 
-function gradient(integral::ScalarIntegral)
-    func(physics, interps, x, t, dt, u, u_n, state_old, state_new, props) = ForwardDiff.gradient(
-        z -> integral.integrand(physics, interps, x, t, dt, z, u_n, state_old, state_new, props),
-        u
-    )
+function gradient(integral::ScalarCellIntegral)
+    func(physics, interps, x, t, dt, u, u_n, state_old, state_new, props) = 
+        ForwardDiff.gradient(
+            z -> integral.integrand(physics, interps, x, t, dt, z, u_n, state_old, state_new, props),
+            u
+        )
     cache = create_field(integral.assembler)
-    return VectorIntegral(integral.assembler, cache, func)
+    return VectorCellIntegral(integral.assembler, cache, func)
+end
+
+function hessian(integral::ScalarCellIntegral)
+    func(physics, interps, x, t, dt, u, u_n, state_old, state_new, props) = 
+        ForwardDiff.hessian(
+            z -> integral.integrand(physics, interps, x, t, dt, z, u_n, state_old, state_new, props),
+            u
+        )
+    cache = create_assembler_cache(asm, AssembledMatrix())
+    return MatrixCellIntegral(integral.assembler, cache, func)
 end
 
-function integrate(integral::ScalarIntegral, U, p)
+function integrate(integral::ScalarCellIntegral, U, p)
     cache, dof = integral.cache, integral.assembler.dof
     func = integral.integrand
     assemble_quadrature_quantity!(cache, dof, nothing, func, U, p)
     return mapreduce(sum, sum, values(cache))
 end
 
-struct MatrixIntegral{A, C <: AbstractVector, I} <: AbstractIntegral{A, C, I}
+struct MatrixCellIntegral{A, C <: AbstractVector, I} <: AbstractIntegral{A, C, I}
     assembler::A
     cache::C
     integrand::I
 end
 
-function MatrixIntegral(asm, integrand)
+function MatrixCellIntegral(asm, integrand)
     cache = create_assembler_cache(asm, AssembledMatrix())
-    return MatrixIntegral(asm, cache, integrand)
+    return MatrixCellIntegral(asm, cache, integrand)
 end
 
-function integrate(integral::MatrixIntegral, U, p)
+function integrate(integral::MatrixCellIntegral, U, p)
     assemble_matrix!(
         integral.cache, 
         integral.assembler.matrix_pattern,
@@ -72,7 +105,7 @@ function integrate(integral::MatrixIntegral, U, p)
     )
 end
 
-function remove_fixed_dofs!(integral::MatrixIntegral)
+function remove_fixed_dofs!(integral::MatrixCellIntegral)
     backend = KA.get_backend(integral)
 
     # TODO change API so we don't need to
@@ -87,18 +120,18 @@ function remove_fixed_dofs!(integral::MatrixIntegral)
     return nothing
 end
 
-struct VectorIntegral{A, C <: AbstractField, I} <: AbstractIntegral{A, C, I}
+struct VectorCellIntegral{A, C <: AbstractField, I} <: AbstractIntegral{A, C, I}
     assembler::A
     cache::C
     integrand::I
 end
 
-function VectorIntegral(asm, integrand)
+function VectorCellIntegral(asm, integrand)
     cache = create_assembler_cache(asm, AssembledVector())
-    return VectorIntegral(asm, cache, integrand)
+    return VectorCellIntegral(asm, cache, integrand)
 end
 
-function integrate(integral::VectorIntegral, U, p)
+function integrate(integral::VectorCellIntegral, U, p)
     assemble_vector!(
         integral.cache, 
         integral.assembler.vector_pattern,
@@ -109,7 +142,7 @@ function integrate(integral::VectorIntegral, U, p)
     return integral.cache
 end
 
-function remove_fixed_dofs!(integral::VectorIntegral)
+function remove_fixed_dofs!(integral::VectorCellIntegral)
     backend = KA.get_backend(integral)
     _adjust_vector_entries_for_constraints!(
         integral.cache, integral.assembler.constraint_storage, backend
diff --git a/src/meshes/Exodus.jl b/src/meshes/Exodus.jl
index 63a842a..57bbeb6 100644
--- a/src/meshes/Exodus.jl
+++ b/src/meshes/Exodus.jl
@@ -202,7 +202,7 @@ function PostProcessor(
       append!(element_var_names, names(var))
     # L2QuadratureField case below
     elseif isa(var.fspace.coords, NamedTuple)
-      max_q_points = mapreduce(num_quadrature_points, maximum, values(var.fspace.ref_fes))
+      max_q_points = map(num_quadrature_points, values(var.fspace.ref_fes)) |> maximum
       temp_names = Symbol[]
       for name in names(var)
         for q in 1:max_q_points
diff --git a/src/solvers/DirectLinearSolver.jl b/src/solvers/DirectLinearSolver.jl
index 1fc8b3d..34a0935 100644
--- a/src/solvers/DirectLinearSolver.jl
+++ b/src/solvers/DirectLinearSolver.jl
@@ -32,3 +32,30 @@ function solve!(solver::DirectLinearSolver, Uu, p)
   map!((x, y) -> x + y, Uu, Uu, solver.ΔUu)
   return nothing
 end
+# struct DirectLinearSolver{
+#   U <: AbstractArray,
+#   A <: MatrixIntegral,
+#   b <: VectorIntegral,
+#   P
+# }
+#   Δu::U
+#   matrix::A
+#   vector::b
+#   precon::P
+#   timer::TimerOutput
+# end
+
+# function DirectLinearSolver(A::MatrixIntegral, b::VectorIntegral)
+#   precon = I
+#   ΔUu = create_unknowns(A)
+#   return DirectLinearSolver(Δu, A, b, precon, TimerOutput())
+# end
+
+# function solve!(solver::DirectLinearSolver, u, p)
+#   R = solver.matrix(u, p)
+#   K = solver.vector(u, p)
+#   # TODO need bcs as well
+#   solver.Δu .= -K \ R
+#   map!((x, y) -> x + y, u, u, solver.Δu)
+#   return nothing 
+# end
diff --git a/test/TestAssemblers.jl b/test/TestAssemblers.jl
index 98336e6..74dab50 100644
--- a/test/TestAssemblers.jl
+++ b/test/TestAssemblers.jl
@@ -26,10 +26,10 @@ function test_sparse_matrix_assembler()
   u = ScalarFunction(V, :u)
   @show asm = SparseMatrixAssembler(u)
   dbcs = DirichletBC[
-    DirichletBC(:u, :sset_1, bc_func),
-    DirichletBC(:u, :sset_2, bc_func),
-    DirichletBC(:u, :sset_3, bc_func),
-    DirichletBC(:u, :sset_4, bc_func),
+    DirichletBC(:u, bc_func; sideset_name = :sset_1),
+    DirichletBC(:u, bc_func; sideset_name = :sset_2),
+    DirichletBC(:u, bc_func; sideset_name = :sset_3),
+    DirichletBC(:u, bc_func; sideset_name = :sset_4),
   ]
   p = create_parameters(mesh, asm, physics, props; dirichlet_bcs=dbcs)
   FiniteElementContainers.initialize!(p)
diff --git a/test/TestBCs.jl b/test/TestBCs.jl
index 5bf7d0c..7177fcf 100644
--- a/test/TestBCs.jl
+++ b/test/TestBCs.jl
@@ -2,13 +2,49 @@ dummy_func_1(x, t) = 5. * t
 dummy_func_2(x, t) = 5. * SVector{2, Float64}(1., 0.)
 
 function test_dirichlet_bc_input()
-  bc = DirichletBC(:my_var, :my_sset, dummy_func_1)
+  # symbol constructors
+  bc = DirichletBC(:my_var, dummy_func_1; block_name = :my_block)
+  @test bc.block_name == :my_block
+  @test bc.nset_name === nothing
+  @test bc.sset_name === nothing
   @test bc.var_name == :my_var
+  @test typeof(bc.func) == typeof(dummy_func_1)
+
+  bc = DirichletBC(:my_var, dummy_func_1; nodeset_name = :my_nset)
+  @test bc.block_name === nothing
+  @test bc.nset_name == :my_nset
+  @test bc.sset_name === nothing
+  @test bc.var_name == :my_var
+  @test typeof(bc.func) == typeof(dummy_func_1)
+
+
+  bc = DirichletBC(:my_var, dummy_func_1; sideset_name = :my_sset)
+  @test bc.block_name === nothing
+  @test bc.nset_name === nothing
   @test bc.sset_name == :my_sset
+  @test bc.var_name == :my_var
   @test typeof(bc.func) == typeof(dummy_func_1)
-  bc = DirichletBC("my_var", "my_sset", dummy_func_1)
+
+  # string constructors
+  bc = DirichletBC("my_var", dummy_func_1; block_name = "my_block")
+  @test bc.block_name == :my_block
+  @test bc.nset_name === nothing
+  @test bc.sset_name === nothing
+  @test bc.var_name == :my_var
+  @test typeof(bc.func) == typeof(dummy_func_1)
+
+  bc = DirichletBC("my_var", dummy_func_1; nodeset_name = "my_nset")
+  @test bc.block_name === nothing
+  @test bc.nset_name == :my_nset
+  @test bc.sset_name === nothing
   @test bc.var_name == :my_var
+  @test typeof(bc.func) == typeof(dummy_func_1)
+
+  bc = DirichletBC("my_var", dummy_func_1; sideset_name = "my_sset")
+  @test bc.block_name === nothing
+  @test bc.nset_name === nothing
   @test bc.sset_name == :my_sset
+  @test bc.var_name == :my_var
   @test typeof(bc.func) == typeof(dummy_func_1)
 end
 
@@ -17,7 +53,7 @@ function test_dirichlet_bc_container_init()
   fspace = FunctionSpace(mesh, H1Field, Lagrange)
   u = VectorFunction(fspace, :displ)
   dof = DofManager(u)
-  bc_in = DirichletBC(:displ_x, :sset_1, dummy_func_1)
+  bc_in = DirichletBC(:displ_x, dummy_func_1; sideset_name = :sset_1)
   bc = FiniteElementContainers.DirichletBCContainer(mesh, dof, bc_in)
   @show bc
 end
diff --git a/test/TestFunctionSpaces.jl b/test/TestFunctionSpaces.jl
index ad0883c..24fd90c 100644
--- a/test/TestFunctionSpaces.jl
+++ b/test/TestFunctionSpaces.jl
@@ -88,6 +88,11 @@ end
 #   # # )
 # end
 
+
+function test_bad_interp_type(mesh)
+  @test_throws MethodError FunctionSpace(mesh, H1Field, :SomethingNotSupported)
+end
+
 function test_fspace_h1_field(mesh)
   @show fspace = FunctionSpace(mesh, H1Field, Lagrange)
 end
@@ -109,6 +114,7 @@ function test_function_spaces()
   # coords = H1Field(coords)
 
   mesh = UnstructuredMesh("mechanics/mechanics.g")
+  test_bad_interp_type(mesh)
   test_fspace_h1_field(mesh)
   test_fspace_l2_element_field(mesh)
   test_fspace_l2_quadrature_field(mesh)
diff --git a/test/TestIntegrals.jl b/test/TestIntegrals.jl
index c12ee6c..4207871 100644
--- a/test/TestIntegrals.jl
+++ b/test/TestIntegrals.jl
@@ -6,7 +6,7 @@ function test_scalar_integral()
     asm = SparseMatrixAssembler(u)
     f(X, _) = 2. * π^2 * sin(π * X[1]) * sin(π * X[2])
     physics = Poisson(f)
-    integral = FiniteElementContainers.ScalarIntegral(asm, energy)
+    integral = FiniteElementContainers.ScalarCellIntegral(asm, energy)
     grad_integral = FiniteElementContainers.gradient(integral)
 end
 
diff --git a/test/mechanics/TestMechanics.jl b/test/mechanics/TestMechanics.jl
index 2a8e231..5fc760f 100644
--- a/test/mechanics/TestMechanics.jl
+++ b/test/mechanics/TestMechanics.jl
@@ -24,10 +24,10 @@ function test_mechanics_dirichlet_only(
   asm = SparseMatrixAssembler(u; use_condensed=use_condensed)
 
   dbcs = DirichletBC[
-    DirichletBC(:displ_x, :sset_3, fixed),
-    DirichletBC(:displ_y, :sset_3, fixed),
-    DirichletBC(:displ_x, :sset_1, fixed),
-    DirichletBC(:displ_y, :sset_1, displace),
+    DirichletBC(:displ_x, fixed; sideset_name = :sset_3),
+    DirichletBC(:displ_y, fixed; sideset_name = :sset_3),
+    DirichletBC(:displ_x, fixed; sideset_name = :sset_1),
+    DirichletBC(:displ_y, displace; sideset_name = :sset_1),
   ]
 
   # pre-setup some scratch arrays
diff --git a/test/poisson/TestPoisson.jl b/test/poisson/TestPoisson.jl
index 54d8c29..da7f050 100644
--- a/test/poisson/TestPoisson.jl
+++ b/test/poisson/TestPoisson.jl
@@ -22,6 +22,7 @@ bc_func_neumann(_, _) = SVector{1, Float64}(1.)
 
 function test_poisson(backend, cond, nlsolver, lsolver)
   test_poisson_dirichlet(backend, cond, nlsolver, lsolver)
+  test_poisson_dirichlet_with_nodesets(backend, cond, nlsolver, lsolver)
   test_poisson_dirichlet_multi_block_quad4_quad4(backend, cond, nlsolver, lsolver)
   test_poisson_dirichlet_multi_block_quad4_tri3(backend, cond, nlsolver, lsolver)
   test_poisson_dirichlet_structured_mesh_quad4(backend, cond, nlsolver, lsolver)
@@ -44,10 +45,10 @@ function test_poisson_dirichlet(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :sset_1, bc_func),
-    DirichletBC(:u, :sset_2, bc_func),
-    DirichletBC(:u, :sset_3, bc_func),
-    DirichletBC(:u, :sset_4, bc_func),
+    DirichletBC(:u, bc_func; sideset_name = :sset_1),
+    DirichletBC(:u, bc_func; sideset_name = :sset_2),
+    DirichletBC(:u, bc_func; sideset_name = :sset_3),
+    DirichletBC(:u, bc_func; sideset_name = :sset_4),
   ]
 
   # setup the parameters
@@ -81,6 +82,55 @@ function test_poisson_dirichlet(
   display(solver.timer)
 end
 
+function test_poisson_dirichlet_with_nodesets(
+  dev, use_condensed,
+  nsolver, lsolver
+)
+  mesh = UnstructuredMesh(mesh_file)
+  V = FunctionSpace(mesh, H1Field, Lagrange) 
+  physics = Poisson(f)
+  props = create_properties(physics)
+  u = ScalarFunction(V, :u)
+  asm = SparseMatrixAssembler(u; use_condensed=use_condensed)
+
+  # setup and update bcs
+  dbcs = DirichletBC[
+    DirichletBC(:u, bc_func; nodeset_name = :nset_1),
+    DirichletBC(:u, bc_func; nodeset_name = :nset_2),
+    DirichletBC(:u, bc_func; nodeset_name = :nset_3),
+    DirichletBC(:u, bc_func; nodeset_name = :nset_4),
+  ]
+
+  # setup the parameters
+  p = create_parameters(mesh, asm, physics, props; dirichlet_bcs=dbcs)
+
+  if dev != cpu
+    p = p |> dev
+    asm = asm |> dev 
+  end
+
+  # setup solver and integrator
+  solver = nsolver(lsolver(asm))
+  integrator = QuasiStaticIntegrator(solver)
+  evolve!(integrator, p)
+
+  if dev != cpu
+    p = p |> cpu
+  end
+
+  U = p.h1_field
+
+  pp = PostProcessor(mesh, output_file, u)
+  write_times(pp, 1, 0.0)
+  write_field(pp, 1, ("u",), U)
+  close(pp)
+
+  if !Sys.iswindows()
+    @test exodiff(output_file, gold_file)
+  end
+  rm(output_file; force=true)
+  display(solver.timer)
+end
 
 function test_poisson_neumann(
   dev, use_condensed,
@@ -95,8 +145,8 @@ function test_poisson_neumann(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :sset_1, bc_func),
-    DirichletBC(:u, :sset_2, bc_func)
+    DirichletBC(:u, bc_func; sideset_name = :sset_1),
+    DirichletBC(:u, bc_func; sideset_name = :sset_2)
   ]
 
   nbcs = NeumannBC[
@@ -154,7 +204,7 @@ function test_poisson_dirichlet_multi_block_quad4_quad4(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :boundary, bc_func)
+    DirichletBC(:u, bc_func; sideset_name = :boundary)
   ]
 
   # setup the parameters
@@ -201,7 +251,7 @@ function test_poisson_dirichlet_multi_block_quad4_tri3(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :boundary, bc_func)
+    DirichletBC(:u, bc_func; sideset_name = :boundary)
   ]
 
   # setup the parameters
@@ -249,10 +299,10 @@ function test_poisson_dirichlet_structured_mesh_quad4(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :bottom, bc_func),
-    DirichletBC(:u, :right, bc_func),
-    DirichletBC(:u, :top, bc_func),
-    DirichletBC(:u, :left, bc_func),
+    DirichletBC(:u, bc_func; sideset_name = :bottom),
+    DirichletBC(:u, bc_func; sideset_name = :right),
+    DirichletBC(:u, bc_func; sideset_name = :top),
+    DirichletBC(:u, bc_func; sideset_name = :left),
   ]
 
   # setup the parameters
@@ -300,10 +350,10 @@ function test_poisson_dirichlet_structured_mesh_tri3(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :bottom, bc_func),
-    DirichletBC(:u, :right, bc_func),
-    DirichletBC(:u, :top, bc_func),
-    DirichletBC(:u, :left, bc_func),
+    DirichletBC(:u, bc_func; sideset_name = :bottom),
+    DirichletBC(:u, bc_func; sideset_name = :right),
+    DirichletBC(:u, bc_func; sideset_name = :top),
+    DirichletBC(:u, bc_func; sideset_name = :left),
   ]
 
   # setup the parameters
@@ -350,8 +400,8 @@ function test_poisson_neumann_structured_mesh_quad4(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :bottom, bc_func),
-    DirichletBC(:u, :right, bc_func)
+    DirichletBC(:u, bc_func; sideset_name = :bottom),
+    DirichletBC(:u, bc_func; sideset_name = :right)
   ]
 
   nbcs = NeumannBC[
@@ -409,8 +459,8 @@ function test_poisson_neumann_structured_mesh_tri3(
 
   # setup and update bcs
   dbcs = DirichletBC[
-    DirichletBC(:u, :bottom, bc_func),
-    DirichletBC(:u, :right, bc_func)
+    DirichletBC(:u, bc_func; sideset_name = :bottom),
+    DirichletBC(:u, bc_func; sideset_name = :right)
   ]
 
   nbcs = NeumannBC[