summaryrefslogtreecommitdiff
path: root/src/physics.c
diff options
context:
space:
mode:
authorPaweł Redman <pawel.redman@gmail.com>2016-04-06 20:17:11 +0200
committerPaweł Redman <pawel.redman@gmail.com>2016-04-06 20:17:11 +0200
commitd44eacf9e805904867dd28524ff2b76540cb6ae1 (patch)
treeb62158b2bbf0d450d565a0c0d781eb35d11f7316 /src/physics.c
parentb3a2f32872391415c0da4e19403da14b4ad66cf8 (diff)
PBC on sides, fixes plane wave source.
Diffstat (limited to 'src/physics.c')
-rw-r--r--src/physics.c84
1 files changed, 54 insertions, 30 deletions
diff --git a/src/physics.c b/src/physics.c
index a50ba1f..0946b06 100644
--- a/src/physics.c
+++ b/src/physics.c
@@ -186,9 +186,9 @@ int phy_sim_create(phy_sim *sim)
assert(sim->rotate_lock);
itc_chan_create(&sim->ctl);
- con_var_set("width", "50");
- con_var_set("height", "50");
- con_var_set("depth", "50");
+ con_var_set("width", "40");
+ con_var_set("height", "40");
+ con_var_set("depth", "70");
con_var_set("scale", "7e+9");
con_var_set("time_delta", "0.1");
con_var_set("pml_x", "0");
@@ -197,6 +197,8 @@ int phy_sim_create(phy_sim *sim)
con_var_set("src_z", "10");
con_var_set("src_angle", "0");
+ con_var_set("src_f", "0.1");
+ con_var_set("src_wtf", "0");
con_var_set("csg_shape", "sphere");
con_var_set("csg_v", NULL);
@@ -231,11 +233,11 @@ static void phy_sim_configure(phy_sim *sim)
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;
-
+ sim->source_z = con_var_get_size("src_z");
+ sim->source_f = con_var_get_float("src_f");
+ sim->source_px = cos(DEG2RAD(con_var_get_float("src_angle")));
+ sim->source_py = sin(DEG2RAD(con_var_get_float("src_angle")));
+ sim->source_wtf = con_var_get_float("src_wtf");
}
int phy_sim_create_fields(phy_sim *sim)
@@ -307,9 +309,9 @@ void phy_sim_destroy_fields(phy_sim *sim)
sim->valid = false;
}
-static float gauss(phy_source_xyplane *s, float t)
+static float phy_sim_eval_source(phy_sim *sim, float t)
{
- return sin(t); //exp(-s->lambda * sq(t - s->time));
+ return sin(t * 2.0f * M_PI * sim->source_f);
}
void phy_sim_add_sources_H(phy_sim *sim)
@@ -318,20 +320,13 @@ void phy_sim_add_sources_H(phy_sim *sim)
size_t x, y, z;
float src_A, src_x, src_y, delay;
- if (!sim->source_xyplane.enabled)
- return;
-
- z = sim->source_xyplane.z;
-
- delay = fi->spacing / (2.0f * LIGHT_SPEED)
- - sim->time_delta / 2.0f + 0.35;
+ z = sim->source_z;
- /*delay = - fi->spacing / (2.0f * LIGHT_SPEED)
- + sim->time_delta / 2.0f;*/
+ delay = sim->source_wtf;
- src_A = gauss(&sim->source_xyplane, sim->time + delay) * 1;
- src_x = -sim->source_xyplane.Py * src_A;
- src_y = sim->source_xyplane.Px * src_A;
+ src_A = phy_sim_eval_source(sim, sim->time + delay);
+ src_x = -sim->source_py * src_A;
+ src_y = sim->source_px * src_A;
for (y = 1; y < fi->dims[1] - 1; y++)
for (x = 1; x < fi->dims[0] - 1; x++) {
@@ -351,14 +346,11 @@ void phy_sim_add_sources_E(phy_sim *sim)
size_t x, y, z;
float src_x, src_y, src_A;
- if (!sim->source_xyplane.enabled)
- return;
+ z = sim->source_z - 1;
- z = sim->source_xyplane.z - 1;
-
- src_A = gauss(&sim->source_xyplane, sim->time) * 1;
- src_x = sim->source_xyplane.Px * src_A;
- src_y = sim->source_xyplane.Py * src_A;
+ src_A = phy_sim_eval_source(sim, sim->time);
+ src_x = sim->source_px * src_A;
+ src_y = sim->source_py * src_A;
for (y = 1; y < fi->dims[1] - 1; y++)
for (x = 1; x < fi->dims[0] - 1; x++) {
@@ -494,6 +486,36 @@ void phy_sim_update_em_E(phy_sim *sim)
}
}
+void phy_sim_wrap_field(phy_sim *sim, float *field)
+{
+ phy_field_info *fi = &sim->field_info;
+ size_t x[3];
+
+#define AT(x,y,z) (field + (x) * fi->xstr + (y) * fi->ystr + \
+ (z) * fi->zstr)
+
+ x[0] = 0;
+ for (x[2] = 0; x[2] < fi->dims[2]; x[2]++)
+ for (x[1] = 0; x[1] < fi->dims[1]; x[1]++)
+ v3_copy(AT(x[0], x[1], x[2]), AT(fi->dims[0] - 2, x[1], x[2]));
+
+ x[0] = fi->dims[0] - 1;
+ for (x[2] = 0; x[2] < fi->dims[2]; x[2]++)
+ for (x[1] = 0; x[1] < fi->dims[1]; x[1]++)
+ v3_copy(AT(fi->dims[0] - 1, x[1], x[2]), AT(1, x[1], x[2]));
+
+ x[1] = 0;
+ for (x[2] = 0; x[2] < fi->dims[2]; x[2]++)
+ for (x[0] = 0; x[0] < fi->dims[0]; x[0]++)
+ v3_copy(AT(x[0], 0, x[2]), AT(x[0], fi->dims[1] - 2, x[2]));
+
+ x[1] = fi->dims[1] - 1;
+ for (x[2] = 0; x[2] < fi->dims[2]; x[2]++)
+ for (x[0] = 0; x[0] < fi->dims[0]; x[0]++)
+ v3_copy(AT(x[0], fi->dims[1] - 1, x[2]), AT(x[0], 1, x[2]));
+#undef AT
+}
+
void phy_sim_step(phy_sim *sim) {
phy_field_em *fields = sim->fields, tmp;
int64_t start_time;
@@ -502,8 +524,10 @@ void phy_sim_step(phy_sim *sim) {
phy_sim_update_aux_H(sim);
phy_sim_update_em_H(sim);
+ phy_sim_wrap_field(sim, fields[2].H);
phy_sim_update_aux_E(sim);
- phy_sim_update_em_E(sim);
+ phy_sim_update_em_E(sim);
+ phy_sim_wrap_field(sim, fields[2].E);
SDL_mutexP(sim->rotate_lock);