找回密码
 欢迎注册
楼主: 笨笨

[求助] 用Mathematica编程求出最大误差函数值时如:10⁻⁸下的(x,a,b,c,d)

  [复制链接]
 楼主| 发表于 2023-9-2 20:45:54 | 显示全部楼层
Jack315 发表于 2023-9-2 10:29
【分享拉马努金 (Ramanujan) 公式推导过程】

首先将超几何函数展开:

QQ截图20230902204445.jpg
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-9-2 20:51:33 | 显示全部楼层
本帖最后由 笨笨 于 2023-9-2 20:54 编辑
Jack315 发表于 2023-9-2 10:47
【补偿函数的构造】
画出误差函数的曲线:

所拟合函数的定义域是0<x<1,而正割函数在\(\left( 0,\frac{\pi}{2}\right)\)图像确实相仿,但它们定义域不同,还要满足我楼上刚发图片上的要求,想想怎么转化之间关系。

点评

将正割函数进行水平和 / 或垂直方向的伸缩应该就行了。  发表于 2023-9-2 21:58
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2023-9-2 21:56:29 | 显示全部楼层

除了拉马努金公式,应该还有其它的函数形式可以达到类似的效果,
只要是包含4个参数,就可以抵消级数展开式中的4个项……
题目中拉马努金公式是必须的吗?

【设想】
除了用一个函数进行补偿外,或许还可能先将拟合函数旋转一个角度,
然后再进行伸缩,使得拟合函数与2F1函数在定义域两端完全重合。
这样误差就可能会减少一(大)部分。这时再用一个函数去补偿或许更容易一些。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-9-2 22:29:44 | 显示全部楼层
本帖最后由 笨笨 于 2023-9-2 23:48 编辑
Jack315 发表于 2023-9-2 21:56
除了拉马努金公式,应该还有其它的函数形式可以达到类似的效果,
只要是包含4个参数,就可以抵消级数展开 ...


按照你的思路走,可否构造一个函数族来且得到什么样的公式?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2023-9-3 10:03:32 | 显示全部楼层
本帖最后由 Jack315 于 2023-9-3 12:12 编辑
笨笨 发表于 2023-9-2 22:29
按照你的思路走,可否构造一个函数族来且得到什么样的公式?


按照设想先选出单纯由初等基本函数构成的拟合函数:
\[\begin{align*}&1+\frac{3x2}{10+\sqrt{4-3x^2}}\tag{1}\\&1+\frac{1}{4}x^2\tag{2}\\&\frac{64-3x^4}{64-16x^2}\tag{3}\\&-1+2e^{\frac{1}{8}x^2}\tag{4}\\&\frac{1}{3}(1+2\cosh(\frac{\sqrt{3}}{2}x))\tag{5}\\&1+\frac{8}{3}\tan^2(\frac{1}{4}\sqrt{\frac{3}{2}}x)\tag{6}\\&\frac{1}{3}\bigg[-7+10\sec\bigg(\frac{1}{2}\sqrt{\frac{3}{5}}x\bigg)\bigg]\tag{7}\\&1+\frac{4}{3}\arcsin^2(\frac{\sqrt{3}}{4}x)\tag{8}\end{align*}\]
(1) 式为拉马努金公式作为参考。
到目前为止构造出的拟合函数,都是单调上升的,
因此最大误差都发生在 \(x=1\) 处。
比较结果拉马努金公式的误差最小,形式也比较优雅。

经过旋转拉伸后,公式表达形式可能会比较复杂。
而由三角函数构成的拟合函数在这一点上或许有一定的优势。
反正 (1) 式经这样处理后,出来的公式特别复杂。

接下来先看看误差能减少多少。
如果再需要补偿的时候,误差函数都是两头为零、中间拱起的形状。
这种情形找补偿函数应该会容易一些。

附件是构造拟合函数的代码,供参考:
超几何函数拟合函数.rar (27.87 KB, 下载次数: 2)

各拟合函数的性能:
拟合函数性能.png

拟合函数的上升速度都比 2F1 函数慢。
如果有一个超过的话,一个线性组全就能使两头端点重合了。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-9-3 20:19:40 | 显示全部楼层
本帖最后由 笨笨 于 2023-9-3 20:35 编辑
Jack315 发表于 2023-9-3 10:03
按照设想先选出单纯由初等基本函数构成的拟合函数:
\[\begin{align*}&1+\frac{3x2}{10+\sqrt{4-3x^2}}\t ...


可以看出拉马努金的公式:\(1+\frac{3 x^2}{10+\sqrt{4-3 x^2}}\)与椭圆周长原级数公式重合度较高,因此只需对拉马努金公式再拟合即可。

所以现在的目标就是找到更好的\(\varphi \left( x \right)\)使得误差越小越好.但还要满足我楼上所提出的要求.

\({}_2{F_1}\left( { - \frac{1}{2}, - \frac{1}{2},1,{x^2}} \right) = 1 + \frac{{3{x^2}}}{{10 + \sqrt {4 - 3{x^2}} }} + \varphi \left( x \right)\)

下面是楼主找到的:
\(Er\left( x \right) \approx {}_2{F_1}\left( { - \frac{1}{2}, - \frac{1}{2},1,{x^2}} \right) - \left( {1 + \frac{{3{x^2}}}{{10 + \sqrt {4 - 3{x^2}} }}} \right) - \frac{3}{{{2^{17}}}}{x^{10}}{\left[ {\left( {\frac{4}{\pi } - \frac{{14}}{{11}}} \right)\frac{{{2^{17}}}}{3}} \right]^{\frac{{{x^2}}}{{{{\left( {1 + 425.524154{x^{ - 2.4565895}}{{\left( {1 - {x^{0.04701004}}} \right)}^{0.94801807}}} \right)}^{0.1537177}}}}}}\)

函数族界最大时的最小值:\(1.383057979\cdots\times10^{-7}\)

你能找到更好的\(\varphi \left( x \right)\)吗?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2023-9-5 08:26:32 | 显示全部楼层
本帖最后由 Jack315 于 2023-9-5 14:14 编辑

问题主要在 x=1 这个奇点及其邻域。
虽然不是剧烈的振荡,但函数的各阶导数却是快速上升。
将这个 2F1 函数的各阶导数画出来就可以看出
要用简单地初等函数获得高精度的拟合有多困难。
下面给个例子。

拟合区间为[0.999, 1.000]:
  1. rng = {0.999, 1.000}
复制代码

拟合函数为二次函数:
  1. f[x_, {a_, b_, c_}] := a + b x + c x^2
复制代码

用于 Table 命令中的计算 x 值的辅助转换函数:
  1. t[n_, {xmin_, xmax_}] := n/2 (xmax - xmin) + xmin
复制代码

利用三个点上两函数值相等来求出 (a,b,c) :
  1. sol = N[Solve[Table[f[t[i, rng], {a, b, c}] ==
  2.      Hypergeometric2F1[-(1/2), -(1/2), 1, (t[i, rng])^2], {i, 0, 2}],
  3.    {a, b, c}], 100]
复制代码

最后画出误差曲线。
  1. Plot[Hypergeometric2F1[-(1/2), -(1/2), 1, x^2] - f[x, {a, b, c}] /. sol, {x, rng[[1]], rng[[2]]},
  2. WorkingPrecision -> 100,
  3. GridLines -> Automatic,
  4. Frame -> True]
复制代码

可以看出精度大致在 1E-8 到 1E-9。

所以理论上,只要拟合区间足够小,
就能达到任意想要的精度,而且是【有限】项的和。
也就是分段拟合能获得任意高的精度。
在实际上,则有可能会受计算精度所限。

而回过头来看楼主找的这个补偿函数。
能达到 1E-7 的水平已经非常优秀。
因为从拟合的角度而言,一个拟合函数要比分段的多个拟合函数好,
只要能达到所需的精度要求就行。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2023-9-5 10:21:16 | 显示全部楼层
本帖最后由 笨笨 于 2023-9-5 10:27 编辑
Jack315 发表于 2023-9-5 08:26
问题主要在 x=1 这个奇点及其邻域。
虽然不是剧烈的振荡,但函数的各阶导数却是快速上升。
将这个 2F1 函数 ...


能否无需分段且对拉马努金公式再拟合达到10^(-8)或10^(-12)?
意思就是说能不能找到比楼主更好的拉马努金公式再拟合函数族?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2023-9-5 10:36:30 | 显示全部楼层
本帖最后由 Jack315 于 2023-9-5 10:39 编辑
笨笨 发表于 2023-9-5 10:21
能否无需分段且对拉马努金公式再拟合达到10^(-8)或10^(-12)?
意思就是说能不能找到比楼主更好的拉马 ...


LZ 构造的补偿函数应该继续发扬光大……真想一睹 LZ 灵感的光亮。
以这个思路下去到 1E-8 以下的精度非常有可能。
精度要求再高,理论上不是没有可能,
但恐怕会遇到计算上的问题,包括:
1. 计算精度
2. 计算时长
让人想起来人家那种挖矿,挖呀挖呀挖……LZ 这个也能挖到矿吗?

点评

谢谢  发表于 2023-9-5 11:44
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2023-9-5 13:48:20 | 显示全部楼层
本帖最后由 Jack315 于 2023-9-5 14:09 编辑

一个函数展开成幂级数时为:
\[\sum_{n=0}^{\infty}c_n(x-x_0)^n\]
在构造补偿函数 \(\varphi(x)\) 时,\(x_0\) 应尽可能选择接近 1 。
这样构造出来的补偿函数,能使 \(x=x_0\) 处的拟合误差最小
(\(c_0=2F1(x_0)\) 时为零),而 \(x=x_0\) 的领域内,误差也会比较小。
由于 \(x=1\) 处是奇点,\(x_0\) 无法选择 1 ,否则这是最佳选择。
基于这一点,LZ 的补偿函数有更新吗?

如果选择 \(x_0=1\) ,函数 2F1 的级数展开式中的系数出现复数。
不是太懂这是什么意思,有人能帮解释一下吗?
这是说函数 2F1 的各阶导数在  \(x=1\) 处依然存在,但不再是实数?
还是说软件在照章办事?

点评

补偿函数暂时没有更新,貌似遇到瓶颈了,很遗憾!  发表于 2023-9-5 23:54
如果按照这种思路走,展开后确实遇见虚数,明显此法不通.  发表于 2023-9-5 23:52
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-10-16 06:13 , Processed in 0.030197 second(s), 18 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表