%vary r
clear
figure(1)
[lvect,gvect,pvect,peak(1),energy(1)]=qswitchsim(100,1.5,250);
[lvect,gvect,pvect,peak(2),energy(2)]=qswitchsim(100,1.75,250);
[lvect,gvect,pvect,peak(3),energy(3)]=qswitchsim(100,2,250);
[lvect,gvect,pvect,peak(4),energy(4)]=qswitchsim(100,2.25,250);
[lvect,gvect,pvect,peak(5),energy(5)]=qswitchsim(100,2.5,250);
[lvect,gvect,pvect,peak(6),energy(6)]=qswitchsim(100,2.75,250);
[lvect,gvect,pvect,peak(7),energy(7)]=qswitchsim(100,3,250);
[lvect,gvect,pvect,peak(8),energy(8)]=qswitchsim(100,3.25,250);
[lvect,gvect,pvect,peak(9),energy(9)]=qswitchsim(100,3.5,250);
[lvect,gvect,pvect,peak(10),energy(10)]=qswitchsim(100,3.75,250);
[lvect,gvect,pvect,peak(11),energy(11)]=qswitchsim(100,4,250);
[lvect,gvect,pvect,peak(12),energy(12)]=qswitchsim(100,4.25,250);
[lvect,gvect,pvect,peak(13),energy(13)]=qswitchsim(100,4.5,250);
[lvect,gvect,pvect,peak(14),energy(14)]=qswitchsim(100,4.75,250);
[lvect,gvect,pvect,peak(15),energy(15)]=qswitchsim(100,5,250);
rs=[1.5 1.75 2 2.25 2.5 2.75 3 3.25 3.5 3.75 4 4.25 4.5 4.75 5];
figure(2)
subplot(2,1,1)
plot(rs, energy)
title('energy vs. r')
ylabel('energy')
subplot(2,1,2)
plot(rs, peak)
title('peak vs. r')
xlabel('Normalized Time')
ylabel('peak')