a=1/16;c=pi^2;
t=1:.01:.5
z=zeros(1,length(t));s1=zeros(1,length(t));s2=zeros(1,length(t));
s3=zeros(1,length(t));s=zeros(1,length(t));t=1:.01:.5;
z=zeros(1,length(t)); s1=zeros(1,length(t));s2=zeros(1,length(t));
s3=zeros(1,length(t));
for p=1:length(t)
for l=1:100
s1(p)=s1(p)+cos(2*pi*l*t(p))*exp(-2*(l*pi)^2*a);
end

 


for m=1:100;
for n=1:100; 
s2(p)=s2(p)+(1/(m*n))*(c0s(2*pi*(m-n)*t(p))*exp(-2*(m-n)^2*c*a)-cos(2*pi*(m+n)*t(p))*exp(-2*(m+n)^2*c*a));
end
end
for m=1:100
for n=1:100
s3(p)=s3(p)+(1/(m*n))*exp(-(m^2+n^2)*c*a)*sin(2*pi*n*t(p))*sin(2*pi*m*t(p));
end
end
s(p)=(s1(p)*4+1)*a+(1/c)*(.5*s2(p)-s3(p));
end
plot(t,s)