最小二乘法求公式中变量最佳值

hanihaoma

现有一物体 长L，宽b,高h
有公式 v(h,b,L，f) 每组不同的h,b,L 可以求出一组不同的v,f为频率。
有v的真实值，和计算值。
现在想通过最小二乘法，求出 最接近真实值的v的计算值所对应的 h, b ,L。
其实就是对这个物体几何大小的优化，因为 物体的长宽高，会影响到v 的计算值，现在想求一组（h,b,L）使得v的计算值，和真实值最接近。
v 的真实值已知。给定一组长宽高，能算出一组v 的计算值。
下面是部分程序 ：
syms h b L
F=1; ro=2330; E=170e9; a1=1.0553; a2=3.7997; b1=3.8018; b2=2.7364;
f1=linspace(5,500,1000); f2=linspace(100,10000,10000);
f_bas=[f1 f1 f1 f1]'; f_haut=[f2 f2 f2 f2]'; omega_bas=2*pi*f_bas;
o_bas=ones(size(omega_bas)); omega_haut=2*pi*f_haut;
o_haut=ones(size(omega_haut)); [n1,nobr_basF]=size(f1); [n2,nobr_hautF]=size(f2);
visco1=0.001*ones(nobr_basF,1); visco2=0.01*ones(nobr_basF,1); visco3=0.05*ones(nobr_basF,1);
visco4=0.1*ones(nobr_basF,1); visco5=1*ones(nobr_basF,1); visco6=5*ones(nobr_basF,1);
visco7=10*ones(nobr_basF,1); visco8=0.001*ones(nobr_hautF,1); visco9=0.01*ones(nobr_hautF,1);
visco10=0.05*ones(nobr_hautF,1); visco11=0.1*ones(nobr_hautF,1); visco12=1*ones(nobr_hautF,1);
visco13=5*ones(nobr_hautF,1); visco14=10*ones(nobr_hautF,1);
b_visco_basF=[visco1;visco2;visco3;visco4]; h_visco_basF=[visco4;visco5;visco6;visco7];
b_visco_hautF=[visco8;visco9;visco10;visco11]; h_visco_hautF=[visco11;visco12;visco13;visco14];

ro_air=1.29; visco_air=1.8e-5;
%                h=113e-6;      
%                b=600e-6;     
%                L=4000e-6;   
               f_vide=1.875^2*h/2/pi/L^2*sqrt(E/12/ro);  
               II=b*h^3/12; m=ro*b*h*L; ml=m/L; k0 = 3*E*II/L^3;
               gamma=a1*o_bas+a2/b*sqrt(2*b_visco_basF./omega_bas/ro_f)-i*(b1/b*sqrt(2*b_visco_basF./omega_bas/ro_f)+2*b2*b_visco_basF./omega_bas/ro_f/b^2);
               AA=(ml/E/II*omega_bas.^2.*(o_bas+pi*ro_f*gamma*b^2/4/ml)).^(1/4);
               W_L=-F*(sin(AA*L).*exp(-AA*L)+cos(AA*L).*exp(-AA*L)-cos(AA*L).*exp(AA*L)+sin(AA*L).*exp(AA*L))./(AA.^3*E*II)./(2*o_bas+exp(-AA*L).*cos(AA*L)+exp(AA*L).*cos(AA*L));
         
               H0=-L^3/3/E/II; H=W_L./F;
               G=20*log(abs(H)); HH=abs(H)/abs(H0); phi=angle(H)-pi; % faire la phase-pi
               g1=-2*pi*ml*abs(H0)*f_vide^2*sin(phi)./(abs(H).*f);
               g2=ml*((o_bas-abs(H0)./abs(H).*cos(phi)).*(f_vide*o_bas./f).^2-o_bas);

               c1=2*a1*b2-a2*b1; c2=a1*b2-a2*b1;
               d=g1*a2-sqrt((2*pi*f.*g2*b1-a2*g1(n,:)).^2+8*pi*f.*g2.*g1*a1*b2);
               ro_f_calcul=(2*pi*f.*g2*c1+a2*d)./(pi^2*b^2*f*a1*c2);
                v=(a1*(2*pi*f.*g2*b1+d).^2)./(pi*c2*(2*pi*f.*g2*c1+a2*d));

程序是截取的，由于比较长，我也没有试，截取的这段能不能运行，估计不能。
给出程序，主要是想让看到的朋友们更清楚我的问题。
v的计算值，是由几个方程，推出来的。不是直接就和长宽高建立的关系。

[ 本帖最后由 ChaChing 于 2009-6-22 23:17 编辑 ]