001: /*
002: * Created on Mar 29, 2004
003: *
004: */
005: package org.jmatlab.linalg;
006:
007: import java.io.Serializable;
008:
009: import org.jmatlab.toolbox.*;
010: import org.jmatlab.util.FormatNumber;
011:
012: /**
013: * @author Ali
014: *
015: */
016: public class DefaultComplexImpl implements IComplex, Serializable {
017:
018: private double real = 0.0;
019: private double imag = 0.0;
020:
021: public DefaultComplexImpl() {
022: this .real = 0.0;
023: this .imag = 0.0;
024: }
025:
026: public DefaultComplexImpl(double d) {
027: this .real = d;
028: }
029:
030: public DefaultComplexImpl(double real, double imag) {
031: this .real = real;
032: this .imag = imag;
033: }
034:
035: public DefaultComplexImpl(IComplex complex) {
036: this .real = complex.getReal();
037: this .imag = complex.getImag();
038: }
039:
040: public double getImag() {
041: return imag;
042: }
043:
044: public void setImag(double d) {
045: this .imag = d;
046: }
047:
048: public double getReal() {
049: return real;
050: }
051:
052: public void setReal(double d) {
053: this .real = d;
054: }
055:
056: public boolean equals(IComplex complex) {
057: if (isNAN() && complex.isNAN()) {
058: return true;
059: } else {
060: return (real == complex.getReal() && imag == complex
061: .getImag());
062: }
063: }
064:
065: public boolean equals(Object obj) {
066: if (obj == null) {
067: return false;
068: } else if (obj instanceof IComplex) {
069: return (equals((IComplex) obj));
070: } else {
071: return false;
072: }
073: }
074:
075: public boolean isNAN() {
076: return (Double.isNaN(real) || Double.isNaN(imag));
077: }
078:
079: public IComplex minus(IComplex complex) {
080: return new DefaultComplexImpl(real - complex.getReal(), imag
081: - complex.getImag());
082: }
083:
084: public IComplex minus(double d) {
085: return new DefaultComplexImpl(real - d, imag);
086: }
087:
088: public IComplex minusReverse(double d) {
089: return new DefaultComplexImpl(d - real, -1 * imag);
090: }
091:
092: private static double copysign(double a, double b) {
093: double abs = Math.abs(a);
094: return ((b < 0) ? -abs : abs);
095: }
096:
097: private static boolean isFinite(double x) {
098: return !(Double.isInfinite(x) || Double.isNaN(x));
099: }
100:
101: private static IComplex divideBy(IComplex x, IComplex y) {
102: double a = x.getReal();
103: double b = x.getImag();
104: double c = y.getReal();
105: double d = y.getImag();
106:
107: double scale = Math.max(Math.abs(c), Math.abs(d));
108: boolean isScaleFinite = isFinite(scale);
109: if (isScaleFinite) {
110: c /= scale;
111: d /= scale;
112: }
113:
114: double den = c * c + d * d;
115: IComplex z = new DefaultComplexImpl((a * c + b * d) / den, (b
116: * c - a * d)
117: / den);
118:
119: if (isScaleFinite) {
120: z.setReal(z.getReal() / scale);
121: z.setImag(z.getImag() / scale);
122: }
123:
124: // Recover infinities and zeros computed as NaN+iNaN.
125: if (Double.isNaN(z.getReal()) && Double.isNaN(z.getImag())) {
126: if (den == 0.0 && (!Double.isNaN(a) || !Double.isNaN(b))) {
127: double s = copysign(Double.POSITIVE_INFINITY, c);
128: z.setReal(s * a);
129: z.setImag(s * b);
130:
131: } else if ((Double.isInfinite(a) || Double.isInfinite(b))
132: && isFinite(c) && isFinite(d)) {
133: a = copysign(Double.isInfinite(a) ? 1.0 : 0.0, a);
134: b = copysign(Double.isInfinite(b) ? 1.0 : 0.0, b);
135: z.setReal(Double.POSITIVE_INFINITY * (a * c + b * d));
136: z.setImag(Double.POSITIVE_INFINITY * (b * c - a * d));
137:
138: } else if (Double.isInfinite(scale) && isFinite(a)
139: && isFinite(b)) {
140: c = copysign(Double.isInfinite(c) ? 1.0 : 0.0, c);
141: d = copysign(Double.isInfinite(d) ? 1.0 : 0.0, d);
142: z.setReal(0.0 * (a * c + b * d));
143: z.setImag(0.0 * (b * c - a * d));
144: }
145: }
146: return z;
147: }
148:
149: public IComplex dividedBy(IComplex complex) {
150: return divideBy(this , complex);
151: }
152:
153: private static IComplex divideBy(IComplex x, double y) {
154: return new DefaultComplexImpl(x.getReal() / y, x.getImag() / y);
155: }
156:
157: public IComplex dividedBy(double d) {
158: return divideBy(this , d);
159: }
160:
161: public IComplex dividedByReverse(double d) {
162: double den, t;
163: IComplex z;
164: if (Math.abs(real) > Math.abs(imag)) {
165: t = imag / real;
166: den = real + imag * t;
167: z = new DefaultComplexImpl(d / den, -d * t / den);
168: } else {
169: t = real / imag;
170: den = imag + real * t;
171: z = new DefaultComplexImpl(d * t / den, -d / den);
172: }
173: return z;
174: }
175:
176: public IComplex plus(IComplex complex) {
177: return new DefaultComplexImpl(real + complex.getReal(), imag
178: + complex.getImag());
179: }
180:
181: public IComplex plus(double d) {
182: return new DefaultComplexImpl(real + d, imag);
183: }
184:
185: private static void timesNaN(IComplex x, IComplex y, IComplex t) {
186: boolean recalc = false;
187: double a = x.getReal();
188: double b = x.getImag();
189: double c = y.getReal();
190: double d = y.getImag();
191:
192: if (Double.isInfinite(a) || Double.isInfinite(b)) {
193: // x is infinite
194: a = copysign(Double.isInfinite(a) ? 1.0 : 0.0, a);
195: b = copysign(Double.isInfinite(b) ? 1.0 : 0.0, b);
196: if (Double.isNaN(c))
197: c = copysign(0.0, c);
198: if (Double.isNaN(d))
199: d = copysign(0.0, d);
200: recalc = true;
201: }
202:
203: if (Double.isInfinite(c) || Double.isInfinite(d)) {
204: // x is infinite
205: a = copysign(Double.isInfinite(c) ? 1.0 : 0.0, c);
206: b = copysign(Double.isInfinite(d) ? 1.0 : 0.0, d);
207: if (Double.isNaN(a))
208: a = copysign(0.0, a);
209: if (Double.isNaN(b))
210: b = copysign(0.0, b);
211: recalc = true;
212: }
213:
214: if (!recalc) {
215: if (Double.isInfinite(a * c) || Double.isInfinite(b * d)
216: || Double.isInfinite(a * d)
217: || Double.isInfinite(b * c)) {
218: // Change all NaNs to 0
219: if (Double.isNaN(a))
220: a = copysign(0.0, a);
221: if (Double.isNaN(b))
222: b = copysign(0.0, b);
223: if (Double.isNaN(c))
224: c = copysign(0.0, c);
225: if (Double.isNaN(d))
226: d = copysign(0.0, d);
227: recalc = true;
228: }
229: }
230:
231: if (recalc) {
232: t.setReal(Double.POSITIVE_INFINITY * (a * c - b * d));
233: t.setImag(Double.POSITIVE_INFINITY * (a * d + b * c));
234: }
235: }
236:
237: public static IComplex times(IComplex x, IComplex y) {
238: IComplex t = new DefaultComplexImpl(x.getReal() * y.getReal()
239: - x.getImag() * y.getImag(), x.getReal() * y.getImag()
240: + x.getImag() * y.getReal());
241: if (Double.isNaN(t.getReal()) && Double.isNaN(t.getImag()))
242: timesNaN(x, y, t);
243: return t;
244: }
245:
246: public IComplex times(IComplex complex) {
247: return times(this , complex);
248: }
249:
250: public IComplex times(double d) {
251: return new DefaultComplexImpl(real * d, imag * d);
252: }
253:
254: public String toString() {
255: if (imag == 0.0) {
256: if (Double.isInfinite(real)) {
257: return "Inf";
258: } else if (Double.isInfinite(-real)) {
259: return "-Inf";
260: } else if (Double.isNaN(real)) {
261: return "NaN";
262: }
263: return FormatNumber.format(real);
264: }
265:
266: String iota = "i";
267:
268: String sign = (imag < 0.0) ? " - " : " + ";
269: String re = FormatNumber.format(real);
270: if (Double.isInfinite(real)) {
271: re = "Inf";
272: } else if (Double.isInfinite(-real)) {
273: re = "-Inf";
274: } else if (Double.isNaN(real)) {
275: re = "NaN";
276: }
277: String im = null;
278: if (Math.abs(imag) == 1) {
279: im = "";
280: } else {
281: im = FormatNumber.format(Math.abs(imag));
282: if (Double.isInfinite(imag)) {
283: im = "Inf";
284: } else if (Double.isInfinite(-imag)) {
285: im = "-Inf";
286: } else if (Double.isNaN(imag)) {
287: im = "NaN";
288: }
289: }
290: return re + sign + im + iota;
291: }
292:
293: public IComplex conjugate() {
294: IComplex c = new DefaultComplexImpl(real, -1 * imag);
295: return c;
296: }
297:
298: public static void main(String[] args) {
299: long start = System.currentTimeMillis();
300: long n = 1000000000;
301: for (int i = 0; i < n; i++) {
302: // Complex c1 = new DefaultComplexImpl(1.2, -9.8);
303: // Complex c2 = new DefaultComplexImpl(-0.2, 0.8);
304: // Complex c = c1.times(c2);
305: double d1 = 2.39876;
306: double d2 = -0.00000001;
307: double d = d1 * d2;
308: }
309: long end = System.currentTimeMillis();
310: long diff = end - start;
311: double t = diff / 1000;
312: System.out.println("Number of multiplications = " + n);
313: System.out.println("Number of mili seconds = " + diff);
314: }
315:
316: public IComplex negative() {
317: IComplex c = null;
318: if (real != 0 && imag != 0) {
319: c = new DefaultComplexImpl(-1 * real, -1 * imag);
320: } else if (real == 0) {
321: c = new DefaultComplexImpl(0, -1 * imag);
322: } else if (imag == 0) {
323: c = new DefaultComplexImpl(-1 * real, 0);
324: }
325: return c;
326: }
327:
328: public IComplex pow(double x) {
329: double absz = abs();
330: IComplex result = null;
331: if (absz == 0.0) {
332: result = this ;
333: } else {
334: double e = Math.pow(absz, x);
335: if (imag == 0) {
336: result = new DefaultComplexImpl(e, 0.0);
337: } else {
338: double a = argument();
339: result = new DefaultComplexImpl(e * Math.cos(x * a), e
340: * Math.sin(x * a));
341: }
342: }
343: return result;
344:
345: }
346:
347: public IComplex pow(IComplex c) {
348: IComplex result = null;
349: //TODO method LOG should be moved to this class
350: result = Elfun.EXP(Elfun.LOG(this ).times(c));
351: return result;
352:
353: }
354:
355: public double abs() {
356: double x = Math.abs(real);
357: double y = Math.abs(imag);
358:
359: if (Double.isInfinite(x) || Double.isInfinite(y))
360: return Double.POSITIVE_INFINITY;
361:
362: if (x + y == 0.0) {
363: return 0.0;
364: } else if (x > y) {
365: y /= x;
366: return x * Math.sqrt(1.0 + y * y);
367: } else {
368: x /= y;
369: return y * Math.sqrt(x * x + 1.0);
370: }
371: }
372:
373: public double argument() {
374: return Math.atan2(imag, real);
375: }
376:
377: public boolean lessThan(IComplex c) {
378: if (abs() < c.abs()) {
379: return true;
380: }
381: return false;
382: }
383:
384: public boolean lessThan(double d) {
385: boolean rtn = false;
386: if (imag == 0.0) {
387: if (real < d)
388: rtn = true;
389: } else {
390: if (abs() < d)
391: rtn = true;
392: }
393: return rtn;
394: }
395:
396: public boolean greaterThan(IComplex c) {
397: if (abs() > c.abs()) {
398: return true;
399: }
400: return false;
401: }
402:
403: public boolean greaterThan(double d) {
404: boolean rtn = false;
405: if (imag == 0.0) {
406: if (real > d)
407: rtn = true;
408: } else {
409: if (abs() > d)
410: rtn = true;
411: }
412: return rtn;
413: }
414:
415: public boolean lessThanEqual(IComplex c) {
416: if (abs() <= c.abs()) {
417: return true;
418: }
419: return false;
420: }
421:
422: public boolean greaterThanEqual(IComplex c) {
423: if (abs() >= c.abs()) {
424: return true;
425: }
426: return false;
427: }
428:
429: public boolean isReal() {
430: return !isComplex();
431: }
432:
433: public boolean isComplex() {
434: if (imag != 0.0)
435: return true;
436: return false;
437: }
438:
439: public IComplex sqr() {
440: return new DefaultComplexImpl(real * real - imag * imag, 2.0
441: * real * imag);
442: }
443:
444: }
|