timour_matlab/mscripts/smoothheatprofile.m

Nx=124;
Ny=Nx;
      Heating = 2.0;

      Rbuoy=0.02;
      Rheat=0.15;
      a1=(Rbuoy-2.0*Rheat)/(4.0*Rbuoy*Rheat);
      a2=1.0/(8.0*Rbuoy*Rheat);
      delX=0.01;
      iC=Nx/2;
      jC=Ny/2;
          for J=1:Ny
          for I=1:Nx
            iD = I-iC;
            jD = J-jC;
            rD = sqrt(iD*iD+jD*jD)*delX;
            if ( rD < (Rheat-Rbuoy) ) 
              Heat(I,J) = Heating;
            elseif (rD < (Rheat+Rbuoy)) 
              Rcool = rD - Rheat - Rbuoy;
              Heat(I,J) = Heating*(a2*Rcool^(2.0)+a1*Rcool);
            else
              Heat(I,J) = 0;
            end
          end
        end
figure
pcolor(Heat(:,:));axis square;colorbar;shading flat
figure
plot(Heat(:,Nx/2))
1	Nx=124;
2	Ny=Nx;
3	Heating = 2.0;
4
5	Rbuoy=0.02;
6	Rheat=0.15;
7	a1=(Rbuoy-2.0Rheat)/(4.0Rbuoy*Rheat);
8	a2=1.0/(8.0RbuoyRheat);
9	delX=0.01;
10	iC=Nx/2;
11	jC=Ny/2;
12	for J=1:Ny
13	for I=1:Nx
14	iD = I-iC;
15	jD = J-jC;
16	rD = sqrt(iDiD+jDjD)*delX;
17	if ( rD < (Rheat-Rbuoy) )
18	Heat(I,J) = Heating;
19	elseif (rD < (Rheat+Rbuoy))
20	Rcool = rD - Rheat - Rbuoy;
21	Heat(I,J) = Heating(a2Rcool^(2.0)+a1*Rcool);
22	else
23	Heat(I,J) = 0;
24	end
25	end
26	end
27	figure
28	pcolor(Heat(:,:));axis square;colorbar;shading flat
29	figure
30	plot(Heat(:,Nx/2))