2011-11-07 6 views
8

Sto eseguendo il porting del codice creato in ottava in pylab. Una delle equazioni portate fornisce risultati drammaticamente differenti in Python di quanto non faccia in ottava.Stessa equazione, risposte diverse da Pylab e Octave

Il modo migliore per spiegare è mostrare i grafici generati da ottava e pylab dalla stessa equazione.

Ecco uno snippet semplificato dell'equazione originale in ottava. In questo piccolo script di test, il risultato di una funzione con phi tenuti a zero è tracciata da ~ (-pi, pi):

clear 
clc 
close all 

L1 = 4.25; % left servo arm length 
L2 = 5.75; % left linkage length 
L3 = 5.75; % right linkage length 
L4 = 4.25; % right servo arm length 
L5 = 11/2; % distance from origin to left servo 
L6 = 11/2; % distance from origin to right servo 

theta_array = [-pi+0.1:0.01:pi-0.1]; 
phi = 0/180*pi; 

for i = 1 : length(theta_array) 

theta = theta_array(i); 

A(i) = -L3*(-((2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1)-2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1))/(2*L3*sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2))-((2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1)-2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1))*(-(L6+L5-cos(phi)*L4-cos(theta)*L1)^2-(sin(phi)*L4-sin(theta)*L1)^2-L3^2+L2^2))/(4*L3*((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)^(3/2)))/sqrt(1-(-(L6+L5-cos(phi)*L4-cos(theta)*L1)^2-(sin(phi)*L4-sin(theta)*L1)^2-L3^2+L2^2)^2/(4*L3^2*((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)))-((cos(theta)*L1)/sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)-((sin(theta)*L1-sin(phi)*L4)*(2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1)-2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1)))/(2*((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)^(3/2)))/sqrt(1-(sin(theta)*L1-sin(phi)*L4)^2/((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)))*sin(acos((-(L6+L5-cos(phi)*L4-cos(theta)*L1)^2-(sin(phi)*L4-sin(theta)*L1)^2-L3^2+L2^2)/(2*L3*sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2)))-asin((sin(theta)*L1-sin(phi)*L4)/sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)^2+(sin(phi)*L4-sin(theta)*L1)^2))); 

end 

plot(theta_array,A) 

La trama ottava risultante è simile al seguente:

Octave result

La stessa equazione è stata copiata e incollata da ottava in python con '^' sostituito con '**', 'acos' sostituito con 'arccos' e 'asin' sostituito con 'arcsin'. La stessa gamma di teta è stata tracciata con la phi tenuti a zero:

from pylab import * 

# physical setup 
L1 = 4.25; # left servo arm length 
L2 = 5.75; # left linkage length 
L3 = 5.75; # right linkage length 
L4 = 4.25; # right servo arm length 
L5 = 11.0/2.0; # distance from origin to left servo 
L6 = 11.0/2.0; # distance from origin to right servo 

theta = arange(-pi+0.1,pi-0.1,0.01); 
phi = 0/180.0*pi 

def func(theta,phi): 

A = -L3*(-((2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1)-2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1))/(2*L3*sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2))-((2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1)-2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1))*(-(L6+L5-cos(phi)*L4-cos(theta)*L1)**2-(sin(phi)*L4-sin(theta)*L1)**2-L3**2+L2**2))/(4*L3*((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2)**(3/2)))/sqrt(1-(-(L6+L5-cos(phi)*L4-cos(theta)*L1)**2-(sin(phi)*L4-sin(theta)*L1)**2-L3**2+L2**2)**2/(4*L3**2*((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2)))-((cos(theta)*L1)/sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin((phi)*L4-sin(theta)*L1)**2)-((sin(theta)*L1-sin(phi)*L4)*(2*sin(theta)*L1*(L6+L5-cos(phi)*L4-cos(theta)*L1)-2*cos(theta)*L1*(sin(phi)*L4-sin(theta)*L1)))/(2*((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2)**(3/2)))/sqrt(1-(sin(theta)*L1-sin(phi)*L4)**2/((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2)))*sin(arccos((-(L6+L5-cos(phi)*L4-cos(theta)*L1)**2-(sin(phi)*L4-sin(theta)*L1)**2-L3**2+L2**2)/(2*L3*sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2)))-arcsin((sin(theta)*L1-sin(phi)*L4)/sqrt((L6+L5-cos(phi)*L4-cos(theta)*L1)**2+(sin(phi)*L4-sin(theta)*L1)**2))) 

return A 

f = figure(); 
a = f.add_subplot(111); 

a.plot(theta,func(theta,phi)) 

ginput(1, timeout=-1); # wait for user to click so we dont lose the plot 

risultato di Python è simile al seguente: Python result

Non posso determinare la causa delle differenze, tutte le idee?

+1

Queste funzioni sono le versioni _simplified della funzione originale? Wow. Qualche possibilità che tu possa staccare pezzi identici da entrambi i pezzi uno alla volta e cercare di trovare qualcosa di più piccolo ancora? :) – sarnold

+0

Data la complessità della funzione, potrebbe trattarsi di errori di precisione e/o arrotondamento in virgola mobile diversi? Hai provato a tracciare parti più piccole della funzione per restringere la causa? –

+0

È stato semplificato nel senso che è stato eliminato tutto il codice estraneo per semplificare il problema per i Gurus di overflow dello stack. –

risposta

12

Prova from __future__ import division per eliminare gli errori derivanti dalla divisione del pavimento.

+0

Huzzah! Grazie! sembra averlo risolto Ci sono altri trucchi matematici che dovrei cercare? –

+0

@Inverse_Jacobian: se questa risposta risolve il problema, è necessario accettarlo (fare clic sul segno di spunta da esso). –