diff --git a/flow.h b/flow.h index 1074d0a..f505aa5 100644 --- a/flow.h +++ b/flow.h @@ -67,7 +67,7 @@ namespace stim { #ifdef __CUDACC__ // for sphere template - __global__ void inside_sphere(const stim::sphere *V, size_t num, size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { + __global__ void inside_sphere(const stim::sphere *V, size_t num, T Z,size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { unsigned ix = blockDim.x * blockIdx.x + threadIdx.x; unsigned iy = blockDim.y * blockIdx.y + threadIdx.y; @@ -78,7 +78,7 @@ namespace stim { stim::vec3 world_pixel; world_pixel[0] = (T)ix * S[1] - x; // translate origin to center of the network world_pixel[1] = (T)iy * S[2] - y; - world_pixel[2] = ((T)z - R[3] / 2) * S[3]; // ???center of box minus half width + world_pixel[2] = ((T)z - Z / 2.0f) * S[3]; // ???center of box minus half width float distance = FLT_MAX; float tmp_distance; @@ -97,7 +97,7 @@ namespace stim { // for cone template - __global__ void inside_cone(const stim::cone *E, size_t num, size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { + __global__ void inside_cone(const stim::cone *E, size_t num, T Z, size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { unsigned ix = blockDim.x * blockIdx.x + threadIdx.x; unsigned iy = blockDim.y * blockIdx.y + threadIdx.y; @@ -107,7 +107,7 @@ namespace stim { stim::vec3 world_pixel; world_pixel[0] = (T)ix * S[1] - x; world_pixel[1] = (T)iy * S[2] - y; - world_pixel[2] = ((T)z - R[3] / 2) * S[3]; + world_pixel[2] = ((T)z - Z / 2.0f) * S[3]; float distance = FLT_MAX; float tmp_distance; @@ -133,7 +133,7 @@ namespace stim { // for source bus template - __global__ void inside_cuboid(const stim::cuboid *B, size_t num, size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { + __global__ void inside_cuboid(const stim::cuboid *B, size_t num, T Z, size_t *R, T *S, unsigned char *ptr, int x, int y, int z) { unsigned ix = blockDim.x * blockIdx.x + threadIdx.x; unsigned iy = blockDim.y * blockIdx.y + threadIdx.y; @@ -143,7 +143,7 @@ namespace stim { stim::vec3 world_pixel; world_pixel[0] = (T)ix * S[1] - x; world_pixel[1] = (T)iy * S[2] - y; - world_pixel[2] = ((T)z - R[3] / 2) * S[3]; + world_pixel[2] = ((T)z - Z / 2.0f) * S[3]; for (unsigned i = 0; i < num; i++) { bool left_outside = false; // flag indicates point is outside the left bound @@ -1042,8 +1042,9 @@ namespace stim { void glSolidCuboid(GLint subdivision, bool manufacture = false, T length = 40.0f, T height = 10.0f) { T width; - stim::vec3 L = bb.A; // get the bottom left corner - stim::vec3 U = bb.B; // get the top right corner + stim::gl_aaboundingbox BB = (*this).boundingbox(); + stim::vec3 L = BB.A; // get the bottom left corner + stim::vec3 U = BB.B; // get the top right corner width = U[2] - L[2] + 10.0f; if (manufacture) @@ -1887,8 +1888,8 @@ namespace stim { // there are cases that the fragment can not satisfy the requirement for width if (width < times * radius || n == 0) { // check feasibility - ratio = 0.0f; // load original lengths - desire_l = tmp_d; + ratio = 1.0f; // load original lengths + desire_l = tmp_d; // loading back-ups origin_l = tmp_l; std::cout << "Warning: current ratio is not feasible, use full original line." << std::endl; @@ -1923,11 +1924,21 @@ namespace stim { height = (T)((desire_l - (1 + 2 * n) * origin_l) / std::pow(2 * n, 2)); } - // degenerated case, compromise + /// degenerated case, use original extending method if (n == 0) { - n = 1; - width = (T)(origin_l) / (2 * n); - height = (T)((desire_l - (1 + 2 * n) * origin_l) / std::pow(2 * n, 2)); + height = (T)((desire_l - origin_l) / 2); + // "up" + tmp_v = tmp_v + stim::vec3(0, coef_up * height, 0); + if (feeder == 1) + inlet[i].V.push_back(tmp_v); + else if (feeder == 0) + outlet[i].V.push_back(tmp_v); + // "left" + tmp_v = tmp_v + stim::vec3(-coef_left * origin_l, 0, 0); + if (feeder == 1) + inlet[i].V.push_back(tmp_v); + else if (feeder == 0) + outlet[i].V.push_back(tmp_v); } // cube-like structure construction @@ -2628,7 +2639,7 @@ namespace stim { /// make full-synthetic binary image stack // prepare for image stack - void preparation(T &Xl, T &Xr, T &Yt, T &Yb, T &Z, bool prototype = false, T length = 40.0f, T height = 10.0f, T radius = 5.0f) { + void preparation(T &Xl, T &Xr, T &Yt, T &Yb, T &Z, bool prototype = false, T length = 40.0f, T height = 10.0f, T radius = 5.0f, T scale = 1.0f) { T max_radius = 0.0f; T top = FLT_MIN; @@ -2658,13 +2669,13 @@ namespace stim { for (unsigned i = 0; i < num_edge; i++) { for (unsigned j = 0; j < E[i].size(); j++) { new_sphere.c = E[i][j]; - new_sphere.r = E[i].r(j); + new_sphere.r = E[i].r(j) * scale; A.push_back(new_sphere); if (j != E[i].size() - 1) { new_cone.c1 = E[i][j]; new_cone.c2 = E[i][j + 1]; - new_cone.r1 = E[i].r(j); - new_cone.r2 = E[i].r(j + 1); + new_cone.r1 = E[i].r(j) * scale; + new_cone.r2 = E[i].r(j + 1) * scale; B.push_back(new_cone); } } @@ -2675,14 +2686,14 @@ namespace stim { for (unsigned i = 0; i < inlet.size(); i++) { for (unsigned j = 1; j < inlet[i].V.size() - 1; j++) { new_sphere.c = inlet[i].V[j]; - new_sphere.r = inlet[i].r; + new_sphere.r = inlet[i].r * scale; A.push_back(new_sphere); } } for (unsigned i = 0; i < outlet.size(); i++) { for (unsigned j = 1; j < outlet[i].V.size() - 1; j++) { new_sphere.c = outlet[i].V[j]; - new_sphere.r = outlet[i].r; + new_sphere.r = outlet[i].r * scale; A.push_back(new_sphere); } } @@ -2691,8 +2702,8 @@ namespace stim { for (unsigned j = 0; j < inlet[i].V.size() - 1; j++) { new_cone.c1 = inlet[i].V[j]; new_cone.c2 = inlet[i].V[j + 1]; - new_cone.r1 = inlet[i].r; - new_cone.r2 = inlet[i].r; + new_cone.r1 = inlet[i].r * scale; + new_cone.r2 = inlet[i].r * scale; B.push_back(new_cone); } } @@ -2700,15 +2711,15 @@ namespace stim { for (unsigned j = 0; j < outlet[i].V.size() - 1; j++) { new_cone.c1 = outlet[i].V[j]; new_cone.c2 = outlet[i].V[j + 1]; - new_cone.r1 = outlet[i].r; - new_cone.r2 = outlet[i].r; + new_cone.r1 = outlet[i].r * scale; + new_cone.r2 = outlet[i].r * scale; B.push_back(new_cone); } } // find out the image stack size - Xl = main_feeder[0][0] - length / 2; // left bound x coordinate - Xr = main_feeder[1][0] + length / 2; // right bound x coordinate + Xl = main_feeder[0][0] - length / 2 - 2 * radius; // left bound x coordinate + Xr = main_feeder[1][0] + length / 2 + 2 * radius; // right bound x coordinate for (unsigned i = 0; i < A.size(); i++) { if (A[i].c[1] > top) @@ -2720,21 +2731,21 @@ namespace stim { bb.B[2] = A[i].c[2]; if (A[i].c[2] < bb.A[2]) bb.A[2] = A[i].c[2]; - if (A[i].r > max_radius) - max_radius = A[i].r; + if (A[i].r * scale > max_radius) + max_radius = A[i].r * scale; } Yt = top + 2 * radius; // top bound y coordinate Yb = bottom - 2 * radius; // bottom bound y coordinate - Z = (bb.B[2] - bb.A[2] + 2 * max_radius); // bounding box width(along z-axis) + Z = (bb.B[2] - bb.A[2] + 4 * radius); // bounding box width(along z-axis) } /// making image stack main function - void make_image_stack(stim::image_stack &I, T dx, T dy, T dz, std::string stackdir, bool prototype = false, T radius = 5.0f) { + void make_image_stack(stim::image_stack &I, T dx, T dy, T dz, std::string stackdir, bool prototype = false, T radius = 5.0f, T scale = 1.0f) { /// preparation for making image stack T X, Xl, Xr, Y, Yt, Yb, Z; - preparation(Xl, Xr, Yt, Yb, Z, prototype); + preparation(Xl, Xr, Yt, Yb, Z, prototype, 40.0f, 10.0f, radius, scale); X = Xr - Xl; // bounding box length(along x-axis) Y = Yt - Yb; // bounding box height(along y-axis) stim::vec3 center = bb.center(); // get the center of bounding box @@ -2744,7 +2755,7 @@ namespace stim { /// make size_x = (int)(X / dx + 1); // set the size of image size_y = (int)(Y / dy + 1); - size_z = (int)(Z / dz + 1); + size_z = (int)(Z / dz + 1); // +3 in order to deal with reminder /// initialize image stack object I.init(1, size_x, size_y, size_z); I.set_dim(dx, dy, dz); @@ -2811,11 +2822,11 @@ namespace stim { dim3 block((unsigned)(size_x / max_thread + 1), (unsigned)(size_y / max_thread + 1)); dim3 thread((unsigned)max_thread, (unsigned)max_thread); - inside_sphere << > > (d_V, A.size(), d_R, d_S, d_ptr, x, y, z); + inside_sphere << > > (d_V, A.size(), Z, d_R, d_S, d_ptr, x, y, z); cudaDeviceSynchronize(); - inside_cone << > > (d_E, B.size(), d_R, d_S, d_ptr, x, y, z); + inside_cone << > > (d_E, B.size(), Z, d_R, d_S, d_ptr, x, y, z); cudaDeviceSynchronize(); - inside_cuboid << > > (d_B, CU.size(), d_R, d_S, d_ptr, x, y, z); + inside_cuboid << > > (d_B, CU.size(), Z, d_R, d_S, d_ptr, x, y, z); HANDLE_ERROR(cudaMemcpy(ptr, d_ptr, num * sizeof(unsigned char), cudaMemcpyDeviceToHost)); diff --git a/main.cu b/main.cu index 79722c0..e0af866 100644 --- a/main.cu +++ b/main.cu @@ -833,7 +833,7 @@ void glut_wheel(int wheel, int direction, int x, int y) { GLdouble posX, posY, posZ; window_to_world(posX, posY, posZ); // get the world coordinates - if (!to_select_pressure && (simulation || build_inlet_outlet)) { // check current pixel position only in simualtion and build_inlet_outlet modes + if (!to_select_pressure && (simulation || build_inlet_outlet || manufacture)) { // check current pixel position only in simualtion and build_inlet_outlet modes bool flag = flow.epsilon_vertex((float)posX, (float)posY, (float)posZ, eps, scale, pressure_index); if (flag && simulation && !glyph_mode) { float tmp_r; @@ -947,7 +947,7 @@ void glut_keyboard(unsigned char key, int x, int y) { case 'm': if (manufacture) { #ifdef __CUDACC__ - flow.make_image_stack(S, dx, dy, dz, stackdir, image_stack); + flow.make_image_stack(S, dx, dy, dz, stackdir, image_stack, default_radius, scale); #else std::cout << "You need to have a gpu to make image stack, sorry." << std::endl; -- libgit2 0.21.4