diff options
| -rw-r--r-- | src/physics.c | 185 | ||||
| -rw-r--r-- | src/physics.h | 7 | ||||
| -rw-r--r-- | src/ui.c | 109 |
3 files changed, 179 insertions, 122 deletions
diff --git a/src/physics.c b/src/physics.c index 49c6552..3981258 100644 --- a/src/physics.c +++ b/src/physics.c @@ -2,7 +2,8 @@ #include "itc.h" #include "physics.h" -void phy_sim_destroy(phy_sim *sim) { +void phy_sim_destroy(phy_sim *sim) +{ size_t i; SDL_DestroyMutex(sim->rotate_lock); @@ -15,7 +16,28 @@ void phy_sim_destroy(phy_sim *sim) { } } -int phy_sim_create(phy_sim *sim) { +void phy_sim_reset(phy_sim *sim) +{ + size_t i; + phy_field_info *fi; + + fi = &sim->field_info; + + for (i = 0; i < 3; i++) { + memset(sim->fields[i].E, 0, fi->size * sizeof(float)); + memset(sim->fields[i].H, 0, fi->size * sizeof(float)); + } + + memset(sim->field_eps, 0, fi->size1 * sizeof(float)); + memset(sim->field_eps, 0, fi->size1 * sizeof(float)); + phy_sim_compute_const_fields(sim); + + sim->time = 0; + sim->running = false; +} + +int phy_sim_create(phy_sim *sim) +{ size_t i; phy_field_info *fi; @@ -23,10 +45,11 @@ int phy_sim_create(phy_sim *sim) { fi = &sim->field_info; - fi->width = 50; - fi->height = 50; + fi->width = 100; + fi->height = 100; fi->depth = 50; - fi->spacing = 1.5e9 / max3(fi->width, fi->height, fi->depth); + fi->spacing = 1.5e10f / max3(fi->width, fi->height, fi->depth); + sim->time_delta = 0.1; fi->xstr = 3; fi->ystr = fi->xstr * fi->width; @@ -39,25 +62,28 @@ int phy_sim_create(phy_sim *sim) { fi->size1 = fi->zstr1 * fi->depth; for (i = 0; i < 3; i++) { - sim->fields[i].E = calloc(fi->size, sizeof(float)); - sim->fields[i].H = calloc(fi->size, sizeof(float)); + sim->fields[i].E = malloc(fi->size * sizeof(float)); + sim->fields[i].H = malloc(fi->size * sizeof(float)); if (!sim->fields[i].E || !sim->fields[i].H) { con_printf("phy_sim_create: out of memory\n"); goto error; } } - sim->field_eps = calloc(fi->size1, sizeof(float)); - sim->field_mu = calloc(fi->size1, sizeof(float)); + sim->field_eps = malloc(fi->size1 * sizeof(float)); + sim->field_mu = malloc(fi->size1 * sizeof(float)); + if (!sim->field_eps || !sim->field_mu) { + con_printf("phy_sim_create: out of memory\n"); + goto error; + } - phy_sim_compute_const_fields(sim); + phy_sim_reset(sim); itc_chan_create(&sim->ctl); sim->running = false; sim->rotate_lock = SDL_CreateMutex(); assert(sim->rotate_lock); - return 0; error: @@ -85,104 +111,89 @@ void phy_sim_compute_const_fields(phy_sim *sim) } } -void phy_sim_step_1(phy_field_info *fi, float dt, float *mul_const, - const float *L0, float *L2, const float *R) +void phy_sim_step_E(phy_field_info *fi, float *mul_const, float dt, + float *En, float *Ep, float *H) { size_t x, y, z; + float dx = fi->spacing; -#define AT(F, x, y, z) ((F) + (z) * fi->zstr + (y) * fi->ystr + (x) * fi->xstr) -#define AT1(F, x, y, z) ((F)[(z) * fi->zstr1 + (y) * fi->ystr1 + \ - (x) * fi->xstr1]) - - // X for (z = 1; z < fi->depth - 1; z++) for (y = 1; y < fi->height - 1; y++) - for (x = 0; x < fi->width; x++) { - float dRzdy, dRydz, dLxdt; - - dRzdy = AT(R, x, y + 1, z)[2] - AT(R, x, y - 1, z)[2]; - dRzdy /= 2 * fi->spacing; - - dRydz = AT(R, x, y, z + 1)[1] - AT(R, x, y, z - 1)[1]; - dRydz /= 2 * fi->spacing; - - dLxdt = AT1(mul_const, x, y, z) * (dRzdy - dRydz); - - AT(L2, x, y, z)[0] = AT(L0, x, y, z)[0] + dLxdt * dt; - } - - // Y - for (z = 1; z < fi->depth - 1; z++) - for (y = 0; y < fi->height; y++) for (x = 1; x < fi->width - 1; x++) { - float dRxdz, dRzdx, dLydt; + size_t offs, offs_nx, offs_ny, offs_nz; + size_t offs1; + float dExdt, dEydt, dEzdt; - dRxdz = AT(R, x, y, z + 1)[0] - AT(R, x, y, z - 1)[0]; - dRxdz /= 2 * fi->spacing; + offs = x * fi->xstr + y * fi->ystr + z * fi->zstr; + offs_nx = offs - fi->xstr; + offs_ny = offs - fi->ystr; + offs_nz = offs - fi->zstr; - dRzdx = AT(R, x + 1, y, z)[2] - AT(R, x - 1, y, z)[2]; - dRzdx /= 2 * fi->spacing; + offs1 = x * fi->xstr1 + y * fi->ystr1 + z * fi->zstr1; - dLydt = AT1(mul_const, x, y, z) * (dRxdz - dRzdx); + dExdt = (H[offs+2] - H[offs_ny+2] - H[offs+1] + H[offs_nz+1]); + dEydt = (H[offs+0] - H[offs_nz+0] - H[offs+2] + H[offs_nx+2]); + dEzdt = (H[offs+1] - H[offs_nx+1] - H[offs+0] + H[offs_ny+0]); - AT(L2, x, y, z)[1] = AT(L0, x, y, z)[1] + dLydt * dt; + En[offs+0] = Ep[offs+0] + dt * mul_const[offs1] * dExdt / dx; + En[offs+1] = Ep[offs+1] + dt * mul_const[offs1] * dEydt / dx; + En[offs+2] = Ep[offs+2] + dt * mul_const[offs1] * dEzdt / dx; } +} - // Z - for (z = 0; z < fi->depth; z++) +void phy_sim_step_H(phy_field_info *fi, float *mul_const, float dt, + float *Hn, float *Hp, float *E) +{ + size_t x, y, z; + float dx = fi->spacing; + + for (z = 1; z < fi->depth - 1; z++) for (y = 1; y < fi->height - 1; y++) for (x = 1; x < fi->width - 1; x++) { - float dRydx, dRxdy, dLzdt; + size_t offs, offs_px, offs_py, offs_pz; + size_t offs1; + float dHxdt, dHydt, dHzdt; - dRydx = AT(R, x + 1, y, z)[1] - AT(R, x - 1, y, z)[1]; - dRydx /= 2 * fi->spacing; + offs = x * fi->xstr + y * fi->ystr + z * fi->zstr; + offs_px = offs + fi->xstr; + offs_py = offs + fi->ystr; + offs_pz = offs + fi->zstr; - dRxdy = AT(R, x, y + 1, z)[0] - AT(R, x, y - 1, z)[0]; - dRxdy /= 2 * fi->spacing; + offs1 = x * fi->xstr1 + y * fi->ystr1 + z * fi->zstr1; - dLzdt = AT1(mul_const, x, y, z) * (dRydx - dRxdy); + dHxdt = (E[offs_py+2] - E[offs+2] - E[offs_pz+1] + E[offs+1]); + dHydt = (E[offs_pz+0] - E[offs+0] - E[offs_px+2] + E[offs+2]); + dHzdt = (E[offs_px+1] - E[offs+1] - E[offs_py+0] + E[offs+0]); - AT(L2, x, y, z)[2] = AT(L0, x, y, z)[2] + dLzdt * dt; + Hn[offs+0] = Hp[offs+0] + dt * mul_const[offs1] * dHxdt / dx; + Hn[offs+1] = Hp[offs+1] + dt * mul_const[offs1] * dHydt / dx; + Hn[offs+2] = Hp[offs+2] + dt * mul_const[offs1] * dHzdt / dx; } -#undef AT -#undef AT1 } void phy_sim_add_sources(phy_sim *sim) { phy_field_info *fi = &sim->field_info; + float *E; size_t x, y, z; - for (z = 0; z < fi->depth; z++) - for (y = 0; y < fi->height; y++) - for (x = 0; x < fi->width; x++) { - float *E, *H, fx, fy, fz, chuj; - - fx = (float)x / fi->width * 2 - 1; - fy = (float)y / fi->height * 2 - 1; - fz = (float)z / fi->depth * 2 - 1; + for (z = 1; z < fi->depth - 1; z++) + for (y = 15; y < fi->height - 15; y++) + for (x = 1; x < fi->width - 1; x++) { + float s; E = sim->fields[2].E + z * fi->zstr + y * fi->ystr + x * fi->xstr; - H = sim->fields[2].H + z * fi->zstr + y * fi->ystr + x * fi->xstr; - chuj = (0.5 - sqrt(fx * fx + fy * fy + fz * fz)) / 0.5; - if (chuj < 0) - chuj = 0; - else - chuj = 1; - - chuj *= sin(sim->time * 3) * 1; - //chuj *= exp(-10.0f * pow(sim->time - 0.5, 2)) * 2; - E[0] += chuj; + s = exp(-600.0f * pow((float)x / fi->width - 0.15f, 2)); + s *= sin(sim->time); + //s *= exp(-2.0f * pow(sim->time - 2.0f, 2)) * 1000; + E[0] += s; } } -void phy_sim_step(phy_sim *sim, float dt) { +void phy_sim_step(phy_sim *sim) { phy_field_info *fi = &sim->field_info; phy_field_em *fields = sim->fields, tmp; - size_t i; - - phy_sim_add_sources(sim); // note: only pointers are rotated, not the entire fields SDL_mutexP(sim->rotate_lock); @@ -192,12 +203,13 @@ void phy_sim_step(phy_sim *sim, float dt) { memcpy(fields + 2, &tmp, sizeof(phy_field_em)); SDL_mutexV(sim->rotate_lock); - phy_sim_step_1(fi, dt, sim->field_eps, - fields[0].E, fields[2].E, fields[1].H); - phy_sim_step_1(fi, dt, sim->field_mu, - fields[0].H, fields[2].H, fields[1].E); + phy_sim_step_H(fi, sim->field_eps, sim->time_delta, + fields[2].H, fields[0].H, fields[1].E); + phy_sim_step_E(fi, sim->field_mu, sim->time_delta, + fields[2].E, fields[0].E, fields[1].H); + phy_sim_add_sources(sim); - sim->time += dt; + sim->time += sim->time_delta; } void phy_run_command(phy_sim *sim, int number, void *data) @@ -212,7 +224,12 @@ void phy_run_command(phy_sim *sim, int number, void *data) break; case PHY_CMD_STEP: - phy_sim_step(sim, 0.1); + sim->running = false; + phy_sim_step(sim); + break; + + case PHY_CMD_RESET: + phy_sim_reset(sim); break; } } @@ -223,15 +240,19 @@ int phy_thread(phy_sim *sim) void *cmd_data; while (1) { +cont_main_loop: if (sim->running) { while (!itc_chan_pop(&sim->ctl, &cmd_num, &cmd_data)) { if (cmd_num == PHY_CMD_QUIT) goto out; phy_run_command(sim, cmd_num, cmd_data); + + if (!sim->running) + goto cont_main_loop; } - phy_sim_step(sim, 0.1); + phy_sim_step(sim); } else { // do it only once per iteration if (!itc_chan_pop_block(&sim->ctl, &cmd_num, diff --git a/src/physics.h b/src/physics.h index 533d8a4..92ec6fd 100644 --- a/src/physics.h +++ b/src/physics.h @@ -21,14 +21,15 @@ enum { PHY_CMD_QUIT, PHY_CMD_PAUSE, PHY_CMD_RESUME, - PHY_CMD_STEP + PHY_CMD_STEP, + PHY_CMD_RESET }; typedef struct { phy_field_info field_info; phy_field_em fields[3]; float *field_eps, *field_mu; //permittivity and permeability - float time; + float time, time_delta; // UI stuff bool running; @@ -39,7 +40,7 @@ typedef struct { void phy_sim_destroy(phy_sim *sim); int phy_sim_create(phy_sim *sim); void phy_sim_compute_const_fields(phy_sim *sim); -void phy_sim_step(phy_sim *sim, float dt); +void phy_sim_step(phy_sim *sim); int phy_thread(phy_sim *sim); @@ -24,7 +24,7 @@ struct { typedef struct { phy_sim *sim; - bool dragging; + bool dragging, selecting; r_xsection_type xsection_type; float xsection_frac; bool select_valid; @@ -143,27 +143,7 @@ void ui_infof(ui_window *uiw, const char *fmt, ...) static void ui_simview_set_select(ui_simview *sv, vec2_t sel_2d) { - if (sel_2d[0] < 0.0f || sel_2d[0] > 1.0f || - sel_2d[1] < 0.0f || sel_2d[1] > 1.0f) { - sv->select_valid = false; - return; - } - - switch (sv->xsection_type) { - case XSECTION_XY: - v3_set(sv->select, sel_2d[0], sel_2d[1], sv->xsection_frac); - break; - - case XSECTION_XZ: - v3_set(sv->select, sel_2d[0], sv->xsection_frac, sel_2d[1]); - break; - case XSECTION_YZ: - v3_set(sv->select, sv->xsection_frac, sel_2d[0], sel_2d[1]); - break; - } - - sv->select_valid = true; } void ui_event_window(SDL_Event *event, ui_window *uiw) @@ -192,22 +172,21 @@ void ui_event_window(SDL_Event *event, ui_window *uiw) break; case SDL_MOUSEBUTTONDOWN: - if (event->button.button == SDL_BUTTON_LEFT && - event->button.y > sv->margin_top && - event->button.y < uiw->h - sv->margin_bottom) - sv->dragging = true; - - if (event->button.button == SDL_BUTTON_RIGHT) { - vec2_t sel_2d; + if (event->button.y > sv->margin_top && + event->button.y < uiw->h - sv->margin_bottom) { + if (event->button.button == SDL_BUTTON_LEFT) + sv->dragging = true; - v2_div_mst2(sel_2d, uiw->mouse, sv->tf_x2s); - ui_simview_set_select(sv, sel_2d); + if (event->button.button == SDL_BUTTON_RIGHT) + sv->selecting = true; } break; case SDL_MOUSEBUTTONUP: if (event->button.button == SDL_BUTTON_LEFT) sv->dragging = false; + else if (event->button.button == SDL_BUTTON_RIGHT) + sv->selecting = false; break; case SDL_MOUSEMOTION: @@ -256,6 +235,14 @@ void ui_event_window(SDL_Event *event, ui_window *uiw) case SDLK_p: itc_chan_push(&sv->sim->ctl, PHY_CMD_RESUME, NULL); break; + + case SDLK_l: + itc_chan_push(&sv->sim->ctl, PHY_CMD_STEP, NULL); + break; + + case SDLK_r: + itc_chan_push(&sv->sim->ctl, PHY_CMD_RESET, NULL); + break; } } } @@ -342,6 +329,38 @@ void ui_draw_simview_selection(ui_window *uiw) phy_sim *sim = sv->sim; vec2_t select_s; + if (sv->selecting) + { + vec2_t select_2d; + + v2_div_mst2(select_2d, uiw->mouse, sv->tf_x2s); + + if (select_2d[0] < 0.0f || select_2d[0] > 1.0f || + select_2d[1] < 0.0f || select_2d[1] > 1.0f) { + sv->select_valid = false; + return; + } + + switch (sv->xsection_type) { + case XSECTION_XY: + v3_set(sv->select, + select_2d[0], select_2d[1], sv->xsection_frac); + break; + + case XSECTION_XZ: + v3_set(sv->select, + select_2d[0], sv->xsection_frac, select_2d[1]); + break; + + case XSECTION_YZ: + v3_set(sv->select, + sv->xsection_frac, select_2d[0], select_2d[1]); + break; + } + + sv->select_valid = true; + } + if (!sv->select_valid) return; @@ -377,6 +396,8 @@ void ui_draw_simview_bars(ui_window *uiw, float dt) ui_simview *sv = &uiw->simview; phy_sim *sim = sv->sim; + phy_field_em *field_em = sim->fields + 1; + phy_field_info *fi = &sim->field_info; // upper @@ -406,16 +427,30 @@ void ui_draw_simview_bars(ui_window *uiw, float dt) if (sv->select_valid) { size_t x, y, z, offs; - sv->margin_bottom += 2 * theme.font_size; + sv->margin_bottom += 3 * theme.font_size; - r_draw_rect(uiw->rw, 0, uiw->h - 3 * theme.font_size, - uiw->w, 2 * theme.font_size, theme.color_main); + r_draw_rect(uiw->rw, 0, uiw->h - 4 * theme.font_size, + uiw->w, 3 * theme.font_size, theme.color_main); - r_draw_text(uiw->rw, 0, uiw->h - theme.font_size, - theme.font_size, - va("x = [%f %f %f]", sv->select[0], - sv->select[1], sv->select[2]), + x = sv->select[0] * fi->width; + y = sv->select[1] * fi->height; + z = sv->select[2] * fi->depth; + offs = z * fi->zstr + y * fi->ystr + x * fi->xstr; + + r_draw_text(uiw->rw, 0, uiw->h - 4 * theme.font_size, + theme.font_size, va("x = [%zu %zu %zu]", x, y, z), theme.color_text, 0); + + + r_draw_text(uiw->rw, 0, uiw->h - 3 * theme.font_size, + theme.font_size, va("E(x) = [%+E %+E %+E]", + field_em->E[offs], field_em->E[offs + 1], + field_em->E[offs + 2]), theme.color_text, 0); + + r_draw_text(uiw->rw, 0, uiw->h - 2 * theme.font_size, + theme.font_size, va("H(x) = [%+E %+E %+E]", + field_em->H[offs], field_em->H[offs + 1], + field_em->H[offs + 2]), theme.color_text, 0); } } |