Closed wjsjtu123 closed 4 months ago
Hi @wjsjtu123,
some possible fixes:
lbm_u
from 0.16f
to 0.1f
or lower - this is the blade tip velocity in LBM units, and it should not be much larger than 0.1f
SUBGRID
is enabledGRAPHICS_U_MAX
to a larger valueKind regards, Moritz
Thank you very much the above problem has been solved and is due to the fact that my propeller tip linear velocity is well above the rendered velocity maximum of 0.25. but I have a new problem.
When performing a propeller rotational motion, I define a rotational motion around the x-axis with a rotational speed of 4145r/min, but in the voxelized motion, instead of following a fixed x-axis, the propeller will rotate periodically off-center, what is the cause of this, and how can I fix it? My code is as follows
void main_setup() { // NACA X-57 lift propeller; required extensions in defines.hpp: FP16C, EQUILIBRIUM_BOUNDARIES, MOVING_BOUNDARIES, SUBGRID, INTERACTIVE_GRAPHICS or GRAPHICS
// ################################################################## define simulation box size, viscosity and volume force ###################################################################
const uint3 lbm_N = resolution(float3(2.0f, 1.0f, 1.0f), 4000u); // input: simulation box aspect ratio and VRAM occupation in MB, output: grid resolution
const float lbm_u = 0.11f;
const float lbm_length = 0.5f*(float)lbm_N.x;
const float si_n =4145.0f; //RPM of propeller
const float si_n_t = 20.0f;
const float si_T = (60.0f/si_n)si_n_t; // 20 revolutions of the propeller
const uint lbm_dt = 30u; // revoxelize rotor every dt time steps
const float si_length=0.984f;
const float si_J = 0.6f; // Advance ratio
const float si_D = 0.57592f; //propeller radius
const float si_u = si_Jsi_D*(si_n/60.0f); //airspeed
const float si_w_u = (si_npif/30.0f)(si_D0.5f); //propeller linear velocity
const float si_nu=1.48E-5f, si_rho=1.225f;
units.set_m_kg_s(lbm_length, lbm_u, 1.0f, si_length, si_u, si_rho);
print_info("lbm_u="+to_string(units.u(si_u)));
LBM lbm(lbm_N, 1u, 1u, 1u, units.nu(si_nu));
// ###################################################################################### define geometry ######################################################################################
Mesh fuselage = read_stl(get_exe_path()+"../stl/propeller_fuselage_semi.stl"); //
Mesh* propeller = read_stl(get_exe_path()+"../stl/propeller_semi.stl"); //
const float scale = lbm_length/fuselage->get_bounding_box_size().x; // scale body and rotors to simulation box size
fuselage->scale(scale);
propeller->scale(scale);
const float3 offset = lbm.center()-fuselage->get_bounding_box_center(); // move body and rotors to simulation box center
fuselage->translate(offset);
propeller->translate(offset);
fuselage->set_center(fuselage->get_bounding_box_center()); // set center of meshes to their bounding box center
propeller->set_center(propeller->get_bounding_box_center());
print_info("diameter="+to_string(units.x(si_D)));
print_info("propeller_max_size="+to_string(propeller->get_max_size()));
const float propeller_radius=0.5funits.x(si_D), propeller_omega=(units.si_u(si_w_u))/(propeller_radius), propeller_domega=propeller_omega(float)lbm_dt;
lbm.voxelize_mesh_on_device(fuselage);
//lbm.voxelize_mesh_on_device(propeller, TYPE_S, propeller->get_center(), float3(0.0f), float3(propeller_omega, 0.0f, 0.0f)); // revoxelize mesh on GPU
const uint Nx=lbm.get_Nx(), Ny=lbm.get_Ny(), Nz=lbm.get_Nz(); parallel_for(lbm.get_N(), [&](ulong n) { uint x=0u, y=0u, z=0u; lbm.coordinates(n, x, y, z);
if(lbm.flags[n]!=TYPE_S) lbm.u.x[n] = lbm_u;
if(x==0u||x==Nx-1u||y==0u||y==Ny-1u||z==0u||z==Nz-1u) lbm.flags[n] = TYPE_E; // all non periodic
}); // ####################################################################### run simulation, export images and data ##########################################################################
lbm.graphics.visualization_modes = VIS_FLAG_SURFACE|VIS_Q_CRITERION;
lbm.run(0u); // initialize simulation
while(lbm.get_t()<=units.t(si_T)) { // main simulation loop
lbm.voxelize_mesh_on_device(propeller, TYPE_S, propeller->get_center(), float3(0.0f), float3(propeller_omega, 0.0f, 0.0f)); // revoxelize mesh on GPU
lbm.run(lbm_dt); // run dt time steps
propeller->rotate(float3x3(float3(1.0f, 0.0f, 0.0f), propeller_domega)); // rotate mesh
#if defined(GRAPHICS) && !defined(INTERACTIVE_GRAPHICS)
if(lbm.graphics.next_frame(units.t(si_T), 10.0f)) {
lbm.graphics.set_camera_free(float3(0.528513f*(float)Nx, 0.102095f*(float)Ny, 1.302283f*(float)Nz), 16.0f, 47.0f, 96.0f);
lbm.graphics.write_frame(get_exe_path()+"export/a/");
lbm.graphics.set_camera_free(float3(0.0f*(float)Nx, -0.114244f*(float)Ny, 0.543265f*(float)Nz), 90.0f+degrees((float)lbm.get_t()/(float)lbm_dtmain_domega), 36.0f, 120.0f);
lbm.graphics.write_frame(get_exe_path()+"export/b/");
lbm.graphics.set_camera_free(float3(0.557719f(float)Nx, -0.503388f*(float)Ny, -0.591976f*(float)Nz), -43.0f, -21.0f, 75.0f);
lbm.graphics.write_frame(get_exe_path()+"export/c/");
lbm.graphics.set_camera_centered(58.0f, 9.0f, 88.0f, 1.648722f);
lbm.graphics.write_frame(get_exe_path()+"export/d/");
lbm.graphics.set_camera_centered(0.0f, 90.0f, 100.0f, 1.100000f);
lbm.graphics.write_frame(get_exe_path()+"export/e/");
lbm.graphics.set_camera_free(float3(0.001612f*(float)Nx, 0.523852f*(float)Ny, 0.992613f*(float)Nz), 90.0f, 37.0f, 94.0f);
lbm.graphics.write_frame(get_exe_path()+"export/f/");
}
#endif // GRAPHICS && !INTERACTIVE_GRAPHICS
}
}
I set the propeller rotational speed to 4145r/min and then the resultant visualization process shows a lot of red noise areas. The code is shown below, I don't know if the code is wrong or not.