@rec-math/math
Advanced tools
| /*! @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"} |
@@ -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"} |
+303
-2
@@ -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 |
+12
-7
| { | ||
| "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 | ||
| // } | ||
| ``` |
+2
-66
@@ -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]; | ||
| }; |
Major refactor
Supply chain riskPackage has recently undergone a major refactor. It may be unstable or indicate significant internal changes. Use caution when updating to versions that include significant changes.
Found 1 instance in 1 package
Minified code
QualityThis package contains minified code. This may be harmless in some cases where minified code is included in packaged libraries, however packages on npm should not minify code.
Found 1 instance in 1 package
83538
114.38%21
75%382
4.66%58
3.57%16
14.29%2
100%1
Infinity%