170 Gpu::DeviceVector<Real> d_xloc(
xloc.size());
171 Gpu::DeviceVector<Real> d_yloc(
yloc.size());
172 Gpu::copy(Gpu::hostToDevice,
xloc.begin(),
xloc.end(), d_xloc.begin());
173 Gpu::copy(Gpu::hostToDevice,
yloc.begin(),
yloc.end(), d_yloc.begin());
175 auto dx = geom.CellSizeArray();
178 const amrex::Box& domain = geom.Domain();
179 int domlo_x = domain.smallEnd(0);
180 int domhi_x = domain.bigEnd(0) + 1;
181 int domlo_y = domain.smallEnd(1);
182 int domhi_y = domain.bigEnd(1) + 1;
183 int domlo_z = domain.smallEnd(2);
184 int domhi_z = domain.bigEnd(2) + 1;
187 mf_vars_simpleAD.setVal(0.0);
189 long unsigned int nturbs =
xloc.size();
191 Gpu::DeviceVector<Real> d_freestream_velocity(nturbs);
192 Gpu::DeviceVector<Real> d_freestream_phi(nturbs);
193 Gpu::DeviceVector<Real> d_disk_cell_count(nturbs);
198 Real* d_freestream_velocity_ptr = d_freestream_velocity.data();
199 Real* d_freestream_phi_ptr = d_freestream_phi.data();
200 Real* d_disk_cell_count_ptr = d_disk_cell_count.data();
207 Gpu::DeviceVector<Real> d_wind_speed(
wind_speed.size());
208 Gpu::DeviceVector<Real> d_thrust_coeff(
thrust_coeff.size());
214 const Real* wind_speed_d = d_wind_speed.dataPtr();
215 const Real* thrust_coeff_d = d_thrust_coeff.dataPtr();
216 const int n_spec_table = d_wind_speed.size();
218 for ( MFIter mfi(cons_in,TilingIfNotGPU()); mfi.isValid(); ++mfi) {
220 const Box& gbx = mfi.growntilebox(1);
221 auto SMark_array = mf_SMark.array(mfi);
222 auto simpleAD_array = mf_vars_simpleAD.array(mfi);
224 ParallelFor(gbx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k) noexcept {
225 int ii = amrex::min(amrex::max(i, domlo_x), domhi_x);
226 int jj = amrex::min(amrex::max(j, domlo_y), domhi_y);
227 int kk = amrex::min(amrex::max(k, domlo_z), domhi_z);
233 int it =
static_cast<int>(SMark_array(ii,jj,kk,1));
236 Real avg_vel = d_freestream_velocity_ptr[it]/(d_disk_cell_count_ptr[it] + 1e-10);
237 Real phi = d_freestream_phi_ptr[it]/(d_disk_cell_count_ptr[it] + 1e-10);
239 Real C_T =
interpolate_1d(wind_speed_d, thrust_coeff_d, avg_vel, n_spec_table);
242 a = 0.5 - 0.5*std::pow(1.0-C_T,0.5);
244 Real Uinfty_dot_nhat = avg_vel*(std::cos(phi)*nx + std::sin(phi)*ny);
246 source_x = -2.0*std::pow(Uinfty_dot_nhat, 2.0)*a*(1.0-a)*dx[1]*dx[2]*std::cos(d_turb_disk_angle)/(dx[0]*dx[1]*dx[2])*std::cos(phi);
247 source_y = -2.0*std::pow(Uinfty_dot_nhat, 2.0)*a*(1.0-a)*dx[1]*dx[2]*std::cos(d_turb_disk_angle)/(dx[0]*dx[1]*dx[2])*std::sin(phi);
250 source_x = -0.5*C_T*std::pow(Uinfty_dot_nhat, 2.0)*dx[1]*dx[2]*std::cos(d_turb_disk_angle)/(dx[0]*dx[1]*dx[2])*std::cos(phi);
251 source_y = -0.5*C_T*std::pow(Uinfty_dot_nhat, 2.0)*dx[1]*dx[2]*std::cos(d_turb_disk_angle)/(dx[0]*dx[1]*dx[2])*std::sin(phi);
255 simpleAD_array(i,j,k,0) = source_x;
256 simpleAD_array(i,j,k,1) = source_y;
void get_turb_disk_angle(amrex::Real &turb_disk_angle)
Definition: ERF_NullWindFarm.H:103
amrex::Real turb_disk_angle
Definition: ERF_SimpleAD.H:47