Annotation of OpenXM_contrib/PHC/Ada/Math_Lib/Numbers/generic_complex_numbers.adb, Revision 1.1
1.1 ! maekawa 1: package body Generic_Complex_Numbers is
! 2:
! 3: TWO : constant number := one+one;
! 4:
! 5: -- CREATORS :
! 6:
! 7: function Create ( i : integer ) return Complex_Number is
! 8:
! 9: res : Complex_Number;
! 10:
! 11: begin
! 12: if i = 0
! 13: then res.re := +zero;
! 14: elsif i = 1
! 15: then res.re := +one;
! 16: else res.re := Create(i);
! 17: end if;
! 18: res.im := +zero;
! 19: return res;
! 20: end Create;
! 21:
! 22: function Create ( f : number ) return Complex_Number is
! 23:
! 24: res : Complex_Number;
! 25:
! 26: begin
! 27: res.re := +f;
! 28: res.im := +zero;
! 29: return res;
! 30: end Create;
! 31:
! 32: function Create ( re,im : number ) return Complex_Number is
! 33:
! 34: res : Complex_Number;
! 35:
! 36: begin
! 37: res.RE := +re;
! 38: res.IM := +im;
! 39: return res;
! 40: end Create;
! 41:
! 42: function Conjugate ( c : Complex_Number ) return Complex_Number is
! 43:
! 44: res : Complex_Number;
! 45:
! 46: begin
! 47: res.RE := +c.RE;
! 48: res.IM := -c.IM;
! 49: return res;
! 50: end Conjugate;
! 51:
! 52: -- COMPARISON and COPYING :
! 53:
! 54: function Equal ( x,y : Complex_Number ) return boolean is
! 55: begin
! 56: return (Equal(x.RE,y.RE) and Equal(x.IM,y.IM));
! 57: end Equal;
! 58:
! 59: procedure Copy ( x : in Complex_Number; y : in out Complex_Number ) is
! 60: begin
! 61: Copy(x.RE,y.RE);
! 62: Copy(x.IM,y.IM);
! 63: end Copy;
! 64:
! 65: function "<" ( x,y : Complex_Number ) return boolean is
! 66:
! 67: avx : number := AbsVal(x);
! 68: avy : number := AbsVal(y);
! 69: res : boolean := (avx < avy);
! 70:
! 71: begin
! 72: Clear(avx);
! 73: Clear(avy);
! 74: return res;
! 75: end "<";
! 76:
! 77: function ">" ( x,y : Complex_Number ) return boolean is
! 78:
! 79: avx : number := AbsVal(x);
! 80: avy : number := AbsVal(y);
! 81: res : boolean := (avx > avy);
! 82:
! 83: begin
! 84: Clear(avx);
! 85: Clear(avy);
! 86: return res;
! 87: end ">";
! 88:
! 89: -- SELECTORS :
! 90:
! 91: function REAL_PART ( x : Complex_Number ) return number is
! 92: begin
! 93: return +x.RE;
! 94: end REAL_PART;
! 95:
! 96: function IMAG_PART ( x : Complex_Number ) return number is
! 97: begin
! 98: return +x.IM;
! 99: end IMAG_PART;
! 100:
! 101: function AbsVal ( x : Complex_Number ) return number is
! 102:
! 103: res : number := AbsVal(x.RE);
! 104: acc : number := AbsVal(x.IM);
! 105:
! 106: begin
! 107: Add(res,acc);
! 108: Clear(acc);
! 109: return res;
! 110: end AbsVal;
! 111:
! 112: function AbsVal ( x : Complex_Number ) return Complex_Number is
! 113:
! 114: abx : number := AbsVal(x);
! 115: res : Complex_Number := Create(abx);
! 116:
! 117: begin
! 118: Clear(abx);
! 119: return res;
! 120: end AbsVal;
! 121:
! 122: -- ARITHMETIC OPERATIONS AS FUNCTIONS :
! 123:
! 124: function "+" ( x : Complex_Number; y : number ) return Complex_Number is
! 125: begin
! 126: return (x.RE+y,+x.IM);
! 127: end "+";
! 128:
! 129: function "-" ( x : Complex_Number; y : number ) return Complex_Number is
! 130: begin
! 131: return (x.RE-y,+x.IM);
! 132: end "-";
! 133:
! 134: function "*" ( x : Complex_Number; y : number ) return Complex_Number is
! 135: begin
! 136: return (x.RE*y,x.IM*y);
! 137: end "*";
! 138:
! 139: function "/" ( x : Complex_Number; y : number ) return Complex_Number is
! 140: begin
! 141: return (x.RE/y,x.IM/y);
! 142: end "/";
! 143:
! 144: function "+" ( x : number; y : Complex_Number ) return Complex_Number is
! 145: begin
! 146: return (x+y.RE,+y.IM);
! 147: end "+";
! 148:
! 149: function "-" ( x : number; y : Complex_Number ) return Complex_Number is
! 150: begin
! 151: return (x-y.RE,-y.IM);
! 152: end "-";
! 153:
! 154: function "*" ( x : number; y : Complex_Number ) return Complex_Number is
! 155: begin
! 156: return (x*y.RE,x*y.IM);
! 157: end "*";
! 158:
! 159: function "/" ( x : number; y : Complex_Number ) return Complex_Number is
! 160:
! 161: res : Complex_Number;
! 162: acc,avyim,avyre : Number;
! 163:
! 164: begin
! 165: if Equal(y.IM,zero)
! 166: then res.RE := x/y.RE;
! 167: res.IM := +zero;
! 168: elsif Equal(y.RE,zero)
! 169: then res.RE := +zero;
! 170: res.IM := x/y.IM; Min(res.IM);
! 171: else avyim := AbsVal(y.IM); avyre := AbsVal(y.RE);
! 172: if avyim < avyre
! 173: then acc := y.IM/y.RE;
! 174: res.RE := x*acc; Mul(acc,y.IM); Add(acc,y.RE); Div(res.RE,acc);
! 175: res.IM := x/acc; Min(res.IM);
! 176: Clear(acc);
! 177: elsif avyim > avyre
! 178: then acc := y.RE/y.IM;
! 179: res.IM := x*acc; Min(res.IM); Mul(acc,y.RE); Add(acc,y.IM);
! 180: res.RE := x/acc;
! 181: Clear(acc);
! 182: elsif Equal(y.IM,y.RE)
! 183: then acc := TWO*y.IM;
! 184: res.RE := x/acc;
! 185: res.IM := -res.RE;
! 186: Clear(acc);
! 187: else -- y.IM = -y.RE then
! 188: acc := TWO*y.IM;
! 189: res.RE := x/acc; Min(res.RE);
! 190: res.IM := x/acc; Min(res.IM);
! 191: Clear(acc);
! 192: end if;
! 193: Clear(avyim); Clear(avyre);
! 194: end if;
! 195: return res;
! 196: end "/";
! 197:
! 198: function "+" ( x : Complex_Number ) return Complex_Number is
! 199: begin
! 200: return (+x.RE,+x.IM);
! 201: end "+";
! 202:
! 203: function "-" ( x : Complex_Number ) return Complex_Number is
! 204: begin
! 205: return (-x.RE,-x.IM);
! 206: end "-";
! 207:
! 208: function "+" ( x,y : Complex_Number ) return Complex_Number is
! 209: begin
! 210: return (x.RE+y.RE,x.IM+y.IM);
! 211: end "+";
! 212:
! 213: function "-" ( x,y : Complex_Number ) return Complex_Number is
! 214: begin
! 215: return (x.RE-y.RE,x.IM-y.IM);
! 216: end "-";
! 217:
! 218: function "*" ( x,y : Complex_Number ) return Complex_Number is
! 219:
! 220: res : Complex_Number;
! 221: acc,avyim,avyre : number;
! 222:
! 223: begin
! 224: if Equal(y.IM,zero)
! 225: then res.RE := x.RE*y.RE;
! 226: res.IM := x.IM*y.RE;
! 227: elsif Equal(y.RE,zero)
! 228: then res.RE := x.IM*y.IM; Min(res.RE);
! 229: res.IM := x.RE*y.IM;
! 230: else avyre := AbsVal(y.RE); avyim := AbsVal(y.IM);
! 231: if avyre < avyim
! 232: then acc := y.RE/y.IM;
! 233: res.RE := x.RE*acc; Sub(res.RE,x.IM); Mul(res.RE,y.IM);
! 234: res.IM := x.IM*acc; Add(res.IM,x.RE); Mul(res.IM,y.IM);
! 235: Clear(acc);
! 236: elsif avyre > avyim
! 237: then acc := y.IM/y.RE;
! 238: res.RE := x.IM*acc; Sub(res.RE,x.RE); Min(res.RE);
! 239: Mul(res.RE,y.RE);
! 240: res.IM := x.RE*acc; Add(res.IM,x.IM); Mul(res.IM,y.RE);
! 241: Clear(acc);
! 242: elsif Equal(y.RE,y.IM)
! 243: then res.RE := x.RE - x.IM; Mul(res.RE,y.RE);
! 244: res.IM := x.IM + x.RE; Mul(res.IM,y.RE);
! 245: else -- y.RE = -y.IM then
! 246: res.RE := x.RE + x.IM; Mul(res.RE,y.RE);
! 247: res.IM := x.IM - x.RE; Mul(res.IM,y.RE);
! 248: end if;
! 249: Clear(avyre); Clear(avyim);
! 250: end if;
! 251: return res;
! 252: end "*";
! 253:
! 254: function "/" ( x,y : Complex_Number ) return Complex_Number is
! 255:
! 256: res : Complex_Number;
! 257: acc,avyre,avyim : number;
! 258:
! 259: begin
! 260: if Equal(y.IM,zero)
! 261: then res.RE := x.RE/y.RE;
! 262: res.IM := x.IM/y.RE;
! 263: elsif Equal(y.RE,zero)
! 264: then res.RE := x.IM/y.IM;
! 265: res.IM := x.RE/y.IM; Min(res.IM);
! 266: else avyre := AbsVal(y.RE); avyim := AbsVal(y.IM);
! 267: if avyre < avyim
! 268: then acc := y.RE/y.IM;
! 269: res.RE := x.RE*acc; Add(res.RE,x.IM);
! 270: res.IM := x.IM*acc; Sub(res.IM,x.RE);
! 271: Mul(acc,y.RE); Add(acc,y.IM);
! 272: Div(res.RE,acc);
! 273: Div(res.IM,acc);
! 274: Clear(acc);
! 275: elsif avyre > avyim
! 276: then acc := y.IM/y.RE;
! 277: res.RE := x.IM*acc; Add(res.RE,x.RE);
! 278: res.IM := x.RE*acc; Sub(res.IM,x.IM); Min(res.IM);
! 279: Mul(acc,y.IM); Add(acc,y.RE);
! 280: Div(res.RE,acc);
! 281: Div(res.IM,acc);
! 282: Clear(acc);
! 283: elsif Equal(y.RE,y.IM)
! 284: then acc := TWO*y.RE;
! 285: res.RE := x.RE + x.IM; Div(res.RE,acc);
! 286: res.IM := x.IM - x.RE; Div(res.IM,acc);
! 287: Clear(acc);
! 288: else -- y.RE = -y.IM then
! 289: acc := TWO*y.RE;
! 290: res.RE := x.RE - x.IM; Div(res.RE,acc);
! 291: res.IM := x.IM + x.RE; Div(res.IM,acc);
! 292: Clear(acc);
! 293: end if;
! 294: Clear(avyre); Clear(avyim);
! 295: end if;
! 296: return res;
! 297: end "/";
! 298:
! 299: function "**" ( x : Complex_Number; m : integer ) return Complex_Number is
! 300:
! 301: res : Complex_Number;
! 302:
! 303: begin
! 304: if m = 0
! 305: then res := Create(1);
! 306: elsif m > 0
! 307: then res := +x; --Copy(x,res); <- Copy CAUSED BUS ERROR!!
! 308: for j in 2..m loop
! 309: Mul(res,x);
! 310: end loop;
! 311: else res := Create(1);
! 312: for j in 1..(-m) loop
! 313: Div(res,x);
! 314: end loop;
! 315: end if;
! 316: return res;
! 317: end "**";
! 318:
! 319: -- ARITHMETIC OPERATIONS AS PROCEDURES :
! 320:
! 321: procedure Add ( x : in out Complex_Number; y : in number ) is
! 322: begin
! 323: Add(x.RE,y);
! 324: end Add;
! 325:
! 326: procedure Sub ( x : in out Complex_Number; y : in number ) is
! 327: begin
! 328: Sub(x.RE,y);
! 329: end Sub;
! 330:
! 331: procedure Mul ( x : in out Complex_Number; y : in number ) is
! 332: begin
! 333: Mul(x.RE,y);
! 334: Mul(x.IM,y);
! 335: end Mul;
! 336:
! 337: procedure Div ( x : in out Complex_Number; y : in number ) is
! 338: begin
! 339: Div(x.RE,y);
! 340: Div(x.IM,y);
! 341: end Div;
! 342:
! 343: procedure Add ( x : in out Complex_Number; y : in Complex_Number ) is
! 344: begin
! 345: Add(x.RE,y.RE);
! 346: Add(x.IM,y.IM);
! 347: end Add;
! 348:
! 349: procedure Sub ( x : in out Complex_Number; y : in Complex_Number ) is
! 350: begin
! 351: Sub(x.RE,y.RE);
! 352: Sub(x.IM,y.IM);
! 353: end Sub;
! 354:
! 355: procedure Min ( x : in out Complex_Number ) is
! 356: begin
! 357: Min(x.RE);
! 358: Min(x.IM);
! 359: end Min;
! 360:
! 361: procedure Mul ( x : in out Complex_Number; y : in Complex_Number ) is
! 362:
! 363: res : Complex_Number;
! 364: acc,avyim,avyre : number;
! 365:
! 366: begin
! 367: if Equal(y.IM,zero)
! 368: then Mul(x.RE,y.RE);
! 369: Mul(x.IM,y.RE);
! 370: elsif Equal(y.RE,zero)
! 371: then res.RE := x.IM*y.IM; Min(res.RE);
! 372: res.IM := x.RE*y.IM;
! 373: Clear(x); x := res;
! 374: else avyre := AbsVal(y.RE); avyim := AbsVal(y.IM);
! 375: if avyre < avyim
! 376: then acc := y.RE/y.IM;
! 377: res.RE := x.RE*acc; Sub(res.RE,x.IM); Mul(res.RE,y.IM);
! 378: res.IM := x.IM*acc; Add(res.IM,x.RE); Mul(res.IM,y.IM);
! 379: Clear(acc);
! 380: elsif avyre > avyim
! 381: then acc := y.IM/y.RE;
! 382: res.RE := x.IM*acc; Sub(res.RE,x.RE); Min(res.RE);
! 383: Mul(res.RE,y.RE);
! 384: res.IM := x.RE*acc; Add(res.IM,x.IM); Mul(res.IM,y.RE);
! 385: Clear(acc);
! 386: elsif Equal(y.RE,y.IM)
! 387: then res.RE := x.RE - x.IM; Mul(res.RE,y.RE);
! 388: res.IM := x.IM + x.RE; Mul(res.IM,y.RE);
! 389: else -- y.RE = -y.IM then
! 390: res.RE := x.RE + x.IM; Mul(res.RE,y.RE);
! 391: res.IM := x.IM - x.RE; Mul(res.IM,y.RE);
! 392: end if;
! 393: Clear(avyre); Clear(avyim);
! 394: Clear(x); x := res;
! 395: end if;
! 396: end Mul;
! 397:
! 398: procedure Div ( x : in out Complex_Number; y : in Complex_Number ) is
! 399:
! 400: res : Complex_Number;
! 401: acc,avyre,avyim : number;
! 402:
! 403: begin
! 404: if Equal(y.IM,zero)
! 405: then Div(x.RE,y.RE);
! 406: Div(x.IM,y.RE);
! 407: elsif Equal(y.IM,zero)
! 408: then res.RE := x.IM/y.IM;
! 409: res.IM := x.RE/y.IM; Min(res.IM);
! 410: Clear(x); x := res;
! 411: else avyre := AbsVal(y.RE); avyim := AbsVal(y.IM);
! 412: if avyre < avyim
! 413: then acc := y.RE/y.IM;
! 414: res.RE := x.RE*acc; Add(res.RE,x.IM);
! 415: res.IM := x.IM*acc; Sub(res.IM,x.RE);
! 416: Mul(acc,y.RE); Add(acc,y.IM);
! 417: Div(res.RE,acc);
! 418: Div(res.IM,acc);
! 419: Clear(acc);
! 420: elsif avyre > avyim
! 421: then acc := y.IM/y.RE;
! 422: res.RE := x.IM*acc; Add(res.RE,x.RE);
! 423: res.IM := x.RE*acc; Sub(res.IM,x.IM); Min(res.IM);
! 424: Mul(acc,y.IM); Add(acc,y.RE);
! 425: Div(res.RE,acc);
! 426: Div(res.IM,acc);
! 427: Clear(acc);
! 428: elsif Equal(y.RE,y.IM)
! 429: then acc := TWO*y.RE;
! 430: res.RE := x.RE + x.IM; Div(res.RE,acc);
! 431: res.IM := x.IM - x.RE; Div(res.IM,acc);
! 432: Clear(acc);
! 433: else -- y.RE = -y.IM then
! 434: acc := TWO*y.RE;
! 435: res.RE := x.RE - x.IM; Div(res.RE,acc);
! 436: res.IM := x.IM + x.RE; Div(res.IM,acc);
! 437: Clear(acc);
! 438: end if;
! 439: Clear(avyre); Clear(avyim);
! 440: Clear(x); x := res;
! 441: end if;
! 442: end Div;
! 443:
! 444: -- DESTRUCTOR :
! 445:
! 446: procedure Clear ( x : in out Complex_Number ) is
! 447: begin
! 448: Clear(x.RE);
! 449: Clear(x.IM);
! 450: end Clear;
! 451:
! 452: end Generic_Complex_Numbers;
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to generic_complex_numbers.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Math_Lib / Numbers