feat(polynomialRoot) (#57)

Implements a simply polynomial root finder function
polynomialRoot, intended to be used for benchmarking
against mathjs.

For this purpose, adds numerous other functions (e.g.
cbrt, arg, cis), refactors sqrt (so that you can
definitely get the complex square root when you want
it), and makes numerous enhancements to the core so
that a template can match after conversions.

Co-authored-by: Glen Whitney <glen@studioinfinity.org>
Reviewed-on: #57

README.md

@@ -25,3 +25,7 @@ Hopefully this shows promise. It is an evolution of the concept first prototyped
 Note that Pocomath allows one implementation to depend just on a specific signature of another function, for efficiency's sake (if for example 'bar(Matrix)' knows it will only call 'foo(Matrix)', it avoids another type-dispatch). That capability is used in sqrt, for example.
 
 Pocomath also lazily reloads operations that depend on the config when that changes, and if an operation has a signature mentioning an undefined type, that signature is ignored until the type is installed, at which point the function lazily redefines itself to use the additional signature.
+
+Pocomath now also allows template operations and template types, also built on top of typed-function (but candidates for integration therein). This is used to make many operations more specific, implement a type-homogeneous Tuple type, and make Complex numbers be type-homogeneous (which it seems like it always should be). One of the cutest consequences of this approach is that with careful definitions of the `Complex<T>` templates, one gets a working quaternion data type absolutely for free as `Complex<Complex<number>>` (and integral quaternions as `Complex<Complex<bigint>>`, etc.)
+
+It also now has a facility to adapt a third-party numeric class as a type in Pocomath, see `src/generic/all.mjs` and `src/generic/Types/adapted.mjs`, which it uses by way of example to incorporate fraction.js Fraction objects into Pocomath.

package.json5

@@ -24,6 +24,7 @@
   },
   dependencies: {
     'bigint-isqrt': '^0.2.1',
+    'fraction.js': '^4.2.0',
     'typed-function': '^3.0.0',
   },
 }

pnpm-lock.yaml

@@ -2,11 +2,13 @@ lockfileVersion: 5.4
 
 specifiers:
   bigint-isqrt: ^0.2.1
+  fraction.js: ^4.2.0
   mocha: ^10.0.0
   typed-function: ^3.0.0
 
 dependencies:
   bigint-isqrt: 0.2.1
+  fraction.js: 4.2.0
   typed-function: 3.0.0
 
 devDependencies:
@@ -192,6 +194,10 @@ packages:
     hasBin: true
     dev: true
 
+  /fraction.js/4.2.0:
+    resolution: {integrity: sha512-MhLuK+2gUcnZe8ZHlaaINnQLl0xRIGRfcGk2yl8xoQAfHrSsL3rYu6FCmBdkdbhc9EPlwyGHewaRsvwRMJtAlA==}
+    dev: false
+
   /fs.realpath/1.0.0:
     resolution: {integrity: sha512-OO0pH2lK6a0hZnAdau5ItzHPI6pUlvI7jMVnxUQRtw4owF2wk8lOSabtGDCTP4Ggrg2MbGnWO9X8K1t4+fGMDw==}
     dev: true

src/bigint/absquare.mjs

@@ -0,0 +1,7 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/bigint.mjs'
+
+/* Absolute value squared */
+export const absquare = {
+    bigint: ({'square(bigint)': sqb}) => Returns('bigint', b => sqb(b))
+}

src/bigint/add.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const add = {'bigint,bigint': () => (a,b) => a+b}
+export const add = {'bigint,bigint': () => Returns('bigint', (a,b) => a+b)}

src/bigint/compare.mjs

@@ -1,5 +1,7 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
 export const compare = {
-   'bigint,bigint': () => (a,b) => a === b ? 0n : (a > b ? 1n : -1n)
+   'bigint,bigint': () => Returns(
+      'boolean', (a,b) => a === b ? 0n : (a > b ? 1n : -1n))
 }

src/bigint/divide.mjs

@@ -1,12 +1,13 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
 export const divide = {
    'bigint,bigint': ({config, 'quotient(bigint,bigint)': quot}) => {
-      if (config.predictable) return quot
-      return (n, d) => {
+      if (config.predictable) return Returns('bigint', (n,d) => quot(n,d))
+      return Returns('bigint|undefined', (n, d) => {
          const q = n/d
          if (q * d == n) return q
          return undefined
-      }
+      })
    }
 }

src/bigint/isZero.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const isZero = {bigint: () => b => b === 0n}
+export const isZero = {bigint: () => Returns('boolean', b => b === 0n)}

src/bigint/multiply.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const multiply = {'bigint,bigint': () => (a,b) => a*b}
+export const multiply = {'bigint,bigint': () => Returns('bigint', (a,b) => a*b)}

src/bigint/native.mjs

@@ -1,11 +1,12 @@
 import gcdType from '../generic/gcdType.mjs'
-import {identity} from '../generic/identity.mjs'
+import {identityType} from '../generic/identity.mjs'
 
 export * from './Types/bigint.mjs'
 
+export {absquare} from './absquare.mjs'
 export {add} from './add.mjs'
 export {compare} from './compare.mjs'
-export const conjugate = {bigint: () => identity}
+export const conjugate = {bigint: identityType('bigint')}
 export {divide} from './divide.mjs'
 export const gcd = gcdType('bigint')
 export {isZero} from './isZero.mjs'

src/bigint/negate.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const negate = {bigint: () => b => -b}
+export const negate = {bigint: () => Returns('bigint', b => -b)}

src/bigint/one.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const one = {bigint: () => () => 1n}
+export const one = {bigint: () => Returns('bigint', () => 1n)}

src/bigint/quotient.mjs

@@ -1,13 +1,14 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-/* Returns the best integer approximation to n/d */
+/* Returns the floor integer approximation to n/d */
 export const quotient = {
-   'bigint,bigint': ({'sign(bigint)': sgn}) => (n, d) => {
+   'bigint,bigint': ({'sign(bigint)': sgn}) => Returns('bigint', (n, d) => {
       const dSgn = sgn(d)
       if (dSgn === 0n) return 0n
       if (sgn(n) === dSgn) return n/d
       const quot = n/d
       if (quot * d == n) return quot
       return quot - 1n
-   }
+   })
 }

src/bigint/roundquotient.mjs

@@ -1,8 +1,9 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
 /* Returns the closest integer approximation to n/d */
 export const roundquotient = {
-   'bigint,bigint': ({'sign(bigint)': sgn}) => (n, d) => {
+   'bigint,bigint': ({'sign(bigint)': sgn}) => Returns('bigint', (n, d) => {
       const dSgn = sgn(d)
       if (dSgn === 0n) return 0n
       const candidate = n/d
@@ -11,5 +12,5 @@ export const roundquotient = {
       if (2n * rem > absd) return candidate + dSgn
       if (-2n * rem >= absd) return candidate - dSgn
       return candidate
-   }
+   })
 }

src/bigint/sign.mjs

@@ -1,9 +1,10 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
 export const sign = {
-   bigint: () => b => {
+   bigint: () => Returns('bigint', b => {
       if (b === 0n) return 0n
       if (b > 0n) return 1n
       return -1n
-   }
+   })
 }

src/bigint/sqrt.mjs

@@ -1,30 +1,39 @@
-export * from './Types/bigint.mjs'
+import Returns from '../core/Returns.mjs'
 import isqrt from 'bigint-isqrt'
+export * from './Types/bigint.mjs'
 
 export const sqrt = {
-   bigint: ({config, complex, 'self(Complex)': complexSqrt}) => {
+   bigint: ({
+      config,
+      'complex(bigint,bigint)': cplx,
+      'negate(bigint)': neg
+   }) => {
       if (config.predictable) {
          // Don't just return the constant isqrt here because the object
          // gets decorated with info that might need to be different
          // for different PocomathInstancss
-         return b => isqrt(b)
+         return Returns('bigint', b => isqrt(b))
       }
-      if (!complexSqrt) {
-         return b => {
+      if (!cplx) {
+         return Returns('bigint|undefined', b => {
             if (b >= 0n) {
                const trial = isqrt(b)
                if (trial * trial === b) return trial
             }
             return undefined
-         }
+         })
       }
-      return b => {
-         if (b >= 0n) {
-            const trial = isqrt(b)
-            if (trial * trial === b) return trial
-            return undefined
+      return Returns('bigint|Complex<bigint>|undefined', b => {
+         if (b === undefined) return undefined
+         let real = true
+         if (b < 0n) {
+            b = neg(b)
+            real = false
          }
-         return complexSqrt(complex(b))
-      }
+         const trial = isqrt(b)
+         if (trial * trial !== b) return undefined
+         if (real) return trial
+         return cplx(0n, trial)
+      })
    }
 }

src/bigint/zero.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/bigint.mjs'
 
-export const zero = {bigint: () => () => 0n}
+export const zero = {bigint: () => Returns('bigint', () => 0n)}

src/complex/Types/Complex.mjs

@@ -1,30 +1,31 @@
+import {Returns, returnTypeOf} from '../../core/Returns.mjs'
 import PocomathInstance from '../../core/PocomathInstance.mjs'
 
-/* Use a plain object with keys re and im for a complex; note the components
- * can be any type (for this proof-of-concept; in reality we'd want to
- * insist on some numeric or scalar supertype).
- */
-function isComplex(z) { 
-   return z && typeof z === 'object' && 're' in z && 'im' in z
-}
-
 const Complex = new PocomathInstance('Complex')
-Complex.installType('Complex', {
-   test: isComplex,
+// Now the template type: Complex numbers are actually always homogeneous
+// in their component types. For an explanation of the meanings of the
+// properties, see ../../tuple/Types/Tuple.mjs
+Complex.installType('Complex<T>', {
+   base: z => z && typeof z === 'object' && 're' in z && 'im' in z,
+   test: testT => z => testT(z.re) && testT(z.im),
+   infer: ({typeOf, joinTypes}) => z => joinTypes([typeOf(z.re), typeOf(z.im)]),
    from: {
-      number: x => ({re: x, im: 0})
-   }
-})
-Complex.installType('GaussianInteger', {
-   test: z => typeof z.re == 'bigint' && typeof z.im == 'bigint',
-   refines: 'Complex',
-   from: {
-      bigint: x => ({re: x, im: 0n})
+      T: t => ({re: t, im: t-t}), // hack: maybe need a way to call zero(T)
+      U: convert => u => {
+         const t = convert(u)
+         return ({re: t, im: t-t})
+      },
+      'Complex<U>': convert => cu => ({re: convert(cu.re), im: convert(cu.im)})
    }
 })
 
 Complex.promoteUnary = {
-   Complex: ({self,complex}) => z => complex(self(z.re), self(z.im))
+   'Complex<T>': ({
+      T,
+      'self(T)': me,
+      complex
+   }) => Returns(
+      `Complex<${returnTypeOf(me)}>`, z => complex(me(z.re), me(z.im)))
 }
 
 export {Complex}

src/complex/abs.mjs

@@ -1,5 +1,20 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const abs = {
-    Complex: ({sqrt, 'absquare(Complex)': absq}) => z => sqrt(absq(z))
+   'Complex<T>': ({
+      sqrt, // Unfortunately no notation yet for the needed signature
+      'absquare(T)': baseabsq,
+      'absquare(Complex<T>)': absq
+   }) => {
+      const midType = returnTypeOf(baseabsq)
+      const sqrtImp = sqrt.fromInstance.resolve('sqrt', midType, sqrt)
+      let retType = returnTypeOf(sqrtImp)
+      if (retType.includes('|')) {
+         // This is a bit of a hack, as it relies on all implementations of
+         // sqrt returning the "typical" return type as the first option
+         retType = retType.split('|',1)[0]
+      }
+      return Returns(retType, z => sqrtImp(absq(z)))
+   }
 }

src/complex/absquare.mjs

@@ -1,5 +1,27 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const absquare = {
-   Complex: ({add, square}) => z => add(square(z.re), square(z.im))
+   'Complex<T>': ({
+      add, // no easy way to write the needed signature; if T is number
+           // it is number,number; but if T is Complex<bigint>, it is just
+           // bigint,bigint. So unfortunately we depend on all of add, and
+           // we extract the needed implementation below.
+      'self(T)': absq
+   }) => {
+      const midType = returnTypeOf(absq)
+      const addImp = add.fromInstance.resolve(
+         'add', `${midType},${midType}`, add)
+      return Returns(
+         returnTypeOf(addImp), z => addImp(absq(z.re), absq(z.im)))
+   }
 }
+
+/* We could imagine notations that Pocomath could support that would simplify
+ * the above, maybe something like
+ * 'Complex<T>': ({
+ *    'self(T): U': absq,
+ *    'add(U,U):V': plus,
+ *    V
+ *  }) => Returns(V, z => plus(absq(z.re), absq(z.im)))
+ */

src/complex/add.mjs

@@ -1,22 +1,10 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const add = {
-   /* Relying on conversions for both complex + number and complex + bigint
-    * leads to an infinite loop when adding a number and a bigint, since they
-    * both convert to Complex.
-    */
-   'Complex,number': ({
-      'self(number,number)': addNum,
-      'complex(number,number)': cplx
-   }) => (z,x) => cplx(addNum(z.re, x), z.im),
-
-   'Complex,bigint': ({
-      'self(bigint,bigint)': addBigInt,
-      'complex(bigint,bigint)': cplx
-   }) => (z,x) => cplx(addBigInt(z.re, x), z.im),
-
-   'Complex,Complex': ({
-      self,
-      complex
-   }) => (w,z) => complex(self(w.re, z.re), self(w.im, z.im))
+   'Complex<T>,Complex<T>': ({
+      T,
+      'self(T,T)': me,
+      'complex(T,T)': cplx
+   }) => Returns(`Complex<${T}>`, (w,z) => cplx(me(w.re, z.re), me(w.im, z.im)))
 }

src/complex/arg.mjs

@@ -0,0 +1,7 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+/* arg is the "argument" or angle theta of z in its form r cis theta */
+export const arg = {
+   'Complex<number>': () => Returns('number', z => Math.atan2(z.im, z.re))
+}

src/complex/associate.mjs

@@ -1,17 +1,18 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 /* Returns true if w is z multiplied by a complex unit */
 export const associate = {
-    'Complex,Complex': ({
-        'multiply(Complex,Complex)': times,
-        'equalTT(Complex,Complex)': eq,
-        zero,
-        one,
-        complex,
-        'negate(Complex)': neg
-    }) => (w,z) => {
+    'Complex<T>,Complex<T>': ({
+        'multiply(Complex<T>,Complex<T>)': times,
+        'equalTT(Complex<T>,Complex<T>)': eq,
+        'zero(T)': zr,
+        'one(T)': uno,
+        'complex(T,T)': cplx,
+        'negate(Complex<T>)': neg
+    }) => Returns('boolean', (w,z) => {
         if (eq(w,z) || eq(w,neg(z))) return true
-        const ti = times(z, complex(zero(z.re), one(z.im)))
+        const ti = times(z, cplx(zr(z.re), uno(z.im)))
         return eq(w,ti) || eq(w,neg(ti))
-    }
+    })
 }

src/complex/cbrtc.mjs

@@ -0,0 +1,28 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+const TAU3 = 2 * Math.PI / 3
+
+/* Complex cube root that returns all three roots as a tuple of complex. */
+/* follows the implementation in mathjs */
+/* Really only works for T = number at the moment because of arg and cbrt */
+export const cbrtc = {
+   'Complex<T>': ({
+      'arg(T)': theta,
+      'divide(T,T)': div,
+      'abs(Complex<T>)': absval,
+      'complex(T)': cplx,
+      'cbrt(T)': cbrtT,
+      'multiply(Complex<T>,Complex<T>)': mult,
+      'cis(T)': cisT,
+      'tuple(...Complex<T>)': tup
+   }) => Returns('Tuple<Complex<T>>', z => {
+      const arg3 = div(theta(z), 3)
+      const r = cplx(cbrtT(absval(z)))
+      return tup(
+         mult(r, cisT(arg3)),
+         mult(r, cisT(arg3 + TAU3)),
+         mult(r, cisT(arg3 - TAU3))
+      )
+   })
+}

src/complex/cis.mjs

@@ -0,0 +1,9 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+/* Returns cosine plus i sin theta */
+export const cis = {
+   'number': ({'complex(number,number)': cplx}) => Returns(
+      'Complex<number>', t => cplx(Math.cos(t), Math.sin(t))
+   )
+}

src/complex/complex.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 export * from '../generic/Types/generic.mjs'
 
@@ -6,11 +7,16 @@ export const complex = {
     * have a numeric/scalar type, e.g. by implementing subtypes in
     * typed-function
     */
-   'undefined': () => u => u,
-   'undefined,any': () => (u, y) => u,
-   'any,undefined': () => (x, u) => u,
-   'undefined,undefined': () => (u, v) => u,
-   'T,T': () => (x, y) => ({re: x, im: y}),
+   'undefined': () => Returns('undefined', u => u),
+   'undefined,any': () => Returns('undefined', (u, y) => u),
+   'any,undefined': () => Returns('undefined', (x, u) => u),
+   'undefined,undefined': () => Returns('undefined', (u, v) => u),
+   'T,T': ({T}) => Returns(`Complex<${T}>`, (x, y) => ({re: x, im: y})),
    /* Take advantage of conversions in typed-function */
-   Complex: () => z => z
+   // 'Complex<T>': () => z => z
+   /* But help out because without templates built in to typed-function,
+    * type inference turns out to be too hard
+    */
+   'T': ({T, 'zero(T)': zr}) => Returns(
+      `Complex<${T}>`, x => ({re: x, im: zr(x)}))
 }

src/complex/conjugate.mjs

@@ -1,6 +1,11 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const conjugate = {
-   Complex: ({negate, complex}) => z => complex(z.re, negate(z.im))
+   'Complex<T>': ({
+      T,
+      'negate(T)': neg,
+      'complex(T,T)': cplx
+   }) => Returns(`Complex<${T}>`, z => cplx(z.re, neg(z.im)))
 }
 

src/complex/equalTT.mjs

@@ -1,19 +1,28 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const equalTT = {
-   'Complex,number': ({
-      'isZero(number)': isZ,
-      'self(number,number)': eqNum
-   }) => (z, x) => eqNum(z.re, x) && isZ(z.im),
+   'Complex<T>,Complex<T>': ({
+      T,
+      'self(T,T)': me
+   }) => Returns('boolean', (w, z) => me(w.re, z.re) && me(w.im, z.im)),
+   // NOTE: Although I do not understand exactly why, with typed-function@3.0's
+   // matching algorithm, the above template must come first to ensure the
+   // most specific match to a template call. I.e, if one of the below
+   // comes first, a call with two complex numbers can match via conversions
+   // with (Complex<Complex<number>>, Complex<number>) (!, hopefully in some
+   // future iteration typed-function will be smart enough to prefer
+   // Complex<T>, Complex<T>. Possibly the problem is in Pocomath's bolted-on
+   // type resolution and the difficulty will go away when features are moved
+   // into typed-function.
+   'Complex<T>,T': ({
+      'isZero(T)': isZ,
+      'self(T,T)': eqReal
+   }) => Returns('boolean', (z, x) => eqReal(z.re, x) && isZ(z.im)),
 
-   'Complex,bigint': ({
-      'isZero(bigint)': isZ,
-      'self(bigint,bigint)': eqBigInt
-   }) => (z, b) => eqBigInt(z.re, b) && isZ(z.im),
+   'T,Complex<T>': ({
+      'isZero(T)': isZ,
+      'self(T,T)': eqReal
+   }) => Returns('boolean', (b, z) => eqReal(z.re, b) && isZ(z.im)),
 
-   'Complex,Complex': ({self}) => (w,z) => self(w.re, z.re) && self(w.im, z.im),
-
-   'GaussianInteger,GaussianInteger': ({
-      'self(bigint,bigint)': eq
-   }) => (a,b) => eq(a.re, b.re) && eq(a.im, b.im)
 }

src/complex/extendToComplex.mjs

@@ -15,4 +15,17 @@ export default async function extendToComplex(pmath) {
          // Guess it wasn't a method available in complex; no worries
       }
    }
+   // Since extension to complex was specifically requested, instantiate
+   // all of the templates so that the associated type conversions will
+   // be available to make function calls work immediately:
+   for (const baseType in pmath.Types) {
+      if (baseType in pmath.Templates || baseType.includes('<')) {
+         continue // don't mess with templates
+      }
+      const ignore = new Set(['undefined', 'any', 'T', 'ground'])
+      if (ignore.has(baseType)) continue
+      // (What we really want is a check for "numeric" types but we don't
+      //  have that concept (yet?)). If we did, we'd instantiate just for those...
+      pmath.instantiateTemplate('Complex', baseType)
+   }
 }

src/complex/gcd.mjs

@@ -1,19 +1,25 @@
 import PocomathInstance from '../core/PocomathInstance.mjs'
-import * as Complex from './Types/Complex.mjs'
+import Returns from '../core/Returns.mjs'
+import * as  Complex from './Types/Complex.mjs'
 import gcdType from '../generic/gcdType.mjs'
 
+const gcdComplexRaw = {}
+Object.assign(gcdComplexRaw, gcdType('Complex<bigint>'))
+Object.assign(gcdComplexRaw, gcdType('Complex<NumInt>'))
 const imps = {
-   gcdGIRaw: gcdType('GaussianInteger'),
+   gcdComplexRaw,
    gcd: { // Only return gcds with positive real part
-      'GaussianInteger,GaussianInteger': ({
-         'gcdGIRaw(GaussianInteger,GaussianInteger)': gcdRaw,
-         'sign(bigint)': sgn,
-         'negate(GaussianInteger)': neg
-      }) => (z,m) => {
+      'Complex<T>,Complex<T>': ({
+         T,
+         'gcdComplexRaw(Complex<T>,Complex<T>)': gcdRaw,
+         'sign(T)': sgn,
+         'one(T)': uno,
+         'negate(Complex<T>)': neg
+      }) => Returns(`Complex<${T}>`, (z,m) => {
          const raw = gcdRaw(z, m)
-         if (sgn(raw.re) === 1n) return raw
+         if (sgn(raw.re) === uno(raw.re)) return raw
          return neg(raw)
-      }
+      })
    }
 }
 

src/complex/invert.mjs

@@ -1,9 +1,16 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const invert = {
-   Complex: ({conjugate, absquare, complex, divide}) => z => {
-      const c = conjugate(z)
-      const d = absquare(z)
-      return complex(divide(c.re, d), divide(c.im, d))
-   }
+   'Complex<T>': ({
+      T,
+      'conjugate(Complex<T>)': conj,
+      'absquare(Complex<T>)': asq,
+      'complex(T,T)': cplx,
+      'divide(T,T)': div
+   }) => Returns(`Complex<${T}>`, z => {
+      const c = conj(z)
+      const d = asq(z)
+      return cplx(div(c.re, d), div(c.im, d))
+   })
 }

src/complex/isReal.mjs

@@ -0,0 +1,7 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+export const isReal = {
+   'Complex<T>': ({'equal(T,T)': eq, 'add(T,T)': plus}) => Returns(
+      'boolean', z => eq(z.re, plus(z.re, z.im)))
+}

src/complex/isZero.mjs

@@ -1,5 +1,7 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const isZero = {
-   Complex: ({self}) => z => self(z.re) && self(z.im)
+   'Complex<T>': ({'self(T)': me}) => Returns(
+      'boolean', z => me(z.re) && me(z.im))
 }

src/complex/multiply.mjs

@@ -1,14 +1,20 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const multiply = {
-   'Complex,Complex': ({
-      'complex(any,any)': cplx,
-      add,
-      subtract,
-      self
-   }) => (w,z) => {
-      return cplx(
-         subtract(self(w.re, z.re), self(w.im, z.im)),
-         add(self(w.re, z.im), self(w.im, z.re)))
-   }
+   'Complex<T>,Complex<T>': ({
+      T,
+      'complex(T,T)': cplx,
+      'add(T,T)': plus,
+      'subtract(T,T)': subt,
+      'self(T,T)': me,
+      'conjugate(T)': conj // makes quaternion multiplication work
+   }) => Returns(
+      `Complex<${T}>`,
+      (w,z) => {
+         const realpart = subt(me(     w.re,  z.re), me(conj(w.im), z.im))
+         const imagpart = plus(me(conj(w.re), z.im), me(     w.im,  z.re))
+         return cplx(realpart, imagpart)
+      }
+   )
 }

src/complex/native.mjs

@@ -3,17 +3,24 @@ export * from './Types/Complex.mjs'
 export {abs} from './abs.mjs'
 export {absquare} from './absquare.mjs'
 export {add} from './add.mjs'
+export {arg} from './arg.mjs'
 export {associate} from './associate.mjs'
+export {cbrtc} from './cbrtc.mjs'
+export {cis} from './cis.mjs'
 export {complex} from './complex.mjs'
 export {conjugate} from './conjugate.mjs'
 export {equalTT} from './equalTT.mjs'
 export {gcd} from './gcd.mjs'
 export {invert} from './invert.mjs'
+export {isReal} from './isReal.mjs'
 export {isZero} from './isZero.mjs'
 export {multiply} from './multiply.mjs'
 export {negate} from './negate.mjs'
+export {polynomialRoot} from './polynomialRoot.mjs'
+export {quaternion} from './quaternion.mjs'
 export {quotient} from './quotient.mjs'
 export {roundquotient} from './roundquotient.mjs'
 export {sqrt} from './sqrt.mjs'
+export {sqrtc} from './sqrtc.mjs'
 export {zero} from './zero.mjs'
 

src/complex/polynomialRoot.mjs

@@ -0,0 +1,118 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+export const polynomialRoot = {
+   'Complex<T>,...Complex<T>': ({
+      T,
+      'tuple(...Complex<T>)': tupCplx,
+      'tuple(...T)': tupReal,
+      'isZero(Complex<T>)': zero,
+      'complex(T)': C,
+      'multiply(Complex<T>,Complex<T>)': mul,
+      'divide(Complex<T>,Complex<T>)': div,
+      'negate(Complex<T>)': neg,
+      'isReal(Complex<T>)': real,
+      'equalTT(Complex<T>,Complex<T>)': eq,
+      'add(Complex<T>,Complex<T>)': plus,
+      'subtract(Complex<T>,Complex<T>)': sub,
+      'sqrtc(Complex<T>)': sqt,
+      'cbrtc(Complex<T>)': cbt
+   }) => Returns(`Tuple<${T}>|Tuple<Complex<${T}>>`, (constant, rest) => {
+      // helper to convert results to appropriate tuple type
+      const typedTup = arr => {
+         if (arr.every(real)) {
+            return tupReal.apply(tupReal, arr.map(z => z.re))
+         }
+         return tupCplx.apply(tupCplx, arr)
+      }
+
+      const coeffs = [constant, ...rest]
+      while (coeffs.length > 0 && zero(coeffs[coeffs.length - 1])) {
+         coeffs.pop()
+      }
+      if (coeffs.length < 2) {
+      }
+      switch (coeffs.length) {
+         case 0: case 1:
+            throw new RangeError(
+               `Polynomial [${constant}, ${rest}] must have at least one`
+                  + 'non-zero non-constant coefficient')
+         case 2: // linear
+            return typedTup([neg(div(coeffs[0], coeffs[1]))])
+         case 3: { // quadratic
+            const [c, b, a] = coeffs
+            const denom = mul(C(2), a)
+            const d1 = mul(b, b)
+            const d2 = mul(C(4), mul(a, c))
+            if (eq(d1, d2)) {
+               return typedTup([div(neg(b), denom)])
+            }
+            let discriminant = sqt(sub(d1, d2))
+            return typedTup([
+               div(sub(discriminant, b), denom),
+               div(sub(neg(discriminant), b), denom)
+            ])
+         }
+         case 4: { // cubic, cf. https://en.wikipedia.org/wiki/Cubic_equation
+            const [d, c, b, a] = coeffs
+            const denom = neg(mul(C(3), a))
+            const asqrd = mul(a, a)
+            const D0_1 = mul(b, b)
+            const bcubed = mul(D0_1, b)
+            const D0_2 = mul(C(3), mul(a, c))
+            const D1_1 = plus(
+               mul(C(2), bcubed), mul(C(27), mul(asqrd, d)))
+            const abc = mul(a, mul(b, c))
+            const D1_2 = mul(C(9), abc)
+            // Check for a triple root
+            if (eq(D0_1, D0_2) && eq(D1_1, D1_2)) {
+               return typedTup([div(b, denom)])
+            }
+            const Delta0 = sub(D0_1, D0_2)
+            const Delta1 = sub(D1_1, D1_2)
+            const csqrd = mul(c, c)
+            const discriminant1 = plus(
+               mul(C(18), mul(abc, d)), mul(D0_1, csqrd))
+            const discriminant2 = plus(
+               mul(C(4), mul(bcubed, d)),
+               plus(
+                  mul(C(4), mul(a, mul(csqrd, c))),
+                  mul(C(27), mul(asqrd, mul(d, d)))))
+            // See if we have a double root
+            if (eq(discriminant1, discriminant2)) {
+               return typedTup([
+                  div(
+                     sub(
+                        mul(C(4), abc),
+                        plus(mul(C(9), mul(asqrd, d)), bcubed)),
+                     mul(a, Delta0)), // simple root
+                  div(
+                     sub(mul(C(9), mul(a, d)), mul(b, c)),
+                     mul(C(2), Delta0)) // double root
+               ])
+            }
+            // OK, we have three distinct roots
+            let Ccubed
+            if (eq(D0_1, D0_2)) {
+               Ccubed = Delta1
+            } else {
+               Ccubed = div(
+                  plus(
+                     Delta1,
+                     sqt(sub(
+                        mul(Delta1, Delta1),
+                        mul(C(4), mul(Delta0, mul(Delta0, Delta0)))))
+                  ),
+                  C(2))
+            }
+            const croots = cbt(Ccubed)
+            return typedTup(cbt(Ccubed).elts.map(
+               C => div(plus(b, plus(C, div(Delta0, C))), denom)))
+         }
+         default:
+            throw new RangeError(
+               'only implemented for cubic or lower-order polynomials, '
+                  + `not ${JSON.stringify(coeffs)}`)
+      }
+   })
+}

src/complex/quaternion.mjs

@@ -0,0 +1,14 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+// Might be nice to have type aliases!
+export const quaternion = {
+    'T,T,T,T': ({
+        T,
+        'complex(T,T)': cplxT,
+        'complex(Complex<T>,Complex<T>)': quat
+    }) => Returns(
+        `Complex<Complex<${T}>>`,
+        (r,i,j,k) => quat(cplxT(r,j), cplxT(i,k))
+    )
+}

src/complex/quotient.mjs

@@ -1,5 +1,9 @@
+import Returns from '../core/Returns.mjs'
 export * from './roundquotient.mjs'
 
 export const quotient = {
-   'Complex,Complex': ({roundquotient}) => (w,z) => roundquotient(w,z)
+   'Complex<T>,Complex<T>': ({
+      T,
+      'roundquotient(Complex<T>,Complex<T>)': rq
+   }) => Returns(`Complex<${T}>`, (w,z) => rq(w,z))
 }

src/complex/roundquotient.mjs

@@ -1,17 +1,19 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const roundquotient = {
-   'Complex,Complex': ({
-      'isZero(Complex)': isZ,
-      conjugate,
-      'multiply(Complex,Complex)': mult,
-      absquare,
-      self,
-      complex
-   }) => (n,d) => {
+   'Complex<T>,Complex<T>': ({
+      T,
+      'isZero(Complex<T>)': isZ,
+      'conjugate(Complex<T>)': conj,
+      'multiply(Complex<T>,Complex<T>)': mult,
+      'absquare(Complex<T>)': asq,
+      'self(T,T)': me,
+      'complex(T,T)': cplx
+   }) => Returns(`Complex<${T}>`, (n,d) => {
       if (isZ(d)) return d
-      const cnum = mult(n, conjugate(d))
-      const dreal = absquare(d)
-      return complex(self(cnum.re, dreal), self(cnum.im, dreal))
-   }
+      const cnum = mult(n, conj(d))
+      const dreal = asq(d)
+      return cplx(me(cnum.re, dreal), me(cnum.im, dreal))
+   })
 }

src/complex/sqrt.mjs

@@ -1,39 +1,34 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
 export * from './Types/Complex.mjs'
 
 export const sqrt = {
-   Complex: ({
+   'Complex<T>': ({
       config,
-      isZero,
-      sign,
-      one,
-      add,
-      complex,
-      multiply,
-      self,
-      divide,
-      'abs(Complex)': abs,
-      subtract
+      'sqrtc(Complex<T>)': predictableSqrt,
+      'isZero(T)': isZ,
    }) => {
+      if (config.checkingDependency) return undefined
+      const complexReturns = returnTypeOf(predictableSqrt)
+      const baseReturns = complexReturns.slice(8, -1); // Complex<WhatWeWant>
       if (config.predictable) {
-         return z => {
-            const reOne = one(z.re)
-            if (isZero(z.im) && sign(z.re) === reOne) return complex(self(z.re))
-            const reTwo = add(reOne, reOne)
-            return complex(
-               multiply(sign(z.im), self(divide(add(abs(z),z.re), reTwo))),
-               self(divide(subtract(abs(z),z.re), reTwo))
-            )
+         return Returns(complexReturns, z => predictableSqrt(z))
+      }
+
+      return Returns(
+         `Complex<${baseReturns}>|${baseReturns}|undefined`,
+         z => {
+            let complexSqrt
+            try {
+               complexSqrt = predictableSqrt(z)
+            } catch (e) {
+               return undefined
+            }
+            if (complexSqrt.re === undefined || complexSqrt.im === undefined) {
+               return undefined
+            }
+            if (isZ(complexSqrt.im)) return complexSqrt.re
+            return complexSqrt
          }
-      }
-      return z => {
-         const reOne = one(z.re)
-         if (isZero(z.im) && sign(z.re) === reOne) return self(z.re)
-         const reTwo = add(reOne, reOne)
-         return complex(
-            multiply(sign(z.im), self(divide(add(abs(z),z.re), reTwo))),
-            self(divide(subtract(abs(z),z.re), reTwo))
-         )
-      }
+      )
    }
 }
-

src/complex/sqrtc.mjs

@@ -0,0 +1,41 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
+export * from './Types/Complex.mjs'
+
+export const sqrtc = {
+   'Complex<T>': ({
+      'isZero(T)': isZ,
+      'sign(T)': sgn,
+      'one(T)': uno,
+      'add(T,T)': plus,
+      'complex(T)': cplxU,
+      'complex(T,T)': cplxB,
+      'multiply(T,T)': mult,
+      'sqrt(T)': sqt,
+      'divide(T,T)': div,
+      'absquare(Complex<T>)': absqC,
+      'subtract(T,T)': sub
+   }) => {
+      if (isZ.checkingDependency) return undefined
+      let baseReturns = returnTypeOf(sqt)
+      if (baseReturns.includes('|')) {
+         // Bit of a hack, because it is relying on other implementations
+         // to list the "typical" value of sqrt first
+         baseReturns = baseReturns.split('|', 1)[0]
+      }
+      return Returns(`Complex<${baseReturns}>`, z => {
+         const reOne = uno(z.re)
+         if (isZ(z.im) && sgn(z.re) === reOne) return cplxU(sqt(z.re))
+         const myabs = sqt(absqC(z))
+         const reTwo = plus(reOne, reOne)
+         const reQuot = div(plus(myabs, z.re), reTwo)
+         const imQuot = div(sub(myabs, z.re), reTwo)
+         if (reQuot === undefined || imQuot === undefined) {
+            throw new TypeError(`Cannot compute sqrt of ${z.re} + {z.im}i`)
+         }
+         return cplxB(
+            mult(sgn(z.im), sqt(div(plus(myabs, z.re), reTwo))),
+            sqt(div(sub(myabs, z.re), reTwo))
+         )
+      })
+   }
+}

src/core/Returns.mjs

@@ -0,0 +1,34 @@
+/* Annotate a function with its return type */
+
+/* Unfortunately JavaScript is missing a way to cleanly clone a function
+ * object, see https://stackoverflow.com/questions/1833588
+ */
+
+const clonedFrom = Symbol('the original function this one was cloned from')
+function cloneFunction(fn) {
+   const behavior = fn[clonedFrom] || fn // don't nest clones
+   const theClone = function () { return behavior.apply(this, arguments) }
+   Object.assign(theClone, fn)
+   theClone[clonedFrom] = body
+   Object.defineProperty(
+      theClone, 'name', {value: fn.name, configurable: true })
+   return theClone
+}
+
+export function Returns(type, fn) {
+   if ('returns' in fn) fn = cloneFunction(fn)
+   fn.returns = type
+   return fn
+}
+
+export function returnTypeOf(fn, signature, pmInstance) {
+   const typeOfReturns = typeof fn.returns
+   if (typeOfReturns === 'undefined') return 'any'
+   if (typeOfReturns === 'string') return fn.returns
+   // not sure if we will need a function to determine the return type,
+   // but allow it for now:
+   return fn.returns(signature, pmInstance)
+}
+
+export default Returns
+

src/core/extractors.mjs

@@ -6,7 +6,7 @@ export function dependencyExtractor(destinationSet) {
    return new Proxy({}, {
       get: (target, property) => {
          destinationSet.add(property)
-         return {}
+         return {checkingDependency: true}
       }
    })
 }

src/core/utils.mjs

@@ -8,6 +8,8 @@ export function subsetOfKeys(set, obj) {
 
 /* Returns a list of the types mentioned in a typed-function signature */
 export function typeListOfSignature(signature) {
+   signature = signature.trim()
+   if (!signature) return []
    return signature.split(',').map(s => s.trim())
 }
 

src/generic/Types/adapted.mjs

@@ -0,0 +1,91 @@
+import PocomathInstance from '../../core/PocomathInstance.mjs'
+import Returns from '../../core/Returns.mjs'
+
+/* creates a PocomathInstance incorporating a new numeric type encapsulated
+ * as a class. (This instance can the be `install()ed` in another to add the
+ * type so created.)
+ *
+ * @param {string} name       The name of the new type
+ * @param {class} Thing       The class implementing the new type
+ * @param {object} overrides  Patches to the auto-generated adaptation
+ */
+export default function adapted(name, Thing, overrides) {
+   const thing = new PocomathInstance('Adapted Thing')
+   const test = overrides.isa || Thing.isa || (x => x instanceof Thing)
+   thing.installType(name, {
+      test,
+      from: overrides.from || {},
+      before: overrides.before || [],
+      refines: overrides.refines || undefined
+   })
+
+   // Build the operations for Thing
+   const operations = {}
+   // first a creator function, with name depending on the name of the thing:
+   const creatorName = overrides.creatorName || name.toLowerCase()
+   const creator = overrides[creatorName]
+      ? overrides[creatorName]['']
+      : Thing[creatorName]
+         ? (Thing[creatorName])
+         : ((...args) => new Thing(...args))
+   const defaultCreatorImps = {
+      '': () => Returns(name, () => creator()),
+      '...any': () => Returns(name, args => creator(...args))
+   }
+   defaultCreatorImps[name] = () => Returns(name, x => x) // x.clone(x)?
+   operations[creatorName] = overrides[creatorName] || defaultCreatorImps
+
+   // We make the default instance, just as a place to check for methods
+   const instance = overrides.instance || creator()
+
+   // Now adapt the methods to typed-function:
+   const unaryOps = {
+      abs: 'abs',
+      ceiling: 'ceil',
+      floor: 'floor',
+      invert: 'inverse',
+      round: 'round',
+      sqrt: 'sqrt',
+      negate: 'neg'
+   }
+   const binaryOps = {
+      add: ['add', name],
+      compare: ['compare', name],
+      divide: ['div', name],
+      equalTT: ['equals', 'boolean'],
+      gcd: ['gcd', name],
+      lcm: ['lcm', name],
+      mod: ['mod', name],
+      multiply: ['mul', name],
+      subtract: ['sub', name]
+   }
+   for (const [mathname, standardname] of Object.entries(unaryOps)) {
+      if (standardname in instance) {
+         operations[mathname] = {}
+         operations[mathname][name] = () => Returns(name, t => t[standardname]())
+      }
+   }
+   operations.zero = {}
+   operations.zero[name] = () => Returns(name, t => creator())
+   operations.one = {}
+   operations.one[name] = () => Returns(name, t => creator(1))
+   operations.conjugate = {}
+   operations.conjugate[name] = () => Returns(name, t => t) // or t.clone() ??
+
+   const binarySignature = `${name},${name}`
+   for (const [mathname, spec] of Object.entries(binaryOps)) {
+      if (spec[0] in instance) {
+         operations[mathname] = {}
+         operations[mathname][binarySignature] = () => Returns(
+            spec[1], (t,u) => t[spec[0]](u))
+      }
+   }
+   if ('operations' in overrides) {
+      Object.assign(operations, overrides.operations)
+   }
+
+   thing.install(operations)
+   return thing
+}
+
+export {adapted}

src/generic/abs.mjs

@@ -0,0 +1,9 @@
+import Returns from '../core/Returns.mjs'
+export const abs = {
+   T: ({
+      T,
+      'smaller(T,T)': lt,
+      'negate(T)': neg,
+      'zero(T)': zr
+   }) => Returns(T, t => (smaller(t, zr(t)) ? neg(t) : t))
+}

src/generic/absquare.mjs

@@ -0,0 +1,9 @@
+import Returns from '../core/Returns.mjs'
+
+export const absquare = {
+   T: ({
+      T,
+      'square(T)': sq,
+      'abs(T)': abval
+   }) => Returns(T, t => sq(abval(t)))
+}

src/generic/all.mjs

@@ -1,2 +1,26 @@
+import {adapted} from './Types/adapted.mjs'
+import Fraction from 'fraction.js/bigfraction.js'
+import Returns from '../core/Returns.mjs'
+
 export * from './arithmetic.mjs'
 export * from './relational.mjs'
+
+export const fraction = adapted('Fraction', Fraction, {
+   before: ['Complex'],
+   from: {number: n => new Fraction(n)},
+   operations: {
+      compare: {
+         'Fraction,Fraction': () => Returns(
+            'Fraction', (f,g) => new Fraction(f.compare(g)))
+      },
+      mod: {
+         'Fraction,Fraction': () => Returns('Fraction', (n,d) => {
+            // patch for "mathematician's modulus"
+            // OK to use full public API of Fraction here
+            const fmod = n.mod(d)
+            if (fmod.s === -1n) return fmod.add(d.abs())
+            return fmod
+         })
+      }
+   }
+})

src/generic/arithmetic.mjs

@@ -2,14 +2,18 @@ import {reducingOperation} from './reducingOperation.mjs'
 
 export * from './Types/generic.mjs'
 
+export {abs} from './abs.mjs'
+export {absquare} from './absquare.mjs'
 export const add = reducingOperation
+export {divide} from './divide.mjs'
 export const gcd = reducingOperation
 export {identity} from './identity.mjs'
 export {lcm} from './lcm.mjs'
 export {mean} from './mean.mjs'
 export {mod} from './mod.mjs'
 export const multiply = reducingOperation
-export {divide} from './divide.mjs'
+export {quotient} from './quotient.mjs'
+export {roundquotient} from './roundquotient.mjs'
 export {sign} from './sign.mjs'
 export {sqrt} from './sqrt.mjs'
 export {square} from './square.mjs'

src/generic/divide.mjs

@@ -1,7 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export const divide = {
    'T,T': ({
+      T,
       'multiply(T,T)': multT,
       'invert(T)': invT
-   }) => (x, y) => multT(x, invT(y))
+   }) => Returns(T, (x, y) => multT(x, invT(y)))
 }
 

src/generic/gcdType.mjs

@@ -1,3 +1,5 @@
+import Returns from '../core/Returns.mjs'
+
 /* Note we do not use a template here so that we can explicitly control
  * which types this is instantiated for, namely the "integer" types, and
  * not simply allow Pocomath to generate instances for any type it encounters.
@@ -7,14 +9,14 @@ export default function(type) {
    const producer = refs => {
       const modder = refs[`mod(${type},${type})`]
       const zeroTester = refs[`isZero(${type})`]
-      return (a,b) => {
+      return Returns(type, (a,b) => {
          while (!zeroTester(b)) {
             const r = modder(a,b)
             a = b
             b = r
          }
          return a
-      }
+      })
    }
    const retval = {}
    retval[`${type},${type}`] = producer

src/generic/identity.mjs

@@ -1,3 +1,11 @@
-export function identity(x) {
-   return x
+import Returns from '../core/Returns.mjs'
+
+export function identityType(type) {
+   return () => Returns(type, x => x)
 }
+
+export function identitySubTypes(type) {
+   return ({T}) => Returns(T, x => x)
+}
+
+export const identity = {T: ({T}) => Returns(T, x => x)}

src/generic/lcm.mjs

@@ -1,10 +1,12 @@
+import Returns from '../core/Returns.mjs'
 import {reducingOperation} from './reducingOperation.mjs'
 
 export const lcm = {
    'T,T': ({
+      T,
       'multiply(T,T)': multT,
       'quotient(T,T)': quotT,
       'gcd(T,T)': gcdT
-   }) => (a,b) => multT(quotT(a, gcdT(a,b)), b)
+   }) => Returns(T, (a,b) => multT(quotT(a, gcdT(a,b)), b))
 }
 Object.assign(lcm, reducingOperation)

src/generic/mean.mjs

@@ -1,3 +1,8 @@
+import Returns from '../core/Returns.mjs'
 export const mean = {
-   '...any': ({add, divide}) => args => divide(add(...args), args.length)
+   '...T': ({
+      T,
+      add,
+      'divide(T,NumInt)': div
+   }) => Returns(T, args => div(add(...args), args.length))
 }

src/generic/mod.mjs

@@ -1,7 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export const mod = {
    'T,T': ({
+      T,
       'subtract(T,T)': subT,
       'multiply(T,T)': multT,
       'quotient(T,T)': quotT
-   }) => (a,m) => subT(a, multT(m, quotT(a,m)))
+   }) => Returns(T, (a,m) => subT(a, multT(m, quotT(a,m))))
 }

src/generic/quotient.mjs

@@ -0,0 +1,9 @@
+import Returns from '../core/Returns.mjs'
+
+export const quotient = {
+    'T,T': ({
+        T,
+        'floor(T)': flr,
+        'divide(T,T)': div
+    }) => Returns(T, (n,d) => flr(div(n,d)))
+}

src/generic/reducingOperation.mjs

@@ -1,12 +1,16 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/generic.mjs'
 
 export const reducingOperation = {
-   'undefined': () => u => u,
-   'undefined,...any': () => (u, rest) => u,
-   'any,undefined': () => (x, u) => u,
-   any: () => x => x,
+   'undefined': () => Returns('undefined', u => u),
+   'undefined,...any': () => Returns('undefined', (u, rest) => u),
+   'any,undefined': () => Returns('undefined', (x, u) => u),
+   'undefined,undefined': () => Returns('undefined', (u,v) => u),
+   T: ({T}) => Returns(T, x => x),
+   // Unfortunately the type language of Pocomath is not (yet?) expressive
+   // enough to properly type the full reduction signature here:
    'any,any,...any': ({
       self
-   }) => (a,b,rest) => [b, ...rest].reduce((x,y) => self(x,y), a)
+   }) => Returns('any', (a,b,rest) => [b, ...rest].reduce((x,y) => self(x,y), a))
 }
 

src/generic/relational.mjs

@@ -1,67 +1,83 @@
+import Returns from '../core/Returns.mjs'
+
 export const compare = {
-   'undefined,undefined': () => () => 0
+   'undefined,undefined': () => Returns('NumInt', () => 0)
 }
 
 export const isZero = {
-   'undefined': () => u => u === 0
+   'undefined': () => Returns('boolean', u => u === 0),
+   T: ({
+      T,
+      'equal(T,T)': eq,
+      'zero(T)': zr
+   }) => Returns('boolean', t => eq(t, zr(t)))
 }
 
 export const equal = {
-   'any,any': ({equalTT, joinTypes, Templates, typeOf}) => (x,y) => {
+   'any,any': ({
+      equalTT,
+      joinTypes,
+      Templates,
+      typeOf
+   }) => Returns('boolean', (x,y) => {
       const resultant = joinTypes([typeOf(x), typeOf(y)], 'convert')
       if (resultant === 'any' || resultant in Templates) {
          return false
       }
       return equalTT(x,y)
-   }
+   })
 }
 
 export const equalTT = {
    'T,T': ({
       'compare(T,T)': cmp,
       'isZero(T)': isZ
-   }) => (x,y) => isZ(cmp(x,y)),
+   }) => Returns('boolean', (x,y) => isZ(cmp(x,y)))
    // If templates were native to typed-function, we should be able to
    // do something like:
    // 'any,any': () => () => false // should only be hit for different types
 }
 
 export const unequal = {
-   'any,any': ({equal}) => (x,y) => !(equal(x,y))
+   'any,any': ({equal}) => Returns('boolean', (x,y) => !(equal(x,y)))
 }
 
 export const larger = {
    'T,T': ({
       'compare(T,T)': cmp,
-      'one(T)' : uno
-   }) => (x,y) => cmp(x,y) === uno(y)
+      'one(T)' : uno,
+      'equalTT(T,T)' : eq
+   }) => Returns('boolean', (x,y) => eq(cmp(x,y), uno(y)))
 }
 
 export const largerEq = {
    'T,T': ({
       'compare(T,T)': cmp,
       'one(T)' : uno,
-      'isZero(T)' : isZ
-   }) => (x,y) => {
+      'isZero(T)' : isZ,
+      'equalTT(T,T)': eq
+   }) => Returns('boolean', (x,y) => {
       const c = cmp(x,y)
-      return isZ(c) || c === uno(y)
-   }
+      return isZ(c) || eq(c, uno(y))
+   })
 }
 
 export const smaller = {
    'T,T': ({
       'compare(T,T)': cmp,
       'one(T)' : uno,
-      'isZero(T)' : isZ
-   }) => (x,y) => {
+      'isZero(T)' : isZ,
+      unequal
+   }) => Returns('boolean', (x,y) => {
       const c = cmp(x,y)
-      return !isZ(c) && c !== uno(y)
-   }
+      return !isZ(c) && unequal(c, uno(y))
+   })
 }
 
 export const smallerEq = {
    'T,T': ({
       'compare(T,T)': cmp,
-      'one(T)' : uno
-   }) => (x,y) => cmp(x,y) !== uno(y)
+      'one(T)' : uno,
+      unequal
+   }) => Returns('boolean', (x,y) => unequal(cmp(x,y), uno(y)))
 }

src/generic/roundquotient.mjs

@@ -0,0 +1,9 @@
+import Returns from '../core/Returns.mjs'
+
+export const roundquotient = {
+   'T,T': ({
+      T,
+      'round(T)': rnd,
+      'divide(T,T)':div
+   }) => Returns(T, (n,d) => rnd(div(n,d)))
+}

src/generic/sign.mjs

@@ -1,3 +1,9 @@
+import Returns from '../core/Returns.mjs'
+
 export const sign = {
-   T: ({'compare(T,T)': cmp, 'zero(T)': Z}) => x => cmp(x, Z(x))
+   T: ({
+      T,
+      'compare(T,T)': cmp,
+      'zero(T)': Z
+   }) => Returns(T, x => cmp(x, Z(x)))
 }

src/generic/sqrt.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/generic.mjs'
 
-export const sqrt = {undefined: () => () => undefined}
+export const sqrt = {undefined: () => Returns('undefined', () => undefined)}

src/generic/square.mjs

@@ -1,3 +1,6 @@
+import {Returns, returnTypeOf} from '../core/Returns.mjs'
+
 export const square = {
-   T: ({'multiply(T,T)': multT}) => x => multT(x,x)
+   T: ({'multiply(T,T)': multT}) => Returns(
+      returnTypeOf(multT), x => multT(x,x))
 }

src/generic/subtract.mjs

@@ -1,3 +1,9 @@
+import Returns from '../core/Returns.mjs'
+
 export const subtract = {
-   'T,T': ({'add(T,T)': addT, 'negate(T)': negT}) => (x,y) => addT(x, negT(y))
+   'T,T': ({
+      T,
+      'add(T,T)': addT,
+      'negate(T)': negT
+   }) => Returns(T, (x,y) => addT(x, negT(y)))
 }

src/number/abs.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
-export const abs = {number: () => n => Math.abs(n)}
+export const abs = {'T:number': ({T}) => Returns(T, n => Math.abs(n))}

src/number/absquare.mjs

@@ -0,0 +1,7 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/number.mjs'
+
+/* Absolute value squared */
+export const absquare = {
+    'T:number': ({T, 'square(T)': sqn}) => Returns(T, n => sqn(n))
+}

src/number/add.mjs

@@ -1,3 +1,8 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
-export const add = {'number,number': () => (m,n) => m+n}
+export const add = {
+    // Note the below assumes that all subtypes of number that will be defined
+    // are closed under addition!
+    'T:number,T': ({T}) => Returns(T, (m,n) => m+n)
+}

src/number/cbrt.mjs

@@ -0,0 +1,19 @@
+import Returns from '../core/Returns.mjs'
+export * from './Types/number.mjs'
+
+/* Returns just the real cube root, following mathjs implementation */
+export const cbrt = {
+   number: ({'negate(number)': neg}) => Returns('number', x => {
+      if (x === 0) return x
+      const negate = x < 0
+      if (negate) x = neg(x)
+      let result = x
+      if (isFinite(x)) {
+         result = Math.exp(Math.log(x) / 3)
+         result = (x / (result * result) + (2 * result)) / 3
+      }
+      if (negate) return neg(result)
+      return result
+   })
+}
+                                                

src/number/compare.mjs

@@ -1,3 +1,5 @@
+import Returns from '../core/Returns.mjs'
+
 /* Lifted from mathjs/src/utils/number.js */
 /**
  * Minimum number added to one that makes the result different than one
@@ -48,5 +50,6 @@ function nearlyEqual (x, y, epsilon) {
 export const compare = {
   'number,number': ({
     config
-  }) => (x,y) => nearlyEqual(x, y, config.epsilon) ? 0 : (x > y ? 1 : -1)
+  }) => Returns(
+    'NumInt', (x,y) => nearlyEqual(x, y, config.epsilon) ? 0 : (x > y ? 1 : -1))
 }

src/number/invert.mjs

@@ -1,3 +1,5 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/number.mjs'
 
-export const invert = {number: () => n => 1/n}
+export const invert = {number: () => Returns('number', n => 1/n)}

src/number/isZero.mjs

@@ -1,3 +1,6 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
-export const isZero = {number: () => n => n === 0}
+export const isZero = {
+    'T:number': () => Returns('boolean', n => n === 0)
+}

src/number/multiply.mjs

@@ -1,3 +1,5 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/number.mjs'
 
-export const multiply = {'number,number': () => (m,n) => m*n}
+export const multiply = {'T:number,T': ({T}) => Returns(T, (m,n) => m*n)}

src/number/native.mjs

@@ -1,12 +1,14 @@
 import gcdType from '../generic/gcdType.mjs'
-import {identity} from '../generic/identity.mjs'
+import {identitySubTypes} from '../generic/identity.mjs'
 
 export * from './Types/number.mjs'
 
 export {abs} from './abs.mjs'
+export {absquare} from './absquare.mjs'
 export {add} from './add.mjs'
+export {cbrt} from './cbrt.mjs'
 export {compare} from './compare.mjs'
-export const conjugate = {number: () => identity}
+export const conjugate = {'T:number': identitySubTypes('number')}
 export const gcd = gcdType('NumInt')
 export {invert} from './invert.mjs'
 export {isZero} from './isZero.mjs'

src/number/negate.mjs

@@ -1,3 +1,6 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
-export const negate = {number: () => n => -n}
+export const negate = {
+   'T:number': ({T}) => Returns(T, n => -n)
+}

src/number/one.mjs

@@ -1,3 +1,5 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/number.mjs'
 
-export const one = {number: () => () => 1}
+export const one = {number: () => Returns('NumInt', () => 1)}

src/number/quotient.mjs

@@ -1,8 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/number.mjs'
 
 export const quotient = {
-   'number,number': () => (n,d) => {
+   'T:number,T': () => Returns('NumInt', (n,d) => {
       if (d === 0) return d
       return Math.floor(n/d)
-   }
+   })
 }

src/number/roundquotient.mjs

@@ -1,8 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/number.mjs'
 
 export const roundquotient = {
-   'number,number': () => (n,d) => {
+   'number,number': () => Returns('NumInt', (n,d) => {
       if (d === 0) return d
       return Math.round(n/d)
-   }
+   })
 }

src/number/sqrt.mjs

@@ -1,14 +1,18 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
 export const sqrt = {
-   number: ({config, complex, 'self(Complex)': complexSqrt}) => {
-      if (config.predictable || !complexSqrt) {
-         return n => isNaN(n) ? NaN : Math.sqrt(n)
+   number: ({
+      config,
+      'complex(number,number)': cplx,
+      'negate(number)': neg}) => {
+         if (config.predictable || !cplx) {
+            return Returns('number', n => isNaN(n) ? NaN : Math.sqrt(n))
+         }
+         return Returns('number|Complex<number>', n => {
+            if (isNaN(n)) return NaN
+            if (n >= 0) return Math.sqrt(n)
+            return cplx(0, Math.sqrt(neg(n)))
+         })
       }
-      return n => {
-         if (isNaN(n)) return NaN
-         if (n >= 0) return Math.sqrt(n)
-         return complexSqrt(complex(n))
-      }
-   }
 }

src/number/zero.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export * from './Types/number.mjs'
 
-export const zero = {number: () => () => 0}
+export const zero = {number: () => Returns('NumInt', () => 0)}

src/ops/choose.mjs

@@ -1,11 +1,13 @@
-/* Note this is not a good algorithm for computing binomial coefficients,
+import Returns from '../core/Returns.mjs'
+
+/* Note this is _not_ a good algorithm for computing binomial coefficients,
  * it's just for demonstration purposes
  */
 export const choose = {
-   'NumInt,NumInt': ({factorial}) => (n,k) => Number(
-      factorial(n) / (factorial(k)*factorial(n-k))),
+   'NumInt,NumInt': ({factorial}) => Returns(
+      'NumInt', (n,k) => Number(factorial(n) / (factorial(k)*factorial(n-k)))),
    'bigint,bigint': ({
       factorial
-   }) => (n,k) => factorial(n) / (factorial(k)*factorial(n-k))
+   }) => Returns('bigint', (n,k) => factorial(n) / (factorial(k)*factorial(n-k)))
 }
         

src/ops/factorial.mjs

@@ -1,8 +1,15 @@
-export function factorial(n) {
+import Returns from '../core/Returns.mjs'
+
+/* Plain functions are OK, too, and they can be decorated with a return type
+ * just like an implementation.
+ */
+const factorial = Returns('bigint', function factorial(n) {
    n = BigInt(n)
    let prod = 1n
    for (let i = n; i > 1n; --i) {
       prod *= i
    }
    return prod
-}
+})
+
+export {factorial}

src/ops/floor.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 import {Complex} from '../complex/Types/Complex.mjs'
 
 /* Note we don't **export** any types here, so that only the options
@@ -5,21 +6,28 @@ import {Complex} from '../complex/Types/Complex.mjs'
  */
 
 export const floor = {
-   bigint: () => x => x,
-   NumInt: () => x => x, // Because Pocomath isn't part of typed-function, or
-   GaussianInteger: () => x => x, // at least have access to the real
-   // typed-function parse, we unfortunately can't coalesce these into one
-   // entry with type `bigint|NumInt|GaussianInteger` because they couldn't
-   // be separately activated then
+   /* Because Pocomath isn't part of typed-function, nor does it have access
+    * to the real typed-function parse, we unfortunately can't coalesce the
+    * first several implementations into one entry with type
+    * `bigint|NumInt|GaussianInteger` because then they couldn't
+    * be separately activated
+    */
+   bigint: () => Returns('bigint', x => x),
+   NumInt: () => Returns('NumInt', x => x),
+   'Complex<bigint>': () => Returns('Complex<bigint>', x => x),
 
-   number: ({'equalTT(number,number)': eq}) => n => {
+   number: ({'equalTT(number,number)': eq}) => Returns('NumInt', n => {
       if (eq(n, Math.round(n))) return Math.round(n)
       return Math.floor(n)
-   },
+   }),
 
-   Complex: Complex.promoteUnary.Complex,
+   'Complex<T>': Complex.promoteUnary['Complex<T>'],
 
    // OK to include a type totally not in Pocomath yet, it'll never be
    // activated.
-   Fraction: ({quotient}) => f => quotient(f.n, f.d),
+   BigNumber: ({
+      'round(BigNumber)': rnd,
+      'equal(BigNumber,BigNumber)': eq
+   }) => Returns('BigNumber', x => eq(x,round(x)) ? round(x) : x.floor())
+
 }

src/tuple/Types/Tuple.mjs

@@ -1,25 +1,24 @@
 /* A template type representing a homeogeneous tuple of elements */
 import PocomathInstance from '../../core/PocomathInstance.mjs'
+import {Returns, returnTypeOf} from '../../core/Returns.mjs'
 
 const Tuple = new PocomathInstance('Tuple')
-// First a base type that will generally not be used directly
-Tuple.installType('Tuple', {
-   test: t => t && typeof t === 'object' && 'elts' in t && Array.isArray(t.elts)
-})
-// Now the template type that is the primary use of this
+
 Tuple.installType('Tuple<T>', {
-   // We are assuming that any 'Type<T>' refines 'Type', so this is
-   // not necessary:
-   // refines: 'Tuple',
-   // But we need there to be a way to determine the type of a tuple:
-   infer: ({typeOf, joinTypes}) => t => joinTypes(t.elts.map(typeOf)),
-   // For the test, we can assume that t is already a base tuple,
-   // and we get the test for T as an input and we have to return
-   // the test for Tuple<T>
+   // A test that "defines" the "base type", which is not really a type
+   // (only fully instantiated types are added to the universe)
+   base: t => t && typeof t === 'object' && 'elts' in t && Array.isArray(t.elts),
+   // The template portion of the test; it takes the test for T as
+   // input and returns the test for an entity _that already passes
+   // the base test_ to be a Tuple<T>:
    test: testT => t => t.elts.every(testT),
-   // These are only invoked for types U such that there is already
-   // a conversion from U to T, and that conversion is passed as an input
-   // and we have to return the conversion to Tuple<T>:
+   // And we need there to be a way to determine the (instantiation)
+   // type of an tuple (that has already passed the base test):
+   infer: ({typeOf, joinTypes}) => t => joinTypes(t.elts.map(typeOf)),
+   // Conversions. Parametrized conversions are only invoked for types
+   // U such that there is already a conversion from U to T, and that
+   // conversion is passed as an input, and we have to return the conversion
+   // function from the indicated template in terms of U to Tuple<T>:
    from: {
       'Tuple<U>': convert => tu => ({elts: tu.elts.map(convert)}),
       // Here since there is no U it's a straight conversion:
@@ -32,48 +31,69 @@ Tuple.installType('Tuple<T>', {
 })
 
 Tuple.promoteUnary = {
-   'Tuple<T>': ({'self(T)': me, tuple}) => t => tuple(...(t.elts.map(me)))
+   'Tuple<T>': ({
+      'self(T)': me,
+      tuple
+   }) => {
+      const compType = me.fromInstance.joinTypes(
+         returnTypeOf(me).split('|'), 'convert')
+      return Returns(
+         `Tuple<${compType}>`, t => tuple(...(t.elts.map(x => me(x)))))
+   }
 }
 
 Tuple.promoteBinaryUnary = {
-   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, 'self(T)': meU, tuple}) => (s,t) => {
-      let i = -1
-      let result = []
-      while (true) {
-         i += 1
-         if (i < s.elts.length) {
-            if (i < t.elts.length) result.push(meB(s.elts[i], t.elts[i]))
-            else result.push(meU(s.elts[i]))
-            continue
+   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, 'self(T)': meU, tuple}) => {
+      const compTypes = returnTypeOf(meB).split('|').concat(
+         returnTypeOf(meU).split('|'))
+      const compType = meU.fromInstance.joinTypes(compTypes, 'convert')
+      return Returns(`Tuple<${compType}>`, (s,t) => {
+         let i = -1
+         let result = []
+         while (true) {
+            i += 1
+            if (i < s.elts.length) {
+               if (i < t.elts.length) result.push(meB(s.elts[i], t.elts[i]))
+               else result.push(meU(s.elts[i]))
+               continue
+            }
+            if (i < t.elts.length) result.push(meU(t.elts[i]))
+            else break
          }
-         if (i < t.elts.length) result.push(meU(t.elts[i]))
-         else break
-      }
-      return tuple(...result)
+         return tuple(...result)
+      })
    }
 }
 
 Tuple.promoteBinary = {
-   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, tuple}) => (s,t) => {
-      const lim = Math.max(s.elts.length, t.elts.length)
-      const result = []
-      for (let i = 0; i < lim; ++i) {
-         result.push(meB(s.elts[i], t.elts[i]))
-      }
-      return tuple(...result)
+   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, tuple}) => {
+      const compType = meB.fromInstance.joinTypes(
+         returnTypeOf(meB).split('|'))
+      return Returns(`Tuple<${compType}>`, (s,t) => {
+         const lim = Math.max(s.elts.length, t.elts.length)
+         const result = []
+         for (let i = 0; i < lim; ++i) {
+            result.push(meB(s.elts[i], t.elts[i]))
+         }
+         return tuple(...result)
+      })
    }
 }
 
 Tuple.promoteBinaryStrict = {
-   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, tuple}) => (s,t) => {
-      if (s.elts.length !== t.elts.length) {
-         throw new RangeError('Tuple length mismatch') // get name of self ??
-      }
-      const result = []
-      for (let i = 0; i < s.elts.length; ++i) {
-         result.push(meB(s.elts[i], t.elts[i]))
-      }
-      return tuple(...result)
+   'Tuple<T>,Tuple<T>': ({'self(T,T)': meB, tuple}) => {
+      const compType = meB.fromInstance.joinTypes(
+         returnTypeOf(meB).split('|'))
+      return Returns(`Tuple<${compType}>`, (s,t) => {
+         if (s.elts.length !== t.elts.length) {
+            throw new RangeError('Tuple length mismatch') // get name of self ??
+         }
+         const result = []
+         for (let i = 0; i < s.elts.length; ++i) {
+            result.push(meB(s.elts[i], t.elts[i]))
+         }
+         return tuple(...result)
+      })
    }
 }
 

src/tuple/equalTT.mjs

@@ -1,11 +1,16 @@
+import Returns from '../core/Returns.mjs'
+
 export * from './Types/Tuple.mjs'
 
 export const equalTT = {
-   'Tuple<T>,Tuple<T>': ({'self(T,T)': me, 'length(Tuple)': len}) => (s,t) => {
+   'Tuple<T>,Tuple<T>': ({
+      'self(T,T)': me,
+      'length(Tuple)': len
+   }) => Returns('boolean', (s,t) => {
       if (len(s) !== len(t)) return false
       for (let i = 0; i < len(s); ++i) {
          if (!me(s.elts[i], t.elts[i])) return false
       }
       return true
-   }
+   })
 }

src/tuple/isZero.mjs

@@ -1,7 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export {Tuple} from './Types/Tuple.mjs'
 
 export const isZero = {
-   'Tuple<T>': ({'self(T)': me}) => t => t.elts.every(e => me(e))
+   'Tuple<T>': ({'self(T)': me}) => Returns(
+      'boolean', t => t.elts.every(e => me(e)))
    // Note we can't just say `every(me)` above since every invokes its
    // callback with more arguments, which then violates typed-function's
    // signature for `me`

src/tuple/length.mjs

@@ -1,3 +1,4 @@
+import Returns from '../core/Returns.mjs'
 export {Tuple} from './Types/Tuple.mjs'
 
-export const length = {Tuple: () => t => t.elts.length}
+export const length = {'Tuple<T>': () => Returns('NumInt', t => t.elts.length)}

src/tuple/tuple.mjs

@@ -1,6 +1,10 @@
+import Returns from '../core/Returns.mjs'
+
 export {Tuple} from './Types/Tuple.mjs'
 
 /* The purpose of the template argument is to ensure that all of the args
  * are convertible to the same type.
  */
-export const tuple = {'...T': () => args => ({elts: args})}
+export const tuple = {
+   '...T': ({T}) => Returns(`Tuple<${T}>`, args => ({elts: args}))
+}

test/_pocomath.mjs

@@ -16,13 +16,60 @@ describe('The default full pocomath instance "math"', () => {
       assert.strictEqual(math.typeOf(-1.5), 'number')
       assert.strictEqual(math.typeOf(-42n), 'bigint')
       assert.strictEqual(math.typeOf(undefined), 'undefined')
-      assert.strictEqual(math.typeOf({re: 15n, im: -2n}), 'GaussianInteger')
-      assert.strictEqual(math.typeOf({re: 6.28, im: 2.72}), 'Complex')
+      assert.strictEqual(math.typeOf({re: 15n, im: -2n}), 'Complex<bigint>')
+      assert.strictEqual(math.typeOf({re: 6.28, im: 2.72}), 'Complex<number>')
+   })
+
+   it('can determine the return types of operations', () => {
+      assert.strictEqual(math.returnTypeOf('negate', 'number'), 'number')
+      assert.strictEqual(math.returnTypeOf('negate', 'NumInt'), 'NumInt')
+      math.negate(math.complex(1.2, 2.8)) // TODO: make this call unnecessary
+      assert.strictEqual(
+         math.returnTypeOf('negate', 'Complex<number>'), 'Complex<number>')
+      assert.strictEqual(math.returnTypeOf('add', 'number,number'), 'number')
+      assert.strictEqual(math.returnTypeOf('add', 'NumInt,NumInt'), 'NumInt')
+      assert.strictEqual(math.returnTypeOf('add', 'NumInt,number'), 'number')
+      assert.strictEqual(math.returnTypeOf('add', 'number,NumInt'), 'number')
+      assert.deepStrictEqual(  // TODO: ditto
+         math.add(3, math.complex(2.5, 1)), math.complex(5.5, 1))
+      assert.strictEqual(
+         math.returnTypeOf('add', 'Complex<number>,NumInt'), 'Complex<number>')
+      // The following is not actually what we want, but the Pocomath type
+      // language isn't powerful enough at this point to capture the true
+      // return type:
+      assert.strictEqual(
+         math.returnTypeOf('add', 'number,NumInt,Complex<number>'), 'any')
+      assert.strictEqual(
+         math.returnTypeOf('chain', 'bigint'), 'Chain<bigint>')
+      assert.strictEqual(
+         math.returnTypeOf('returnTypeOf', 'string,string'), 'string')
+      assert.strictEqual(
+         math.returnTypeOf('conjugate', 'bigint'), 'bigint')
+      assert.strictEqual(
+         math.returnTypeOf('gcd', 'bigint,bigint'), 'bigint')
+      math.identity(math.fraction(3,5)) // TODO: ditto
+      assert.strictEqual(math.returnTypeOf('identity', 'Fraction'), 'Fraction')
+      assert.strictEqual(
+         math.returnTypeOf('quotient', 'bigint,bigint'), 'bigint')
+      math.abs(math.complex(2,1)) //TODO: ditto
+      assert.strictEqual(
+         math.returnTypeOf('abs','Complex<NumInt>'), 'number')
+      math.multiply(math.quaternion(1,1,1,1), math.quaternion(1,-1,1,-1)) // dit
+      const quatType = math.returnTypeOf(
+         'quaternion', 'NumInt,NumInt,NumInt,NumInt')
+      assert.strictEqual(quatType, 'Complex<Complex<NumInt>>')
+      assert.strictEqual(
+         math.returnTypeOf('multiply', quatType + ',' + quatType), quatType)
+      assert.strictEqual(math.returnTypeOf('isZero', 'NumInt'), 'boolean')
+      assert.strictEqual(
+         math.returnTypeOf('roundquotient', 'NumInt,number'), 'NumInt')
+      assert.strictEqual(
+         math.returnTypeOf('factorial', 'NumInt'), 'bigint')
    })
 
    it('can subtract numbers', () => {
       assert.strictEqual(math.subtract(12, 5), 7)
-      //assert.strictEqual(math.subtract(3n, 1.5), 1.5)
+      assert.throws(() => math.subtract(3n, 1.5), 'TypeError')
    })
 
    it('can add numbers', () => {
@@ -105,11 +152,9 @@ describe('The default full pocomath instance "math"', () => {
 
    it('calculates multi-way gcds and lcms', () => {
       assert.strictEqual(math.gcd(30,105,42), 3)
-      assert.ok(
-         math.associate(
-            math.lcm(
-               math.complex(2n,1n), math.complex(1n,1n), math.complex(0n,1n)),
-            math.complex(1n,3n)))
+      const gaussianLCM = math.lcm(
+         math.complex(2n,1n), math.complex(1n,1n), math.complex(0n,1n))
+      assert.strictEqual(math.associate(gaussianLCM, math.complex(1n,3n)), true)
    })
 
 })

test/complex/_all.mjs

@@ -39,16 +39,30 @@ describe('complex', () => {
    })
 
    it('checks for equality', () => {
-      assert.ok(math.equal(math.complex(3,0), 3))
-      assert.ok(math.equal(math.complex(3,2), math.complex(3, 2)))
+      assert.ok(math.equal(math.complex(3, 0), 3))
+      assert.ok(math.equal(math.complex(3, 2), math.complex(3, 2)))
       assert.ok(!(math.equal(math.complex(45n, 3n), math.complex(45n, -3n))))
       assert.ok(!(math.equal(math.complex(45n, 3n), 45n)))
    })
 
+   it('tests for reality', () => {
+      assert.ok(math.isReal(math.complex(3, 0)))
+      assert.ok(!(math.isReal(math.complex(3, 2))))
+   })
+
    it('computes gcd', () => {
       assert.deepStrictEqual(
          math.gcd(math.complex(53n, 56n), math.complex(47n, -13n)),
          math.complex(4n, 5n))
+      // And now works for NumInt, too!
+      assert.deepStrictEqual(
+         math.gcd(math.complex(53,56), math.complex(47, -13)),
+         math.complex(4, 5))
+      // But properly fails for general complex
+      assert.throws(
+         () => math.gcd(math.complex(5.3,5.6), math.complex(4.7, -1.3)),
+         TypeError
+      )
    })
 
    it('computes floor', () => {
@@ -59,4 +73,24 @@ describe('complex', () => {
       assert.strictEqual(math.floor(gi), gi) // literally a no-op
    })
 
+   it('performs rudimentary quaternion calculations', () => {
+      const q0 = math.quaternion(1, 0, 1, 0)
+      const q1 = math.quaternion(1, 0.5, 0.5, 0.75)
+      assert.deepStrictEqual(
+         q1,
+         math.complex(math.complex(1, 0.5), math.complex(0.5, 0.75)))
+      assert.deepStrictEqual(
+         math.add(q0,q1),
+         math.quaternion(2, 0.5, 1.5, 0.75))
+      assert.deepStrictEqual(
+         math.multiply(q0, q1),
+         math.quaternion(0.5, 1.25, 1.5, 0.25))
+      assert.deepStrictEqual(
+         math.multiply(q0, math.quaternion(2, 1, 0.1, 0.1)),
+         math.quaternion(1.9, 1.1, 2.1, -0.9))
+      math.absquare(math.complex(1.25, 2.5)) //HACK: need absquare(Complex<number>)
+      assert.strictEqual(math.abs(q0), Math.sqrt(2))
+      assert.strictEqual(math.abs(q1), Math.sqrt(33)/4)
+   })
+
 })

test/complex/_polynomialRoot.mjs

@@ -0,0 +1,63 @@
+import assert from 'assert'
+import * as approx from '../../tools/approx.mjs'
+import math from '../../src/pocomath.mjs'
+
+describe('polynomialRoot', () => {
+   it('should solve a linear equation with real coefficients', function () {
+      assert.deepEqual(math.polynomialRoot(6, 3), math.tuple(-2))
+      assert.deepEqual(
+         math.polynomialRoot(math.complex(-3, 2), 2),
+         math.tuple(math.complex(1.5, -1)))
+      assert.deepEqual(
+         math.polynomialRoot(math.complex(3, 1), math.complex(-1, -1)),
+         math.tuple(math.complex(2, -1)))
+   })
+   // Should be safe now to capture the functions:
+   const complex = math.complex
+   const pRoot = math.polynomialRoot
+   const tup = math.tuple
+   it('should solve a quadratic equation with a double root', function () {
+      assert.deepEqual(pRoot(4, 4, 1), tup(-2))
+      assert.deepEqual(
+         pRoot(complex(0, 2), complex(2, 2), 1), tup(complex(-1, -1)))
+   })
+   it('should solve a quadratic with two distinct roots', function () {
+      assert.deepEqual(pRoot(-3, 2, 1), tup(1, -3))
+      assert.deepEqual(pRoot(-2, 0, 1), tup(math.sqrt(2), -math.sqrt(2)))
+      assert.deepEqual(
+         pRoot(4, 2, 1),
+         tup(complex(-1, math.sqrt(3)), complex(-1, -math.sqrt(3))))
+      assert.deepEqual(
+         pRoot(complex(3, 1), -3, 1), tup(complex(1, 1), complex(2, -1)))
+   })
+   it('should solve a cubic with a triple root', function () {
+      assert.deepEqual(pRoot(8, 12, 6, 1), tup(-2))
+      assert.deepEqual(
+         pRoot(complex(-2, 11), complex(9, -12), complex(-6, 3), 1),
+         tup(complex(2, -1)))
+   })
+   it('should solve a cubic with one simple and one double root', function () {
+      assert.deepEqual(pRoot(4, 0, -3, 1), tup(-1, 2))
+      assert.deepEqual(
+         pRoot(complex(9, 9), complex(15, 6), complex(7, 1), 1),
+         tup(complex(-1, -1), -3))
+      assert.deepEqual(
+         pRoot(complex(0, 6), complex(6, 8), complex(5, 2), 1),
+         tup(-3, complex(-1, -1)))
+      assert.deepEqual(
+         pRoot(complex(2, 6), complex(8, 6), complex(5, 1), 1),
+         tup(complex(-3, 1), complex(-1, -1)))
+   })
+   it('should solve a cubic with three distinct roots', function () {
+      approx.deepEqual(pRoot(6, 11, 6, 1), tup(-3, -1, -2))
+      approx.deepEqual(
+         pRoot(-1, -2, 0, 1),
+         tup(-1, (1 + math.sqrt(5)) / 2, (1 - math.sqrt(5)) / 2))
+      approx.deepEqual(
+         pRoot(1, 1, 1, 1),
+         tup(-1, complex(0, -1), complex(0, 1)))
+      approx.deepEqual(
+         pRoot(complex(0, -10), complex(8, 12), complex(-6, -3), 1),
+         tup(complex(1, 1), complex(3, 1), complex(2, 1)))
+   })
+})

test/core/_utils.mjs

@@ -0,0 +1,8 @@
+import assert from 'assert'
+import * as utils from '../../src/core/utils.mjs'
+
+describe('typeListOfSignature', () => {
+   it('returns an empty list for the empty signature', () => {
+      assert.deepStrictEqual(utils.typeListOfSignature(''), [])
+   })
+})

test/custom.mjs

@@ -3,12 +3,14 @@ import math from '../src/pocomath.mjs'
 import PocomathInstance from '../src/core/PocomathInstance.mjs'
 import * as numbers from '../src/number/all.mjs'
 import * as numberAdd from '../src/number/add.mjs'
+import * as numberZero from '../src/number/zero.mjs'
 import {add as genericAdd} from '../src/generic/arithmetic.mjs'
 import * as complex from '../src/complex/all.mjs'
 import * as complexAdd from '../src/complex/add.mjs'
 import * as complexNegate from '../src/complex/negate.mjs'
 import * as complexComplex from '../src/complex/complex.mjs'
 import * as bigintAdd from '../src/bigint/add.mjs'
+import * as bigintZero from '../src/bigint/zero.mjs'
 import * as concreteSubtract from '../src/generic/subtract.concrete.mjs'
 import * as genericSubtract from '../src/generic/subtract.mjs'
 import extendToComplex from '../src/complex/extendToComplex.mjs'
@@ -54,9 +56,9 @@ describe('A custom instance', () => {
          math.complex(-5, -1))
       // And now floor has been activated for Complex as well, since the type
       // is present
-      assert.deepStrictEqual(
-         pm.floor(math.complex(1.9, 0)),
-         math.complex(1))
+      const fracComplex = math.complex(1.9, 0)
+      const intComplex = math.complex(1)
+      assert.deepStrictEqual(pm.floor(fracComplex), intComplex)
       // And the chain functions refresh themselves:
       assert.deepStrictEqual(
          pm.chain(5).add(pm.chain(0).complex(7).value).value, math.complex(5,7))
@@ -65,10 +67,11 @@ describe('A custom instance', () => {
    it("can defer definition of (even used) types", () => {
       const dt = new PocomathInstance('Deferred Types')
       dt.install(numberAdd)
+      dt.install(numberZero) // for promoting numbers to complex, to fill in im
       dt.install({times: {
          'number,number': () => (m,n) => m*n,
-         'Complex,Complex': ({complex}) => (w,z) => {
-            return complex(w.re*z.re - w.im*z.im, w.re*z.im + w.im*z.re)
+         'Complex<T>,Complex<T>': ({'complex(T,T)': cplx}) => (w,z) => {
+            return cplx(w.re*z.re - w.im*z.im, w.re*z.im + w.im*z.re)
          }
       }})
       // complex type not present but should still be able to add numbers:
@@ -82,6 +85,7 @@ describe('A custom instance', () => {
    it("can selectively import in cute ways", async function () {
       const cherry = new PocomathInstance('cherry')
       cherry.install(numberAdd)
+      cherry.install(numberZero) // for complex promotion
       await extendToComplex(cherry)
       cherry.install({add: genericAdd})
       /* Now we have an instance that supports addition for number and complex
@@ -124,19 +128,15 @@ describe('A custom instance', () => {
       inst.install(complexAdd)
       inst.install(complexComplex)
       inst.install(bigintAdd)
+      inst.install(bigintZero) // for complex promotion
       assert.strictEqual(
          inst.typeMerge(6n, inst.complex(3n, 2n)),
-         'Merge to GaussianInteger')
+         'Merge to Complex<bigint>')
       assert.strictEqual(
          inst.typeMerge(3, inst.complex(4.5,2.1)),
-         'Merge to Complex')
-      // The following is the current behavior, since 3 converts to 3+0i
-      // which is technically the same Complex type as 3n+0ni.
-      // This should clear up when Complex is templatized
-      assert.strictEqual(inst.typeMerge(3, inst.complex(3n)), 'Merge to Complex')
-      // But types that truly cannot be merged should throw a TypeError
-      // Should add a variation of this with a more usual type once there is
-      // one not interconvertible with others...
+         'Merge to Complex<number>')
+      assert.throws(
+         () => inst.typeMerge(3, inst.complex(3n)), TypeError)
       inst.install(genericSubtract)
       assert.throws(() => inst.typeMerge(3, undefined), TypeError)
    })

test/generic/_all.mjs

@@ -0,0 +1,20 @@
+import assert from 'assert'
+import math from '../../src/pocomath.mjs'
+
+describe('generic', () => {
+   it('calculates mean', () => {
+      assert.strictEqual(math.mean(1,2.5,3.25,4.75), 2.875)
+      assert.strictEqual(
+         math.returnTypeOf('mean', 'number,number,number,number'),
+         'number'
+      )
+   })
+   it('compares things', () => {
+      assert.strictEqual(math.larger(7n, 3n), true)
+      assert.strictEqual(
+         math.returnTypeOf('larger', 'bigint,bigint'), 'boolean')
+      assert.strictEqual(math.smallerEq(7.2, 3), false)
+      assert.strictEqual(
+         math.returnTypeOf('smallerEq', 'number,NumInt'), 'boolean')
+   })
+})

test/generic/fraction.mjs

@@ -0,0 +1,101 @@
+import assert from 'assert'
+import math from '../../src/pocomath.mjs'
+import Fraction from 'fraction.js/bigfraction.js'
+
+describe('fraction', () => {
+   const half = new Fraction(1/2)
+   const tf = new Fraction(3, 4)
+   let zero // will fill in during a test
+   const one = new Fraction(1)
+
+   it('supports typeOf', () => {
+      assert.strictEqual(math.typeOf(half), 'Fraction')
+   })
+
+   it('can be built', () => {
+      zero = math.fraction()
+      assert.deepStrictEqual(zero, new Fraction(0))
+      assert.deepStrictEqual(math.fraction(1/2), half)
+      assert.strictEqual(math.fraction(half), half) // maybe it should be a clone?
+      assert.strictEqual(math.fraction(9, 16).valueOf(), 9/16)
+      assert.strictEqual(math.fraction(9n, 16n).valueOf(), 9/16)
+   })
+
+   it('has abs and sign', () => {
+      assert.deepStrictEqual(math.abs(math.fraction('-1/2')), half)
+      assert.deepStrictEqual(math.sign(math.negate(tf)), math.negate(one))
+   })
+
+   it('can add and multiply', () => {
+      assert.strictEqual(math.add(half, 1).valueOf(), 1.5)
+      assert.strictEqual(math.multiply(2, half).valueOf(), 1)
+   })
+
+   it('can subtract and divide', () => {
+      assert.strictEqual(math.subtract(half,tf).valueOf(), -0.25)
+      assert.strictEqual(math.divide(tf,half).valueOf(), 1.5)
+   })
+
+   it('computes mod', () => {
+      assert.strictEqual(math.mod(tf, half).valueOf(), 0.25)
+      assert.strictEqual(math.mod(tf, math.negate(half)).valueOf(), 0.25)
+      assert.strictEqual(math.mod(math.negate(tf), half).valueOf(), 0.25)
+      assert.strictEqual(
+         math.mod(math.negate(tf), math.negate(half)).valueOf(),
+         0.25)
+      assert.deepStrictEqual(
+         math.mod(math.fraction(-1, 3), half),
+         math.fraction(1, 6))
+   })
+
+   it('supports conjugate', () => {
+      assert.strictEqual(math.conjugate(half), half)
+   })
+
+   it('can compare fractions', () => {
+      assert.deepStrictEqual(math.compare(tf, half), one)
+      assert.strictEqual(math.equal(half, math.fraction("2/4")), true)
+      assert.strictEqual(math.smaller(half, tf), true)
+      assert.strictEqual(math.larger(half, tf), false)
+      assert.strictEqual(math.smallerEq(tf, math.fraction(0.75)), true)
+      assert.strictEqual(math.largerEq(tf, half), true)
+      assert.strictEqual(math.unequal(half, tf), true)
+      assert.strictEqual(math.isZero(math.zero(tf)), true)
+      assert.strictEqual(math.isZero(half), false)
+   })
+
+   it('computes gcd and lcm', () => {
+      assert.strictEqual(math.gcd(half,tf).valueOf(), 0.25)
+      assert.strictEqual(math.lcm(half,tf).valueOf(), 1.5)
+   })
+
+   it('computes additive and multiplicative inverses', () => {
+      assert.strictEqual(math.negate(half).valueOf(), -0.5)
+      assert.deepStrictEqual(math.invert(tf), math.fraction('4/3'))
+   })
+
+   it('computes integer parts and quotients', () => {
+      assert.deepStrictEqual(math.floor(tf), zero)
+      assert.deepStrictEqual(math.round(tf), one)
+      assert.deepStrictEqual(math.ceiling(half), one)
+      assert.deepStrictEqual(math.quotient(tf, half), one)
+      assert.deepStrictEqual(
+         math.roundquotient(math.fraction(7/8), half),
+         math.multiply(2,math.one(tf)))
+   })
+
+   it('has no sqrt (although that should be patched)', () => {
+      assert.throws(() => math.sqrt(math.fraction(9/16)), TypeError)
+   })
+
+   it('but it can square', () => {
+      assert.deepStrictEqual(math.square(tf), math.fraction(9/16))
+   })
+
+   it('knows the types of its operations', () => {
+      assert.deepStrictEqual(
+         math.returnTypeOf('ceiling', 'Fraction'), 'Fraction')
+      assert.deepStrictEqual(
+         math.returnTypeOf('multiply', 'Fraction,Fraction'), 'Fraction')
+   })
+})

test/tuple/_native.mjs

@@ -8,14 +8,12 @@ describe('tuple', () => {
 
    it('does not allow unification by converting consecutive arguments', () => {
       assert.throws(() => math.tuple(3, 5.2, 2n), /TypeError.*unif/)
-      // Hence, the order matters in a slightly unfortunate way,
-      // but I think being a little ragged in these edge cases is OK:
       assert.throws(
          () => math.tuple(3, 2n, math.complex(5.2)),
          /TypeError.*unif/)
-      assert.deepStrictEqual(
-         math.tuple(3, math.complex(2n), 5.2),
-         {elts: [math.complex(3), math.complex(2n), math.complex(5.2)]})
+      assert.throws(
+         () => math.tuple(3, math.complex(2n), 5.2),
+         /TypeError.*unif/)
    })
 
    it('can be tested for zero and equality', () => {
@@ -56,6 +54,9 @@ describe('tuple', () => {
       assert.deepStrictEqual(
          math.subtract(math.tuple(3n,4n,5n), math.tuple(2n,1n,0n)),
          math.tuple(1n,3n,5n))
+      assert.deepStrictEqual(
+         math.returnTypeOf('subtract', 'Tuple<bigint>,Tuple<bigint>'),
+         'Tuple<bigint>')
       assert.throws(
          () => math.subtract(math.tuple(5,6), math.tuple(7)),
          /RangeError/)
@@ -106,9 +107,16 @@ describe('tuple', () => {
    })
 
    it('supports sqrt', () => {
+      const mixedTuple = math.tuple(2, math.complex(0,2), 1.5)
       assert.deepStrictEqual(
-         math.sqrt(math.tuple(4,-4,2.25)),
-         math.tuple(2, math.complex(0,2), 1.5))
+         mixedTuple,
+         math.tuple(math.complex(2), math.complex(0,2), math.complex(1.5)))
+      assert.strictEqual(
+         math.returnTypeOf('tuple', 'NumInt, Complex<NumInt>, number'),
+         'Tuple<Complex<number>>')
+      assert.deepStrictEqual(math.sqrt(math.tuple(4,-4,2.25)), mixedTuple)
+      assert.strictEqual(
+         math.returnTypeOf('sqrt', 'Tuple<NumInt>'), 'Tuple<Complex<number>>')
    })
 
 })

tools/approx.mjs

@@ -0,0 +1,46 @@
+import assert from 'assert'
+
+export const epsilon = 1e-12
+
+const isNumber = entity => (typeof entity === 'number')
+
+export function equal(a, b) {
+   if (isNumber(a) && isNumber(b)) {
+      if (a === b) return true
+      if (isNaN(a)) return assert.strictEqual(a.toString(), b.toString())
+      const message = `${a} ~= ${b} (to ${epsilon})`
+      if (a === 0) return assert.ok(Math.abs(b) < epsilon, message)
+      if (b === 0) return assert.ok(Math.abs(a) < epsilon, message)
+      const diff = Math.abs(a - b)
+      const maxDiff = Math.abs(epsilon * Math.max(Math.abs(a), Math.abs(b)))
+      return assert.ok(diff <= maxDiff, message)
+   }
+   return assert.strictEqual(a, b)
+}
+
+export function deepEqual(a, b) {
+   if (Array.isArray(a) && Array.isArray(b)) {
+      const alen = a.length
+      assert.strictEqual(alen, b.length, `${a} ~= ${b}`)
+      for (let i = 0; i < alen; ++i) deepEqual(a[i], b[i])
+      return true
+   }
+   if (typeof a === 'object' && typeof b === 'object') {
+      for (const prop in a) {
+         if (a.hasOwnProperty(prop)) {
+            assert.ok(
+               b.hasOwnProperty(prop), `a[${prop}] = ${a[prop]} ~= ${b[prop]}`)
+            deepEqual(a[prop], b[prop])
+         }
+      }
+
+      for (const prop in b) {
+         if (b.hasOwnProperty(prop)) {
+            assert.ok(
+               a.hasOwnProperty(prop), `${a[prop]} ~= ${b[prop]} = b[${prop}]`)
+         }
+      }
+      return true
+   }
+   return equal(a, b)
+}