From e24b81a206123f6e0bdda81c28791eb849f92c02 Mon Sep 17 00:00:00 2001 From: Glen Whitney Date: Fri, 12 Dec 2025 23:29:49 -0800 Subject: [PATCH] feat: add cbrt for complex numbers and quaternions, returning three values --- src/complex/__test__/arithmetic.spec.js | 11 +++++++ src/complex/__test__/type.spec.js | 7 +++++ src/complex/arithmetic.js | 41 +++++++++++++++++++++++++ src/complex/type.js | 34 +++++++++++--------- 4 files changed, 79 insertions(+), 14 deletions(-) diff --git a/src/complex/__test__/arithmetic.spec.js b/src/complex/__test__/arithmetic.spec.js index 6194370..7c8c4c7 100644 --- a/src/complex/__test__/arithmetic.spec.js +++ b/src/complex/__test__/arithmetic.spec.js @@ -132,4 +132,15 @@ describe('complex arithmetic operations', () => { assert.deepStrictEqual( sub(cplx(z,z), cplx(10,20)), cplx(cplx(-7, 4), cplx(-17, 4))) }) + it('cube roots complex numbers and quaternions', () => { + assert(math.equal(math.cbrt(cplx(0, 8))[0], cplx(Math.sqrt(3), 1))) + assert.deepStrictEqual( + math.cbrt(cplx(2, 3)), + [cplx(1.4518566183526649, 0.49340353410400467), + cplx(-1.1532283040274218, 1.0106429470939742), + cplx(-0.29862831432524256, -1.5040464811979786)]) + const quat = cplx(cplx(1, 2), cplx(2, 1)) + const root = math.cbrt(quat)[0] + assert(math.equal(math.multiply(root, root, root), quat)) + }) }) diff --git a/src/complex/__test__/type.spec.js b/src/complex/__test__/type.spec.js index fa19704..dc4fdaa 100644 --- a/src/complex/__test__/type.spec.js +++ b/src/complex/__test__/type.spec.js @@ -19,6 +19,13 @@ describe('complex type operations', () => { assert.strictEqual(math.arg(cplx(1, Math.sqrt(3))), Math.PI/3) assert.strictEqual(math.arg(cplx(true, true)), Math.PI/4) }) + it('calculates the real parts of complex numbers and quaternions', () => { + assert.strictEqual(math.re(cplx(2, -3)), 2) + assert.strictEqual(math.re(cplx(cplx(0.5), cplx(8, -7))), 0.5) + }) + it('calculates the argument of a quaternion', () => { + assert(math.equal(math.arg(cplx(cplx(1, 1), cplx(1, 1))), Math.acos(0.5))) + }) it('detects associates of a complex number', () => { const z = cplx(3, 4) const assoc = math.associate diff --git a/src/complex/arithmetic.js b/src/complex/arithmetic.js index 045c0f9..3a77a11 100644 --- a/src/complex/arithmetic.js +++ b/src/complex/arithmetic.js @@ -4,6 +4,7 @@ import {ResolutionError} from '#core/helpers.js' import {Returns, ReturnType, ReturnTyping} from '#core/Type.js' import {match} from '#core/TypePatterns.js' import {ReturnsAs} from '#generic/helpers.js' +import {NumberT} from '#number/NumberT.js' const {conservative, full, free} = ReturnTyping @@ -141,4 +142,44 @@ export const sqrt = match(Complex, (math, C, strategy) => { return ReturnsAs(prune, z => prune(sqrtImp(z))) }) +const TAU3 = 2 * Math.PI / 3 + +export const cbrt = match(Complex, (math, C) => { + const arg = math.arg.resolve(C) + const divArg = math.divide.resolve([ReturnType(arg), NumberT]) + const absC = math.abs.resolve(C) + const cbrtR = math.cbrt.resolve(ReturnType(absC), conservative) + const im = math.im.resolve(C) + const absIm = math.abs.resolve(ReturnType(im)) + const divIm = math.divide.resolve([ReturnType(im), ReturnType(absIm)]) + // TODO: replace with nanomath cos and sin when available + // const cos = math.cos.resolve(ReturnType(divArg)) + // const CosType = ReturnType(cos) // and similarly for sin + const cos = Math.cos + const CosType = NumberT + const sin = Math.sin + const SinType = NumberT + const mulIm = math.multiply.resolve([ReturnType(divIm), SinType]) + const addUp = math.add.resolve([CosType, ReturnType(mulIm)]) + const addAngle = math.add.resolve([ReturnType(divArg), math.typeOf(TAU3)]) + const subAngle = math.subtract.resolve( + [ReturnType(divArg), math.typeOf(TAU3)]) + const scale = math.multiply.resolve([ReturnType(cbrtR), ReturnType(addUp)]) + const Cout = ReturnType(scale) + const vec = math.vector.resolve([Cout, Cout, Cout]) + return ReturnsAs(vec, z => { + const arg3 = divArg(arg(z), 3) + const mag = absC(z) + const r = cbrtR(mag) + const imz = im(z) + const unit = divIm(imz, absIm(imz)) + // At this point, z = mag*(cos(arg) + unit * sin(arg)) + // so one cube root is r*(cos(arg3) + unit * sin(arg3)) + // and we get the other two by adjusting arg3 by +- TAU3 + const cus = theta => scale(r, addUp(cos(theta), mulIm(unit, sin(theta)))) + return vec( + cus(arg3), cus(addAngle(arg3, TAU3)), cus(subAngle(arg3, TAU3))) + }) +}) + export const subtract = promoteBinary('subtract') diff --git a/src/complex/type.js b/src/complex/type.js index b3586e4..a36094f 100644 --- a/src/complex/type.js +++ b/src/complex/type.js @@ -31,20 +31,26 @@ export const complex = [ }) ] -export const arg = // [ // enable when we have atan2 in mathjs -// match(Complex, (math, C) => { -// const re = math.re.resolve(C) -// const R = ReturnType(re) -// const im = math.im.resolve(C) -// const abs = math.abs.resolve(C) -// const atan2 = math.atan2.resolve([R, R], conservative) -// return Returns(R, z => atan2(abs(im(z)), re(z))) -// }), // note always between 0 and tau/2; need to use in conjunction -// // with a complex unit function that gives you the proper -// // imaginary unit, ±i in the simple complex case, to restore the -// // full circle of values for the direction of a complex number - match(Complex(NumberT), Returns(NumberT, z => Math.atan2(z.im, z.re))) -//] +export const arg = match(Complex, (math, C) => { + const re = math.re.resolve(C) + const atan2 = Math.atan2 + const ArgType = NumberT + if (ReturnType(re) === C.Component) { + // "plain" Complex (not nested, i.e. not quaternions) + // FIXME: use math.atan2 once available + // const atan2 = math.atan2.resolve( + // [C.Component, C.Component], conservative) + // const argType = ReturnType(atan2) + return Returns(ArgType, z => atan2(z.im, z.re)) + } + const im = math.im.resolve(C) + const abs = math.abs.resolve(ReturnType(im)) + // FIXME: use math.atan2 once available + // const atan2 = math.atan2.resolve( + // [ReturnType(abs), ReturnType(re)], conservative) + // const ArgType = ReturnType(atan2) + return Returns(ArgType, z => atan2(abs(im(z)), re(z))) +}) /* Returns true if w is z multiplied by a complex unit */ export const associate = match(