diff --git a/src/generate_data/Perturbations.jl b/src/generate_data/Perturbations.jl
index 85cc96a..df0c95d 100644
--- a/src/generate_data/Perturbations.jl
+++ b/src/generate_data/Perturbations.jl
@@ -494,8 +494,11 @@ function add_surrogate_perturbation!(
         load_multiplier_range[1]
     Pnew = PSY.get_impedance_active_power(l_new) * multiplier
     Qnew = PSY.get_impedance_reactive_power(l_new) * multiplier
-    push!(psid_perturbations, PSID.LoadChange(time, l_new, :P_ref_impedance, Pnew))
-    push!(psid_perturbations, PSID.LoadChange(time, l_new, :Q_ref_impedance, Qnew))
+    psid_perturbation_1 = PSID.LoadChange(time, l_new, :P_ref_impedance, Pnew)
+    psid_perturbation_2 = PSID.LoadChange(time, l_new, :Q_ref_impedance, Qnew)
+    @info "Perturbation", psid_perturbation_1, psid_perturbation_2
+    push!(psid_perturbations, psid_perturbation_1)
+    push!(psid_perturbations, psid_perturbation_2)
 end
 
 ###############################################################################
@@ -533,25 +536,19 @@ function add_surrogate_perturbation!(
             sys_aux,
         ),
     )
-    println(PSY.get_name.(static_injectors))
     if length(static_injectors) === 0
         @error "Trying to change a dynamic injector but a dynamic injector not found in system"
         return
     end
     s = rand(static_injectors)
     s_new = PSY.get_component(typeof(s), sys, PSY.get_name(s))
-    println(PSY.get_active_power(s_new))
     multiplier =
         rand() * (P_multiplier_range[2] - P_multiplier_range[1]) + P_multiplier_range[1]
-    println(multiplier)
     Pnew = PSY.get_P_ref(PSY.get_dynamic_injector(s_new)) * multiplier
     #Qnew = PSY.get_impedance_reactive_power(l_new) * multiplier
-    println(
-        PSID.ControlReferenceChange(time, PSY.get_dynamic_injector(s_new), :P_ref, Pnew),
-    )
-    push!(
-        psid_perturbations,
-        PSID.ControlReferenceChange(time, PSY.get_dynamic_injector(s_new), :P_ref, Pnew),
-    )
+    psid_perturbation =
+        PSID.ControlReferenceChange(time, PSY.get_dynamic_injector(s_new), :P_ref, Pnew)
+    @info "Perturbation", psid_perturbation
+    push!(psid_perturbations, psid_perturbation)
     #push!(psid_perturbations, PSID.LoadChange(time, d_new, :Q_ref_impedance, Qnew))
 end
diff --git a/src/generate_data/datasets/generate.jl b/src/generate_data/datasets/generate.jl
index 0f7574d..9a2aff4 100644
--- a/src/generate_data/datasets/generate.jl
+++ b/src/generate_data/datasets/generate.jl
@@ -147,6 +147,7 @@ function generate_surrogate_data(
     ########################################################################################
     sys = deepcopy(sys_main)
     update_operating_point!(sys, operating_points[1], sys_aux)
+    PSID.Simulation!(PSID.MassMatrixModel, sys, pwd(), (0.0, 0.0))  #Run power flow and re-init devices before defining perturbation
     psid_perturbations = PSID.Perturbation[]
     for p_single in perturbations[1]
         add_surrogate_perturbation!(sys, psid_perturbations, p_single, sys_aux)
@@ -192,6 +193,7 @@ function generate_surrogate_data(
         for (ix_p, p) in enumerate(perturbations)
             sys = deepcopy(sys_main)
             update_operating_point!(sys, o, sys_aux)
+            PSID.Simulation!(PSID.MassMatrixModel, sys, pwd(), (0.0, 0.0))  #Run power flow and re-init devices before defining perturbation
             psid_perturbations = PSID.Perturbation[]
             for p_single in p
                 add_surrogate_perturbation!(sys, psid_perturbations, p_single, sys_aux)