<- previous    index    next ->

Lecture 11a, More Complex Arithmetic


If you have code that works for double, extending to complex
may be easy.

Python evaluating a polynomial on "real" coefficients

test_polyfit.py polyfit in Python, needs numpy from numpy import array from numpy import polyfit from numpy import polyval import pylab def f(x): return 5.0 + 4.0*x + 3.0*x*x + 2.0*x*x*x + 1.0*x*x*x*x print "test_polyfit.py a x,y,5) " print "fit 5.0 + 4.0*x + 3.0*x*x + 2.0*x*x*x + 1.0*x*x*x*x" xx = [0.0 for i in range(5)] yy = [0.0 for i in range(5)] for i in range(5): xx[i] = 0.1*(i+1) yy[i] = f(xx[i]) # function we will fit print "i=", print i, print " ,xx=", print xx[i], print " ,yy=", print yy[i] print "p=polyfit(xx,yy,4) for 5 points" p=polyfit(xx,yy,4) # 5 input values, 4th order print "polyfit coefficients" print p print "backwards? largest power first, expected 5, 4, 3, 2, 1" print " " print "polyval values of fit" yy1=polyval(p,xx) print yy1 print "should be values above"

output

test_polyfit.py a x,y,5) fit 5.0 + 4.0*x + 3.0*x*x + 2.0*x*x*x + 1.0*x*x*x*x i= 0 ,xx= 0.1 ,yy= 5.4321 i= 1 ,xx= 0.2 ,yy= 5.9376 i= 2 ,xx= 0.3 ,yy= 6.5321 i= 3 ,xx= 0.4 ,yy= 7.2336 i= 4 ,xx= 0.5 ,yy= 8.0625 p=polyfit(xx,yy,4) for 5 points polyfit coefficients [ 1. 2. 3. 4. 5.] backwards? largest power first, expected 5, 4, 3, 2, 1 polyval values of fit [ 5.4321 5.9376 6.5321 7.2336 8.0625] should be values above

Python evaluating a polynomial on "complex" coefficients

test_cxpolyfit.py complex polyfit in Python, needs numpy from numpy import array from numpy import polyfit from numpy import polyval import pylab def f(x): # changed to some complex coefficients return (5.0+1.0j) + (4.0+1.1j)*x + (3.0+1.2j)*x*x + 2.0*x*x*x + 1.0*x*x*x*x print "test_cxpolyfit.py a x,y,5) using complex numbers " print "fit (5.0+1.0j) + (4.0+1.1j)*x + (3.0+1.2j)*x*x + 2.0*x*x*x + 1.0*x*x*x*x" xx = [0.0 for i in range(5)] yy = [0.0 for i in range(5)] for i in range(5): xx[i] = 0.1*(i+1) yy[i] = f(xx[i]) # now has complex coefficients print "i=", print i, print " ,xx=", print xx[i], print " ,yy=", print yy[i] print "p=polyfit(xx,yy,4) for 5 points" p=polyfit(xx,yy,4) # 5 input values, 4th order print "polyfit coefficients" print p print "backwards? largest power first, expected about 5, 4, 3, 2, 1" print " " print "polyval values of fit" yy1=polyval(p,xx) print yy1 print "should be values above"

output

test_cxpolyfit.py a x,y,5) using complex numbers fit (5.0+1.0j) + (4.0+1.1j)*x + (3.0+1.2j)*x*x + 2.0*x*x*x + 1.0*x*x*x*x i= 0 ,xx= 0.1 ,yy= (5.4321+1.122j) i= 1 ,xx= 0.2 ,yy= (5.9376+1.268j) i= 2 ,xx= 0.3 ,yy= (6.5321+1.438j) i= 3 ,xx= 0.4 ,yy= (7.2336+1.632j) i= 4 ,xx= 0.5 ,yy= (8.0625+1.85j) p=polyfit(xx,yy,4) for 5 points polyfit coefficients [ 1. +1.45680467e-12j 2. -1.77472663e-12j 3. +1.20000000e+00j 4. +1.10000000e+00j 5. +1.00000000e+00j] backwards? largest power first, expected about 5, 4, 3, 2, 1 polyval values of fit [ 5.4321+1.122j 5.9376+1.268j 6.5321+1.438j 7.2336+1.632j 8.0625+1.85j ] should be values above

No difference in polyfit or polyval, they accept real or complex.

similar files to above

test_polyfit.py test_polyfit_py.out test_cxpolyfit.py test_cxpolyfit_py.out

"C" language polyval on double

double polyval(int order, double p[], double x) // return y = p(x) { // using Horner's rule double y; int i; y = p[order]*x; // p one larger than order for(i=order-1; i>0; i--) y = (p[i]+y)*x; y = p[0]+y; return y; } // end polyval

"C" language polyval on complex

easy editing a*b becomes cxmul(a,b) a+b becomes cxadd(a,b) double becomes complex // #include "complex.h" compile with complex.c complex cxpolyval(int order, complex p[], complex x) // return y = p(x) { // using Horner's rule complex y; int i; y = cxmul(p[order],x); // p one larger than order for(i=order-1; i>0; i--) y = cxmul(cxadd(p[i],y),x); y = cxadd(p[0],y); return y; } // end cxpolyval

Java polyval on real

// polyval.java with simple test class polyval { polyval(){} double polyval(int order, double p[], double x) // return y = p(x) { // using Horner's rule double y; y = p[order]*x; // p one larger than order for(int i=order-1; i>0; i--) y = (p[i]+y)*x; y = p[0]+y; return y; } // end polyval public static void main (String[] args) { double y; double x = 2.0; double p[] = {1.0, 2.0, 3.0}; polyval c = new polyval(); y = c.polyval(2, p, x); System.out.println("polyval.java on real"); System.out.println("p = {"+p[0]+","+p[1]+","+p[2]+"}"); System.out.println("y=p(x) x="+x+", y="+y); } // end main } // end polyval.java

polyval.java real output

polyval.java on real p = {1.0,2.0,3.0} y=p(x) x=2.0, y=17.0

Java polyval on complex

easy editing a*b becomes a.multiply(b) a+b becomes a.add(b) double becomes Complex // cxpolyval.java with simple test on Complex // needs class Complex.class that defines type and functions class cxpolyval { cxpolyval(){} Complex cxpolyval(int order, Complex p[], Complex x) // return y = p(x) { // using Horner's rule Complex y; y = p[order].multiply(x); // p one larger than order for(int i=order-1; i>0; i--) y = (p[i].add(y)).multiply(x); y = p[0].add(y); return y; } // end cxpolyval public static void main (String[] args) { Complex y; Complex x = new Complex(2.0, 0.1); Complex p[] = {new Complex(1.0,0.01), new Complex(2.0,0.02), new Complex(3.0)}; cxpolyval c = new cxpolyval(); y = c.cxpolyval(2, p, x); System.out.println("cxpolyval.java on Complex"); System.out.println("p = {"+p[0]+","+p[1]+","+p[2]+"}"); System.out.println("y=p(x) x="+x+", y="+y); } // end main } // end cxpolyval.java

cxpolyval.java complex output

cxpolyval.java on Complex p = {(1.0,0.01),(2.0,0.02),(3.0,0.0)} y=p(x) x=(2.0,0.1), y=(16.968,1.4500000000000002) The notation (1.0,2.0) comes from Fortran where the syntax was 1.0 real, 2.0 imaginary.

many Java files built for Complex

Complex.java ComplexMatrix.java ComplexPolynomial.java

Even least square fit can be performed on complex values and the polynomial may have powers of all variables

* Y_actual X1 X2 * -------- ----- ----- * 32.5+i21.2 1.0+i1.0 2.5+i2.0 * 7.2+i5.3 2.0+i1.1 2.5 * 6.9+i1.1 3.0+i1.2 2.7 * 22.4+i0.4 2.2+i2.2 2.1+i2.1 * 10.4+i2.3 1.5+i3.4 2.0 * 11.3+i4.0 1.6+i1.0 2.0 * ... ... ... need at least as many as c's * * Find coefficients c0, c1, c2, ... such that * Y_approximate = c0 + c1*X1 + c2*X2 + c3*X1^2 + c4*X1*X2 + c5*X2^2 + * c6*X1^3 + c7*X1^2*X2 + c8*X1*X2^2 + c9*X2^3 + ... * and such that the sum of (Y_actual - Y_approximate) squared is minimized. * The same simultaneous equations are built, be sure to use a simultaneous * equation solver where all arithmetic was changed to complex. some examples with multiple variables to higher powers (Sorry, modifying to complex is left as an exercise to the student.) least_square_fit_2d.c least_square_fit_3d.c least_square_fit_4d.c least_square_fit_4d.java least_square_fit_4d_java.out
    <- previous    index    next ->

Other links

Go to top