diff --git a/sbp/src/euler.rs b/sbp/src/euler.rs index b86efd9..840d301 100644 --- a/sbp/src/euler.rs +++ b/sbp/src/euler.rs @@ -54,49 +54,20 @@ impl System { #[allow(clippy::many_single_char_names)] pub fn init_with_vortex(&mut self, x0: f32, y0: f32) { // Should parametrise such that we have radius, drop in pressure at center, etc - let rstar = 1.0; - let eps = 3.0; - #[allow(non_snake_case)] - let M = 0.5; + let vortex_parameters = VortexParameters { + x0, + y0, + rstar: 1.0, + eps: 3.0, + mach: 0.5, + }; - let p_inf = 1.0 / (GAMMA * M * M); - let t = 0.0; - - let nx = self.grid.nx(); - let ny = self.grid.ny(); - - for j in 0..ny { - for i in 0..nx { - let x = self.grid.x[(j, i)]; - let y = self.grid.y[(j, i)]; - - let dx = (x - x0) - t; - let dy = y - y0; - let f = (1.0 - (dx * dx + dy * dy)) / (rstar * rstar); - - use std::f32::consts::PI; - let u = - 1.0 - eps * dy / (2.0 * PI * p_inf.sqrt() * rstar * rstar) * (f / 2.0).exp(); - let v = - 0.0 + eps * dx / (2.0 * PI * p_inf.sqrt() * rstar * rstar) * (f / 2.0).exp(); - let rho = f32::powf( - 1.0 - eps * eps * (GAMMA - 1.0) * M * M - / (8.0 * PI * PI * p_inf * rstar * rstar) - * f.exp(), - 1.0 / (GAMMA - 1.0), - ); - assert!(rho > 0.0); - let p = p_inf * rho.powf(GAMMA); - assert!(p > 0.0); - let e = p / (GAMMA - 1.0) + rho * (u * u + v * v) / 2.0; - assert!(e > 0.0); - - self.sys.0[(0, j, i)] = rho; - self.sys.0[(1, j, i)] = rho * u; - self.sys.0[(2, j, i)] = rho * v; - self.sys.0[(3, j, i)] = e; - } - } + self.sys.0.vortex( + self.grid.x.view(), + self.grid.y.view(), + 0.0, + vortex_parameters, + ) } pub fn field(&self) -> &Field { @@ -242,6 +213,56 @@ impl Field { fn west_mut(&mut self) -> ArrayViewMut2 { self.slice_mut(s![.., .., 0]) } + + fn vortex( + &mut self, + x: ArrayView2, + y: ArrayView2, + t: f32, + vortex_param: VortexParameters, + ) { + assert_eq!(x.shape(), y.shape()); + assert_eq!(x.shape()[1], self.nx()); + assert_eq!(x.shape()[0], self.ny()); + + let (rho, rhou, rhov, e) = self.components_mut(); + + let eps = vortex_param.eps; + let m = vortex_param.mach; + let rstar = vortex_param.rstar; + let p_inf = 1.0 / (GAMMA * m * m); + azip!((rho in rho, + rhou in rhou, + rhov in rhov, + e in e, + x in x, + y in y) + { + use std::f32::consts::PI; + let dx = (x - vortex_param.x0) - t; + let dy = y - vortex_param.y0; + let f = (1.0 - (dx*dx + dy*dy))/(rstar*rstar); + + *rho = f32::powf(1.0 - eps*eps*(GAMMA - 1.0)*m*m/(8.0*PI*PI*p_inf*rstar*rstar)*f.exp(), 1.0/(GAMMA - 1.0)); + assert!(*rho > 0.0); + let p = f32::powf(*rho, GAMMA)*p_inf; + assert!(p > 0.0); + let u = 1.0 - eps*dy/(2.0*PI*p_inf.sqrt()*rstar*rstar)*(f/2.0).exp(); + let v = eps*dx/(2.0*PI*p_inf.sqrt()*rstar*rstar)*(f/2.0).exp(); + *rhou = *rho*u; + *rhov = *rho*v; + *e = p/(GAMMA - 1.0) + *rho*(u*u + v*v)/2.0; + }); + } +} + +#[derive(Copy, Clone)] +pub struct VortexParameters { + x0: f32, + y0: f32, + rstar: f32, + eps: f32, + mach: f32, } fn pressure(gamma: f32, rho: f32, rhou: f32, rhov: f32, e: f32) -> f32 { diff --git a/webfront/Cargo.toml b/webfront/Cargo.toml index 91b5984..0e70c41 100644 --- a/webfront/Cargo.toml +++ b/webfront/Cargo.toml @@ -9,8 +9,8 @@ crate-type = ["cdylib"] path = "lib.rs" [dependencies] -wasm-bindgen = "0.2.54" -console_error_panic_hook = { version = "0.1.6" } -wee_alloc = { version = "0.4.5" } +wasm-bindgen = "0.2.58" +console_error_panic_hook = "0.1.6" +wee_alloc = "0.4.5" sbp = { path = "../sbp" } ndarray = "0.13.0" diff --git a/webfront/lib.rs b/webfront/lib.rs index edeba8f..13f3a9a 100644 --- a/webfront/lib.rs +++ b/webfront/lib.rs @@ -8,7 +8,6 @@ static ALLOC: wee_alloc::WeeAlloc = wee_alloc::WeeAlloc::INIT; #[wasm_bindgen] pub fn set_panic_hook() { - #[cfg(feature = "console_error_panic_hook")] console_error_panic_hook::set_once(); } diff --git a/webfront/make_wasm.py b/webfront/make_wasm.py index c1235c7..c94cb84 100755 --- a/webfront/make_wasm.py +++ b/webfront/make_wasm.py @@ -24,7 +24,7 @@ if __name__ == "__main__": publish.mkdir(exist_ok=True) target_triple = "wasm32-unknown-unknown" - command = ["env", "RUSTFLAGS=-Clto=thin", "cargo", "build", "--target", target_triple] + command = ["cargo", "build", "--target", target_triple] if args.release: command.append("--release")