Skip to content

Commit

Permalink
Refactored Simplex significantly. No longer implementing SortedSet, a…
Browse files Browse the repository at this point in the history
…nd now based on a Case Class.
  • Loading branch information
Mikael Vejdemo-Johansson committed May 24, 2024
1 parent 88d3675 commit 5869658
Show file tree
Hide file tree
Showing 12 changed files with 158 additions and 278 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -62,14 +62,15 @@ object FiniteMetricSpace {
val metricSpace: FiniteMetricSpace[VertexT]
) extends PartialFunction[Simplex[VertexT], Double] {
def isDefinedAt(spx: Simplex[VertexT]): Boolean =
spx.forall(v => metricSpace.contains(v))
spx.vertices.forall(v => metricSpace.contains(v))

def apply(spx: Simplex[VertexT]): Double =
if (spx.size <= 1)
if (spx.dim <= 0)
0.0
else
spx
.flatMap(v => spx.filter(_ > v).map(w => metricSpace.distance(v, w)))
.vertices
.flatMap(v => spx.vertices.filter(_ > v).map(w => metricSpace.distance(v, w)))
.max
}
}
Expand Down
6 changes: 3 additions & 3 deletions src/main/scala/org/appliedtopology/tda4j/Homology.scala
Original file line number Diff line number Diff line change
Expand Up @@ -182,11 +182,11 @@ class SimplicialHomologyByDimensionContext[VertexT: Ordering, CoefficientT: Frac
// secondly, we can read off homology completely from a minimal spanning tree
val kruskal = new Kruskal[Simplex[VertexT]](
cycles.keys.toSeq,
{ (x: Simplex[VertexT], y: Simplex[VertexT]) => stream.filtrationValue(x.union(y)) }
{ (x: Simplex[VertexT], y: Simplex[VertexT]) => stream.filtrationValue(Simplex.from(x.vertices ++ y.vertices)) }
)(using stream.filtrationOrdering)

kruskal.mstIterator.foreach { (src, tgt) =>
val edge: Simplex[VertexT] = src.union(tgt)
val edge: Simplex[VertexT] = Simplex.from(src.vertices ++ tgt.vertices)

// the edge src -- tgt will connect src to tgt thus removing one of the cycles
val dEdge = edge.boundary
Expand All @@ -199,7 +199,7 @@ class SimplicialHomologyByDimensionContext[VertexT: Ordering, CoefficientT: Frac
}

kruskal.cyclesIterator.foreach { (src, tgt) =>
val edge: Simplex[VertexT] = src.union(tgt)
val edge: Simplex[VertexT] = Simplex.from(src.vertices ++ tgt.vertices)

// the edge src -- tgt will connect src to tgt thus closing a loop
val dEdge = edge.boundary
Expand Down
30 changes: 16 additions & 14 deletions src/main/scala/org/appliedtopology/tda4j/RipserStream.scala
Original file line number Diff line number Diff line change
Expand Up @@ -40,9 +40,9 @@ class SimplexIndexing(val vertexCount: Int) {
}

def cofacetIterator(simplex: Simplex[Int]): Iterator[Int] =
cofacetIterator(apply(simplex), simplex.size, true)
cofacetIterator(apply(simplex), simplex.vertices.size, true)
def topCofacetIterator(simplex: Simplex[Int]): Iterator[Int] =
cofacetIterator(apply(simplex), simplex.size, false)
cofacetIterator(apply(simplex), simplex.vertices.size, false)
def cofacetIterator(
index: Int,
size: Int,
Expand All @@ -60,7 +60,7 @@ class SimplexIndexing(val vertexCount: Int) {
) { (s, iB, iA, k, j) =>
if (j < 0) {
None // end iteration when we're done
} else if (s.contains(j)) {
} else if (s.vertices.contains(j)) {
if (!allCofacets)
None
else
Expand All @@ -78,7 +78,7 @@ class SimplexIndexing(val vertexCount: Int) {
.map((os: Option[Int]) => os.get)

def facetIterator(index: Int, size: Int): Iterator[Int] =
Iterator.unfold((apply(index, size).toSeq.sorted, index, 0, size - 1)) { (s, iB, iA, k) =>
Iterator.unfold((apply(index, size).vertices.toSeq.sorted, index, 0, size - 1)) { (s:Seq[Int], iB:Int, iA:Int, k:Int) =>
if (k < 0) None
else {
val j = s(k)
Expand All @@ -89,8 +89,8 @@ class SimplexIndexing(val vertexCount: Int) {
}

def apply(simplex: Simplex[Int]): Int =
simplex.toSeq.sorted.reverse.zipWithIndex.map { (v, i) =>
binomial(v, simplex.size - i)
simplex.vertices.toSeq.sorted.reverse.zipWithIndex.map { (v, i) =>
binomial(v, simplex.vertices.size - i)
}.sum
}

Expand Down Expand Up @@ -181,11 +181,11 @@ class RipserStreamSparse(
fV <- simplexCache.map(_._1).iterator
previousSimplex <- simplexCache.filter(_._1 < fV).iterator.map(_._2)
nextVertex <- metricSpace.elements
.filter(!previousSimplex.contains(_))
.filter(nV => previousSimplex.map(oV => metricSpace.distance(oV, nV)).max <= fV)
simplex: Simplex[Int] <- List(previousSimplex + nextVertex)
if zeroApparentCofacet(si(simplex), simplex.size).isEmpty
if zeroApparentFacet(si(simplex), simplex.size).isEmpty
.filter(!previousSimplex.vertices.contains(_))
.filter(nV => previousSimplex.vertices.map(oV => metricSpace.distance(oV, nV)).max <= fV)
simplex: Simplex[Int] = Simplex.from(previousSimplex.vertices.appended(nextVertex))
if zeroApparentCofacet(si(simplex), simplex.vertices.size).isEmpty
if zeroApparentFacet(si(simplex), simplex.vertices.size).isEmpty
// also check if this is cleared?
yield simplex
}.iterator
Expand Down Expand Up @@ -294,10 +294,12 @@ class RipserStreamOf[VertexT: Ordering](
RipserStream(intMetricSpace, maxFiltrationValue, maxDimension)

override def iterator: Iterator[Simplex[VertexT]] =
rs.iterator.map(s => s.map(v => vertices(v)))
rs.iterator.map(s => Simplex.from(s.vertices.map(v => vertices(v))))

override def filtrationValue: PartialFunction[Simplex[VertexT], Double] =
rs.filtrationValue.compose(s => s.map(v => vertices.indexOf(v)))
override def filtrationValue: PartialFunction[Simplex[VertexT], Double] = { spx =>
val indices : Seq[Int] = spx.vertices.map(vertices.indexOf)
rs.filtrationValue(Simplex.from(indices))
}
}

class SymmetricRipserCliqueFinder[KeyT](
Expand Down
152 changes: 45 additions & 107 deletions src/main/scala/org/appliedtopology/tda4j/Simplex.scala
Original file line number Diff line number Diff line change
Expand Up @@ -4,130 +4,68 @@ import scala.collection.{mutable, SortedIterableFactory, SortedSetFactoryDefault
import scala.collection.immutable.{Set, SortedMap, SortedSet, SortedSetOps, TreeSet}
import scala.math.Ordering.IntOrdering
import scala.math.Ordering.Double.IeeeOrdering
import math.Ordering.Implicits.sortedSetOrdering

given simplexOrdering[VertexT: Ordering]: Ordering[Simplex[VertexT]] = sortedSetOrdering[Simplex, VertexT]

/** Class representing an abstract simplex. Abstract simplices are sets (of totally ordered vertices) and thus are
* represented by a type that implements the interface of a `SortedSet` fully, and also inherits from `Cell` so that
* the class has a `boundary` method.
*
* With this `SortedSet` interface in place, there is scope for using `Simplex` for combinatorial algebraic topology
* tasks unrelated to persistence, and possibly for writing persistence algorithms more smoothly.
*
* You should never have reason to use the constructor directly (...and if you do, you should make sure to give the
* internal `SortedSet` yourself) - instead use the factory method in the companion object. In code this means that
* instead of `new Simplex[T](a,b,c)` you would write `Simplex[T](a,b,c)`.
*
* @param vertexSet
* Vertices of the simplex
* @param ordering
* Ordering of the vertex type
* @tparam VertexT
* Vertex type
*/
class Simplex[VertexT](protected val vertexSet: SortedSet[VertexT])( //vertexSet variable defined here
using val ordering: Ordering[VertexT] // ordering definied here
) extends Cell[Simplex[VertexT]]
with SortedSet[VertexT]
with SortedSetOps[VertexT, Simplex, Simplex[VertexT]]
with SortedSetFactoryDefaults[VertexT, Simplex, Set] {
self => // gives methods access to the object that's calling it in the first place
import math.Ordering.Implicits.seqOrdering

given simplexOrdering[VertexT: Ordering]: Ordering[Simplex[VertexT]] =
Ordering.by{(spx:Simplex[VertexT]) => spx.vertices}(seqOrdering[Seq, VertexT])

/** Class representing an abstract simplex. Abstract simplices are given by sets (of totally ordered vertices)
* and inherit from `Cell` so that the class has a `boundary` and a `dim` method.
*
* You should never have reason to use the constructor directly (...and if you do, you should make sure to give the
* internal `SortedSet` yourself) - instead use the factory method in the companion object. In code this means that
* instead of `new Simplex[T](a,b,c)` you would write `Simplex[T](a,b,c)`.
*
* @param vertices
* Vertices of the simplex
* @param ordering
* Ordering of the vertex type
* @tparam VertexT
* Vertex type
*/
case class Simplex[VertexT : Ordering] private (vertices : VertexT*) extends Cell[Simplex[VertexT]] {
override def equals(obj: Any): Boolean = obj match {
case other : Simplex[VertexT] => vertices == other.vertices
case _ => super.equals(obj)
}

override def dim: Int = vertices.size-1

override def toString(): String =
vertexSet.mkString(s"∆(", ",", ")")

// ***** Simplex specific operations
vertices.mkString(s"∆(", ",", ")")

/** Implementation of the face enumeration for the simplicial boundary map: exclude one element at a time, and use
* alternating signs.
*
* @return
* A sequence of boundary cells.
* @todo
* Change `List` to `Chain` once we have an implementation of `Chain`
*/
override def boundary[CoefficientT](using
fr: Fractional[CoefficientT]
): Chain[Simplex[VertexT], CoefficientT] =
if (dim == 0) Chain()
fr: Fractional[CoefficientT]
): Chain[Simplex[VertexT], CoefficientT] =
if (dim <= 0) Chain()
else
Chain.from(
self
vertices
.to(Seq)
.map(vtx => self - vtx)
.zipWithIndex
.map((vtx,i) => Simplex.from(vertices.drop(i)))
.zip(Iterator.unfold(fr.one)(s => Some((s, fr.negate(s)))))
)

def dim: Int = size - 1

// ***** Overriding for inheriting and extending standard library constructions
override def className = "Simplex"

override def iterator: Iterator[VertexT] = vertexSet.iterator

override def excl(elem: VertexT): Simplex[VertexT] =
new Simplex[VertexT](
vertexSet.excl(elem)
)

override def incl(elem: VertexT): Simplex[VertexT] =
new Simplex[VertexT](
vertexSet.incl(elem)
)

override def contains(elem: VertexT): Boolean = vertexSet.contains(elem)
def union(other: Simplex[VertexT]) =
new Simplex(vertices.union(other.vertices))
}

override def rangeImpl(
from: Option[VertexT],
until: Option[VertexT]
): Simplex[VertexT] =
new Simplex[VertexT](vertexSet.rangeImpl(from, until))
/** Simplex companion object with factory methods
*/
object Simplex {
def apply[VertexT: Ordering](vertices: VertexT*) =
new Simplex[VertexT](Seq.from(vertices.sorted) : _*)

override def iteratorFrom(start: VertexT): Iterator[VertexT] =
vertexSet.iteratorFrom(start)
def empty[VertexT: Ordering]: Simplex[VertexT] =
new Simplex[VertexT]()

override def sortedIterableFactory: SortedIterableFactory[Simplex] =
Simplex
def from[VertexT: Ordering](source: IterableOnce[VertexT]): Simplex[VertexT] =
new Simplex(Seq.from(source.iterator.toSeq.sorted) : _*)
}

/** Convenience method for defining simplices
*
* The character ◊, used to avoid hogging `s`, is typed as Alt+Shift+V on Mac GB layout, and has Unicode code 0x25CA.
*
* The character ∆ is typed as Alt+J on Mac GB layout, and has unicode code 0x0394.
*
* The character σ has Windows Alt-code 229, and unicode code 0x03c3.
*/
def s[T: Ordering](ts: T*): Simplex[T] = Simplex.from(ts)
def [T: Ordering](ts: T*): Simplex[T] = Simplex.from(ts)
def [T: Ordering](ts: T*): Simplex[T] = Simplex.from(ts)
def σ[T: Ordering](ts: T*): Simplex[T] = Simplex.from(ts)

/** Simplex companion object with factory methods
*/
object Simplex extends SortedIterableFactory[Simplex] {
override def apply[VertexT: Ordering](vertices: VertexT*) =
new Simplex[VertexT](TreeSet[VertexT](vertices: _*))
override def empty[VertexT: Ordering]: Simplex[VertexT] =
new Simplex[VertexT](TreeSet[VertexT]())

override def from[VertexT: Ordering](
source: IterableOnce[VertexT]
): Simplex[VertexT] =
(newBuilder[VertexT] ++= source).result()

override def newBuilder[VertexT: Ordering]: mutable.Builder[VertexT, Simplex[VertexT]] =
new mutable.ReusableBuilder[VertexT, Simplex[VertexT]] {
private val elems = mutable.Set[VertexT]()

override def clear(): Unit = elems.clear()

def addOne(elem: VertexT): this.type =
elems += elem; this

override def result(): Simplex[VertexT] =
new Simplex[VertexT](TreeSet[VertexT](elems.to(Seq): _*))
}

}
23 changes: 11 additions & 12 deletions src/main/scala/org/appliedtopology/tda4j/SimplexStream.scala
Original file line number Diff line number Diff line change
Expand Up @@ -142,12 +142,11 @@ class ExplicitStreamBuilder[VertexT: Ordering, FiltrationT](using
}

override def result(): ExplicitStream[VertexT, FiltrationT] = {
def lt(
fs1: (FiltrationT, Simplex[VertexT]),
fs2: (FiltrationT, Simplex[VertexT])
): Boolean =
(fs1._1, fs1._2.to(Seq)) < (fs2._1, fs2._2.to(Seq))
simplices.sortInPlaceWith(lt)
given filtrationOrdering : Ordering[(FiltrationT, Simplex[VertexT])] =
Ordering
.by[(FiltrationT,Simplex[VertexT]),FiltrationT]{(f:FiltrationT,s:Simplex[VertexT]) => f}(ordering)
.orElseBy{(f:FiltrationT,s:Simplex[VertexT]) => s}(simplexOrdering[VertexT])
simplices.sortInPlace

new ExplicitStream(filtrationValues.toMap, simplices.map((_, s) => s).toSeq)(using filterable)
}
Expand All @@ -170,21 +169,21 @@ class FilteredSimplexOrdering[VertexT, FiltrationT](
case (x, y) if filtration.filtrationValue.isDefinedAt(x) && filtration.filtrationValue.isDefinedAt(y) =>
filtrationOrdering.compare(filtration.filtrationValue(x), filtration.filtrationValue(y)) match {
case 0 =>
if (Ordering.Int.compare(x.size, y.size) == 0)
if (Ordering.Int.compare(x.dim, y.dim) == 0)
Ordering.Implicits
.seqOrdering[Seq, VertexT](vertexOrdering)
.compare(x.to(Seq), y.to(Seq))
.compare(x.vertices, y.vertices)
else
Ordering.Int.compare(x.size, y.size)
Ordering.Int.compare(x.dim, y.dim)
case cmp if cmp != 0 => cmp
}
case (x, y) => // at least one does not have a filtration value defined; just go by dimension and lexicographic
if (Ordering.Int.compare(x.size, y.size) == 0)
if (Ordering.Int.compare(x.dim, y.dim) == 0)
Ordering.Implicits
.seqOrdering[Seq, VertexT](vertexOrdering)
.compare(x.to(Seq), y.to(Seq))
.compare(x.vertices, y.vertices)
else
Ordering.Int.compare(x.size, y.size)
Ordering.Int.compare(x.dim, y.dim)
}
}

Expand Down
10 changes: 5 additions & 5 deletions src/main/scala/org/appliedtopology/tda4j/SymmetryGroup.scala
Original file line number Diff line number Diff line change
Expand Up @@ -56,14 +56,14 @@ trait SymmetryGroup[KeyT, VertexT: Ordering]() {
def orbitSeq(
simplex: Simplex[VertexT]
): Set[Simplex[VertexT]] =
keys.map(k => simplex.map(apply(k))).toSet
keys.map(k => Simplex.from(simplex.vertices.map(apply(k)))).toSet

def orbitPar(
simplex: Simplex[VertexT]
): Set[Simplex[VertexT]] = {
val futures: Iterable[Future[Simplex[VertexT]]] =
for (k <- keys) yield Future {
simplex.map(apply(k))
Simplex.from(simplex.vertices.map(apply(k)))
}

val allfutures = Future.sequence(futures)
Expand All @@ -82,7 +82,7 @@ trait SymmetryGroup[KeyT, VertexT: Ordering]() {
def representative(
simplex: Simplex[VertexT]
): Simplex[VertexT] =
keys.map(k => simplex.map(apply(k))).min
keys.map(k => Simplex.from(simplex.vertices.map(apply(k)))).min

/** Check if `simplex` is the canonical representative of its own orbit.
*
Expand Down Expand Up @@ -490,7 +490,7 @@ class HyperCubeSymmetryGeneratorsBitSet(val bitlength: Int) extends HyperCubeSym
if (representatives.contains(simplex)) {
simplex == representatives(simplex)
} else {
if (generators.forall(g => simplex <= simplex.map(s => g(s)))) {
if (generators.forall(g => simplex <= Simplex.from(simplex.vertices.map(s => g(s))))) {
// simplex is a pseudo-minimum
// time to check the entire orbit
representatives(simplex) = super.representative(simplex)
Expand Down Expand Up @@ -524,7 +524,7 @@ class HyperCubeSymmetryGenerators(val bitlength: Int) extends HyperCubeSymmetry(
if (representatives.contains(simplex)) {
simplex == representatives(simplex)
} else {
if (generators.par.forall(g => simplex <= simplex.map(s => g(s)))) {
if (generators.par.forall(g => simplex <= Simplex.from(simplex.vertices.map(s => g(s))))) {
// simplex is a pseudo-minimum
// time to check the entire orbit
representatives(simplex) = super.representative(simplex)
Expand Down
Loading

0 comments on commit 5869658

Please sign in to comment.