diff --git a/.github/workflows/style.yml b/.github/workflows/style.yml
index 2694b4cc7..8da70be16 100644
--- a/.github/workflows/style.yml
+++ b/.github/workflows/style.yml
@@ -20,13 +20,13 @@ jobs:
       - name: Install clang-format
         run: |
           wget -O- https://apt.llvm.org/llvm-snapshot.gpg.key | sudo apt-key add -
-          sudo add-apt-repository 'deb http://apt.llvm.org/jammy/ llvm-toolchain-jammy-18 main'
+          sudo add-apt-repository 'deb http://apt.llvm.org/jammy/ llvm-toolchain-jammy-19 main'
           sudo apt-get update -q
-          sudo apt-get install -y clang-format-18
+          sudo apt-get install -y clang-format-19
 
       - name: Check style
         run: |
-          ./scripts/format-source --clang-format clang-format-18
+          ./scripts/format-source --clang-format clang-format-19
           if [[ `git status -s | wc -l` -ne 0 ]]; then
             echo 'Error: Commit is not formatted!'
             echo 'Run '\`'./scripts/format-source'\`' in the source root directory'
diff --git a/examples/cpw/cpw.jl b/examples/cpw/cpw.jl
index e2ae1f56d..dd5db420f 100644
--- a/examples/cpw/cpw.jl
+++ b/examples/cpw/cpw.jl
@@ -74,7 +74,7 @@ function generate_cpw_data(; num_processors::Integer=1)
     data_wave_uniform = data_wave_uniform[:, 2:end]
 
     # Wrap phases
-    for p = 1:(size(data_lumped_adaptive, 2) ÷ 2)
+    for p = 1:(size(data_lumped_adaptive, 2)÷2)
         idx = (data_lumped_adaptive[:, 2 * p] .< 0.0)
         data_lumped_adaptive[idx, 2 * p] = data_lumped_adaptive[idx, 2 * p] .+ 180.0
 
diff --git a/examples/cpw/mesh/mesh.jl b/examples/cpw/mesh/mesh.jl
index 449418892..6c9766a7f 100644
--- a/examples/cpw/mesh/mesh.jl
+++ b/examples/cpw/mesh/mesh.jl
@@ -42,19 +42,19 @@ Generate a mesh for the coplanar waveguide with wave ports using Gmsh
 """
 function generate_cpw_wave_mesh(;
     filename::AbstractString,
-    refinement::Integer       = 0,
-    order::Integer            = 1,
+    refinement::Integer        = 0,
+    order::Integer             = 1,
     trace_width_μm::Real      = 30.0,
     gap_width_μm::Real        = 18.0,
     separation_width_μm::Real = 200.0,
     ground_width_μm::Real     = 800.0,
     substrate_height_μm::Real = 500.0,
     metal_height_μm::Real     = 0.0,
-    remove_metal_vol::Bool    = true,
+    remove_metal_vol::Bool     = true,
     length_μm::Real           = 4000.0,
-    coax_ports::Bool          = false,
-    verbose::Integer          = 5,
-    gui::Bool                 = false
+    coax_ports::Bool           = false,
+    verbose::Integer           = 5,
+    gui::Bool                  = false
 )
     @assert refinement >= 0
     @assert order > 0
@@ -224,14 +224,13 @@ function generate_cpw_wave_mesh(;
     kernel.synchronize()
 
     # Add physical groups
-    metal_domains =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(x -> x[1] == 3 && x[2] in metal, geom_dimtags)]
-                )
+    metal_domains = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 3 && x[2] in metal, geom_dimtags)]
             )
         )
+    )
 
     si_domain = last.(geom_map[findfirst(x -> x == (3, substrate), geom_dimtags)])
     @assert length(si_domain) == 1
@@ -244,8 +243,7 @@ function generate_cpw_wave_mesh(;
 
     if length(metal_domains) > 0 && remove_metal_vol
         remove_dimtags = [(3, x) for x in metal_domains]
-        for tag in
-            last.(
+        for tag in last.(
             filter(
                 x -> x[1] == 2,
                 gmsh.model.getBoundary(
@@ -289,17 +287,13 @@ function generate_cpw_wave_mesh(;
     end1_group = gmsh.model.addPhysicalGroup(2, end1, -1, "end1")
     end2_group = gmsh.model.addPhysicalGroup(2, end2, -1, "end2")
 
-    farfield =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in domain_boundary,
-                        geom_dimtags
-                    )]
-                )
+    farfield = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in domain_boundary, geom_dimtags)]
             )
         )
+    )
     filter!(
         x -> !(
             x in port1 || x in port2 || x in port3 || x in port4 || x in end1 || x in end2
@@ -309,43 +303,34 @@ function generate_cpw_wave_mesh(;
 
     farfield_group = gmsh.model.addPhysicalGroup(2, farfield, -1, "farfield")
 
-    trace =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in metal_boundary,
-                        geom_dimtags
-                    )]
-                )
+    trace = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in metal_boundary, geom_dimtags)]
             )
         )
-    gap =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in [n1, n2, n3, n4],
-                        geom_dimtags
-                    )]
-                )
+    )
+    gap = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in [n1, n2, n3, n4], geom_dimtags)]
             )
         )
+    )
 
     trace_group = gmsh.model.addPhysicalGroup(2, trace, -1, "trace")
     gap_group = gmsh.model.addPhysicalGroup(2, gap, -1, "gap")
 
-    trace_top =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in metal_boundary_top,
-                        geom_dimtags
-                    )]
-                )
+    trace_top = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(
+                    x -> x[1] == 2 && x[2] in metal_boundary_top,
+                    geom_dimtags
+                )]
             )
         )
+    )
     filter!(
         x -> !(
             x in port1 || x in port2 || x in port3 || x in port4 || x in end1 || x in end2
@@ -362,20 +347,18 @@ function generate_cpw_wave_mesh(;
     gmsh.option.setNumber("Mesh.MeshSizeFromCurvature", 0)
     gmsh.option.setNumber("Mesh.MeshSizeExtendFromBoundary", 0)
 
-    gap_points =
-        last.(
-            filter(
-                x -> x[1] == 0,
-                gmsh.model.getBoundary([(2, z) for z in gap], false, true, true)
-            )
+    gap_points = last.(
+        filter(
+            x -> x[1] == 0,
+            gmsh.model.getBoundary([(2, z) for z in gap], false, true, true)
         )
-    gap_curves =
-        last.(
-            filter(
-                x -> x[1] == 1,
-                gmsh.model.getBoundary([(2, z) for z in gap], false, false, false)
-            )
+    )
+    gap_curves = last.(
+        filter(
+            x -> x[1] == 1,
+            gmsh.model.getBoundary([(2, z) for z in gap], false, false, false)
         )
+    )
 
     gmsh.model.mesh.field.add("Distance", 1)
     gmsh.model.mesh.field.setNumbers(1, "PointsList", gap_points)
@@ -478,18 +461,18 @@ Generate a mesh for the coplanar waveguide with lumped ports using Gmsh
 """
 function generate_cpw_lumped_mesh(;
     filename::AbstractString,
-    refinement::Integer       = 0,
-    order::Integer            = 1,
+    refinement::Integer        = 0,
+    order::Integer             = 1,
     trace_width_μm::Real      = 30.0,
     gap_width_μm::Real        = 18.0,
     separation_width_μm::Real = 200.0,
     ground_width_μm::Real     = 800.0,
     substrate_height_μm::Real = 500.0,
     metal_height_μm::Real     = 0.0,
-    remove_metal_vol::Bool    = true,
+    remove_metal_vol::Bool     = true,
     length_μm::Real           = 4000.0,
-    verbose::Integer          = 5,
-    gui::Bool                 = false
+    verbose::Integer           = 5,
+    gui::Bool                  = false
 )
     @assert refinement >= 0
     @assert order > 0
@@ -594,14 +577,13 @@ function generate_cpw_lumped_mesh(;
     kernel.synchronize()
 
     # Add physical groups
-    metal_domains =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(x -> x[1] == 3 && x[2] in metal, geom_dimtags)]
-                )
+    metal_domains = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 3 && x[2] in metal, geom_dimtags)]
             )
         )
+    )
 
     si_domain = last.(geom_map[findfirst(x -> x == (3, substrate), geom_dimtags)])
     @assert length(si_domain) == 1
@@ -614,8 +596,7 @@ function generate_cpw_lumped_mesh(;
 
     if length(metal_domains) > 0 && remove_metal_vol
         remove_dimtags = [(3, x) for x in metal_domains]
-        for tag in
-            last.(
+        for tag in last.(
             filter(
                 x -> x[1] == 2,
                 gmsh.model.getBoundary(
@@ -641,17 +622,13 @@ function generate_cpw_lumped_mesh(;
     si_domain_group = gmsh.model.addPhysicalGroup(3, [si_domain], -1, "si")
     metal_domain_group = gmsh.model.addPhysicalGroup(3, metal_domains, -1, "metal")
 
-    farfield =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in domain_boundary,
-                        geom_dimtags
-                    )]
-                )
+    farfield = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in domain_boundary, geom_dimtags)]
             )
         )
+    )
 
     farfield_group = gmsh.model.addPhysicalGroup(2, farfield, -1, "farfield")
 
@@ -673,28 +650,20 @@ function generate_cpw_lumped_mesh(;
     port3b_group = gmsh.model.addPhysicalGroup(2, port3b, -1, "port3b")
     port4b_group = gmsh.model.addPhysicalGroup(2, port4b, -1, "port4b")
 
-    trace =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in metal_boundary,
-                        geom_dimtags
-                    )]
-                )
+    trace = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in metal_boundary, geom_dimtags)]
             )
         )
-    gap =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in [n1, n2, n3, n4],
-                        geom_dimtags
-                    )]
-                )
+    )
+    gap = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(x -> x[1] == 2 && x[2] in [n1, n2, n3, n4], geom_dimtags)]
             )
         )
+    )
     filter!(
         x -> !(
             x in port1a ||
@@ -712,17 +681,16 @@ function generate_cpw_lumped_mesh(;
     trace_group = gmsh.model.addPhysicalGroup(2, trace, -1, "trace")
     gap_group = gmsh.model.addPhysicalGroup(2, gap, -1, "gap")
 
-    trace_top =
-        last.(
-            collect(
-                Iterators.flatten(
-                    geom_map[findall(
-                        x -> x[1] == 2 && x[2] in metal_boundary_top,
-                        geom_dimtags
-                    )]
-                )
+    trace_top = last.(
+        collect(
+            Iterators.flatten(
+                geom_map[findall(
+                    x -> x[1] == 2 && x[2] in metal_boundary_top,
+                    geom_dimtags
+                )]
             )
         )
+    )
 
     trace_top_group = gmsh.model.addPhysicalGroup(2, trace_top, -1, "trace2")
 
@@ -733,20 +701,18 @@ function generate_cpw_lumped_mesh(;
     gmsh.option.setNumber("Mesh.MeshSizeFromCurvature", 0)
     gmsh.option.setNumber("Mesh.MeshSizeExtendFromBoundary", 0)
 
-    gap_points =
-        last.(
-            filter(
-                x -> x[1] == 0,
-                gmsh.model.getBoundary([(2, z) for z in gap], false, true, true)
-            )
+    gap_points = last.(
+        filter(
+            x -> x[1] == 0,
+            gmsh.model.getBoundary([(2, z) for z in gap], false, true, true)
         )
-    gap_curves =
-        last.(
-            filter(
-                x -> x[1] == 1,
-                gmsh.model.getBoundary([(2, z) for z in gap], false, false, false)
-            )
+    )
+    gap_curves = last.(
+        filter(
+            x -> x[1] == 1,
+            gmsh.model.getBoundary([(2, z) for z in gap], false, false, false)
         )
+    )
 
     gmsh.model.mesh.field.add("Distance", 1)
     gmsh.model.mesh.field.setNumbers(1, "PointsList", gap_points)
diff --git a/examples/cylinder/cavity.jl b/examples/cylinder/cavity.jl
index aa7614944..15482acb7 100644
--- a/examples/cylinder/cavity.jl
+++ b/examples/cylinder/cavity.jl
@@ -236,19 +236,19 @@ function generate_cavity_convergence_data(;
         0,
         1,
         0;
-        ϵᵣ   = 2.08,
-        μᵣ   = 1.0,
-        a_cm = radius,
-        d_cm = aspect_ratio * 2 * radius
+        ϵᵣ = 2.08,
+        μᵣ = 1.0,
+        a_cm  = radius,
+        d_cm  = aspect_ratio * 2 * radius
     )
     f_TE_111_true, ~ = frequency_transverse(
         1,
         1,
         1;
-        ϵᵣ   = 2.08,
-        μᵣ   = 1.0,
-        a_cm = radius,
-        d_cm = aspect_ratio * 2 * radius
+        ϵᵣ = 2.08,
+        μᵣ = 1.0,
+        a_cm  = radius,
+        d_cm  = aspect_ratio * 2 * radius
     )
 
     dof = Vector{Vector{Int}}()
diff --git a/examples/cylinder/mesh/mesh.jl b/examples/cylinder/mesh/mesh.jl
index 59d986e52..5ab6b27cb 100644
--- a/examples/cylinder/mesh/mesh.jl
+++ b/examples/cylinder/mesh/mesh.jl
@@ -106,42 +106,39 @@ function generate_cylindrical_mesh(;
     # Add physical groups
     cylinder_group = gmsh.model.addPhysicalGroup(3, cylinder, -1, "cylinder")
 
-    @show top =
-        last.(
-            gmsh.model.getEntitiesInBoundingBox(
-                -1.1 * radius,
-                -1.1 * radius,
-                0.5 * height,
-                1.1 * radius,
-                1.1 * radius,
-                1.1 * height,
-                2
-            )
+    @show top = last.(
+        gmsh.model.getEntitiesInBoundingBox(
+            -1.1 * radius,
+            -1.1 * radius,
+            0.5 * height,
+            1.1 * radius,
+            1.1 * radius,
+            1.1 * height,
+            2
         )
-    @show bottom =
-        last.(
-            gmsh.model.getEntitiesInBoundingBox(
-                -1.1 * radius,
-                -1.1 * radius,
-                -0.1 * height,
-                1.1 * radius,
-                1.1 * radius,
-                0.5 * height,
-                2
-            )
+    )
+    @show bottom = last.(
+        gmsh.model.getEntitiesInBoundingBox(
+            -1.1 * radius,
+            -1.1 * radius,
+            -0.1 * height,
+            1.1 * radius,
+            1.1 * radius,
+            0.5 * height,
+            2
         )
-    @show symmetry =
-        last.(
-            gmsh.model.getEntitiesInBoundingBox(
-                -1.1 * radius,
-                -0.1 * radius,
-                -0.1 * height,
-                1.1 * radius,
-                0.1 * radius,
-                1.1 * height,
-                2
-            )
+    )
+    @show symmetry = last.(
+        gmsh.model.getEntitiesInBoundingBox(
+            -1.1 * radius,
+            -0.1 * radius,
+            -0.1 * height,
+            1.1 * radius,
+            0.1 * radius,
+            1.1 * height,
+            2
         )
+    )
     exterior = []
     for domain in cylinder
         _, domain_boundaries = gmsh.model.getAdjacencies(3, domain)
@@ -174,25 +171,24 @@ function generate_cylindrical_mesh(;
     gmsh.option.setNumber("Mesh.Algorithm3D", 10)
 
     if (mesh_type == 2)
-        base_boundaries =
-            last.(
-                gmsh.model.getEntitiesInBoundingBox(
-                    -1.1 * radius,
-                    -1.1 * radius,
-                    -0.1 * height,
-                    1.1 * radius,
-                    1.1 * radius,
-                    0.1 * height,
-                    2
-                )
+        base_boundaries = last.(
+            gmsh.model.getEntitiesInBoundingBox(
+                -1.1 * radius,
+                -1.1 * radius,
+                -0.1 * height,
+                1.1 * radius,
+                1.1 * radius,
+                0.1 * height,
+                2
             )
+        )
         for boundary in base_boundaries
             gmsh.model.mesh.setRecombine(2, boundary)
         end
     end
 
     gmsh.model.mesh.generate(3)
-    for i = 0:(refinement - 1)
+    for i = 0:(refinement-1)
         gmsh.model.mesh.refine()
     end
     gmsh.model.mesh.setOrder(order)
diff --git a/examples/rings/mesh/mesh.jl b/examples/rings/mesh/mesh.jl
index 71c16337f..4693cea11 100644
--- a/examples/rings/mesh/mesh.jl
+++ b/examples/rings/mesh/mesh.jl
@@ -38,7 +38,7 @@ function generate_ring_mesh(;
     outer_radius                       = 100.0,
     rot_center::AbstractVector{<:Real} = [0.0, 0.0, 0.0],
     rot_axis::AbstractVector{<:Real}   = [0.0, 0.0, 1.0],
-    rot_θ::Real                        = π / 2,
+    rot_θ::Real                       = π / 2,
     verbose::Integer                   = 5,
     gui::Bool                          = false
 )
@@ -185,15 +185,14 @@ function generate_ring_mesh(;
     gmsh.model.mesh.field.setNumber(1, "DistMax", 6.0 * outer_radius)
     gmsh.model.mesh.field.setNumber(1, "SizeMax", l_farfield)
 
-    mesh_curves =
-        last.(
-            gmsh.model.getBoundary(
-                [(2, x) for x in [inner_ring, outer_ring, inner_terminal, outer_terminal]],
-                true,
-                false,
-                false
-            )
+    mesh_curves = last.(
+        gmsh.model.getBoundary(
+            [(2, x) for x in [inner_ring, outer_ring, inner_terminal, outer_terminal]],
+            true,
+            false,
+            false
         )
+    )
 
     gmsh.model.mesh.field.add("Distance", 2)
     gmsh.model.mesh.field.setNumbers(2, "CurvesList", mesh_curves)
diff --git a/palace/fem/errorindicator.cpp b/palace/fem/errorindicator.cpp
index 507ab696a..a1316b661 100644
--- a/palace/fem/errorindicator.cpp
+++ b/palace/fem/errorindicator.cpp
@@ -38,10 +38,9 @@ void ErrorIndicator::AddIndicator(const Vector &indicator)
   const int Dn = n;
   const auto *DI = indicator.Read();
   auto *DL = local.ReadWrite();
-  mfem::forall_switch(use_dev, N,
-                      [=] MFEM_HOST_DEVICE(int i) {
-                        DL[i] = std::sqrt((DL[i] * DL[i] * Dn + DI[i] * DI[i]) / (Dn + 1));
-                      });
+  mfem::forall_switch(
+      use_dev, N, [=] MFEM_HOST_DEVICE(int i)
+      { DL[i] = std::sqrt((DL[i] * DL[i] * Dn + DI[i] * DI[i]) / (Dn + 1)); });
 
   // More samples have been added, update for the running average.
   n += 1;
diff --git a/palace/utils/geodata.cpp b/palace/utils/geodata.cpp
index b3073de45..65f535e83 100644
--- a/palace/utils/geodata.cpp
+++ b/palace/utils/geodata.cpp
@@ -1362,7 +1362,8 @@ BoundingBox BoundingBallFromPointCloud(MPI_Comm comm,
       // Find the next point as the vertex furthest from the initial axis.
       const Eigen::Vector3d n_1 = (*p_2 - *p_1).normalized();
       auto p_3 = std::max_element(vertices.begin(), vertices.end(),
-                                  [&](const auto &x, const auto &y) {
+                                  [&](const auto &x, const auto &y)
+                                  {
                                     return PerpendicularDistance({n_1}, *p_1, x) <
                                            PerpendicularDistance({n_1}, *p_1, y);
                                   });