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

#include <ERF_Fitch.H>

Inheritance diagram for Fitch:
Collaboration diagram for Fitch:

Public Member Functions

 Fitch ()
 
virtual ~Fitch ()=default
 
void advance (const amrex::Geometry &geom, const amrex::Real &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &mf_vars_fitch, amrex::MultiFab &U_old, amrex::MultiFab &V_old, amrex::MultiFab &W_old, const amrex::MultiFab &mf_Nturb, const amrex::MultiFab &mf_SMark, const amrex::Real &time) override
 
void source_terms_cellcentered (const amrex::Geometry &geom, const amrex::MultiFab &cons_in, amrex::MultiFab &mf_vars_ewp, const amrex::MultiFab &U_old, const amrex::MultiFab &V_old, const amrex::MultiFab &W_old, const amrex::MultiFab &mf_Nturb)
 
void update (const amrex::Real &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &U_old, amrex::MultiFab &V_old, const amrex::MultiFab &mf_vars_fitch)
 
- 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::Real > xloc
 
amrex::Vector< amrex::Real > yloc
 
amrex::Real hub_height
 
amrex::Real rotor_rad
 
amrex::Real thrust_coeff_standing
 
amrex::Real nominal_power
 
amrex::Vector< amrex::Real > wind_speed
 
amrex::Vector< amrex::Real > thrust_coeff
 
amrex::Vector< amrex::Real > power
 
- Protected Attributes inherited from NullWindFarm
amrex::Vector< amrex::Real > m_xloc
 
amrex::Vector< amrex::Real > m_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::Real > m_wind_speed
 
amrex::Vector< amrex::Real > m_thrust_coeff
 
amrex::Vector< amrex::Real > m_power
 
amrex::Vector< amrex::Real > m_bld_rad_loc
 
amrex::Vector< amrex::Real > m_bld_twist
 
amrex::Vector< amrex::Real > m_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::Real > m_velocity
 
amrex::Vector< amrex::Real > m_C_P
 
amrex::Vector< amrex::Real > m_C_T
 
amrex::Vector< amrex::Real > m_rotor_RPM
 
amrex::Vector< amrex::Real > m_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

◆ Fitch()

Fitch::Fitch ( )
inline
13 {}

◆ ~Fitch()

virtual Fitch::~Fitch ( )
virtualdefault

Member Function Documentation

◆ advance()

void Fitch::advance ( const amrex::Geometry &  geom,
const amrex::Real &  dt_advance,
amrex::MultiFab &  cons_in,
amrex::MultiFab &  mf_vars_fitch,
amrex::MultiFab &  U_old,
amrex::MultiFab &  V_old,
amrex::MultiFab &  W_old,
const amrex::MultiFab &  mf_Nturb,
const amrex::MultiFab &  mf_SMark,
const amrex::Real &  time 
)
overridevirtual

Implements NullWindFarm.

58 {
59  AMREX_ALWAYS_ASSERT(W_old.nComp() > 0);
60  AMREX_ALWAYS_ASSERT(mf_SMark.nComp() > 0);
61  AMREX_ALWAYS_ASSERT(time > -1.0);
62  source_terms_cellcentered(geom, cons_in, mf_vars_fitch, U_old, V_old, W_old, mf_Nturb);
63  update(dt_advance, cons_in, U_old, V_old, mf_vars_fitch);
64 }
void update(const amrex::Real &dt_advance, amrex::MultiFab &cons_in, amrex::MultiFab &U_old, amrex::MultiFab &V_old, const amrex::MultiFab &mf_vars_fitch)
Definition: ERF_AdvanceFitch.cpp:68
void source_terms_cellcentered(const amrex::Geometry &geom, const amrex::MultiFab &cons_in, amrex::MultiFab &mf_vars_ewp, const amrex::MultiFab &U_old, const amrex::MultiFab &V_old, const amrex::MultiFab &W_old, const amrex::MultiFab &mf_Nturb)
Definition: ERF_AdvanceFitch.cpp:102

◆ source_terms_cellcentered()

void Fitch::source_terms_cellcentered ( const amrex::Geometry &  geom,
const amrex::MultiFab &  cons_in,
amrex::MultiFab &  mf_vars_ewp,
const amrex::MultiFab &  U_old,
const amrex::MultiFab &  V_old,
const amrex::MultiFab &  W_old,
const amrex::MultiFab &  mf_Nturb 
)
109 {
110 
113 
114  auto dx = geom.CellSizeArray();
115 
116  // Domain valid box
117  const amrex::Box& domain = geom.Domain();
118  int domlo_z = domain.smallEnd(2);
119  int domhi_z = domain.bigEnd(2) + 1;
120 
121 
122  //Real sum = 0.0;
123  //Real *sum_area = &sum;
124 
125  // The order of variables are - Vabs dVabsdt, dudt, dvdt, dTKEdt
126  mf_vars_fitch.setVal(0.0);
127  Real d_hub_height = hub_height;
128  Real d_rotor_rad = rotor_rad;
129  Gpu::DeviceVector<Real> d_wind_speed(wind_speed.size());
130  Gpu::DeviceVector<Real> d_thrust_coeff(thrust_coeff.size());
131 
132  // Copy data from host vectors to device vectors
133  Gpu::copy(Gpu::hostToDevice, wind_speed.begin(), wind_speed.end(), d_wind_speed.begin());
134  Gpu::copy(Gpu::hostToDevice, thrust_coeff.begin(), thrust_coeff.end(), d_thrust_coeff.begin());
135 
136  for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
137 
138  const Box& gbx = mfi.growntilebox(1);
139  auto fitch_array = mf_vars_fitch.array(mfi);
140  auto Nturb_array = mf_Nturb.array(mfi);
141  auto u_vel = U_old.array(mfi);
142  auto v_vel = V_old.array(mfi);
143  auto w_vel = W_old.array(mfi);
144 
145  const Real* wind_speed_d = d_wind_speed.dataPtr();
146  const Real* thrust_coeff_d = d_thrust_coeff.dataPtr();
147  const int n_spec_table = d_wind_speed.size();
148 
149  ParallelFor(gbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
150  int kk = amrex::min(amrex::max(k, domlo_z), domhi_z);
151 
152  Real z_k = kk*dx[2];
153  Real z_kp1 = (kk+1)*dx[2];
154 
155  Real A_ijk = compute_Aijk(z_k, z_kp1, d_hub_height, d_rotor_rad);
156 
157  // Compute Fitch source terms
158 
159  Real Vabs = std::pow(u_vel(i,j,k)*u_vel(i,j,k) +
160  v_vel(i,j,k)*v_vel(i,j,k) +
161  w_vel(i,j,kk)*w_vel(i,j,kk), 0.5);
162 
163  Real C_T = interpolate_1d(wind_speed_d, thrust_coeff_d, Vabs, n_spec_table);
164  Real C_TKE = 0.0;
165 
166  fitch_array(i,j,k,0) = Vabs;
167  fitch_array(i,j,k,1) = -0.5*Nturb_array(i,j,k)/(dx[0]*dx[1])*C_T*Vabs*Vabs*A_ijk/(z_kp1 - z_k);
168  fitch_array(i,j,k,2) = u_vel(i,j,k)/Vabs*fitch_array(i,j,k,1);
169  fitch_array(i,j,k,3) = v_vel(i,j,k)/Vabs*fitch_array(i,j,k,1);
170  fitch_array(i,j,k,4) = 0.5*Nturb_array(i,j,k)/(dx[0]*dx[1])*C_TKE*std::pow(Vabs,3)*A_ijk/(z_kp1 - z_k);
171 
172  //amrex::Gpu::Atomic::Add(sum_area, A_ijk);
173  });
174  }
175  //std::cout << "Checking sum here...." <<"\n";
176  //printf("%0.15g, %0.15g\n", *sum_area , PI*R*R);
177  //exit(0);
178 }
AMREX_FORCE_INLINE AMREX_GPU_DEVICE Real compute_Aijk(const Real z_k, const Real z_kp1, const Real hub_height, const Real rotor_rad)
Definition: ERF_AdvanceFitch.cpp:25
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:12
amrex::Vector< amrex::Real > thrust_coeff
Definition: ERF_Fitch.H:44
amrex::Vector< amrex::Real > power
Definition: ERF_Fitch.H:44
amrex::Real rotor_rad
Definition: ERF_Fitch.H:43
amrex::Real thrust_coeff_standing
Definition: ERF_Fitch.H:43
amrex::Real hub_height
Definition: ERF_Fitch.H:43
amrex::Vector< amrex::Real > wind_speed
Definition: ERF_Fitch.H:44
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:

◆ update()

void Fitch::update ( const amrex::Real &  dt_advance,
amrex::MultiFab &  cons_in,
amrex::MultiFab &  U_old,
amrex::MultiFab &  V_old,
const amrex::MultiFab &  mf_vars_fitch 
)
72 {
73 
74  for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
75 
76  Box bx = mfi.tilebox();
77  Box tbx = mfi.nodaltilebox(0);
78  Box tby = mfi.nodaltilebox(1);
79 
80  auto cons_array = cons_in.array(mfi);
81  auto fitch_array = mf_vars_fitch.array(mfi);
82  auto u_vel = U_old.array(mfi);
83  auto v_vel = V_old.array(mfi);
84 
85  ParallelFor(tbx, tby, bx,
86  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
87  {
88  u_vel(i,j,k) = u_vel(i,j,k) + (fitch_array(i-1,j,k,2) + fitch_array(i,j,k,2))/2.0*dt_advance;
89  },
90  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
91  {
92  v_vel(i,j,k) = v_vel(i,j,k) + (fitch_array(i,j-1,k,3) + fitch_array(i,j,k,3))/2.0*dt_advance;
93  },
94  [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
95  {
96  cons_array(i,j,k,RhoKE_comp) = cons_array(i,j,k,RhoKE_comp) + fitch_array(i,j,k,4)*dt_advance;
97  });
98  }
99 }
#define RhoKE_comp
Definition: ERF_IndexDefines.H:38

Member Data Documentation

◆ hub_height

amrex::Real Fitch::hub_height
protected

◆ nominal_power

amrex::Real Fitch::nominal_power
protected

◆ power

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

◆ rotor_rad

amrex::Real Fitch::rotor_rad
protected

◆ thrust_coeff

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

◆ thrust_coeff_standing

amrex::Real Fitch::thrust_coeff_standing
protected

◆ wind_speed

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

◆ xloc

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

◆ yloc

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

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