001// --- BEGIN LICENSE BLOCK --- 002/* 003 * Copyright (c) 2009, Mikio L. Braun 004 * All rights reserved. 005 * 006 * Redistribution and use in source and binary forms, with or without 007 * modification, are permitted provided that the following conditions are 008 * met: 009 * 010 * * Redistributions of source code must retain the above copyright 011 * notice, this list of conditions and the following disclaimer. 012 * 013 * * Redistributions in binary form must reproduce the above 014 * copyright notice, this list of conditions and the following 015 * disclaimer in the documentation and/or other materials provided 016 * with the distribution. 017 * 018 * * Neither the name of the Technische Universität Berlin nor the 019 * names of its contributors may be used to endorse or promote 020 * products derived from this software without specific prior 021 * written permission. 022 * 023 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 024 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 025 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 026 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT 027 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 028 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT 029 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, 030 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 031 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 032 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE 033 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 034 */ 035// --- END LICENSE BLOCK --- 036 037package org.jblas; 038 039/** 040 * A complex value with double precision. 041 * 042 * @author Mikio L. Braun 043 */ 044public class ComplexDouble { 045 046 private double r, i; 047 public static final ComplexDouble UNIT = new ComplexDouble(1.0, 0.0); 048 public static final ComplexDouble I = new ComplexDouble(0.0, 1.0); 049 public static final ComplexDouble NEG_UNIT = new ComplexDouble(-1.0, 0.0); 050 public static final ComplexDouble NEG_I = new ComplexDouble(0.0, -1.0); 051 public static final ComplexDouble ZERO = new ComplexDouble(0.0); 052 053 public ComplexDouble(double real, double imag) { 054 r = real; 055 i = imag; 056 } 057 058 public ComplexDouble(double real) { 059 this(real, 0.0); 060 } 061 062 public String toString() { 063 if (i >= 0) { 064 return r + " + " + i + "i"; 065 } else { 066 return r + " - " + (-i) + "i"; 067 } 068 } 069 070 public ComplexDouble set(double real, double imag) { 071 r = real; 072 i = imag; 073 return this; 074 } 075 076 public double real() { 077 return r; 078 } 079 080 public double imag() { 081 return i; 082 } 083 084 public ComplexDouble dup() { 085 return new ComplexDouble(r, i); 086 } 087 088 public ComplexDouble copy(ComplexDouble other) { 089 r = other.r; 090 i = other.i; 091 return this; 092 } 093 094 /** 095 * Add two complex numbers in-place 096 * 097 * @param c other complex number 098 * @param result complex number where result is stored 099 * @return same as result 100 */ 101 public ComplexDouble addi(ComplexDouble c, ComplexDouble result) { 102 if (this == result) { 103 r += c.r; 104 i += c.i; 105 } else { 106 result.r = r + c.r; 107 result.i = i + c.i; 108 } 109 return result; 110 } 111 112 /** 113 * Add two complex numbers in-place storing the result in this. 114 * 115 * @param c other complex number 116 * @return resulting complex number 117 */ 118 public ComplexDouble addi(ComplexDouble c) { 119 return addi(c, this); 120 } 121 122 /** 123 * Add two complex numbers. 124 * 125 * @param c other complex number 126 * @return new complex number with result 127 */ 128 public ComplexDouble add(ComplexDouble c) { 129 return dup().addi(c); 130 } 131 132 /** 133 * Add a real number to a complex number in-place. 134 * 135 * @param a real number to add 136 * @param result complex number to hold result 137 * @return same as result 138 */ 139 public ComplexDouble addi(double a, ComplexDouble result) { 140 if (this == result) { 141 r += a; 142 } else { 143 result.r = r + a; 144 result.i = i; 145 } 146 return result; 147 } 148 149 /** 150 * Add a real number to complex number in-place, storing the result in this. 151 * 152 * @param c real number to add 153 * @return resulting complex number 154 */ 155 public ComplexDouble addi(double c) { 156 return addi(c, this); 157 } 158 159 /** 160 * Add a real number to a complex number. 161 * 162 * @param c real number to add 163 * @return new complex number with result 164 */ 165 public ComplexDouble add(double c) { 166 return dup().addi(c); 167 } 168 169 /** 170 * Subtract two complex numbers, in-place 171 * 172 * @param c complex number to subtract 173 * @param result resulting complex number 174 * @return same as result 175 */ 176 public ComplexDouble subi(ComplexDouble c, ComplexDouble result) { 177 if (this == result) { 178 r -= c.r; 179 i -= c.i; 180 } else { 181 result.r = r - c.r; 182 result.i = i - c.i; 183 } 184 return this; 185 } 186 187 public ComplexDouble subi(ComplexDouble c) { 188 return subi(c, this); 189 } 190 191 /** 192 * Subtract two complex numbers 193 * 194 * @param c complex number to subtract 195 * @return new complex number with result 196 */ 197 public ComplexDouble sub(ComplexDouble c) { 198 return dup().subi(c); 199 } 200 201 public ComplexDouble subi(double a, ComplexDouble result) { 202 if (this == result) { 203 r -= a; 204 } else { 205 result.r = r - a; 206 result.i = i; 207 } 208 return result; 209 } 210 211 public ComplexDouble subi(double a) { 212 return subi(a, this); 213 } 214 215 public ComplexDouble sub(double r) { 216 return dup().subi(r); 217 } 218 219 /** 220 * Multiply two complex numbers, in-place 221 * 222 * @param c other complex number 223 * @param result complex number where product is stored 224 * @return same as result 225 */ 226 public ComplexDouble muli(ComplexDouble c, ComplexDouble result) { 227 double newR = r * c.r - i * c.i; 228 double newI = r * c.i + i * c.r; 229 result.r = newR; 230 result.i = newI; 231 return result; 232 } 233 234 public ComplexDouble muli(ComplexDouble c) { 235 return muli(c, this); 236 } 237 238 /** 239 * Multiply two complex numbers 240 * 241 * @param c other complex number 242 * @return new complex number with product of this and c 243 */ 244 public ComplexDouble mul(ComplexDouble c) { 245 return dup().muli(c); 246 } 247 248 public ComplexDouble mul(double v) { 249 return dup().muli(v); 250 } 251 252 public ComplexDouble muli(double v, ComplexDouble result) { 253 if (this == result) { 254 r *= v; 255 i *= v; 256 } else { 257 result.r = r * v; 258 result.i = i * v; 259 } 260 return this; 261 } 262 263 public ComplexDouble muli(double v) { 264 return muli(v, this); 265 } 266 267 /** 268 * Divide two complex numbers 269 * 270 * @param c complex number to divide this by 271 * @return new complex number with quotient of this and c 272 */ 273 public ComplexDouble div(ComplexDouble c) { 274 return dup().divi(c); 275 } 276 277 /** 278 * Divide two complex numbers, in-place 279 * 280 * @param c complex number to divide this by 281 * @param result complex number to hold result 282 * @return same as result 283 */ 284 public ComplexDouble divi(ComplexDouble c, ComplexDouble result) { 285 double d = c.r * c.r + c.i * c.i; 286 double newR = (r * c.r + i * c.i) / d; 287 double newI = (i * c.r - r * c.i) / d; 288 result.r = newR; 289 result.i = newI; 290 return result; 291 } 292 293 public ComplexDouble divi(ComplexDouble c) { 294 return divi(c, this); 295 } 296 297 public ComplexDouble divi(double v, ComplexDouble result) { 298 if (this == result) { 299 r /= v; 300 i /= v; 301 } else { 302 result.r = r / v; 303 result.i = i / v; 304 } 305 return this; 306 } 307 308 public ComplexDouble divi(double v) { 309 return divi(v, this); 310 } 311 312 public ComplexDouble div(double v) { 313 return dup().divi(v); 314 } 315 316 /** 317 * Return the absolute value 318 * 319 * @return the result (length of the vector in 2d plane) 320 */ 321 public double abs() { 322 return (double) Math.sqrt(r * r + i * i); 323 } 324 325 /** 326 * Returns the argument of a complex number. 327 * 328 * @return the result (angle in radians of the vector in 2d plane) 329 */ 330 public double arg() { 331 return (double) Math.atan2(i, r); 332 } 333 334 public ComplexDouble invi() { 335 double d = r * r + i * i; 336 r = r / d; 337 i = -i / d; 338 return this; 339 } 340 341 public ComplexDouble inv() { 342 return dup().invi(); 343 } 344 345 public ComplexDouble neg() { 346 return dup().negi(); 347 } 348 349 public ComplexDouble negi() { 350 r = -r; 351 i = -i; 352 return this; 353 } 354 355 public ComplexDouble conji() { 356 i = -i; 357 return this; 358 } 359 360 public ComplexDouble conj() { 361 return dup().conji(); 362 } 363 364 public ComplexDouble sqrt() { 365 double a = abs(); 366 double s2 = (double) Math.sqrt(2); 367 double p = (double) Math.sqrt(a + r) / s2; 368 double sgn = Math.signum(i); 369 if (sgn == 0.0) { 370 sgn = 1.0; 371 } 372 double q = (double) Math.sqrt(a - r) / s2 * Math.signum(sgn); 373 return new ComplexDouble(p, q); 374 } 375 376 /** 377 * Comparing two ComplexDouble values. 378 * 379 * @param o object to compare this against 380 * @return true if both numbers have the same value 381 */ 382 public boolean equals(Object o) { 383 if (!(o instanceof ComplexDouble)) { 384 return false; 385 } 386 ComplexDouble c = (ComplexDouble) o; 387 388 return r == c.r && i == c.i; 389 } 390 391 public int hashCode() { 392 return Double.valueOf(r).hashCode() ^ Double.valueOf(i).hashCode(); 393 } 394 395 public boolean eq(ComplexDouble c) { 396 return Math.abs(r - c.r) + Math.abs(i - c.i) < (double) 1e-6; 397 } 398 399 public boolean ne(ComplexDouble c) { 400 return !eq(c); 401 } 402 403 public boolean isZero() { 404 return r == 0.0 && i == 0.0; 405 } 406 407 public boolean isReal() { 408 return i == 0.0; 409 } 410 411 public boolean isImag() { 412 return r == 0.0; 413 } 414}