Complex numbers in JavaScript Our goal for today is to perform the dot product of two complex vectors in JavaScript. To make it a bit more challenging we are going to make use of typed arrays for efficient data storage. JavaScript typed arrays are array-like objects and provide a mechanism for accessing raw binary data. As you may already know, Array objects grow and shrink dynamically and can have any JavaScript value. JavaScript engines perform optimizations so that these arrays are fast. However, as web applications become more and more powerful, adding features such as audio and video manipulation, access to raw data using WebSockets, and so forth, it has become clear that there are times when it would be helpful for JavaScript code to be able to quickly and easily manipulate raw binary data in typed arrays. –Mozilla Developer Network The basis of our implementation will be the zdotc routine from the Reference BLAS package by Netlib. We'll go through that a bit later, but let's start with a crash course in complex arithmetic. A complex number is a number that can be expressed in the form , where and are real numbers, and is a solution of the equation . –Wikipedia So, a complex number consists of two parts, one real part and one imaginary part . We can also imagine it as a coordinate in a 2 dimensional space where is the x axis and is the y axis. Enough chit chat, let's write it up!
class Complex16 extends Float64Array { constructor(...args) { super(...args); if (this.length !== 2) { throw new Error('invalid complex number'); } } // setters set re(value) { this[0] = value; } set im(value) { this[1] = value; } // getters get re() { return this[0]; } get im() { return this[1]; } // for visualization purposes only toString() { const { im, re } = this; let real = ''; if (re !== 0 || im === 0) { real = re; } let op = ''; if (re !== 0 && im > 0) { op = '+'; } let imaginary = ''; if (im !== 0) { if (im === 1) { imaginary = `i`; } else if (im === -1) { imaginary = `-i` } else { imaginary = `${im}i`; } } return `${real}${op}${imaginary}`; } }
For readability purposes I'm subclassing Float64Array, which let's us borrow some of its internal methods. We are ready to create our first complex number! Hint: Try changing a number and press ⇧↩︎ on your keyboard.
{ const z = new Complex16(2); z.re = 4; z.im = 3; return tex.block` \begin{aligned} z &= ${z} \\ \operatorname{Re}(z) &= ${z.re} \\ \operatorname{Im}(z) &= ${z.im} \\ \end{aligned} `; }
Addition To reach our final goal, we need to build a toolbox of functions that operate on complex numbers. Our first member is complex addition and it is performed by adding the real and imaginary parts of two numbers separately.
function add(z1, z2) { return new Complex16([ z1.re + z2.re, z1.im + z2.im ]); }
{ const z1 = new Complex16([2, 1]); const z2 = new Complex16([3, 1]); return tex.block` \begin{aligned} z_1 &= ${z1} \\ z_2 &= ${z2} \\ z_1+z_2 &= ${add(z1, z2)} \end{aligned} `; }
Multiplication The second tool we are going to need is complex multiplication. What we need to keep in mind here is that the square of equals . Otherwise it's just normal multiplication.
function mul(z1, z2) { return new Complex16([ z1.re * z2.re - z1.im * z2.im, z1.re * z2.im + z1.im * z2.re ]); }
{ const z1 = new Complex16([2, -1]); const z2 = new Complex16([3, 1]); return tex.block` \begin{aligned} z_1 &= ${z1} \\ z_2 &= ${z2} \\ z_1z_2 &= (${new Complex16([z1.re * z2.re, 0])}) + (${new Complex16([0, z1.re * z2.im])}) + (${new Complex16([0, z1.im * z2.re])}) + (${new Complex16([0, z1.im * z2.im])}^2) \\ &= ${mul(z1, z2)} \\ \end{aligned} `; }
Conjugation Now we just need one more piece of arithmetic, the complex conjugate. In mathematics, the complex conjugate of a complex number is the number with an equal real part and an imaginary part equal in magnitude but opposite in sign. – Wikipedia Well, that's not even a challenge. We just have to negate the imaginary part.
function conjg(z) { return new Complex16([z.re, -z.im]); }
{ const z = new Complex16([1, 2]); return tex.block` \begin{aligned} z &= ${z} \\ \overline{z} &= ${conjg(z)} \end{aligned} `; }
Up until now we have been working with a single complex value, Complex16, for our calculations. To be able to compute the dot product of two complex vectors, we need to introduce a complex vector data structure. If we were lazy we could just save a bunch of Complex16 instances in a JavaScript Array and loop through it. But we like it fast, right? Regular JavaScript arrays are not. It turns out that an efficient way of packing complex vectors is to store one complex number after the other, like this: That is good news for us, it means we can still use Float64Array as our base class, just in a slightly different way.
class Complex16Array extends Float64Array { constructor(...args) { super(...args); if (this.length % 2 !== 0) { throw new Error('you supplied an uneven amount of elements'); } } // get the ith complex number from the vector get(i) { const pos = i * 2; this.checkBounds(pos); // see https://developer.mozilla.org/en-US/docs/Web/JavaScript/Reference/Global_Objects/TypedArray/subarray return new Complex16(this.subarray(pos, pos + 2)); } // set the ith complex number in the vector set(i, z) { const pos = i * 2; this.checkBounds(pos); this[pos] = z.re; this[pos + 1] = z.im; } // helper to check if accessed index is out of bounds checkBounds(pos) { if (pos < 0 || pos > this.length - 2) { throw new Error('index out of bounds'); } } // helper to get count of complex numbers instead of actual length of the array get size() { return this.length / 2; } // again, just for visualization purposes... toString() { const res = []; for (let i = 0; i < this.size; i++) { res.push(this.get(i)); } return `(${res.join(', ')})`; } }
Here is our complex vector in action:
{ // allocate new complex vector with 4 memory slots const u = new Complex16Array(4); u.set(0, new Complex16([1, 2])); u.set(1, new Complex16([3, 4])); const out = []; for (let i = 0; i < u.size; i++) { out.push(`u_${i} &= ${u.get(i)} \\\\`); } return tex.block` \begin{aligned} \vec u &= ${u} \\ ${out.join('')} \end{aligned} `; }
Dot product The complex dot product is (usually) defined as the sum of the th conjugated element of multiplied by the th element of . I mentioned the zdotc routine in the beginning of the article, so we're going to take a look at it now. For those of you that have never heard of BLAS, it stands for Basic Linear Algebra Subprograms and they are the de facto standard low level routines for linear algebra libraries. They are written in Fortran. To make it a bit shorter I took the liberty of ignoring increments (incx and incy) and removing comments. You can see the full source code here. COMPLEX*16 FUNCTION zdotc(N,ZX,ZY) * * .. Scalar Arguments .. INTEGER N * .. * .. Array Arguments .. COMPLEX*16 ZX(*),ZY(*) * * .. Local Scalars .. COMPLEX*16 ZTEMP INTEGER I * .. * .. Intrinsic Functions .. INTRINSIC dconjg * .. ztemp = (0.0d0,0.0d0) zdotc = (0.0d0,0.0d0) IF (n.LE.0) RETURN DO i = 1,n ztemp = ztemp + dconjg(zx(i))*zy(i) END DO zdotc = ztemp RETURN END Not too hard to translate!
function zdotc(n, zx, zy) { let ztemp = new Complex16([0, 0]); if (n < 0) { return ztemp; } for (let i = 0; i < n; i++) { ztemp = add(ztemp, mul(conjg(zx.get(i)), zy.get(i))); } return ztemp; }
{ const n = 3; const u = new Complex16Array([1, 2, 3, 4, 5, 6]); const v = new Complex16Array([4, 3, 2, 1, 0, -1]); return tex.block` \begin{aligned} \vec{u} &= ${u} \\ \vec{v} &= ${v} \\ \vec{u}\cdot\vec{v} &= \sum_{i = 0}^n{\overline{u_i} v_i} \\ &= ${zdotc(n, u, v)} \end{aligned} `; }
If you have feedback, questions, or maybe even suggestions for a future post you'll find my contact info in my GitHub profile.