🚀 Socket Launch Week Day 5:Introducing Repository Access Permissions and Custom Roles.Learn more
Sign In

@rec-math/math

Package Overview
Dependencies
Maintainers
1
Versions
4
Alerts
File Explorer

Advanced tools

Socket logo

Install Socket

Detect and block malicious and high-risk dependencies

Install

@rec-math/math - npm Package Compare versions

Comparing version
1.2.0
to
1.3.0
+9
dist/rec-math-numerical.min.js
/*! @rec-math/math-numerical v1.3.0 2022-06-28 23:54:20
*
* This file includes only the RecMath.numerical module.
*
* https://github.com/rec-math/ts-math#readme
* Copyright pbuk (https://github.com/pb-uk) MIT license.
*/
!function(t){"use strict";const e={epsilon:16*Number.EPSILON,maxDepth:1/0},r=(t,r,n,i,a={})=>{const o=Object.assign(Object.assign({},e),a),{epsilon:c,maxDepth:h}=o,p={steps:0,errorEstimate:0,depth:1},[u,l,b]=((t,e,r)=>{if(!isFinite(r))return isFinite(e)?[r=>{const s=1/(1-r);return t(e+r*s)*s*s},0,1]:[e=>{const r=e*e,s=1/(1-r);return t(e*s)*(1+r)*s*s},-1,1];if(!isFinite(e))return[e=>t(r-(1-e)/e)/(e*e),0,1];return[t,e,r]})(r,n,i);let[f,g]=s(t,u,l,b,1,1,0,Object.assign({},p));const d=Math.abs(c*f/(b-l));return[f,g]=s(t,u,l,b,1,h,d,p),[f,Object.assign(Object.assign({},p),{errorEstimate:g})]},s=(t,e,r,n,i,a,o,c)=>{const[h,p]=t(e,r,n),u=Math.abs(p-h);if(i>=a)return++c.steps,[h,u];if(u<=o*(n-r))return++c.steps,[h,u];const l=(r+n)/2;if(r>=l||l>=n){c.isUnreliable=!0;const t=isNaN(u)?0:u;return isNaN(h)||Math.abs(h)===1/0?[0,t]:[h,t]}++i>c.depth&&(c.depth=i);const[b,f]=s(t,e,r,l,i,a,o,c),[g,d]=s(t,e,l,n,i,a,o,c);return[b+g,f+d]},n=(t,e,r)=>{const s=(r-e)/2,n=(e+r)/2,i=.20778495500789848*s,a=.4058451513773972*s,o=.5860872354676911*s,c=.7415311855993945*s,h=.8648644233597691*s,p=.9491079123427585*s,u=.9914553711208126*s,l=.4179591836734694*s,b=.3818300505051189*s,f=.27970539148927664*s,g=.1294849661688697*s,d=.20948214108472782*s,m=.20443294007529889*s,E=.19035057806478542*s,O=.1690047266392679*s,j=.14065325971552592*s,M=.10479001032225019*s,N=.06309209262997856*s,_=.022935322010529224*s,w=t(n),y=t(n+i),R=t(n-i),v=t(n+a),x=t(n-a),F=t(n+o),A=t(n-o),D=t(n+c),P=t(n-c),q=t(n+h),z=t(n-h),I=t(n+p),L=t(n-p);return[d*w+m*(y+R)+E*(v+x)+O*(F+A)+j*(D+P)+M*(q+z)+N*(I+L)+_*(t(n+u)+t(n-u)),l*w+b*(v+x)+f*(D+P)+g*(I+L)]},i="integration range must be an array of at least two endpoints";var a=Object.freeze({__proto__:null,quad:(t,e,s={})=>{if(!Array.isArray(e))throw new RangeError(i);if(2===e.length){const[i,a]=e;return r(n,t,i,a,s)}if(e.length<2)throw new RangeError(i);let a=0;const o={steps:0,errorEstimate:0,points:[],depth:0};for(let i=0;i<e.length-1;++i){const c=r(n,t,e[i],e[i+1],s);o.points.push(c),a+=c[0],o.steps+=c[1].steps,o.errorEstimate+=c[1].errorEstimate,o.depth=Math.max(o.depth,c[1].depth)}return[a,o]}});t.numerical=a,t.version="1.3.0",Object.defineProperty(t,"__esModule",{value:!0})}(this.RecMath=this.RecMath||{});
//# sourceMappingURL=rec-math-numerical.min.js.map
{"version":3,"file":"rec-math-numerical.min.js","sources":["../src/version.ts","../src/numerical/quad/adaptive-quadrature.ts","../src/numerical/quad/gauss-kronrod-g7k15.ts","../src/numerical/quad/index.ts"],"sourcesContent":["// rec-math/src/version.ts\n\nexport const version = '1.3.0';\n","// rec-math/src/numerical/quad/adaptive-quadrature.ts\n\nimport type {\n IntegrandCallback,\n IntegrationStep,\n QuadratureInfo,\n QuadratureOptions,\n} from '.';\n\n/** Defaults for integration. */\nconst defaults = {\n /**\n * Target error estimate for a step: if this is smaller it can lead to\n * accumulation of roundoff errors.\n */\n epsilon: Number.EPSILON * 16,\n /**\n * Maximum depth \\\\( d_{max} \\\\)for integration: the maximum number of steps\n * can be \\\\( 2^{d_{max}} \\\\).\n * */\n maxDepth: Infinity,\n};\n\n/**\n * Perform a substitution with an appropriate change of variables to deal with\n * infinite ranges.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns An array with any necessary substitution.\n */\nconst changeOfVariables = (\n f: IntegrandCallback,\n a: number,\n b: number,\n): [f: IntegrandCallback, a: number, b: number] => {\n if (!isFinite(b)) {\n if (!isFinite(a)) {\n // Change variables to integrate between - and + infinity.\n const _f = (t: number): number => {\n const tSquared = t * t;\n const oneOverOneMinusTSquared = 1 / (1 - tSquared);\n return (\n f(t * oneOverOneMinusTSquared) *\n (1 + tSquared) *\n oneOverOneMinusTSquared *\n oneOverOneMinusTSquared\n );\n };\n return [_f, -1, 1];\n }\n // Change variables to integrate up to infinity.\n const _f = (t: number): number => {\n const oneOverOneMinusT = 1 / (1 - t);\n return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;\n };\n return [_f, 0, 1];\n }\n\n if (!isFinite(a)) {\n // Change variables to integrate up from negative infinity.\n const _f = (t: number): number => {\n return f(b - (1 - t) / t) / (t * t);\n };\n return [_f, 0, 1];\n }\n return [f, a, b];\n};\n\nexport const integrate = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number,\n b: number,\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Establish settings.\n const settings = { ...defaults, ...options };\n const { epsilon, maxDepth } = settings;\n\n const info = { steps: 0, errorEstimate: 0, depth: 1 };\n\n // Allow for a change of variables to deal with infinite limits.\n const [_f, _a, _b] = changeOfVariables(f, a, b);\n\n // Get estimate so we can work out the acceptable global error.\n // Use a depth of 1 to calculate a 15 point Kronrod quadrature.\n let [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n 1, // Maximum depth.\n 0, // Acceptable error per unit step.\n { ...info }, // Statistics.\n );\n\n // Now calculate using the target global error.\n const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));\n [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n\n return [result, { ...info, errorEstimate }];\n};\n\nconst integrate_part = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number, // Lower limit.\n b: number, // Upper limit.\n depth: number, // New depth.\n maxDepth: number, // Maximum depth.\n acceptableUnitError: number, // Acceptable error per unit step.\n info: QuadratureInfo, // Statistics.\n): [a: number, b: number] => {\n // Initialize things.\n const [currentEstimate, poorEstimate] = integrationStep(\n f, // Integrand.\n a, // Lower limit.\n b, // Upper limit.\n );\n\n const errorEstimate = Math.abs(poorEstimate - currentEstimate);\n\n if (depth >= maxDepth) {\n // Reached the maximum allowable depth so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const acceptableError = acceptableUnitError * (b - a);\n if (errorEstimate <= acceptableError) {\n // Error is acceptable for the size of step so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const mid = (a + b) / 2;\n if (a >= mid || mid >= b) {\n // We can't make this step any smaller: looks like a discontinuity.\n info.isUnreliable = true;\n const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;\n if (isNaN(currentEstimate) || Math.abs(currentEstimate) === Infinity) {\n return [0, safeErrorEstimate];\n }\n return [currentEstimate, safeErrorEstimate];\n }\n\n // Recurse deeper.\n ++depth;\n if (depth > info.depth) {\n info.depth = depth;\n }\n const [leftEstimate, leftErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n a, // Lower limit.\n mid, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n const [rightEstimate, rightErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n mid, // Lower limit.\n b, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n return [leftEstimate + rightEstimate, leftErrorEstimate + rightErrorEstimate];\n};\n","// rec-math/src/numerical/quad/gauss-kronrod-g7k15.ts\n\nimport type { IntegrandCallback, IntegrationStep } from '.';\n\n// Gauss-Kronrod constants (G7, K15) on [-1, 1].\n// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/\n\n// node_g0k0 = 0;\n// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;\nconst node_g1k2 = 0.4058451513773972;\n// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;\nconst node_g2k4 = 0.7415311855993945;\n// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;\nconst node_g3k6 = 0.9491079123427585;\n\n// node_k1 = 2.077849550078984676006894037732449e-01 exact;\nconst node_k1 = 0.20778495500789848;\n// node_k3 = 5.860872354676911302941448382587296e-01 exact;\nconst node_k3 = 0.5860872354676911;\n// node_k5 = 8.648644233597690727897127886409262e-01 exact;\nconst node_k5 = 0.8648644233597691;\n// node_k7 = 9.914553711208126392068546975263285e-01 exact;\nconst node_k7 = 0.9914553711208126;\n\n// weight_g0 = 4.179591836734693877551020408163265e-01 exact;\nconst weight_g0 = 0.4179591836734694;\n// weight_g1 = 3.818300505051189449503697754889751e-01 exact;\nconst weight_g1 = 0.3818300505051189;\n// weight_g2 = 2.797053914892766679014677714237796e-01 exact;\nconst weight_g2 = 0.27970539148927664;\n// weight_g3 = 1.294849661688696932706114326790820e-01 exact;\nconst weight_g3 = 0.1294849661688697;\n\n// weight_k0 = 2.094821410847278280129991748917143e-01 exact;\nconst weight_k0 = 0.20948214108472782;\n// weight_k1 = 2.044329400752988924141619992346491e-01 exact;\nconst weight_k1 = 0.20443294007529889;\n// weight_k2 = 1.903505780647854099132564024210137e-01 exact;\nconst weight_k2 = 0.19035057806478542;\n// weight_k3 = 1.690047266392679028265834265985503e-01 exact;\nconst weight_k3 = 0.1690047266392679;\n// weight_k4 = 1.406532597155259187451895905102379e-01 exact;\nconst weight_k4 = 0.14065325971552592;\n// weight_k5 = 1.047900103222501838398763225415180e-01 exact;\nconst weight_k5 = 0.10479001032225019;\n// weight_k6 = 6.309209262997855329070066318920429e-02 exact;\nconst weight_k6 = 0.06309209262997856;\n// weight_k7 = 2.293532201052922496373200805896959e-02 exact;\nconst weight_k7 = 0.022935322010529224;\n\n/**\n * Perform a single integration step.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns [current best estimate, poor estimate]\n */\nexport const integrationStep: IntegrationStep = (\n f: IntegrandCallback, // Integrand.\n a: number, // Lower limit.\n b: number, // Upper limit.\n): [a: number, b: number] => {\n // Gauss-Kronrod (G7, K15).\n const scale = (b - a) / 2;\n const offset = (a + b) / 2;\n\n // const scaled_g0k0 = 0;\n const scaled_k1 = node_k1 * scale;\n const scaled_g1k2 = node_g1k2 * scale;\n const scaled_k3 = node_k3 * scale;\n const scaled_g2k4 = node_g2k4 * scale;\n const scaled_k5 = node_k5 * scale;\n const scaled_g3k6 = node_g3k6 * scale;\n const scaled_k7 = node_k7 * scale;\n\n const scaled_weight_g0 = weight_g0 * scale;\n const scaled_weight_g1 = weight_g1 * scale;\n const scaled_weight_g2 = weight_g2 * scale;\n const scaled_weight_g3 = weight_g3 * scale;\n\n const scaled_weight_k0 = weight_k0 * scale;\n const scaled_weight_k1 = weight_k1 * scale;\n const scaled_weight_k2 = weight_k2 * scale;\n const scaled_weight_k3 = weight_k3 * scale;\n const scaled_weight_k4 = weight_k4 * scale;\n const scaled_weight_k5 = weight_k5 * scale;\n const scaled_weight_k6 = weight_k6 * scale;\n const scaled_weight_k7 = weight_k7 * scale;\n\n const f_g0k0 = f(offset);\n\n const f_k1_h = f(offset + scaled_k1);\n const f_k1_l = f(offset - scaled_k1);\n const f_g1k2_h = f(offset + scaled_g1k2);\n const f_g1k2_l = f(offset - scaled_g1k2);\n\n const f_k3_h = f(offset + scaled_k3);\n const f_k3_l = f(offset - scaled_k3);\n const f_g2k4_h = f(offset + scaled_g2k4);\n const f_g2k4_l = f(offset - scaled_g2k4);\n\n const f_k5_h = f(offset + scaled_k5);\n const f_k5_l = f(offset - scaled_k5);\n const f_g3k6_h = f(offset + scaled_g3k6);\n const f_g3k6_l = f(offset - scaled_g3k6);\n\n const f_k7_h = f(offset + scaled_k7);\n const f_k7_l = f(offset - scaled_k7);\n\n const poorEstimate =\n scaled_weight_g0 * f_g0k0 + // g0\n scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1\n scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2\n scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3\n\n const currentEstimate =\n scaled_weight_k0 * f_g0k0 + // k0\n scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1\n scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2\n scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3\n scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4\n scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5\n scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6\n scaled_weight_k7 * (f_k7_h + f_k7_l); // k7\n\n return [currentEstimate, poorEstimate];\n};\n","// rec-math/src/numerical/quad/index.ts\n\nimport { integrate } from './adaptive-quadrature.js';\nimport { integrationStep } from './gauss-kronrod-g7k15.js';\n\nexport type IntegrandCallback = (a: number) => number;\n\nexport interface QuadratureInfo {\n /** Number of steps _used_ in the integration. */\n steps: number;\n /** Estimate of the global error. */\n errorEstimate: number;\n /**\n * Set to true if there was a problem (e.g. could not reach required accuracy\n * with a step at machine precision).\n */\n isUnreliable?: true;\n /**\n * If the range has multiple (i.e. more than 2) points, contains the results\n * for intermediate ranges.\n */\n points?: [r: number, i: QuadratureInfo][];\n /** The maximum depth used. */\n depth: number;\n}\n\nexport interface QuadratureOptions {\n /** Used to control global error. */\n epsilon?: number;\n /** Maximum depth to use for adaptive step length. */\n maxDepth?: number;\n}\n\nexport type IntegrationStep = (\n /** The integrand. */\n f: (a: number) => number,\n a: number, // Lower limit.\n b: number, // Upper limit.\n) => [a: number, b: number];\n\nconst rangeErrorMessage =\n 'integration range must be an array of at least two endpoints';\n\n/**\n * Numerically compute a definite integral.\n *\n * @param f Callback returning value of integrand.\n * @param range A range of at least 2 endpoints.\n * @param options Options for the computation.\n *\n * @returns The results of the computation.\n */\nexport const quad = (\n f: IntegrandCallback,\n range: number[],\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Interpret the range.\n if (!Array.isArray(range)) {\n throw new RangeError(rangeErrorMessage);\n }\n\n if (range.length === 2) {\n // Integrate over a single range.\n const [a, b] = range;\n return integrate(integrationStep, f, a, b, options);\n }\n\n if (range.length < 2) {\n // Can't integrate at a point!\n throw new RangeError(rangeErrorMessage);\n }\n\n // Integrate over multiple ranges.\n let result = 0;\n const points: [r: number, i: QuadratureInfo][] = [];\n const info = {\n steps: 0,\n errorEstimate: 0,\n points,\n depth: 0,\n };\n\n for (let i = 0; i < range.length - 1; ++i) {\n const single = integrate(\n integrationStep,\n f,\n range[i],\n range[i + 1],\n options,\n );\n\n info.points.push(single);\n // Update the result.\n result += single[0];\n // Update the cumulative statistics.\n info.steps += single[1].steps;\n info.errorEstimate += single[1].errorEstimate;\n info.depth = Math.max(info.depth, single[1].depth);\n }\n\n return [result, info];\n};\n"],"names":["defaults","epsilon","Number","EPSILON","maxDepth","Infinity","integrate","integrationStep","f","a","b","options","settings","Object","assign","info","steps","errorEstimate","depth","_f","_a","_b","isFinite","t","oneOverOneMinusT","tSquared","oneOverOneMinusTSquared","changeOfVariables","result","integrate_part","acceptableUnitError","Math","abs","currentEstimate","poorEstimate","mid","isUnreliable","safeErrorEstimate","isNaN","leftEstimate","leftErrorEstimate","rightEstimate","rightErrorEstimate","scale","offset","scaled_k1","scaled_g1k2","scaled_k3","scaled_g2k4","scaled_k5","scaled_g3k6","scaled_k7","scaled_weight_g0","scaled_weight_g1","scaled_weight_g2","scaled_weight_g3","scaled_weight_k0","scaled_weight_k1","scaled_weight_k2","scaled_weight_k3","scaled_weight_k4","scaled_weight_k5","scaled_weight_k6","scaled_weight_k7","f_g0k0","f_k1_h","f_k1_l","f_g1k2_h","f_g1k2_l","f_k3_h","f_k3_l","f_g2k4_h","f_g2k4_l","f_k5_h","f_k5_l","f_g3k6_h","f_g3k6_l","rangeErrorMessage","range","Array","isArray","RangeError","length","points","i","single","push","max"],"mappings":";;;;;;;0BAEO,MCQDA,EAAW,CAKfC,QAA0B,GAAjBC,OAAOC,QAKhBC,SAAUC,KAkDCC,EAAY,CACvBC,EACAC,EACAC,EACAC,EACAC,EAA6B,MAG7B,MAAMC,EAAgBC,OAAAC,OAAAD,OAAAC,OAAA,GAAAd,GAAaW,IAC7BV,QAAEA,EAAOG,SAAEA,GAAaQ,EAExBG,EAAO,CAAEC,MAAO,EAAGC,cAAe,EAAGC,MAAO,IAG3CC,EAAIC,EAAIC,GApDS,EACxBb,EACAC,EACAC,KAEA,IAAKY,SAASZ,GACZ,OAAKY,SAASb,GAmBP,CAJKc,IACV,MAAMC,EAAmB,GAAK,EAAID,GAClC,OAAOf,EAAEC,EAAIc,EAAIC,GAAoBA,EAAmBA,GAE9C,EAAG,GAPN,CAVKD,IACV,MAAME,EAAWF,EAAIA,EACfG,EAA0B,GAAK,EAAID,GACzC,OACEjB,EAAEe,EAAIG,IACL,EAAID,GACLC,EACAA,IAGS,EAAG,GAUpB,IAAKJ,SAASb,GAKZ,MAAO,CAHKc,GACHf,EAAEE,GAAK,EAAIa,GAAKA,IAAMA,EAAIA,GAEvB,EAAG,GAEjB,MAAO,CAACf,EAAGC,EAAGC,IAiBOiB,CAAkBnB,EAAGC,EAAGC,GAI7C,IAAKkB,EAAQX,GAAiBY,EAC5BtB,EACAY,EACAC,EACAC,EACA,EACA,EACA,EAACR,OAAAC,OAAA,GACIC,IAIP,MAAMe,EAAsBC,KAAKC,IAAK/B,EAAU2B,GAAWP,EAAKD,IAYhE,OAXCQ,EAAQX,GAAiBY,EACxBtB,EACAY,EACAC,EACAC,EACA,EACAjB,EACA0B,EACAf,GAGK,CAACa,EAAMf,OAAAC,OAAAD,OAAAC,OAAA,GAAOC,GAAM,CAAAE,oBAGvBY,EAAiB,CACrBtB,EACAC,EACAC,EACAC,EACAQ,EACAd,EACA0B,EACAf,KAGA,MAAOkB,EAAiBC,GAAgB3B,EACtCC,EACAC,EACAC,GAGIO,EAAgBc,KAAKC,IAAIE,EAAeD,GAE9C,GAAIf,GAASd,EAGX,QADEW,EAAKC,MACA,CAACiB,EAAiBhB,GAI3B,GAAIA,GADoBa,GAAuBpB,EAAID,GAIjD,QADEM,EAAKC,MACA,CAACiB,EAAiBhB,GAG3B,MAAMkB,GAAO1B,EAAIC,GAAK,EACtB,GAAID,GAAK0B,GAAOA,GAAOzB,EAAG,CAExBK,EAAKqB,cAAe,EACpB,MAAMC,EAAoBC,MAAMrB,GAAiB,EAAIA,EACrD,OAAIqB,MAAML,IAAoBF,KAAKC,IAAIC,KAAqB5B,IACnD,CAAC,EAAGgC,GAEN,CAACJ,EAAiBI,KAIzBnB,EACUH,EAAKG,QACfH,EAAKG,MAAQA,GAEf,MAAOqB,EAAcC,GAAqBX,EACxCtB,EACAC,EACAC,EACA0B,EACAjB,EACAd,EACA0B,EACAf,IAEK0B,EAAeC,GAAsBb,EAC1CtB,EACAC,EACA2B,EACAzB,EACAQ,EACAd,EACA0B,EACAf,GAEF,MAAO,CAACwB,EAAeE,EAAeD,EAAoBE,IC7H/CnC,EAAmC,CAC9CC,EACAC,EACAC,KAGA,MAAMiC,GAASjC,EAAID,GAAK,EAClBmC,GAAUnC,EAAIC,GAAK,EAGnBmC,EApDQ,mBAoDcF,EACtBG,EA5DU,kBA4DgBH,EAC1BI,EApDQ,kBAoDcJ,EACtBK,EA5DU,kBA4DgBL,EAC1BM,EApDQ,kBAoDcN,EACtBO,EA5DU,kBA4DgBP,EAC1BQ,EApDQ,kBAoDcR,EAEtBS,EAnDU,kBAmDqBT,EAC/BU,EAlDU,kBAkDqBV,EAC/BW,EAjDU,mBAiDqBX,EAC/BY,EAhDU,kBAgDqBZ,EAE/Ba,EA/CU,mBA+CqBb,EAC/Bc,EA9CU,mBA8CqBd,EAC/Be,EA7CU,mBA6CqBf,EAC/BgB,EA5CU,kBA4CqBhB,EAC/BiB,EA3CU,mBA2CqBjB,EAC/BkB,EA1CU,mBA0CqBlB,EAC/BmB,EAzCU,mBAyCqBnB,EAC/BoB,EAxCU,oBAwCqBpB,EAE/BqB,EAASxD,EAAEoC,GAEXqB,EAASzD,EAAEoC,EAASC,GACpBqB,EAAS1D,EAAEoC,EAASC,GACpBsB,EAAW3D,EAAEoC,EAASE,GACtBsB,EAAW5D,EAAEoC,EAASE,GAEtBuB,EAAS7D,EAAEoC,EAASG,GACpBuB,EAAS9D,EAAEoC,EAASG,GACpBwB,EAAW/D,EAAEoC,EAASI,GACtBwB,EAAWhE,EAAEoC,EAASI,GAEtByB,EAASjE,EAAEoC,EAASK,GACpByB,EAASlE,EAAEoC,EAASK,GACpB0B,EAAWnE,EAAEoC,EAASM,GACtB0B,EAAWpE,EAAEoC,EAASM,GAqB5B,MAAO,CATLM,EAAmBQ,EACnBP,GAAoBQ,EAASC,GAC7BR,GAAoBS,EAAWC,GAC/BT,GAAoBU,EAASC,GAC7BV,GAAoBW,EAAWC,GAC/BX,GAAoBY,EAASC,GAC7BZ,GAAoBa,EAAWC,GAC/Bb,GAjBavD,EAAEoC,EAASO,GACX3C,EAAEoC,EAASO,IAGxBC,EAAmBY,EACnBX,GAAoBc,EAAWC,GAC/Bd,GAAoBiB,EAAWC,GAC/BjB,GAAoBoB,EAAWC,KC1E7BC,EACJ,wGAWkB,CAClBrE,EACAsE,EACAnE,EAA6B,MAG7B,IAAKoE,MAAMC,QAAQF,GACjB,MAAM,IAAIG,WAAWJ,GAGvB,GAAqB,IAAjBC,EAAMI,OAAc,CAEtB,MAAOzE,EAAGC,GAAKoE,EACf,OAAOxE,EAAUC,EAAiBC,EAAGC,EAAGC,EAAGC,GAG7C,GAAImE,EAAMI,OAAS,EAEjB,MAAM,IAAID,WAAWJ,GAIvB,IAAIjD,EAAS,EACb,MACMb,EAAO,CACXC,MAAO,EACPC,cAAe,EACfkE,OAJ+C,GAK/CjE,MAAO,GAGT,IAAK,IAAIkE,EAAI,EAAGA,EAAIN,EAAMI,OAAS,IAAKE,EAAG,CACzC,MAAMC,EAAS/E,EACbC,EACAC,EACAsE,EAAMM,GACNN,EAAMM,EAAI,GACVzE,GAGFI,EAAKoE,OAAOG,KAAKD,GAEjBzD,GAAUyD,EAAO,GAEjBtE,EAAKC,OAASqE,EAAO,GAAGrE,MACxBD,EAAKE,eAAiBoE,EAAO,GAAGpE,cAChCF,EAAKG,MAAQa,KAAKwD,IAAIxE,EAAKG,MAAOmE,EAAO,GAAGnE,OAG9C,MAAO,CAACU,EAAQb,8BHnGK"}
{"version":3,"file":"index.js","sources":["../src/version.ts","../src/numerical/quad/adaptive-quadrature.ts","../src/numerical/quad/gauss-kronrod-g7k15.ts","../src/numerical/quad/index.ts","../src/numerical/index.ts"],"sourcesContent":["// rec-math/src/version.ts\n\nexport const version = '1.3.0';\n","// rec-math/src/numerical/quad/adaptive-quadrature.ts\n\nimport type {\n IntegrandCallback,\n IntegrationStep,\n QuadratureInfo,\n QuadratureOptions,\n} from '.';\n\n/** Defaults for integration. */\nconst defaults = {\n /**\n * Target error estimate for a step: if this is smaller it can lead to\n * accumulation of roundoff errors.\n */\n epsilon: Number.EPSILON * 16,\n /**\n * Maximum depth \\\\( d_{max} \\\\)for integration: the maximum number of steps\n * can be \\\\( 2^{d_{max}} \\\\).\n * */\n maxDepth: Infinity,\n};\n\n/**\n * Perform a substitution with an appropriate change of variables to deal with\n * infinite ranges.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns An array with any necessary substitution.\n */\nconst changeOfVariables = (\n f: IntegrandCallback,\n a: number,\n b: number,\n): [f: IntegrandCallback, a: number, b: number] => {\n if (!isFinite(b)) {\n if (!isFinite(a)) {\n // Change variables to integrate between - and + infinity.\n const _f = (t: number): number => {\n const tSquared = t * t;\n const oneOverOneMinusTSquared = 1 / (1 - tSquared);\n return (\n f(t * oneOverOneMinusTSquared) *\n (1 + tSquared) *\n oneOverOneMinusTSquared *\n oneOverOneMinusTSquared\n );\n };\n return [_f, -1, 1];\n }\n // Change variables to integrate up to infinity.\n const _f = (t: number): number => {\n const oneOverOneMinusT = 1 / (1 - t);\n return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;\n };\n return [_f, 0, 1];\n }\n\n if (!isFinite(a)) {\n // Change variables to integrate up from negative infinity.\n const _f = (t: number): number => {\n return f(b - (1 - t) / t) / (t * t);\n };\n return [_f, 0, 1];\n }\n return [f, a, b];\n};\n\nexport const integrate = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number,\n b: number,\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Establish settings.\n const settings = { ...defaults, ...options };\n const { epsilon, maxDepth } = settings;\n\n const info = { steps: 0, errorEstimate: 0, depth: 1 };\n\n // Allow for a change of variables to deal with infinite limits.\n const [_f, _a, _b] = changeOfVariables(f, a, b);\n\n // Get estimate so we can work out the acceptable global error.\n // Use a depth of 1 to calculate a 15 point Kronrod quadrature.\n let [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n 1, // Maximum depth.\n 0, // Acceptable error per unit step.\n { ...info }, // Statistics.\n );\n\n // Now calculate using the target global error.\n const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));\n [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n\n return [result, { ...info, errorEstimate }];\n};\n\nconst integrate_part = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number, // Lower limit.\n b: number, // Upper limit.\n depth: number, // New depth.\n maxDepth: number, // Maximum depth.\n acceptableUnitError: number, // Acceptable error per unit step.\n info: QuadratureInfo, // Statistics.\n): [a: number, b: number] => {\n // Initialize things.\n const [currentEstimate, poorEstimate] = integrationStep(\n f, // Integrand.\n a, // Lower limit.\n b, // Upper limit.\n );\n\n const errorEstimate = Math.abs(poorEstimate - currentEstimate);\n\n if (depth >= maxDepth) {\n // Reached the maximum allowable depth so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const acceptableError = acceptableUnitError * (b - a);\n if (errorEstimate <= acceptableError) {\n // Error is acceptable for the size of step so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const mid = (a + b) / 2;\n if (a >= mid || mid >= b) {\n // We can't make this step any smaller: looks like a discontinuity.\n info.isUnreliable = true;\n const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;\n if (isNaN(currentEstimate) || Math.abs(currentEstimate) === Infinity) {\n return [0, safeErrorEstimate];\n }\n return [currentEstimate, safeErrorEstimate];\n }\n\n // Recurse deeper.\n ++depth;\n if (depth > info.depth) {\n info.depth = depth;\n }\n const [leftEstimate, leftErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n a, // Lower limit.\n mid, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n const [rightEstimate, rightErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n mid, // Lower limit.\n b, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n return [leftEstimate + rightEstimate, leftErrorEstimate + rightErrorEstimate];\n};\n","// rec-math/src/numerical/quad/gauss-kronrod-g7k15.ts\n\nimport type { IntegrandCallback, IntegrationStep } from '.';\n\n// Gauss-Kronrod constants (G7, K15) on [-1, 1].\n// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/\n\n// node_g0k0 = 0;\n// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;\nconst node_g1k2 = 0.4058451513773972;\n// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;\nconst node_g2k4 = 0.7415311855993945;\n// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;\nconst node_g3k6 = 0.9491079123427585;\n\n// node_k1 = 2.077849550078984676006894037732449e-01 exact;\nconst node_k1 = 0.20778495500789848;\n// node_k3 = 5.860872354676911302941448382587296e-01 exact;\nconst node_k3 = 0.5860872354676911;\n// node_k5 = 8.648644233597690727897127886409262e-01 exact;\nconst node_k5 = 0.8648644233597691;\n// node_k7 = 9.914553711208126392068546975263285e-01 exact;\nconst node_k7 = 0.9914553711208126;\n\n// weight_g0 = 4.179591836734693877551020408163265e-01 exact;\nconst weight_g0 = 0.4179591836734694;\n// weight_g1 = 3.818300505051189449503697754889751e-01 exact;\nconst weight_g1 = 0.3818300505051189;\n// weight_g2 = 2.797053914892766679014677714237796e-01 exact;\nconst weight_g2 = 0.27970539148927664;\n// weight_g3 = 1.294849661688696932706114326790820e-01 exact;\nconst weight_g3 = 0.1294849661688697;\n\n// weight_k0 = 2.094821410847278280129991748917143e-01 exact;\nconst weight_k0 = 0.20948214108472782;\n// weight_k1 = 2.044329400752988924141619992346491e-01 exact;\nconst weight_k1 = 0.20443294007529889;\n// weight_k2 = 1.903505780647854099132564024210137e-01 exact;\nconst weight_k2 = 0.19035057806478542;\n// weight_k3 = 1.690047266392679028265834265985503e-01 exact;\nconst weight_k3 = 0.1690047266392679;\n// weight_k4 = 1.406532597155259187451895905102379e-01 exact;\nconst weight_k4 = 0.14065325971552592;\n// weight_k5 = 1.047900103222501838398763225415180e-01 exact;\nconst weight_k5 = 0.10479001032225019;\n// weight_k6 = 6.309209262997855329070066318920429e-02 exact;\nconst weight_k6 = 0.06309209262997856;\n// weight_k7 = 2.293532201052922496373200805896959e-02 exact;\nconst weight_k7 = 0.022935322010529224;\n\n/**\n * Perform a single integration step.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns [current best estimate, poor estimate]\n */\nexport const integrationStep: IntegrationStep = (\n f: IntegrandCallback, // Integrand.\n a: number, // Lower limit.\n b: number, // Upper limit.\n): [a: number, b: number] => {\n // Gauss-Kronrod (G7, K15).\n const scale = (b - a) / 2;\n const offset = (a + b) / 2;\n\n // const scaled_g0k0 = 0;\n const scaled_k1 = node_k1 * scale;\n const scaled_g1k2 = node_g1k2 * scale;\n const scaled_k3 = node_k3 * scale;\n const scaled_g2k4 = node_g2k4 * scale;\n const scaled_k5 = node_k5 * scale;\n const scaled_g3k6 = node_g3k6 * scale;\n const scaled_k7 = node_k7 * scale;\n\n const scaled_weight_g0 = weight_g0 * scale;\n const scaled_weight_g1 = weight_g1 * scale;\n const scaled_weight_g2 = weight_g2 * scale;\n const scaled_weight_g3 = weight_g3 * scale;\n\n const scaled_weight_k0 = weight_k0 * scale;\n const scaled_weight_k1 = weight_k1 * scale;\n const scaled_weight_k2 = weight_k2 * scale;\n const scaled_weight_k3 = weight_k3 * scale;\n const scaled_weight_k4 = weight_k4 * scale;\n const scaled_weight_k5 = weight_k5 * scale;\n const scaled_weight_k6 = weight_k6 * scale;\n const scaled_weight_k7 = weight_k7 * scale;\n\n const f_g0k0 = f(offset);\n\n const f_k1_h = f(offset + scaled_k1);\n const f_k1_l = f(offset - scaled_k1);\n const f_g1k2_h = f(offset + scaled_g1k2);\n const f_g1k2_l = f(offset - scaled_g1k2);\n\n const f_k3_h = f(offset + scaled_k3);\n const f_k3_l = f(offset - scaled_k3);\n const f_g2k4_h = f(offset + scaled_g2k4);\n const f_g2k4_l = f(offset - scaled_g2k4);\n\n const f_k5_h = f(offset + scaled_k5);\n const f_k5_l = f(offset - scaled_k5);\n const f_g3k6_h = f(offset + scaled_g3k6);\n const f_g3k6_l = f(offset - scaled_g3k6);\n\n const f_k7_h = f(offset + scaled_k7);\n const f_k7_l = f(offset - scaled_k7);\n\n const poorEstimate =\n scaled_weight_g0 * f_g0k0 + // g0\n scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1\n scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2\n scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3\n\n const currentEstimate =\n scaled_weight_k0 * f_g0k0 + // k0\n scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1\n scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2\n scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3\n scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4\n scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5\n scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6\n scaled_weight_k7 * (f_k7_h + f_k7_l); // k7\n\n return [currentEstimate, poorEstimate];\n};\n","// rec-math/src/numerical/quad/index.ts\n\nimport { integrate } from './adaptive-quadrature.js';\nimport { integrationStep } from './gauss-kronrod-g7k15.js';\n\nexport type IntegrandCallback = (a: number) => number;\n\nexport interface QuadratureInfo {\n /** Number of steps _used_ in the integration. */\n steps: number;\n /** Estimate of the global error. */\n errorEstimate: number;\n /**\n * Set to true if there was a problem (e.g. could not reach required accuracy\n * with a step at machine precision).\n */\n isUnreliable?: true;\n /**\n * If the range has multiple (i.e. more than 2) points, contains the results\n * for intermediate ranges.\n */\n points?: [r: number, i: QuadratureInfo][];\n /** The maximum depth used. */\n depth: number;\n}\n\nexport interface QuadratureOptions {\n /** Used to control global error. */\n epsilon?: number;\n /** Maximum depth to use for adaptive step length. */\n maxDepth?: number;\n}\n\nexport type IntegrationStep = (\n /** The integrand. */\n f: (a: number) => number,\n a: number, // Lower limit.\n b: number, // Upper limit.\n) => [a: number, b: number];\n\nconst rangeErrorMessage =\n 'integration range must be an array of at least two endpoints';\n\n/**\n * Numerically compute a definite integral.\n *\n * @param f Callback returning value of integrand.\n * @param range A range of at least 2 endpoints.\n * @param options Options for the computation.\n *\n * @returns The results of the computation.\n */\nexport const quad = (\n f: IntegrandCallback,\n range: number[],\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Interpret the range.\n if (!Array.isArray(range)) {\n throw new RangeError(rangeErrorMessage);\n }\n\n if (range.length === 2) {\n // Integrate over a single range.\n const [a, b] = range;\n return integrate(integrationStep, f, a, b, options);\n }\n\n if (range.length < 2) {\n // Can't integrate at a point!\n throw new RangeError(rangeErrorMessage);\n }\n\n // Integrate over multiple ranges.\n let result = 0;\n const points: [r: number, i: QuadratureInfo][] = [];\n const info = {\n steps: 0,\n errorEstimate: 0,\n points,\n depth: 0,\n };\n\n for (let i = 0; i < range.length - 1; ++i) {\n const single = integrate(\n integrationStep,\n f,\n range[i],\n range[i + 1],\n options,\n );\n\n info.points.push(single);\n // Update the result.\n result += single[0];\n // Update the cumulative statistics.\n info.steps += single[1].steps;\n info.errorEstimate += single[1].errorEstimate;\n info.depth = Math.max(info.depth, single[1].depth);\n }\n\n return [result, info];\n};\n","// rec-math/src/numerical/index.ts\n\nexport { quad } from './quad/index.js';\n"],"names":[],"mappings":";;;;;AAAA;AAEO,MAAM,OAAO,GAAG;;ACFvB;AASA;AACA,MAAM,QAAQ,GAAG;AACf;;;AAGG;AACH,IAAA,OAAO,EAAE,MAAM,CAAC,OAAO,GAAG,EAAE;AAC5B;;;AAGK;AACL,IAAA,QAAQ,EAAE,QAAQ;CACnB,CAAC;AAEF;;;;;;;;AAQG;AACH,MAAM,iBAAiB,GAAG,CACxB,CAAoB,EACpB,CAAS,EACT,CAAS,KACuC;AAChD,IAAA,IAAI,CAAC,QAAQ,CAAC,CAAC,CAAC,EAAE;AAChB,QAAA,IAAI,CAAC,QAAQ,CAAC,CAAC,CAAC,EAAE;;AAEhB,YAAA,MAAM,EAAE,GAAG,CAAC,CAAS,KAAY;AAC/B,gBAAA,MAAM,QAAQ,GAAG,CAAC,GAAG,CAAC,CAAC;gBACvB,MAAM,uBAAuB,GAAG,CAAC,IAAI,CAAC,GAAG,QAAQ,CAAC,CAAC;AACnD,gBAAA,QACE,CAAC,CAAC,CAAC,GAAG,uBAAuB,CAAC;qBAC7B,CAAC,GAAG,QAAQ,CAAC;oBACd,uBAAuB;AACvB,oBAAA,uBAAuB,EACvB;AACJ,aAAC,CAAC;YACF,OAAO,CAAC,EAAE,EAAE,CAAC,CAAC,EAAE,CAAC,CAAC,CAAC;AACpB,SAAA;;AAED,QAAA,MAAM,EAAE,GAAG,CAAC,CAAS,KAAY;YAC/B,MAAM,gBAAgB,GAAG,CAAC,IAAI,CAAC,GAAG,CAAC,CAAC,CAAC;AACrC,YAAA,OAAO,CAAC,CAAC,CAAC,GAAG,CAAC,GAAG,gBAAgB,CAAC,GAAG,gBAAgB,GAAG,gBAAgB,CAAC;AAC3E,SAAC,CAAC;AACF,QAAA,OAAO,CAAC,EAAE,EAAE,CAAC,EAAE,CAAC,CAAC,CAAC;AACnB,KAAA;AAED,IAAA,IAAI,CAAC,QAAQ,CAAC,CAAC,CAAC,EAAE;;AAEhB,QAAA,MAAM,EAAE,GAAG,CAAC,CAAS,KAAY;AAC/B,YAAA,OAAO,CAAC,CAAC,CAAC,GAAG,CAAC,CAAC,GAAG,CAAC,IAAI,CAAC,CAAC,IAAI,CAAC,GAAG,CAAC,CAAC,CAAC;AACtC,SAAC,CAAC;AACF,QAAA,OAAO,CAAC,EAAE,EAAE,CAAC,EAAE,CAAC,CAAC,CAAC;AACnB,KAAA;AACD,IAAA,OAAO,CAAC,CAAC,EAAE,CAAC,EAAE,CAAC,CAAC,CAAC;AACnB,CAAC,CAAC;AAEK,MAAM,SAAS,GAAG,CACvB,eAAgC,EAChC,CAAoB,EACpB,CAAS,EACT,CAAS,EACT,OAA6B,GAAA,EAAE,KACG;;AAElC,IAAA,MAAM,QAAQ,GAAQ,MAAA,CAAA,MAAA,CAAA,MAAA,CAAA,MAAA,CAAA,EAAA,EAAA,QAAQ,CAAK,EAAA,OAAO,CAAE,CAAC;AAC7C,IAAA,MAAM,EAAE,OAAO,EAAE,QAAQ,EAAE,GAAG,QAAQ,CAAC;AAEvC,IAAA,MAAM,IAAI,GAAG,EAAE,KAAK,EAAE,CAAC,EAAE,aAAa,EAAE,CAAC,EAAE,KAAK,EAAE,CAAC,EAAE,CAAC;;AAGtD,IAAA,MAAM,CAAC,EAAE,EAAE,EAAE,EAAE,EAAE,CAAC,GAAG,iBAAiB,CAAC,CAAC,EAAE,CAAC,EAAE,CAAC,CAAC,CAAC;;;AAIhD,IAAA,IAAI,CAAC,MAAM,EAAE,aAAa,CAAC,GAAG,cAAc,CAC1C,eAAe,EACf,EAAE;AACF,IAAA,EAAE;AACF,IAAA,EAAE;AACF,IAAA,CAAC;AACD,IAAA,CAAC;IACD,CAAC,EAAA,MAAA,CAAA,MAAA,CAAA,EAAA,EACI,IAAI,CAAA,CACV,CAAC;;AAGF,IAAA,MAAM,mBAAmB,GAAG,IAAI,CAAC,GAAG,CAAC,CAAC,OAAO,GAAG,MAAM,KAAK,EAAE,GAAG,EAAE,CAAC,CAAC,CAAC;IACrE,CAAC,MAAM,EAAE,aAAa,CAAC,GAAG,cAAc,CACtC,eAAe,EACf,EAAE;AACF,IAAA,EAAE;AACF,IAAA,EAAE;AACF,IAAA,CAAC;AACD,IAAA,QAAQ;AACR,IAAA,mBAAmB;AACnB,IAAA,IAAI,CACL,CAAC;AAEF,IAAA,OAAO,CAAC,MAAM,EAAA,MAAA,CAAA,MAAA,CAAA,MAAA,CAAA,MAAA,CAAA,EAAA,EAAO,IAAI,CAAE,EAAA,EAAA,aAAa,IAAG,CAAC;AAC9C,CAAC,CAAC;AAEF,MAAM,cAAc,GAAG,CACrB,eAAgC,EAChC,CAAoB,EACpB,CAAS;AACT,CAAS;AACT,KAAa;AACb,QAAgB;AAChB,mBAA2B;AAC3B,IAAoB,KACM;;IAE1B,MAAM,CAAC,eAAe,EAAE,YAAY,CAAC,GAAG,eAAe,CACrD,CAAC;AACD,IAAA,CAAC;AACD,IAAA,CAAC,CACF,CAAC;IAEF,MAAM,aAAa,GAAG,IAAI,CAAC,GAAG,CAAC,YAAY,GAAG,eAAe,CAAC,CAAC;IAE/D,IAAI,KAAK,IAAI,QAAQ,EAAE;;QAErB,EAAE,IAAI,CAAC,KAAK,CAAC;AACb,QAAA,OAAO,CAAC,eAAe,EAAE,aAAa,CAAC,CAAC;AACzC,KAAA;IAED,MAAM,eAAe,GAAG,mBAAmB,IAAI,CAAC,GAAG,CAAC,CAAC,CAAC;IACtD,IAAI,aAAa,IAAI,eAAe,EAAE;;QAEpC,EAAE,IAAI,CAAC,KAAK,CAAC;AACb,QAAA,OAAO,CAAC,eAAe,EAAE,aAAa,CAAC,CAAC;AACzC,KAAA;IAED,MAAM,GAAG,GAAG,CAAC,CAAC,GAAG,CAAC,IAAI,CAAC,CAAC;AACxB,IAAA,IAAI,CAAC,IAAI,GAAG,IAAI,GAAG,IAAI,CAAC,EAAE;;AAExB,QAAA,IAAI,CAAC,YAAY,GAAG,IAAI,CAAC;AACzB,QAAA,MAAM,iBAAiB,GAAG,KAAK,CAAC,aAAa,CAAC,GAAG,CAAC,GAAG,aAAa,CAAC;AACnE,QAAA,IAAI,KAAK,CAAC,eAAe,CAAC,IAAI,IAAI,CAAC,GAAG,CAAC,eAAe,CAAC,KAAK,QAAQ,EAAE;AACpE,YAAA,OAAO,CAAC,CAAC,EAAE,iBAAiB,CAAC,CAAC;AAC/B,SAAA;AACD,QAAA,OAAO,CAAC,eAAe,EAAE,iBAAiB,CAAC,CAAC;AAC7C,KAAA;;AAGD,IAAA,EAAE,KAAK,CAAC;AACR,IAAA,IAAI,KAAK,GAAG,IAAI,CAAC,KAAK,EAAE;AACtB,QAAA,IAAI,CAAC,KAAK,GAAG,KAAK,CAAC;AACpB,KAAA;AACD,IAAA,MAAM,CAAC,YAAY,EAAE,iBAAiB,CAAC,GAAG,cAAc,CACtD,eAAe,EACf,CAAC;AACD,IAAA,CAAC;AACD,IAAA,GAAG;AACH,IAAA,KAAK;AACL,IAAA,QAAQ;AACR,IAAA,mBAAmB;AACnB,IAAA,IAAI,CACL,CAAC;AACF,IAAA,MAAM,CAAC,aAAa,EAAE,kBAAkB,CAAC,GAAG,cAAc,CACxD,eAAe,EACf,CAAC;AACD,IAAA,GAAG;AACH,IAAA,CAAC;AACD,IAAA,KAAK;AACL,IAAA,QAAQ;AACR,IAAA,mBAAmB;AACnB,IAAA,IAAI,CACL,CAAC;IACF,OAAO,CAAC,YAAY,GAAG,aAAa,EAAE,iBAAiB,GAAG,kBAAkB,CAAC,CAAC;AAChF,CAAC;;ACxLD;AAIA;AACA;AAEA;AACA;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AACrC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AACrC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AAErC;AACA,MAAM,OAAO,GAAG,mBAAmB,CAAC;AACpC;AACA,MAAM,OAAO,GAAG,kBAAkB,CAAC;AACnC;AACA,MAAM,OAAO,GAAG,kBAAkB,CAAC;AACnC;AACA,MAAM,OAAO,GAAG,kBAAkB,CAAC;AAEnC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AACrC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AACrC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AAErC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,kBAAkB,CAAC;AACrC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,mBAAmB,CAAC;AACtC;AACA,MAAM,SAAS,GAAG,oBAAoB,CAAC;AAEvC;;;;;;;AAOG;AACI,MAAM,eAAe,GAAoB,CAC9C,CAAoB;AACpB,CAAS;AACT,CAAS,KACiB;;IAE1B,MAAM,KAAK,GAAG,CAAC,CAAC,GAAG,CAAC,IAAI,CAAC,CAAC;IAC1B,MAAM,MAAM,GAAG,CAAC,CAAC,GAAG,CAAC,IAAI,CAAC,CAAC;;AAG3B,IAAA,MAAM,SAAS,GAAG,OAAO,GAAG,KAAK,CAAC;AAClC,IAAA,MAAM,WAAW,GAAG,SAAS,GAAG,KAAK,CAAC;AACtC,IAAA,MAAM,SAAS,GAAG,OAAO,GAAG,KAAK,CAAC;AAClC,IAAA,MAAM,WAAW,GAAG,SAAS,GAAG,KAAK,CAAC;AACtC,IAAA,MAAM,SAAS,GAAG,OAAO,GAAG,KAAK,CAAC;AAClC,IAAA,MAAM,WAAW,GAAG,SAAS,GAAG,KAAK,CAAC;AACtC,IAAA,MAAM,SAAS,GAAG,OAAO,GAAG,KAAK,CAAC;AAElC,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAE3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAC3C,IAAA,MAAM,gBAAgB,GAAG,SAAS,GAAG,KAAK,CAAC;AAE3C,IAAA,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,CAAC,CAAC;IAEzB,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IACzC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IAEzC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IACzC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IAEzC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IACzC,MAAM,QAAQ,GAAG,CAAC,CAAC,MAAM,GAAG,WAAW,CAAC,CAAC;IAEzC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;IACrC,MAAM,MAAM,GAAG,CAAC,CAAC,MAAM,GAAG,SAAS,CAAC,CAAC;AAErC,IAAA,MAAM,YAAY,GAChB,gBAAgB,GAAG,MAAM;AACzB,QAAA,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC;AACxC,QAAA,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC;QACxC,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC,CAAC;AAE3C,IAAA,MAAM,eAAe,GACnB,gBAAgB,GAAG,MAAM;AACzB,QAAA,gBAAgB,IAAI,MAAM,GAAG,MAAM,CAAC;AACpC,QAAA,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC;AACxC,QAAA,gBAAgB,IAAI,MAAM,GAAG,MAAM,CAAC;AACpC,QAAA,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC;AACxC,QAAA,gBAAgB,IAAI,MAAM,GAAG,MAAM,CAAC;AACpC,QAAA,gBAAgB,IAAI,QAAQ,GAAG,QAAQ,CAAC;QACxC,gBAAgB,IAAI,MAAM,GAAG,MAAM,CAAC,CAAC;AAEvC,IAAA,OAAO,CAAC,eAAe,EAAE,YAAY,CAAC,CAAC;AACzC,CAAC;;AC/HD;AAwCA,MAAM,iBAAiB,GACrB,8DAA8D,CAAC;AAEjE;;;;;;;;AAQG;AACI,MAAM,IAAI,GAAG,CAClB,CAAoB,EACpB,KAAe,EACf,OAAA,GAA6B,EAAE,KACG;;AAElC,IAAA,IAAI,CAAC,KAAK,CAAC,OAAO,CAAC,KAAK,CAAC,EAAE;AACzB,QAAA,MAAM,IAAI,UAAU,CAAC,iBAAiB,CAAC,CAAC;AACzC,KAAA;AAED,IAAA,IAAI,KAAK,CAAC,MAAM,KAAK,CAAC,EAAE;;AAEtB,QAAA,MAAM,CAAC,CAAC,EAAE,CAAC,CAAC,GAAG,KAAK,CAAC;AACrB,QAAA,OAAO,SAAS,CAAC,eAAe,EAAE,CAAC,EAAE,CAAC,EAAE,CAAC,EAAE,OAAO,CAAC,CAAC;AACrD,KAAA;AAED,IAAA,IAAI,KAAK,CAAC,MAAM,GAAG,CAAC,EAAE;;AAEpB,QAAA,MAAM,IAAI,UAAU,CAAC,iBAAiB,CAAC,CAAC;AACzC,KAAA;;IAGD,IAAI,MAAM,GAAG,CAAC,CAAC;IACf,MAAM,MAAM,GAAqC,EAAE,CAAC;AACpD,IAAA,MAAM,IAAI,GAAG;AACX,QAAA,KAAK,EAAE,CAAC;AACR,QAAA,aAAa,EAAE,CAAC;QAChB,MAAM;AACN,QAAA,KAAK,EAAE,CAAC;KACT,CAAC;AAEF,IAAA,KAAK,IAAI,CAAC,GAAG,CAAC,EAAE,CAAC,GAAG,KAAK,CAAC,MAAM,GAAG,CAAC,EAAE,EAAE,CAAC,EAAE;QACzC,MAAM,MAAM,GAAG,SAAS,CACtB,eAAe,EACf,CAAC,EACD,KAAK,CAAC,CAAC,CAAC,EACR,KAAK,CAAC,CAAC,GAAG,CAAC,CAAC,EACZ,OAAO,CACR,CAAC;AAEF,QAAA,IAAI,CAAC,MAAM,CAAC,IAAI,CAAC,MAAM,CAAC,CAAC;;AAEzB,QAAA,MAAM,IAAI,MAAM,CAAC,CAAC,CAAC,CAAC;;QAEpB,IAAI,CAAC,KAAK,IAAI,MAAM,CAAC,CAAC,CAAC,CAAC,KAAK,CAAC;QAC9B,IAAI,CAAC,aAAa,IAAI,MAAM,CAAC,CAAC,CAAC,CAAC,aAAa,CAAC;AAC9C,QAAA,IAAI,CAAC,KAAK,GAAG,IAAI,CAAC,GAAG,CAAC,IAAI,CAAC,KAAK,EAAE,MAAM,CAAC,CAAC,CAAC,CAAC,KAAK,CAAC,CAAC;AACpD,KAAA;AAED,IAAA,OAAO,CAAC,MAAM,EAAE,IAAI,CAAC,CAAC;AACxB,CAAC;;ACtGD;;;;;;;;;"}
export { quad } from './quad/index.js';
//# sourceMappingURL=index.d.ts.map
{"version":3,"file":"index.d.ts","sourceRoot":"","sources":["../../src/numerical/index.ts"],"names":[],"mappings":"AAEA,OAAO,EAAE,IAAI,EAAE,MAAM,iBAAiB,CAAC"}
import type { IntegrandCallback, IntegrationStep, QuadratureInfo, QuadratureOptions } from '.';
export declare const integrate: (integrationStep: IntegrationStep, f: IntegrandCallback, a: number, b: number, options?: QuadratureOptions) => [r: number, i: QuadratureInfo];
//# sourceMappingURL=adaptive-quadrature.d.ts.map
{"version":3,"file":"adaptive-quadrature.d.ts","sourceRoot":"","sources":["../../../src/numerical/quad/adaptive-quadrature.ts"],"names":[],"mappings":"AAEA,OAAO,KAAK,EACV,iBAAiB,EACjB,eAAe,EACf,cAAc,EACd,iBAAiB,EAClB,MAAM,GAAG,CAAC;AA+DX,eAAO,MAAM,SAAS,oBACH,eAAe,KAC7B,iBAAiB,KACjB,MAAM,KACN,MAAM,YACA,iBAAiB,KACzB,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CAqC/B,CAAC"}
import type { IntegrationStep } from '.';
/**
* Perform a single integration step.
*
* @param f Integrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns [current best estimate, poor estimate]
*/
export declare const integrationStep: IntegrationStep;
//# sourceMappingURL=gauss-kronrod-g7k15.d.ts.map
{"version":3,"file":"gauss-kronrod-g7k15.d.ts","sourceRoot":"","sources":["../../../src/numerical/quad/gauss-kronrod-g7k15.ts"],"names":[],"mappings":"AAEA,OAAO,KAAK,EAAqB,eAAe,EAAE,MAAM,GAAG,CAAC;AAgD5D;;;;;;;GAOG;AACH,eAAO,MAAM,eAAe,EAAE,eAqE7B,CAAC"}
export declare type IntegrandCallback = (a: number) => number;
export interface QuadratureInfo {
/** Number of steps _used_ in the integration. */
steps: number;
/** Estimate of the global error. */
errorEstimate: number;
/**
* Set to true if there was a problem (e.g. could not reach required accuracy
* with a step at machine precision).
*/
isUnreliable?: true;
/**
* If the range has multiple (i.e. more than 2) points, contains the results
* for intermediate ranges.
*/
points?: [r: number, i: QuadratureInfo][];
/** The maximum depth used. */
depth: number;
}
export interface QuadratureOptions {
/** Used to control global error. */
epsilon?: number;
/** Maximum depth to use for adaptive step length. */
maxDepth?: number;
}
export declare type IntegrationStep = (
/** The integrand. */
f: (a: number) => number, a: number, // Lower limit.
b: number) => [a: number, b: number];
/**
* Numerically compute a definite integral.
*
* @param f Callback returning value of integrand.
* @param range A range of at least 2 endpoints.
* @param options Options for the computation.
*
* @returns The results of the computation.
*/
export declare const quad: (f: IntegrandCallback, range: number[], options?: QuadratureOptions) => [r: number, i: QuadratureInfo];
//# sourceMappingURL=index.d.ts.map
{"version":3,"file":"index.d.ts","sourceRoot":"","sources":["../../../src/numerical/quad/index.ts"],"names":[],"mappings":"AAKA,oBAAY,iBAAiB,GAAG,CAAC,CAAC,EAAE,MAAM,KAAK,MAAM,CAAC;AAEtD,MAAM,WAAW,cAAc;IAC7B,iDAAiD;IACjD,KAAK,EAAE,MAAM,CAAC;IACd,oCAAoC;IACpC,aAAa,EAAE,MAAM,CAAC;IACtB;;;OAGG;IACH,YAAY,CAAC,EAAE,IAAI,CAAC;IACpB;;;OAGG;IACH,MAAM,CAAC,EAAE,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CAAC,EAAE,CAAC;IAC1C,8BAA8B;IAC9B,KAAK,EAAE,MAAM,CAAC;CACf;AAED,MAAM,WAAW,iBAAiB;IAChC,oCAAoC;IACpC,OAAO,CAAC,EAAE,MAAM,CAAC;IACjB,qDAAqD;IACrD,QAAQ,CAAC,EAAE,MAAM,CAAC;CACnB;AAED,oBAAY,eAAe,GAAG;AAC5B,qBAAqB;AACrB,CAAC,EAAE,CAAC,CAAC,EAAE,MAAM,KAAK,MAAM,EACxB,CAAC,EAAE,MAAM,EAAE,eAAe;AAC1B,CAAC,EAAE,MAAM,KACN,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,MAAM,CAAC,CAAC;AAK5B;;;;;;;;GAQG;AACH,eAAO,MAAM,IAAI,MACZ,iBAAiB,SACb,MAAM,EAAE,YACN,iBAAiB,KACzB,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CA8C/B,CAAC"}
export declare const version = "1.3.0";
//# sourceMappingURL=version.d.ts.map
{"version":3,"file":"version.d.ts","sourceRoot":"","sources":["../src/version.ts"],"names":[],"mappings":"AAEA,eAAO,MAAM,OAAO,UAAU,CAAC"}
+2
-2

@@ -1,6 +0,6 @@

/*! @rec-math/math v1.2.0 2022-06-14 11:50:09
/*! @rec-math/math v1.3.0 2022-06-28 23:54:20
* https://github.com/rec-math/ts-math#readme
* Copyright pbuk (https://github.com/pb-uk) MIT license.
*/
var RecMath=function(t){"use strict";const e={epsilon:16*Number.EPSILON,maxDepth:Number.POSITIVE_INFINITY},r=(t,r,n,i,o={})=>{const a=Object.assign(Object.assign({},e),o),{epsilon:u,maxDepth:c}=a,p={steps:0,errorEstimate:0,depth:1},[h,l,b]=((t,e,r)=>{if(!isFinite(r))return isFinite(e)?[r=>{const s=1/(1-r);return t(e+r*s)*s*s},0,1]:[e=>{const r=e*e,s=1/(1-r);return t(e*s)*(1+r)*s*s},-1,1];if(!isFinite(e))return[e=>t(r-(1-e)/e)/(e*e),0,1];return[t,e,r]})(r,n,i);let[f,g]=s(t,h,l,b,1,1,0,Object.assign({},p));const m=Math.abs(u*f/(b-l));return[f,g]=s(t,h,l,b,1,c,m,p),[f,Object.assign(Object.assign({},p),{errorEstimate:g})]},s=(t,e,r,n,i,o,a,u)=>{const c=t(e,r,n);let p=c[0];const h=c[1];let l=Math.abs(h-p);if(i>=o)return++u.steps,[p,l];if(l<=a*(n-r))return++u.steps,[p,l];const b=(r+n)/2;if(r>=b||b>=n){u.isUnreliable=!0;const t=isNaN(l)?0:l;return p===Number.POSITIVE_INFINITY?[0,t]:[p,t]}++i>u.depth&&(u.depth=i);const f=s(t,e,r,b,i,o,a,u),g=s(t,e,b,n,i,o,a,u);return p=f[0]+g[0],l=f[1]+g[1],[p,l]},n=(t,e,r)=>{const s=(r-e)/2,n=(e+r)/2,i=.20778495500789848*s,o=.4058451513773972*s,a=.5860872354676911*s,u=.7415311855993945*s,c=.8648644233597691*s,p=.9491079123427585*s,h=.9914553711208126*s,l=.4179591836734694*s,b=.3818300505051189*s,f=.27970539148927664*s,g=.1294849661688697*s,m=.20948214108472782*s,d=.20443294007529889*s,I=.19035057806478542*s,E=.1690047266392679*s,N=.14065325971552592*s,O=.10479001032225019*s,_=.06309209262997856*s,j=.022935322010529224*s,w=t(n),F=t(n+i),M=t(n-i),v=t(n+o),y=t(n-o),P=t(n+a),T=t(n-a),x=t(n+u),R=t(n-u),S=t(n+c),A=t(n-c),D=t(n+p),V=t(n-p);return[m*w+d*(F+M)+I*(v+y)+E*(P+T)+N*(x+R)+O*(S+A)+_*(D+V)+j*(t(n+h)+t(n-h)),l*w+b*(v+y)+f*(x+R)+g*(D+V)]},i="integration range must be an array of at least two endpoints";var o=Object.freeze({__proto__:null,quad:(t,e,s={})=>{if(!Array.isArray(e))throw new RangeError(i);if(2===e.length){const[i,o]=e;return r(n,t,i,o,s)}if(e.length<2)throw new RangeError(i);let o=0;const a={steps:0,errorEstimate:0,points:[],depth:0};for(let i=0;i<e.length-1;++i){const u=r(n,t,e[i],e[i+1],s);a.points.push(u),o+=u[0],a.steps+=u[1].steps,a.errorEstimate+=u[1].errorEstimate,a.depth=Math.max(a.depth,u[1].depth)}return[o,a]}});return t.integrate=o,t.version="1.2.0",Object.defineProperty(t,"__esModule",{value:!0}),t}({});
!function(t){"use strict";const e={epsilon:16*Number.EPSILON,maxDepth:1/0},r=(t,r,n,i,a={})=>{const o=Object.assign(Object.assign({},e),a),{epsilon:c,maxDepth:h}=o,p={steps:0,errorEstimate:0,depth:1},[u,l,b]=((t,e,r)=>{if(!isFinite(r))return isFinite(e)?[r=>{const s=1/(1-r);return t(e+r*s)*s*s},0,1]:[e=>{const r=e*e,s=1/(1-r);return t(e*s)*(1+r)*s*s},-1,1];if(!isFinite(e))return[e=>t(r-(1-e)/e)/(e*e),0,1];return[t,e,r]})(r,n,i);let[f,g]=s(t,u,l,b,1,1,0,Object.assign({},p));const d=Math.abs(c*f/(b-l));return[f,g]=s(t,u,l,b,1,h,d,p),[f,Object.assign(Object.assign({},p),{errorEstimate:g})]},s=(t,e,r,n,i,a,o,c)=>{const[h,p]=t(e,r,n),u=Math.abs(p-h);if(i>=a)return++c.steps,[h,u];if(u<=o*(n-r))return++c.steps,[h,u];const l=(r+n)/2;if(r>=l||l>=n){c.isUnreliable=!0;const t=isNaN(u)?0:u;return isNaN(h)||Math.abs(h)===1/0?[0,t]:[h,t]}++i>c.depth&&(c.depth=i);const[b,f]=s(t,e,r,l,i,a,o,c),[g,d]=s(t,e,l,n,i,a,o,c);return[b+g,f+d]},n=(t,e,r)=>{const s=(r-e)/2,n=(e+r)/2,i=.20778495500789848*s,a=.4058451513773972*s,o=.5860872354676911*s,c=.7415311855993945*s,h=.8648644233597691*s,p=.9491079123427585*s,u=.9914553711208126*s,l=.4179591836734694*s,b=.3818300505051189*s,f=.27970539148927664*s,g=.1294849661688697*s,d=.20948214108472782*s,m=.20443294007529889*s,E=.19035057806478542*s,O=.1690047266392679*s,j=.14065325971552592*s,M=.10479001032225019*s,N=.06309209262997856*s,_=.022935322010529224*s,w=t(n),y=t(n+i),R=t(n-i),v=t(n+a),x=t(n-a),F=t(n+o),A=t(n-o),D=t(n+c),P=t(n-c),q=t(n+h),z=t(n-h),I=t(n+p),L=t(n-p);return[d*w+m*(y+R)+E*(v+x)+O*(F+A)+j*(D+P)+M*(q+z)+N*(I+L)+_*(t(n+u)+t(n-u)),l*w+b*(v+x)+f*(D+P)+g*(I+L)]},i="integration range must be an array of at least two endpoints";var a=Object.freeze({__proto__:null,quad:(t,e,s={})=>{if(!Array.isArray(e))throw new RangeError(i);if(2===e.length){const[i,a]=e;return r(n,t,i,a,s)}if(e.length<2)throw new RangeError(i);let a=0;const o={steps:0,errorEstimate:0,points:[],depth:0};for(let i=0;i<e.length-1;++i){const c=r(n,t,e[i],e[i+1],s);o.points.push(c),a+=c[0],o.steps+=c[1].steps,o.errorEstimate+=c[1].errorEstimate,o.depth=Math.max(o.depth,c[1].depth)}return[a,o]}});t.numerical=a,t.version="1.3.0",Object.defineProperty(t,"__esModule",{value:!0})}(this.RecMath=this.RecMath||{});
//# sourceMappingURL=rec-math.min.js.map

@@ -1,1 +0,1 @@

{"version":3,"file":"rec-math.min.js","sources":["../esm/integrate/quad/adaptive-quadrature.js","../esm/integrate/quad/gauss-kronrod-g7k15.js","../esm/integrate/quad.js","../esm/index.js"],"sourcesContent":["// Export the API.\r\nexport { quadrature };\r\nconst defaults = {\r\n epsilon: Number.EPSILON * 16,\r\n maxDepth: Number.POSITIVE_INFINITY,\r\n};\r\n/**\r\n * Perform a substitution with an appropriate change of variables to deal with\r\n * infinite ranges.\r\n *\r\n * @param f Intgrand callback.\r\n * @param a Lower limit.\r\n * @param b Upper limit.\r\n * @returns An array with any necessary substitution.\r\n */\r\nconst changeOfVariables = (f, a, b) => {\r\n if (!isFinite(b)) {\r\n if (!isFinite(a)) {\r\n // Change variables to integrate between - and + infinity.\r\n const _f = (t) => {\r\n const tSquared = t * t;\r\n const oneOverOneMinusTSquared = 1 / (1 - tSquared);\r\n return (f(t * oneOverOneMinusTSquared) *\r\n (1 + tSquared) *\r\n oneOverOneMinusTSquared *\r\n oneOverOneMinusTSquared);\r\n };\r\n return [_f, -1, 1];\r\n }\r\n // Change variables to integrate up to infinity.\r\n const _f = (t) => {\r\n const oneOverOneMinusT = 1 / (1 - t);\r\n return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;\r\n };\r\n return [_f, 0, 1];\r\n }\r\n if (!isFinite(a)) {\r\n // Change variables to integrate up from negative infinity.\r\n const _f = (t) => {\r\n return f(b - (1 - t) / t) / (t * t);\r\n };\r\n return [_f, 0, 1];\r\n }\r\n return [f, a, b];\r\n};\r\nconst quadrature = (integrationStep, f, a, b, options = {}) => {\r\n // Establish settings.\r\n const settings = Object.assign(Object.assign({}, defaults), options);\r\n const { epsilon, maxDepth } = settings;\r\n const info = { steps: 0, errorEstimate: 0, depth: 1 };\r\n // Allow for a change of variables to deal with infinite limits.\r\n const [_f, _a, _b] = changeOfVariables(f, a, b);\r\n // Get estimate so we can work out the acceptable global error.\r\n // Use a depth of 1 to calculate a 15 point Kronrod quadrature.\r\n let [result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.\r\n _a, // Lower limit.\r\n _b, // Upper limit.\r\n 1, // New depth.\r\n 1, // Maximum depth.\r\n 0, Object.assign({}, info));\r\n // Now calculate using the target global error.\r\n const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));\r\n [result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.\r\n _a, // Lower limit.\r\n _b, // Upper limit.\r\n 1, // New depth.\r\n maxDepth, // Maximum depth.\r\n acceptableUnitError, // Acceptable error per unit step.\r\n info);\r\n return [result, Object.assign(Object.assign({}, info), { errorEstimate })];\r\n};\r\nconst integrate_part = (integrationStep, f, a, // Lower limit.\r\nb, // Upper limit.\r\ndepth, // New depth.\r\nmaxDepth, // Maximum depth.\r\nacceptableUnitError, // Acceptable error per unit step.\r\ninfo) => {\r\n // Initialize things.\r\n const estimates = integrationStep(f, // Integrand.\r\n a, // Lower limit.\r\n b);\r\n let currentEstimate = estimates[0];\r\n const poorEstimate = estimates[1];\r\n let errorEstimate = Math.abs(poorEstimate - currentEstimate);\r\n if (depth >= maxDepth) {\r\n // Reached the maximum allowable depth so return the partial sum.\r\n ++info.steps;\r\n return [currentEstimate, errorEstimate];\r\n }\r\n const acceptableError = acceptableUnitError * (b - a);\r\n if (errorEstimate <= acceptableError) {\r\n // Error is acceptable for the size of step so return the partial sum.\r\n ++info.steps;\r\n return [currentEstimate, errorEstimate];\r\n }\r\n const mid = (a + b) / 2;\r\n if (a >= mid || mid >= b) {\r\n // We can't make this step any smaller: looks like a discontinuity.\r\n info.isUnreliable = true;\r\n const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;\r\n if (currentEstimate === Number.POSITIVE_INFINITY) {\r\n return [0, safeErrorEstimate];\r\n }\r\n return [currentEstimate, safeErrorEstimate];\r\n }\r\n // Recurse deeper.\r\n ++depth;\r\n if (depth > info.depth) {\r\n info.depth = depth;\r\n }\r\n const leftResult = integrate_part(integrationStep, f, // Integrand.\r\n a, // Lower limit.\r\n mid, // Upper limit.\r\n depth, // New depth.\r\n maxDepth, // Maximum depth.\r\n acceptableUnitError, // Acceptable error per unit step.\r\n info);\r\n const rightResult = integrate_part(integrationStep, f, // Integrand.\r\n mid, // Lower limit.\r\n b, // Upper limit.\r\n depth, // New depth.\r\n maxDepth, // Maximum depth.\r\n acceptableUnitError, // Acceptable error per unit step.\r\n info);\r\n currentEstimate = leftResult[0] + rightResult[0];\r\n errorEstimate = leftResult[1] + rightResult[1];\r\n return [currentEstimate, errorEstimate];\r\n};\r\n","// Export the API.\r\nexport { integrationStep };\r\n// Gauss-Kronrod constants (G7, K15) on [-1, 1].\r\n// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/\r\n// node_g0k0 = 0;\r\n// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;\r\nconst node_g1k2 = 0.4058451513773972;\r\n// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;\r\nconst node_g2k4 = 0.7415311855993945;\r\n// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;\r\nconst node_g3k6 = 0.9491079123427585;\r\n// node_k1 = 2.077849550078984676006894037732449e-01 exact;\r\nconst node_k1 = 0.20778495500789848;\r\n// node_k3 = 5.860872354676911302941448382587296e-01 exact;\r\nconst node_k3 = 0.5860872354676911;\r\n// node_k5 = 8.648644233597690727897127886409262e-01 exact;\r\nconst node_k5 = 0.8648644233597691;\r\n// node_k7 = 9.914553711208126392068546975263285e-01 exact;\r\nconst node_k7 = 0.9914553711208126;\r\n// weight_g0 = 4.179591836734693877551020408163265e-01 exact;\r\nconst weight_g0 = 0.4179591836734694;\r\n// weight_g1 = 3.818300505051189449503697754889751e-01 exact;\r\nconst weight_g1 = 0.3818300505051189;\r\n// weight_g2 = 2.797053914892766679014677714237796e-01 exact;\r\nconst weight_g2 = 0.27970539148927664;\r\n// weight_g3 = 1.294849661688696932706114326790820e-01 exact;\r\nconst weight_g3 = 0.1294849661688697;\r\n// weight_k0 = 2.094821410847278280129991748917143e-01 exact;\r\nconst weight_k0 = 0.20948214108472782;\r\n// weight_k1 = 2.044329400752988924141619992346491e-01 exact;\r\nconst weight_k1 = 0.20443294007529889;\r\n// weight_k2 = 1.903505780647854099132564024210137e-01 exact;\r\nconst weight_k2 = 0.19035057806478542;\r\n// weight_k3 = 1.690047266392679028265834265985503e-01 exact;\r\nconst weight_k3 = 0.1690047266392679;\r\n// weight_k4 = 1.406532597155259187451895905102379e-01 exact;\r\nconst weight_k4 = 0.14065325971552592;\r\n// weight_k5 = 1.047900103222501838398763225415180e-01 exact;\r\nconst weight_k5 = 0.10479001032225019;\r\n// weight_k6 = 6.309209262997855329070066318920429e-02 exact;\r\nconst weight_k6 = 0.06309209262997856;\r\n// weight_k7 = 2.293532201052922496373200805896959e-02 exact;\r\nconst weight_k7 = 0.022935322010529224;\r\n/**\r\n * Perform a single integration step.\r\n *\r\n * @param f Integrand callback.\r\n * @param a Lower limit.\r\n * @param b Upper limit.\r\n * @returns [current best estimate, poor estimate]\r\n */\r\nconst integrationStep = (f, // Integrand.\r\na, // Lower limit.\r\nb) => {\r\n // Gauss-Kronrod (G7, K15).\r\n const scale = (b - a) / 2;\r\n const offset = (a + b) / 2;\r\n // const scaled_g0k0 = 0;\r\n const scaled_k1 = node_k1 * scale;\r\n const scaled_g1k2 = node_g1k2 * scale;\r\n const scaled_k3 = node_k3 * scale;\r\n const scaled_g2k4 = node_g2k4 * scale;\r\n const scaled_k5 = node_k5 * scale;\r\n const scaled_g3k6 = node_g3k6 * scale;\r\n const scaled_k7 = node_k7 * scale;\r\n const scaled_weight_g0 = weight_g0 * scale;\r\n const scaled_weight_g1 = weight_g1 * scale;\r\n const scaled_weight_g2 = weight_g2 * scale;\r\n const scaled_weight_g3 = weight_g3 * scale;\r\n const scaled_weight_k0 = weight_k0 * scale;\r\n const scaled_weight_k1 = weight_k1 * scale;\r\n const scaled_weight_k2 = weight_k2 * scale;\r\n const scaled_weight_k3 = weight_k3 * scale;\r\n const scaled_weight_k4 = weight_k4 * scale;\r\n const scaled_weight_k5 = weight_k5 * scale;\r\n const scaled_weight_k6 = weight_k6 * scale;\r\n const scaled_weight_k7 = weight_k7 * scale;\r\n const f_g0k0 = f(offset);\r\n const f_k1_h = f(offset + scaled_k1);\r\n const f_k1_l = f(offset - scaled_k1);\r\n const f_g1k2_h = f(offset + scaled_g1k2);\r\n const f_g1k2_l = f(offset - scaled_g1k2);\r\n const f_k3_h = f(offset + scaled_k3);\r\n const f_k3_l = f(offset - scaled_k3);\r\n const f_g2k4_h = f(offset + scaled_g2k4);\r\n const f_g2k4_l = f(offset - scaled_g2k4);\r\n const f_k5_h = f(offset + scaled_k5);\r\n const f_k5_l = f(offset - scaled_k5);\r\n const f_g3k6_h = f(offset + scaled_g3k6);\r\n const f_g3k6_l = f(offset - scaled_g3k6);\r\n const f_k7_h = f(offset + scaled_k7);\r\n const f_k7_l = f(offset - scaled_k7);\r\n const poorEstimate = scaled_weight_g0 * f_g0k0 + // g0\r\n scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1\r\n scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2\r\n scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3\r\n const currentEstimate = scaled_weight_k0 * f_g0k0 + // k0\r\n scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1\r\n scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2\r\n scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3\r\n scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4\r\n scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5\r\n scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6\r\n scaled_weight_k7 * (f_k7_h + f_k7_l); // k7\r\n return [currentEstimate, poorEstimate];\r\n};\r\n","import { quadrature } from './quad/adaptive-quadrature.js';\r\nimport { integrationStep } from './quad/gauss-kronrod-g7k15.js';\r\nconst rangeErrorMessage = 'integration range must be an array of at least two endpoints';\r\n/**\r\n * Numerically compute a definite integral.\r\n *\r\n * @param f Callback returning value of integrand.\r\n * @param range A range of at least 2 endpoints.\r\n * @param options Options for the computation.\r\n *\r\n * @returns The results of the computation.\r\n */\r\nexport const quad = (f, range, options = {}) => {\r\n // Interpret the range.\r\n if (!Array.isArray(range)) {\r\n throw new RangeError(rangeErrorMessage);\r\n }\r\n if (range.length === 2) {\r\n // Integrate over a single range.\r\n const [a, b] = range;\r\n return quadrature(integrationStep, f, a, b, options);\r\n }\r\n if (range.length < 2) {\r\n // Can't integrate at a point!\r\n throw new RangeError(rangeErrorMessage);\r\n }\r\n // Integrate over multiple ranges.\r\n let result = 0;\r\n const points = [];\r\n const info = {\r\n steps: 0,\r\n errorEstimate: 0,\r\n points,\r\n depth: 0,\r\n };\r\n for (let i = 0; i < range.length - 1; ++i) {\r\n const single = quadrature(integrationStep, f, range[i], range[i + 1], options);\r\n info.points.push(single);\r\n // Update the result.\r\n result += single[0];\r\n // Update the cumulative statistics.\r\n info.steps += single[1].steps;\r\n info.errorEstimate += single[1].errorEstimate;\r\n info.depth = Math.max(info.depth, single[1].depth);\r\n }\r\n return [result, info];\r\n};\r\n","export const version = '1.2.0';\r\nexport * as integrate from './integrate.js';\r\n"],"names":["defaults","epsilon","Number","EPSILON","maxDepth","POSITIVE_INFINITY","quadrature","integrationStep","f","a","b","options","settings","Object","assign","info","steps","errorEstimate","depth","_f","_a","_b","isFinite","t","oneOverOneMinusT","tSquared","oneOverOneMinusTSquared","changeOfVariables","result","integrate_part","acceptableUnitError","Math","abs","estimates","currentEstimate","poorEstimate","mid","isUnreliable","safeErrorEstimate","isNaN","leftResult","rightResult","scale","offset","scaled_k1","scaled_g1k2","scaled_k3","scaled_g2k4","scaled_k5","scaled_g3k6","scaled_k7","scaled_weight_g0","scaled_weight_g1","scaled_weight_g2","scaled_weight_g3","scaled_weight_k0","scaled_weight_k1","scaled_weight_k2","scaled_weight_k3","scaled_weight_k4","scaled_weight_k5","scaled_weight_k6","scaled_weight_k7","f_g0k0","f_k1_h","f_k1_l","f_g1k2_h","f_g1k2_l","f_k3_h","f_k3_l","f_g2k4_h","f_g2k4_l","f_k5_h","f_k5_l","f_g3k6_h","f_g3k6_l","rangeErrorMessage","range","Array","isArray","RangeError","length","points","i","single","push","max"],"mappings":";;;;qCAEA,MAAMA,EAAW,CACbC,QAA0B,GAAjBC,OAAOC,QAChBC,SAAUF,OAAOG,mBAyCfC,EAAa,CAACC,EAAiBC,EAAGC,EAAGC,EAAGC,EAAU,MAEpD,MAAMC,EAAWC,OAAOC,OAAOD,OAAOC,OAAO,GAAId,GAAWW,IACtDV,QAAEA,EAAOG,SAAEA,GAAaQ,EACxBG,EAAO,CAAEC,MAAO,EAAGC,cAAe,EAAGC,MAAO,IAE3CC,EAAIC,EAAIC,GApCO,EAACb,EAAGC,EAAGC,KAC7B,IAAKY,SAASZ,GACV,OAAKY,SAASb,GAiBP,CAJKc,IACR,MAAMC,EAAmB,GAAK,EAAID,GAClC,OAAOf,EAAEC,EAAIc,EAAIC,GAAoBA,EAAmBA,GAEhD,EAAG,GAPJ,CARKD,IACR,MAAME,EAAWF,EAAIA,EACfG,EAA0B,GAAK,EAAID,GACzC,OAAQjB,EAAEe,EAAIG,IACT,EAAID,GACLC,EACAA,IAEK,EAAG,GASxB,IAAKJ,SAASb,GAKV,MAAO,CAHKc,GACDf,EAAEE,GAAK,EAAIa,GAAKA,IAAMA,EAAIA,GAEzB,EAAG,GAEnB,MAAO,CAACf,EAAGC,EAAGC,IAQOiB,CAAkBnB,EAAGC,EAAGC,GAG7C,IAAKkB,EAAQX,GAAiBY,EAAetB,EAAiBY,EAC9DC,EACAC,EACA,EACA,EACA,EAAGR,OAAOC,OAAO,GAAIC,IAErB,MAAMe,EAAsBC,KAAKC,IAAK/B,EAAU2B,GAAWP,EAAKD,IAQhE,OAPCQ,EAAQX,GAAiBY,EAAetB,EAAiBY,EAC1DC,EACAC,EACA,EACAjB,EACA0B,EACAf,GACO,CAACa,EAAQf,OAAOC,OAAOD,OAAOC,OAAO,GAAIC,GAAO,CAAEE,oBAEvDY,EAAiB,CAACtB,EAAiBC,EAAGC,EAC5CC,EACAQ,EACAd,EACA0B,EACAf,KAEI,MAAMkB,EAAY1B,EAAgBC,EAClCC,EACAC,GACA,IAAIwB,EAAkBD,EAAU,GAChC,MAAME,EAAeF,EAAU,GAC/B,IAAIhB,EAAgBc,KAAKC,IAAIG,EAAeD,GAC5C,GAAIhB,GAASd,EAGT,QADEW,EAAKC,MACA,CAACkB,EAAiBjB,GAG7B,GAAIA,GADoBa,GAAuBpB,EAAID,GAI/C,QADEM,EAAKC,MACA,CAACkB,EAAiBjB,GAE7B,MAAMmB,GAAO3B,EAAIC,GAAK,EACtB,GAAID,GAAK2B,GAAOA,GAAO1B,EAAG,CAEtBK,EAAKsB,cAAe,EACpB,MAAMC,EAAoBC,MAAMtB,GAAiB,EAAIA,EACrD,OAAIiB,IAAoBhC,OAAOG,kBACpB,CAAC,EAAGiC,GAER,CAACJ,EAAiBI,KAG3BpB,EACUH,EAAKG,QACbH,EAAKG,MAAQA,GAEjB,MAAMsB,EAAaX,EAAetB,EAAiBC,EACnDC,EACA2B,EACAlB,EACAd,EACA0B,EACAf,GACM0B,EAAcZ,EAAetB,EAAiBC,EACpD4B,EACA1B,EACAQ,EACAd,EACA0B,EACAf,GAGA,OAFAmB,EAAkBM,EAAW,GAAKC,EAAY,GAC9CxB,EAAgBuB,EAAW,GAAKC,EAAY,GACrC,CAACP,EAAiBjB,IC3EvBV,EAAkB,CAACC,EACzBC,EACAC,KAEI,MAAMgC,GAAShC,EAAID,GAAK,EAClBkC,GAAUlC,EAAIC,GAAK,EAEnBkC,EA9CM,mBA8CgBF,EACtBG,EArDQ,kBAqDkBH,EAC1BI,EA9CM,kBA8CgBJ,EACtBK,EArDQ,kBAqDkBL,EAC1BM,EA9CM,kBA8CgBN,EACtBO,EArDQ,kBAqDkBP,EAC1BQ,EA9CM,kBA8CgBR,EACtBS,EA7CQ,kBA6CuBT,EAC/BU,EA5CQ,kBA4CuBV,EAC/BW,EA3CQ,mBA2CuBX,EAC/BY,EA1CQ,kBA0CuBZ,EAC/Ba,EAzCQ,mBAyCuBb,EAC/Bc,EAxCQ,mBAwCuBd,EAC/Be,EAvCQ,mBAuCuBf,EAC/BgB,EAtCQ,kBAsCuBhB,EAC/BiB,EArCQ,mBAqCuBjB,EAC/BkB,EApCQ,mBAoCuBlB,EAC/BmB,EAnCQ,mBAmCuBnB,EAC/BoB,EAlCQ,oBAkCuBpB,EAC/BqB,EAASvD,EAAEmC,GACXqB,EAASxD,EAAEmC,EAASC,GACpBqB,EAASzD,EAAEmC,EAASC,GACpBsB,EAAW1D,EAAEmC,EAASE,GACtBsB,EAAW3D,EAAEmC,EAASE,GACtBuB,EAAS5D,EAAEmC,EAASG,GACpBuB,EAAS7D,EAAEmC,EAASG,GACpBwB,EAAW9D,EAAEmC,EAASI,GACtBwB,EAAW/D,EAAEmC,EAASI,GACtByB,EAAShE,EAAEmC,EAASK,GACpByB,EAASjE,EAAEmC,EAASK,GACpB0B,EAAWlE,EAAEmC,EAASM,GACtB0B,EAAWnE,EAAEmC,EAASM,GAe5B,MAAO,CARiBM,EAAmBQ,EACvCP,GAAoBQ,EAASC,GAC7BR,GAAoBS,EAAWC,GAC/BT,GAAoBU,EAASC,GAC7BV,GAAoBW,EAAWC,GAC/BX,GAAoBY,EAASC,GAC7BZ,GAAoBa,EAAWC,GAC/Bb,GAbWtD,EAAEmC,EAASO,GACX1C,EAAEmC,EAASO,IACLC,EAAmBY,EACpCX,GAAoBc,EAAWC,GAC/Bd,GAAoBiB,EAAWC,GAC/BjB,GAAoBoB,EAAWC,KC7FjCC,EAAoB,wGAUN,CAACpE,EAAGqE,EAAOlE,EAAU,MAErC,IAAKmE,MAAMC,QAAQF,GACf,MAAM,IAAIG,WAAWJ,GAEzB,GAAqB,IAAjBC,EAAMI,OAAc,CAEpB,MAAOxE,EAAGC,GAAKmE,EACf,OAAOvE,EAAWC,EAAiBC,EAAGC,EAAGC,EAAGC,GAEhD,GAAIkE,EAAMI,OAAS,EAEf,MAAM,IAAID,WAAWJ,GAGzB,IAAIhD,EAAS,EACb,MACMb,EAAO,CACTC,MAAO,EACPC,cAAe,EACfiE,OAJW,GAKXhE,MAAO,GAEX,IAAK,IAAIiE,EAAI,EAAGA,EAAIN,EAAMI,OAAS,IAAKE,EAAG,CACvC,MAAMC,EAAS9E,EAAWC,EAAiBC,EAAGqE,EAAMM,GAAIN,EAAMM,EAAI,GAAIxE,GACtEI,EAAKmE,OAAOG,KAAKD,GAEjBxD,GAAUwD,EAAO,GAEjBrE,EAAKC,OAASoE,EAAO,GAAGpE,MACxBD,EAAKE,eAAiBmE,EAAO,GAAGnE,cAChCF,EAAKG,MAAQa,KAAKuD,IAAIvE,EAAKG,MAAOkE,EAAO,GAAGlE,OAEhD,MAAO,CAACU,EAAQb,qCC7CG"}
{"version":3,"file":"rec-math.min.js","sources":["../src/version.ts","../src/numerical/quad/adaptive-quadrature.ts","../src/numerical/quad/gauss-kronrod-g7k15.ts","../src/numerical/quad/index.ts"],"sourcesContent":["// rec-math/src/version.ts\n\nexport const version = '1.3.0';\n","// rec-math/src/numerical/quad/adaptive-quadrature.ts\n\nimport type {\n IntegrandCallback,\n IntegrationStep,\n QuadratureInfo,\n QuadratureOptions,\n} from '.';\n\n/** Defaults for integration. */\nconst defaults = {\n /**\n * Target error estimate for a step: if this is smaller it can lead to\n * accumulation of roundoff errors.\n */\n epsilon: Number.EPSILON * 16,\n /**\n * Maximum depth \\\\( d_{max} \\\\)for integration: the maximum number of steps\n * can be \\\\( 2^{d_{max}} \\\\).\n * */\n maxDepth: Infinity,\n};\n\n/**\n * Perform a substitution with an appropriate change of variables to deal with\n * infinite ranges.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns An array with any necessary substitution.\n */\nconst changeOfVariables = (\n f: IntegrandCallback,\n a: number,\n b: number,\n): [f: IntegrandCallback, a: number, b: number] => {\n if (!isFinite(b)) {\n if (!isFinite(a)) {\n // Change variables to integrate between - and + infinity.\n const _f = (t: number): number => {\n const tSquared = t * t;\n const oneOverOneMinusTSquared = 1 / (1 - tSquared);\n return (\n f(t * oneOverOneMinusTSquared) *\n (1 + tSquared) *\n oneOverOneMinusTSquared *\n oneOverOneMinusTSquared\n );\n };\n return [_f, -1, 1];\n }\n // Change variables to integrate up to infinity.\n const _f = (t: number): number => {\n const oneOverOneMinusT = 1 / (1 - t);\n return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;\n };\n return [_f, 0, 1];\n }\n\n if (!isFinite(a)) {\n // Change variables to integrate up from negative infinity.\n const _f = (t: number): number => {\n return f(b - (1 - t) / t) / (t * t);\n };\n return [_f, 0, 1];\n }\n return [f, a, b];\n};\n\nexport const integrate = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number,\n b: number,\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Establish settings.\n const settings = { ...defaults, ...options };\n const { epsilon, maxDepth } = settings;\n\n const info = { steps: 0, errorEstimate: 0, depth: 1 };\n\n // Allow for a change of variables to deal with infinite limits.\n const [_f, _a, _b] = changeOfVariables(f, a, b);\n\n // Get estimate so we can work out the acceptable global error.\n // Use a depth of 1 to calculate a 15 point Kronrod quadrature.\n let [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n 1, // Maximum depth.\n 0, // Acceptable error per unit step.\n { ...info }, // Statistics.\n );\n\n // Now calculate using the target global error.\n const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));\n [result, errorEstimate] = integrate_part(\n integrationStep,\n _f, // Integrand.\n _a, // Lower limit.\n _b, // Upper limit.\n 1, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n\n return [result, { ...info, errorEstimate }];\n};\n\nconst integrate_part = (\n integrationStep: IntegrationStep,\n f: IntegrandCallback,\n a: number, // Lower limit.\n b: number, // Upper limit.\n depth: number, // New depth.\n maxDepth: number, // Maximum depth.\n acceptableUnitError: number, // Acceptable error per unit step.\n info: QuadratureInfo, // Statistics.\n): [a: number, b: number] => {\n // Initialize things.\n const [currentEstimate, poorEstimate] = integrationStep(\n f, // Integrand.\n a, // Lower limit.\n b, // Upper limit.\n );\n\n const errorEstimate = Math.abs(poorEstimate - currentEstimate);\n\n if (depth >= maxDepth) {\n // Reached the maximum allowable depth so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const acceptableError = acceptableUnitError * (b - a);\n if (errorEstimate <= acceptableError) {\n // Error is acceptable for the size of step so return the partial sum.\n ++info.steps;\n return [currentEstimate, errorEstimate];\n }\n\n const mid = (a + b) / 2;\n if (a >= mid || mid >= b) {\n // We can't make this step any smaller: looks like a discontinuity.\n info.isUnreliable = true;\n const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;\n if (isNaN(currentEstimate) || Math.abs(currentEstimate) === Infinity) {\n return [0, safeErrorEstimate];\n }\n return [currentEstimate, safeErrorEstimate];\n }\n\n // Recurse deeper.\n ++depth;\n if (depth > info.depth) {\n info.depth = depth;\n }\n const [leftEstimate, leftErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n a, // Lower limit.\n mid, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n const [rightEstimate, rightErrorEstimate] = integrate_part(\n integrationStep,\n f, // Integrand.\n mid, // Lower limit.\n b, // Upper limit.\n depth, // New depth.\n maxDepth, // Maximum depth.\n acceptableUnitError, // Acceptable error per unit step.\n info, // Statistics.\n );\n return [leftEstimate + rightEstimate, leftErrorEstimate + rightErrorEstimate];\n};\n","// rec-math/src/numerical/quad/gauss-kronrod-g7k15.ts\n\nimport type { IntegrandCallback, IntegrationStep } from '.';\n\n// Gauss-Kronrod constants (G7, K15) on [-1, 1].\n// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/\n\n// node_g0k0 = 0;\n// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;\nconst node_g1k2 = 0.4058451513773972;\n// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;\nconst node_g2k4 = 0.7415311855993945;\n// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;\nconst node_g3k6 = 0.9491079123427585;\n\n// node_k1 = 2.077849550078984676006894037732449e-01 exact;\nconst node_k1 = 0.20778495500789848;\n// node_k3 = 5.860872354676911302941448382587296e-01 exact;\nconst node_k3 = 0.5860872354676911;\n// node_k5 = 8.648644233597690727897127886409262e-01 exact;\nconst node_k5 = 0.8648644233597691;\n// node_k7 = 9.914553711208126392068546975263285e-01 exact;\nconst node_k7 = 0.9914553711208126;\n\n// weight_g0 = 4.179591836734693877551020408163265e-01 exact;\nconst weight_g0 = 0.4179591836734694;\n// weight_g1 = 3.818300505051189449503697754889751e-01 exact;\nconst weight_g1 = 0.3818300505051189;\n// weight_g2 = 2.797053914892766679014677714237796e-01 exact;\nconst weight_g2 = 0.27970539148927664;\n// weight_g3 = 1.294849661688696932706114326790820e-01 exact;\nconst weight_g3 = 0.1294849661688697;\n\n// weight_k0 = 2.094821410847278280129991748917143e-01 exact;\nconst weight_k0 = 0.20948214108472782;\n// weight_k1 = 2.044329400752988924141619992346491e-01 exact;\nconst weight_k1 = 0.20443294007529889;\n// weight_k2 = 1.903505780647854099132564024210137e-01 exact;\nconst weight_k2 = 0.19035057806478542;\n// weight_k3 = 1.690047266392679028265834265985503e-01 exact;\nconst weight_k3 = 0.1690047266392679;\n// weight_k4 = 1.406532597155259187451895905102379e-01 exact;\nconst weight_k4 = 0.14065325971552592;\n// weight_k5 = 1.047900103222501838398763225415180e-01 exact;\nconst weight_k5 = 0.10479001032225019;\n// weight_k6 = 6.309209262997855329070066318920429e-02 exact;\nconst weight_k6 = 0.06309209262997856;\n// weight_k7 = 2.293532201052922496373200805896959e-02 exact;\nconst weight_k7 = 0.022935322010529224;\n\n/**\n * Perform a single integration step.\n *\n * @param f Integrand callback.\n * @param a Lower limit.\n * @param b Upper limit.\n * @returns [current best estimate, poor estimate]\n */\nexport const integrationStep: IntegrationStep = (\n f: IntegrandCallback, // Integrand.\n a: number, // Lower limit.\n b: number, // Upper limit.\n): [a: number, b: number] => {\n // Gauss-Kronrod (G7, K15).\n const scale = (b - a) / 2;\n const offset = (a + b) / 2;\n\n // const scaled_g0k0 = 0;\n const scaled_k1 = node_k1 * scale;\n const scaled_g1k2 = node_g1k2 * scale;\n const scaled_k3 = node_k3 * scale;\n const scaled_g2k4 = node_g2k4 * scale;\n const scaled_k5 = node_k5 * scale;\n const scaled_g3k6 = node_g3k6 * scale;\n const scaled_k7 = node_k7 * scale;\n\n const scaled_weight_g0 = weight_g0 * scale;\n const scaled_weight_g1 = weight_g1 * scale;\n const scaled_weight_g2 = weight_g2 * scale;\n const scaled_weight_g3 = weight_g3 * scale;\n\n const scaled_weight_k0 = weight_k0 * scale;\n const scaled_weight_k1 = weight_k1 * scale;\n const scaled_weight_k2 = weight_k2 * scale;\n const scaled_weight_k3 = weight_k3 * scale;\n const scaled_weight_k4 = weight_k4 * scale;\n const scaled_weight_k5 = weight_k5 * scale;\n const scaled_weight_k6 = weight_k6 * scale;\n const scaled_weight_k7 = weight_k7 * scale;\n\n const f_g0k0 = f(offset);\n\n const f_k1_h = f(offset + scaled_k1);\n const f_k1_l = f(offset - scaled_k1);\n const f_g1k2_h = f(offset + scaled_g1k2);\n const f_g1k2_l = f(offset - scaled_g1k2);\n\n const f_k3_h = f(offset + scaled_k3);\n const f_k3_l = f(offset - scaled_k3);\n const f_g2k4_h = f(offset + scaled_g2k4);\n const f_g2k4_l = f(offset - scaled_g2k4);\n\n const f_k5_h = f(offset + scaled_k5);\n const f_k5_l = f(offset - scaled_k5);\n const f_g3k6_h = f(offset + scaled_g3k6);\n const f_g3k6_l = f(offset - scaled_g3k6);\n\n const f_k7_h = f(offset + scaled_k7);\n const f_k7_l = f(offset - scaled_k7);\n\n const poorEstimate =\n scaled_weight_g0 * f_g0k0 + // g0\n scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1\n scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2\n scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3\n\n const currentEstimate =\n scaled_weight_k0 * f_g0k0 + // k0\n scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1\n scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2\n scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3\n scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4\n scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5\n scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6\n scaled_weight_k7 * (f_k7_h + f_k7_l); // k7\n\n return [currentEstimate, poorEstimate];\n};\n","// rec-math/src/numerical/quad/index.ts\n\nimport { integrate } from './adaptive-quadrature.js';\nimport { integrationStep } from './gauss-kronrod-g7k15.js';\n\nexport type IntegrandCallback = (a: number) => number;\n\nexport interface QuadratureInfo {\n /** Number of steps _used_ in the integration. */\n steps: number;\n /** Estimate of the global error. */\n errorEstimate: number;\n /**\n * Set to true if there was a problem (e.g. could not reach required accuracy\n * with a step at machine precision).\n */\n isUnreliable?: true;\n /**\n * If the range has multiple (i.e. more than 2) points, contains the results\n * for intermediate ranges.\n */\n points?: [r: number, i: QuadratureInfo][];\n /** The maximum depth used. */\n depth: number;\n}\n\nexport interface QuadratureOptions {\n /** Used to control global error. */\n epsilon?: number;\n /** Maximum depth to use for adaptive step length. */\n maxDepth?: number;\n}\n\nexport type IntegrationStep = (\n /** The integrand. */\n f: (a: number) => number,\n a: number, // Lower limit.\n b: number, // Upper limit.\n) => [a: number, b: number];\n\nconst rangeErrorMessage =\n 'integration range must be an array of at least two endpoints';\n\n/**\n * Numerically compute a definite integral.\n *\n * @param f Callback returning value of integrand.\n * @param range A range of at least 2 endpoints.\n * @param options Options for the computation.\n *\n * @returns The results of the computation.\n */\nexport const quad = (\n f: IntegrandCallback,\n range: number[],\n options: QuadratureOptions = {},\n): [r: number, i: QuadratureInfo] => {\n // Interpret the range.\n if (!Array.isArray(range)) {\n throw new RangeError(rangeErrorMessage);\n }\n\n if (range.length === 2) {\n // Integrate over a single range.\n const [a, b] = range;\n return integrate(integrationStep, f, a, b, options);\n }\n\n if (range.length < 2) {\n // Can't integrate at a point!\n throw new RangeError(rangeErrorMessage);\n }\n\n // Integrate over multiple ranges.\n let result = 0;\n const points: [r: number, i: QuadratureInfo][] = [];\n const info = {\n steps: 0,\n errorEstimate: 0,\n points,\n depth: 0,\n };\n\n for (let i = 0; i < range.length - 1; ++i) {\n const single = integrate(\n integrationStep,\n f,\n range[i],\n range[i + 1],\n options,\n );\n\n info.points.push(single);\n // Update the result.\n result += single[0];\n // Update the cumulative statistics.\n info.steps += single[1].steps;\n info.errorEstimate += single[1].errorEstimate;\n info.depth = Math.max(info.depth, single[1].depth);\n }\n\n return [result, info];\n};\n"],"names":["defaults","epsilon","Number","EPSILON","maxDepth","Infinity","integrate","integrationStep","f","a","b","options","settings","Object","assign","info","steps","errorEstimate","depth","_f","_a","_b","isFinite","t","oneOverOneMinusT","tSquared","oneOverOneMinusTSquared","changeOfVariables","result","integrate_part","acceptableUnitError","Math","abs","currentEstimate","poorEstimate","mid","isUnreliable","safeErrorEstimate","isNaN","leftEstimate","leftErrorEstimate","rightEstimate","rightErrorEstimate","scale","offset","scaled_k1","scaled_g1k2","scaled_k3","scaled_g2k4","scaled_k5","scaled_g3k6","scaled_k7","scaled_weight_g0","scaled_weight_g1","scaled_weight_g2","scaled_weight_g3","scaled_weight_k0","scaled_weight_k1","scaled_weight_k2","scaled_weight_k3","scaled_weight_k4","scaled_weight_k5","scaled_weight_k6","scaled_weight_k7","f_g0k0","f_k1_h","f_k1_l","f_g1k2_h","f_g1k2_l","f_k3_h","f_k3_l","f_g2k4_h","f_g2k4_l","f_k5_h","f_k5_l","f_g3k6_h","f_g3k6_l","rangeErrorMessage","range","Array","isArray","RangeError","length","points","i","single","push","max"],"mappings":";;;;0BAEO,MCQDA,EAAW,CAKfC,QAA0B,GAAjBC,OAAOC,QAKhBC,SAAUC,KAkDCC,EAAY,CACvBC,EACAC,EACAC,EACAC,EACAC,EAA6B,MAG7B,MAAMC,EAAgBC,OAAAC,OAAAD,OAAAC,OAAA,GAAAd,GAAaW,IAC7BV,QAAEA,EAAOG,SAAEA,GAAaQ,EAExBG,EAAO,CAAEC,MAAO,EAAGC,cAAe,EAAGC,MAAO,IAG3CC,EAAIC,EAAIC,GApDS,EACxBb,EACAC,EACAC,KAEA,IAAKY,SAASZ,GACZ,OAAKY,SAASb,GAmBP,CAJKc,IACV,MAAMC,EAAmB,GAAK,EAAID,GAClC,OAAOf,EAAEC,EAAIc,EAAIC,GAAoBA,EAAmBA,GAE9C,EAAG,GAPN,CAVKD,IACV,MAAME,EAAWF,EAAIA,EACfG,EAA0B,GAAK,EAAID,GACzC,OACEjB,EAAEe,EAAIG,IACL,EAAID,GACLC,EACAA,IAGS,EAAG,GAUpB,IAAKJ,SAASb,GAKZ,MAAO,CAHKc,GACHf,EAAEE,GAAK,EAAIa,GAAKA,IAAMA,EAAIA,GAEvB,EAAG,GAEjB,MAAO,CAACf,EAAGC,EAAGC,IAiBOiB,CAAkBnB,EAAGC,EAAGC,GAI7C,IAAKkB,EAAQX,GAAiBY,EAC5BtB,EACAY,EACAC,EACAC,EACA,EACA,EACA,EAACR,OAAAC,OAAA,GACIC,IAIP,MAAMe,EAAsBC,KAAKC,IAAK/B,EAAU2B,GAAWP,EAAKD,IAYhE,OAXCQ,EAAQX,GAAiBY,EACxBtB,EACAY,EACAC,EACAC,EACA,EACAjB,EACA0B,EACAf,GAGK,CAACa,EAAMf,OAAAC,OAAAD,OAAAC,OAAA,GAAOC,GAAM,CAAAE,oBAGvBY,EAAiB,CACrBtB,EACAC,EACAC,EACAC,EACAQ,EACAd,EACA0B,EACAf,KAGA,MAAOkB,EAAiBC,GAAgB3B,EACtCC,EACAC,EACAC,GAGIO,EAAgBc,KAAKC,IAAIE,EAAeD,GAE9C,GAAIf,GAASd,EAGX,QADEW,EAAKC,MACA,CAACiB,EAAiBhB,GAI3B,GAAIA,GADoBa,GAAuBpB,EAAID,GAIjD,QADEM,EAAKC,MACA,CAACiB,EAAiBhB,GAG3B,MAAMkB,GAAO1B,EAAIC,GAAK,EACtB,GAAID,GAAK0B,GAAOA,GAAOzB,EAAG,CAExBK,EAAKqB,cAAe,EACpB,MAAMC,EAAoBC,MAAMrB,GAAiB,EAAIA,EACrD,OAAIqB,MAAML,IAAoBF,KAAKC,IAAIC,KAAqB5B,IACnD,CAAC,EAAGgC,GAEN,CAACJ,EAAiBI,KAIzBnB,EACUH,EAAKG,QACfH,EAAKG,MAAQA,GAEf,MAAOqB,EAAcC,GAAqBX,EACxCtB,EACAC,EACAC,EACA0B,EACAjB,EACAd,EACA0B,EACAf,IAEK0B,EAAeC,GAAsBb,EAC1CtB,EACAC,EACA2B,EACAzB,EACAQ,EACAd,EACA0B,EACAf,GAEF,MAAO,CAACwB,EAAeE,EAAeD,EAAoBE,IC7H/CnC,EAAmC,CAC9CC,EACAC,EACAC,KAGA,MAAMiC,GAASjC,EAAID,GAAK,EAClBmC,GAAUnC,EAAIC,GAAK,EAGnBmC,EApDQ,mBAoDcF,EACtBG,EA5DU,kBA4DgBH,EAC1BI,EApDQ,kBAoDcJ,EACtBK,EA5DU,kBA4DgBL,EAC1BM,EApDQ,kBAoDcN,EACtBO,EA5DU,kBA4DgBP,EAC1BQ,EApDQ,kBAoDcR,EAEtBS,EAnDU,kBAmDqBT,EAC/BU,EAlDU,kBAkDqBV,EAC/BW,EAjDU,mBAiDqBX,EAC/BY,EAhDU,kBAgDqBZ,EAE/Ba,EA/CU,mBA+CqBb,EAC/Bc,EA9CU,mBA8CqBd,EAC/Be,EA7CU,mBA6CqBf,EAC/BgB,EA5CU,kBA4CqBhB,EAC/BiB,EA3CU,mBA2CqBjB,EAC/BkB,EA1CU,mBA0CqBlB,EAC/BmB,EAzCU,mBAyCqBnB,EAC/BoB,EAxCU,oBAwCqBpB,EAE/BqB,EAASxD,EAAEoC,GAEXqB,EAASzD,EAAEoC,EAASC,GACpBqB,EAAS1D,EAAEoC,EAASC,GACpBsB,EAAW3D,EAAEoC,EAASE,GACtBsB,EAAW5D,EAAEoC,EAASE,GAEtBuB,EAAS7D,EAAEoC,EAASG,GACpBuB,EAAS9D,EAAEoC,EAASG,GACpBwB,EAAW/D,EAAEoC,EAASI,GACtBwB,EAAWhE,EAAEoC,EAASI,GAEtByB,EAASjE,EAAEoC,EAASK,GACpByB,EAASlE,EAAEoC,EAASK,GACpB0B,EAAWnE,EAAEoC,EAASM,GACtB0B,EAAWpE,EAAEoC,EAASM,GAqB5B,MAAO,CATLM,EAAmBQ,EACnBP,GAAoBQ,EAASC,GAC7BR,GAAoBS,EAAWC,GAC/BT,GAAoBU,EAASC,GAC7BV,GAAoBW,EAAWC,GAC/BX,GAAoBY,EAASC,GAC7BZ,GAAoBa,EAAWC,GAC/Bb,GAjBavD,EAAEoC,EAASO,GACX3C,EAAEoC,EAASO,IAGxBC,EAAmBY,EACnBX,GAAoBc,EAAWC,GAC/Bd,GAAoBiB,EAAWC,GAC/BjB,GAAoBoB,EAAWC,KC1E7BC,EACJ,wGAWkB,CAClBrE,EACAsE,EACAnE,EAA6B,MAG7B,IAAKoE,MAAMC,QAAQF,GACjB,MAAM,IAAIG,WAAWJ,GAGvB,GAAqB,IAAjBC,EAAMI,OAAc,CAEtB,MAAOzE,EAAGC,GAAKoE,EACf,OAAOxE,EAAUC,EAAiBC,EAAGC,EAAGC,EAAGC,GAG7C,GAAImE,EAAMI,OAAS,EAEjB,MAAM,IAAID,WAAWJ,GAIvB,IAAIjD,EAAS,EACb,MACMb,EAAO,CACXC,MAAO,EACPC,cAAe,EACfkE,OAJ+C,GAK/CjE,MAAO,GAGT,IAAK,IAAIkE,EAAI,EAAGA,EAAIN,EAAMI,OAAS,IAAKE,EAAG,CACzC,MAAMC,EAAS/E,EACbC,EACAC,EACAsE,EAAMM,GACNN,EAAMM,EAAI,GACVzE,GAGFI,EAAKoE,OAAOG,KAAKD,GAEjBzD,GAAUyD,EAAO,GAEjBtE,EAAKC,OAASqE,EAAO,GAAGrE,MACxBD,EAAKE,eAAiBoE,EAAO,GAAGpE,cAChCF,EAAKG,MAAQa,KAAKwD,IAAIxE,EAAKG,MAAOmE,EAAO,GAAGnE,OAG9C,MAAO,CAACU,EAAQb,8BHnGK"}

@@ -1,2 +0,303 @@

export const version = '1.2.0';
export * as integrate from './integrate.js';
/*! @rec-math/math v1.3.0 2022-06-28 23:54:20
* https://github.com/rec-math/ts-math#readme
* Copyright pbuk (https://github.com/pb-uk) MIT license.
*/
// rec-math/src/version.ts
const version = '1.3.0';
// rec-math/src/numerical/quad/adaptive-quadrature.ts
/** Defaults for integration. */
const defaults = {
/**
* Target error estimate for a step: if this is smaller it can lead to
* accumulation of roundoff errors.
*/
epsilon: Number.EPSILON * 16,
/**
* Maximum depth \\( d_{max} \\)for integration: the maximum number of steps
* can be \\( 2^{d_{max}} \\).
* */
maxDepth: Infinity,
};
/**
* Perform a substitution with an appropriate change of variables to deal with
* infinite ranges.
*
* @param f Integrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns An array with any necessary substitution.
*/
const changeOfVariables = (f, a, b) => {
if (!isFinite(b)) {
if (!isFinite(a)) {
// Change variables to integrate between - and + infinity.
const _f = (t) => {
const tSquared = t * t;
const oneOverOneMinusTSquared = 1 / (1 - tSquared);
return (f(t * oneOverOneMinusTSquared) *
(1 + tSquared) *
oneOverOneMinusTSquared *
oneOverOneMinusTSquared);
};
return [_f, -1, 1];
}
// Change variables to integrate up to infinity.
const _f = (t) => {
const oneOverOneMinusT = 1 / (1 - t);
return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;
};
return [_f, 0, 1];
}
if (!isFinite(a)) {
// Change variables to integrate up from negative infinity.
const _f = (t) => {
return f(b - (1 - t) / t) / (t * t);
};
return [_f, 0, 1];
}
return [f, a, b];
};
const integrate = (integrationStep, f, a, b, options = {}) => {
// Establish settings.
const settings = Object.assign(Object.assign({}, defaults), options);
const { epsilon, maxDepth } = settings;
const info = { steps: 0, errorEstimate: 0, depth: 1 };
// Allow for a change of variables to deal with infinite limits.
const [_f, _a, _b] = changeOfVariables(f, a, b);
// Get estimate so we can work out the acceptable global error.
// Use a depth of 1 to calculate a 15 point Kronrod quadrature.
let [result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.
_a, // Lower limit.
_b, // Upper limit.
1, // New depth.
1, // Maximum depth.
0, Object.assign({}, info));
// Now calculate using the target global error.
const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));
[result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.
_a, // Lower limit.
_b, // Upper limit.
1, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
return [result, Object.assign(Object.assign({}, info), { errorEstimate })];
};
const integrate_part = (integrationStep, f, a, // Lower limit.
b, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info) => {
// Initialize things.
const [currentEstimate, poorEstimate] = integrationStep(f, // Integrand.
a, // Lower limit.
b);
const errorEstimate = Math.abs(poorEstimate - currentEstimate);
if (depth >= maxDepth) {
// Reached the maximum allowable depth so return the partial sum.
++info.steps;
return [currentEstimate, errorEstimate];
}
const acceptableError = acceptableUnitError * (b - a);
if (errorEstimate <= acceptableError) {
// Error is acceptable for the size of step so return the partial sum.
++info.steps;
return [currentEstimate, errorEstimate];
}
const mid = (a + b) / 2;
if (a >= mid || mid >= b) {
// We can't make this step any smaller: looks like a discontinuity.
info.isUnreliable = true;
const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;
if (isNaN(currentEstimate) || Math.abs(currentEstimate) === Infinity) {
return [0, safeErrorEstimate];
}
return [currentEstimate, safeErrorEstimate];
}
// Recurse deeper.
++depth;
if (depth > info.depth) {
info.depth = depth;
}
const [leftEstimate, leftErrorEstimate] = integrate_part(integrationStep, f, // Integrand.
a, // Lower limit.
mid, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
const [rightEstimate, rightErrorEstimate] = integrate_part(integrationStep, f, // Integrand.
mid, // Lower limit.
b, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
return [leftEstimate + rightEstimate, leftErrorEstimate + rightErrorEstimate];
};
// rec-math/src/numerical/quad/gauss-kronrod-g7k15.ts
// Gauss-Kronrod constants (G7, K15) on [-1, 1].
// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/
// node_g0k0 = 0;
// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;
const node_g1k2 = 0.4058451513773972;
// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;
const node_g2k4 = 0.7415311855993945;
// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;
const node_g3k6 = 0.9491079123427585;
// node_k1 = 2.077849550078984676006894037732449e-01 exact;
const node_k1 = 0.20778495500789848;
// node_k3 = 5.860872354676911302941448382587296e-01 exact;
const node_k3 = 0.5860872354676911;
// node_k5 = 8.648644233597690727897127886409262e-01 exact;
const node_k5 = 0.8648644233597691;
// node_k7 = 9.914553711208126392068546975263285e-01 exact;
const node_k7 = 0.9914553711208126;
// weight_g0 = 4.179591836734693877551020408163265e-01 exact;
const weight_g0 = 0.4179591836734694;
// weight_g1 = 3.818300505051189449503697754889751e-01 exact;
const weight_g1 = 0.3818300505051189;
// weight_g2 = 2.797053914892766679014677714237796e-01 exact;
const weight_g2 = 0.27970539148927664;
// weight_g3 = 1.294849661688696932706114326790820e-01 exact;
const weight_g3 = 0.1294849661688697;
// weight_k0 = 2.094821410847278280129991748917143e-01 exact;
const weight_k0 = 0.20948214108472782;
// weight_k1 = 2.044329400752988924141619992346491e-01 exact;
const weight_k1 = 0.20443294007529889;
// weight_k2 = 1.903505780647854099132564024210137e-01 exact;
const weight_k2 = 0.19035057806478542;
// weight_k3 = 1.690047266392679028265834265985503e-01 exact;
const weight_k3 = 0.1690047266392679;
// weight_k4 = 1.406532597155259187451895905102379e-01 exact;
const weight_k4 = 0.14065325971552592;
// weight_k5 = 1.047900103222501838398763225415180e-01 exact;
const weight_k5 = 0.10479001032225019;
// weight_k6 = 6.309209262997855329070066318920429e-02 exact;
const weight_k6 = 0.06309209262997856;
// weight_k7 = 2.293532201052922496373200805896959e-02 exact;
const weight_k7 = 0.022935322010529224;
/**
* Perform a single integration step.
*
* @param f Integrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns [current best estimate, poor estimate]
*/
const integrationStep = (f, // Integrand.
a, // Lower limit.
b) => {
// Gauss-Kronrod (G7, K15).
const scale = (b - a) / 2;
const offset = (a + b) / 2;
// const scaled_g0k0 = 0;
const scaled_k1 = node_k1 * scale;
const scaled_g1k2 = node_g1k2 * scale;
const scaled_k3 = node_k3 * scale;
const scaled_g2k4 = node_g2k4 * scale;
const scaled_k5 = node_k5 * scale;
const scaled_g3k6 = node_g3k6 * scale;
const scaled_k7 = node_k7 * scale;
const scaled_weight_g0 = weight_g0 * scale;
const scaled_weight_g1 = weight_g1 * scale;
const scaled_weight_g2 = weight_g2 * scale;
const scaled_weight_g3 = weight_g3 * scale;
const scaled_weight_k0 = weight_k0 * scale;
const scaled_weight_k1 = weight_k1 * scale;
const scaled_weight_k2 = weight_k2 * scale;
const scaled_weight_k3 = weight_k3 * scale;
const scaled_weight_k4 = weight_k4 * scale;
const scaled_weight_k5 = weight_k5 * scale;
const scaled_weight_k6 = weight_k6 * scale;
const scaled_weight_k7 = weight_k7 * scale;
const f_g0k0 = f(offset);
const f_k1_h = f(offset + scaled_k1);
const f_k1_l = f(offset - scaled_k1);
const f_g1k2_h = f(offset + scaled_g1k2);
const f_g1k2_l = f(offset - scaled_g1k2);
const f_k3_h = f(offset + scaled_k3);
const f_k3_l = f(offset - scaled_k3);
const f_g2k4_h = f(offset + scaled_g2k4);
const f_g2k4_l = f(offset - scaled_g2k4);
const f_k5_h = f(offset + scaled_k5);
const f_k5_l = f(offset - scaled_k5);
const f_g3k6_h = f(offset + scaled_g3k6);
const f_g3k6_l = f(offset - scaled_g3k6);
const f_k7_h = f(offset + scaled_k7);
const f_k7_l = f(offset - scaled_k7);
const poorEstimate = scaled_weight_g0 * f_g0k0 + // g0
scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1
scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2
scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3
const currentEstimate = scaled_weight_k0 * f_g0k0 + // k0
scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1
scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2
scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3
scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4
scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5
scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6
scaled_weight_k7 * (f_k7_h + f_k7_l); // k7
return [currentEstimate, poorEstimate];
};
// rec-math/src/numerical/quad/index.ts
const rangeErrorMessage = 'integration range must be an array of at least two endpoints';
/**
* Numerically compute a definite integral.
*
* @param f Callback returning value of integrand.
* @param range A range of at least 2 endpoints.
* @param options Options for the computation.
*
* @returns The results of the computation.
*/
const quad = (f, range, options = {}) => {
// Interpret the range.
if (!Array.isArray(range)) {
throw new RangeError(rangeErrorMessage);
}
if (range.length === 2) {
// Integrate over a single range.
const [a, b] = range;
return integrate(integrationStep, f, a, b, options);
}
if (range.length < 2) {
// Can't integrate at a point!
throw new RangeError(rangeErrorMessage);
}
// Integrate over multiple ranges.
let result = 0;
const points = [];
const info = {
steps: 0,
errorEstimate: 0,
points,
depth: 0,
};
for (let i = 0; i < range.length - 1; ++i) {
const single = integrate(integrationStep, f, range[i], range[i + 1], options);
info.points.push(single);
// Update the result.
result += single[0];
// Update the cumulative statistics.
info.steps += single[1].steps;
info.errorEstimate += single[1].errorEstimate;
info.depth = Math.max(info.depth, single[1].depth);
}
return [result, info];
};
// rec-math/src/numerical/index.ts
var index = /*#__PURE__*/Object.freeze({
__proto__: null,
quad: quad
});
export { index as numerical, version };
//# sourceMappingURL=index.js.map
{
"name": "@rec-math/math",
"version": "1.2.0",
"version": "1.3.0",
"description": "Mathematics for the browser (and TypeScript, JavaScript).",
"browser": "dist/rec-math.min.js",
"module": "esm/index.js",
"type": "module",
"sideEffects": false,
"types": "types/index.d.ts",
"files": [

@@ -15,6 +11,13 @@ "dist",

],
"type": "module",
"browser": "dist/rec-math.min.js",
"module": "esm/index.js",
"exports": {
"import": "./esm/index.js"
},
"types": "types",
"scripts": {
"build": "rimraf dist esm types && npm run lint && npm run build:esm && rollup -c && npm run test:dist",
"build:esm": "tsc --project tsconfig.types.json && tsc --project tsconfig.build.json",
"build": "rimraf dist esm types && npm run lint && rollup -c && tsc --project tsconfig.types.json && npm run test:dist",
"ci:build": "npm run lint && npm run test:unit && npm run build",
"docs": "typedoc",
"lint": "eslint . && prettier . --check",

@@ -38,2 +41,3 @@ "lint:fix": "eslint . --fix && prettier . --write",

"devDependencies": {
"@rollup/plugin-typescript": "^8.3.3",
"@types/chai": "^4.3.1",

@@ -52,4 +56,5 @@ "@types/mocha": "^9.1.1",

"ts-node": "^10.8.1",
"typedoc": "^0.23.2",
"typescript": "^4.7.3"
}
}
+14
-12

@@ -26,8 +26,8 @@ # RecMath mathematics module.

```Javascript
const [result, info] = RecMath.integrate.quad(
(x) => Math.exp(x), // A function to integrate.
[0, Number.POSITIVE_INFINITY] // A range to integrate over.
const [result, info] = RecMath.numerical.quad(
(x) => Math.exp(-x), // A function to integrate.
[0, Infinity], // A range to integrate over.
);
console log(result); // 1
console log(info);
console.log(result); // 1
console.log(info);
// { steps: 14, errorEstimate: 3.384539692172424e-16, depth: 7 }

@@ -39,19 +39,21 @@ ```

```Javascript
const [, { points }] = RecMath.integrate.quad(
const [result, { steps, points }] = RecMath.numerical.quad(
// Normal distribution.
(t) => Math.exp(-0.5 * t * t) / Math.sqrt(2 * Math.PI),
[Number.NEGATIVE_INFINITY, -3, -2, -1, 1, 2, 3, Number.POSITIVE_INFINITY]
[-Infinity, -3, -2, -1, 1, 2, 3, Infinity],
);
// 99.7%, 96% and 68% confidence intervals.
console.log(result, steps); // { result: 1, steps: 32 }
// 68%, 96% and 99.7% confidence intervals.
const oneSigma = points[3][0];
const twoSigma = oneSigma + points[2][0] + points[4][0];
const threeSigma = 1 - points[0][0] - points[6][0];
const twoSigma = threeSigma - points[1][0] - points[5][0];
const oneSigma = twoSigma - points[2][0] - points[4][0];
console.log({ oneSigma, twoSigma, threeSigma });
// {
// oneSigma: 0.6826894921370859,
// twoSigma: 0.9544997361036416,
// oneSigma: 0.682689492137086,
// twoSigma: 0.9544997361036417,
// threeSigma: 0.9973002039367398
// }
```

@@ -1,67 +0,3 @@

declare module "integrate/quad/adaptive-quadrature" {
import type { IntegrandCallback, IntegrationStep, QuadratureInfo, QuadratureOptions } from "integrate/quad";
export { quadrature };
const quadrature: (integrationStep: IntegrationStep, f: IntegrandCallback, a: number, b: number, options?: QuadratureOptions) => [r: number, i: QuadratureInfo];
}
declare module "integrate/quad/gauss-kronrod-g7k15" {
import type { IntegrationStep } from "integrate/quad";
export { integrationStep };
/**
* Perform a single integration step.
*
* @param f Integrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns [current best estimate, poor estimate]
*/
const integrationStep: IntegrationStep;
}
declare module "integrate/quad" {
export type IntegrandCallback = (a: number) => number;
export interface QuadratureInfo {
/** Number of steps _used_ in the integration. */
steps: number;
/** Estimate of the global error. */
errorEstimate: number;
/**
* Set to true if there was a problem (e.g. could not reach required accuracy
* with a step at machine precision).
*/
isUnreliable?: true;
/**
* If the range has multiple (i.e. more than 2) points, contains the results
* for intermediate ranges.
*/
points?: [r: number, i: QuadratureInfo][];
/** The maximum depth used. */
depth: number;
}
export interface QuadratureOptions {
/** Used to control global error. */
epsilon?: number;
/** Maximum depth to use for adaptive step length. */
maxDepth?: number;
}
export type IntegrationStep = (
/** The integrand. */
f: (a: number) => number, a: number, // Lower limit.
b: number) => [a: number, b: number];
/**
* Numerically compute a definite integral.
*
* @param f Callback returning value of integrand.
* @param range A range of at least 2 endpoints.
* @param options Options for the computation.
*
* @returns The results of the computation.
*/
export const quad: (f: IntegrandCallback, range: number[], options?: QuadratureOptions) => [r: number, i: QuadratureInfo];
}
declare module "integrate" {
export { quad } from "integrate/quad";
}
declare module "index" {
export const version = "1.2.0";
export * as integrate from "integrate";
}
export { version } from './version.js';
export * as numerical from './numerical/index.js';
//# sourceMappingURL=index.d.ts.map

@@ -1,1 +0,1 @@

{"version":3,"file":"index.d.ts","sourceRoot":"","sources":["../src/integrate/quad/adaptive-quadrature.ts","../src/integrate/quad/gauss-kronrod-g7k15.ts","../src/integrate/quad.ts","../src/integrate.ts","../src/index.ts"],"names":[],"mappings":";IAAA,OAAO,KAAK,EACV,iBAAiB,EACjB,eAAe,EACf,cAAc,EACd,iBAAiB,EAClB,uBAAgB;IAGjB,OAAO,EAAE,UAAU,EAAE,CAAC;IAsDtB,MAAM,UAAU,oBACG,eAAe,KAC7B,iBAAiB,KACjB,MAAM,KACN,MAAM,YACA,iBAAiB,KACzB,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CAqC/B,CAAC;;;ICzGF,OAAO,KAAK,EAAqB,eAAe,EAAE,uBAAgB;IAGlE,OAAO,EAAE,eAAe,EAAE,CAAC;IAgD3B;;;;;;;OAOG;IACH,MAAM,eAAe,EAAE,eAqEtB,CAAC;;;IC7HF,MAAM,MAAM,iBAAiB,GAAG,CAAC,CAAC,EAAE,MAAM,KAAK,MAAM,CAAC;IAEtD,MAAM,WAAW,cAAc;QAC7B,iDAAiD;QACjD,KAAK,EAAE,MAAM,CAAC;QACd,oCAAoC;QACpC,aAAa,EAAE,MAAM,CAAC;QACtB;;;WAGG;QACH,YAAY,CAAC,EAAE,IAAI,CAAC;QACpB;;;WAGG;QACH,MAAM,CAAC,EAAE,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CAAC,EAAE,CAAC;QAC1C,8BAA8B;QAC9B,KAAK,EAAE,MAAM,CAAC;KACf;IAED,MAAM,WAAW,iBAAiB;QAChC,oCAAoC;QACpC,OAAO,CAAC,EAAE,MAAM,CAAC;QACjB,qDAAqD;QACrD,QAAQ,CAAC,EAAE,MAAM,CAAC;KACnB;IAED,MAAM,MAAM,eAAe,GAAG;IAC5B,qBAAqB;IACrB,CAAC,EAAE,CAAC,CAAC,EAAE,MAAM,KAAK,MAAM,EACxB,CAAC,EAAE,MAAM,EAAE,eAAe;IAC1B,CAAC,EAAE,MAAM,KACN,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,MAAM,CAAC,CAAC;IAK5B;;;;;;;;OAQG;IACH,MAAM,CAAC,MAAM,IAAI,MACZ,iBAAiB,SACb,MAAM,EAAE,YACN,iBAAiB,KACzB,CAAC,CAAC,EAAE,MAAM,EAAE,CAAC,EAAE,cAAc,CA8C/B,CAAC;;;ICnGF,OAAO,EAAE,IAAI,EAAE,uBAA4B;;;ICD3C,MAAM,CAAC,MAAM,OAAO,UAAU,CAAC;IAE/B,OAAO,KAAK,SAAS,kBAAuB"}
{"version":3,"file":"index.d.ts","sourceRoot":"","sources":["../src/index.ts"],"names":[],"mappings":"AAEA,OAAO,EAAE,OAAO,EAAE,MAAM,cAAc,CAAC;AAEvC,OAAO,KAAK,SAAS,MAAM,sBAAsB,CAAC"}
// Export the default method.
export { quad } from './integrate/quad.js';
import { quadrature } from './quad/adaptive-quadrature.js';
import { integrationStep } from './quad/gauss-kronrod-g7k15.js';
const rangeErrorMessage = 'integration range must be an array of at least two endpoints';
/**
* Numerically compute a definite integral.
*
* @param f Callback returning value of integrand.
* @param range A range of at least 2 endpoints.
* @param options Options for the computation.
*
* @returns The results of the computation.
*/
export const quad = (f, range, options = {}) => {
// Interpret the range.
if (!Array.isArray(range)) {
throw new RangeError(rangeErrorMessage);
}
if (range.length === 2) {
// Integrate over a single range.
const [a, b] = range;
return quadrature(integrationStep, f, a, b, options);
}
if (range.length < 2) {
// Can't integrate at a point!
throw new RangeError(rangeErrorMessage);
}
// Integrate over multiple ranges.
let result = 0;
const points = [];
const info = {
steps: 0,
errorEstimate: 0,
points,
depth: 0,
};
for (let i = 0; i < range.length - 1; ++i) {
const single = quadrature(integrationStep, f, range[i], range[i + 1], options);
info.points.push(single);
// Update the result.
result += single[0];
// Update the cumulative statistics.
info.steps += single[1].steps;
info.errorEstimate += single[1].errorEstimate;
info.depth = Math.max(info.depth, single[1].depth);
}
return [result, info];
};
// Export the API.
export { quadrature };
const defaults = {
epsilon: Number.EPSILON * 16,
maxDepth: Number.POSITIVE_INFINITY,
};
/**
* Perform a substitution with an appropriate change of variables to deal with
* infinite ranges.
*
* @param f Intgrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns An array with any necessary substitution.
*/
const changeOfVariables = (f, a, b) => {
if (!isFinite(b)) {
if (!isFinite(a)) {
// Change variables to integrate between - and + infinity.
const _f = (t) => {
const tSquared = t * t;
const oneOverOneMinusTSquared = 1 / (1 - tSquared);
return (f(t * oneOverOneMinusTSquared) *
(1 + tSquared) *
oneOverOneMinusTSquared *
oneOverOneMinusTSquared);
};
return [_f, -1, 1];
}
// Change variables to integrate up to infinity.
const _f = (t) => {
const oneOverOneMinusT = 1 / (1 - t);
return f(a + t * oneOverOneMinusT) * oneOverOneMinusT * oneOverOneMinusT;
};
return [_f, 0, 1];
}
if (!isFinite(a)) {
// Change variables to integrate up from negative infinity.
const _f = (t) => {
return f(b - (1 - t) / t) / (t * t);
};
return [_f, 0, 1];
}
return [f, a, b];
};
const quadrature = (integrationStep, f, a, b, options = {}) => {
// Establish settings.
const settings = Object.assign(Object.assign({}, defaults), options);
const { epsilon, maxDepth } = settings;
const info = { steps: 0, errorEstimate: 0, depth: 1 };
// Allow for a change of variables to deal with infinite limits.
const [_f, _a, _b] = changeOfVariables(f, a, b);
// Get estimate so we can work out the acceptable global error.
// Use a depth of 1 to calculate a 15 point Kronrod quadrature.
let [result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.
_a, // Lower limit.
_b, // Upper limit.
1, // New depth.
1, // Maximum depth.
0, Object.assign({}, info));
// Now calculate using the target global error.
const acceptableUnitError = Math.abs((epsilon * result) / (_b - _a));
[result, errorEstimate] = integrate_part(integrationStep, _f, // Integrand.
_a, // Lower limit.
_b, // Upper limit.
1, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
return [result, Object.assign(Object.assign({}, info), { errorEstimate })];
};
const integrate_part = (integrationStep, f, a, // Lower limit.
b, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info) => {
// Initialize things.
const estimates = integrationStep(f, // Integrand.
a, // Lower limit.
b);
let currentEstimate = estimates[0];
const poorEstimate = estimates[1];
let errorEstimate = Math.abs(poorEstimate - currentEstimate);
if (depth >= maxDepth) {
// Reached the maximum allowable depth so return the partial sum.
++info.steps;
return [currentEstimate, errorEstimate];
}
const acceptableError = acceptableUnitError * (b - a);
if (errorEstimate <= acceptableError) {
// Error is acceptable for the size of step so return the partial sum.
++info.steps;
return [currentEstimate, errorEstimate];
}
const mid = (a + b) / 2;
if (a >= mid || mid >= b) {
// We can't make this step any smaller: looks like a discontinuity.
info.isUnreliable = true;
const safeErrorEstimate = isNaN(errorEstimate) ? 0 : errorEstimate;
if (currentEstimate === Number.POSITIVE_INFINITY) {
return [0, safeErrorEstimate];
}
return [currentEstimate, safeErrorEstimate];
}
// Recurse deeper.
++depth;
if (depth > info.depth) {
info.depth = depth;
}
const leftResult = integrate_part(integrationStep, f, // Integrand.
a, // Lower limit.
mid, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
const rightResult = integrate_part(integrationStep, f, // Integrand.
mid, // Lower limit.
b, // Upper limit.
depth, // New depth.
maxDepth, // Maximum depth.
acceptableUnitError, // Acceptable error per unit step.
info);
currentEstimate = leftResult[0] + rightResult[0];
errorEstimate = leftResult[1] + rightResult[1];
return [currentEstimate, errorEstimate];
};
// Export the API.
export { integrationStep };
// Gauss-Kronrod constants (G7, K15) on [-1, 1].
// https://www.advanpix.com/2011/11/07/gauss-kronrod-quadrature-nodes-weights/
// node_g0k0 = 0;
// node_g1k2 = 4.058451513773971669066064120769615e-01 exact;
const node_g1k2 = 0.4058451513773972;
// node_g2k4 = 7.415311855993944398638647732807884e-01 exact;
const node_g2k4 = 0.7415311855993945;
// node_g3k6 = 9.491079123427585245261896840478513e-01 exact;
const node_g3k6 = 0.9491079123427585;
// node_k1 = 2.077849550078984676006894037732449e-01 exact;
const node_k1 = 0.20778495500789848;
// node_k3 = 5.860872354676911302941448382587296e-01 exact;
const node_k3 = 0.5860872354676911;
// node_k5 = 8.648644233597690727897127886409262e-01 exact;
const node_k5 = 0.8648644233597691;
// node_k7 = 9.914553711208126392068546975263285e-01 exact;
const node_k7 = 0.9914553711208126;
// weight_g0 = 4.179591836734693877551020408163265e-01 exact;
const weight_g0 = 0.4179591836734694;
// weight_g1 = 3.818300505051189449503697754889751e-01 exact;
const weight_g1 = 0.3818300505051189;
// weight_g2 = 2.797053914892766679014677714237796e-01 exact;
const weight_g2 = 0.27970539148927664;
// weight_g3 = 1.294849661688696932706114326790820e-01 exact;
const weight_g3 = 0.1294849661688697;
// weight_k0 = 2.094821410847278280129991748917143e-01 exact;
const weight_k0 = 0.20948214108472782;
// weight_k1 = 2.044329400752988924141619992346491e-01 exact;
const weight_k1 = 0.20443294007529889;
// weight_k2 = 1.903505780647854099132564024210137e-01 exact;
const weight_k2 = 0.19035057806478542;
// weight_k3 = 1.690047266392679028265834265985503e-01 exact;
const weight_k3 = 0.1690047266392679;
// weight_k4 = 1.406532597155259187451895905102379e-01 exact;
const weight_k4 = 0.14065325971552592;
// weight_k5 = 1.047900103222501838398763225415180e-01 exact;
const weight_k5 = 0.10479001032225019;
// weight_k6 = 6.309209262997855329070066318920429e-02 exact;
const weight_k6 = 0.06309209262997856;
// weight_k7 = 2.293532201052922496373200805896959e-02 exact;
const weight_k7 = 0.022935322010529224;
/**
* Perform a single integration step.
*
* @param f Integrand callback.
* @param a Lower limit.
* @param b Upper limit.
* @returns [current best estimate, poor estimate]
*/
const integrationStep = (f, // Integrand.
a, // Lower limit.
b) => {
// Gauss-Kronrod (G7, K15).
const scale = (b - a) / 2;
const offset = (a + b) / 2;
// const scaled_g0k0 = 0;
const scaled_k1 = node_k1 * scale;
const scaled_g1k2 = node_g1k2 * scale;
const scaled_k3 = node_k3 * scale;
const scaled_g2k4 = node_g2k4 * scale;
const scaled_k5 = node_k5 * scale;
const scaled_g3k6 = node_g3k6 * scale;
const scaled_k7 = node_k7 * scale;
const scaled_weight_g0 = weight_g0 * scale;
const scaled_weight_g1 = weight_g1 * scale;
const scaled_weight_g2 = weight_g2 * scale;
const scaled_weight_g3 = weight_g3 * scale;
const scaled_weight_k0 = weight_k0 * scale;
const scaled_weight_k1 = weight_k1 * scale;
const scaled_weight_k2 = weight_k2 * scale;
const scaled_weight_k3 = weight_k3 * scale;
const scaled_weight_k4 = weight_k4 * scale;
const scaled_weight_k5 = weight_k5 * scale;
const scaled_weight_k6 = weight_k6 * scale;
const scaled_weight_k7 = weight_k7 * scale;
const f_g0k0 = f(offset);
const f_k1_h = f(offset + scaled_k1);
const f_k1_l = f(offset - scaled_k1);
const f_g1k2_h = f(offset + scaled_g1k2);
const f_g1k2_l = f(offset - scaled_g1k2);
const f_k3_h = f(offset + scaled_k3);
const f_k3_l = f(offset - scaled_k3);
const f_g2k4_h = f(offset + scaled_g2k4);
const f_g2k4_l = f(offset - scaled_g2k4);
const f_k5_h = f(offset + scaled_k5);
const f_k5_l = f(offset - scaled_k5);
const f_g3k6_h = f(offset + scaled_g3k6);
const f_g3k6_l = f(offset - scaled_g3k6);
const f_k7_h = f(offset + scaled_k7);
const f_k7_l = f(offset - scaled_k7);
const poorEstimate = scaled_weight_g0 * f_g0k0 + // g0
scaled_weight_g1 * (f_g1k2_h + f_g1k2_l) + // g1
scaled_weight_g2 * (f_g2k4_h + f_g2k4_l) + // g2
scaled_weight_g3 * (f_g3k6_h + f_g3k6_l); // g3
const currentEstimate = scaled_weight_k0 * f_g0k0 + // k0
scaled_weight_k1 * (f_k1_h + f_k1_l) + // k1
scaled_weight_k2 * (f_g1k2_h + f_g1k2_l) + // k2
scaled_weight_k3 * (f_k3_h + f_k3_l) + // k3
scaled_weight_k4 * (f_g2k4_h + f_g2k4_l) + // k4
scaled_weight_k5 * (f_k5_h + f_k5_l) + // k5
scaled_weight_k6 * (f_g3k6_h + f_g3k6_l) + // k6
scaled_weight_k7 * (f_k7_h + f_k7_l); // k7
return [currentEstimate, poorEstimate];
};