24小时热门版块排行榜    

查看: 1671  |  回复: 6

科学啄木鸟

铁虫 (小有名气)

[求助] 用matlab写得一个程序,运行不了,不知如何改 已有1人参与

clc;clear
syms M1 M2 M3 %1Cl 2-K 3-Cd
z=[1,1,2];
load dataCl_K_Cd_25;
B0_12=0.04661;B1_12=0.22341;C_phi_12=-0.00044;
B0_13=0.01624;B1_13=0.43945;C_phi_13=-0.00109;
I=1/2*(M1.*z(1)^2+M2.*z(2)^2+M3.*z(3)^2);
Z=z(1)*M1+z(2)*M2+z(3)*M3;
A_phi=0.3920000000;
g=@(x)('2*(1-(1+x).*exp(-x))./x.^2');
g_der=@(x)('-2*(1-(1+x+x.^2/2).*exp(-x))./(x.^2)');
G=g(2*I.^0.5);
G_der=g_der(2*I.^0.5);
x_23=6*z(2)*z(3)*A_phi.*(I.^0.5);x_22=6*z(2)*z(2)*A_phi.*(I.^0.5);x_33=6*z(3)*z(3)*A_phi.*(I.^0.5);
J=@(x)('x.*(4+4.581*x.^(-0.7237).*exp(-0.0120*x.^(0.528))).^(-1)');
J_der=@(x)('(4+4.581*x.^(-0.7237).*exp(-0.0120*x.^(0.528))).^(-1)+(4+4.581*x.^(-0.7237).*exp(-0.0120*x.^(0.528))).^(-2).*(4.581*x.*exp(-0.0120*x.^(0.528)).*(0.7237*x.^(-0.7237-1)+0.0120*0.528*x.^(0.528-1).*x.^(-0.7237)))');
theta_23=(z(2)*z(3)./(4.*I)).*(J(x_23)-J(x_22)/2-J(x_23)/2);
theta_der_23=-(theta_23./I)+(z(2)*z(3)./(8*I.^2)).*(x_23.*J_der(x_23)-x_22.*J_der(x_22)/2-x_33.*J_der(x_33)/2);
syms theta23 psi
Phi_phi_23=theta23+theta_23+I.*theta_der_23;
Phi_23=theta23+theta_23;
Phi_der_23=theta_der_23;
B_phi_12=B0_12+B1_12*exp(-2.0*I.^0.5);
B_phi_13=B0_13+B1_13*exp(-2.0*I.^0.5);
B_12=B0_12+B1_12.*G;
B_13=B0_13+B1_13.*G;
B_der_12=B1_12*G_der./I;
B_der_13=B1_13*G_der./I;
C_12=C_phi_12./(2*(z(1)*z(2))^0.5);C_13=C_phi_13./(2*(z(1)*z(3))^0.5);
F=-A_phi*((I.^0.5)./(1+1.2*I.^0.5)+2*log(1+1.2*I.^0.5)/1.2)+M1.*M2.*B_der_12+M1.*M3.*B_der_13+M2.*M3.*Phi_der_23;
lnR1=z(2)^2*F+M1.*(2*B_12+Z*C_12)+M3.*(2*Phi_23+M1.*psi)+z(2)*(M1.*M2.*C_12+M1.*M3.*C_13);
lnR2=z(1)^2*F+M2.*(2*B_12+Z*C_12)+M3.*(2*B_13+Z*C_13)+z(1)*(M1.*M2.*C_12+M1.*M3.*C_13)+M2.*M3.*psi;
f1=1/2*(lnR1+lnR2);
x1=diff(f1,'theta23');
x2=diff(f1,'psi');
x0=1/2*(z(1)^2*F+M2.*(2*B_12+Z.*C_12)+M3.*(2*B_13+Z.*C_13)+z(1)*(M1.*M2.*C_12+M1.*M3.*C_13)+z(2)^2*F+M1.*(2*B_12+Z.*C_12)+2*M3.*theta_23+z(2)*(M1.*M2.*C_12+M1.*M3.*C_13));

x1=eval(x1);
x2=eval(x2);
X=[ones(5,1),x1,x2];
[b,bint,r,rint,stats]=regress(Y,X)




错误使用  .*
矩阵维度必须一致。

出错 Cl_Na_Sr_250000 (line 17)
theta_23=(z(2)*z(3)./(4.*I)).*(J(x_23)-J(x_22)/2-J(x_23)/2);



有没有哪位matlab大神帮看一下,问题出在哪,如何改
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

koobeexia

至尊木虫 (文坛精英)

【答案】应助回帖

感谢参与,应助指数 +1
需要转置后再运行,试一下

发自小木虫Android客户端
2楼2018-03-06 20:42:22
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

科学啄木鸟

铁虫 (小有名气)

引用回帖:
2楼: Originally posted by koobeexia at 2018-03-06 20:42:22
需要转置后再运行,试一下

哪里需要转置下
3楼2018-03-06 20:49:48
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

独孤神宇

版主 (知名作家)


月只蓝: 金币+1, 鼓励交流 2018-04-04 07:44:20
太过混乱了,又是符号函数,又是匿名函数,参数传递写法不对吧。可以考虑用嵌套函数进行参数传递。P .S. 你这是做参数回归,可以考虑用1stOpt

发自小木虫Android客户端
数值计算
4楼2018-03-06 20:57:29
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

科学啄木鸟

铁虫 (小有名气)

引用回帖:
4楼: Originally posted by 独孤神宇 at 2018-03-06 20:57:29
太过混乱了,又是符号函数,又是匿名函数,参数传递写法不对吧。可以考虑用嵌套函数进行参数传递。P .S. 你这是做参数回归,可以考虑用1stOpt

这是以前师兄他们编写的,以前他们用能拟合出来,但我在matlab 里运行,就报错,捉摸了挺久,没觉得哪里有问题,报错说矩阵维数不一致,我检查了挺多遍,觉得没错啊
5楼2018-03-06 21:04:16
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

rominglv

铜虫 (小有名气)

做拟合不知道用其他工具吗?如minitab或者jmp

发自小木虫IOS客户端
6楼2018-04-04 07:08:23
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

cherish02you

新虫 (初入文坛)


其实可以直接把问题贴出来,说不定有其他的解决办法,你这个参数太多了,看着眼花缭乱

发自小木虫Android客户端
7楼2018-04-16 21:39:59
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 科学啄木鸟 的主题更新
信息提示
请填处理意见