diff --git a/src/main/scala/org/appliedtopology/tda4j/PrintingHelper.scala b/src/main/scala/org/appliedtopology/tda4j/PrintingHelper.scala
new file mode 100644
index 00000000..1c0203e1
--- /dev/null
+++ b/src/main/scala/org/appliedtopology/tda4j/PrintingHelper.scala
@@ -0,0 +1,67 @@
+package org.appliedtopology.tda4j
+package unicode
+
+def applyTranslation(trMap: Map[Char,String]): (String => String) = s =>
+  s.flatMap(trMap.orElse(_.toString))
+
+def applyCharTranslation(trMap: Map[Char,Char]): (String => String) = s =>
+  s.map(trMap.orElse{c => c})
+
+val superScriptMap = Map(
+  '0' -> '⁰',
+  '1' -> '¹',
+  '2' -> '²',
+  '3' -> '³',
+  '4' -> '⁴',
+  '5' -> '⁵',
+  '6' -> '⁶',
+  '7' -> '⁷',
+  '8' -> '⁸',
+  '9' -> '⁹',
+  '+' -> '⁺',
+  '-' -> '⁻',
+  '=' -> '⁼',
+  '(' -> '⁽',
+  ')' -> '⁾',
+  'i' -> 'ⁱ',
+  'n' -> 'ⁿ',
+  'h' -> 'ʰ',
+)
+
+val subScriptMap = Map(
+  'i' -> 'ᵢ',
+  'r' -> 'ᵣ',
+  'u' -> 'ᵤ',
+  'v' -> 'ᵥ',
+  '0' -> '₀',
+  '1' -> '₁',
+  '2' -> '₂',
+  '3' -> '₃',
+  '4' -> '₄',
+  '5' -> '₅',
+  '6' -> '₆',
+  '7' -> '₇',
+  '8' -> '₈',
+  '9' -> '₉',
+  '+' -> '₊',
+  '-' -> '₋',
+  '=' -> '₌',
+  '(' -> '₍',
+  ')' -> '₎',
+  'a' -> 'ₐ',
+  'e' -> 'ₔ',
+  'o' -> 'ₒ',
+  'x' -> 'ₓ',
+  'h' -> 'ₕ',
+  'k' -> 'ₖ',
+  'l' -> 'ₗ',
+  'm' -> 'ₘ',
+  'n' -> 'ₙ',
+  'p' -> 'ₚ',
+  's' -> 'ₛ',
+  't' -> 'ₜ',
+  'j' -> 'ⱼ',
+)
+
+def unicodeSuperScript = applyCharTranslation(superScriptMap)
+def unicodeSubScript = applyCharTranslation(subScriptMap)
\ No newline at end of file
diff --git a/src/main/scala/org/appliedtopology/tda4j/SimplicialSet.scala b/src/main/scala/org/appliedtopology/tda4j/SimplicialSet.scala
index d42be613..8374c473 100644
--- a/src/main/scala/org/appliedtopology/tda4j/SimplicialSet.scala
+++ b/src/main/scala/org/appliedtopology/tda4j/SimplicialSet.scala
@@ -1,11 +1,13 @@
 package org.appliedtopology.tda4j
 
+import org.appliedtopology.tda4j.unicode
 import scala.annotation.tailrec
 
 trait SimplicialSetElement {
   def base: SimplicialSetElement
   def dim: Int
   def degeneracies: List[Int]
+  override def toString: String = unicode.unicodeSuperScript(s"∆$dim")
 }
 
 def simplicialGenerator(generatorDim: Int): SimplicialSetElement =
@@ -14,8 +16,29 @@ def simplicialGenerator(generatorDim: Int): SimplicialSetElement =
     override def dim: Int = generatorDim
     override def degeneracies: List[Int] = List.empty
 
+case class SimplicialWrapper[T:OrderedCell](wrapped : T) extends SimplicialSetElement {
+  override def base: SimplicialSetElement = this
+  override def dim: Int = wrapped.dim
+  override def degeneracies: List[Int] = List.empty
+}
+
+
 case class DegenerateSimplicialSetElement private (base: SimplicialSetElement, degeneracies: List[Int]) extends SimplicialSetElement:
   override def dim: Int = base.dim + degeneracies.size
+  override def toString: String =
+    degeneracies
+      .map{d => "s" + unicode.unicodeSubScript(d.toString)}
+      .mkString("", "  ", " ") + s"${base.toString}"
+  // Handle cases where the base is also a simplicial set element
+  def join(): DegenerateSimplicialSetElement = base match {
+    case oldbase : DegenerateSimplicialSetElement =>
+      val newbase : SimplicialSetElement = oldbase.join()
+      new DegenerateSimplicialSetElement(
+        newbase.base,
+        DegenerateSimplicialSetElement.normalizeDegeneracies(newbase.degeneracies, degeneracies)
+      )
+    case _ => this
+  }
 
 object DegenerateSimplicialSetElement {
   @tailrec
@@ -32,7 +55,7 @@ object DegenerateSimplicialSetElement {
 }
 
 trait SimplicialSetLike {
-  def generators : List[SimplicialSetElement]
+  def generators : Seq[SimplicialSetElement]
   def face(index: Int): PartialFunction[SimplicialSetElement,SimplicialSetElement]
   def contains(sse: SimplicialSetElement): Boolean
 
@@ -43,15 +66,23 @@ trait SimplicialSetLike {
         n.to(0, -1)
           .combinations(n - g.dim)
           .map { deg => DegenerateSimplicialSetElement(g.base, (deg ++ g.degeneracies).toList) }
-      }
+      }.toList
 }
 
-case class SimplicialSet private (generators: List[SimplicialSetElement],
+/**
+ * Class contract: generators are ascending in dimension. This is assumed by things like n_skeleton
+ * @param generators
+ * @param faceMapping
+ */
+case class SimplicialSet (generators: LazyList[SimplicialSetElement],
                          faceMapping: PartialFunction[SimplicialSetElement, List[SimplicialSetElement]])
   extends SimplicialSetLike {
   override def face(index: Int): PartialFunction[SimplicialSetElement, SimplicialSetElement] = {
     case sse if (0 <= index) && (index <= sse.dim) => {
-      val split = sse.degeneracies.groupBy { j => (index < j, index == j || index == j + 1, index > j + 1) }
+      val split = sse
+        .degeneracies
+        .groupBy { j => (index < j, index == j || index == j + 1, index > j + 1) }
+        .orElse { _ => List.empty }
       inline def INDEX_SMALL = (true, false, false)
       inline def INDEX_MATCH = (false, true, false)
       inline def INDEX_LARGE = (false, false, true)
@@ -62,20 +93,52 @@ case class SimplicialSet private (generators: List[SimplicialSetElement],
 
       val newDegeneracies: List[Int] =
         split(INDEX_SMALL).map(_ - 1) ++ // degeneracy indices decrease by one each time they commute with smaller face index
-          split(INDEX_MATCH).tail ++ // one of the up to two matching degeneracies go away
+          split(INDEX_MATCH).drop(1) ++ // one of the up to two matching degeneracies go away
           split(INDEX_LARGE) // degeneracies unchanged once the face index is large enough
 
       newFace match {
         case None =>
           DegenerateSimplicialSetElement(sse.base, newDegeneracies)
         case Some(i) =>
-          DegenerateSimplicialSetElement(faceMapping.apply(sse.base)(i), newDegeneracies)
+          DegenerateSimplicialSetElement(faceMapping.apply(sse.base)(i), newDegeneracies).join()
       }
     }
   }
 
+  given sseOrdering : Ordering[SimplicialSetElement] =
+    Ordering.by{(sse : SimplicialSetElement) => sse.dim}.orElseBy(_.hashCode())
+  given normalizedHomologyCell(using seOrd : Ordering[SimplicialSetElement]) : OrderedCell[SimplicialSetElement] with {
+    extension (t: SimplicialSetElement) {
+      override def boundary[CoefficientT](using fr: Fractional[CoefficientT]): Chain[SimplicialSetElement, CoefficientT] =
+        if(t.dim == 0) Chain() else {
+        val pmOne = List(fr.one, fr.negate(fr.one))
+        if (t.degeneracies.isEmpty) Chain.from(
+          (0 to t.base.dim)
+            .map { i => face(i)(t) }
+            .zipWithIndex
+            .filter { (s,i) => s.degeneracies.isEmpty }
+            .map { (s, i) => s -> pmOne(i%2) }
+          )
+        else Chain()
+      }
+      override def dim: Int = t.base.dim + t.degeneracies.size
+    }
+
+    override def compare(x: SimplicialSetElement, y: SimplicialSetElement): Int =
+      seOrd.compare(x,y)
+  }
+
   override def contains(sse: SimplicialSetElement): Boolean =
     generators.exists{g => g.base == sse.base}
+
+  def f_vector(maxDim : Int = 10): Seq[Int] =
+    (0 to maxDim)
+      .map(generators
+        .takeWhile(_.dim <= maxDim)
+        .groupBy{g => g.dim}
+        .orElse{ _ => LazyList.empty }
+      )
+      .map(_.size)
 }
 
 
@@ -84,11 +147,33 @@ object SimplicialSet {
             faceMapping: PartialFunction[SimplicialSetElement, List[SimplicialSetElement]]): SimplicialSet = {
     assert(generators.forall { g => faceMapping.isDefinedAt(g) })
     assert(generators.forall { g => if (g.dim > 0) faceMapping.apply(g).size == g.dim + 1 else true })
-    new SimplicialSet(generators, faceMapping)
+    new SimplicialSet(LazyList.from(generators), faceMapping)
   }
 }
 
 
+/** Notes on boundaries of degenerate elements
+ *
+ * Suppose w = s_j_ z
+ * Suppose k > j+1, i < j
+ * Then by the simplicial set laws:
+ * d_k_ s_j_ z = s_j d_k-1_ z
+ * d_i_ s_j_ z = s_j-1_ d_i_ z
+ * d_j_ s_j_ z = d_j+1_ s_j_ z = z
+ *
+ * So from all this follows:
+ * ∂ s_j_ z = ∑ (-1)^n^ d_i_ s_j_ z =
+ *   s_j_ (d_0_ - d_1_ +  ... ± d_j-1_) z + (z - z) ± s_j_ (d_j+1_ - d_j+2_ + ... ± d_n_) z =
+ *   s_j_ ∂ z (??? not sure it all lines up right here?)
+ *
+ *
+ *   Goerss-Jardine III:2, esp. Theorem 2.1:
+ *
+ *   Normalized Chain Complex NA is isomorphic to the quotient A/DA chain complex clearing out the
+ *   degenerate simplices.
+ */
+
+
 /** Normalization is with respect to simplicial identities:
  *
  * Face(i)Face(j) = Face(j-1)Face(i) if i < j
@@ -109,11 +194,36 @@ object SimplicialSet {
  * Constructions
  */
 
+import math.Ordering.Implicits.sortedSetOrdering
+class Singular[VertexT:Ordering](val underlying : Seq[Simplex[VertexT]]) extends SimplicialSetLike {
+  val generators : Seq[SimplicialWrapper[Simplex[VertexT]]] = underlying.sortBy(_.dim).map { spx => SimplicialWrapper(spx) }
+  val faceMapping : Map[SimplicialSetElement, List[SimplicialSetElement]] = Map.from(
+    generators.map { spx => spx ->
+      (0 to spx.dim).map { i =>
+        SimplicialWrapper(Simplex(spx.wrapped.vertices.drop(i)))
+      }.toList
+    }
+  )
+  val simplicialSet = SimplicialSet(LazyList.from(generators), faceMapping)
+
+  override def face(index: Int): PartialFunction[SimplicialSetElement, SimplicialSetElement] =
+    simplicialSet.face(index)
+  override def contains(sse: SimplicialSetElement): Boolean =
+    simplicialSet.contains(sse)
+}
+
+
+
 def nSkeleton(n:Int)(ss: SimplicialSet) : SimplicialSet =
-  SimplicialSet(ss.generators.filter(_.dim <= n), ss.faceMapping)
+  SimplicialSet(ss.generators.takeWhile(_.dim <= n), ss.faceMapping)
+
+
+def alternateLazyLists[A](left: LazyList[A], right: LazyList[A]): LazyList[A] =
+  if(left.isEmpty) right
+  else left.head #:: alternateLazyLists(right, left.tail)
 
 case class DisjointUnion(left: SimplicialSet, right: SimplicialSet) extends SimplicialSetLike {
-  override def generators: List[SimplicialSetElement] = left.generators ++ right.generators
+  override def generators: LazyList[SimplicialSetElement] = alternateLazyLists(left.generators, right.generators)
   override def face(index: Int): PartialFunction[SimplicialSetElement, SimplicialSetElement] = {
     case sse if(left.face(index).isDefinedAt(sse)) => left.face(index)(sse)
     case sse if(right.face(index).isDefinedAt(sse)) => right.face(index)(sse)
@@ -124,7 +234,7 @@ case class DisjointUnion(left: SimplicialSet, right: SimplicialSet) extends Simp
 
 
 case class SubSimplicialSet(ss: SimplicialSet, ambient: SimplicialSet, inclusion: PartialFunction[SimplicialSetElement,SimplicialSetElement]) extends SimplicialSetLike:
-  override def generators: List[SimplicialSetElement] = ss.generators
+  override def generators: LazyList[SimplicialSetElement] = ss.generators
   override def face(index: Int): PartialFunction[SimplicialSetElement, SimplicialSetElement] = ss.face(index)
   override def contains(sse: SimplicialSetElement): Boolean = ss.contains(sse)
 
diff --git a/src/test/scala/org/appliedtopology/tda4j/SimplicialSetSpec.scala b/src/test/scala/org/appliedtopology/tda4j/SimplicialSetSpec.scala
new file mode 100644
index 00000000..d095e0d6
--- /dev/null
+++ b/src/test/scala/org/appliedtopology/tda4j/SimplicialSetSpec.scala
@@ -0,0 +1,33 @@
+package org.appliedtopology.tda4j
+
+import org.specs2.mutable
+
+class SimplicialSetSpec extends mutable.Specification {
+
+  "Homology of 4-sphere" >> {
+    val sphere4 = sphere(4)
+    import sphere4.given
+    given chc: CellularHomologyContext[SimplicialSetElement, Double, Double]()
+    import chc.{*, given}
+
+    val sphere4stream = new CellStream[SimplicialSetElement, Double] {
+      override def filtrationValue = _ => 0.0
+
+      override def iterator = sphere4.generators.iterator
+
+      override def filtrationOrdering = sseOrdering
+
+      override def smallest = Double.NegativeInfinity
+
+      override def largest = Double.PositiveInfinity
+    }
+
+    val ph = persistentHomology(sphere4stream)
+    val dgm = ph.diagramAt(1.0)
+
+    dgm must containTheSameElementsAs(List(
+      (0, 0.0, Double.PositiveInfinity),
+      (4, 0.0, Double.PositiveInfinity)
+    ))
+  }
+}