diff options
Diffstat (limited to 'src/physics.c')
-rw-r--r-- | src/physics.c | 190 |
1 files changed, 142 insertions, 48 deletions
diff --git a/src/physics.c b/src/physics.c index fedc4c4..c7cc477 100644 --- a/src/physics.c +++ b/src/physics.c @@ -55,9 +55,48 @@ static void phy_calc_uec_E(float *uec, float dt, float fc0, float fc1, uec[3] *= 1.0f / eps; } -static float remap(float x, float a0, float b0, float a1, float b1) + +static float phy_sim_csg_subsample(phy_csg_step *stack, float *x) +{ + phy_csg_step *step; + float rv = 1.0f; + + for (step = stack; step; step = step->next) { + switch (step->shape) { + case CSG_SPHERE: + if (sq(x[0] - step->x[0]) + sq(x[1] - step->x[1]) + + sq(x[2] - step->x[2]) > sq(step->d[0])) + goto next_step; + break; + } + + rv = step->value; + + next_step: + ; + } + + return rv; +} + +static float phy_sim_csg_sample(phy_csg_step *stack, size_t *x) { - return (x - a0) * (b1 - a1) / (b0 - a0) + a1; + vec3_t xf = {x[0], x[1], x[2]}; + size_t dx, dy, dz; + float total; + + for (dz = 0; dz < 2; dz++) + for (dy = 0; dy < 2; dy++) + for (dx = 0; dx < 2; dx++) { + vec3_t xf; + + v3_set(xf, x[0] + 0.5f * dx, x[1] + 0.5f * dy, + x[2] + 0.5f * dz); + + total += phy_sim_csg_subsample(stack, xf); + } + + return total / 8.0f; } static void phy_sim_zero_aux(phy_sim *sim) @@ -68,6 +107,8 @@ static void phy_sim_zero_aux(phy_sim *sim) memset(sim->aux, 0, fi->size1 * sizeof(phy_field_aux_point)); + sim->eps_min = sim->eps_max = 1.0f; + for (x[2] = 1; x[2] < fi->dims[2] - 1; x[2]++) for (x[1] = 1; x[1] < fi->dims[1] - 1; x[1]++) for (x[0] = 1; x[0] < fi->dims[0] - 1; x[0]++) { @@ -77,23 +118,14 @@ static void phy_sim_zero_aux(phy_sim *sim) aux = sim->aux + x[2] * fi->zstr1 + x[1] * fi->ystr1 + x[0] * fi->xstr1; - aux->eps = 1; - aux->mu = 1; - - if (0){ - float cx, cy, cz, cr, d; + aux->eps = phy_sim_csg_sample(sim->csg_stack, x); - cx = 25, cy = 25, cz = 65; - cr = 10; + if (aux->eps > sim->eps_max) + sim->eps_max = aux->eps; + if (aux->eps < sim->eps_min) + sim->eps_min = aux->eps; - d = sqrt(sq(x[0] - cx) + sq(x[1] - cy) + sq(x[2] - cz)); - - if (d < cr) { - aux->eps = remap(d, cr, cr * 0.8, 1, 2); - if (aux->eps > 2) - aux->eps = 2; - } - } + aux->mu = 1; v3_set(aux->int_curl_E, 0, 0, 0); v3_set(aux->int_curl_H, 0, 0, 0); @@ -122,6 +154,10 @@ static void phy_sim_zero(phy_sim *sim) size_t i; phy_field_info *fi; + sim->pml_widths[0] = con_var_get_size("pml_x"); + sim->pml_widths[1] = con_var_get_size("pml_y"); + sim->pml_widths[2] = con_var_get_size("pml_z"); + fi = &sim->field_info; for (i = 0; i < 3; i++) { @@ -149,6 +185,18 @@ int phy_sim_create(phy_sim *sim) con_var_set("depth", "50"); con_var_set("scale", "7e+9"); con_var_set("time_delta", "0.1"); + con_var_set("pml_x", "0"); + con_var_set("pml_y", "0"); + con_var_set("pml_z", "0"); + + con_var_set("src_z", "10"); + con_var_set("src_angle", "0"); + + con_var_set("csg_x", NULL); + con_var_set("csg_y", NULL); + con_var_set("csg_z", NULL); + con_var_set("csg_r", NULL); + con_var_set("csg_v", NULL); phy_sim_create_fields(sim); @@ -157,16 +205,28 @@ int phy_sim_create(phy_sim *sim) void phy_sim_destroy(phy_sim *sim) { - size_t i; - + phy_sim_csg_clear(sim); + phy_sim_destroy_fields(sim); SDL_DestroyMutex(sim->rotate_lock); - itc_chan_destroy(&sim->ctl); +} + +static void phy_sim_configure(phy_sim *sim) +{ + phy_field_info *fi = &sim->field_info; + + fi->spacing = con_var_get_float("scale") / + min3(fi->dims[0], fi->dims[1], fi->dims[2]); + sim->time_delta = con_var_get_float("time_delta"); + + con_var_set("courant", va("%f", 3 * LIGHT_SPEED * sim->time_delta / + fi->spacing)); + + sim->source_xyplane.z = con_var_get_size("src_z"); + sim->source_xyplane.Px = cos(DEG2RAD(con_var_get_float("src_angle"))); + sim->source_xyplane.Py = sin(DEG2RAD(con_var_get_float("src_angle"))); + sim->source_xyplane.enabled = true; - for (i = 0; i < 3; i++) { - free(sim->fields[i].E); - free(sim->fields[i].H); - } } int phy_sim_create_fields(phy_sim *sim) @@ -179,29 +239,9 @@ int phy_sim_create_fields(phy_sim *sim) fi->dims[0] = con_var_get_size("width"); fi->dims[1] = con_var_get_size("height"); fi->dims[2] = con_var_get_size("depth"); - fi->spacing = con_var_get_float("scale") / - min3(fi->dims[0], fi->dims[1], fi->dims[2]); - sim->time_delta = con_var_get_float("time_delta"); - con_var_set("courant", va("%f", 3 * LIGHT_SPEED * sim->time_delta / - fi->spacing)); - - #if 0 - sim->pml_widths[0] = 0; - sim->pml_widths[1] = 0; - sim->pml_widths[2] = fi->dims[2] / 5 * 1; - #else - sim->pml_widths[0] = 0; - sim->pml_widths[1] = 0; - sim->pml_widths[2] = 0; - #endif - - sim->source_xyplane.z = fi->dims[2] / 2;//sim->pml_widths[2] + 5; - sim->source_xyplane.Px = 0; - sim->source_xyplane.Py = 1; - sim->source_xyplane.lambda = 20; - sim->source_xyplane.time = 2; - sim->source_xyplane.enabled = true; + printf("phy_sim_create_fields: %zux%zux%zu\n", + fi->dims[0], fi->dims[1], fi->dims[2]); fi->xstr = 3; fi->ystr = fi->xstr * fi->dims[0]; @@ -228,6 +268,7 @@ int phy_sim_create_fields(phy_sim *sim) goto error; } + phy_sim_configure(sim); phy_sim_zero(sim); sim->running = false; @@ -470,6 +511,44 @@ void phy_sim_step(phy_sim *sim) { SDL_mutexV(sim->rotate_lock); } +void phy_sim_csg_push(phy_sim *sim) +{ + char *shape_str; + phy_csg_shape shape; + phy_csg_step *step; + + // TODO shape_str = con_var_get("csg_shape"); + + shape = CSG_SPHERE; + + step = malloc(sizeof(phy_csg_step)); + assert(step); + + step->next = sim->csg_stack; + sim->csg_stack = step; + + step->value = con_var_get_float("csg_v"); + + switch (shape) { + case CSG_SPHERE: + step->x[0] = con_var_get_float("csg_x"); + step->x[1] = con_var_get_float("csg_y"); + step->x[2] = con_var_get_float("csg_z"); + step->d[0] = con_var_get_float("csg_r"); + break; + } +} + +void phy_sim_csg_clear(phy_sim *sim) +{ + phy_csg_step *step, *next; + + for (step = sim->csg_stack; step; step = next) { + next = step->next; + free(step); + } +} + void phy_sim_debug(phy_sim *sim, size_t *coords) { phy_field_info *fi = &sim->field_info; @@ -538,10 +617,25 @@ void phy_run_command(phy_sim *sim, int number, void *data) case PHY_CMD_RESET: SDL_LockMutex(sim->rotate_lock); - phy_sim_destroy_fields(sim); - phy_sim_create_fields(sim); + if (sim->field_info.dims[0] == con_var_get_size("width") && + sim->field_info.dims[1] == con_var_get_size("height") && + sim->field_info.dims[2] == con_var_get_size("depth")) { + phy_sim_configure(sim); + phy_sim_zero(sim); + } else { + phy_sim_destroy_fields(sim); + phy_sim_create_fields(sim); + } SDL_UnlockMutex(sim->rotate_lock); break; + + case PHY_CMD_CSG_PUSH: + phy_sim_csg_push(sim); + break; + + case PHY_CMD_CSG_CLEAR: + phy_sim_csg_clear(sim); + break; } } |