«Монография УДК ББК Авторский коллектив Прохоров В.Т., Осина Т.М., Жихарев А.П., Михайлов, А.Б., Михайлова И.Д. Рецензенты: Доктор технических наук, профессор А. Жаворонков (г.Шахты) Доктор ...»
Винили- 2,30 0,0842 0,0951 0,1003 0,1042 0,1131 0,1272 0, скожа-Т + пенополиуретан + трикотаж (ПА волокно) Ткань ка- 2,45 0,0791 0,0963 0,1001 0,1104 0,1234 0,1263 0, проновая + пенополиуретан + трикотаж (ПА волокно) Байка х/б 0,75 0,0304 0,0343 0,0365 0,0392 0,0433 0,0445 0, Кожа хро- 1,60 0,0612 0,0707 0,0742 0,0782 0,082 0,0852 0, мовая для верха обуви Топография теплового потока у мужчин и женщин Часть тела Голова 51,8±2,0 62,9±2,8 100,5±2,9 97,6±1,9 79,0±2,2 105,6±3, (лоб) Туловище 46,5±1,3 47,9±2,2 84,1±1,4 67,0±2,0 32,7±1,2 29,3±0, Плечо 41,5±1,4 41,2±1,7 72,6±2,1 56,1±1,7 40,6±1,9 33,8±1, Кисть 51,7±1,9 57,9±2,9 70,2±3,8 62,6±4,0 80,8±2,7 75,2±8, Бедро 40,8±1,6 40,4±1,8 70,5±1,9 62,2±2,9 39,1±1,7 46,2±2, Голень 49,4±1,6 50,5±2,2 90,2±2,7 68,3±4,9 51,1±2,1 65,1±1, Стопа 37,3±2,7 54,3±3,4 57,6±3,5 66,6±4,2 53,0±1,4 64,0±3, Средневзвешенный теп- 45,2±1,0 47,6±1,7 78,8±1,6 66,8±1,8 45,1±1,6 49,1±0, ловой поток, Вт/м 1 Программа расчета зависимости температуры внутриобувного пространства от времени для многослойного плоского пакета материалов restart:
"Введите число слоев":n1:= d:=array(1..n1):s:=array(1..n1):a:=array(1..n1):
"Введите толщину слоев пакета":d[1]:= :d[2]:= :
"Введите коэффициенты теплопроводности слоев": s[1]:= :s[2]:= :
"Введите коэффициенты температуропроводности слоев": a[1]:= :a[2]:= :
"Введите плотность теплового потока":POT:= :
"Введите температуру окружающей среды в градусах K" : TC:= :
"Введите начальную температуру обуви в градусах K": TP:= :
F:=(z,l1,l2,l3,h1,h2,h3,k1,k2,k3)-(-h1/sqrt(k1) *sin(z*l1/sqrt(k1))*cos(z*(l2-l1)/sqrt(k2))h2/sqrt(k2)*cos(z*l1/sqrt(k1))*sin(z*(l2-l1)/sqrt(k2)))* sin(z*(l3- l2)/sqrt(k3))+h3/sqrt(k3)*cos(z*(l3-l2)/ sqrt(k3))*(cos(z*l1/sqrt(k1))*cos(z*(l2-l1)/sqrt(k2))- h1/sqrt(k1)*sqrt(k2)/h2*sin(z*l1/sqrt(k1))*sin(z*(l2-l1)/sqrt(k2))):
h:=array(1..3):k:=array(1..3):l:=array(0..3):
l[0]:=0:l[1]:=d[1]:l[2]:=sum(d[n2],n2=2..n1-1): l[3]:=d[n1]+l[2]:
h[1]:=s[1]:h[2]:=l[2]/sum(d[n3]/s[n3],n3=2..n1-1): h[3]:=s[n1]:
k[1]:=a[1]:k[2]:=l[2]*h[2]/sum(d[n4]*s[n4]/a[n4], n4=2..n1k[3]:=a[n1]:
z:=array(1..5):
i:=1:x:=0.01:
p:= array(1..15):p[1]:=0.01:
while i13 do if evalf(F(x,l[1],l[2],l[3],h[1],h[2],h[3],k[1],k[2],k[3])*F(x+0.1,l[1],l[2],l[3],h [1],h[2],h[3],k[1],k[2],k[3]))0 then x:=x+0.1 else p[i+1]:=x:
p[i+2]:=x+0.1: i:=i+2: x:=x+0.1: fi; od;
for n to 5 do z[n]:=fsolve(F(z,l[1],l[2],l[3],h[1],h[2],h[3],k[1],k[2],k[3])=0,z=p[2*np[2*n+2]):od:
z[1]:=z[1]:z[2]:=z[2]:z[3]:=z[3]:z[4]:=z[4]:z[5]:=z[5]:
A:=array(1..5,1..3):
for i to 5 do A[i,1]:=0:A[i,2]:=(tan(z[i]*l[1]/sqrt(k[2]))h[1]*sqrt(k[2])/(h[2]*sqrt(k[1]))*tan(z[i]*l[1]/sqrt(k[1])))/(1+h[1]*sqrt(k[ ])/(h[2]*sqrt(k[1]))*tan(z[i]*l[1]/sqrt(k[1]))*tan(z[i]*l[1]/sqrt(k[2]))):A[i,3] :=-1/tan(z[i]*l[3]/sqrt(k[3])) od:
M:=array(1..5,1..3):
M[1,1]:=1:M[2,1]:=1:M[3,1]:=1:M[4,1]:=1:M[5,1]:=1:
for n to 5 do for m from 2 to 3 do M[n,m]:=product((cos(z[n]*l[r-1]/k[r-1]^(1/2))+A[n,r-1]* sin(z[n]*l[rk[r-1]^(1/2)))/(cos(z[n]*l[r-1]/ k[r]^(1/2))+A[n,r]*sin(z[n]*l[r- 1]/k[r]^(1/2))),r=2..m)od od;
H:=array(1..3): C:=array(1..5):
for m to 3 do H[m]:=-POT/h[m] od:
B:=array(1..5):
for s to 3 do B[s]:=POT*(l[3]/h[3]+sum(l[k]*(1/h[k]-1/h[k+1]),k=s..2))od:
g:=array(1..5,1..5):
for s from 1 to 5 do for q from 1 to 5 do g[q,s]:=evalf(sum(h[d]/k[d]*M[q,d]*M[s,d]*int((cos(z[q]*v/sqrt(k[d]))+A[ q,d]*sin(z[q]*v/sqrt(k[d])))*(cos(z[s]*v/sqrt(k[d]))+A[s,d]*sin(z[s]*v/sqrt( k[d]))),v=l[d-1].. l[d]),d=1..3)) od od;
print(g);
for y to 5 do C[y]:=sum(h[p]/k[p]*M[y,p]*int((TP-TC-H[p]*v-B[p])* (cos(z[y]*v/sqrt(k[p]))+ A[y,p]*sin(z[y]*v/sqrt(k[p]))),v=l[pl[p]),p=1..3)/g[y,y] od ;
print(C);
x:=0:
T1:=sum(C[l]*M[l,1]*exp(-z[l]^2*t)*(cos(z[l]*x/ sqrt(k[1]))+A[l,1]*sin(z[l]*x/sqrt(k[1]))),l=1..5)+TC-273+H[1]*x+B[1];
plot([T1,21],t=0.01..8,T=0..50,color=[black,black]);
2 Программа расчета зависимости температуры внутриобувного пространства от времени для многослойного цилиндрического пакета материалов restart:
"Введите число слоев":n1:= :
d:=array(1..n1):s:=array(1..n1):a:=array(1..n1):
"Введите внутренний радиус пакета":r0:= :
"Введите толщину слоев пакета": d[1]:= :d[2]:= :
"Введите коэффициенты теплопроводности слоев":
s[1]:=0.05:s[2]:=0.163:s[3]:=0.04:s[4]:=0.038:s[5]:=0.07:s[6]:=0.06:
"Введите коэффициенты температуропроводности слоев": a[1]:= :a[2]:= :
"Введите плотность теплового потока":POT:= :
"Введите температуру окружающей среды в градусах K":TC:= :
"Введите начальную температуру обуви в градусах K": TP:= :
"Введите коэффициент теплоотдачи":al:= :
F:=(z,l1,l2,l3,h1,h2,h3,k1,k2,k3,R0,R1,R2,R3,al)- h[2]*(BesselJ(1, z*R[2]/k[2]^(1/2))*(h[1]*(BesselJ(1, z*R[1]/k[1]^(1/2))*BesselY(1, z*R[0]/k[1]^(1/2))-BesselJ(1, z*R[0]/k[1]^(1/2))*BesselY(1, z*R[1]/ k[1]^(1/2))) k[2]^(1/2)*BesselY(0, z*R[1]/k[2]^(1/2))k[1]^(1/2)*(BesselJ(0, z*R[1]/k[1]^(1/2))*BesselY(1, z*R[0]/k[1]^(1/2))BesselJ(1, z*R[0]/k[1]^(1/2))* BesselY(0, z*R[1]/k[1]^(1/2)))*h[2]*BesselY(1, z*R[1] /k[2]^(1/2)))+(h[1]*(BesselJ(1, z*R[1]/k[1]^(1/2))* BesselY(1, z*R[0]/k[1]^(1/2))BesselJ(1, z*R[0]/ k[1]^(1/2))*BesselY(1, z*R[1]/k[1]^(1/2))) *k[2]^(1/2)* BesselJ(0, z*R[1]/k[2]^(1/2))+k[1]^(1/2)*(BesselJ(0, z*R[1]/k[1]^(1/2))*BesselY(1, z*R[0]/k[1]^(1/2))-BesselJ(1, z*R[0]/k[1]^(1/2))*BesselY(0, z*R[1]/ k[1]^(1/2)))*h[2]*BesselJ(1, z*R[1]/k[2]^(1/2)))* BesselY(1, z*R[2]/k[2]^(1/2)))*(BesselJ(0, z*R[2]/k[3]^(1/2))*(al*BesselY(0, z*R[3]/k[3]^(1/2))-h[3]*z*BesselY(1, z*R[3]/k[3]^(1/2))/k[3]^(1/2))+ (h[3]*z*BesselJ(1, z*R[3]/k[3]^(1/2))/k[3]^(1/2)-al*BesselJ(0, z*R[3]/ k[3]^(1/2)))*BesselY(0, z*R[2]/k[3]^(1/2)))/k[2]^(1/2)-(h[3]*(BesselJ(1, z*R[2]/k[3]^(1/2))*(al*BesselY(0, z*R[3]/k[3]^(1/2))-h[3]*z*BesselY(1, z*R[3]/k[3]^(1/2))/k[3]^(1/2))+ (h[3]*z*BesselJ(1, z*R[3]/k[3]^(1/2))/k[3]^(1/2)-al*BesselJ(0, z*R[3]/k[3]^(1/2)))*BesselY(1, z*R[2]/ k[3]^(1/2)))*(BesselJ(0, z*R[2]/k[2]^(1/2))* (h[1]* (BesselJ(1, z*R[1]/k[1]^(1/2))*BesselY(1, z*R[0]/ k[1]^(1/2))-BesselJ(1, z*R[0]/k[1]^(1/2))*BesselY(1, z*R[1]/ k[1]^(1/2)))* k[2]^(1/2)*BesselY(0, z*R[1]/ k[2]^(1/2))k[1]^(1/2)*(BesselJ(0, z*R[1]/k[1]^(1/2))* BesselY(1, z*R[0]/k[1]^(1/2))BesselJ(1, z*R[0]/k[1]^(1/2))*BesselY(0, z*R[1]/ k[1]^(1/2)))* h[2]*BesselY(1, z*R[1]/k[2]^(1/2)))+(-h[1]*(BesselJ(1, z*R[1]/ k[1]^(1/2))*BesselY(1, z*R[0]/k[1]^(1/2))-BesselJ(1, z*R[0]/k[1]^(1/2))*BesselY(1, z*R[1]/ k[1]^(1/2)))*k[2]^(1/2)*BesselJ(0, z*R[1]/k[2]^(1/2))+ k[1]^(1/2)* (BesselJ(0, z*R[1]/k[1]^(1/2))*BesselY(1, z*R[0]/k[1]^(1/2))-BesselJ(1, z*R[0]/ k[1]^(1/2))* BesselY(0, z*R[1]/k[1]^(1/2)))*h[2]*BesselJ(1, z*R[1]/ k[2]^(1/2)))*BesselY(0, z*R[2]/k[2]^(1/2)))/k[3]^(1/2)):
R:=array(0..3):R[0]:=r0: R[1]:=R[0]+d[1]:
R[2]:=R[1]+sum(d[n2],n2=2..n1-1):R[3]:=d[n1]+R[2]:
h[1]:=s[1]:h[2]:=(R[2]-R[1])/sum(d[n3]/s[n3],n3=2..n1-1):h[3]:=s[n1]:
k[1]:=a[1]:k[2]:=(R[2]-R[1])*h[2]/ sum(d[n4]*s[n4]/a[n4],n4=2..n1k[3]:=a[n1]:
z:=array(1..5):
i:=1:x:=0.01:
p:= array(1..15):p[1]:=0.01:
while i13 do if evalf(F(x,l[1],l[2],l[3],h[1],h[2],h[3],k[1],k[2],k[3],R[0],R[1],R[2],R[3],al) *F(x+0.1,l[1],l[2],l[3],h[1],h[2],h[3],k[1],k[2],k[3],R[0],R[1],R[2],R[3],al)) 0 then x:=x+0.1 else p[i+1]:=x: p[i+2]:=x+0.1: i:=i+2: x:=x+0.1: fi: od:
for n to 5 do z[n]:=fsolve(F(z,l[1],l[2],l[3],h[1],h[2],h[3],k[1],k[2],k[3],R[0],R[1],R[2],R [3],al)=0,z=p[2*n-1]..p[2*n+2]):od:
z[1]:=z[1]:z[2]:=z[2]:z[3]:=z[3]:z[4]:=z[4]:z[5]:=z[5]:
y:=array(1..5,1..3):
for n to 5 do y[n,1]:=-BesselJ(1,z[n]*R[0]/k[1]^(1/2))/ BesselY(1,z[n]*R[0]/k[1]^(1/2)):
y[n,3]:=(h[3]/k[3]^(1/2)*z[n]*BesselJ(1,z[n]*R[3]/k[3]^(1/2))al*BesselJ(0,z[n]*R[3]/k[3]^(1/2)))/ (al*BesselY(0,z[n]*R[3]/k[3]^(1/2))h[3]/ k[3]^(1/2)*z[n]*BesselY(1,z[n]*R[3]/k[3]^(1/2))):
y[n,2]:= (-h[1]*(BesselJ(1, z[n]*R[1]/k[1]^(1/2))-BesselJ(1, z[n]*R[0]/k[1]^(1/2))*BesselY(1, z[n]*R[1]/ k[1]^(1/2))/BesselY(1, z[n]*R[0]/k[1]^(1/2)))* k[2]^(1/2)*BesselJ(0, z[n]*R[1]/k[2]^(1/2))+ k[1]^(1/2)* (BesselJ(0, z[n]*R[1]/k[1]^(1/2))-BesselJ(1, z[n]*R[0]/ k[1]^(1/2))*BesselY(0, z[n]*R[1]/k[1]^(1/2))/BesselY(1, z[n]*R[0]/k[1]^(1/2)))*h[2]*BesselJ(1, z[n]*R[1]/ k[2]^(1/2)))/(h[1]*(BesselJ(1, z[n]*R[1]/k[1]^(1/2))-BesselJ(1, z[n]*R[0]/k[1]^(1/2))*BesselY(1, z[n]*R[1]/ k[1]^(1/2))/BesselY(1, z[n]*R[0]/k[1]^(1/2)))* k[2]^(1/2)*BesselY(0, z[n]*R[1]/k[2]^(1/2))k[1]^(1/2)* (BesselJ(0, z[n]*R[1]/k[1]^(1/2))-BesselJ(1, z[n]*R[0]/ k[1]^(1/2))*BesselY(0, z[n]*R[1]/k[1]^(1/2))/BesselY(1, z[n]*R[0]/k[1]^(1/2)))*h[2]*BesselY(1, z[n]*R[1]/ k[2]^(1/2)))od:
M:=array(1..5,1..3):
M[1,1]:=1:M[2,1]:=1:M[3,1]:=1:M[4,1]:=1:M[5,1]:=1:
for j to 5 do for m from 2 to 3 do M[j,m]:=product((BesselJ(0,z[j]*R[r-1]/k[r-1]^(1/2))+ y[j,rBesselY(0,z[j]*R[r-1]/k[r-1]^(1/2)))/ (BesselJ(0,z[j]*R[rk[r]^(1/2))+y[j,r]*BesselY(0,z[j]* R[r-1]/k[r]^(1/2))),r=2..m)od od:
g:=array(1..5,1..5):
for q to 5 do for w from 1 to 5 do g[w,q]:=sum(h[b]/k[b]*M[w,b]*M[q,b]*int(v*(BesselJ(0,z[w]*v/sqrt(k[b]) )+y[w,b]*BesselY(0,z[w]*v/sqrt(k[b])))*(BesselJ(0,z[q]*v/sqrt(k[b]))+y[q, b]*BesselY(0,z[q]*v/sqrt(k[b]))),v=R[b-1]..R[b]),b=1..3)od od:
C:=array(1..5):
A:=array(1..3):B:=array(1..3):
for i to 3 do A[i]:=-POT*R[0]/h[i]:
B[i]:=POT*R[0]*(1/(al*R[3])+ln(R[i])/h[i]+sum(ln(R[k]/R[kh[k],k=i+1..3)) od:
r:=R[0]:
for s to 5 do C[s]:=sum(h[p]/k[p]*M[s,p]*int(v*(TP-TC-A[p]*ln(v)- B[p])*(BesselJ(0,z[s]*v/sqrt(k[p]))+y[s,p]*BesselY(0,z[s]*v/sqrt(k[p]))),v =R[p-1]..R[p]),p=1..3)/g[s,s]od:
T1:=sum(C[l]*exp(-z[l]^2*t)*(BesselJ(0,z[l]*r /sqrt(k[1]))+y[l,1]*BesselY(0,z[l]*r/sqrt(k[1]))),l=1..5)+TCA[1]*ln(r)+B[1]:
plot([T1], t=0.01..5,T=-10..40,color=[black]);
3 Программа расчета зависимости температуры внутриобувного пространства от времени для многослойного шарового пакета материалов restart:
"Введите число слоев":n1:= :
d:=array(1..n1):s:=array(1..n1):k:=array(1..n1):
"Введите внутренний радиус пакета":r0:= :
"Введите толщину слоев пакета": d[1]:= :d[2]:= :
"Введите коэффициенты теплопроводности слоев": s[1]:= :s[2]:= :
"Введите коэффициенты температуропроводности слоев": k[1]:= :k[2]:= :
"Введите плотность теплового потока":POT:= :
"Введите температуру окружающей среды в градусах K":TC:= :
"Введите начальную температуру обуви в градусах K":TP:= :
"Введите коэффициент теплоотдачи":al:= :
F:=(z,h1,h2,h3,a1,a2,a3,R0,R1,R2,R3,al)-(-h2*(cos(z*R2/a2^(1/2))* (cos(z*R1/a1^(1/2))*sin(z*R1/a2^(1/2))*(h2*(z*cos(z*R0/a1^(1/2))/a1^(1/ 2)-sin(z*R0/a1^(1/2))/R0)/R1-h1*(z*cos(z*R0/a1^(1/2))/ a1^(1/2)sin(z*R0/a1^(1/2))/R0)/R1+h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+ cos(z*R0/a1^(1/2))/R0)*z/a1^(1/2))-sin(z*R1/a1^(1/2))* sin(z*R1/a2^(1/2))* (h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1- h2*(z*sin(z*R0/ a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1+h1*z*(z*cos(z*R0/a1^(1/ ))/a1^(1/2)-sin(z*R0/a1^(1/2))/R0)/a1^(1/2))-h2*z*cos(z*R1/a2^(1/2))* (cos(z*R1/ a1^(1/2))*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/a1^(1/2)) /R0)+(z* sin(z*R0/ a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*sin(z*R1/a1^(1/2)))/a2^(1/2))+ (cos(z*R1/a1^(1/2))*cos(z*R1/a2^(1/2))*(- h2*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/ a1^(1/2))/R0)/R1+h1*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)/R1- h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*z/a1^(1/2))+si n(z*R1/a1^(1/2))*cos(z*R1/a2^(1/2))*(h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2) +cos(z*R0/a1^(1/2))/R0)/R1- h2*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1+ h1*z*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/a1^(1/2))/R0)/a1^(1/2))- h2*z* sin(z*R1/a2^(1/2))*(cos(z*R1/a1^(1/2))*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z* R0/a1^(1/2))/R0)+(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0) *sin(z*R1/a1^(1/2)))/a2^(1/2))*sin(z*R2/a2^(1/2)))/R2+h2*(- z*sin(z*R2/a2^(1/2))* (cos(z*R1/a1^(1/2))*sin(z*R1/a2^(1/2))*(h2*(z*cos(z*R0/a1^(1/2))/a1^(1/ 2)-sin(z*R0/a1^(1/2))/R0)/R1-h1*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/ a1^(1/2))/R0)/R1+h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/ R0)*z/a1^(1/2))- sin(z*R1/a1^(1/2))*sin(z*R1/a2^(1/2))*(h1*(z*sin(z*R0/a1^(1/2))/ a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1-h2*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+ cos(z*R0/a1^(1/2))/R0)/R1+h1*z*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/ a1^(1/2))/R0)/a1^(1/2))- h2*z*cos(z*R1/a2^(1/2))*(cos(z*R1/a1^(1/2))*(z* cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)+(z*sin(z*R0/a1^(1/2))/ a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*sin(z*R1/a1^(1/2)))/a2^(1/2))/a2^(1/2)+( cos(z*R1/a1^(1/2))*cos(z*R1/a2^(1/2))*(- h2*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/ a1^(1/2))/R0)/R1+h1*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)/R1- h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*z/a1^(1/2))+si n(z*R1/a1^(1/2))*cos(z*R1/a2^(1/2))*(h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2) +cos(z*R0/a1^(1/2))/R0)/R1- h2*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1+ h1*z*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/a1^(1/2))/R0)/a1^(1/2))- h2*z*sin(z*R1/a2^(1/2))*(cos(z*R1/a1^(1/2))*(z*cos(z*R0/a1^(1/2))/a1^( 1/2)- sin(z*R0/a1^(1/2))/R0)+(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/ ))/R0)*sin(z*R1/a1^(1/2)))/a2^(1/2))*z*cos(z*R2/a2^(1/2))/a2^(1/2)))*(co s(z*R2/a3^(1/2))*(- h3*sin(z*R3/a3^(1/2))/R3+h3*z*cos(z*R3/a3^(1/2))/a3^(1/2)+ al*sin(z* R3/a3^(1/2)))+(h3*cos(z*R3/a3^(1/2))/R3+h3*z*sin(z*R3/a3^(1/2))/a3^( /2)- al*cos(z*R3/a3^(1/2)))*sin(z*R2/a3^(1/2)))+(h3*(cos(z*R2/a3^(1/2))*(- h3*sin(z* R3/a3^(1/2))/R3+h3*z*cos(z*R3/a3^(1/2))/a3^(1/2)+al*sin(z*R3/a3^(1/2)) )+(h3*cos(z*R3/a3^(1/2))/R3+h3*z*sin(z*R3/a3^(1/2))/a3^(1/2)- al*cos(z*R3/ a3^(1/2)))*sin(z*R2/a3^(1/2)))/R2-h3*(- z*sin(z*R2/a3^(1/2))*(-h3*sin(z*R /a3^(1/2))/R3+h3*z*cos(z*R3/a3^(1/2))/a3^(1/2)+al*sin(z*R3/a3^(1/2)))/a 3^(1/2)+(h3*cos(z*R3/a3^(1/2))/R3+h3*z*sin(z*R3/a3^(1/2))/a3^(1/2)- al*cos(z*R3/ a3^(1/2)))*z*cos(z*R2/a3^(1/2))/a3^(1/2)))*(cos(z*R2/a2^(1/2))*(cos(z*R 1/a1^(1/2))*sin(z*R1/a2^(1/2))*(h2*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/ R0)/R1-h1*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)/R1+ h1* (z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*z/a1^(1/2))- sin(z*R1/ a1^(1/2))*sin(z*R1/a2^(1/2))*(h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z* R0/a1^(1/2))/R0)/R1- h2*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1+ h1*z*(z*cos(z*R0/a1^(1/2))/a1^(1/2)-sin(z*R0/a1^(1/2))/R0)/a1^(1/2))- h2*z*cos(z*R1/a2^(1/2))*(cos(z*R1/a1^(1/2))*(z*cos(z*R0/a1^(1/2))/a1^( 1/2)- sin(z*R0/a1^(1/2))/R0)+(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/ ))/R0)*sin(z*R1/a1^(1/2)))/a2^(1/2))+(cos(z*R1/a1^(1/2))*cos(z*R1/a2^( /2))*(-h2*(z* cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)/R1+h1*(z* cos(z*R0/ a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)/R1-h1*(z*sin(z*R0/a1^(1/2)) /a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*z/a1^(1/2))+sin(z*R1/a1^(1/2))*cos(z* R1/a2^(1/2))*(h1*(z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/ R1-h2*(z* sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)/R1+h1*z*(z*cos(z*R 0/a1^(1/2))/a1^(1/2)-sin(z*R0/a1^(1/2))/R0)/a1^(1/2))- h2*z*sin(z*R1/a2^(1/2))* (cos(z*R1/a1^(1/2))*(z*cos(z*R0/a1^(1/2))/a1^(1/2)- sin(z*R0/a1^(1/2))/R0)+ (z*sin(z*R0/a1^(1/2))/a1^(1/2)+cos(z*R0/a1^(1/2))/R0)*sin(z*R1/a1^(1/2) ))/a2^(1/2))*sin(z*R2/a2^(1/2))):
R:=array(0..3):h:=array(1..3):a:=array(1..3):
R[0]:=r0: R[1]:=R[0]+d[1]:R[2]:=R[1]+sum(d[n2],n2=2..n1-1):
R[3]:=d[n1]+R[2]:
h[1]:=s[1]:h[2]:=(R[2]-R[1])/sum(d[n3]/s[n3],n3=2..n1-1):h[3]:=s[n1]:
a[1]:=k[1]:a[2]:=(R[2]-R[1])*h[2]/sum(d[n4]*s[n4]/k[n4],n4=2..n1-1):
a[3]:=k[n1]:
z:=array(1..5):
i:=1:x:=0.01:
p:= array(1..15):p[1]:=0.01:
while i13 do if evalf(F(x,h[1],h[2],h[3],a[1],a[2],a[3],R[0],R[1],R[2],R[3],al)*F(x+0.1,h[1],h[2],h[3],a[1],a[2],a[3],R[0],R[1],R[2],R[3],al))0 then x:=x+0.1 else p[i+1]:=x: p[i+2]:=x+0.1: i:=i+2: x:=x+0.1: fi: od:
for n to 5 do z[n]:=fsolve(F(z,h[1],h[2],h[3],a[1],a[2],a[3],R[0],R[1],R[2],R[3],al)=0,z=p [2*n-1]..p[2*n+2]);od:
z[1]:=z[1]:z[2]:=z[2]:z[3]:=z[3]:z[4]:=z[4]:z[5]:=z[5]:
y:=array(1..5,1..3):
for n to 5 do y[n,1]:=(z[n]*sin(z[n]*R[0]/a[1]^(1/2))/a[1]^(1/2)+cos(z[n]*R[0]/a[1]^(1/ ))/R[0])/(z[n]*cos(z[n]*R[0]/a[1]^(1/2))/a[1]^(1/2)sin(z[n]*R[0]/a[1]^(1/2))/R[0]):
y[n,3]:=(h[3]*cos(z[n]*R[3]/a[3]^(1/2))/R[3]+h[3]*z[n]*sin(z[n]*R[3]/a[3] ^(1/2))/a[3]^(1/2)-al*cos(z[n]*R[3]/a[3]^(1/2)))/(- h[3]*sin(z[n]*R[3]/a[3]^(1/2))/R[3]+h[3]*z[n]*cos(z[n]*R[3]/a[3]^(1/2))/a [3]^(1/2)+al*sin(z[n]*R[3]/a[3]^(1/2))):
y[n,2]:=(cos(z[n]*R[1]/a[1]^(1/2))*cos(z[n]*R[1]/a[2]^(1/2))*(-h[2]/R[1]+ h[1]/R[1]-h[1]*y[n,1]*z[n]/a[1]^(1/2))+sin(z[n]*R[1] /a[1]^(1/2))* cos(z[n]* R[1]/a[2]^(1/2))*(h[1]*y[n,1]/R[1]- h[2]*y[n,1]/R[1]+h[1]*z[n]/a[1]^(1/2))- h[2]*z[n]*sin(z[n]*R[1]/a[2]^(1/2))*(cos(z[n]*R[1]/a[1]^(1/2))+y[n,1]*sin (z[n]*R[1]/a[1]^(1/2)))/a[2]^(1/2))/(cos(z[n]*R[1]/a[1]^(1/2))*sin(z[n]*R[ ]/a[2]^(1/2))*(h[2]/R[1]-h[1]/R[1]+h[1]*y[n,1]*z[n]/a[1]^(1/2))- sin(z[n]*R[1]/a[1]^(1/2))* sin(z[n]*R[1]/a[2]^(1/2))*(h[1]*y[n,1]/R[1]- h[2]*y[n,1]/R[1]+h[1]* z[n] /a[1]^(1/2))- h[2]*z[n]*cos(z[n]*R[1]/a[2]^(1/2))*(cos(z[n]*R[1]/a[1]^(1/2))+ y[n,1]*sin(z[n]*R[1]/a[1]^(1/2)))/a[2]^(1/2)):od:
M:=array(1..5,1..3):
M[1,1]:=1:M[2,1]:=1:M[3,1]:=1:M[4,1]:=1:M[5,1]:=1:
for n to 5 do for m from 2 to 3 do M[n,m]:=product((cos(z[n]*R[r-1]/a[r-1]^(1/2))+y[n,r-1]*sin(z[n]*R[ra[r-1]^(1/2)))/(cos(z[n]*R[r-1]/a[r]^(1/2))+y[n,r]*sin(z[n]*R[r- 1]/a[r]^(1/2))),r=2..m) od od:
g:=array(1..5,1..5):
for s from 1 to 5 do for q from 1 to 5 do g[q,s]:=evalf(sum(h[d]/a[d]*M[q,d]*M[s,d]*int((cos(z[q]*v/sqrt(a[d]))+y[q,d]*sin(z[q]*v/sqrt(a[d])))*(cos(z[s]*v/sqrt(a[d]))+y[s,d]*sin(z[s]*v/sqrt(a[ d]))),v=R[d-1]..R[d]),d=1..3)) od od:
A:=array(1..3):
C:=array(1..5):
for m to 3 do A[m]:=POT*R[0]^2/h[m] od:
B:=array(1..5):
for s to 3 do B[s]:=POT*R[0]^2*(1/(al*R[3]^2)-1/(R[3]*h[3])+sum(1/R[3-k]*(1/h[3k+1]-1/h[3-k]),k=1..3-s))od:
for l to 5 do C[l]:=sum(h[p]/a[p]*M[l,p]*int((TP-TC-A[p]/v-B[p])*v*(cos(z[l]*v/ sqrt(a[p]))+y[l,p]*sin(z[l]*v/sqrt(a[p]))),v=R[p-1]..R[p]),p=1..3)/g[l,l] od:
r:=R[0]:
T1:=sum(C[j]*M[j,1]*exp(-z[j]^2*t)*(cos(z[j]*r/sqrt(a[1]))+y[j,1]* sin(z[j]*r/sqrt(a[1]))),j=1..5)/r+TC-273+A[1]/r+B[1]:
plot(T1, t=0.01..5,T=-10..40,color=black);
«