Unsteady FLLC

amr-wind/wind_energy/actuator/FLLC.H

@@ -46,6 +46,28 @@ struct FLLCData
     VecList nonuniform_lift;             // non-uniform lift
 };
 
+struct FLLCUnsteadyData
+{
+    // constants
+    amrex::Real epsilon;
+    amrex::Real fllc_start_time{0.0};
+    int n_times_history{100};
+
+    bool initialized{false};
+    RealList times;
+    RealList dt;
+
+    // computed values
+    // a vector of a vector of a vector:
+    //  - listed in time (n_times_history)
+    //  - listed in space (npts)
+    //  - coordinate directions
+    amrex::Vector<VecList> induced_velocity_eps_LES;
+    amrex::Vector<VecList> induced_velocity_eps_opt;
+    amrex::Vector<VecList> delta_velocity;
+    amrex::Vector<VecList> force_coeff;
+};
+
 /**
  * \brief Function to capture common parsing requirements for the filtered
  * lifting line correction
@@ -55,6 +77,8 @@ struct FLLCData
  */
 void fllc_parse(const utils::ActParser& pp, FLLCData& data);
 
+void fllc_unsteady_parse(const utils::ActParser& pp, FLLCUnsteadyData& data);
+
 /**
  * \brief Initialize FLLC data structure. This should be called at the end of
  * the first ComputeForceOp to ensure the data is fully populated
@@ -66,6 +90,11 @@ void fllc_parse(const utils::ActParser& pp, FLLCData& data);
 void fllc_init(
     FLLCData& data, const ComponentView& view, const amrex::Real eps_chord);
 
+void fllc_unsteady_init(
+    FLLCUnsteadyData& data,
+    const ComponentView& view,
+    const amrex::Real eps_chord);
+
 } // namespace amr_wind::actuator
 
 #endif /* FLLC */

amr-wind/wind_energy/actuator/FLLC.cpp

@@ -129,4 +129,53 @@ void fllc_parse(const utils::ActParser& pp, FLLCData& data)
     }
 }
 
+void fllc_unsteady_init(
+    FLLCUnsteadyData& data,
+    const ComponentView& view,
+    const amrex::Real eps_chord)
+{
+
+    // Do init work here
+
+    // Set up time-related lists
+    data.times.resize(data.n_times_history);
+    data.dt.resize(data.n_times_history - 1);
+    data.times.assign(data.n_times_history, 0.0);
+    data.dt.assign(data.n_times_history - 1, 0.0);
+
+    // Set up velocity lists according to n_times_history, npts
+    data.induced_velocity_eps_LES.resize(data.n_times_history);
+    data.induced_velocity_eps_opt.resize(data.n_times_history);
+    data.delta_velocity.resize(data.n_times_history);
+    data.force_coeff.resize(data.n_times_history);
+
+    const int npts = static_cast<int>(view.pos.size());
+    for (int n = 0; n < npts; ++n) {
+        data.induced_velocity_eps_LES[n].resize(data.n_times_history);
+        data.induced_velocity_eps_opt[n].resize(data.n_times_history);
+        data.delta_velocity[n].resize(data.n_times_history);
+        data.force_coeff[n].resize(data.n_times_history);
+
+        data.induced_velocity_eps_LES[n].assign(npts, vs::Vector::zero());
+        data.induced_velocity_eps_opt[n].assign(npts, vs::Vector::zero());
+        data.delta_velocity[n].assign(npts, vs::Vector::zero());
+        data.force_coeff[n].assign(npts, vs::Vector::zero());
+    }
+
+    data.initialized = true;
+}
+
+void fllc_unsteady_parse(const utils::ActParser& pp, FLLCUnsteadyData& data)
+{
+    pp.query("epsilon", data.epsilon);
+    pp.query("fllc_start_time", data.fllc_start_time);
+    pp.query("fllc_unsteady_number_timesteps", data.n_times_history);
+
+    if (!pp.contains("epsilon") || !pp.contains("epsilon_chord")) {
+        amrex::Abort(
+            "Actuators using the filtered lifting line correction (FLLC) "
+            "require specification 'epsilon' and 'epsilon_chord'");
+    }
+}
+
 } // namespace amr_wind::actuator

amr-wind/wind_energy/actuator/FLLCOp.H

@@ -293,6 +293,46 @@ struct FLLCOp
     }
 };
 
+struct FLLCUnsteadyOp
+{
+    void operator()(
+        ComponentView& data,
+        FLLCUnsteadyData& fllc,
+        const amrex::Real current_time)
+    {
+        if (!fllc.initialized) {
+            return;
+        }
+        unsteady_correction(data, fllc, current_time);
+    }
+
+    static void unsteady_correction(
+        ComponentView& data,
+        FLLCUnsteadyData& fllc,
+        const amrex::Real current_time)
+    {
+        // Update time list
+        // Shift elements toward back
+        for (int n = fllc.times.size() - 1; n > 0; ++n) {
+            fllc.times[n] = fllc.times[n - 1];
+        }
+        // Replace new first element
+        fllc.times[0] = current_time;
+        // Recalculate dt
+        for (int n = 0; n < fllc.dt.size(); ++n) {
+            fllc.dt[n] = fllc.times[n + 1] - fllc.times[n];
+        }
+
+        // TEST COMMENT
+
+        // Store current velocities to lists
+
+        // Compute integrals
+
+        // Update velocities with correction
+    }
+};
+
 } // namespace amr_wind::actuator
 
 #endif /* FLLCOP_H */

amr-wind/wind_energy/actuator/wing/ActuatorWing.H

@@ -9,6 +9,7 @@
 namespace amr_wind::actuator {
 
 struct FLLCData;
+struct FLLCUnsteadyData;
 
 /** Base data representation of a fixed-wing.
  *
@@ -73,6 +74,7 @@ struct WingBaseData
 
     //! Filtered Lifting Line Correction data
     std::unique_ptr<FLLCData> fllc;
+    std::unique_ptr<FLLCUnsteadyData> fllc_unsteady;
 
     //! Data structure required for passing data to the FLLC
     ComponentView component_view;

amr-wind/wind_energy/actuator/wing/fixed_wing_ops.H

@@ -31,6 +31,9 @@ struct ReadInputsOp<FixedWing, ActSrcLine>
         bool use_fllc = false;
         pp.query("fllc", use_fllc);
 
+        bool use_fllc_unsteady = false;
+        pp.query("fllc_unsteady", use_fllc_unsteady);
+
         // If spanwise components of Gaussian should be ignored
         pp.query("disable_spanwise_gaussian", wdata.gauss_2D);
         // If spanwise components of Gaussian and forces should be normalized
@@ -89,6 +92,11 @@ struct ReadInputsOp<FixedWing, ActSrcLine>
             fllc_parse(pp, *(wdata.fllc));
         }
 
+        if (use_fllc_unsteady) {
+            wdata.fllc_unsteady = std::make_unique<FLLCUnsteadyData>();
+            fllc_unsteady_parse(pp, *(wdata.fllc_unsteady));
+        }
+
         if (!pp.contains("epsilon") && !pp.contains("epsilon_chord")) {
             amrex::Abort(
                 "Actuator fixed wing requires specification of one or both "

amr-wind/wind_energy/actuator/wing/wing_ops.H

@@ -109,6 +109,11 @@ struct UpdateVelOp<
         if (meta.fllc) {
             FLLCOp()(meta.component_view, *(meta.fllc));
         }
+        if (meta.fllc_unsteady) {
+            FLLCUnsteadyOp()(
+                meta.component_view, *(meta.fllc_unsteady),
+                data.sim().time().current_time());
+        }
     }
 };
 
@@ -249,6 +254,17 @@ struct ComputeForceOp<
                     wdata.epsilon_chord[0]);
             }
         }
+        if (wdata.fllc_unsteady) {
+            if (!(wdata.fllc_unsteady->initialized) &&
+                (time.current_time() > wdata.fllc_unsteady->fllc_start_time)) {
+                wdata.component_view =
+                    amr_wind::actuator::wing::make_component_view<ActTrait>(
+                        data);
+                fllc_unsteady_init(
+                    *(wdata.fllc_unsteady.get()), wdata.component_view,
+                    wdata.epsilon_chord[0]);
+            }
+        }
     }
 };