Jack315 发表于 2023-8-26 11:25:27

为了表示对大师的敬意,也为了尊重 LZ 出的题,
还是希望能得到一个答案。先来个临时结果:


对这个临时结果做个简单说明,细节部分看文档。
在 28# 提了个思路:
【解题思路】
假设给一组 (a,b,c,d) 的值,利用 (EQ:4) 定义的
Err 函数就能求出 x 在其定义域 范围内的最大值。
将这个过程包装到一个新的函数中:
这个新的函数在文档中被称为模型的”评价函数“。
在数学(优化)的术语中,称之为响应函数。

在文档中虽然误差函数都一样,但使用了两个不同的响应函数。
一个是求误差函数的绝对值的最大值;
另一个是求误差函数平方的积分。由于积分运算比较费时,
这个响应函数又改成了离散化后的误差平方和的平均值。
两者精度接近,但计算速度提高很多。

第一个响应函数运行了很长时间,CPU 也开始冒烟了,
但仍然还没有结果出来。原因应该是响应曲面不光滑,
可能会是比较粗糙:参看 20# 的图。
NMinimize 命令是求全局最优解的,因而导致了计算量剧增。

第二个响应函数瞬间就出结果了:(a->0, b->0, c->0, d->0) 。
有一定的补偿效果,看起来至少不是不靠谱了。
由于在 NMinimize 命令的设置里,只使用了参数大于零的空间 (1/16) 。
而搜索结果指向这部分空间的边界,说明这部分空间里没有极值,
或者说最优解至少应该有一个参数必须小于零。

补偿函数的计算只是在 x=0 这个点上有问题,后续代码再作修改。

笨笨 发表于 2023-8-26 13:57:23

本帖最后由 笨笨 于 2023-8-26 14:28 编辑

Jack315 发表于 2023-8-26 11:25
为了表示对大师的敬意,也为了尊重 LZ 出的题,
还是希望能得到一个答案。先来个临时结果:



本楼作为补充说明用:

昨晚整理的,这是四参数拟合的效果图看看有没有可靠的信息。



另外补充一下,一个小细节转换:

Jack315 发表于 2023-8-26 16:50:01

本帖最后由 Jack315 于 2023-8-26 17:36 编辑

能分享下这个公式的由来吗?
\[\begin{matrix}H(x,a,b,c,d)=\frac{3x^{10}}{2^{17}}(1+\frac{[\frac{79}{48}+(\mu-1-\frac{79}{48})x]x^2}{^d}),&\mu=\frac{2^{17}}{3}(\frac{4}{\pi}-\frac{14}{11})\end{matrix}\tag{EQ:3}\]

Jack315 发表于 2023-8-27 08:08:49

本帖最后由 Jack315 于 2023-8-27 08:11 编辑

这个题有可能会走进死胡同了……

先用下列代码定义相关的函数:
ClearAll["Global`*"]
x0 = x -> 0.997;
param = {a -> 0.53, b -> 2.75, c -> 0.82, d -> 3.38};
target := Hypergeometric2F1[-(1/2), -(1/2), 1, x^2]
ramanujan := 1 + (3 x^2)/(10 + Sqrt)
compensation := ((3 x^10)/2^17
(1 + ((79/48 + (\ - 1 - 79/48) x) x^2) /
(1 + x^a (1 - x^b)^c)^d)) /. \ -> 2^17/3 (4/\ - 14/11)
err := target - ramanujan - compensation
然后逐一画出与各参数对应的误差函数的图来:
Plot3D[{err /. param, 0},
{x, 0, 1}, {varA, 0, 2 a /. param},
WorkingPrecision -> 30,
AxesLabel -> Automatic,
PlotLegends -> {"误差", "零平面"}]

Plot3D[{err /. param, 0},
{x, 0, 1}, {varB, 0, 2 b /. param},
WorkingPrecision -> 30,
AxesLabel -> Automatic,
PlotLegends -> {"误差", "零平面"}]

Plot3D[{err /. param, 0},
{x, 0, 1}, {varC, 0, 2 c /. param},
WorkingPrecision -> 30,
AxesLabel -> Automatic,
PlotLegends -> {"误差", "零平面"}]

Plot3D[{err /. param, 0},
{x, 0, 1}, {varD, 0, 2 d /. param},
WorkingPrecision -> 30,
AxesLabel -> Automatic,
PlotLegends -> {"误差", "零平面"}]

这些图说明:
1. 误差函数的过零点(解)就在 param 这个点附近;
2. 在过零点附近,误差函数的曲面如陡崖峭壁——斜率非常大。

这意味着四个参数要精确定位到零点,其精度要求非常高。
是非常非常非常的高。而导致这种情况的原因
很可能就是补偿函数的形式导致了误差函数不是单调的。
或许应该重新寻找合适的拟合函数形式。
包括大师提供的这个拟合函数。
当然也可以继续保留——下一个贴说明原因。

Jack315 发表于 2023-8-27 08:41:26

本帖最后由 Jack315 于 2023-8-27 08:56 编辑

【寻找拟合函数的迭代过程】
任意给定一函数 \(\begin{matrix} y=F_0(x),&x\in\end{matrix}\),
将这个函数看作是与目标函数 \(y=0\) 的误差函数。
目标是寻找一组基本初等函数 \(\begin{matrix}y_i=F_i(x),&i=1,2,...,n\end{matrix}\) 的组合,
以使叠加上这组函数后,误差函数与横坐标之间的面积最小。

将这些基本初等函数加上横、竖两个方向进行伸缩和移位的参数,
并与这些函数某些特定的组合函数作为备选函数,按下列步骤进行迭代:
1. 从备选函数中选出与误差函数形状最接近的一个函数。
2. 优化这个函数的参数。
3. 叠加上这个函数后形成新的误差函数。
4. 判断拟合精度是否满足要求,
不满足要求就重复步骤 1~3 的迭代过程。

这是一个设想,直觉上应该可行。
因为基本初等函数大部分都有比较好的特性,
只要每次叠加上选出的函数后,
误差函数与坐标轴之间的面积都能(大幅)减少就行。

所以实际过程是将每个备选函数的参数先进行优化,
与误差函数之间所围面积最小的就是被选中的函数。

在选择备选函数及其相应参数的选择上,
应避免导致计算困难的形式,以便进行参数的优化。
大师提供的拟合函数也符合这些要求。

笨笨 发表于 2023-8-27 09:03:47

本帖最后由 笨笨 于 2023-8-27 09:08 编辑

Jack315 发表于 2023-8-27 08:41
【寻找拟合函数的迭代过程】
任意给定一函数 \(\begin{matrix} y=F_0(x),&x\in\end{matrix}\),
将 ...

用简单的初等函数构造拟合函数貌似行不通?

不过楼上的elim老师曾经说过:

现在可设法用“初等超越函数”\(\varphi\)拟合楼上的 \(Er\), 使 \(\|E_r-\varphi\|_{\infty}\le 10^{-12}\)

作为阶段性的有理论合实际意义的结果.

附:论坛处于维护中,还是无法上传图片。

Jack315 发表于 2023-8-27 09:36:17

【拟合函数相关定义】

备选初等函数清单,不包括其组合:
1) 常数函数:\(f(x)=c\)
2) 幂函数:\(f(x)=x^a\)
3) 指数函数:\(f(x)=a^x\)
4) 对数函数:\(f(x)=\log_a(x)\)
5) 三角函数:\(f(x)=\sin(x),\ \cos(x),\ \tan(x),\ \cot(x),\ \sec(x),\ \csc(x)\)
6) 反三角函数:\(f(x)=\sin^{-1}(x),\ \cos^{-1}(x),\ \tan^{-1}(x),\ \cot^{-1}(x),\ \sec^{-1}(x),\ \csc^{-1}(x)\)

进行组合的运算:
(有限次的)加、减、乘、除、乘方、开方。

笨笨 发表于 2023-8-28 11:02:32

本帖最后由 笨笨 于 2023-8-28 11:13 编辑

Jack315 发表于 2023-8-26 16:50
能分享下这个公式的由来吗?
\[\begin{matrix}H(x,a,b,c,d)=\frac{3x^{10}}{2^{17}}(1+\frac{[\frac{79}{48 ...
原构造思想来源于楼上“elim”老师的贴子,楼主只是作小幅度改动。不过大致构造思想如下:



另外楼主这几天重新构造一个,效果不错。

\(Er\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}} }} - \frac{3}{{{2^{17}}}}{x^{10}}{\mu ^{\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}\)


Jack315 发表于 2023-8-28 12:21:58

本帖最后由 Jack315 于 2023-8-28 12:24 编辑

笨笨 发表于 2023-8-28 11:02
原构造思想来源于楼上“elim”老师的贴子,楼主只是作小幅度改动。不过大致构造思想如下:




感谢 LZ 分享!

这是补偿函数一步到位的思路。
无论最终用哪一类函数族作为 \(\varphi(x)\),都可能会导致下列问题:
如果误差函数比较粗糙,就是坑坑洼洼的那种,
如 \(f(x,y)=\sin(x^2+y^2)\) 这样的曲面,
软件命令就可能直接无法找到全局最优解,除非运气好。
如果是用光滑函数进行迭代,可能软件命令就能直接搞定。

为了解决找到可能的全局最优解这个问题,或可考虑先进行一个初步的搜索。
期望这个搜索能找到至少在平均意义上拟合效果最佳的区域。
这个就是 elim 在 21# 提出的思想。
现在正在尝试实现这个想法,对一步到位的这类函数族作最后一次努力。

在找到这个最佳区域后,再用软件的命令在这个区域内找出局部最优解。
希望这次能有令人满意的结果。

笨笨 发表于 2023-8-28 12:42:08

本帖最后由 笨笨 于 2023-8-28 12:53 编辑

Jack315 发表于 2023-8-28 12:21
感谢 LZ 分享!

这是补偿函数一步到位的思路。


33楼elim的方法有望突破10^(-12).

static/image/hrline/1.gif

设拉马努金拟合误差函数为\(\psi_0,\;\psi_n\)是函数\(Er_{n-1}:=\psi_0-\cdots-\psi_{n-1}\)
的拟合, 使得 \(\lVert Er_k\rVert /\lVert Er_{k-1}\rVert\le 10^{-2},\:\lVert f\rVert
:=\max|f|((0,1))\)
则对某\(n\le 6\) 就有 \(\lVert Er\rVert \le 10^{-12}\).

以上方法可称为辗转拟合法.

static/image/hrline/1.gif

相信这里的探讨的理论阐述会相当繁琐,或许可以从优化参数序列的(包络)极限来理解?总之很妙很有趣!

现在可设法用初等超越函数\(\varphi_2\)拟合楼上的 \(Er\), 使 \(\|E_r-\varphi_2\|_{\infty}\le 10^{-12}\)
作为阶段性的有理论合实际意义的结果.
页: 1 2 3 4 [5] 6 7 8 9 10 11 12 13
查看完整版本: 用Mathematica编程求出最大误差函数值时如:10⁻⁸下的(x,a,b,c,d)