基于自适应高斯-切比雪夫求积的奇异点附近振荡函数积分的误差可控分区策略
题目描述：
考虑计算有限区间 \([-1, 1]\) 上具有可积奇异点（例如端点奇异性、可去奇异性）且被积函数在奇点附近呈现振荡行为的积分，例如
\[I = \int_{-1}^{1} \frac{\cos(\omega x)}{\sqrt{1-x^2}} \, dx \]
其中 \(\omega\) 是较大的频率参数。被积函数在端点 \(x = \pm 1\) 处具有 \(1/\sqrt{1-x^2}\) 奇异性，同时由于 \(\cos(\omega x)\) 因子，整个被积函数在区间内是振荡的。直接应用标准高斯-切比雪夫求积（其权函数正好是 \(1/\sqrt{1-x^2}\)）虽然能处理奇异性，但对高频振荡可能导致需要大量节点才能达到所需精度。本题目要求设计一种自适应分区策略，在奇异点附近自动加密节点，而在平滑区域使用较少节点，并结合高斯-切比雪夫求积，使得整体误差可控且计算高效。
解题过程循序渐进讲解：
1. 问题分析与核心难点
- 被积函数形式为 \(f(x) = g(x) / \sqrt{1-x^2}\)，其中 \(g(x) = \cos(\omega x)\) 是振荡部分。
- 高斯-切比雪夫（第一类）求积公式的节点和权重正好对应权函数 \(w(x)=1/\sqrt{1-x^2}\)，即
\[ \int_{-1}^{1} \frac{g(x)}{\sqrt{1-x^2}} \, dx \approx \sum_{i=1}^{n} w_i \, g(x_i) \]
其中节点 \(x_i = \cos\left( \frac{(2i-1)\pi}{2n} \right)\)，权重 \(w_i = \frac{\pi}{n}\)。该公式对多项式 \(g(x)\) 能达到最高代数精度。
- 但当 \(g(x)\) 是高振荡函数时，多项式逼近需要很高阶数，导致 \(n\) 必须很大，计算量增加。
- 关键观察：振荡导致函数变化剧烈，但在奇异点附近（靠近 \(x=\pm 1\)），函数值可能因奇异性而增长，同时振荡频率在变换后的变量中可能被压缩或拉伸，使得局部逼近更加困难。因此，需要在奇异点附近采用更密集的采样（即更小的子区间）来捕捉振荡细节。
2. 自适应分区策略设计
自适应积分的基本思想：将整个区间分成若干子区间，在每个子区间上应用高斯-切比雪夫求积（需作变量变换，因为标准公式只适用于 \([-1,1]\) 和权函数 \(1/\sqrt{1-x^2}\)），并根据局部误差估计决定是否进一步细分。
步骤 2.1：区间划分与局部变量变换
- 将 \([-1,1]\) 划分为子区间 \([a,b]\)。对每个子区间，作线性变换 \(x = \frac{b-a}{2}t + \frac{a+b}{2}\)，将积分变为标准形式：
\[ I_{[a,b]} = \int_{a}^{b} \frac{g(x)}{\sqrt{1-x^2}} \, dx = \int_{-1}^{1} \frac{g\left( \frac{b-a}{2}t + \frac{a+b}{2} \right)}{\sqrt{1 - \left( \frac{b-a}{2}t + \frac{a+b}{2} \right)^2}} \cdot \frac{b-a}{2} \, dt \]
- 记新被积函数为 \(G(t) = \frac{b-a}{2} \cdot \frac{g(\cdots)}{\sqrt{1-(\cdots)^2}}\)，此时积分区间为 \([-1,1]\)，但权函数已不是简单的 \(1/\sqrt{1-t^2}\)，因此不能直接应用标准高斯-切比雪夫公式。
- 然而，我们可以将 \(G(t)\) 视为一个整体，在子区间 \([-1,1]\) 上使用普通的高斯-勒让德求积（权函数为1）。但这样会失去对奇异性的专门处理。为了利用高斯-切比雪夫公式处理奇异性，我们换一种方式：将原积分拆分为奇异部分和振荡部分，但更直接的做法是在每个子区间上仍使用高斯-切比雪夫公式，但需调整节点和权重。
- 实际上，对于每个子区间 \([a,b]\)，我们可以构造该区间上带权函数 \(1/\sqrt{1-x^2}\) 的高斯-切比雪夫公式，这需要计算该权函数在 \([a,b]\) 上的正交多项式（复杂）。为了避免这个复杂性，通常采用全局自适应策略：先将整个 \([-1,1]\) 用标准高斯-切比雪夫公式计算，然后根据误差估计决定在哪些参数域（即 \(t\) 域）细分，而不是在物理 \(x\) 域细分。
步骤 2.2：基于参数域的自适应细分
标准高斯-切比雪夫求积的节点是 \(t_i = \cos\theta_i\) 其中 \(\theta_i = (2i-1)\pi/(2n)\)。在参数 \(\theta\) 域（即反余弦变换后的角度域）上，节点是等距分布的。振荡函数 \(g(\cos\theta)\) 在 \(\theta\) 域中可能呈现出非均匀的振荡行为（因为 \(x=\cos\theta\) 是非线性变换）。
- 在 \(\theta\) 域进行自适应细分更为自然，因为权函数 \(1/\sqrt{1-x^2}\) 在 \(\theta\) 域变成常数：当 \(x=\cos\theta\)，有 \(dx = -\sin\theta d\theta\)，则
\[ \int_{-1}^{1} \frac{g(x)}{\sqrt{1-x^2}}dx = \int_{0}^{\pi} g(\cos\theta) \, d\theta \]
这是一个无权重积分，被积函数为 \(g(\cos\theta)\)，区间为 \([0,\pi]\)。
- 因此，原积分转化为在 \(\theta\) 区间 \([0,\pi]\) 上对 \(g(\cos\theta)\) 的积分，此时奇异性消失，但振荡频率可能变化（因为 \(g(\cos\theta)\) 的频率与 \(\omega\) 和 \(\cos\theta\) 有关）。
- 在 \(\theta\) 域应用自适应积分策略（例如自适应辛普森或自适应高斯-克朗罗德）更为直接，因为积分是无权重的。但本题要求使用高斯-切比雪夫求积，我们可以这样做：
- 将 \(\theta\) 区间划分为子区间 \([\theta_{k-1}, \theta_k]\)。
- 在每个子区间上，用高斯-勒让德求积计算 \(\int_{\theta_{k-1}}^{\theta_k} g(\cos\theta) d\theta\)，因为此时被积函数平滑（但可能振荡）。
- 但这样与高斯-切比雪夫公式的原形式脱节。为了保持联系，我们注意到：标准 \(n\) 点高斯-切比雪夫公式等价于在 \(\theta\) 域上用中点公式计算，即
\[ \int_{0}^{\pi} g(\cos\theta) d\theta \approx \frac{\pi}{n} \sum_{i=1}^{n} g\left( \cos\left( \frac{2i-1}{2n}\pi \right) \right) \]
这是矩形法，精度不高。为了得到高精度，需要增加 $n$ 或采用更高阶的复合求积。
- 因此，自适应策略可以设计为：在 \(\theta\) 域上，将 \([0,\pi]\) 划分为若干子区间，在每个子区间上用高斯-勒让德求积计算积分，并利用误差估计控制细分。这样虽然名义上不是“高斯-切比雪夫求积”，但实质是等价且更灵活的。
然而，题目要求基于高斯-切比雪夫求积，所以我们可以将每个子区间 \([\theta_{k-1}, \theta_k]\) 再通过变量变换映射回 \(x\) 域的子区间 \([a,b]\)，然后在 \([a,b]\) 上使用带适当权函数的高斯-切比雪夫型公式。但这样实现复杂。
3. 实用的自适应高斯-切比雪夫求积策略
实际中常用的方法是：在 \(x\) 域上，将积分写成
\[I = \int_{-1}^{1} \frac{g(x)}{\sqrt{1-x^2}} dx = \sum_{k} \int_{a_k}^{b_k} \frac{g(x)}{\sqrt{1-x^2}} dx \]
对每个子区间 \([a_k,b_k]\)，用变换 \(x = \cos\theta\) 将其映射为 \(\theta\) 区间 \([\alpha, \beta]\)，其中 \(\alpha = \arccos b_k, \beta = \arccos a_k\)（注意反余弦单调递减）。于是
\[\int_{a_k}^{b_k} \frac{g(x)}{\sqrt{1-x^2}} dx = \int_{\alpha}^{\beta} g(\cos\theta) d\theta \]
这样，每个子积分转化为无权重积分。然后在 \([\alpha, \beta]\) 上应用自适应高斯-勒让德求积（因为被积函数 \(g(\cos\theta)\) 可能振荡，但无奇异性）。
- 误差估计：在子区间 \([\alpha, \beta]\) 上，计算两个不同精度的积分近似，比如用 \(m\) 点和 \(m+1\) 点高斯-勒让德公式，比较其差值作为误差估计。
- 细分准则：如果误差估计大于给定容差，则将 \([\alpha, \beta]\) 对半分，并递归计算。
- 由于 \(x = \cos\theta\) 变换将 \(x\) 端点附近映射为 \(\theta\) 的较大区间（例如 \(x\) 靠近1时，\(\theta\) 靠近0，但区间长度较大），因此振荡在 \(\theta\) 域可能被拉伸，使得自适应细分在振荡剧烈处自然发生。
4. 算法步骤总结
- 输入：函数 \(g(x)\)，频率 \(\omega\)，误差容差 \(\epsilon\)，最大递归深度 \(max\_depth\)。
- 初始化积分总和 \(I \leftarrow 0\)，将初始 \(\theta\) 区间 \([0,\pi]\) 压入栈。
- 当栈非空时：
a. 弹出当前区间 \([\alpha, \beta]\)。
b. 计算该区间上的积分近似 \(Q_1\)（例如用低阶高斯-勒让德，如5点）和 \(Q_2\)（用高阶，如10点）。
c. 估计误差 \(err = |Q_2 - Q_1|\)。
d. 如果 \(err < \epsilon \times (\beta-\alpha)/\pi\)（局部容差按区间长度比例分配），则接受 \(Q_2\)，累加到 \(I\)。
e. 否则，如果递归深度未超限，将区间对半分为 \([\alpha, \gamma]\) 和 \([\gamma, \beta]\) 其中 \(\gamma=(\alpha+\beta)/2\)，分别压入栈，深度加1。 - 输出 \(I\) 作为积分近似值。
5. 误差控制与振荡处理
- 此方法本质是在 \(\theta\) 域做自适应高斯-勒让德积分，避免了奇异性，且能自动在 \(g(\cos\theta)\) 变化剧烈处（对应 \(x\) 靠近端点或振荡密集处）细分。
- 对于高频 \(\omega\)，\(g(\cos\theta)=\cos(\omega \cos\theta)\) 在 \(\theta\) 域的振荡频率不是常数，但自适应细分能保证局部近似精度。
- 由于在 \(\theta\) 域积分无权重，标准高斯-勒让德求积的高代数精度可有效处理振荡，比直接用复合高斯-切比雪夫公式更灵活且节点数更少。
6. 示例计算
以 \(g(x)=\cos(10x)\) 为例，在 \([-1,1]\) 上积分 \(I=\int_{-1}^{1} \frac{\cos(10x)}{\sqrt{1-x^2}} dx\)。
- 变换为 \(\theta\) 域：\(I=\int_{0}^{\pi} \cos(10\cos\theta) d\theta\)。
- 应用上述自适应策略，在 \(\theta\) 域细分，比如在 \(\theta\) 靠近 \(0\) 和 \(\pi\) 处，\(\cos\theta\) 变化快，可能需要更多子区间。
- 最终得到满足误差要求的近似积分值。
这种方法结合了高斯-切比雪夫求积消除奇异性的思想（通过变换到 \(\theta\) 域），又利用自适应高斯-勒让德求积灵活处理振荡，实现了在奇异点附近自动加密采样，且整体误差可控。