mulimoen
/
SummationByParts
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

do EW/NS boundaries in same loop

This commit is contained in:
Magnus Ulimoen 2019-09-03 20:26:07 +02:00
parent 25f177b20e
commit 3513cc496a
1 changed files with 11 additions and 17 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
src/maxwell.rs
View File

@ -106,33 +106,30 @@ impl System {
let hinv = 1.0 / (SBP::h()[0] / (nx - 1) as f32); let hinv = 1.0 / (SBP::h()[0] / (nx - 1) as f32);
// East boundary
for j in 0..ny { for j in 0..ny {
// East boundary
let tau = -1.0; let tau = -1.0;
let g = (y.0[(j, 0)], y.1[(j, 0)], y.2[(j, 0)]); let g = (y.0[(j, 0)], y.1[(j, 0)], y.2[(j, 0)]);
let v = (y.0[(j, nx - 1)], y.1[(j, nx - 1)], y.2[(j, nx - 1)]); let v = (y.0[(j, nx - 1)], y.1[(j, nx - 1)], y.2[(j, nx - 1)]);
// A+ = (0, 0, 0; 0, 1/2, -1/2; 0, -1/2, 1/2); // A+ = (0, 0, 0; 0, 1/2, -1/2; 0, -1/2, 1/2);
k[i].0[(j, nx - 1)] += 0.0; // k[i].0[(j, nx - 1)] += 0.0;
k[i].1[(j, nx - 1)] += tau * hinv * (0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2)); k[i].1[(j, nx - 1)] += tau * hinv * (0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2));
k[i].2[(j, nx - 1)] += tau * hinv * (-0.5 * (v.1 - g.1) + 0.5 * (v.2 - g.2)); k[i].2[(j, nx - 1)] += tau * hinv * (-0.5 * (v.1 - g.1) + 0.5 * (v.2 - g.2));
}
// West boundary
for j in 0..ny {
let tau = 1.0;
let g = (y.0[(j, nx - 1)], y.1[(j, nx - 1)], y.2[(j, nx - 1)]);
let v = (y.0[(j, 0)], y.1[(j, 0)], y.2[(j, 0)]);
// West boundary
let tau = 1.0;
let (v, g) = (g, v);
// A- = (0, 0, 0; 0, -1/2, -1/2; 0, -1/2, -1/2); // A- = (0, 0, 0; 0, -1/2, -1/2; 0, -1/2, -1/2);
k[i].0[(j, 0)] += 0.0; // k[i].0[(j, 0)] += 0.0;
k[i].1[(j, 0)] += tau * hinv * (-0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2)); k[i].1[(j, 0)] += tau * hinv * (-0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2));
k[i].2[(j, 0)] += tau * hinv * (-0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2)); k[i].2[(j, 0)] += tau * hinv * (-0.5 * (v.1 - g.1) - 0.5 * (v.2 - g.2));
} }
let hinv = 1.0 / (SBP::h()[0] / (ny - 1) as f32); let hinv = 1.0 / (SBP::h()[0] / (ny - 1) as f32);
// North boundary
for j in 0..nx { for j in 0..nx {
// North boundary
let tau = -1.0; let tau = -1.0;
let g = (y.0[(0, j)], y.1[(0, j)], y.2[(0, j)]); let g = (y.0[(0, j)], y.1[(0, j)], y.2[(0, j)]);
let v = (y.0[(ny - 1, j)], y.1[(ny - 1, j)], y.2[(ny - 1, j)]); let v = (y.0[(ny - 1, j)], y.1[(ny - 1, j)], y.2[(ny - 1, j)]);
@ -140,19 +137,16 @@ impl System {
// B+ = (1/2, 1/2, 0; 1/2, 1/2, 0; 0, 0, 0) // B+ = (1/2, 1/2, 0; 1/2, 1/2, 0; 0, 0, 0)
k[i].0[(ny - 1, j)] += tau * hinv * (0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1)); k[i].0[(ny - 1, j)] += tau * hinv * (0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1));
k[i].1[(ny - 1, j)] += tau * hinv * (0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1)); k[i].1[(ny - 1, j)] += tau * hinv * (0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1));
k[i].2[(ny - 1, j)] += 0.0; // k[i].2[(ny - 1, j)] += 0.0;
}
// South boundary // South boundary
for j in 0..nx {
let tau = 1.0; let tau = 1.0;
let g = (y.0[(ny - 1, j)], y.1[(ny - 1, j)], y.2[(ny - 1, j)]); let (g, v) = (v, g);
let v = (y.0[(0, j)], y.1[(0, j)], y.2[(0, j)]);
// B- = (-1/2, 1/2, 0; 1/2, -1/2, 0; 0, 0, 0); // B- = (-1/2, 1/2, 0; 1/2, -1/2, 0; 0, 0, 0);
k[i].0[(0, j)] += tau * hinv * (-0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1)); k[i].0[(0, j)] += tau * hinv * (-0.5 * (v.0 - g.0) + 0.5 * (v.1 - g.1));
k[i].1[(0, j)] += tau * hinv * (0.5 * (v.0 - g.0) - 0.5 * (v.1 - g.1)); k[i].1[(0, j)] += tau * hinv * (0.5 * (v.0 - g.0) - 0.5 * (v.1 - g.1));
k[i].2[(0, j)] += 0.0; // k[i].2[(0, j)] += 0.0;
} }
} }