作者shyfang (水深火热)
看板NCTU-STAT94G
标题解非线性方程式
时间Fri Aug 14 12:38:00 2009
%in Matlab
% 想要解出x1,x2的方程式, a1,a2,a3 和y1,y2是data
% a1*x1 + a2*x1*x2 = y1
% a3*x2 + a2*x1*x2 = y2
%
% i.e 使f1=0 f2=0 的(x1,x2)
% f1 = a1*x1 + a2*x1*x2 - y1
% f2 = a3*x2 + a2*x1*x2 - y2
%
%
function f = fun(x,a,y)
f=zeros(2,1);
f(1)=a(1)*x(1)+a(2)*x(1)*x(2)-y(1);
f(2)=a(3)*x(2)+a(2)*x(1)*x(2)-y(2);
end
%通常是把function那一段另外存一个*.m档，放在同一个目录下
y=[ 2.1,0.5];
a=[200199000000,58546439234,1];
%解出0根
h=fsolve(@(x) fun(x,a,y),[0;0]);
# in R
# 虽然用Matlab解很方便，但是matlab 2007非常占资源
# 在R里面要用optim (找最小值的function)
# function改成找||f||的。使得norm最小就行了
# 想解的只有x，是x的函式
fun=function(x){
f=rep(0,length(x));
f[1]=a[1]*x[1]+a[2]*x[1]*x[2]-y[1]
f[2]=a[3]*x[2]+a[2]*x[1]*x[2]-y[2]
normf=sum(f^2) # f 的 norm^2
return(normf)
}
y=c( 2.1,0.5);
a=c(200199000000,58546439234,1);
# 把 y&a 的值 assign 好了之後就可以使用
opt<-optim(rep(0,nunknown),fun)
# 因为function不能丢其它的变数，要小心资料的污染。
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 140.114.36.142
1F：推 catherinejoy:有个 stanford 作的 cvx 可在matlab 上使用 也很不赖 11/04 23:19