c++ - Precisión creciente de la solución de la ecuación trascendental

math geometry

Question

Tengo una cinemática específica como parte de una máquina más compleja y necesito calcular algunos parámetros físicos que son muy difíciles (más como imposibles) de medir con la precisión adecuada con los instrumentos que tengo a mi disposición.

[cinemática]

A primera vista, es un brazo simple de 1 grado de libertad (negro) que puede girar alrededor del eje x . Tiene un peso para obligarlo a subir siempre hasta que llegue al punto final mecánico (ángulo a0 ) o algún tubo (azul) con radio r0 . El centro de rotación del brazo está en y0 . El tubo se puede mover a cualquier altura y(t) .

[uso]

Esto se utiliza para medir el radio de un tubo para su posterior procesamiento. El radio se puede calcular (por goniometría básica) que conduce a la ecuación en la parte inferior de la imagen. Las constantes a0,y0,z0 son muy difíciles de medir (está dentro de maquinaria compleja), por lo que la precisión de medición para distancias es mín. 0.1 mm y ángulo 0.1 deg e incluso eso es cuestionable.

[calibración]

Así que decidí intentar calcular estos parámetros a partir del conjunto de mediciones realizadas por la propia máquina (autocalibración). Entonces tengo un tubo de calibración con radio conocido r0 . Todos los parámetros verdes se pueden manejar como constantes. Ahora coloco el tubo a lo largo del eje y para cubrir tantos ángulos de brazo como pueda. Lamentablemente, el rango es de solo unos 20 degrees (para la configuración actual de la máquina) recordando la medida a(t) para el preajuste y(t) ... como conjunto de datos de n puntos. Esto me da un sistema de n ecuaciones trascendentes. A partir de esto, intento / adivino "todas" las posibilidades de a0,y0,z0 recordando la mejor solución (más cercana a r0 )

[aproximación de a0, y0, z0]

la aproximación se basa en esta clase de mina:

//---------------------------------------------------------------------------
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
            }
        }
    };
//---------------------------------------------------------------------------

Busca el rango completo de una sola variable por algún paso inicial y luego encuentra el punto de desviación mínima. Después de eso, cambie el rango y el paso para cerrar el área de este punto y aumente recursivamente la precisión.

La solución en sí se ve así:

// (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

Esto lleva a una solución cercana a los valores medidos, pero dentro de la simulación el resultado aún no es lo suficientemente preciso. Es de 0.1 mm a 0.5 mm dependiendo del número de puntos y el rango del ángulo. Si z0 correctamente z0 e ignoro su aproximación, entonces la precisión aumenta significativamente dejando y0 sin error (en simulación) y a0 con error alrededor de 0.3 grados

Q1, ¿cómo puedo mejorar aún más la precisión de la solución?

No puedo aumentar el rango angular. El número de puntos es mejor alrededor de 100 mayor es la mejor precisión, pero por encima de 150 el resultado es inestable (para algunos radios está completamente apagado). No tengo ni idea de por qué. El número de recurrencias por encima de 6 no tiene mucho efecto.

¿Podría ayudar a ponderar las desviaciones según la distancia angular de 0 degree ? Pero lamentablemente a(t) rango a(t) no necesariamente incluye 0 degrees

la precisión deseada es 0.01 mm para y0,z0 y 0.01 degree para a0

Q2 ¿hay algo que me haya perdido?

Como aproximaciones erróneamente anidadas o alguna simplificación matemática o enfoque diferente

[notas]

El ángulo debe estar en forma de a(t)+a0 porque es medido por IRC con reinicio de SW ( 16000 steps/round ). Se restablece cuando en una posición 0 no cuento las vibraciones y la excentricidad del tubo de calibración que ya se han solucionado y mi primer objetivo es hacer que esto funcione en simulación sin ellas. El tubo y(t) se puede colocar a voluntad y la medición de a a(t) se puede hacer a voluntad.

En este momento, el proceso de calibración escanea los puntos a lo largo del eje y (movimiento desde a0 hacia abajo). El cálculo con 6 recursiones toma alrededor de 35 segundos (así que sea paciente). 5 recursiones toman alrededor de 22 segundos

[edit1] aquí cómo se realiza la simulación

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;

[edit2] algunos valores

Me acabo de dar cuenta de que solo tenía 4 recursiones en el código de simulación para que coincida con la precisión del IRC de entrada, entonces debe haber 6 recursiones. Después de cambiarlo (también en la edición anterior) aquí hay algunos resultados

                | 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|

Por lo tanto, la precisión con z0 medida está casi en el rango deseado, pero con z0 desconocido, el error sigue siendo ~10 veces mayor de lo necesario. El aumento de la precisión de la simulación no tiene ningún efecto por encima de 6 recursiones y tampoco tiene sentido porque los datos de entrada reales tampoco serán más precisos.

Aquí los puntos simulados / medidos para probar con las configuraciones simuladas anteriores:

 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 

[edit3] actualización de progreso

algunas aclaraciones para @Ben

cómo funciona

la ecuación coloreada debajo de la primera imagen le da el radio r0 está formado por 2 triángulos unidos de 90 degree (trigonometría básica)

cosas rojas:

• y(t) es la posición del motor y se conoce
• a(t) es el estado IRC también conocido

Cosa verde:

• a0,y0,z0 son dimensiones mecánicas y son conocidas pero no precisas, por lo que mido muchas a(t) para diferentes posiciones de y(t) con el tubo de calibración conocido r0 y calculo a0,y0,z0 con mayor precisión a partir de él

mejora de precisión adicional

Realmente logré hacerlo más preciso midiendo y1=y0+z0*cos(a0) partir del movimiento de calibración especial con una precisión de alrededor de 0.03 mm y mejor. Es la altura de la intersección entre el brazo en posición a0 y eje de movimiento del tubo y . Se mide e interpola desde la situación en la que el brazo se pone en contacto por primera vez cuando el tubo sale de arriba hacia abajo, pero la posición real debe calcularse nuevamente con el radio utilizado y a0 ... porque el punto de contacto no está en este eje ... r0=0.0 ). Esto también elimina un bucle de aproximación de la calibración porque y1,a0,z0 son dependientes y pueden calcularse entre sí. Además, eliminar el doble alias de la medición de IRC debido a la forma discontinua de medición y las posiciones a a(t),y(t) ayudaron mucho a aumentar la precisión y la estabilidad de cálculo (en una máquina real). No puedo evaluar con precisión la precisión en este momento porque al analizar muchos ciclos medidos encontré algunos problemas mecánicos en la máquina, así que espero hasta que se repare. De todos modos, la precisión de calibración vs. simulación para r0=80.03 mm con contabilidad de ambos enfoques y _irc_calib_n=30 es ahora:

    ;      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

Cuanto mayor sea la calibración r0 menor r0 la precisión (debido a a(t) rango a a(t) más limitado) al calcular todo a0,y0,(y1),z1 nada se mide directamente ni se conoce. Esto ya es aceptable, pero como escribí antes, es necesario verificar la máquina cuando esté lista. Solo para completar aquí es cómo se ven las mediciones simuladas ahora:

[edit4] vea Cómo funciona la búsqueda de aproximación

Answer 1

Si entiendo esto correctamente, está tratando de inferir (pero no medir ) el radio r0 del tubo a partir de las mediciones para y y a.

Aplicando la propagación de error habitual a su fórmula para r0, se obtiene (una estimación de) el error de la r0 resultante. En el límite de ángulos pequeños (aplicable aquí, ya que a (t) está limitado a 20 grados), esto da aproximadamente (usando la aproximación de ángulo pequeño para las funciones trigonométricas)

dr0 ^ 2 ~ = dy ^ 2 + z0 ^ 2 (pi * da / 180) ^ 2

Por lo tanto, en el caso de r0 mucho más pequeño que z0, el error relativo en r0 siempre es mucho mayor que los errores relativos de y y z0 * sin (a). Esto ya está claro en su gráfico: las cantidades medidas dependen solo débilmente de r0.

En otras palabras, esta no es una forma inteligente de determinar el radio r0. No hay mucho que pueda hacer sobre esta limitación fundamental (excepto que puede aumentar el rango del ángulo a). Presumiblemente, hacer muchas mediciones (el método habitual para vencer el ruido / errores) probablemente no ayudará, porque estas mediciones no son independientes entre sí debido al funcionamiento interno de su máquina. Entonces, la única ayuda serían mediciones más precisas.

Para analizar la situación, recomiendo hacer gráficos / figuras de, digamos, el r0 inferido como función de y o de y como función de a para r0 fijo.