/*
 * This file:
 *   http://anggtwu.net/MAXIMA/laurent1.mac.html
 *   http://anggtwu.net/MAXIMA/laurent1.mac
 *          (find-angg "MAXIMA/laurent1.mac")
 * Author: Eduardo Ochs <eduardoochs@gmail.com>
 * See: (find-es "maxima" "laurent")
 *      http://anggtwu.net/eev-maxima.html
 *
 * Dot notation for Laurent polynomials, version 1.
 * Idea: 456.78 is
 *
 *   4*10^2 + 5*10^1 + 6*10^0 + 7*10^-1 + 8*10^-2,
 *
 * and we write it (roughly) as:
 *
 *   [ 4  5  6  .  7  8 ]
 *
 * We can generalize the 10 to x.
 * With the functions of this file we have:
 *
 *   (%i2) lpx(4*x^2 + 5*x^1 + 6*x^0 + 7*x^-1 + 8*x^-2);
 *   (%o2)           [ 4  5  6  .  7  8 ]
 *
 * The "[ 4 5 6 . 7 8 ]" is a horizontal matrix that has the string
 * "." in its fourth position; the "." separates the "integer part"
 * from the "fractional part".
 *
 * These Laurent polynomials are also great for teaching trigonomotric
 * identities:
 *
 *   If we define          E = e^(ix),
 *                   icos(x) = E + E^-1,
 *            and    isin(x) = E - E^-1,
 *   then we have     cos(x) =    1/2 * icos(x)
 *            and     sin(x) = 1/(2i) * isin(x).
 *
 * Using Laurent polynomials on E instead of on x we have:
 *
 *   (%i1) lpE(icos(x));
 *   (%o1)                           [ 1  0  .  1 ]
 *   (%i2) lpE(isin(x));
 *   (%o2)                          [ 1  0  .  - 1 ]
 *   (%i3) lpE(isin(3*x));
 *   (%o3)                    [ 1  0  0  0  .  0  0  - 1 ]
 *   (%i4) lpE(isin(3*x) * icos(x));
 *   (%o4)                [ 1  0  1  0  0  .  0  - 1  0  - 1 ]
 *   (%i5) lpE(isin(4*x) + isin(2*x));
 *   (%o5)                [ 1  0  1  0  0  .  0  - 1  0  - 1 ]
 *   (%i6) lpe( sin(4*x) +  sin(2*x));
 *   (%o6)                         sin(4 x) + sin(2 x)
 *   (%i7) 
 *
 * See: (find-es "maxima" "laurent")
 * (defun e () (interactive) (find-angg "MAXIMA/laurent1.mac"))
*/

lpdegrees(lp,var) := block(
  [revlp,origposdeg,orignegdeg,posdeg,negdeg],
  revlp      : subst([var=var^-1], lp),
  origposdeg :  hipow(   lp, var),
  orignegdeg : -hipow(revlp, var),
  posdeg     : max(origposdeg, 0),
  negdeg     : min(orignegdeg, 0),
  [posdeg,negdeg])$

lpcoeffs0(lp,var,hilo) := makelist(ratcoef(lp,var,k), k, hilo[1],hilo[2], -1)$
lpcoeffs (lp,var)      := apply('lpcoeffs0, [lp,var, lpdegrees(lp,var)]);

lpdot    (lp,var)      := block([hi,lo,posdigits,negdigits],
  [hi,lo]    : lpdegrees(lp,var),
  posdigits  : lpcoeffs0(lp,var,[hi,0]),
  negdigits  : lpcoeffs0(lp,var,[-1,lo]),
  matrix(append(posdigits, ["."], negdigits))
  )$

Exponentialize(f) := subst([x=1/%i,%e=E], expand(exponentialize(f)))$

lpx(lp) := lpdot(expand(lp),x)$
lpE(f)  := lpdot(Exponentialize(f),E)$
lpe(f)  := expand(demoivre(expand(exponentialize(f))))$

/*
* (eepitch-maxima)
* (eepitch-kill)
* (eepitch-maxima)
load("laurent1.mac");

lpx(4*x^2 + 5*x^1 + 6*x^0 + 7*x^-1 + 8*x^-2);

lpx(x^4);
lpx(x^-4);

g1 : 4*x^2 + 5*x^1 + 6*x^0 + 7*x^-1 + 8*x^-2;
g1 : 4*x^2 + 5*x^1 + 6*x^0 + 7*x^-1 + 8*x^-2 + 9*x^-3;
lpx(g1);
lpx(g1 * x);
lpx(g1 * x^2);
lpx(g1 *     x     );
lpx(g1 *  10*x     );
lpx(g1 * (10*x + 1));

h : cos(x)^3;
Exponentialize(h);
lpE(h);
lpE(h * 8);

icos(x) := 2    * cos(x);
isin(x) := 2*%i * sin(x);

lpE( cos(x)   * 2);
lpE( cos(x)^2 * 4);
lpE( cos(x)^3 * 8);
lpE( sin(x)     *  2*%i);
lpE( sin(x)^2   * (2*%i)^2);
lpE( sin(x)^3   * (2*%i)^3);
lpE( sin(2*x)   * (2*%i));
lpE( sin(2*x)^2 * (2*%i)^2);

linenum : 0;
lpE(icos(x));
lpE(isin(x));
lpE(isin(3*x));
lpE(isin(3*x) * icos(x));
lpE(isin(4*x) + isin(2*x));
lpe( sin(4*x) +  sin(2*x));

*/