From 0e540b48cdad4e682b9516dcb393a8532c73793d Mon Sep 17 00:00:00 2001 From: Maciej Waruszewski Date: Tue, 28 Nov 2023 12:50:38 -0700 Subject: [PATCH] Exact inverse in 3d --- dynamics/spam/src/hamiltonians/variableset.h | 30 +++++++++++++++++--- 1 file changed, 26 insertions(+), 4 deletions(-) diff --git a/dynamics/spam/src/hamiltonians/variableset.h b/dynamics/spam/src/hamiltonians/variableset.h index 9bd9f6df..15b5ecf5 100644 --- a/dynamics/spam/src/hamiltonians/variableset.h +++ b/dynamics/spam/src/hamiltonians/variableset.h @@ -227,10 +227,6 @@ template class VariableSetBase : public T { "The number of crm cells in the horizontal " "has to be odd when using the couple_wind_exact_inverse option"); } - if (ndims > 1) { - throw std::runtime_error( - "couple_wind_exact_inverse not implemented in 3d"); - } } varset.thermo = thermodynamics; varset.reference_state = refstate; @@ -823,6 +819,32 @@ void convert_coupler_to_dynamics_wind(const VariableSetBase &varset, i + pis, n); } }); + if (ndims > 1) { + parallel_for( + "Coupler to Dynamics State Primal V", + SimpleBounds<3>(primal_topology.ni, primal_topology.n_cells_x, + primal_topology.nens), + YAKL_LAMBDA(int k, int i, int n) { + real x0 = 0; + for (int j = 0; j < primal_topology.n_cells_y; ++j) { + x0 += (j % 2 == 0 ? 1 : -1) * dm_vvel(k, j, i, n); + } + prog_vars.fields_arr[VVAR].data(1, k + pks, pjs, i + pis, n) = x0; + prog_vars.fields_arr[VVAR].data(1, k + pks, pjs, i + pis, n) *= + primal_geometry.get_area_10entity(1, k + pks, pjs, i + pis, + n); + + for (int j = 1; j < primal_topology.n_cells_y; ++j) { + x0 = 2 * dm_vvel(k, j - 1, i, n) - x0; + prog_vars.fields_arr[VVAR].data(1, k + pks, j + pjs, i + pis, + n) = x0; + prog_vars.fields_arr[VVAR].data(1, k + pks, j + pjs, i + pis, + n) *= + primal_geometry.get_area_10entity(1, k + pks, j + pjs, + i + pis, n); + } + }); + } parallel_for( "Coupler to Dynamics State Primal W", SimpleBounds<3>(primal_topology.n_cells_y, primal_topology.n_cells_x,