2016-06-30 12 views
0

Ich versuche eine Funktion von Matlab zu Python zu konvertieren. Die Matlab-Funktion ist:Gleiche Gleichung gibt unterschiedliche Werte in Matlab und Numpy?

function [f,df_dr1,df_dr2,g,dg_dr1,dg_dr2] = f_eval_2eq(r1,r2,r3,z1,z2,z3,n1,n2,n3) 

f  = (r1)./sqrt(r1.^2 + z1.^2)... 
- (n2/n1)*(r2-r1)./sqrt((r2-r1).^2 + z2.^2); 

df_dr1 = 1./sqrt(r1.^2 + z1.^2)... 
- r1.^2./(r1.^2 + z1.^2).^(3/2)... 
+ (n2/n1)./sqrt(z2.^2 + (r1-r2).^2)... 
- (n2/n1).*(r1-r2).*(2*r1-2*r2)./(2*((r1-r2).^2 + z2.^2).^(3/2)); 

df_dr2 = (n2/n1).*(r1-r2).*(2*r1-2*r2)./(2*((r1-r2).^2 + z2.^2).^(3/2))... 
- (n2/n1)./sqrt(z2.^2 + (r1-r2).^2); 


g  = (r2-r1)./sqrt((r2-r1).^2 + z2.^2)... 
- (n3/n2)*(r3-r2)./sqrt((r3-r2).^2 + z3.^2); 

dg_dr1 = (r1-r2).*(2*r1-2*r2)./(2*((r1-r2).^2 + z2.^2).^(3/2))... 
- 1./sqrt(z2.^2 + (r1-r2).^2); 

dg_dr2 = 1./sqrt((r1-r2).^2 + z2.^2)... 
+ (n3/n2)./sqrt(z3.^2 + (r2-r3).^2)... 
- (r1-r2).*(2*r1-2*r2)./(2*((r1-r2).^2 + z2.^2).^(3/2))... 
- (n3/n2).*(r2-r3).*(2*r2-2*r3)./(2*((r2-r3).^2 + z3.^2).^(3/2)); 

end 

%test code 
K>> a=[1,2,3];b=a+1;c=b+1;d=a;e=b;f=c;g=1;h=2;i=3; 
K>> [f,df_dr1,df_dr2,g,dg_dr1,dg_dr2] = f_eval_2eq(a,b,c,d,e,f,g,h,i) 

Die Python-Funktion, die ich geschrieben habe, ist:

def f_eval_2eq(r1,r2,r3,z1,z2,z3,n1,n2,n3): 
    #evaluate gradients 
    #n_ are scalars 
    f = (r1)/np.sqrt(r1**2 + z1**2) \ 
     - (n2/n1)*(r2-r1)/np.sqrt((r2-r1)**2 + z2**2); 

    df_dr1 = 1/np.sqrt(r1**2 + z1**2) \ 
      - r1**2/((r1**2 + z1**2)**(3/2)) \ 
      + (n2/n1)/np.sqrt(z2**2 + (r1-r2)**2) \ 
      - (n2/n1)*(r1-r2)*(2*r1-2*r2)/(2*((r1-r2)**2 + z2**2)**(3/2)); 

    df_dr2 = (n2/n1)*(r1-r2)*(2*r1-2*r2)/(2*((r1-r2)**2 + z2**2)**(3/2)) \ 
      - (n2/n1)/np.sqrt(z2**2 + (r1-r2)**2); 


    g  = (r2-r1)/np.sqrt((r2-r1)**2 + z2**2) \ 
      - (n3/n2)*(r3-r2)/np.sqrt((r3-r2)**2 + z3**2); 

    dg_dr1 = (r1-r2)*(2*r1-2*r2)/(2*((r1-r2)**2 + z2**2)**(3/2)) \ 
      - 1/np.sqrt(z2**2 + (r1-r2)**2); 

    dg_dr2 = 1/np.sqrt((r1-r2)**2 + z2**2) \ 
      + (n3/n2)/np.sqrt(z3**2 + (r2-r3)**2) \ 
      - (r1-r2)*(2*r1-2*r2)/(2*((r1-r2)**2 + z2**2)**(3/2)) \ 
      - (n3/n2)*(r2-r3)*(2*r2-2*r3)/(2*((r2-r3)**2 + z3**2)**(3/2)); 

    return (f,df_dr1,df_dr2,g,dg_dr1,dg_dr2) 

#test code 
A=np.array([1,2,3]) 
B=A+1 
C=B+1 
D=A 
E=B 
F=C 
G=1 
H=2 
I=3 
[f,df_dr1,df_dr2,g,dg_dr1,dg_dr2] =f_eval_2eq(A,B,C,D,E,F,G,H,I) 
print ('f= '+str(f) +'\n'+'df_dr1= '+str(df_dr1) +'\n' +'df_dr2='+str(df_dr2) +'\n'+'g= '+str(g) +'\n'+'dg_dr1= '+str(dg_dr1) +'\n'+'dg_dr2= '+str(dg_dr2) +'\n') 

Der Ausgang für f ist die gleiche in beiden, aber alle anderen Werte sind unterschiedlich, und ich kann nicht herausfinden, warum? ??

Jede Hilfe wird geschätzt.

+4

Python 2.x oder 3.x? Sie haben '3/2' anstelle von' 3./2' überall – Suever

Antwort

2

In Python 2.x, wenn Sie teilen zwei ganze Zahlen (wie 2 und 3) das Ergebnis als eine ganze Zahl als auch gegossen wird:

x = 3/2 
# 1 

type(x) 
# <type 'int'> 

Sie müssen explizit entweder den Zähler oder Nenner spezifizieren Sei ein Gleitkommawert anstelle einer Ganzzahl, die einen Dezimalpunkt verwendet, und dies ermöglicht auch, dass die Ausgabe ein Gleitkommawert ist.

y = 3./2 
# 1.5 

type(y) 
# <type 'float'> 

Alternativ kann, wie durch @rayryeng vorgeschlagen, können Sie das folgende an der Spitze des Codes legen Sie das Verhalten, das Sie erwarten, zu bekommen.

from __future__ import division 
+0

, dass ... klingt wie es könnte das Problem sein. Würde man '[r1, r2, r3, z1, z2, z3, n1, n2, n3] = np.float _ ([r1, r2, r3, z1, z2, z3, n1, n2, n3])' fixieren ? – JoshuaF

+0

@ user3769033 Fügen Sie einfach einen expliziten '.' nach all Ihren Ganzzahlen in Ihren Code ein. – Suever

+0

Nur hinzugefügt, dass. Jetzt sind f und g korrekt, aber die anderen Werte sind noch unterschiedlich? – JoshuaF

2

Sie können auch

from __future__ import division 

zum Anfang der Datei hinzufügen, wenn Sie Python verwenden 2, um die Python 3 Verhalten zu erhalten, das heißt immer float Division verwenden.

Verwandte Themen