DIRECTORY
Complex,
FS,
IO,
Polynomial,
RealFns,
Roots,
Rope;

RootsImpl: PROGRAM
IMPORTS Complex, FS, IO, Polynomial, RealFns
EXPORTS Roots =

BEGIN

Vec: TYPE = Complex.Vec;
RootArray: TYPE = Roots.RootArray;
LinRoots: TYPE = Roots.LinRoots;
CubeRoots: TYPE = Roots.CubeRoots;
QuarticRoots: TYPE = Roots.QuarticRoots;
testCount: NAT = 17;
shortTestCount: NAT = 5;
TestVector: TYPE = ARRAY [1..testCount] OF REAL;
ShortTestVector: TYPE = ARRAY [1..shortTestCount] OF REAL;
Inconsistent: PUBLIC SIGNAL = CODE;
root3over2: REAL _ RealFns.SqRt[3.0]/2.0;
coeffs: TestVector _ [-25.0, -17.0, -12.0, -5.6, -2.2, -1.0, -0.5, -0.001, 0.0, 0.001, 0.5, 1.0, 2.2, 5.6, 12.0, 17.0, 25.0];
shortCoeffs: ShortTestVector _ [-17.0, -1.0, 0.0, 1.0, 17.0];


ComplexLinearFormula: PUBLIC PROC [a, b: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
IF a = 0 THEN {
IF b = 0 THEN {roots[1] _ [0.0, 0.0]; rootCount _ 0; RETURN}
ELSE roots[1] _ [0.0, 0.0]; rootCount _ 0; RETURN}	-- inconsistent case
ELSE {
IF b = 0  THEN {roots[1] _ [0.0, 0.0]; rootCount _ 1; RETURN}
ELSE {roots[1] _ [-b/a, 0.0]; rootCount _ 1; RETURN};
};
};

ComplexQuadratic: PUBLIC PROC [a, b, c: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
IF a = 0 THEN {
linRoots: RootArray;
[linRoots, rootCount] _ ComplexLinearFormula[b, c];
RETURN[linRoots, rootCount];
};
IF c = 0 THEN {
linRoots: RootArray;
roots[1] _ [0.0, 0.0];
[linRoots, rootCount] _ ComplexLinearFormula[a, b];
roots[2] _ linRoots[1];
rootCount _ rootCount + 1;
RETURN};
RETURN QuadraticAux[a, b, c, (b*b)-(4.0*a*c)];
};	-- end of ComplexQuadraticFormula

QuadraticAux: PRIVATE PROC [a, b, c, discriminant: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
IF discriminant < 0 THEN {
twoA: REAL _ 2.0*a;
RETURN QuadraticConjugate[[-b/twoA, RealFns.SqRt[-discriminant]/twoA]];
};
IF discriminant = 0 THEN {
twoA: REAL _ 2.0*a;
roots[1] _ roots[2] _ [-b/twoA, 0.0];
rootCount _ 2;
RETURN;
};
IF b < 0 THEN RETURN QuadraticReal[a, RealFns.SqRt[discriminant] - b, c]
ELSE RETURN QuadraticReal[a, -RealFns.SqRt[discriminant] - b, c];
};	-- end of QuadraticAux

QuadraticReal: PRIVATE PROC [a, term, c: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
roots[1] _ [term/(2.0*a), 0.0];
roots[2] _ [(2.0*c)/term, 0.0];
rootCount _ 2;
RETURN;
};	-- end of QuadraticReal

QuadraticConjugate: PRIVATE PROC [z: Vec] RETURNS [roots: RootArray, rootCount: NAT] = {
roots[1] _ z;
roots[2] _ [z.x, -z.y];
rootCount _ 2;
};	-- end of QuadraticConjugate

ComplexCubic: PUBLIC PROC [a, b, c, d: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
bb, bbb, cubeTerm: REAL;
quadRoots: RootArray;
quadCount: NAT;
IF a = 0 THEN RETURN ComplexQuadratic[b, c, d];
bb _ b*b;
bbb _ bb*b;
cubeTerm _ (bb - 3.0*a*c);
cubeTerm _ cubeTerm*cubeTerm*cubeTerm;
[quadRoots, quadCount] _ ComplexQuadratic[1.0, 2.0*bbb + 9.0*a*(3.0*a*d-b*c), cubeTerm];
RETURN CubicAux[a,
					b,
					quadRoots,
					quadCount];
};	-- end of ComplexCubic

CubicAux: PRIVATE PROC [a, b: REAL, quadRoots: RootArray, quadCount: NAT] RETURNS [roots: RootArray, rootCount: NAT] = {
IF quadRoots[1].y = 0 THEN	-- quadRoots real
RETURN CubicConjugate[a, b, CubeRoot[quadRoots[1].x], CubeRoot[quadRoots[2].x]]
ELSE RETURN CubicReal[a, b, RealFns.SqRt[(quadRoots[1].x*quadRoots[1].x) +
													(quadRoots[1].y*quadRoots[1].y)],
									 RealFns.ArcTan[quadRoots[1].y, quadRoots[1].x]];
};	-- end of CubicAux

CubeRoot: PUBLIC PROCEDURE [z: REAL] RETURNS [root: REAL] = {
IF z < 0 THEN {
posZ: REAL _ -z;
root _ RealFns.Exp[RealFns.Ln[posZ]/3.0];
root _ -(2.0*root + posZ/(root*root))/3.0;
}
ELSE {
root _ RealFns.Exp[RealFns.Ln[z]/3.0];
root _ (2.0*root + z/(root*root))/3.0;
};
};

CubicConjugate: PRIVATE PROC [a, b, r, s: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
threeA: REAL _ 3.0*a;
RETURN CubicConjugateAux[(r+s-b)/threeA, (-(r+s)/2.0 - b)/threeA, (r-s)*root3over2/threeA];
};	-- end of CubicConjugate

CubicConjugateAux: PRIVATE PROC [realRoot, real, imaginary: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
roots[1] _ [realRoot, 0.0];
roots[2] _ [real, imaginary];
roots[3] _ [real, -imaginary];
rootCount _ 3;
};	-- end of CubicConjugateAux

CubicReal: PRIVATE PROC [a, b, rho, theta: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
RETURN CubicRealAux[a, b, 2.0*CubeRoot[rho], RealFns.Cos[theta/3.0]*-0.5,
RealFns.Sin[theta/3.0]*-root3over2];
};

CubicRealAux: PROCEDURE [a, b, rd, cd, sd: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
threeA: REAL _ 3.0*a;
roots[1] _ [((-2.0*rd*cd) - b)/threeA, 0.0];
roots[2] _ [(rd*(cd+sd) - b)/threeA, 0.0];
roots[3] _ [(rd*(cd-sd) - b)/threeA, 0.0];
rootCount _ 3;
}; 

ComplexQuartic: PROCEDURE [a, b, c, d, e: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
cubicCount: NAT;
cubicRoots: RootArray;
IF a=0 THEN RETURN ComplexCubic[b, c, d, e];
IF e=0 THEN {
[cubicRoots, cubicCount] _ ComplexCubic[a, b, c, d];
roots[1] _ [0.0, 0.0];
FOR i: NAT IN [1..cubicCount] DO
roots[i+1] _ cubicRoots[i];
ENDLOOP;
rootCount _ cubicCount+1;
RETURN;
};
[cubicRoots, cubicCount] _ ComplexCubic[1.0, -c, b*d-4.0*a*e, 4.0*a*c*e-(a*d*d+b*b*e)];
RETURN QuarticAux[a, b, c, d, e, cubicRoots[1].x];
};	-- end of ComplexQuartic

QuarticAux: PROCEDURE [a, b, c, d, e, s: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
RETURN QuarticSplit[a, b, RealFns.SqRt[ABS[b*b - 4.0*a*(c-s)]], s,
		    RealFns.SqRt[ABS[s*s-4.0*a*e]], b*s - 2.0*a*d];
};	-- end of QuarticAux

QuarticSplit: PROCEDURE [a, b, r1, s, r2, bsMinus2ad: REAL] RETURNS [roots: RootArray, rootCount: NAT] = {
quadRoots: RootArray;
quadCount: NAT;
IF (r1*r2*bsMinus2ad) < 0 THEN {
[roots, rootCount] _ ComplexQuadratic[2.0*a, b - r1, s+r2];
[quadRoots, quadCount] _ ComplexQuadratic[2.0*a, b+r1, s-r2];
roots[rootCount+1] _ quadRoots[1];
roots[rootCount+2] _ quadRoots[2];
rootCount _ rootCount + quadCount;
RETURN;
}
ELSE {
[roots, rootCount] _ ComplexQuadratic[2.0*a, b - r1, s-r2];
[quadRoots, quadCount] _ ComplexQuadratic[2.0*a, b+r1, s+r2];
roots[rootCount+1] _ quadRoots[1];
roots[rootCount+2] _ quadRoots[2];
rootCount _ rootCount + quadCount;
RETURN;
};
};	-- end of QuarticSplit

RealQuartic: PUBLIC PROC [a, b, c, d, e: REAL] RETURNS [roots: ARRAY[1..4] OF REAL, rootCount: NAT] = {
complexRoots: RootArray;
complexCount: NAT;
recipMaxMag: REAL;
recipMaxMag _ 1.0/MaxFive[a, b, c, d, e];
[complexRoots, complexCount] _ ComplexQuartic[a*recipMaxMag, b*recipMaxMag, c*recipMaxMag, d*recipMaxMag, e*recipMaxMag];
rootCount _ 0;
FOR i: NAT IN [1..complexCount] DO
IF complexRoots[i].y = 0.0 THEN {
rootCount _ rootCount + 1;
roots[rootCount] _ complexRoots[i].x;
};
ENDLOOP;
};	-- end of RealQuartic

MaxFive: PROCEDURE [a, b, c, d, e: REAL] RETURNS [max: REAL] = {
aMag, bMag, cMag, dMag, eMag: REAL;
aMag _ ABS[a]; bMag _ ABS[b]; cMag _ ABS[c]; dMag _ ABS[d]; eMag _ ABS[e];
max _ aMag;
IF bMag>max THEN max _ bMag;
IF cMag>max THEN max _ cMag;
IF dMag>max THEN max _ dMag;
IF eMag>max THEN max _ eMag;
};

EvalQuartic: PROCEDURE [r: Vec, a, b, c, d, e: REAL] RETURNS [result: Vec] = {
OPEN Complex;
t: Vec;
t _ ComplexAdd[ComplexMul[[a, 0.0], r], [b, 0.0]];	-- t _ ar+b;
t _ ComplexAdd[ComplexMul[t, r], [c, 0.0]];	-- t _ (ar+b)r+c;
t _ ComplexAdd[ComplexMul[t, r], [d, 0.0]];	-- t _ ((ar+b)r+c)r+d;
t _ ComplexAdd[ComplexMul[t, r], [e, 0.0]];	-- t _ (((ar+b)r+c)r+d)r+e;
result _ t;
};

EvalQuarticAtReal: PROCEDURE [r: REAL, a, b, c, d, e: REAL] RETURNS [result: REAL] = {
OPEN Complex;
t: REAL;
t _ (((a*r+b)*r+c)*r+d)*r+e;
result _ t;
};

EvalCubic: PROCEDURE [r: Vec, a, b, c, d: REAL] RETURNS [result: Vec] = {
OPEN Complex;
t: Vec;
t _ ComplexAdd[ComplexMul[[a, 0.0], r], [b, 0.0]];	-- t _ ar+b;
t _ ComplexAdd[ComplexMul[t, r], [c, 0.0]];	-- t _ (ar+b)r+c;
t _ ComplexAdd[ComplexMul[t, r], [d, 0.0]];	-- t _ ((ar+b)r+c)r+d;
result _ t;
};

EvalQuad: PROCEDURE [r: Vec, a, b, c: REAL] RETURNS [result: Vec] = {
OPEN Complex;
t: Vec;
t _ ComplexAdd[ComplexMul[[a, 0.0], r], [b, 0.0]];	-- t _ ar+b;
t _ ComplexAdd[ComplexMul[t, r], [c, 0.0]];	-- t _ (ar+b)r+c;
result _ t;
};

TestQuad: PUBLIC PROC [filename: Rope.ROPE] = {
file: IO.STREAM;
roots: RootArray;
rootCount: NAT;
eval: Vec;
file _ FS.StreamOpen[filename, $create];
FOR i: NAT IN [1..testCount] DO
FOR j: NAT IN [1..testCount] DO
[roots, rootCount] _ ComplexQuadratic[coeffs[i], coeffs[j], 1.0];
file.PutF["poly: [%g, %g, %g] ", [real[coeffs[i]]], [real[coeffs[j]]], [real[1.0]]];
file.PutF["roots: [(%g, %g), (%g, %g)]\n", [real[roots[1].x]], [real[roots[1].y]],
[real[roots[2].x]], [real[roots[2].y]]]; 
FOR k: NAT IN [1..rootCount] DO
eval _ EvalQuad[roots[k], coeffs[i], coeffs[j], 1.0];
file.PutF["		root %g evals to (%g, %g).\n", [integer[k]], [real[eval.x]], [real[eval.y]]];
ENDLOOP;
ENDLOOP;
ENDLOOP;
file.Close[];
};

TestCubic: PUBLIC PROC [filename: Rope.ROPE] = {
file: IO.STREAM;
roots: RootArray;
rootCount: NAT;
eval: Vec;
file _ FS.StreamOpen[filename, $create];
FOR i: NAT IN [1..shortTestCount] DO
FOR j: NAT IN [1..shortTestCount] DO
FOR l: NAT IN [1..shortTestCount] DO
[roots, rootCount] _ ComplexCubic[shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], 1.0];
file.PutF["poly: [%g, %g, %g, %g] ", [real[shortCoeffs[i]]], [real[shortCoeffs[j]]], [real[shortCoeffs[l]]], [real[1.0]]];
file.PutF["roots: ["];
file.PutF["[%g, %g]", [real[roots[1].x]], [real[roots[1].y]]];
FOR k: NAT IN [2..rootCount] DO
file.PutF[", [%g, %g]", [real[roots[k].x]], [real[roots[k].y]]];
ENDLOOP;
file.PutF["]\n"];
FOR k: NAT IN [1..rootCount] DO
eval _ EvalCubic[roots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], 1.0];
file.PutF["		root %g evals to (%g, %g).\n", [integer[k]], [real[eval.x]], [real[eval.y]]];
ENDLOOP;
ENDLOOP;
ENDLOOP;
ENDLOOP;
file.Close[];
};

TestQuartic: PUBLIC PROC [filename: Rope.ROPE] = {
file: IO.STREAM;
roots: RootArray;
rootCount: NAT;
eval: Vec;
file _ FS.StreamOpen[filename, $create];
FOR i: NAT IN [1..shortTestCount] DO
FOR j: NAT IN [1..shortTestCount] DO
FOR l: NAT IN [1..shortTestCount] DO
FOR m: NAT IN [1..shortTestCount] DO
[roots, rootCount] _ ComplexQuartic[shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["poly: [%g, %g, %g, %g, %g] ", [real[shortCoeffs[i]]], [real[shortCoeffs[j]]], [real[shortCoeffs[l]]], [real[shortCoeffs[m]]], [real[1.0]]];
file.PutF["roots: ["];
file.PutF["[%g, %g]", [real[roots[1].x]], [real[roots[1].y]]];
FOR k: NAT IN [2..rootCount] DO
file.PutF[", [%g, %g]", [real[roots[k].x]], [real[roots[k].y]]];
ENDLOOP;
file.PutF["]\n"];
FOR k: NAT IN [1..rootCount] DO
eval _ EvalQuartic[roots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["		root %g evals to (%g, %g).\n", [integer[k]], [real[eval.x]], [real[eval.y]]];
ENDLOOP;
ENDLOOP;
ENDLOOP;
ENDLOOP;
ENDLOOP;
file.Close[];
};

globalPoly: Polynomial.Ref;

CompareQuartic: PUBLIC PROC [filename: Rope.ROPE] = {
file: IO.STREAM;
roots: ARRAY[1..4] OF REAL;
rootCount: NAT;
realRoots: REF Polynomial.RealRootRec;
problem: BOOL;
eval: REAL;
file _ FS.StreamOpen[filename, $create];
FOR i: NAT IN [1..shortTestCount] DO
FOR j: NAT IN [1..shortTestCount] DO
FOR l: NAT IN [1..shortTestCount] DO
FOR m: NAT IN [1..shortTestCount] DO
[roots, rootCount] _ RealQuartic[shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
globalPoly[4] _ shortCoeffs[i];
globalPoly[3] _ shortCoeffs[j];
globalPoly[2] _ shortCoeffs[l];
globalPoly[1] _ shortCoeffs[m];
globalPoly[0] _ 1.0;
realRoots _ Polynomial.RealRoots[globalPoly];
IF realRoots.nRoots # rootCount THEN {
file.PutF["Different number of roots ** ** ** **\n"];
file.PutF["poly: [%g, %g, %g, %g, %g] ", [real[shortCoeffs[i]]], [real[shortCoeffs[j]]], [real[shortCoeffs[l]]], [real[shortCoeffs[m]]], [real[1.0]]];
file.PutF["Closed form:\n"];
file.PutF["roots: (%g)[", [integer[rootCount]] ];
IF rootCount > 0 THEN {
file.PutF["%g", [real[roots[1]]] ];
FOR k: NAT IN [2..rootCount] DO
file.PutF[", %g", [real[roots[k]]] ];
ENDLOOP;
};
file.PutF["]\n"];
FOR k: NAT IN [1..rootCount] DO
eval _ EvalQuarticAtReal[roots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["		root %g evals to %g.\n", [integer[k]], [real[eval]] ];
ENDLOOP;
file.PutF["Iterative form:\n"];
file.PutF["roots: (%g)[", [integer[realRoots.nRoots]] ];
IF realRoots.nRoots > 0 THEN {
file.PutF["%g", [real[realRoots[0]]] ];
FOR k: NAT IN [1..realRoots.nRoots-1] DO
file.PutF[", %g", [real[realRoots[k]]] ];
ENDLOOP;
};
file.PutF["]\n"];
FOR k: NAT IN [0..realRoots.nRoots) DO
eval _ EvalQuarticAtReal[realRoots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["		root %g evals to %g.\n", [integer[k]], [real[eval]] ];
ENDLOOP;
}
ELSE {
problem _ FALSE;
FOR k: NAT IN [1..rootCount] DO
eval _ EvalQuarticAtReal[roots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
IF eval > 1.0e-3 THEN problem _ TRUE;
ENDLOOP;
FOR k: NAT IN [0..realRoots.nRoots) DO
eval _ EvalQuarticAtReal[realRoots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
IF eval > 1.0e-3 THEN problem _ TRUE;
ENDLOOP;
IF problem THEN {
file.PutF["EVAL TOO LARGE\n"];
file.PutF["poly: [%g, %g, %g, %g, %g] ", [real[shortCoeffs[i]]], [real[shortCoeffs[j]]], [real[shortCoeffs[l]]], [real[shortCoeffs[m]]], [real[1.0]]];
file.PutF["Closed form:\n"];
file.PutF["roots: (%g)[", [integer[rootCount]] ];
file.PutF["%g", [real[roots[1]]] ];
FOR k: NAT IN [2..rootCount] DO
file.PutF[", %g", [real[roots[k]]] ];
ENDLOOP;
file.PutF["]\n"];
FOR k: NAT IN [1..rootCount] DO
eval _ EvalQuarticAtReal[roots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["		root %g evals to %g.\n", [integer[k]], [real[eval]] ];
ENDLOOP;
file.PutF["Iterative form:\n"];
file.PutF["roots: (%g)[", [integer[realRoots.nRoots]] ];
file.PutF["%g", [real[realRoots[0]]] ];
FOR k: NAT IN [1..realRoots.nRoots-1] DO
file.PutF[", %g", [real[realRoots[k]]] ];
ENDLOOP;
file.PutF["]\n"];
FOR k: NAT IN [0..realRoots.nRoots) DO
eval _ EvalQuarticAtReal[realRoots[k], shortCoeffs[i], shortCoeffs[j], shortCoeffs[l], shortCoeffs[m], 1.0];
file.PutF["		root %g evals to %g.\n", [integer[k]], [real[eval]] ];
ENDLOOP;
};
};
ENDLOOP;
ENDLOOP;
ENDLOOP;
ENDLOOP;
file.Close[];
};

ComplexMul: PROC [r, s: Vec] RETURNS [z: Vec] = {
z _ Complex.Mul[r, s];
};

ComplexAdd: PROC [r, s: Vec] RETURNS [z: Vec] = {
z _ Complex.Add[r, s];
};

Init: PROC = {
globalPoly _ Polynomial.Quartic[[0,0,0,0,0]];
};

Init[];

END.
	 File: RootsImpl.mesa
Copyright c 1984 by Xerox Corporation.  All rights reserved.
Last edited by Bier on December 31, 1984 12:58:55 pm PST
Contents: The quadratic formula, cubic formula and quartic formula as translated from LISP from 
"LISP" by Patrick Winston and Berthold Horn.  All functions find complex roots.

GLOBALS
The solution to the equation "ax + b = 0".
The solution to the equation "ax2+bx+c=0".  If a=0 this is just a linear equation.  If
c = 0, one root is zero and we solve a linear equation.  Otherwise, we use the
quadratic formula in either the form [-b+-(b2-4ac)1/2]/2a or (2c)/[-b-+(b2-4ac)1/2]
depending on ...
Finds the roots of equations of the form "ax^3+bx^2+cx+d".
r and s are cube roots
Evaluate the quartic polynomial ar4 + br3 + cr2 + dr + e using Horner's rule:
(((ar+b)r+c)r+d)r+e.
Evaluate the quartic polynomial ar4 + br3 + cr2 + dr + e using Horner's rule:
(((ar+b)r+c)r+d)r+e.
Evaluate the quartic polynomial ar3 + br2 + cr + d using Horner's rule:
((ar+b)r+c)r+d.
Evaluate the quartic polynomial ar2 + br + c using Horner's rule:
(ar+b)r+c.
Letting coefficients a, and b take on any of the values in the test vector, we use our quad finder to find all of the complex roots of each polynomial.  For each root, we evaluate the polynomial at that root.  The result should come out zero.  For now, write out a file entry for each root
Letting coefficients a, b, and c take on any of the values in the test vector, we use our quad finder to find all of the complex roots of each polynomial.  For each root, we evaluate the polynomial at that root.  The result should come out zero.  For now, write out a file entry for each root.
Letting coefficients a, b, and c take on any of the values in the test vector, we use our quad finder to find all of the complex roots of each polynomial.  For each root, we evaluate the polynomial at that root.  The result should come out zero.  For now, write out a file entry for each root.
Letting coefficients a, b, and c take on any of the values in the test vector, we use our two quad finders to find all of the real roots of each polynomial.  If the finders disagree on the number of roots, we print out the complete situation.  For each root, we evaluate the polynomial at that root.  If either finder is off by more than e-3, we print out the situation.
ΚΗ˜Ihead1š
œ™Jšœ
Οm
œ1™<Iprocš
œ8™8Lš
ϡ
™°
L™šœ
Οk	˜
L˜L˜L˜L˜Lš
œ˜Lš
œ˜Lš
œ˜—šœž˜Lšžœ%˜,Lšžœ˜—Lšœ
ž˜˜Lšœžœ˜Lšœžœ˜"Lšœ
žœ˜ Jšœžœ˜"Jšœžœ˜(Jšœžœ˜Jš
œ˜Lš	œžœžœžœ
žœ
˜0Lš	œžœžœžœ
žœ
˜:—
head2™Lšœžœ
žœžœ
˜#Lšœžœ˜)L˜}Lš
œ=˜=Lš
œ
˜
—šΟnœžœ
žœžœžœžœ˜]L™*šžœžœ˜Lšžœžœ(žœ
˜<Lšžœ'žœΟc˜G—šžœ˜Lšžœžœ(žœ
˜=Lšžœ)žœ˜5L˜—Lš
œ˜—šŸœžœ
žœžœžœžœ˜\Ilead1šœ Οu
ϱ
‘
œ‘œ‘
œ‘œ™Ššžœžœ˜Lš
œ˜Lš
œ3˜3Lšžœ˜L˜—šžœžœ˜L˜Lš
œ˜Lš
œ3˜3Lš
œ˜Lš
œ˜Lšžœ˜—Lšžœ(˜.Lšœ  œ˜%—šŸœžœ
žœžœžœžœ˜gšžœžœ˜Lšœžœ	˜LšžœA˜GL˜—šžœžœ˜Lšœžœ	˜Lš
œ%˜%Lš
œ˜Lšžœ
˜L˜—Lšžœžœ
žœ4˜HLšžœ
žœ6˜ALšœ œ
˜—šŸ
œžœ
žœžœžœžœ˜]Lš
œ˜L˜Lš
œ˜Lšžœ
˜Lšœ ˜—š
Ÿœžœ
žœ
žœžœ˜XLš
œ
˜
L˜Lš
œ˜Lšœ ˜ —šŸœžœ
žœžœžœžœ˜[L™:Lšœžœ
˜Lš
œ˜Lšœžœ
˜Lšžœžœ
žœ˜/L˜	L˜L˜L˜&Lš
œX˜Xšžœ˜Lš
œ(˜(—Lšœ ˜—šŸœžœ
žœžœ#žœžœžœ˜xšžœžœ
 ˜,LšžœI˜OLšžœ
žœ©
˜΄
—Lšœ ˜—šŸœžœ
ž	œžœžœžœ˜=šžœžœ˜Jšœžœ˜J˜)J˜*J˜
—šžœ˜J˜&J˜&J˜—J˜—šŸœžœ
žœžœžœžœ˜^J™Jšœžœ	˜JšžœU˜[Jšœ ˜—šŸœžœ
žœžœžœžœ˜pJš
œ˜Jš
œ˜Jš
œ˜Jš
œ˜Jšœ ˜—šŸ	œžœ
žœžœžœžœ˜_Jšžœh˜nJ˜—š
Ÿœž	œžœžœžœ˜_Jšœžœ	˜Jš
œ,˜,Jš
œ*˜*Jš
œ*˜*Jš
œ˜Jš
œ˜—š
Ÿœž	œžœžœžœ˜^Jšœžœ
˜Jš
œ˜Jšžœžœ
žœ˜,šžœžœ˜
J˜4Jš
œ˜šžœžœ
žœž˜ Jš
œ˜—Jšžœ
˜Jš
œ˜Jšžœ
˜J˜—J˜WJšžœ,˜2Jšœ ˜—š
Ÿ
œž	œžœžœžœ˜]Jšžœr˜xJšœ ˜—š
Ÿœž	œžœžœžœ˜jJš
œ˜Jšœžœ
˜šžœžœ˜ Jš
œ;˜;Jš
œ=˜=Jš
œ"˜"Jš
œ"˜"Jš
œ"˜"Jšžœ
˜J˜
—šžœ˜Jš
œ;˜;Jš
œ=˜=Jš
œ"˜"Jš
œ"˜"Jš
œ"˜"Jšžœ
˜J˜—Jšœ œ
˜—šŸœžœ
žœžœžœ	žœžœ
žœ
žœ˜gJš
œ˜Jšœžœ
˜Jšœ
žœ
˜Jš
œ)˜)Jš
œy˜yJš
œ˜šžœžœ
žœž˜"šžœžœ˜!Jš
œ˜Jš
œ%˜%J˜——Jšžœ
˜Jšœ ˜—š
Ÿœž	œžœžœžœ˜@Jšœžœ
˜#J˜JJš
œ˜Jšžœ
žœ˜Jšžœ
žœ˜Jšžœ
žœ˜Jšžœ
žœ˜Jš
œ˜—šŸœž	œžœžœ˜NJšœ"‘
œ‘
œ‘
œ3™bJšžœ	˜
Jš
œ˜Jšœ3 ˜?Jšœ, ˜=Jšœ, ˜BJšœ, ˜GJš
œ˜Jš
œ˜—šŸœž	œžœžœžœ
žœ˜VJšœ"‘
œ‘
œ‘
œ3™bJšžœ	˜
Jšœžœ
˜Jš
œ˜Jš
œ˜Jš
œ˜—šŸ	œž	œžœžœ˜IJšœ"‘
œ‘
œ.™WJšžœ	˜
Jš
œ˜Jšœ3 ˜?Jšœ, ˜=Jšœ, ˜BJš
œ˜Jš
œ˜—šŸœž	œžœžœ˜EJšœ"‘
œ)™LJšžœ	˜
Jš
œ˜Jšœ3 ˜?Jšœ, ˜=Jš
œ˜Jš
œ˜—šŸœžœ
žœžœ˜/Jšœžœ
žœ
˜J˜Jšœžœ
˜Jš
œ
˜
Jš
œ(˜(J™‘šžœžœ
žœž˜šžœžœ
žœž˜J˜AJš
œT˜TJš
œ|˜|šžœžœ
žœž˜Jš
œ5˜5Jš
œZ˜Z—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜Jš
œ
˜
J˜—šŸ	œžœ
žœžœ˜0Jšœžœ
žœ
˜J˜Jšœžœ
˜Jš
œ
˜
Jš
œ(˜(J™₯šžœžœ
žœž˜$šžœžœ
žœž˜$šžœžœ
žœž˜$J˜WJš
œz˜zJš
œ˜Jš
œ>˜>šžœžœ
žœž˜Jš
œ@˜@—Jšžœ
˜Jš
œ˜šžœžœ
žœž˜Jš
œP˜PJš
œZ˜Z—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜Jš
œ
˜
J˜—šŸœžœ
žœžœ˜2Jšœžœ
žœ
˜J˜Jšœžœ
˜Jš
œ
˜
Jš
œ(˜(J™₯šžœžœ
žœž˜$šžœžœ
žœž˜$šžœžœ
žœž˜$šžœžœ
žœž˜$J˜iJš
œ–
˜–
Jš
œ˜Jš
œ>˜>šžœžœ
žœž˜Jš
œ@˜@—Jšžœ
˜Jš
œ˜šžœžœ
žœž˜Jš
œb˜bJš
œZ˜Z—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜Jš
œ
˜
J˜—Jš
œ˜šŸœžœ
žœžœ˜5Jšœžœ
žœ
˜Jšœžœžœ
žœ
˜Jšœžœ
˜Jšœžœ˜&Jšœ	žœ
˜Jš
œ˜Jš
œ(˜(J™ςšžœžœ
žœž˜$šžœžœ
žœž˜$šžœžœ
žœž˜$šžœžœ
žœž˜$JšœΟbœF˜fJš
œ˜Jš
œ˜Jš
œ˜Jš
œ˜Jš
œ˜Jšœ’œ
˜-šžœžœ˜&Jš
œ5˜5Jš
œ–
˜–
Jš
œ˜Jš
œ1˜1š
œ˜Jš
œ#˜#šžœžœ
žœž˜Jš
œ%˜%—Jšžœ
˜J˜—Jš
œ˜šžœžœ
žœž˜Jš
œh˜hJš
œC˜C—Jšžœ
˜Jš
œ˜Jš
œ8˜8š
œ˜Jš
œ'˜'šžœžœ
žœž˜(Jš
œ)˜)—Jšžœ
˜J˜—Jš
œ˜šžœžœ
žœž˜&Jš
œl˜lJš
œC˜C—Jšžœ
˜J˜
—šžœ˜Jšœ
žœ
˜šžœžœ
žœž˜Jš
œh˜hJšžœžœžœ
˜%—Jšžœ
˜šžœžœ
žœž˜&Jš
œl˜lJšžœžœžœ
˜%—Jšžœ
˜šžœ	žœ˜Jš
œ˜Jš
œ–
˜–
Jš
œ˜Jš
œ1˜1Jš
œ#˜#šžœžœ
žœž˜Jš
œ%˜%—Jšžœ
˜Jš
œ˜šžœžœ
žœž˜Jš
œh˜hJš
œC˜C—Jšžœ
˜Jš
œ˜Jš
œ8˜8Jš
œ'˜'šžœžœ
žœž˜(Jš
œ)˜)—Jšžœ
˜Jš
œ˜šžœžœ
žœž˜&Jš
œl˜lJš
œC˜C—Jšžœ
˜J˜—J˜——Jšžœ
˜—Jšžœ
˜—Jšžœ
˜—Jšžœ
˜Jš
œ
˜
J˜—šŸ
œžœ
žœ
˜1J˜J˜—šŸ
œžœ
žœ
˜1J˜J˜—šŸœžœ˜Jš
œ-˜-Jš
œ˜—š
œ˜J˜—Jšžœ
˜—…—8<R#