function hcaculate
[filename,pathname]=uigetfile('*.bmp');
file=[pathname,filename];
I=imread(file);
I1=rgb2gray(I);
Height=double(I1);
Hmax=double(max(max(I1)));
[m1 m2 m3]=size(I1);
xlabel('width','FontSize',12),ylabel('length','FontSize',12),zlabel('Height','FontSize',12);
Ds=(5.11*5)/(m1*m2);%area of a single grid
load= 47.2675;
fC=0.18;%friction coefficient
v=0.3;%sliding velocity
Hdist=7.88928291806752E-02;%heat distribution coefficient
Hcond=0.14;%heat conductivity
Time=1;%testing duration
e=zeros(m1,m2);%strain
f=zeros(m1,m2);%load
dw=zeros(m1,m2);%friction work
for n=1:150
for i=1:m1
for j=1:m2
if(Height(i,j)>=Hmax-n)
e(i,j)=(Height(i,j)-(Hmax-n))/Height(i,j);
if(e(i,j)<=0.02)
f(i,j)=e(i,j)*Ds*695.16;
else
f(i,j)=1.7129*e(i,j).^4-1.4222*e(i,j).^3+0.1745*e(i,j).^2+0.1541*e(i,j)+0.0107;
end
end
end
end
sf(n)=sum(sum(f));
sload(n)=abs(sf(n)-load);
end
n=find(sload==min(sload))
e=zeros(m1,m2);
f=zeros(m1,m2);
dw=zeros(m1,m2);
for i=1:m1
for j=1:m2
if(Height(i,j)>=Hmax-n)
e(i,j)=(Height(i,j)-(Hmax-n))/Height(i,j);
if(e(i,j)<=0.02)
f(i,j)=e(i,j)*Ds*695.16;
else
f(i,j)=1.7129*e(i,j).^4-1.4222*e(i,j).^3+0.1745*e(i,j).^2+0.1541*e(i,j)+0.0107;
end
end
dw(i, j) = f(i, j) * fC * v;
pt(i, j) = 1000 * dw(i, j) * Hdist * (1 + Time* 60 * 0.001) / (4 * Hcond * ((Ds/ 3.14) .^ 0.5));
end
end
for k=1:Time*30
for i=2:m1-1
for j=21:m2-1
ptt(i,j)=(pt(i-1,j)+pt(i+1,j)+pt(i,j-1)+pt(i,j+1))*0.25;
end
end
for i=2:m1-1
for j=21:m2-1
if(abs(ptt(i, j) - pt(i, j)) > 0.001)
pt(i, j) = ptt(i, j);
end
end
end
end
figure;
surf(ptt),shading interp,colorbar;
axis([1 m2 1 m1]), axis('ij'), shading('interp'),colorbar,view(30,80);
xlabel('width','FontSize',12),ylabel('length','FontSize',12),zlabel('Temperature/C','FontSize',12);
