%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Numerische Blochgleichungen mit Relaxation, n = (n11 n22 n12 n21).
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function Bloch

figure(1);					clf;
V=[0:.1:2*pi];

subplot(221),plot3(sin(V),cos(V),0*V,'b','MarkerSize',10,'LineWidth',1.2)
set(gca,'xLim',[-1 1],'yLim',[-1 1],'zLim',[-1 1],'FontSize',18,'LineWidth',1.2)
hold on
subplot(221),plot3(0*V,sin(V),cos(V),'b','MarkerSize',10,'LineWidth',1.2)
set(gca,'xLim',[-1 1],'yLim',[-1 1],'zLim',[-1 1],'FontSize',18,'LineWidth',1.2)
subplot(221),plot3(cos(V),0*V,sin(V),'b','MarkerSize',10,'LineWidth',1.2)
set(gca,'xLim',[-1 1],'yLim',[-1 1],'zLim',[-1 1],'FontSize',18,'LineWidth',1.2)

dt  =      1e-9;
Ga  =  2*pi*0e6;
rho0=[1 0 0 0]';	
De  =  2*pi*6e6;
Om  =  2*pi*6e6;
Ta  =     pi/De*.15;

for jj=1:2
	T =Ta*jj;
	t =0;
	t1=pi/2/Om;
	t2=pi/2/Om+T;
	t3=pi/2/Om+T+pi/Om;
	t4=pi/2/Om+T+pi/Om+T;
	t5=pi/2/Om+T+pi/Om+T+pi/2/Om;
	
	while t<t5
		d=0 *(t<t1) + De*(t1<t)*(t<t2) + 0 *(t2<t)*(t<t3) + De*(t3<t)*(t<t4) + 0 *(t4<t)*(t<t5-2*dt);
		o=Om*(t<t1) + 0 *(t1<t)*(t<t2) + Om*(t2<t)*(t<t3) + 0 *(t3<t)*(t<t4) + Om*(t4<t)*(t<t5-2*dt);

		M=[	0			Ga			i/2*o		-i/2*o
			0			-Ga			-i/2*o		i/2*o
			i/2*o		-i/2*o		-Ga/2-i*d	0
			-i/2*o		i/2*o		0			-Ga/2+i*d];
	
		t=t+dt;
		rho=expm(M*dt)*rho0;
	
		SubPlot(t,dt,rho0,rho,jj);
		rho0=rho;
	end
end
hold off
end


function SubPlot(t,dt,rho0,rho,jj)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Plots.

	col='rg';
	n0(1)=2*imag(rho0(3));
	n0(2)=2*real(rho0(3));
	n0(3)=rho0(2)-rho0(1);
	n(1) =2*imag(rho(3));
	n(2) =2*real(rho(3));
	n(3) =rho(2)-rho(1);
	
	subplot(221),plot3([n0(1) n(1)],[n0(2) n(2)],[n0(3) n(3)],col(jj),'MarkerSize',10,'LineWidth',1.8)
	set(gca,'xLim',[-1 1],'yLim',[-1 1],'zLim',[-1 1],'FontSize',18,'LineWidth',1.2)
	title(['|\rho| =  ',num2str(norm(n))])
	xlabel('\fontsize{18}Im {\it\rho_{12}}')
	ylabel('\fontsize{18}Re {\it\rho_{12}}')
	zlabel('\fontsize{18}2{\it\rho_{22}}-1')
	subplot(222),plot([t-dt t]*1e6,[n0(3) n(3)],col(jj),'MarkerSize',10,'LineWidth',1.4)
	set(gca,'xLim',[0 t]*1e6,'yLim',[-1 1],'FontSize',18,'LineWidth',1.2)
	xlabel('\fontsize{18}{\itt}  (\mus)')
	ylabel('\fontsize{18}2{\it\rho_{22}}-1')
	hold on
	subplot(223),plot([t-dt t]*1e6,[n0(1) n(1)],col(jj),'MarkerSize',10,'LineWidth',1.4)
	set(gca,'xLim',[0 t]*1e6,'yLim',[-1 1],'FontSize',18,'LineWidth',1.2)
	xlabel('\fontsize{18}{\itt}  (\mus)')
	ylabel('\fontsize{18}Im {\it\rho_{12}}')
	hold on
	subplot(224),plot([t-dt t]*1e6,[n0(2) n(2)],col(jj),'MarkerSize',10,'LineWidth',1.4)
	set(gca,'xLim',[0 t]*1e6,'yLim',[-1 1],'FontSize',18,'LineWidth',1.2)
	xlabel('\fontsize{18}{\itt}  (\mus)')
	ylabel('\fontsize{18}Re {\it\rho_{12}}')
	hold on
	pause(.01);
end