h2 = {{0, -Ω/2}, {-Ω/2, k v - Δ}};
ρ = {{ρ11[t], ρ12[t]}, {ρ21[t], ρ22[t]}};
ρprime = -I (h2.ρ - ρ.h2) + {{1/2 γ ρ22[t], -γ ρ12[t]}, {-γ ρ21[t], -1/2 γ ρ22[t]}};
replace3 = {Δ -> 1/(16*10^-9), γ -> 1/(16*10^-9), Ω -> 10^8, k -> ( 2 [Pi])/(500*10^-9)};

var = Flatten@ρ;

{eq, ic} = {D[var, t] == Flatten@[Rho]prime // Thread,
var == {1, 0, 0, 0} /. t -> 0 // Thread} /. replace3;

tvar = LaplaceTransform[var, t, s];
tsol = tvar /.
First@Solve[LaplaceTransform[eq, t, s] /. Rule @@@ ic, tvar] // Simplify;

(sol = InverseLaplaceTransform[tsol, s, t]) // AbsoluteTiming;

ContourPlot[v == 10^3 - 700000 t Re[sol[[4]]], {t, 0, 10}, {v, 0, 10^3}]