DIRECTORY Complex, FileIO, IO, Polynomial, RealFns, Roots, Rope; RootsImpl: PROGRAM IMPORTS Complex, FileIO, 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; }; 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 Abs: PROCEDURE [z: REAL] RETURNS [REAL] = { RETURN[IF z < 0 THEN -z ELSE z]; }; 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 _ FileIO.Open[filename, overwrite]; 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 _ FileIO.Open[filename, overwrite]; 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 _ FileIO.Open[filename, overwrite]; 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 _ FileIO.Open[filename, overwrite]; 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 Author: Eric Bier on June 10, 1983 6:56 pm Last edited by Bier on June 15, 1983 3:49 pm 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šΟc™Iprocšœ*™*Lšœ,™,Lšœ°™°L™šœΟk ˜ L˜L˜L˜L˜ Lšœ˜Lšœ˜Lšœ˜—šœ ž˜Lšžœ)˜0Lšžœ˜—Lšœž˜˜Lšœžœ˜Lšœ žœ˜"Lšœ žœ˜ Jšœ žœ˜"Jšœžœ˜(Jšœ žœ˜Jšœ˜Lš œ žœžœžœžœ˜0Lš œžœžœžœžœ˜:—head2™Lšœžœžœžœ˜#Lšœ žœ˜)L˜}Lšœ=˜=Lšœ˜—š Οnœžœžœžœžœžœ˜]L™*šžœžœ˜Lšžœžœ(žœ˜˜>šžœžœžœž˜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˜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ŒRl