My attempt at http://www.azspcs.net/Contest/Factorials/

Factorials.scala

import collection.mutable
import collection.JavaConversions._
import org.mapdb.DBMaker   // Great disk based map -  https://github.com/jankotek/MapDB

/**
 * My attempt at http://www.azspcs.net/Contest/Factorials/
 * Brief explanation:
 * Brute force BFS (ignoring non-positive numbers and duplicates) upto length `fbf`
 * Higher the `fbf`, the better quality result. I managed to get to `fbf`=8 before running out of memory
 * Beyond `fbf`, try to multiply your way to n! Since the last few steps involve only multiplication,
 * We can safely ignore all sequences that do not have the right number of prime factors of n!
 * (it is trivial to prime factorize n! since we know it has no prime factors >n)
 * This simple heuristic is enough to get you a score of 20+
 *
 * Hypothetical Bests:
 * (record: 17->12(13),19->13(14),20->14(15),21->14(16),22->14(16),23->15(17),24->15(17),25->15(17), 26->15(18), 27->16(19), 28->16(19), 29->17(21), 30->17(21), 31->18(24), 32->18(22) , 33->18(24), 34->18(24), 35->19(27), 36->19(27), 37->20(28))
 *
 * Some resources:
 * OEIS: http://oeis.org/A216999 http://oeis.org/A001035 http://oeis.org/A027423 http://oeis.org/A173419 http://oeis.org/A217032
 * Forums: http://tech.groups.yahoo.com/group/AlZimmermannsProgrammingContests/messages http://dxdy.ru/topic67743-315.html http://e-science.ru/forum/index.php?showtopic=36028&st=380
 */
object Factorials extends App {
  val primes: Seq[BigInt] = Seq(2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37)

  def factorial(n: Int): BigInt = if (n < 1) 1 else n * factorial(n - 1)

  def multiplicity(n: BigInt)(p: BigInt): Int = if (n % p == 0) 1 + multiplicity(n / p)(p) else 0

  def multiply(primeFactorization: Seq[Int]) = primes zip primeFactorization map {i => i._1 pow i._2} product

  type SLP = Seq[BigInt]    // a straight line program

  def factors(n: BigInt, c: Int = fbf) = {           // select c factorizations of n (defaults to fbf)
    val f = primes map multiplicity(n) map {1+}      // factoring n is same as 2-partitioning prime-factors of n
    val p = f.scanLeft(1){_ * _}                     // optimized way of running-product of f
    val t = p.last                                   // total t factorizations possible (t/2 symmetric)

    def factor(k: Int) = multiply({0 until f.length} map {i => k / p(i) % f(i)})      //kth factor of n

    for (i <- Math.max(t/2 - c, 1) to t/2) yield (factor(i), factor(t - i - 1))    // n = factor(i) * factor(t - i - 1)
  }

  def bfs(depth: Int) = {                               // bfs till certain depth
    var s = Set(Seq(BigInt(1)))
    val known = mutable.Map(BigInt(1) -> s)             // known contains known shortest slps for each key

    val operations = Array[(BigInt, BigInt) => BigInt](_ * _, _ + _, _ - _)

    def operate(slp: SLP) = {
      val seen = mutable.Set[BigInt]()
      for {
        i <- 0 until slp.length
        j <- 0 to i
        op <- operations
        v = op(slp(i), slp(j))
        if v>0
        if !(known contains v)
        if seen.add(v)
      } yield slp :+ v
    }

    for (i <- 1 to depth) {
      val (t, time, mem) = profile {s flatMap operate}
      known ++= t groupBy {_.last}
      val (width, total) = (s.size, known.size)
      println(f"BFS {depth = $i%2d, width = $width%7d, total = $total%9d, memory = $mem%4d MB, time = $time%5.0fs}")
      s = t
    }
    known
  }

  val cache: mutable.Map[BigInt, SLP] = DBMaker.newTempHashMap[BigInt, SLP]()

  def memo(n: BigInt) = cache getOrElseUpdate (n, (factors(n) map dp minBy {_.size}) :+ n)

  def dp(f: (BigInt, BigInt)) = merge(getKnown(f._1), getKnown(f._2))

  def getKnown(n: BigInt): Set[SLP] = known getOrElse (n, Set(memo(n)))

  def merge(a: Set[SLP], b: Set[SLP]) = {for (i <- a; j <- b) yield (i ++ j).distinct} minBy {_.size}

  def solve(N: Int)  = {
    val (s, time, mem) = profile { memo(factorial(N)) }
    val (output, len) = (s.mkString(", "), s.size-1)
    println(f"$N%3d {length = $len%2d, cache = ${cache.size}%7d, memory = $mem%5d MB, time = $time%10.0fs}: $output%s")
    output
  }

  import System.{currentTimeMillis => cTime}
  println(f"Running with memory = ${Runtime.getRuntime.maxMemory>>30}GB")
  def profile[R](code: => R, t: Long = cTime) = (code, (cTime - t)/1000.0, Runtime.getRuntime.totalMemory>>20)

  // start here
  val fbf = 10000      // select upto this many factorizations aka factorization-branching-factor
  val known = bfs(depth = 8)
  var output = {for {i <- 2 to 37; result = solve(i); if i>=13} yield result} mkString (";")
  println(output)
}