scipy.optimize.最小二乘法和matlab lsqnonlin差分

Question

我想找出这些脚本中给出的方程的最小值。它看起来非常混乱（但不需要对方程进行深刻的理解——我想）。在def的末尾是要最小化的表达式：

 vys1=-Qd-40*sqrt(5)*sqrt((ch+cm)*ep*kb*Na*T)*w1
 vys2=fi0-fib-Q0/cq 
 vys3=fib-fid+Qd/cq1 
 vysf= np.array([vys1,vys2,vys3])
 return vysf

我在matlab中使用lsqnonlin来比较结果。Matlab的结果似乎非常准确。结果为（fi0、fib、fid）

^{pr2}$

注意，这个脚本会检查公式中的错误（如果它们在python和matlab中是相同的） 对于[fi0,fib,fid]=[-0.120, -0.0750 ,-0.011]，结果是相同的[vys1,vys2,vys3]-

python [0.00069376  0.05500097 -0.06179421]
matlab [0.0006937598,0.05500096 -0.06179421]

在least_squares中是否有改进结果的选项？谢谢你的帮助（抱歉误会英语）

Python

import scipy as sc
import numpy as np
from math import sinh
import matplotlib as plt
from numpy import exp, sqrt
from scipy.optimize import leastsq,least_squares

def q(par,ep,Na,kb,T,e,gamaal,gamasi,gamax,k1,k2,k3,k4,cq,cq1,ch,cm):
    fi0,fib,fid=np.array([par[0],par[1],par[2]])
    AlOH= gamaal*k1*exp(e*fi0/(T*kb))/(ch + k1*exp(e*fi0/(T*kb)))
    AlOH2= ch*gamaal/(ch + k1*exp(e*fi0/(T*kb)))
    SiO= gamasi*k2*exp(e*fi0/(T*kb))/(ch + k2*exp(e*fi0/(T*kb)))
    SiOH= ch*gamasi/(ch + k2*exp(e*fi0/(T*kb)))
    X= gamax*k3*k4*exp(e*fib/(T*kb))/(ch*k4 + cm*k3 + k3*k4*exp(e*fib/  (T*kb)))
    XH= ch*gamax*k4/(ch*k4 + cm*k3 + k3*k4*exp(e*fib/(T*kb)))
    Xm= cm*gamax*k3/(ch*k4 + cm*k3 + k3*k4*exp(e*fib/(T*kb)))
    Q0=e*(0.5*(AlOH2+SiOH-AlOH-SiO)-gamax)
    Qb=e*(XH+Xm)
    Qd=-Q0-Qb
    w1=sc.sinh(0.5*e*fid/kb/T)
    vys1=-Qd-40*sqrt(5)*sqrt((ch+cm)*ep*kb*Na*T)*w1
    vys2=fi0-fib-Q0/cq 
    vys3=fib-fid+Qd/cq1 
    vysf= np.array([vys1,vys2,vys3])
    return vysf

kb=1.38E-23;T=300;e=1.6e-19;Na=6.022e23;gamaal=1e16;gamasi=1e16
gamax=1e18;k1=1e-4;k2=1e5;k3=1e-4;k4=1e-4;cq=1.6;cq1=0.2
cm=1e-3;ep=80*8.8e-12
ch1=np.array([1e-3,1e-5,1e-7,1e-10])

# Check the equations, if they are same
x0=np.array([-0.120,    -0.0750 ,-0.011])
val=q(x0,ep,Na,kb,T,e,gamaal,gamasi,gamax,k1,k2,k3,k4,cq,cq1,ch1[0],cm)
print(val)
w1=least_squares(q,x0, args=(kb,ep,Na,T,e,gamaal,gamasi,gamax,k1,k2,k3,
                             k4,cq,cq1,ch1[0],cm))
print(w1['x'])

matlab语言

function[F1,poten,fval]=test()
kb=1.38E-23;T=300;e=1.6e-19;Na=6.022e23;gamaal=1e16;gamasi=1e16;gamax=1e18;
k1=1e-4;k2=1e5;k3=1e-4;k4=1e-4;cq=1.6;cq1=0.2;ch=[1e-3];cm=1e-3;ep=80*8.8e-   12;
% Test if equation are same 
x0=[-0.120, -0.0750 ,-0.011];
F1=rovnica(x0,ch) ;
[poten,fval]= lsqnonlin(@(c) rovnica(c,ch(1)),x0);
function[F]=rovnica(c,ch) 
fi0=c(1);
fib=c(2);
fid=c(3);
aloh=exp(1).^(e.*fi0.*kb.^(-1).*T.^(-1)).*gamaal.*k1.*(ch+exp(1).^(e.* ...
fi0.*kb.^(-1).*T.^(-1)).*k1).^(-1);
aloh2=ch.*gamaal.*(ch+exp(1).^(e.*fi0.*kb.^(-1).*T.^(-1)).*k1).^(-1);
sioh=ch.*gamasi.*(ch+exp(1).^(e.*fi0.*kb.^(-1).*T.^(-1)).*k2).^(-1);
sio=exp(1).^(e.*fi0.*kb.^(-1).*T.^(-1)).*gamasi.*k2.*(ch+exp(1).^(e.* ...
fi0.*kb.^(-1).*T.^(-1)).*k2).^(-1);
Xm=cm.*gamax.*k3.*(cm.*k3+ch.*k4+exp(1).^(e.*fib.*kb.^(-1).*T.^(-1)) ...
.*k3.*k4).^(-1);
 XH=ch.*gamax.*k4.*(cm.*k3+ch.*k4+exp(1).^(e.*fib.*kb.^(-1).*T.^(-1)) ...
.*k3.*k4).^(-1);
Q0=e*(0.5*(aloh2+sioh-aloh-sio)-gamax);
Qb=e*(XH+Xm);
Qd=-Q0-Qb;
F=[-Qd+(-40).*5.^(1/2).*((ch+cm).*ep.*kb.*Na.*T).^(1/2).*sinh((1/2).*e.* ...
fid.*kb.^(-1).*T.^(-1));...
 fi0-fib-Q0/cq;...
(fib-fid+Qd/cq1)];
 end
end

Answer 1

这一行有个错误：

w1=least_squares(q,x0, args=(kb,ep,Na,T,e,gamaal,gamasi,gamax,k1,k2,k3,
                             k4,cq,cq1,ch1[0],cm))

参数kb在错误的位置。q的签名是：

^{pr2}$

参数kb介于Na和{}之间。如果修复了least_squares调用中的args参数：

w1 = least_squares(q, x0, args=(ep, Na, kb, T, e, gamaal, gamasi, gamax,
                                k1, k2, k3, k4, cq, cq1, ch1[0], cm))

然后Python脚本的输出是

[ 0.00069376  0.05500097 -0.06179421]
[-0.13253313 -0.03253254 -0.02131043]

这与Matlab的输出结果一致。在