use Pochhammer.frink

/** Lebesgue elliptic integral of the first kind.  This is like the simplest
    series I found.  In fact, it's the only series expansion that contains only
    *one* infinite loop and not *two* infinite loops nested LOL.

    https://functions.wolfram.com/EllipticIntegrals/EllipticF/06/01/03/
    https://functions.wolfram.com/EllipticIntegrals/EllipticF/06/01/03/0003/

    Typeset equation for this series:
    https://functions.wolfram.com/EllipticIntegrals/EllipticF/06/01/03/0003/MainEq1.L.gif

    When m is close to 1, this requires working with like 100 digits of precision    and/or exact rational arguments or it blows up.
*/
ellipticIntegralF[z, m] :=
{
   sum = z

   KLOOP:                   // THINK ABOUT:
   for k = 1 to 80 step 2   // When arguments are real, even k adds 0 to sum
   {
      zk = z^k
      sumq = 0
      for q = 1 to k-1
      {
         sj = 0
         for j = 0 to q
         {
            si = 0
            for ii = 0 to j
               si = si + binomial[j,ii] (2 ii - j)^(k-1) zk
            
            sj = sj + binomial[q,j] m^j (2-m)^(q-j) 2^(k-j-q-1) si
         }
         sumq = sumq + (-1)^q/(2 q + 1) binomial[k-1,q] sj
         println["sumq is $sumq"]
      }

//      println["i=$i k=$k, Thing = " + (i^(k-1))/(k!)]
      sk = i^(k-1)/k! generalizedBinomial[k-1/2, k-1] sumq

      println["k=$k, sk = $sk, sum = $sum"]
      sum = sum + sk
   }

   return sum
}