ComplexesImpl.mesa
Last Edited by: Arnon, July 19, 1985 2:54:46 pm PDT
DIRECTORY
Rope,
IO,
Basics,
Ieee,
Real,
RealFns,
Vector,
Vector2,
Atom,
Convert,
AlgebraClasses,
MathConstructors,
Bools,
Ints,
Reals,
Complexes;
ComplexesImpl: CEDAR PROGRAM
IMPORTS IO, Convert, Rope, Real, RealFns, AlgebraClasses, MathConstructors, Ints, Reals
EXPORTS Complexes
= BEGIN OPEN Complexes, Convert, AC: AlgebraClasses;
Types and Constants
ComplexesError: PUBLIC SIGNAL [reason: ATOM ← $Unspecified] = CODE;
bitsPerWord: CARDINAL = Basics.bitsPerWord;
CARD: TYPE = LONG CARDINAL;
ROPE: TYPE = Rope.ROPE;
STREAM: TYPE = IO.STREAM;
Object: TYPE = AC.Object;
Method: TYPE = AC.Method;
Structure Operations
PrintName: PUBLIC AC.ToRopeOp = {
RETURN["Complexes"];
};
ShortPrintName: PUBLIC AC.ToRopeOp = {
RETURN["C"];
};
Characteristic: PUBLIC AC.StructureRankOp = {
RETURN[ 0 ]
};
I/O and Conversion
Recast: PUBLIC AC.BinaryOp = {
args are a StructureElement and a Structure
IF AC.StructureEqual[firstArg.class, Complexes] THEN RETURN[firstArg];
IF Reals.CanRecast[firstArg, Reals.Reals] THEN {
real: REAL ← Reals.ToREAL[Reals.Recast[firstArg, Reals.Reals] ];
RETURN[FromPairREAL[real, 0.0] ];
};
RETURN[NIL];
};
CanRecast: PUBLIC AC.BinaryPredicate = {
firstArgStructure: Object ← IF firstArg.flavor = StructureElement THEN firstArg.class ELSE IF firstArg.flavor = Structure THEN firstArg ELSE ERROR;
SELECT TRUE FROM
AC.StructureEqual[firstArgStructure, Complexes] => RETURN[TRUE];
Reals.CanRecast[firstArg, Reals.Reals] => RETURN[TRUE];
ENDCASE;
RETURN[FALSE];
};
LegalFirstChar: PUBLIC AC.LegalFirstCharOp = {
SELECT char FROM
= '( => RETURN[TRUE];
ENDCASE;
RETURN[Reals.LegalFirstChar[char, structure] ];
};
Read: PUBLIC AC.ReadOp ~ {
puncChar, imagChar: CHAR;
var: Rope.ROPE;
real, imag: REAL;
negative: BOOL ← FALSE;
ReadVLFail: PUBLIC ERROR [subclass: ATOM ← $Unspecified] = CODE;
{
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[ ];
IF puncChar # '( THEN GO TO Error;
[]← in.SkipWhitespace[];
real ← in.GetReal[];
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[];
IF puncChar = '- THEN negative ← TRUE ELSE IF puncChar # '+ THEN GO TO Error;
[]← in.SkipWhitespace[];
imagChar ← in.PeekChar[];
IF imagChar # '. AND NOT imagChar IN ['0..'9] THEN imag ← 1.0 ELSE imag ← in.GetReal[];
IF negative THEN imag ← - imag;
var ← in.GetID[];
IF NOT Rope.Equal[var,"i"] AND NOT Rope.Equal[var,"I"] THEN GO TO Error;
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[];
IF puncChar # ') THEN GO TO Error;
RETURN[ FromPairREAL[real, imag] ];
EXITS
Error => {
in.SetIndex[0];
RETURN[Recast[ Reals.Read[in, Reals.Reals], Complexes] ] };
};
};
FromRope: PUBLIC AC.FromRopeOp = {
stream: IO.STREAM ← IO.RIS[in];
RETURN[ Read[stream, structure] ];
};
ToRope: PUBLIC AC.ToRopeOp = {
data: ComplexData ← NARROW[in.data];
out ← "(";
out ← Rope.Concat[ out, RopeFromReal[data.x] ];
IF data.y < 0 THEN out ← Rope.Concat[ out, " - " ] ELSE out ← Rope.Concat[ out, " + " ];
out ← Rope.Concat[ out, RopeFromReal[ABS[data.y]] ];
out ← Rope.Concat[ out, " i )" ];
};
Write: PUBLIC AC.WriteOp = {
IO.PutRope[ stream, ToRope[in] ]
};
ToExpr: PUBLIC AC.ToExprOp = {
data: ComplexData ← NARROW[in.data];
RETURN[MathConstructors.MakeComplex[
MathConstructors.MakeReal[data.x],
MathConstructors.MakeReal[data.y]] ]
};
FromPairReal: PUBLIC AC.BinaryOp = {
RETURN[ FromPairREAL[Reals.ToREAL[firstArg], Reals.ToREAL[secondArg] ] ];
};
RealArgsDesired: PUBLIC AC.UnaryToListOp ~ {
RETURN[ LIST[Reals.Reals] ];
};
FromPairREAL: PUBLIC PROC [realPart, imagPart: REAL] RETURNS [out: Complex] = {
RETURN[ NEW[AC.ObjectRec ← [class: Complexes, flavor: StructureElement, data: NEW[Vector2.VEC ← [realPart, imagPart] ] ] ] ];
};
ToPairREAL: PUBLIC PROC [in: Complex] RETURNS [realPart, imagPart: REAL] = {
data: ComplexData ← NARROW[in.data];
RETURN[realPart: data.x, imagPart: data.y];
};
Comparison
Equal: PUBLIC AC.BinaryPredicate = {
RETURN[ AlmostEqual[firstArg, secondArg] ]
};
Arithmetic
Zero: PUBLIC AC.NullaryOp = {
RETURN[ ComplexZero ]
};
One: PUBLIC AC.NullaryOp = {
RETURN[ ComplexOne ]
};
Add: PUBLIC AC.BinaryOp ~ {
firstData: ComplexData ← NARROW[firstArg.data];
secondData: ComplexData ← NARROW[secondArg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [firstData.x + secondData.x, firstData.y + secondData.y] ]
] ] ];
};
Negate: PUBLIC AC.UnaryOp ~ {
data: ComplexData ← NARROW[arg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [-data.x, -data.y] ]
] ] ];
};
Subtract: PUBLIC AC.BinaryOp ~ {
firstData: ComplexData ← NARROW[firstArg.data];
secondData: ComplexData ← NARROW[secondArg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [firstData.x - secondData.x, firstData.y - secondData.y] ]
] ] ];
};
Multiply: PUBLIC AC.BinaryOp ~ {
firstData: ComplexData ← NARROW[firstArg.data];
secondData: ComplexData ← NARROW[secondArg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [
(firstData.x*secondData.x - firstData.y*secondData.y),
(firstData.x*secondData.y + firstData.y*secondData.x)
] ]
] ] ];
};
Conjugate: PUBLIC AC.UnaryOp ~ {
data: ComplexData ← NARROW[arg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [data.x, -data.y] ]
] ] ];
};
Modulus: PUBLIC AC.UnaryOp = {
data: ComplexData ← NARROW[arg.data];
RETURN[Reals.FromREAL[Real.SqRt[data.x*data.x+data.y*data.y] ] ]
}; -- same as Vector.Mag
ModulusSquared: PUBLIC PROC [a: Complex] RETURNS [REAL] ~ {
data: ComplexData ← NARROW[a.data];
RETURN[data.x*data.x + data.y*data.y];
};
ObjectAndIntDesired: PUBLIC AC.UnaryToListOp ~ {
RETURN[ LIST[arg, Ints.Ints] ]; -- arg assumed to be a Structure
};
Power: PUBLIC AC.BinaryOp ~ { -- this simple algorithm is Structure independent
power: INT ← Ints.ToINT[secondArg];
structure: Object ← firstArg.class;
one: Object ← AC.ApplyLkpNoRecastObject[$one, structure, LIST[structure] ];
productMethod: Method ← AC.LookupMethodInStructure[$product, structure];
IF power < 0 THEN {
invertMethod: Method ← AC.LookupMethodInStructure[$invert, structure];
temp: Object;
IF invertMethod = NIL THEN ERROR;
temp ← Power[firstArg, Ints.FromINT[ABS[power] ] ];
RETURN[AC.ApplyNoLkpNoRecastObject[invertMethod, LIST[temp] ] ];
};
IF power = 0 THEN RETURN[one];
result ← firstArg;
FOR i:INT IN [2..power] DO
result ← AC.ApplyNoLkpNoRecastObject[productMethod, LIST[firstArg, result] ];
ENDLOOP;
};
Invert: PUBLIC AC.UnaryOp ~ {
m: REAL ← ModulusSquared[arg];
data: ComplexData ← NARROW[arg.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [data.x / m , -data.y / m ] ]
] ] ];
};
Divide: PUBLIC AC.BinaryOp ~ {
RETURN[ Multiply[firstArg, Invert[secondArg]] ];
};
Paren: PUBLIC AC.UnaryOp ~ {
RETURN[NEW[AC.ObjectRec ← [
flavor: StructureElement,
class: Complexes,
data: arg.data
] ] ];
};
Exponent: PROCEDURE [x: REAL] RETURNS [INTEGER] = INLINE
BEGIN
fl: Ieee.SingleReal ← LOOPHOLE[x];
RETURN[fl.exp];
END;
AlmostEqual: PUBLIC PROCEDURE [a: Complex, b: Complex, mag:[-126..0] ← -20]
RETURNS [BOOLEAN] =
BEGIN
sumSqrAbs: REAL ← ModulusSquared[a] + ModulusSquared[b];
sqrAbsDif: REAL ← ModulusSquared[Subtract[a, b]];
RETURN [Exponent[sumSqrAbs]+mag+mag-1>Exponent[sqrAbsDif]];
END;
FromPolar: PUBLIC PROCEDURE [r: REAL, radians: REAL] RETURNS [Complex] ~ {
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [ r*RealFns.Cos[radians] , r*RealFns.Sin[radians] ] ]
] ] ];
};
Arg: PUBLIC PROCEDURE [a: Complex] RETURNS [REAL] = {
data: ComplexData ← NARROW[a.data];
RETURN[RealFns.ArcTan[data.y, data.x]]
};
Exp: PUBLIC PROCEDURE [a: Complex] RETURNS [Complex] = {
data: ComplexData ← NARROW[a.data];
RETURN[FromPolar[RealFns.Exp[data.x], data.y]]};
complex exponential function
Ln: PUBLIC PROCEDURE [a: Complex] RETURNS [Complex] = {
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [ RealFns.Ln[Reals.ToREAL[Modulus[a] ] ], Arg[a] ] ]
] ] ];
};
Sqr: PUBLIC PROCEDURE [a: Complex] RETURNS [Complex] = {
data: ComplexData ← NARROW[a.data];
RETURN[NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [(data.x*data.x - data.y*data.y) , (2*data.x*data.y) ] ]
] ] ];
};
SqRt: PUBLIC PROCEDURE [a: Complex] RETURNS [Complex] = {
Find the complex square root, using Newton's iteration
data: ComplexData ← NARROW[a.data];
z: Complex ← (IF data.y>=0 THEN FromPairREAL[1, 1] ELSE FromPairREAL[1, -1]);
zData: ComplexData ← NARROW[z.data];
oldz:Complex;
IF data.y = 0 THEN
BEGIN
IF data.x>0 THEN zData.y𡤀 -- real square root
ELSE IF data.x<0 THEN zData.x𡤀 -- pure imaginary
ELSE RETURN[FromPairREAL[0, 0]]
END;
FOR I:NAT IN [0..50) DO
oldz ← z;
z ← Add[z, Divide[a, z]];
zData ← NARROW[z.data];
z ← NEW[AC.ObjectRec ← [
class: Complexes,
flavor: StructureElement,
data: NEW[Vector2.VEC ← [zData.x / 0.5 , zData.y / 0.5 ] ]
] ];
IF AlmostEqual[z, oldz, -20] THEN EXIT;
ENDLOOP;
RETURN[z];
};
Standard Desired Arg Structures
ComplexesDesired: PUBLIC AC.UnaryToListOp ~ {
Name Complexes explicitly, instead of using AC.DefaultDesiredArgStructures, so that if a Complexes method found by lookup from a subclasss, then will recast its args correctly (i.e. to Complexes)
RETURN[ LIST[Complexes] ];
};
Start Code
ComplexClass: AC.StructureClass ← NEW[AC.StructureClassRec ← [
characteristic: ClassCharacteristic,
isElementOf: AC.defaultElementOfProc,
completeField: TRUE,
realField: FALSE,
realClosedField: FALSE,
algebraicallyClosedField: TRUE,
] ];
ComplexClass: Object ← AC.MakeClass["ComplexClass", NIL, NIL];
Complexes: PUBLIC Object ← AC.MakeStructure["Complexes", ComplexClass, NIL];
ComplexOne: Complex ← FromPairREAL[1.0, 0.0]; -- do after Complexes set
ComplexZero: Complex ← FromPairREAL[0.0, 0.0];
categoryMethod: Method ← AC.MakeMethod[Value, FALSE, NEW[AC.Category ← field], NIL, "category"];
groundStructureMethod: Method ← AC.MakeMethod[Value, FALSE, NIL, NIL, "groundStructure"];
shortPrintNameMethod: Method ← AC.MakeMethod[ToRopeOp, FALSE, NEW[AC.ToRopeOp ← ShortPrintName], NIL, "shortPrintName"];
recastMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Recast], NIL, "recast"];
canRecastMethod: Method ← AC.MakeMethod[BinaryPredicate, TRUE, NEW[AC.BinaryPredicate ← CanRecast], NIL, "canRecast"];
legalFirstCharMethod: Method ← AC.MakeMethod[LegalFirstCharOp, FALSE, NEW[AC.LegalFirstCharOp ← LegalFirstChar], NIL, "legalFirstChar"];
readMethod: Method ← AC.MakeMethod[ReadOp, FALSE, NEW[AC.ReadOp ← Read], NIL, "read"];
fromRopeMethod: Method ← AC.MakeMethod[FromRopeOp, TRUE, NEW[AC.FromRopeOp ← FromRope], NIL, "fromRope"];
toRopeMethod: Method ← AC.MakeMethod[ToRopeOp, FALSE, NEW[AC.ToRopeOp ← ToRope], NIL, "toRope"];
toExprMethod: Method ← AC.MakeMethod[ToExprOp, FALSE, NEW[AC.ToExprOp ← ToExpr], NEW[AC.UnaryToListOp ← ComplexesDesired], "toExpr"];
complexMethod: Method ← AC.MakeMethod[BinaryOp, FALSE, NEW[AC.BinaryOp ← FromPairReal], NEW[AC.UnaryToListOp ← RealArgsDesired], "complex"];
zeroMethod: Method ← AC.MakeMethod[NullaryOp, FALSE, NEW[AC.NullaryOp ← Zero], NIL, "zero"];
oneMethod: Method ← AC.MakeMethod[NullaryOp, FALSE, NEW[AC.NullaryOp ← One], NIL, "one"];
parenMethod: Method ← AC.MakeMethod[UnaryOp, FALSE, NEW[AC.UnaryOp ← Paren], NIL, "paren"];
sumMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Add], NEW[AC.UnaryToListOp ← ComplexesDesired], "sum"];
negationMethod: Method ← AC.MakeMethod[UnaryOp, TRUE, NEW[AC.UnaryOp ← Negate], NEW[AC.UnaryToListOp ← ComplexesDesired], "negation"];
differenceMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Subtract], NEW[AC.UnaryToListOp ← ComplexesDesired], "difference"];
productMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Multiply], NEW[AC.UnaryToListOp ← ComplexesDesired], "product"];
commutativeMethod: Method ← AC.MakeMethod[Value, FALSE, NIL, NIL, "commutative"];
powerMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Power], NEW[AC.UnaryToListOp ← ObjectAndIntDesired], "power"];
conjugateMethod: Method ← AC.MakeMethod[UnaryOp, TRUE, NEW[AC.UnaryOp ← Conjugate], NEW[AC.UnaryToListOp ← ComplexesDesired], "conjugate"];
modulusMethod: Method ← AC.MakeMethod[UnaryOp, TRUE, NEW[AC.UnaryOp ← Modulus], NEW[AC.UnaryToListOp ← ComplexesDesired], "modulus"];
invertMethod: Method ← AC.MakeMethod[UnaryOp, TRUE, NEW[AC.UnaryOp ← Invert], NEW[AC.UnaryToListOp ← ComplexesDesired], "invert"];
fractionMethod: Method ← AC.MakeMethod[BinaryOp, TRUE, NEW[AC.BinaryOp ← Divide], NEW[AC.UnaryToListOp ← ComplexesDesired], "fraction"];
equalMethod: Method ← AC.MakeMethod[BinaryPredicate, TRUE, NEW[AC.BinaryPredicate ← Equal], NEW[AC.UnaryToListOp ← ComplexesDesired], "equals"];
AC.AddMethodToClass[$category, categoryMethod, ComplexClass];
AC.AddMethodToClass[$groundStructure, categoryMethod, ComplexClass];
AC.AddMethodToClass[$shortPrintName, shortPrintNameMethod, ComplexClass];
AC.AddMethodToClass[$recast, recastMethod, ComplexClass];
AC.AddMethodToClass[$canRecast, canRecastMethod, ComplexClass];
AC.AddMethodToClass[$legalFirstChar, legalFirstCharMethod, ComplexClass];
AC.AddMethodToClass[$read, readMethod, ComplexClass];
AC.AddMethodToClass[$fromRope, fromRopeMethod, ComplexClass];
AC.AddMethodToClass[$toRope, toRopeMethod, ComplexClass];
AC.AddMethodToClass[$toExpr, toExprMethod, ComplexClass];
AC.AddMethodToClass[$complex, complexMethod, ComplexClass];
AC.AddMethodToClass[$zero, zeroMethod, ComplexClass];
AC.AddMethodToClass[$one, oneMethod, ComplexClass];
AC.AddMethodToClass[$paren, parenMethod, ComplexClass];
AC.AddMethodToClass[$sum, sumMethod, ComplexClass];
AC.AddMethodToClass[$negation, negationMethod, ComplexClass];
AC.AddMethodToClass[$difference, differenceMethod, ComplexClass];
AC.AddMethodToClass[$product, productMethod, ComplexClass];
AC.AddMethodToClass[$commutative, commutativeMethod, ComplexClass];
AC.AddMethodToClass[$pow, powerMethod, ComplexClass];
AC.AddMethodToClass[$conjugate, conjugateMethod, ComplexClass];
AC.AddMethodToClass[$modulus, modulusMethod, ComplexClass];
AC.AddMethodToClass[$invert, invertMethod, ComplexClass];
AC.AddMethodToClass[$fraction, fractionMethod, ComplexClass];
AC.AddMethodToClass[$eqFormula, equalMethod, ComplexClass];
AC.InstallStructure[Complexes];
AC.SetSuperClass[Reals.Reals, Complexes];
END.