用Javascript实现了函数积分的算法。
从计算结果可以看出,算法精度受随机数发生器的影响较大。
//Monte carro求给定函数的定积分 //给定一个函数 function f1(x, xmin, xmax) { if (x >= xmin && x <= xmax) return x * x * Math.sqrt(1 + x * x * x) - 1; } function f2(x, xmin, xmax) { if (x >= xmin && x <= xmax) return x * Math.log(1 + x); } function f3(x, xmin, xmax) { if (x >= xmin && x <= xmax) return Math.sqrt(1 - x * x); } function f4(x, xmin, xmax) { if (x >= xmin && x <= xmax) return Math.cos(x) - 0.9; } //取得函数的上下极限——其精确度影响很大---由于用于积分,故若y>0,底限0;若y<0,则高限0;若y穿越x轴,则不修改默认极限 function getLimit(func, xmin, xmax) { //边界极限 var t1 = func(xmin, xmin, xmax), t2 = func(xmax, xmin, xmax); var limit = {y0: Math.min(t1, t2), y1: Math.max(t1, t2)}; var delta = 1e-3, lastX, lastY, y; for(var x = xmin, i = 0; x <= xmax; x += delta, i++) { y = func(x, xmin, xmax); if (y < limit.y0) limit.y0 = y; if (y > limit.y1) limit.y1 = y; if (i == 0) { lastX = x; lastY = y; } else { var tx = (lastX + x) / 2; var ty = func(tx, xmin, xmax); if (Math.abs(lastY) > Math.abs(y)) { if (Math.abs(lastY) < Math.abs(ty)) { lastX = tx; lastY = ty; if (ty < 0 && limit.y0 > ty) limit.y0 = ty; else if (ty >= 0 && limit.y1 <= ty) limit.y1 = ty; } } if (Math.abs(lastY) < Math.abs(y)) { if (Math.abs(y) < Math.abs(ty)) { lastX = tx; lastY = ty; if (ty < 0 && limit.y0 > ty) limit.y0 = ty; else if (ty >= 0 && limit.y1 <= ty) limit.y1 = ty; } } } } if (limit.y0 > 0) limit.y0 = 0; if (limit.y1 < 0) limit.y1 = 0; return limit; } //积分 function cumulate(func, N, xmin, xmax) { var limit = getLimit(func, xmin, xmax); //console.log('Limit=', limit.y0, limit.y1); var x, yr, y; var cntA = 0, cntB = 0; for(var i = 0; i < N; i++) { x = xmin + (xmax - xmin) * Math.random(); yr = limit.y0 + (limit.y1 - limit.y0) * Math.random(); y = func(x, xmin, xmax); if ( y >= 0 && yr >= 0 && yr <= y) cntA++; else if (y <= 0 && yr <= 0 && yr >= y) cntB++; } //console.log('counts=', cntA, cntB); return (cntA - cntB) * (xmax - xmin) * (limit.y1 - limit.y0) / N; } //测试 var t1, t2; var rst; var N = 1e7; var fs = [f1, f2, f3, f4]; var ans = [-0.593682861167513, 0.25, Math.PI/4, -0.058529015192103493347497678369701]; for(var i = 0; i < 4; i++) { t1 = new Date(); rst = cumulate(fs[i], N, 0, 1); t2 = new Date(); console.log('[', i, ']=', ans[i].toFixed(15), rst.toFixed(15), t2 - t1, 'ms'); }