算法 {(自适应)辛普森积分}

算法 {(自适应)辛普森积分}
@LOC: 2

辛普森積分

定義

若函數$F(x)$ 在$[l,r]$上可積, 令$mid=(l+r)/2$, $G(x)$為經過$(l,f(l)),(mid,f(mid)),(r,f(r))$三點的二次方程(三點確定一條拋物線 (如果他们共线 则该抛物线会退化为直线)), 則用${\int }_l^{r}G(x)dx$來近似${\int }_l^{r}F(x)dx$;
. 总之, 是用(抛物线/直线)来替代该函数图像;

算法

自適應辛普森積分

錯誤優化

如果超時了, 可以將double_epsilon設置小一點;
如果出錯了, 估計是發生了@LINK: @LOC_0, 此時將空間旋轉一個角度;

筆記

我們將F函數 設置為分段函數 比如是[l0,r0] ... [ln,rn] 然後對每個子區間[li,ri]進行積分 累加起來就是答案;
那麼 你可能認為: 反正就是求這些子區間的面積的和, 那麼對於沒有定義的點 設置為0, 然後直接對[l0, rn]進行積分 更方便呀;
. 這個思路 從積分角度 本質是完全正確的 沒有問題, 但是 由於@MARK:@LOC_1這種情況 (這是完全由於自適應算法本身的問題 與你的思路無關) 因此對於算法而言 是不可行的, 但思路沒問題;

@DELI;

比如你要对$f(x)$在[l,r]上进行积分, 请确保 函数在该区间是可积的;
从算法角度, 请确保 所有的调用f(x0) 他是有定义的 值是有界的;
比如 求若干个圆的并集, 你需要先將他倆投影到X軸上, 然後用區間合併 會得到[l0, r0] [l1,r1] ... [ln.rn], 不可以直接對[l0, rn]進行積分 因為他裡面存在未定義的點;

@DELI;

@MARK: @LOC_0;
想象一下 第一次计算[l, r](他是三个点l,m,r) 第二次计算[l, ml, m, mr, r] 他俩的面积一样(误差极小) 即l,ml,m,mr,r这5个点 都在同一个抛物线上(得到的都是同一个抛物线P), 于是我们认为答案确定了;
但是 函数本身图像 并不是抛物线P, 即函数在(l, ml) (ml,m) (m,mr) (mr,r)的里面 有很多拐点, 此时这个算法是出错了;

但是 此时也说明 函数本身 他的阶数非常高, 一般不会遇到函数阶数特别高的情况;
因此 上述这种情况, 即这5个点 是在同一个抛物线上 一般不会出现, 假如题目给定的图像 都是阶数比较小的;

@MARK: @LOC_1;
但是有一种情况是常见的: l,ml,m,mr,r這5個點 是共線的 但是圖像本身不是一個直線 (此時辛普森積分的拋物線 退化為直線來擬合), 此時也會出錯;
. 举个例子, 对于圆的并的面积, 圆心为(0,0) (2,0) (4,0) (6,0) 半径都为1, 此时l=-1, ml=1, m=3, mr= 5, r=7 他们是共线的, 得到面积均为0, 就出错了;
. 此时的做法是: 将空间旋转一个角度后 (即将这些圆进行旋转) 这个操作 不会影响圆的并集的面积, 此时的f(x) 除了l|r端点 中间不会有零点;

代碼

class __SimpsonIntegration{
public:
    Tools::Double (* Ptr_func)( Tools::Double);

    __SimpsonIntegration(){}
    inline Tools::Double __Get_simpson( Tools::Double _l, Tools::Double _r){ return (Ptr_func(_l) + Ptr_func( (_l+_r)/2) * 4 + Ptr_func(_r)) / 6 * (_r - _l);}
    Tools::Double __Dfs( Tools::Double _l, Tools::Double _r, Tools::Double _cur){
        Tools::Double mid = ( _l + _r) / 2;
        Tools::Double nex_l = __Get_simpson( _l, mid), nex_r = __Get_simpson( mid, _r);
        if( Double_cmp( nex_l + nex_r, _cur) == 0){ return nex_l + nex_r;}
        auto ANS = __Dfs( _l, mid, nex_l);
        ANS += __Dfs( mid, _r, nex_r);
        return ANS;
    }
    Tools::Double Get_integration( Tools::Double _l, Tools::Double _r){
        return __Dfs( _l, _r, __Get_simpson( _l, _r));
    }
}; // class __SimpsonIntegration
__SimpsonIntegration SimpsonIntegration;
@TODO( SimpsonIntegration.Ptr_func = ?);

int main(){
	求若干個圖像S的並集的面積;
	for( s : S){
		s.Set_rotate( {0,0}, 0.3); // 旋轉 避免出現`@LINK:@LOC_1`這種情況;
	}
	F(x): S的並集與直線`X=x`的交集的長度;
	SimpsonIntegration.Ptr_func = F;
	Double ANS = 0;
	for( [l,r] : 函数F的各个定义域){ // 即F在`[l,r]`上有定义;
		ANS += SimpsonIntegration.Get_integration( preLef, preMaxRig);
	}	
}

例題

模板(圓的面積並): @LINK: https://editor.csdn.net/md/?not_checkout=1&articleId=132989973;