From cec0cc945f525b81bfffd7f87b8f8e7a5f54f001 Mon Sep 17 00:00:00 2001 From: "Timothy A. Smith" <2timothy18@gmail.com> Date: Fri, 15 Nov 2019 21:05:55 -0600 Subject: [PATCH 1/4] refactor out metric term creation into fixture --- test/data_for_test.py | 61 +++++++++++++++++++++++------------- test/test_eigensystem.py | 2 +- test/test_flux_window_ops.py | 2 +- 3 files changed, 42 insertions(+), 23 deletions(-) diff --git a/test/data_for_test.py b/test/data_for_test.py index 2f12e9c..1b1db78 100644 --- a/test/data_for_test.py +++ b/test/data_for_test.py @@ -31,7 +31,7 @@ class WindowData: ndim = ref.gas.ndim dirs = ref.gas.dirs - def __init__(self, states_str, direction): + def __init__(self, states_str, direction, window_metrics): self.direction = self.dirs[direction] self.dir_internal = self.direction-1 @@ -47,30 +47,49 @@ class WindowData: self.fluxes = ref.pointwise_fluxes( self.states)[:,:,self.dir_internal].T.copy(order="F") - self.metrics = np.array([np.identity(self.ndim) for i in range(6)], - dtype=np.float64, order="F") - self.jacobians = np.repeat(1.0, 6) + self.metrics, self.jacobians = window_metrics + + +@pytest.fixture(scope="session") +def metrics_generator(): + ndim = ref.gas.ndim + + def window_metrics(metric_str): + if (metric_str == "uniform"): + metrics = np.array([np.identity(ndim) for i in range(6)], + dtype=np.float64, + order="F") + jacobians = np.repeat(1.0, 6) + else: + raise ValueError("Metric string {} not supported".format(metric_str)) + return metrics, jacobians + + return window_metrics @pytest.fixture(scope="session", params=[ - ("1 1 1 1 5.5,1 1 1 1 5.5", "x"), - ("1 1 1 1 5.5,1 1 1 1 5.5", "y"), - ("1 1 1 1 5.5,1 1 1 1 5.5", "z"), - ("2 4 4 4 20,1 1 1 1 5.5", "x"), - ("2 4 4 4 20,1 1 1 1 5.5", "y"), - ("2 4 4 4 20,1 1 1 1 5.5", "z"), - ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "x"), - ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "y"), - ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "z"), - ("2 4 8 12 64,1 1 2 3 11", "x"), - ("2 8 12 4 64,1 2 3 1 11", "y"), - ("2 12 4 8 64,1 3 1 2 11", "z"), - ("1 -1 -2 -3 11,2 -4 -8 -12 64", "x"), - ("1 -2 -3 -1 11,2 -8 -12 -4 64", "y"), - ("1 -3 -1 -2 11,2 -12 -4 -8 64", "z") + ("1 1 1 1 5.5,1 1 1 1 5.5", "x", "uniform"), + ("1 1 1 1 5.5,1 1 1 1 5.5", "y", "uniform"), + ("1 1 1 1 5.5,1 1 1 1 5.5", "z", "uniform"), + ("2 4 4 4 20,1 1 1 1 5.5", "x", "uniform"), + ("2 4 4 4 20,1 1 1 1 5.5", "y", "uniform"), + ("2 4 4 4 20,1 1 1 1 5.5", "z", "uniform"), + ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "x", "uniform"), + ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "y", "uniform"), + ("1 -1 -1 -1 5.5,2 -4 -4 -4 20", "z", "uniform"), + ("2 4 8 12 64,1 1 2 3 11", "x", "uniform"), + ("2 8 12 4 64,1 2 3 1 11", "y", "uniform"), + ("2 12 4 8 64,1 3 1 2 11", "z", "uniform"), + ("1 -1 -2 -3 11,2 -4 -8 -12 64", "x", "uniform"), + ("1 -2 -3 -1 11,2 -8 -12 -4 64", "y", "uniform"), + ("1 -3 -1 -2 11,2 -12 -4 -8 64", "z", "uniform") ]) -def window_data(request): - return WindowData(*request.param) +def window_data(request, metrics_generator): + states_str = request.param[0] + dir_str = request.param[1] + metric_str = request.param[2] + + return WindowData(states_str, dir_str, metrics_generator(metric_str)) # vim: foldmethod=marker diff --git a/test/test_eigensystem.py b/test/test_eigensystem.py index e8f3a9c..b8b0364 100644 --- a/test/test_eigensystem.py +++ b/test/test_eigensystem.py @@ -35,7 +35,7 @@ import pytest import utilities as u import weno_reference_implementation as ref -from data_for_test import window_data +from data_for_test import window_data, metrics_generator class PairData: diff --git a/test/test_flux_window_ops.py b/test/test_flux_window_ops.py index b48e98a..a32d5bb 100644 --- a/test/test_flux_window_ops.py +++ b/test/test_flux_window_ops.py @@ -35,7 +35,7 @@ import pytest import utilities as u import weno_reference_implementation as ref -from data_for_test import window_data +from data_for_test import window_data, metrics_generator class WindowResults: -- GitLab From 1e864da060155d519b0826ec6a492e3b6e0873f9 Mon Sep 17 00:00:00 2001 From: "Timothy A. Smith" <2timothy18@gmail.com> Date: Fri, 15 Nov 2019 21:50:04 -0600 Subject: [PATCH 2/4] update reference implementation of frozen metric computation to match paper --- test/data_for_test.py | 25 +++++++++++++++++++++++++ test/test_eigensystem.py | 25 +------------------------ test/test_flux_window_ops.py | 29 ++++++++++++++--------------- 3 files changed, 40 insertions(+), 39 deletions(-) diff --git a/test/data_for_test.py b/test/data_for_test.py index 1b1db78..9e62317 100644 --- a/test/data_for_test.py +++ b/test/data_for_test.py @@ -92,4 +92,29 @@ def window_data(request, metrics_generator): return WindowData(states_str, dir_str, metrics_generator(metric_str)) +class PairData: + def __init__(self, data): + self.nvars = data.nvars + self.ndim = data.ndim + self.direction = data.direction + self.dir_internal = data.dir_internal + + self.state_pair = data.states[:,2:4] + + # FIXME: these should be generalized fluxes + # FIXME: make a clear distinction between fluxes in physical and + # generalized coordinates + self.flux_pair = data.fluxes[:,2:4] + + weighted_metric_slice = data.metrics[2:4]/data.jacobians[2:4,None,None] + self.frozen_metrics = np.mean(weighted_metric_slice, axis=0) + self.frozen_jacobian = np.mean(data.jacobians[2:4]) + self.combined_frozen_metrics = np.sum(self.frozen_metrics**2, axis=1) + + +@pytest.fixture(scope="session") +def pair_data(window_data): + return PairData(window_data) + + # vim: foldmethod=marker diff --git a/test/test_eigensystem.py b/test/test_eigensystem.py index b8b0364..c077e68 100644 --- a/test/test_eigensystem.py +++ b/test/test_eigensystem.py @@ -35,30 +35,7 @@ import pytest import utilities as u import weno_reference_implementation as ref -from data_for_test import window_data, metrics_generator - - -class PairData: - def __init__(self, data): - self.nvars = data.nvars - self.ndim = data.ndim - self.direction = data.direction - self.dir_internal = data.dir_internal - - self.state_pair = data.states[:,2:4] - - # FIXME: these should be generalized fluxes - # FIXME: make a clear distinction between fluxes in physical and - # generalized coordinates - self.flux_pair = data.fluxes[:,2:4] - - # FIXME: is this the right computation? - self.frozen_metrics = np.mean(data.metrics[2:4], axis=0) - - -@pytest.fixture(scope="session") -def pair_data(window_data): - return PairData(window_data) +from data_for_test import pair_data, window_data, metrics_generator def test_roe_uniform_grid_ideal_gas(queue, pair_data): diff --git a/test/test_flux_window_ops.py b/test/test_flux_window_ops.py index a32d5bb..8387bda 100644 --- a/test/test_flux_window_ops.py +++ b/test/test_flux_window_ops.py @@ -35,11 +35,13 @@ import pytest import utilities as u import weno_reference_implementation as ref -from data_for_test import window_data, metrics_generator +from data_for_test import PairData, window_data, metrics_generator class WindowResults: def __init__(self, queue, data): + pair = PairData(data) + self.nvars = data.nvars self.ndim = data.ndim self.direction = data.direction @@ -47,17 +49,14 @@ class WindowResults: self.states = data.states self.fluxes = data.fluxes - # FIXME: duplicate of code in PairData - self.state_pair = data.states[:,2:4] + self.state_pair = pair.state_pair self.metrics = data.metrics self.jacobians = data.jacobians - # FIXME: should be computed directly from the metrics and jacobians - # FIXME: partial duplicate of code in PairData - self.frozen_metrics = np.mean(self.metrics[2:4], axis=0) - self.frozen_jacobian = np.mean(self.jacobians[2:4], axis=0) - self.combined_frozen_metrics = 1.0 + self.frozen_metrics = pair.frozen_metrics + self.frozen_jacobian = pair.frozen_jacobian + self.combined_frozen_metric = pair.combined_frozen_metrics[self.dir_internal] self.lam_pointwise = ref.lambda_pointwise( self.states, self.metrics, self.dir_internal) @@ -74,9 +73,9 @@ class WindowResults: self.oscillation_neg = ref.oscillation(self.char_fluxes_neg[:,::-1]) self.weno_weights_pos = ref.weno_weights( - self.oscillation_pos, self.combined_frozen_metrics) + self.oscillation_pos, self.combined_frozen_metric) self.weno_weights_neg = ref.weno_weights( - self.oscillation_neg, self.combined_frozen_metrics) + self.oscillation_neg, self.combined_frozen_metric) self.consistent = ref.consistent_part(self.fluxes) self.dissipation_pos = ref.dissipation_part( @@ -137,7 +136,7 @@ def test_weno_flux(queue, window_results): generalized_fluxes=data.fluxes, characteristic_fluxes_pos=data.char_fluxes_pos, characteristic_fluxes_neg=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metrics, + combined_frozen_metrics=data.combined_frozen_metric, R=data.R, flux=flux_dev) @@ -169,7 +168,7 @@ def test_dissipation_part_pos(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_pos, - combined_frozen_metrics=data.combined_frozen_metrics, + combined_frozen_metrics=data.combined_frozen_metric, R=data.R, dissipation_pos=dissipation_dev) @@ -186,7 +185,7 @@ def test_dissipation_part_neg(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metrics, + combined_frozen_metrics=data.combined_frozen_metric, R=data.R, dissipation_neg=dissipation_dev) @@ -203,7 +202,7 @@ def test_weno_weights_pos(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_pos, - combined_frozen_metrics=data.combined_frozen_metrics, + combined_frozen_metrics=data.combined_frozen_metric, w=weights_dev) sum_weights = np.sum(weights_dev.get(), axis=1) @@ -222,7 +221,7 @@ def test_weno_weights_neg(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metrics, + combined_frozen_metrics=data.combined_frozen_metric, w=weights_dev) sum_weights = np.sum(weights_dev.get(), axis=1) -- GitLab From 84a6e5289aa4a8830d7298bd21809f9e9b5747bd Mon Sep 17 00:00:00 2001 From: "Timothy A. Smith" <2timothy18@gmail.com> Date: Fri, 15 Nov 2019 21:55:40 -0600 Subject: [PATCH 3/4] add note on provenance of algorithm --- test/weno_reference_implementation.py | 1 + 1 file changed, 1 insertion(+) diff --git a/test/weno_reference_implementation.py b/test/weno_reference_implementation.py index 179402a..a7b4095 100644 --- a/test/weno_reference_implementation.py +++ b/test/weno_reference_implementation.py @@ -31,6 +31,7 @@ import loopy as lp # noqa import utilities as u import ideal_gas as gas +### WENO algorithm taken from Nonomura, et al., Comput. Fluids 107 (2015) 242-255 def pointwise_fluxes(states): return np.array([gas.flux(s) for s in states.T]) -- GitLab From 91b7427684cd9826d89fe8c524f402c75673e639 Mon Sep 17 00:00:00 2001 From: "Timothy A. Smith" <2timothy18@gmail.com> Date: Fri, 15 Nov 2019 22:17:39 -0600 Subject: [PATCH 4/4] rename combined_frozen_metrics in some locations to reflect that there's only one --- WENO.F90 | 44 ++++++++++++++++++------------------ test/test_flux_window_ops.py | 10 ++++---- 2 files changed, 27 insertions(+), 27 deletions(-) diff --git a/WENO.F90 b/WENO.F90 index d798cdb..5a7a6f4 100644 --- a/WENO.F90 +++ b/WENO.F90 @@ -509,14 +509,14 @@ subroutine split_characteristic_fluxes(nvars, & end subroutine subroutine weno_flux(nvars, generalized_fluxes, characteristic_fluxes_pos, & - characteristic_fluxes_neg, combined_frozen_metrics, R, flux) + characteristic_fluxes_neg, combined_frozen_metric, R, flux) implicit none integer, intent(in) :: nvars real*8, intent(in) :: generalized_fluxes(nvars, -2:3) real*8, intent(in) :: characteristic_fluxes_pos(nvars, -2:3) real*8, intent(in) :: characteristic_fluxes_neg(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(in) :: R(nvars, nvars) real*8, intent(out) :: flux(nvars) @@ -527,9 +527,9 @@ subroutine weno_flux(nvars, generalized_fluxes, characteristic_fluxes_pos, & call consistent_part(nvars, generalized_fluxes, consistent) call dissipation_part_pos(nvars, characteristic_fluxes_pos, & - combined_frozen_metrics, R, dissipation_pos) + combined_frozen_metric, R, dissipation_pos) call dissipation_part_neg(nvars, characteristic_fluxes_neg, & - combined_frozen_metrics, R, dissipation_neg) + combined_frozen_metric, R, dissipation_neg) do v=1,nvars flux(v) = consistent(v) + dissipation_pos(v) + dissipation_neg(v) @@ -554,36 +554,36 @@ subroutine consistent_part(nvars, generalized_fluxes, consistent) end subroutine subroutine dissipation_part_pos(nvars, characteristic_fluxes, & - combined_frozen_metrics, R, dissipation_pos) + combined_frozen_metric, R, dissipation_pos) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(in) :: R(nvars, nvars) real*8, intent(out) :: dissipation_pos(nvars) real*8 combined_fluxes(nvars) - call weno_combination_pos(nvars, characteristic_fluxes, combined_frozen_metrics, combined_fluxes) + call weno_combination_pos(nvars, characteristic_fluxes, combined_frozen_metric, combined_fluxes) call mult_mat_vec(nvars, nvars, -1.0d0/60, R, combined_fluxes, dissipation_pos) end subroutine subroutine weno_combination_pos(nvars, characteristic_fluxes, & - combined_frozen_metrics, combined_fluxes) + combined_frozen_metric, combined_fluxes) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(out) :: combined_fluxes(nvars) real*8 w(nvars, 3) real*8 flux_differences(nvars, 3) integer v - call weno_weights_pos(nvars, characteristic_fluxes, combined_frozen_metrics, w) + call weno_weights_pos(nvars, characteristic_fluxes, combined_frozen_metric, w) call flux_differences_pos(nvars, characteristic_fluxes, flux_differences) do v=1,nvars @@ -593,12 +593,12 @@ subroutine weno_combination_pos(nvars, characteristic_fluxes, & end do end subroutine -subroutine weno_weights_pos(nvars, characteristic_fluxes, combined_frozen_metrics, w) +subroutine weno_weights_pos(nvars, characteristic_fluxes, combined_frozen_metric, w) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(out) :: w(nvars, 3) real*8 :: C(3) = (/0.1d0, 0.6d0, 0.3d0/) @@ -607,7 +607,7 @@ subroutine weno_weights_pos(nvars, characteristic_fluxes, combined_frozen_metric integer p, i, j real*8 sum_alpha(1) - eps = 1.0d-6!*combined_frozen_metrics + eps = 1.0d-6!*combined_frozen_metric p = 2 call oscillation_pos(nvars, characteristic_fluxes, IS) @@ -671,36 +671,36 @@ subroutine flux_differences_pos(nvars, characteristic_fluxes, flux_differences) end subroutine subroutine dissipation_part_neg(nvars, characteristic_fluxes, & - combined_frozen_metrics, R, dissipation_neg) + combined_frozen_metric, R, dissipation_neg) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(in) :: R(nvars, nvars) real*8, intent(out) :: dissipation_neg(nvars) real*8 combined_fluxes(nvars) - call weno_combination_neg(nvars, characteristic_fluxes, combined_frozen_metrics, combined_fluxes) + call weno_combination_neg(nvars, characteristic_fluxes, combined_frozen_metric, combined_fluxes) call mult_mat_vec(nvars, nvars, 1.0d0/60, R, combined_fluxes, dissipation_neg) end subroutine subroutine weno_combination_neg(nvars, characteristic_fluxes, & - combined_frozen_metrics, combined_fluxes) + combined_frozen_metric, combined_fluxes) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(out) :: combined_fluxes(nvars) real*8 w(nvars, 3) real*8 flux_differences(nvars, 3) integer v - call weno_weights_neg(nvars, characteristic_fluxes, combined_frozen_metrics, w) + call weno_weights_neg(nvars, characteristic_fluxes, combined_frozen_metric, w) call flux_differences_neg(nvars, characteristic_fluxes, flux_differences) do v=1,nvars @@ -710,12 +710,12 @@ subroutine weno_combination_neg(nvars, characteristic_fluxes, & end do end subroutine -subroutine weno_weights_neg(nvars, characteristic_fluxes, combined_frozen_metrics, w) +subroutine weno_weights_neg(nvars, characteristic_fluxes, combined_frozen_metric, w) implicit none integer, intent(in) :: nvars real*8, intent(in) :: characteristic_fluxes(nvars, -2:3) - real*8, intent(in) :: combined_frozen_metrics + real*8, intent(in) :: combined_frozen_metric real*8, intent(out) :: w(nvars, 3) real*8 :: C(3) = (/0.1d0, 0.6d0, 0.3d0/) @@ -724,7 +724,7 @@ subroutine weno_weights_neg(nvars, characteristic_fluxes, combined_frozen_metric integer p, i, j real*8 sum_alpha(1) - eps = 1.0d-6!*combined_frozen_metrics + eps = 1.0d-6!*combined_frozen_metric p = 2 call oscillation_neg(nvars, characteristic_fluxes, IS) diff --git a/test/test_flux_window_ops.py b/test/test_flux_window_ops.py index 8387bda..0f4742a 100644 --- a/test/test_flux_window_ops.py +++ b/test/test_flux_window_ops.py @@ -136,7 +136,7 @@ def test_weno_flux(queue, window_results): generalized_fluxes=data.fluxes, characteristic_fluxes_pos=data.char_fluxes_pos, characteristic_fluxes_neg=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metric, + combined_frozen_metric=data.combined_frozen_metric, R=data.R, flux=flux_dev) @@ -168,7 +168,7 @@ def test_dissipation_part_pos(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_pos, - combined_frozen_metrics=data.combined_frozen_metric, + combined_frozen_metric=data.combined_frozen_metric, R=data.R, dissipation_pos=dissipation_dev) @@ -185,7 +185,7 @@ def test_dissipation_part_neg(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metric, + combined_frozen_metric=data.combined_frozen_metric, R=data.R, dissipation_neg=dissipation_dev) @@ -202,7 +202,7 @@ def test_weno_weights_pos(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_pos, - combined_frozen_metrics=data.combined_frozen_metric, + combined_frozen_metric=data.combined_frozen_metric, w=weights_dev) sum_weights = np.sum(weights_dev.get(), axis=1) @@ -221,7 +221,7 @@ def test_weno_weights_neg(queue, window_results): prg(queue, nvars=data.nvars, characteristic_fluxes=data.char_fluxes_neg, - combined_frozen_metrics=data.combined_frozen_metric, + combined_frozen_metric=data.combined_frozen_metric, w=weights_dev) sum_weights = np.sum(weights_dev.get(), axis=1) -- GitLab