超越方程解的精度越来越高

我有一个特定的运动学作为一个更复杂的机器的一部分，需要计算一些非常困难 （更不可能）的物理参数，以适当的准确度来测量我所用的仪器

[运动学]

首先看它是一个简单的1自由度arm（黑色），可以围绕x轴旋转。 它有一个重量，迫使它总是上升，直到它碰到机械终点（angular度a0 ）或半径r0一些pipe子（蓝色）。 arm的旋转中心位于y0 。 pipe可以移动到任何y(t)高度。

[用法]

这是用来测量pipe的半径进一步处理。 可以计算半径（通过基本测angular法），这导致在图像底部的等式。 常数a0,y0,z0很难测量（它是在复杂的机械中），所以距离的测量精度最小0.1 mm ，angular度0.1 deg ，甚至是有问题的。

[校准]

所以我决定尝试从机器本身完成的一系列测量（自动校准）中计算这些参数。 所以我有已知半径r0校准pipe。 所有的绿色参数都可以作为常量处理。 现在我把pipe子沿着y轴放置，以尽可能多地覆盖arm的angular度。 可悲的是，范围只有20 degrees （对于当前的机器设置）记住测量预设的y(t) …作为n点数据集的a(t) 。 这给了我n超越方程的系统。 从这里我试着猜测a0,y0,z0记得最好的解决scheme（最接近r0 ）的“所有”可能性，

[近似于a0，y0，z0]

近似值是基于这个类的：

 //--------------------------------------------------------------------------- class approx { public: double a,aa,a0,a1,da,*e,e0; int i,n; bool done,stop; approx() { a=0.0; aa=0.0; a0=0.0; a1=1.0; da=0.1; e=NULL; e0=NULL; i=0; n=5; done=true; } approx(approx& a) { *this=a; } ~approx() {} approx* operator = (const approx *a) { *this=*a; return this; } //approx* operator = (const approx &a) { ...copy... return this; } void init(double _a0,double _a1,double _da,int _n,double *_e) { if (_a0<=_a1) { a0=_a0; a1=_a1; } else { a0=_a1; a1=_a0; } da=fabs(_da); n =_n ; e =_e ; e0=-1.0; i=0; a=a0; aa=a0; done=false; stop=false; } void step() { if ((e0<0.0)||(e0>*e)) { e0=*e; aa=a; } // better solution if (stop) // increase accuracy { i++; if (i>=n) { done=true; a=aa; return; } // final solution a0=aa-fabs(da); a1=aa+fabs(da); a=a0; da*=0.1; a0+=da; a1-=da; stop=false; } else{ a+=da; if (a>a1) { a=a1; stop=true; } // next point } } }; //--------------------------------------------------------------------------- 

它通过一些初始步骤search单variables的全部范围，然后find最小偏差点。 之后改变范围，并逐步closures这一点的区域，并recursion地提高精度。

解决scheme本身如下所示：

 // (global) input data #define _irc_calib_n 100 #define _irc_approx_n 5 int irc_calib_ix; // number of measured points double irc_calib_y[_irc_calib_n]; // y(t) double irc_calib_a[_irc_calib_n]; // a(t) double irc_calib_r; // calibration tube radius + arm radius // approximation int ix=0; double e,a,deg=M_PI/180.0; approx aa,ay,az; // min max step recursions ErrorOfSolutionVariable for (aa.init(-90.0*deg,+90.0*deg,10.0*deg,_irc_approx_n,&e);!aa.done;aa.step()) for (ay.init( 0.0 ,200.0 ,10.0 ,_irc_approx_n,&e);!ay.done;ay.step()) for (az.init( 50.0 ,400.0 ,10.0 ,_irc_approx_n,&e);!az.done;az.step()) { for (e=0.0,ix=0;ix<_irc_calib_n;ix++) // test all measured points (e is cumulative error) { a=irc_calib_a[ix]+aa.a; if (a> pi) a-=pi2; if (a<-pi) a+=pi2; if (fabs(a)>0.5*pi) { e=100.0; break; } // ignore too far angles e+=fabs(+(cos(a)*(irc_calib_y[ix]-ay.a)) -(sin(a)*(az.a)) -(irc_calib_r)); } } // here aa.a,ay.a,az.a holds the result 

这导致解决scheme接近测量值，但在模拟内部结果仍不够准确。 根据点数和angular度范围，它从0.1毫米到0.5毫米。 如果我正确地测量z0并且忽略它的近似值，那么精确度被显着地提高，使得y0没有误差（在模拟中）和a0在0.3度左右

Q1如何进一步提高解决scheme的准确性？

我不能增加angular度范围。 点的数量最好在100左右，精度越高越好，但是在150以上的结果是不稳定的（对于某些半径，它是完全closures的）。 绝对不知道为什么。 6以上的recursion数并没有太大的影响

根据从0 degreeangular度距离可以帮助加权偏差？ 但令人遗憾的a(t)范围不一定包括0 degrees

对于y0,z0所需的精确度是0.01 mm y0,z0对于a0 0.01 degree

Q2有没有我错过的东西？

像错误的嵌套近似或一些math简化或不同的方法

[笔记]

angular度必须以a(t)+a0forms表示，因为它是通过SW复位（ 16000 steps/round ）的IRC测量的。 它在a0位置时被复位，我不计算振动和校准pipe的偏心率，他们已经被处理了，我的第一个目标是在没有它们的情况下进行仿真。 pipey(t)可以自由定位， a(t)测量可以随意完成。

现在校准过程沿y轴扫描点（从a0向下移动）。 计算6recursion需要大约35秒（所以耐心）。 5recursion需要大约22秒

[edit1]这里是如何完成模拟的

 approx aa; double e; for (aa.init(-90.0*deg,+90.0*deg,10.0*deg,6,&e);!aa.done;aa.step()) e=fabs(+(cos(aa.a)*(y(t)-y0)) -(sin(aa.a)*(z0)) -(irc_calib_r)); if (aa.a<a0) aa.a=a0; 

一些价值

刚刚意识到我只有4recursion模拟代码来匹配inputIRC精度，那么必须有6recursion。 改变它之后（也在以前的编辑）这里有一些结果

  | a0[deg]| y0[mm] | z0[mm] | simulated | -7.4510|191.2590|225.9000| z0 known | -7.4441|191.1433|225.9000| z0 unknown | -7.6340|191.8074|225.4971| 

所以测量z0的准确度几乎在所需的范围内，但z0未知，误差仍然是需要的大~10倍。 提高模拟精度对6recursion没有影响，也是没有意义的，因为实际input数据也不会更精确。

这里用以上仿真设置进行testing的模拟/测量点：

  ix a [deg] y [mm] 0 -0.2475 +105.7231 1 -0.4500 +104.9231 2 -0.6525 +104.1231 3 -0.8550 +103.3231 4 -1.0575 +102.5231 5 -1.2600 +101.7231 6 -1.4625 +100.9231 7 -1.6650 +100.1231 8 -1.8675 +99.3231 9 -2.0700 +98.5231 10 -2.2725 +97.7231 11 -2.4750 +96.9231 12 -2.6775 +96.1231 13 -2.8575 +95.3077 14 -3.0600 +94.5154 15 -3.2625 +93.7231 16 -3.4650 +92.9308 17 -3.6675 +92.1385 18 -3.8700 +91.3462 19 -4.0725 +90.5538 20 -4.2750 +89.7615 21 -4.4877 +88.9692 22 -4.6575 +88.1769 23 -4.8825 +87.3615 24 -5.0850 +86.5154 25 -5.2650 +85.7000 26 -5.4675 +84.9077 27 -5.6700 +84.1154 28 -5.8725 +83.3231 29 -6.0750 +82.5308 30 -6.2775 +81.7000 31 -6.5025 +80.8462 32 -6.6825 +80.0462 33 -6.8850 +79.2538 34 -7.0875 +78.4615 35 -7.2900 +77.6538 36 -7.5159 +76.7692 37 -7.6725 +75.9769 38 -7.8750 +75.1846 39 -8.1049 +74.3692 40 -8.2800 +73.5000 41 -8.4825 +72.7077 42 -8.6850 +71.9154 43 -8.9100 +71.0308 44 -9.0900 +70.2231 45 -9.2925 +69.4308 46 -9.5175 +68.5462 47 -9.6975 +67.7462 48 -9.9000 +66.9462 49 -10.1025 +66.0615 50 -10.3148 +65.2692 51 -10.4850 +64.3769 52 -10.6875 +63.5846 53 -10.9125 +62.7462 54 -11.0925 +61.9077 55 -11.2950 +61.0846 56 -11.4975 +60.2231 57 -11.7000 +59.3923 58 -11.9025 +58.5308 59 -12.1288 +57.6692 60 -12.3075 +56.8385 61 -12.5100 +55.9462 62 -12.7125 +55.1538 63 -12.9150 +54.2615 64 -13.1175 +53.4000 65 -13.2975 +52.5769 66 -13.5000 +51.6846 67 -13.7025 +50.7923 68 -13.9050 +50.0000 69 -14.1075 +49.1077 70 -14.3100 +48.2154 71 -14.5350 +47.3615 72 -14.7150 +46.5308 73 -14.9175 +45.6385 74 -15.1200 +44.7462 75 -15.3225 +43.8538 76 -15.5250 +42.9615 77 -15.7490 +42.0692 78 -15.9075 +41.2769 79 -16.1100 +40.3846 80 -16.3125 +39.4923 81 -16.5150 +38.6000 82 -16.7175 +37.7077 83 -16.9200 +36.8154 84 -17.1225 +35.9231 85 -17.3250 +34.9308 86 -17.5275 +34.0385 87 -17.7300 +33.1462 88 -17.9325 +32.2538 89 -18.1350 +31.3615 90 -18.3405 +30.4692 91 -18.5175 +29.4769 92 -18.7200 +28.5846 93 -18.9225 +27.6923 94 -19.1250 +26.8000 95 -19.3275 +25.8077 96 -19.5300 +24.9154 97 -19.7325 +23.9231 98 -19.9350 +23.0308 99 -20.1375 +22.1385 

进度更新

对@Ben做了一些澄清

怎么运行的

在第一个图像下面的有色方程给出半径r0它由2个连接的90 degree三angular形（基本三angular）组成，

红色的东西：

• y(t)是电机的位置，它是已知的
• a(t)也是IRC状态

绿色的东西：

• a0,y0,z0是机械尺寸，已知但不精确，所以我用已知的校准piper0测量y(t)不同位置的许多a(t) ，并以更高的精度计算a0,y0,z0

进一步提高准确度

实际上，通过测量y1=y0+z0*cos(a0) ，我们设法使其更精确，精确度大约0.03 mm ，更好。 它是a0位置臂与y轴运动轴交点的高度。 它是从臂上下第一次接触时的情况开始测量和插补的，但是实际位置必须用已用半径和a0 …重新计算，因为接触点不在该轴上…（除非r0=0.0 ）。 由于y1,a0,z0是相关的并且可以相互计算y1,a0,z0所以这也消除了校准的一个近似回路。 由于不连续的测量方式和a(t),y(t)位置a(t),y(t)从消除IRC的测量中去除双重混叠也有助于提高精度和计算稳定性（在真机上）。 我现在无法可靠地评估精确度，因为通过分析许多测量周期，我发现机器上存在一些机械问题，所以我等待修复。 无论如何， r0=80.03 mm的校准与模拟精确度（使用会计两种方法）和_irc_calib_n=30现在是：

  ; computed simulated |delta| a0= -6.915840 ; -6.916710 +0.000870 deg y0=+186.009765 ;+186.012822 +0.003057 mm y1=+158.342452 ;+158.342187 +0.000264 mm z0=+228.102470 ;+228.100000 +0.002470 mm 

校准越大，准确度越低（由于a(t)范围受到更多的限制），这是通过计算所有的a0,y0,(y1),z1没有直接测量或已知的。 这已经是可以接受的，但是正如我以前写的，需要在机器准备好时检查。 在这里完成模拟测量如何看起来像现在：

请参阅近似search如何工作

Answer 1

如果我正确地理解了这一点，你试图从y和a的测量结果推断（但不 测量 ）pipe的半径r0。

将通常的错误传播应用到r0的公式中，就可以得到（估计）所得到的r0的误差。 在小angular度（在此适用，由于a（t）被限制为20度）的极限，这大致给出了（使用三angular函数的小angular度近似）

dr0 ^ 2 = dy ^ 2 + z0 ^ 2（pi * da / 180）^ 2

因此，在r0远小于z0的情况下，r0上的相对误差总是远大于y和z0 * sin（a）的相对误差。 这从图中已经很清楚了：测量的数量只依赖于r0。

换句话说，这不是确定半径r0的巧妙方法。 对于这个基本限制，你可以做的事情不多（除了你可以增加angular度范围a）。 做许多测量（通常的方法来减less噪音/错误）大概不会有帮助，因为这些测量不是由于你的机器的内部工作而相互独立的。 所以，唯一的帮助就是更精确的测量。

为了分析这种情况，我build议把例如推断的r0作为y或y函数的graphics/graphics作为固定r0的函数。