diff --git a/Project.toml b/Project.toml
index 19326adbe..5e859db13 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "NQCDynamics"
 uuid = "36248dfb-79eb-4f4d-ab9c-e29ea5f33e14"
 authors = ["James <james.gardner1421@gmail.com>"]
-version = "0.13.5"
+version = "0.13.6"
 
 [deps]
 AdvancedHMC = "0bf59076-c3b1-5ca4-86bd-e02cd72cde3d"
diff --git a/src/DynamicsOutputs.jl b/src/DynamicsOutputs.jl
index c147374dd..01e964d52 100644
--- a/src/DynamicsOutputs.jl
+++ b/src/DynamicsOutputs.jl
@@ -172,6 +172,12 @@ Output the total adiabatic population at each timestep during the trajectory.
 OutputTotalAdiabaticPopulation(sol, i) = sum.(Estimators.adiabatic_population.(sol.prob.p, sol.u))
 export OutputTotalAdiabaticPopulation
 
+"""
+Output the first point of each trajectory in DynamicsVariables format. (Useful when using distributions for initial conditions.)
+"""
+OutputInitial(sol, i) = first(sol.u)
+export OutputInitial
+
 """
 Output the end point of each trajectory.
 """