diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml
index 21a2c237..3ab301ac 100644
--- a/.github/workflows/build.yml
+++ b/.github/workflows/build.yml
@@ -23,7 +23,7 @@ jobs:
             c-compiler: 'clang'
             cxx-compiler: 'clang++'
       fail-fast: false
-    name: ${{ matrix.c-compiler}} Python ${{ matrix.py-version }}
+    name: Linux ${{ matrix.c-compiler}} Python ${{ matrix.py-version }}
 
 
     runs-on: ubuntu-latest
@@ -45,7 +45,7 @@ jobs:
       uses: actions/cache/restore@v3
       with:
         path: /home/runner/.cache/ccache
-        key: ccache-${{ matrix.c-compiler }}-${{ matrix.py-version }}
+        key: ccache-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}
 
     - uses: actions/setup-python@v4
       if: always()
@@ -77,7 +77,7 @@ jobs:
         GH_TOKEN: ${{ github.token }}
       run: |
         gh extension install actions/gh-actions-cache
-        gh actions-cache delete ccache-${{ matrix.c-compiler}}-${{ matrix.py-version }} --confirm
+        gh actions-cache delete ccache-${{ runner.os }}-${{ matrix.c-compiler}}-${{ matrix.py-version }} --confirm
       continue-on-error: true
 
     - name: Push ccache cache
@@ -85,7 +85,7 @@ jobs:
       uses: actions/cache/save@v3
       with:
         path: /home/runner/.cache/ccache
-        key: ccache-${{ matrix.c-compiler }}-${{ matrix.py-version }}
+        key: ccache-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}
 
     - name: Ccache report
       if: always()
@@ -93,27 +93,27 @@ jobs:
 
     - name: Run tests
       if: always()
-      run: python3 -m pytest --junit-xml=test-results-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml tests/
+      run: python3 -m pytest --junit-xml=test-results-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml tests/
 
     - name: Run flake8
       if: always()
       run: |
         make flake8
-        mv flake8.xml flake8-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
+        mv flake8.xml flake8-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
 
     - name: Archive test results
       uses: actions/upload-artifact@v3
       if: always()
       with:
         name: Test results
-        path: test-results-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
+        path: test-results-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
 
     - name: Report test results
       uses: dorny/test-reporter@v1
       if: always()
       with:
-        name: Test report (${{ matrix.c-compiler }}, Python ${{ matrix.py-version }})
-        path: test-results-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
+        name: Test report (${{ runner.os }}, ${{ matrix.c-compiler }}, Python ${{ matrix.py-version }})
+        path: test-results-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
         reporter: java-junit
         fail-on-error: true
 
@@ -121,7 +121,101 @@ jobs:
       uses: dorny/test-reporter@v1
       if: always()
       with:
-        name: Flake8 report (${{ matrix.c-compiler }}, Python ${{ matrix.py-version }})
-        path: flake8-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
+        name: Flake8 report (${{ runner.os }}, ${{ matrix.c-compiler }}, Python ${{ matrix.py-version }})
+        path: flake8-${{ runner.os }}-${{ matrix.c-compiler }}-${{ matrix.py-version }}.xml
         reporter: java-junit
         fail-on-error: false
+
+  mac:
+
+    strategy:
+      matrix:
+        py-version: ['3.10']
+      fail-fast: false
+    name: Mac clang Python ${{ matrix.py-version }}
+
+    runs-on: macos-latest
+    timeout-minutes: 180
+    env:
+      MPIEXEC_FLAGS: "--allow-run-as-root --oversubscribe"
+      PYNUCLEUS_BUILD_PARALLELISM: 3
+      SUITESPARSE_INCLUDE_DIR: /usr/local/Cellar/brewsci-suite-sparse/4.5.5_2/include
+      SUITESPARSE_LIBRARY_DIR: /usr/local/Cellar/brewsci-suite-sparse/4.5.5_2/lib
+      LDFLAGS: -L/usr/local/opt/brewsci-metis/lib -L/usr/local/opt/brewsci-parmetis/lib
+      CPPFLAGS: -I/usr/local/opt/brewsci-metis/include -I/usr/local/opt/brewsci-parmetis/include
+
+    steps:
+    - name: Check out repo
+      uses: actions/checkout@v3
+
+    - name: Pull ccache cache
+      if: always()
+      id: ccache-restore
+      uses: actions/cache/restore@v3
+      with:
+        path: /Users/runner/Library/Caches/ccache
+        key: ccache-${{ runner.os }}-${{ matrix.py-version }}
+
+    - name: Setup GNU Fortran
+      uses: modflowpy/install-gfortran-action@v1
+
+    - uses: actions/setup-python@v4
+      with:
+        python-version: ${{ matrix.py-version }}
+
+    - name: Install packages
+      run: |
+        brew install open-mpi ccache
+        brew tap brewsci/num
+        brew install brewsci-metis brewsci-parmetis brewsci-suite-sparse
+
+    - name: Get ccache config dir
+      if: always()
+      run: ccache --show-config
+
+    - name: Install Python dependencies
+      run: make prereq && make prereq-extra
+
+    - name: Install
+      run: PIP_INSTALL_FLAGS=-vvv make install
+
+    - name: Remove ccache cache
+      if: ${{ steps.ccache-restore.outputs.cache-hit }}
+      shell: bash
+      env:
+        GH_TOKEN: ${{ github.token }}
+      run: |
+        gh extension install actions/gh-actions-cache
+        gh actions-cache delete ccache-${{ runner.os }}-${{ matrix.py-version }} --confirm
+      continue-on-error: true
+
+    - name: Push ccache cache
+      if: always()
+      uses: actions/cache/save@v3
+      with:
+        path: /Users/runner/Library/Caches/ccache
+        key: ccache-${{ runner.os }}-${{ matrix.py-version }}
+
+    - name: Ccache report
+      if: always()
+      run: ccache -s
+
+    - name: Run tests
+      if: always()
+      run: python3 -m pytest --junit-xml=test-results-${{ runner.os }}-${{ matrix.py-version }}.xml tests/
+
+    - name: Archive test results
+      uses: actions/upload-artifact@v3
+      if: always()
+      with:
+        name: Test results
+        path: test-results-${{ runner.os }}-${{ matrix.py-version }}.xml
+
+    - name: Report test results
+      uses: dorny/test-reporter@v1
+      if: always()
+      with:
+        name: Test report (${{ runner.os }}, Python ${{ matrix.py-version }})
+        path: test-results-${{ runner.os }}-${{ matrix.py-version }}.xml
+        reporter: java-junit
+        fail-on-error: true
diff --git a/.github/workflows/mac.yml b/.github/workflows/mac.yml
deleted file mode 100644
index 8ce421d6..00000000
--- a/.github/workflows/mac.yml
+++ /dev/null
@@ -1,105 +0,0 @@
-
-name: CI MacOS
-
-on:
-  push:
-    branches: [ "master" ]
-  pull_request:
-    branches: [ "master" ]
-  workflow_dispatch:
-
-jobs:
-
-  mac:
-
-    strategy:
-      matrix:
-        py-version: ['3.10']
-      fail-fast: false
-    name: Python ${{ matrix.py-version }}
-
-    runs-on: macos-latest
-    timeout-minutes: 180
-    env:
-      MPIEXEC_FLAGS: "--allow-run-as-root --oversubscribe"
-      PYNUCLEUS_BUILD_PARALLELISM: 3
-      SUITESPARSE_INCLUDE_DIR: /usr/local/Cellar/brewsci-suite-sparse/4.5.5_2/include
-      SUITESPARSE_LIBRARY_DIR: /usr/local/Cellar/brewsci-suite-sparse/4.5.5_2/lib
-      LDFLAGS: -L/usr/local/opt/brewsci-metis/lib -L/usr/local/opt/brewsci-parmetis/lib
-      CPPFLAGS: -I/usr/local/opt/brewsci-metis/include -I/usr/local/opt/brewsci-parmetis/include
-
-    steps:
-    - name: Check out repo
-      uses: actions/checkout@v3
-
-    - name: Pull ccache cache
-      if: always()
-      id: ccache-restore
-      uses: actions/cache/restore@v3
-      with:
-        path: /Users/runner/Library/Caches/ccache
-        key: ccache-mac-${{ matrix.py-version }}
-
-    - name: Setup GNU Fortran
-      uses: modflowpy/install-gfortran-action@v1
-
-    - uses: actions/setup-python@v4
-      with:
-        python-version: ${{ matrix.py-version }}
-
-    - name: Install packages
-      run: |
-        brew install open-mpi ccache
-        brew tap brewsci/num
-        brew install brewsci-metis brewsci-parmetis brewsci-suite-sparse
-
-    - name: Get ccache config dir
-      if: always()
-      run: ccache --show-config
-
-    - name: Install Python dependencies
-      run: make prereq && make prereq-extra
-
-    - name: Install
-      run: PIP_INSTALL_FLAGS=-vvv make install
-
-    - name: Remove ccache cache
-      if: ${{ steps.ccache-restore.outputs.cache-hit }}
-      shell: bash
-      env:
-        GH_TOKEN: ${{ github.token }}
-      run: |
-        gh extension install actions/gh-actions-cache
-        gh actions-cache delete ccache-mac-${{ matrix.py-version }} --confirm
-      continue-on-error: true
-
-    - name: Push ccache cache
-      if: always()
-      uses: actions/cache/save@v3
-      with:
-        path: /Users/runner/Library/Caches/ccache
-        key: ccache-mac-${{ matrix.py-version }}
-
-    - name: Ccache report
-      if: always()
-      run: ccache -s
-
-    - name: Run tests
-      if: always()
-      run: python3 -m pytest --junit-xml=test-results-Mac-${{ matrix.py-version }}.xml tests/
-
-    - name: Archive test results
-      uses: actions/upload-artifact@v3
-      if: always()
-      with:
-        name: Test results
-        path: test-results-${{ matrix.py-version }}.xml
-
-    - name: Report test results
-      uses: dorny/test-reporter@v1
-      if: always()
-      with:
-        name: Test report (Mac, Python ${{ matrix.py-version }})
-        path: test-results-Mac-${{ matrix.py-version }}.xml
-        reporter: java-junit
-        fail-on-error: true
diff --git a/fem/PyNucleus_fem/mesh.py b/fem/PyNucleus_fem/mesh.py
index 0411c5d0..d1fa56ee 100644
--- a/fem/PyNucleus_fem/mesh.py
+++ b/fem/PyNucleus_fem/mesh.py
@@ -108,13 +108,20 @@ def pacman(h=0.1, **kwargs):
     return mesh
 
 
-def uniformSquare(N, M=None, ax=0, ay=0, bx=1, by=1, crossed=False, preserveLinesHorizontal=[], preserveLinesVertical=[]):
-    if M is None:
-        M = max(int(np.around((by-ay)/(bx-ax)))*N, 2)
-    assert N >= 2
-    assert M >= 2
-    xVals = np.linspace(ax, bx, N)
-    yVals = np.linspace(ay, by, M)
+def uniformSquare(N=2, M=None, ax=0, ay=0, bx=1, by=1, crossed=False, preserveLinesHorizontal=[], preserveLinesVertical=[], xVals=None, yVals=None):
+    if xVals is None:
+        assert N is not None
+        assert N >= 2
+        xVals = np.linspace(ax, bx, N)
+    else:
+        N = xVals.shape[0]
+    if yVals is None:
+        if M is None:
+            M = max(int(np.around((by-ay)/(bx-ax)))*N, 2)
+        assert M >= 2
+        yVals = np.linspace(ay, by, M)
+    else:
+        M = yVals.shape[0]
     x, y = np.meshgrid(xVals, yVals)
     for yVal in preserveLinesHorizontal:
         assert (yVals-yVal).min() < 1e-10
@@ -310,70 +317,54 @@ def squareWithInteractions(ax, ay, bx, by,
     if not uniform:
         from . meshConstruction import (circularSegment,
                                         line,
+                                        polygon,
                                         transformationRestriction)
         if h is None:
             h = horizon
         bottomLeft = np.array([ax, ay])
+        bottomRight = np.array([bx, ay])
+        topRight = np.array([bx, by])
+        topLeft = np.array([ax, by])
+
+        horizontalOffset = np.array([horizon, 0.])
+        verticalOffset = np.array([0., horizon])
+
         center = np.array([(ax+bx)/2, (ay+by)/2])
 
         numPointsPerUnitLength = int(np.ceil(1/h))
 
         assert len(preserveLinesVertical) == 0 or len(preserveLinesHorizontal) == 0
-        if len(preserveLinesVertical)+len(preserveLinesHorizontal) > 0:
-            preserve = preserveLinesVertical+preserveLinesHorizontal
-
-            c1 = circularSegment(bottomLeft, horizon, np.pi, 3/2*np.pi, numPointsPerUnitLength)
-
-            x1 = preserve[0]
-            c2 = line((ax, ay), (x1, ay))
-            for k in range(len(preserve)-1):
-                x1 = preserve[k]
-                x2 = preserve[k+1]
-                c2 = c2+line((x1, ay), (x2, ay))
-            x2 = preserve[-1]
-            c2 = c2+line((x2, ay), (bx, ay))
-            c1 = c1 + (c2+(0., -horizon))
-        else:
-            c1 = circularSegment(bottomLeft, horizon, np.pi, 3/2*np.pi, numPointsPerUnitLength)
-            c2 = line((ax, ay), (bx, ay))
-            c1 = c1 + (c2+(0., -horizon))
-        c3 = line((ax, ay), (ax, ay-horizon))
-        c4 = line((ax, ay), (ax-horizon, ay))
-        c = c1+c2+c3+c4
-
-        frame = (c + (c*(center, np.pi/2)) + (c*(center, np.pi)) + (c*(center, -np.pi/2)))
-
-        if len(preserveLinesVertical) > 0:
-            d = line((0, ay-horizon), (0, ay))
-            x1 = preserve[0]
-            d = d + line((0, ay), (0, x1))
-            for k in range(len(preserve)-1):
-                x1 = preserve[k]
-                x2 = preserve[k+1]
-                d += line((0, x1), (0, x2))
-            x2 = preserve[-1]
-            d = d + line((0, x2), (0, by))
-            d = d + line((0, by), (0, by+horizon))
-            for x in preserveLinesVertical:
-                assert ax <= x <= bx
-                frame += (d+(x, 0.))
-        if len(preserveLinesHorizontal) > 0:
-            d = line((ax-horizon, 0), (ax, 0))+line((ax, 0), (bx, 0))+line((bx, 0), (bx+horizon, 0))
-
-            d = line((ax-horizon, 0), (ax, 0))
-            x1 = preserve[0]
-            d = d + line((ax, 0), (x1, 0))
-            for k in range(len(preserve)-1):
-                x1 = preserve[k]
-                x2 = preserve[k+1]
-                d += line((x1, 0), (x2, 0))
-            x2 = preserve[-1]
-            d = d + line((x2, 0), (bx, 0))
-            d = d + line((bx, 0), (bx+horizon, 0))
-
-            for y in preserveLinesHorizontal:
-                assert ay <= y <= by
-                frame += (d+(0, y))
+
+        lineHorizontal = polygon([(0., 0.)] + [(p-ax, 0.) for p in preserveLinesVertical] + [(bx-ax, 0.)], doClose=False)
+        lineVertical = polygon([(0., 0.)] + [(0., p-ay) for p in preserveLinesHorizontal] + [(0., by-ay)], doClose=False)
+
+        d1 = (circularSegment(bottomLeft, horizon, np.pi, 1.5*np.pi, numPointsPerUnitLength) +
+              line(bottomLeft, bottomLeft-horizontalOffset) +
+              line(bottomLeft, bottomLeft-verticalOffset) +
+              (lineHorizontal+bottomLeft)+
+              (lineHorizontal+(bottomLeft-verticalOffset)))
+
+        d2 = (circularSegment(bottomRight, horizon, -0.5*np.pi, 0., numPointsPerUnitLength) +
+              line(bottomRight, bottomRight+horizontalOffset) +
+              line(bottomRight, bottomRight-verticalOffset) +
+              (lineVertical+(bottomRight+horizontalOffset)) +
+              (lineVertical+bottomRight))
+
+        d3 = (circularSegment(topRight, horizon, 0, 0.5*np.pi, numPointsPerUnitLength) +
+              line(topRight, topRight+horizontalOffset) +
+              line(topRight, topRight+verticalOffset) +
+              (lineHorizontal+topLeft) +
+              (lineHorizontal+(topLeft+verticalOffset)))
+
+        d4 = (circularSegment(topLeft, horizon, 0.5*np.pi, np.pi, numPointsPerUnitLength) +
+              line(topLeft, topLeft-horizontalOffset) +
+              line(topLeft, topLeft+verticalOffset) +
+              (lineVertical+bottomLeft) +
+              (lineVertical+(bottomLeft-horizontalOffset)))
+
+        frame = d1 + d2 + d3 + d4
+
+        frame.holes.append(center)
 
         if innerRadius > 0:
             frame += transformationRestriction(circularSegment(center, innerRadius, 0, 2*np.pi, numPointsPerUnitLength),
@@ -381,25 +372,30 @@ def squareWithInteractions(ax, ay, bx, by,
                                                center+(innerRadius, innerRadius))
             mesh = frame.mesh(max_volume=h**2, min_angle=30, **kwargs)
         else:
-            frame.holes.append(center)
             mesh = frame.mesh(max_volume=h**2, min_angle=30, **kwargs)
 
             eps = 1e-10
-            N1 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 0]-ax) < eps,
-                                np.logical_and(mesh.vertices_as_array[:, 1] >= ay-eps,
-                                               mesh.vertices_as_array[:, 1] <= by+eps)).sum()
-            N2 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 0]-bx) < eps,
-                                np.logical_and(mesh.vertices_as_array[:, 1] >= ay-eps,
-                                               mesh.vertices_as_array[:, 1] <= by+eps)).sum()
-            M1 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 1]-ay) < eps,
+            idx1 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 0]-ax) < eps,
+                                  np.logical_and(mesh.vertices_as_array[:, 1] >= ay-eps,
+                                                 mesh.vertices_as_array[:, 1] <= by+eps))
+            idx2 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 0]-bx) < eps,
+                                  np.logical_and(mesh.vertices_as_array[:, 1] >= ay-eps,
+                                                 mesh.vertices_as_array[:, 1] <= by+eps))
+            yVals1 = np.sort(mesh.vertices_as_array[idx1, 1])
+            yVals2 = np.sort(mesh.vertices_as_array[idx2, 1])
+            assert np.allclose(yVals1, yVals2), (yVals1, yVals2)
+
+            idx3 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 1]-ay) < eps,
                                 np.logical_and(mesh.vertices_as_array[:, 0] >= ax-eps,
-                                               mesh.vertices_as_array[:, 0] <= bx+eps)).sum()
-            M2 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 1]-by) < eps,
+                                               mesh.vertices_as_array[:, 0] <= bx+eps))
+            idx4 = np.logical_and(np.absolute(mesh.vertices_as_array[:, 1]-by) < eps,
                                 np.logical_and(mesh.vertices_as_array[:, 0] >= ax-eps,
-                                               mesh.vertices_as_array[:, 0] <= bx+eps)).sum()
-            assert N1 == N2, (N1, N2)
-            assert M1 == M2, (M1, M2)
-            mesh2 = uniformSquare(N=N1, M=M1, ax=ax, ay=ay, bx=bx, by=by)
+                                               mesh.vertices_as_array[:, 0] <= bx+eps))
+            xVals3 = np.sort(mesh.vertices_as_array[idx3, 0])
+            xVals4 = np.sort(mesh.vertices_as_array[idx4, 0])
+            assert np.allclose(xVals3, xVals4), (xVals3, xVals4)
+            mesh2 = uniformSquare(ax=ax, ay=ay, bx=bx, by=by, xVals=xVals3, yVals=yVals1)
+            mesh2.plot()
             mesh = snapMeshes(mesh, mesh2)
 
         location = uninitialized((mesh.num_vertices), dtype=INDEX)
diff --git a/fem/PyNucleus_fem/meshConstruction.py b/fem/PyNucleus_fem/meshConstruction.py
index 76b9dada..5869b7ca 100644
--- a/fem/PyNucleus_fem/meshConstruction.py
+++ b/fem/PyNucleus_fem/meshConstruction.py
@@ -23,7 +23,7 @@ def __init__(self, points, facets, holes=[]):
         self.meshTransformations = []
 
     def __add__(self, other):
-        if isinstance(other, tuple):
+        if isinstance(other, (tuple, np.ndarray)):
             newPoints = [(other[0]+p[0], other[1]+p[1]) for p in self.points]
             newHoles = [(other[0]+p[0], other[1]+p[1]) for p in self.holes]
             newSegment = segment(newPoints, self.facets, newHoles)