ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
GeneralAD Class Reference

#include <ERF_GeneralAD.H>

Inheritance diagram for GeneralAD:
Collaboration diagram for GeneralAD:

Public Member Functions

 GeneralAD ()
 
virtual ~GeneralAD ()=default
 
void advance (const amrex::Geometry &geom, const double &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &mf_vars_windfarm, amrex::MultiFab &U_old, amrex::MultiFab &V_old, amrex::MultiFab &W_old, const amrex::MultiFab &mf_Nturb, const amrex::MultiFab &mf_SMark, const double &time) override
 
void compute_freestream_velocity (const amrex::MultiFab &cons_in, const amrex::MultiFab &U_old, const amrex::MultiFab &V_old, const amrex::MultiFab &mf_SMark)
 
void source_terms_cellcentered (const amrex::Geometry &geom, const amrex::MultiFab &cons_in, const amrex::MultiFab &mf_Smark, amrex::MultiFab &mf_vars_generalAD)
 
void update (const double &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &U_old, amrex::MultiFab &V_old, amrex::MultiFab &W_old, const amrex::MultiFab &mf_vars)
 
void compute_power_output (const double &time)
 
- Public Member Functions inherited from NullWindFarm
 NullWindFarm ()
 
virtual ~NullWindFarm ()=default
 
virtual void set_turb_spec (const amrex::Real &rotor_rad, const amrex::Real &hub_height, const amrex::Real &thrust_coeff_standing, const amrex::Vector< amrex::Real > &wind_speed, const amrex::Vector< amrex::Real > &thrust_coeff, const amrex::Vector< amrex::Real > &power)
 
virtual void set_turb_loc (const amrex::Vector< amrex::Real > &xloc, const amrex::Vector< amrex::Real > &yloc)
 
virtual void set_turb_disk_angle (const amrex::Real &turb_disk_angle)
 
virtual void set_blade_spec (const amrex::Vector< amrex::Real > &bld_rad_loc, const amrex::Vector< amrex::Real > &bld_twist, const amrex::Vector< amrex::Real > &bld_chord)
 
virtual void set_blade_airfoil_spec (const amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_aoa, const amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cl, const amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cd)
 
virtual void set_turb_spec_extra (const amrex::Vector< amrex::Real > &velocity, const amrex::Vector< amrex::Real > &C_P, const amrex::Vector< amrex::Real > &C_T, const amrex::Vector< amrex::Real > &rotor_RPM, const amrex::Vector< amrex::Real > &blade_pitch)
 
void get_turb_spec (amrex::Real &rotor_rad, amrex::Real &hub_height, amrex::Real &thrust_coeff_standing, amrex::Vector< amrex::Real > &wind_speed, amrex::Vector< amrex::Real > &thrust_coeff, amrex::Vector< amrex::Real > &power)
 
void get_turb_loc (amrex::Vector< amrex::Real > &xloc, amrex::Vector< amrex::Real > &yloc)
 
void get_turb_disk_angle (amrex::Real &turb_disk_angle)
 
void get_blade_spec (amrex::Vector< amrex::Real > &bld_rad_loc, amrex::Vector< amrex::Real > &bld_twist, amrex::Vector< amrex::Real > &bld_chord)
 
void get_blade_airfoil_spec (amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_aoa, amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cl, amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cd)
 
void get_turb_spec_extra (amrex::Vector< amrex::Real > &velocity, amrex::Vector< amrex::Real > &C_P, amrex::Vector< amrex::Real > &C_T, amrex::Vector< amrex::Real > &rotor_RPM, amrex::Vector< amrex::Real > &blade_pitch)
 

Protected Attributes

amrex::Vector< amrex::Realxloc
 
amrex::Vector< amrex::Realyloc
 
amrex::Real turb_disk_angle
 
amrex::Real hub_height
 
amrex::Real rotor_rad
 
amrex::Real thrust_coeff_standing
 
amrex::Real nominal_power
 
amrex::Vector< amrex::Realwind_speed
 
amrex::Vector< amrex::Realthrust_coeff
 
amrex::Vector< amrex::Realpower
 
amrex::Vector< amrex::Realfreestream_velocity
 
amrex::Vector< amrex::Realfreestream_phi
 
amrex::Vector< amrex::Realdisk_cell_count
 
amrex::Vector< amrex::Realbld_rad_loc
 
amrex::Vector< amrex::Realbld_twist
 
amrex::Vector< amrex::Realbld_chord
 
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_aoa
 
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_Cl
 
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_Cd
 
amrex::Vector< amrex::Realvelocity
 
amrex::Vector< amrex::RealC_P
 
amrex::Vector< amrex::RealC_T
 
amrex::Vector< amrex::Realrotor_RPM
 
amrex::Vector< amrex::Realblade_pitch
 
- Protected Attributes inherited from NullWindFarm
amrex::Vector< amrex::Realm_xloc
 
amrex::Vector< amrex::Realm_yloc
 
amrex::Real m_turb_disk_angle
 
amrex::Real m_hub_height
 
amrex::Real m_rotor_rad
 
amrex::Real m_thrust_coeff_standing
 
amrex::Real m_nominal_power
 
amrex::Vector< amrex::Realm_wind_speed
 
amrex::Vector< amrex::Realm_thrust_coeff
 
amrex::Vector< amrex::Realm_power
 
amrex::Vector< amrex::Realm_bld_rad_loc
 
amrex::Vector< amrex::Realm_bld_twist
 
amrex::Vector< amrex::Realm_bld_chord
 
amrex::Vector< amrex::Vector< amrex::Real > > m_bld_airfoil_aoa
 
amrex::Vector< amrex::Vector< amrex::Real > > m_bld_airfoil_Cl
 
amrex::Vector< amrex::Vector< amrex::Real > > m_bld_airfoil_Cd
 
amrex::Vector< amrex::Realm_velocity
 
amrex::Vector< amrex::Realm_C_P
 
amrex::Vector< amrex::Realm_C_T
 
amrex::Vector< amrex::Realm_rotor_RPM
 
amrex::Vector< amrex::Realm_blade_pitch
 

Additional Inherited Members

- Static Public Member Functions inherited from NullWindFarm
static AMREX_GPU_DEVICE bool find_if_marked (amrex::Real x1, amrex::Real x2, amrex::Real y1, amrex::Real y2, amrex::Real x0, amrex::Real y0, amrex::Real nx, amrex::Real ny, amrex::Real d_hub_height, amrex::Real d_rotor_rad, amrex::Real z)
 

Constructor & Destructor Documentation

◆ GeneralAD()

GeneralAD::GeneralAD ( )
inline
12 {}

◆ ~GeneralAD()

virtual GeneralAD::~GeneralAD ( )
virtualdefault

Member Function Documentation

◆ advance()

void GeneralAD::advance ( const amrex::Geometry &  geom,
const double &  dt_advance,
amrex::MultiFab &  cons_in,
amrex::MultiFab &  mf_vars_windfarm,
amrex::MultiFab &  U_old,
amrex::MultiFab &  V_old,
amrex::MultiFab &  W_old,
const amrex::MultiFab &  mf_Nturb,
const amrex::MultiFab &  mf_SMark,
const double &  time 
)
overridevirtual

Implements NullWindFarm.

20 {
21  AMREX_ALWAYS_ASSERT(W_old.nComp() > 0);
22  AMREX_ALWAYS_ASSERT(mf_Nturb.nComp() > 0);
23  AMREX_ALWAYS_ASSERT(mf_vars_generalAD.nComp() > 0);
24  AMREX_ALWAYS_ASSERT(time > -one);
25  compute_freestream_velocity(cons_in, U_old, V_old, mf_SMark);
26  source_terms_cellcentered(geom, cons_in, mf_SMark, mf_vars_generalAD);
27  update(dt_advance, cons_in, U_old, V_old, W_old, mf_vars_generalAD);
29 }
AMREX_ALWAYS_ASSERT(bx.length()[2]==khi+1)
constexpr amrex::Real one
Definition: ERF_NumericalConstants.H:30
void update(const double &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &U_old, amrex::MultiFab &V_old, amrex::MultiFab &W_old, const amrex::MultiFab &mf_vars)
Definition: ERF_AdvanceGeneralAD.cpp:61
void source_terms_cellcentered(const amrex::Geometry &geom, const amrex::MultiFab &cons_in, const amrex::MultiFab &mf_Smark, amrex::MultiFab &mf_vars_generalAD)
Definition: ERF_AdvanceGeneralAD.cpp:323
void compute_freestream_velocity(const amrex::MultiFab &cons_in, const amrex::MultiFab &U_old, const amrex::MultiFab &V_old, const amrex::MultiFab &mf_SMark)
Definition: ERF_AdvanceGeneralAD.cpp:97
void compute_power_output(const double &time)
Definition: ERF_AdvanceGeneralAD.cpp:32
Here is the call graph for this function:

◆ compute_freestream_velocity()

void GeneralAD::compute_freestream_velocity ( const amrex::MultiFab &  cons_in,
const amrex::MultiFab &  U_old,
const amrex::MultiFab &  V_old,
const amrex::MultiFab &  mf_SMark 
)
101 {
103  freestream_velocity.clear();
104  freestream_phi.clear();
105  disk_cell_count.clear();
106  freestream_velocity.resize(xloc.size(),zero);
107  freestream_phi.resize(xloc.size(),zero);
108  disk_cell_count.resize(xloc.size(),zero);
109 
110  Gpu::DeviceVector<Real> d_freestream_velocity(xloc.size());
111  Gpu::DeviceVector<Real> d_freestream_phi(yloc.size());
112  Gpu::DeviceVector<Real> d_disk_cell_count(yloc.size());
113  Gpu::copy(Gpu::hostToDevice, freestream_velocity.begin(), freestream_velocity.end(), d_freestream_velocity.begin());
114  Gpu::copy(Gpu::hostToDevice, freestream_phi.begin(), freestream_phi.end(), d_freestream_phi.begin());
115  Gpu::copy(Gpu::hostToDevice, disk_cell_count.begin(), disk_cell_count.end(), d_disk_cell_count.begin());
116 
117  Real* d_freestream_velocity_ptr = d_freestream_velocity.data();
118  Real* d_freestream_phi_ptr = d_freestream_phi.data();
119  Real* d_disk_cell_count_ptr = d_disk_cell_count.data();
120 
121 
122  for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
123 
124  auto SMark_array = mf_SMark.array(mfi);
125  auto u_vel = U_old.array(mfi);
126  auto v_vel = V_old.array(mfi);
127  // NOTE: this reduction is driven by the cell-centered SMark, so it must
128  // run over the cell-centered tilebox. Using nodaltilebox(0) here
129  // would include the plane at bigEnd(0)+1, which is a ghost cell of
130  // this box and a valid cell of its x-neighbor, and so would count
131  // that plane twice in the sums below.
132  Box tbx = mfi.tilebox();
133 
134  ParallelFor(tbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
135 
136  if(SMark_array(i,j,k,0) != -one) {
137  int turb_index = static_cast<int>(SMark_array(i,j,k,0));
138  Real phi = std::atan2(v_vel(i,j,k),u_vel(i,j,k)); // Wind direction w.r.t the x-direction
139  Gpu::Atomic::Add(&d_freestream_velocity_ptr[turb_index],std::pow(u_vel(i,j,k)*u_vel(i,j,k) + v_vel(i,j,k)*v_vel(i,j,k),myhalf));
140  Gpu::Atomic::Add(&d_disk_cell_count_ptr[turb_index],one);
141  Gpu::Atomic::Add(&d_freestream_phi_ptr[turb_index],phi);
142  }
143  });
144  }
145 
146  // Copy back to host
147  Gpu::copy(Gpu::deviceToHost, d_freestream_velocity.begin(), d_freestream_velocity.end(), freestream_velocity.begin());
148  Gpu::copy(Gpu::deviceToHost, d_freestream_phi.begin(), d_freestream_phi.end(), freestream_phi.begin());
149  Gpu::copy(Gpu::deviceToHost, d_disk_cell_count.begin(), d_disk_cell_count.end(), disk_cell_count.begin());
150 
151  // Reduce the data on every processor
152  amrex::ParallelAllReduce::Sum(freestream_velocity.data(),
153  static_cast<int>(freestream_velocity.size()),
154  amrex::ParallelContext::CommunicatorAll());
155 
156  amrex::ParallelAllReduce::Sum(freestream_phi.data(),
157  static_cast<int>(freestream_phi.size()),
158  amrex::ParallelContext::CommunicatorAll());
159 
160 
161  amrex::ParallelAllReduce::Sum(disk_cell_count.data(),
162  static_cast<int>(disk_cell_count.size()),
163  amrex::ParallelContext::CommunicatorAll());
164 
166 
167 
168  /*if (ParallelDescriptor::IOProcessor()){
169  for(int it=0; it<xloc.size(); it++){
170  std::cout << "turbine index, freestream velocity is " << it << " " << freestream_velocity[it] << " " <<
171  disk_cell_count[it] << " " <<
172  freestream_velocity[it]/(disk_cell_count[it] + 1e-10) << " " <<
173  freestream_phi[it]/(disk_cell_count[it] + 1e-10) << "\n";
174  }
175  }*/
176 }
ParallelFor(fab_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
constexpr amrex::Real zero
Definition: ERF_NumericalConstants.H:29
constexpr amrex::Real myhalf
Definition: ERF_NumericalConstants.H:34
amrex::Real Real
Definition: ERF_ShocInterface.H:19
amrex::Vector< amrex::Real > freestream_velocity
Definition: ERF_GeneralAD.H:51
amrex::Vector< amrex::Real > xloc
Definition: ERF_GeneralAD.H:47
amrex::Vector< amrex::Real > disk_cell_count
Definition: ERF_GeneralAD.H:51
amrex::Vector< amrex::Real > freestream_phi
Definition: ERF_GeneralAD.H:51
amrex::Vector< amrex::Real > yloc
Definition: ERF_GeneralAD.H:47
void get_turb_loc(amrex::Vector< amrex::Real > &xloc, amrex::Vector< amrex::Real > &yloc)
Definition: ERF_NullWindFarm.H:96
Here is the call graph for this function:

◆ compute_power_output()

void GeneralAD::compute_power_output ( const double &  time)
33 {
37 
38  const int n_spec_table = static_cast<int>(wind_speed.size());
39  // Compute power based on the look-up table
40 
41  if (ParallelDescriptor::IOProcessor()){
42  static std::ofstream file("power_output_GeneralAD.txt", std::ios::app);
43  // Check if the file opened successfully
44  if (!file.is_open()) {
45  std::cerr << "Error opening file!" << std::endl;
46  Abort("Could not open file to write power output in ERF_AdvanceSimpleAD.cpp");
47  }
48  Real total_power = zero;
49  for(int it=0; it<static_cast<int>(xloc.size()); it++){
50  Real avg_vel = freestream_velocity[it]/(disk_cell_count[it] + Real(1e-10));
51  Real turb_power = interpolate_1d(wind_speed.data(), power.data(), avg_vel, n_spec_table);
52  total_power = total_power + turb_power;
53  }
54  file << time << " " << total_power << "\n";
55  file.flush();
56  }
57 }
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real interpolate_1d(const amrex::Real *alpha, const amrex::Real *beta, const amrex::Real alpha_interp, const int alpha_size)
Definition: ERF_Interpolation_1D.H:14
amrex::Vector< amrex::Real > wind_speed
Definition: ERF_GeneralAD.H:50
amrex::Vector< amrex::Real > power
Definition: ERF_GeneralAD.H:50
amrex::Vector< amrex::Real > thrust_coeff
Definition: ERF_GeneralAD.H:50
amrex::Real rotor_rad
Definition: ERF_GeneralAD.H:49
amrex::Real hub_height
Definition: ERF_GeneralAD.H:49
amrex::Real thrust_coeff_standing
Definition: ERF_GeneralAD.H:49
void get_turb_spec(amrex::Real &rotor_rad, amrex::Real &hub_height, amrex::Real &thrust_coeff_standing, amrex::Vector< amrex::Real > &wind_speed, amrex::Vector< amrex::Real > &thrust_coeff, amrex::Vector< amrex::Real > &power)
Definition: ERF_NullWindFarm.H:84
Here is the call graph for this function:

◆ source_terms_cellcentered()

void GeneralAD::source_terms_cellcentered ( const amrex::Geometry &  geom,
const amrex::MultiFab &  cons_in,
const amrex::MultiFab &  mf_Smark,
amrex::MultiFab &  mf_vars_generalAD 
)
327 {
328 
330 
333 
335 
337 
339 
340  Real d_hub_height = hub_height;
341  Real d_rotor_rad = rotor_rad;
342 
343  Gpu::DeviceVector<Real> d_xloc(xloc.size());
344  Gpu::DeviceVector<Real> d_yloc(yloc.size());
345  Gpu::copy(Gpu::hostToDevice, xloc.begin(), xloc.end(), d_xloc.begin());
346  Gpu::copy(Gpu::hostToDevice, yloc.begin(), yloc.end(), d_yloc.begin());
347 
348  auto dx = geom.CellSizeArray();
349 
350  // Domain valid box
351  const amrex::Box& domain = geom.Domain();
352  auto ProbLoArr = geom.ProbLoArray();
353  int domlo_x = domain.smallEnd(0);
354  int domhi_x = domain.bigEnd(0) + 1;
355  int domlo_y = domain.smallEnd(1);
356  int domhi_y = domain.bigEnd(1) + 1;
357  int domlo_z = domain.smallEnd(2);
358  int domhi_z = domain.bigEnd(2) + 1;
359 
360  // The order of variables are - Vabs dVabsdt, dudt, dvdt, dTKEdt
361  mf_vars_generalAD.setVal(0.0);
362 
363  long unsigned int nturbs = static_cast<long unsigned int>(xloc.size());
364 
365  // This is the angle phi in Fig. 10 in Mirocha et. al. 2014
366  // set_turb_disk angle in ERF_InitWindFarm.cpp sets this phi as
367  // the turb_disk_angle
369  Real d_turb_disk_angle = turb_disk_angle;
370 
371  Gpu::DeviceVector<Real> d_freestream_velocity(nturbs);
372  Gpu::DeviceVector<Real> d_disk_cell_count(nturbs);
373  Gpu::copy(Gpu::hostToDevice, freestream_velocity.begin(), freestream_velocity.end(), d_freestream_velocity.begin());
374  Gpu::copy(Gpu::hostToDevice, disk_cell_count.begin(), disk_cell_count.end(), d_disk_cell_count.begin());
375 
376  Real* d_xloc_ptr = d_xloc.data();
377  Real* d_yloc_ptr = d_yloc.data();
378  Real* d_freestream_velocity_ptr = d_freestream_velocity.data();
379  Real* d_disk_cell_count_ptr = d_disk_cell_count.data();
380 
381  int n_bld_sections = static_cast<int>(bld_rad_loc.size());
382 
383  Gpu::DeviceVector<Real> d_bld_rad_loc(n_bld_sections);
384  Gpu::DeviceVector<Real> d_bld_twist(n_bld_sections);
385  Gpu::DeviceVector<Real> d_bld_chord(n_bld_sections);
386 
387  Gpu::copy(Gpu::hostToDevice, bld_rad_loc.begin(), bld_rad_loc.end(), d_bld_rad_loc.begin());
388  Gpu::copy(Gpu::hostToDevice, bld_twist.begin(), bld_twist.end(), d_bld_twist.begin());
389  Gpu::copy(Gpu::hostToDevice, bld_chord.begin(), bld_chord.end(), d_bld_chord.begin());
390 
391  Real* bld_rad_loc_ptr = d_bld_rad_loc.data();
392  Real* bld_twist_ptr = d_bld_twist.data();
393  Real* bld_chord_ptr = d_bld_chord.data();
394 
395  Vector<Gpu::DeviceVector<Real>> d_bld_airfoil_aoa(n_bld_sections);
396  Vector<Gpu::DeviceVector<Real>> d_bld_airfoil_Cl(n_bld_sections);
397  Vector<Gpu::DeviceVector<Real>> d_bld_airfoil_Cd(n_bld_sections);
398 
399  Vector<int> h_n_pts_airfoil(n_bld_sections);
400  Gpu::DeviceVector<int> d_n_pts_airfoil(n_bld_sections);
401 
402  for (int i = 0; i < n_bld_sections; ++i) {
403  h_n_pts_airfoil[i] = static_cast<int>(bld_airfoil_aoa[i].size());
404 
405  d_bld_airfoil_aoa[i].resize(bld_airfoil_aoa[i].size());
406  d_bld_airfoil_Cl[i].resize(bld_airfoil_Cl[i].size());
407  d_bld_airfoil_Cd[i].resize(bld_airfoil_Cd[i].size());
408 
409  Gpu::copy(Gpu::hostToDevice,
410  bld_airfoil_aoa[i].begin(), bld_airfoil_aoa[i].end(),
411  d_bld_airfoil_aoa[i].begin());
412 
413  Gpu::copy(Gpu::hostToDevice,
414  bld_airfoil_Cl[i].begin(), bld_airfoil_Cl[i].end(),
415  d_bld_airfoil_Cl[i].begin());
416 
417  Gpu::copy(Gpu::hostToDevice,
418  bld_airfoil_Cd[i].begin(), bld_airfoil_Cd[i].end(),
419  d_bld_airfoil_Cd[i].begin());
420  }
421 
422  Gpu::copy(Gpu::hostToDevice,
423  h_n_pts_airfoil.begin(), h_n_pts_airfoil.end(),
424  d_n_pts_airfoil.begin());
425 
426  int* d_n_pts_airfoil_ptr = d_n_pts_airfoil.data();
427 
428  Vector<Real*> hp_bld_airfoil_aoa, hp_bld_airfoil_Cl, hp_bld_airfoil_Cd;
429  for (auto & v :d_bld_airfoil_aoa) {
430  hp_bld_airfoil_aoa.push_back(v.data());
431  }
432  for (auto & v :d_bld_airfoil_Cl) {
433  hp_bld_airfoil_Cl.push_back(v.data());
434  }
435  for (auto & v :d_bld_airfoil_Cd) {
436  hp_bld_airfoil_Cd.push_back(v.data());
437  }
438 
439  Gpu::AsyncArray<Real*> aoa(hp_bld_airfoil_aoa.data(), n_bld_sections);
440  Gpu::AsyncArray<Real*> Cl(hp_bld_airfoil_Cl.data(), n_bld_sections);
441  Gpu::AsyncArray<Real*> Cd(hp_bld_airfoil_Cd.data(), n_bld_sections);
442 
443  auto d_bld_airfoil_aoa_ptr = aoa.data();
444  auto d_bld_airfoil_Cl_ptr = Cl.data();
445  auto d_bld_airfoil_Cd_ptr = Cd.data();
446 
447  int n_spec_extra = static_cast<int>(velocity.size());
448 
449  Gpu::DeviceVector<Real> d_velocity(n_spec_extra);
450  Gpu::DeviceVector<Real> d_rotor_RPM(n_spec_extra);
451  Gpu::DeviceVector<Real> d_blade_pitch(n_spec_extra);
452 
453  Gpu::copy(Gpu::hostToDevice, velocity.begin(), velocity.end(), d_velocity.begin());
454  Gpu::copy(Gpu::hostToDevice, rotor_RPM.begin(), rotor_RPM.end(), d_rotor_RPM.begin());
455  Gpu::copy(Gpu::hostToDevice, blade_pitch.begin(), blade_pitch.end(), d_blade_pitch.begin());
456 
457  auto d_velocity_ptr = d_velocity.data();
458  auto d_rotor_RPM_ptr = d_rotor_RPM.data();
459  auto d_blade_pitch_ptr = d_blade_pitch.data();
460 
461  for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
462 
463  const Box& gbx = mfi.growntilebox(1);
464  auto SMark_array = mf_SMark.array(mfi);
465  auto generalAD_array = mf_vars_generalAD.array(mfi);
466 
467  ParallelFor(gbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
468  int ii = amrex::min(amrex::max(i, domlo_x), domhi_x);
469  int jj = amrex::min(amrex::max(j, domlo_y), domhi_y);
470  int kk = amrex::min(amrex::max(k, domlo_z), domhi_z);
471 
472  Real x = ProbLoArr[0] + (ii+myhalf)*dx[0];
473  Real y = ProbLoArr[1] + (jj+myhalf)*dx[1];
474  Real z = ProbLoArr[2] + (kk+myhalf)*dx[2];
475  // ?? Density needed here
476 
477  int check_int = 0;
478 
479  Real source_x = zero, source_y = zero, source_z = zero;
480  std::array<Real,2> Fn_and_Ft;
481 
482  for(long unsigned int it=0;it<nturbs;it++) {
483  Real avg_vel = d_freestream_velocity_ptr[it]/(d_disk_cell_count_ptr[it] + Real(1e-10));
484  Real phi = d_turb_disk_angle;
485 
486  // This if check makes sure it is a point on the actuator disk
487  if(SMark_array(ii,jj,kk,1) == static_cast<double>(it)) {
488  check_int++;
489 
490  // Find radial distance of the point and the zeta angle
491  Real rad = std::pow( (x-d_xloc_ptr[it])*(x-d_xloc_ptr[it]) +
492  (y-d_yloc_ptr[it])*(y-d_yloc_ptr[it]) +
493  (z-d_hub_height)*(z-d_hub_height), myhalf );
494 
495  int index = find_rad_loc_index(rad, bld_rad_loc_ptr, n_bld_sections);
496 
497  // This if check makes sure it is a point with radial distance
498  // between the hub radius and the rotor radius.
499  // ?? hub radius needed here
500  if(rad >= two and rad <= d_rotor_rad) {
501  //AMREX_ASSERT( (z-d_hub_height) <= rad );
502  // Consider the vector that joines the point and the turbine center.
503  // Dot it on to the vector that joins the turbine center and along
504  // the plane of the disk. See fig. 10 in Mirocha et. al. Real(2014.)
505 
506  Real vec_proj = (x-d_xloc_ptr[it])*(std::sin(phi)) +
507  (y-d_yloc_ptr[it])*(-std::cos(phi));
508 
509 
510  Real zeta = std::atan2(z-d_hub_height, vec_proj);
511  //printf("zeta val is %0.15g\n", zeta*Real(180.0)/PI);
512  Fn_and_Ft = compute_source_terms_Fn_Ft(rad, avg_vel,
513  bld_rad_loc_ptr,
514  bld_twist_ptr,
515  bld_chord_ptr,
516  n_bld_sections,
517  d_bld_airfoil_aoa_ptr[index],
518  d_bld_airfoil_Cl_ptr[index],
519  d_bld_airfoil_Cd_ptr[index],
520  d_n_pts_airfoil_ptr[index],
521  d_velocity_ptr,
522  d_rotor_RPM_ptr,
523  d_blade_pitch_ptr,
524  n_spec_extra);
525 
526  Real Fn = three*Fn_and_Ft[0];
527  Real Ft = three*Fn_and_Ft[1];
528  // Compute the source terms - pass in radial distance, free stream velocity
529 
530  Real Fx = Fn*std::cos(phi) + Ft*std::sin(zeta)*std::sin(phi);
531  Real Fy = Fn*std::sin(phi) - Ft*std::sin(zeta)*std::cos(phi);
532  Real Fz = -Ft*std::cos(zeta);
533 
534  //Real dn = (
535 
536  source_x = -Fx/(two*PI*rad*dx[0])*one/std::pow(two*PI,myhalf);
537  source_y = -Fy/(two*PI*rad*dx[0])*one/std::pow(two*PI,myhalf);
538  source_z = -Fz/(two*PI*rad*dx[0])*one/std::pow(two*PI,myhalf);
539 
540 
541  //printf("Val source_x, is %0.15g, %0.15g, %0.15g %0.15g %0.15g %0.15g\n", rad, Fn, Ft, source_x, source_y, source_z);
542  }
543  }
544  }
545 
546  if(check_int > 1){
547  amrex::Error("Actuator disks are overlapping. Visualize actuator_disks.vtk "
548  "and check the windturbine locations input file. Exiting..");
549  }
550 
551  generalAD_array(i,j,k,0) = source_x;
552  generalAD_array(i,j,k,1) = source_y;
553  generalAD_array(i,j,k,2) = source_z;
554  });
555  }
556 }
AMREX_FORCE_INLINE AMREX_GPU_DEVICE int find_rad_loc_index(const Real rad, const Real *bld_rad_loc, const int n_bld_sections)
Definition: ERF_AdvanceGeneralAD.cpp:180
AMREX_FORCE_INLINE AMREX_GPU_DEVICE std::array< Real, 2 > compute_source_terms_Fn_Ft(const Real rad, const Real avg_vel, const Real *bld_rad_loc, const Real *bld_twist, const Real *bld_chord, int n_bld_sections, const Real *bld_airfoil_aoa, const Real *bld_airfoil_Cl, const Real *bld_airfoil_Cd, const int n_pts_airfoil, const Real *velocity, const Real *rotor_RPM, const Real *blade_pitch, const int n_spec_extra)
Definition: ERF_AdvanceGeneralAD.cpp:213
const Real dx
Definition: ERF_InitCustomPert_ABL.H:44
constexpr amrex::Real three
Definition: ERF_NumericalConstants.H:32
constexpr amrex::Real two
Definition: ERF_NumericalConstants.H:31
constexpr amrex::Real PI
Definition: ERF_NumericalConstants.H:39
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_Cl
Definition: ERF_GeneralAD.H:53
amrex::Vector< amrex::Real > velocity
Definition: ERF_GeneralAD.H:54
amrex::Vector< amrex::Real > blade_pitch
Definition: ERF_GeneralAD.H:54
amrex::Vector< amrex::Real > C_P
Definition: ERF_GeneralAD.H:54
amrex::Vector< amrex::Real > bld_twist
Definition: ERF_GeneralAD.H:52
amrex::Vector< amrex::Real > bld_chord
Definition: ERF_GeneralAD.H:52
amrex::Vector< amrex::Real > bld_rad_loc
Definition: ERF_GeneralAD.H:52
amrex::Real turb_disk_angle
Definition: ERF_GeneralAD.H:48
amrex::Vector< amrex::Real > C_T
Definition: ERF_GeneralAD.H:54
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_aoa
Definition: ERF_GeneralAD.H:53
amrex::Vector< amrex::Real > rotor_RPM
Definition: ERF_GeneralAD.H:54
amrex::Vector< amrex::Vector< amrex::Real > > bld_airfoil_Cd
Definition: ERF_GeneralAD.H:53
void get_turb_spec_extra(amrex::Vector< amrex::Real > &velocity, amrex::Vector< amrex::Real > &C_P, amrex::Vector< amrex::Real > &C_T, amrex::Vector< amrex::Real > &rotor_RPM, amrex::Vector< amrex::Real > &blade_pitch)
Definition: ERF_NullWindFarm.H:126
void get_blade_airfoil_spec(amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_aoa, amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cl, amrex::Vector< amrex::Vector< amrex::Real >> &bld_airfoil_Cd)
Definition: ERF_NullWindFarm.H:117
void get_turb_disk_angle(amrex::Real &turb_disk_angle)
Definition: ERF_NullWindFarm.H:103
void get_blade_spec(amrex::Vector< amrex::Real > &bld_rad_loc, amrex::Vector< amrex::Real > &bld_twist, amrex::Vector< amrex::Real > &bld_chord)
Definition: ERF_NullWindFarm.H:108
Here is the call graph for this function:

◆ update()

void GeneralAD::update ( const double &  dt_advance,
amrex::MultiFab &  cons_in,
amrex::MultiFab &  U_old,
amrex::MultiFab &  V_old,
amrex::MultiFab &  W_old,
const amrex::MultiFab &  mf_vars 
)
67 {
68 
69  for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
70 
71  Box tbx = mfi.nodaltilebox(0);
72  Box tby = mfi.nodaltilebox(1);
73  Box tbz = mfi.nodaltilebox(2);
74 
75  auto generalAD_array = mf_vars_generalAD.array(mfi);
76  auto u_vel = U_old.array(mfi);
77  auto v_vel = V_old.array(mfi);
78  auto w_vel = W_old.array(mfi);
79 
80  ParallelFor(tbx, tby, tbz,
81  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
82  {
83  u_vel(i,j,k) = u_vel(i,j,k) + (generalAD_array(i-1,j,k,0) + generalAD_array(i,j,k,0))/two*dt_advance;
84  },
85  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
86  {
87  v_vel(i,j,k) = v_vel(i,j,k) + (generalAD_array(i,j-1,k,1) + generalAD_array(i,j,k,1))/two*dt_advance;
88  },
89  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
90  {
91  w_vel(i,j,k) = w_vel(i,j,k) + (generalAD_array(i,j,k-1,2) + generalAD_array(i,j,k,2))/two*dt_advance;
92  });
93 
94  }
95 }
Here is the call graph for this function:

Member Data Documentation

◆ blade_pitch

amrex::Vector<amrex::Real> GeneralAD::blade_pitch
protected

◆ bld_airfoil_aoa

amrex::Vector<amrex::Vector<amrex::Real> > GeneralAD::bld_airfoil_aoa
protected

◆ bld_airfoil_Cd

amrex::Vector<amrex::Vector<amrex::Real> > GeneralAD::bld_airfoil_Cd
protected

◆ bld_airfoil_Cl

amrex::Vector<amrex::Vector<amrex::Real> > GeneralAD::bld_airfoil_Cl
protected

◆ bld_chord

amrex::Vector<amrex::Real> GeneralAD::bld_chord
protected

◆ bld_rad_loc

amrex::Vector<amrex::Real> GeneralAD::bld_rad_loc
protected

◆ bld_twist

amrex::Vector<amrex::Real> GeneralAD::bld_twist
protected

◆ C_P

amrex::Vector<amrex::Real> GeneralAD::C_P
protected

◆ C_T

amrex::Vector<amrex::Real> GeneralAD::C_T
protected

◆ disk_cell_count

amrex::Vector<amrex::Real> GeneralAD::disk_cell_count
protected

◆ freestream_phi

amrex::Vector<amrex::Real> GeneralAD::freestream_phi
protected

◆ freestream_velocity

amrex::Vector<amrex::Real> GeneralAD::freestream_velocity
protected

◆ hub_height

amrex::Real GeneralAD::hub_height
protected

◆ nominal_power

amrex::Real GeneralAD::nominal_power
protected

◆ power

amrex::Vector<amrex::Real> GeneralAD::power
protected

◆ rotor_rad

amrex::Real GeneralAD::rotor_rad
protected

◆ rotor_RPM

amrex::Vector<amrex::Real> GeneralAD::rotor_RPM
protected

◆ thrust_coeff

amrex::Vector<amrex::Real> GeneralAD::thrust_coeff
protected

◆ thrust_coeff_standing

amrex::Real GeneralAD::thrust_coeff_standing
protected

◆ turb_disk_angle

amrex::Real GeneralAD::turb_disk_angle
protected

◆ velocity

amrex::Vector<amrex::Real> GeneralAD::velocity
protected

◆ wind_speed

amrex::Vector<amrex::Real> GeneralAD::wind_speed
protected

◆ xloc

amrex::Vector<amrex::Real> GeneralAD::xloc
protected

◆ yloc

amrex::Vector<amrex::Real> GeneralAD::yloc
protected

The documentation for this class was generated from the following files: