σz=PauliMatrix[3] ;
σx=PauliMatrix[1];
σy=PauliMatrix[2];

H=((ℏ ω0)/2)σz;
MatrixForm[H]
(* Ec. de Schrodinger es I (d/dt)ψ = H ψ  o I ℏ (d/dt)ψ - H ψ =0 *)

ψ={ψ1[t],ψ2[t]};

mat=I ℏ D[ψ,t]-H.ψ

eqs=Table[el==0,{el,mat}]

eqsn=eqs/.{ℏ->1,ω0->1}

ci={ψ1[0]==1,ψ2[0]==0}

deqsn=Join[eqsn,ci]

sol=DSolve[deqsn,ψ,t]

(* Sx=<ψ|((σx)/2)|ψ> *)

Sx=Conjugate[ψ].σx.ψ/.sol[[1]]
Sy=Conjugate[ψ].σy.ψ/.sol[[1]]
Sz=Conjugate[ψ].σz.ψ/.sol[[1]]

Simplify[Sz,Assumptions->{t>0}]

H=((ℏ ω0)/2)σz;
MatrixForm[H]
ψ={ψ1[t],ψ2[t]};
mat=I ℏ D[ψ,t]-H.ψ;
eqs=Table[el==0,{el,mat}];
eqsn=eqs/.{ℏ->1,ω0->1};
ci={ψ1[0]==(1/(Sqrt[2])),ψ2[0]==(1/(Sqrt[2]))};
deqsn=Join[eqsn,ci];
sol=DSolve[deqsn,ψ,t];
Sx=Simplify[Conjugate[ψ].σx.ψ/.sol[[1]],Assumptions->{t>0}]
Sy=Simplify[Conjugate[ψ].σy.ψ/.sol[[1]],Assumptions->{t>0}]
Sz=Simplify[Conjugate[ψ].σz.ψ/.sol[[1]],Assumptions->{t>0}]

ParametricPlot3D[{Sx,Sy,Sz},{t,0,8}]

H=-((ℏ ω0)/2)σz+((ℏ ω1)/2)(Cos[ω t]σx-Sin[ω t]σy);
MatrixForm[H]
(* Ec. de Schrodinger es I (d/dt)ψ = H ψ  o I ℏ (d/dt)ψ - H ψ =0 *)

ψ={ψ1[t],ψ2[t]};

mat=I ℏ D[ψ,t]-H.ψ

eqs=Table[el==0,{el,mat}]

eqsn=eqs/.{ℏ->1,ω0->1,ω->0.5,ω1->0.2}

ci={ψ1[0]==1,ψ2[0]==0}

deqsn=Join[eqsn,ci]

sol=NDSolve[deqsn,ψ,{t,0,20}]

(* Sx=<ψ|((σx)/2)|ψ> *)

Sx=Conjugate[ψ].σx.ψ/.sol[[1]]
Sy=Conjugate[ψ].σy.ψ/.sol[[1]]
Sz=Conjugate[ψ].σz.ψ/.sol[[1]]

Sx=Simplify[Conjugate[ψ].σx.ψ/.sol[[1]],Assumptions->{t>0}]
Sy=Simplify[Conjugate[ψ].σy.ψ/.sol[[1]],Assumptions->{t>0}]
Sz=Simplify[Conjugate[ψ].σz.ψ/.sol[[1]],Assumptions->{t>0}]

ParametricPlot3D[{Sx,Sy,Sz},{t,0,18}]

P=((Abs[ψ[[2]]/.sol[[1]]])^2)

Plot[P,{t,0,20}]

H=-((ℏ ω0)/2)σz+((ℏ ω1)/2)(Cos[ω t]σx-Sin[ω t]σy);
mat=I ℏ D[ψ,t]-H.ψ;
eqs=Table[el==0,{el,mat}];
eqsn=eqs/.{ℏ->1,ω0->1,ω->1,ω1->0.2}
ci={ψ1[0]==1,ψ2[0]==0};
deqsn=Join[eqsn,ci];
sol=NDSolve[deqsn,ψ,{t,0,100}];
P=((Abs[ψ2[t]/.sol[[1]]])^2);
Plot[P,{t,0,50}]

H=-((ℏ ω0)/2)σz+((ℏ ω1)/2)(Cos[ω t]σx-Sin[ω t]σy);
mat=I ℏ D[ψ,t]-H.ψ;
eqs=Table[el==0,{el,mat}];
eqsn=eqs/.{ℏ->1,ω0->1,ω->1,ω1->0.8}
ci={ψ1[0]==1,ψ2[0]==0};
deqsn=Join[eqsn,ci];
sol=NDSolve[deqsn,ψ,{t,0,100}];
P=((Abs[ψ2[t]/.sol[[1]]])^2);
Plot[P,{t,0,50}]

conds={ℏ->1,ω0->1,ω->1,ω1->0.8,t0->20}

Plot[{Cos[ω t](1/(Sqrt[2π]))Exp[-(1/2)(((((t-t0)/((π/ω1)))))^2)]/.conds,Exp[-(1/2)(((((t-t0)/((π/ω1)))))^2)]/.conds},{t,0,40},PlotRange->All]

H=-((ℏ ω0)/2)σz+((ℏ ω1)/2)(Cos[ω t]σx-Sin[ω t]σy)(1/(Sqrt[2π]))Exp[-(1/2)(((((t-t0)/((π/ω1)))))^2)];
mat=I ℏ D[ψ,t]-H.ψ;
eqs=Table[el==0,{el,mat}];
eqsn=eqs/.conds
ci={ψ1[0]==1,ψ2[0]==0};
deqsn=Join[eqsn,ci];
sol=NDSolve[deqsn,ψ,{t,0,100}];
P=((Abs[ψ2[t]/.sol[[1]]])^2);
Plot[P,{t,0,50}]

ψ/.sol[[1]]/.{t->0}//MatrixForm

ψ/.sol[[1]]/.{t->50}//MatrixForm