24小时热门版块排行榜    

查看: 282  |  回复: 0

yangjian911x

铁虫 (初入文坛)

[求助] 求大神帮忙解决一个matlab算例问题

下面是我写的一串代码:
k0=0;
for h=0.01:0.01:50;
%////////cengshu
p=101325;%wendu
t=300;
a=100;
N=50;
l=100;
g=h/a;
px=124.5;%xianweimidu
pa=110;%qiningjiaomidu
p1     =         16.8367680681458;
p2     =         351.719716799583;
p3      =        170.039755829837;
p4     =         23.1916198325517;
p5     =         1269.41233622451;
p6      =        -4.76831022571228;
p7      =        -124.095938410379;
p8             =0.158164420679919;
p9      =        3.48022582551722;
kgo=1.6241*10^(-3)+8.4798*10^(-5)*t+2.8587*10^(-9)*t^2-3.7581*10^(-11)*t^3+1.6705*10^(-14)*t^4;
kg=kgo/(1+sqrt(2)*1.75*1.38*10^(-23)*t/(pi*0.366^(2)*10^(-18)*p*l*10^(-9)));
%
a2=0.8;
s=(324.3/pa+5.03)*10^5;%比表面积
kng=60.22*10^(5)*p*t^(-0.5)/(0.25*s*(1-pa/2200)^(-1)*pa+4.01*10^(9)*p*t^(-1));
d=12*(1-(1-pa/2200))/((2+a2^2)*pa*s);
D=(-pi*sqrt(1-a2^2)*d^2/(pa*s)+((pi*sqrt(1-a2^2)*d^2/(pa*s))^2+(-pi*d/pa/s)^3)^0.5)^(1/3)+(-pi*sqrt(1-a2^2)*d^2/(pa*s)-((pi*sqrt(1-a2^2)*d^2/(pa*s))^2+(-pi*d/pa/s)^3)^0.5)^(1/3);
a1=d/D;
C1=5.2*10^(-6);
kSO=0.75264+0.0031286*t-0.0000045242*t^2+3.5253*10^(-9)*t^3;
kS0=0.75264+0.0031286*300-0.0000045242*300^2+3.5253*10^(-9)*300^3;
kns=C1*pa^1.5*kSO/kS0;%气凝胶固体热导率
b=1-kng/kns;
n=50;
ka=(pi*a1^2*a2^2/(4*(1-b))+(1-a1^2)-pi*a1^2*(1-a2^2)*(b+log(1-b))/(2*b^2)+pi*(1/sqrt(1-a2^2)-a1)*(log((1-b*a1*a2)/(1-b*a1))/(b*a1)-(1-a2))/b)*kng;
%//气凝胶热导率
l1=0.508742684953864;
l2=0.0783770550695127;
l3=-0.004081649272934;
l4=0.000177301411935905;
l5=2.59695488944144*10^(-7);
l6=-3.97089287852354*10^(-7);
l7=1.22847546022621*10^(-8);
l8=-7.92885898906499*10^(-11);
q1     =         596118.307345212;
q2     =         -113950.80374852;
q3      =        9006.52284616573;
q4      =        -372.348258004001;
q5      =        8.38526805820327;
q6      =        -0.0955903338053084;
q7      =        0.000414336358839436;%气凝胶
ks=l1+l2*t^0.5+l3*t+l4*t^1.5+l5*t^2+l6*t^2.5+l7*t^3+l8*t^3.5;
kc=0.8*ka+0.2*ks;
kcg=2*kc*kg/(kc+kg);
%复合材料单元体密度
KR1=(p1+p3*t^0.5+p5*t+p7*t^1.5+p9*t^2)/(1+p2*t^0.5+p4*t+p6*t^1.5+p8*t^2);
KR2=(q1+q2*t^0.5+q3*t+q4*t^1.5+q5*t^2+q6*t^2.5+q7*t^3)/pa;
KR=px*KR1+0.8*0.95*pa*KR2;
kr=(16*5.67*10^(-8)*t^3)/(3*KR);%fushe
ke=kg+4*g*(sqrt(3)-g)*(kc-kg)/3;%danceng
if k0<=0;
for i=1:N-1;
A1=rand(1)*((1.5*sqrt(3)*a^2-6*a*h+2*sqrt(3)*h^2)-(1.5*sqrt(3)*a^2-12*a*h+8*sqrt(3)*h^2))+1.5*sqrt(3)*a^2-12*a*h+8*sqrt(3)*h^2;
A2=rand(1)*((6*a*h-2*sqrt(3)*h^2)-(4*sqrt(3)*h^2))+4*sqrt(3)*h^2;
A3=1.5*sqrt(3)*a^2-(A1+A2);
K=(A1*kg+A2*kc+A3*kcg)/(1.5*sqrt(3)*a^2);
k0=k0+1/K;
hold on;
end;
end;
%///////
k=N/(k0+1/ke)+kr;
hold on;
plot(h,k,'b-.*');
end
想求各位大神,能不能保证代码中(红色部分)随机取样对应于不同的h只运行一次。

[ Last edited by jjdg on 2012-11-17 at 23:55 ]
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 yangjian911x 的主题更新
信息提示
请填处理意见