Skip to content

Commit 1c8640b

Browse files
mbkuhnmbkuhn
authored
Upgrade field norm capabilities - vector magnitude, different norm types (#1704)
* first commit, existing tests passing * use the default better * introduce unit test for different norm types * remove print statement * set up test for vector magnitude functionality * add to documentation * fix docs warning
1 parent 8c1a6a3 commit 1c8640b

File tree

6 files changed

+341
-48
lines changed

6 files changed

+341
-48
lines changed

‎amr-wind/utilities/sampling/FieldNorms.H‎

Lines changed: 12 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -34,9 +34,12 @@ public:
3434

3535
void post_regrid_actions() override {}
3636

37-
//! calculate the L2 norm of a given field and component
38-
static amrex::Real
39-
l2_norm(const amr_wind::Field& field, const int comp, const bool use_mask);
37+
static amrex::Real get_norm(
38+
const amr_wind::Field& field,
39+
const int comp,
40+
const int ncomp,
41+
const int norm_type,
42+
const bool use_mask);
4043

4144
const amrex::Vector<std::string>& var_names() const { return m_var_names; }
4245

@@ -77,8 +80,14 @@ private:
7780
//! precision in ASCII output
7881
int m_precision{12};
7982

83+
//! Type of norm
84+
int m_norm_type{2};
85+
8086
//! Mask redundant grids (finest available grid only)
8187
bool m_use_mask{true};
88+
89+
//! Consider vector magnitude, not separate vector components
90+
bool m_use_vector_magnitude{false};
8291
};
8392

8493
} // namespace amr_wind::field_norms

‎amr-wind/utilities/sampling/FieldNorms.cpp‎

Lines changed: 121 additions & 42 deletions
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@
55
#include <utility>
66
#include "AMReX_ParmParse.H"
77
#include "amr-wind/utilities/IOManager.H"
8+
#include "amr-wind/utilities/constants.H"
89

910
namespace amr_wind::field_norms {
1011

@@ -22,20 +23,42 @@ void FieldNorms::initialize()
2223
amrex::ParmParse pp(m_label);
2324
populate_output_parameters(pp);
2425
pp.query("mask_redundant_grids", m_use_mask);
26+
pp.query("use_vector_magnitude", m_use_vector_magnitude);
27+
std::string norm_type(std::to_string(m_norm_type));
28+
pp.query("norm_type", norm_type);
29+
if (norm_type == "2") {
30+
m_norm_type = 2;
31+
} else if (norm_type == "1") {
32+
m_norm_type = 1;
33+
} else if (amrex::toLower(norm_type) == "infinity") {
34+
m_norm_type = -1;
35+
} else {
36+
amrex::Abort(
37+
"FieldNorms: norm_type not recognized. Valid options are 1, 2, "
38+
"and infinity.");
39+
}
2540
}
2641

2742
const auto& io_mng = m_sim.io_manager();
2843
for (const auto& fld : io_mng.plot_fields()) {
29-
ioutils::add_var_names(m_var_names, fld->name(), fld->num_comp());
44+
if (m_use_vector_magnitude) {
45+
m_var_names.push_back(fld->name());
46+
} else {
47+
ioutils::add_var_names(m_var_names, fld->name(), fld->num_comp());
48+
}
3049
}
3150

3251
m_fnorms.resize(m_var_names.size(), 0.0);
3352

3453
prepare_ascii_file();
3554
}
3655

37-
amrex::Real FieldNorms::l2_norm(
38-
const amr_wind::Field& field, const int comp, const bool use_mask)
56+
amrex::Real FieldNorms::get_norm(
57+
const amr_wind::Field& field,
58+
const int comp,
59+
const int ncomp,
60+
const int norm_type,
61+
const bool use_mask)
3962
{
4063
amrex::Real nrm = 0.0;
4164

@@ -81,55 +104,111 @@ amrex::Real FieldNorms::l2_norm(
81104
level_mask.setVal(1.);
82105
}
83106

84-
nrm += amrex::ReduceSum(
85-
level_mask, field(lev), nghost,
86-
[=] AMREX_GPU_HOST_DEVICE(
87-
amrex::Box const& bx,
88-
amrex::Array4<amrex::Real const> const& mask_arr,
89-
amrex::Array4<amrex::Real const> const& field_arr)
90-
-> amrex::Real {
91-
amrex::Real nrm_fab = 0.0;
92-
93-
const auto& fbx =
94-
it_sum == 1
95-
? amrex::surroundingNodes(bx, node_dir)
96-
: (it_sum > 1 ? amrex::surroundingNodes(bx) : bx);
97-
98-
amrex::Loop(fbx, [=, &nrm_fab](int i, int j, int k) noexcept {
99-
const amrex::IntVect iv{i, j, k};
100-
amrex::IntVect iv_cc = iv;
101-
// adjust volume for different data locations
102-
amrex::Real data_vol = cell_vol;
103-
for (int nn = 0; nn < AMREX_SPACEDIM; ++nn) {
104-
const bool at_node_dir =
105-
index_type[nn] == 1 && (iv[nn] == fbx.bigEnd(nn) ||
106-
iv[nn] == fbx.smallEnd(nn));
107-
data_vol *= at_node_dir ? 0.5 : 1.0;
108-
// limit mask array indices to cell-centered
109-
iv_cc[nn] = amrex::min(
110-
bx.bigEnd(nn),
111-
amrex::max(bx.smallEnd(nn), iv_cc[nn]));
112-
}
113-
nrm_fab += data_vol * field_arr(iv, comp) *
114-
field_arr(iv, comp) * mask_arr(iv_cc);
107+
if (norm_type > 0) {
108+
nrm += amrex::ReduceSum(
109+
level_mask, field(lev), nghost,
110+
[=] AMREX_GPU_HOST_DEVICE(
111+
amrex::Box const& bx,
112+
amrex::Array4<amrex::Real const> const& mask_arr,
113+
amrex::Array4<amrex::Real const> const& field_arr)
114+
-> amrex::Real {
115+
amrex::Real nrm_fab = 0.0;
116+
117+
const auto& fbx =
118+
it_sum == 1
119+
? amrex::surroundingNodes(bx, node_dir)
120+
: (it_sum > 1 ? amrex::surroundingNodes(bx) : bx);
121+
122+
amrex::Loop(
123+
fbx, [=, &nrm_fab](int i, int j, int k) noexcept {
124+
const amrex::IntVect iv{i, j, k};
125+
amrex::IntVect iv_cc = iv;
126+
// adjust volume for different data locations
127+
amrex::Real data_vol = cell_vol;
128+
for (int nn = 0; nn < AMREX_SPACEDIM; ++nn) {
129+
const bool at_node_dir =
130+
index_type[nn] == 1 &&
131+
(iv[nn] == fbx.bigEnd(nn) ||
132+
iv[nn] == fbx.smallEnd(nn));
133+
data_vol *= at_node_dir ? 0.5 : 1.0;
134+
// limit mask array indices to cell-centered
135+
iv_cc[nn] = amrex::min(
136+
bx.bigEnd(nn),
137+
amrex::max(bx.smallEnd(nn), iv_cc[nn]));
138+
}
139+
amrex::Real fval =
140+
std::abs(field_arr(i, j, k, comp));
141+
// Calculate magnitude if requested
142+
for (int n = 1; n < ncomp; ++n) {
143+
fval *= fval;
144+
fval += field_arr(i, j, k, n) *
145+
field_arr(i, j, k, n);
146+
fval = std::sqrt(fval);
147+
}
148+
fval *= norm_type == 2 ? fval : 1.;
149+
nrm_fab += data_vol * fval * mask_arr(iv_cc);
150+
});
151+
return nrm_fab;
152+
});
153+
} else {
154+
// Can directly use the field box for the infinity norm
155+
const amrex::Real nrm_lev = amrex::ReduceMax(
156+
field(lev), level_mask, nghost,
157+
[=] AMREX_GPU_HOST_DEVICE(
158+
amrex::Box const& bx,
159+
amrex::Array4<amrex::Real const> const& field_arr,
160+
amrex::Array4<amrex::Real const> const& mask_arr)
161+
-> amrex::Real {
162+
amrex::Real nrm_fab = 0.0;
163+
164+
amrex::Loop(
165+
bx, [=, &nrm_fab](int i, int j, int k) noexcept {
166+
amrex::Real fval =
167+
std::abs(field_arr(i, j, k, comp));
168+
// Calculate magnitude if requested
169+
for (int n = 1; n < ncomp; ++n) {
170+
fval *= fval;
171+
fval += field_arr(i, j, k, n) *
172+
field_arr(i, j, k, n);
173+
fval = std::sqrt(fval);
174+
}
175+
if (std::abs(mask_arr(i, j, k) - 1.) <
176+
constants::TIGHT_TOL) {
177+
nrm_fab = std::max(nrm_fab, fval);
178+
}
179+
});
180+
return nrm_fab;
115181
});
116-
return nrm_fab;
117-
});
182+
nrm = amrex::max(nrm_lev, nrm);
183+
}
118184
}
119185

120-
amrex::ParallelDescriptor::ReduceRealSum(nrm);
121-
122-
const amrex::Real total_volume = geom[0].ProbDomain().volume();
186+
if (norm_type > 0) {
187+
amrex::ParallelDescriptor::ReduceRealSum(nrm);
188+
const amrex::Real total_volume = geom[0].ProbDomain().volume();
189+
nrm /= total_volume;
190+
if (norm_type == 2) {
191+
nrm = std::sqrt(nrm);
192+
}
193+
} else {
194+
amrex::ParallelDescriptor::ReduceRealMax(nrm);
195+
}
123196

124-
return std::sqrt(nrm / total_volume);
197+
return nrm;
125198
}
126199

127200
void FieldNorms::process_field_norms()
128201
{
129202
int ind = 0;
130203
for (const auto& fld : m_sim.io_manager().plot_fields()) {
131-
for (int comp = 0; comp < fld->num_comp(); ++comp) {
132-
m_fnorms[ind++] = l2_norm(*fld, comp, m_use_mask);
204+
if (m_use_vector_magnitude) {
205+
m_fnorms[ind++] =
206+
get_norm(*fld, 0, fld->num_comp(), m_norm_type, m_use_mask);
207+
} else {
208+
for (int comp = 0; comp < fld->num_comp(); ++comp) {
209+
m_fnorms[ind++] =
210+
get_norm(*fld, comp, 1, m_norm_type, m_use_mask);
211+
}
133212
}
134213
}
135214
}

‎docs/sphinx/user/inputs.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -101,6 +101,7 @@ documentation provided here might not work with older releases.
101101
inputs_Averaging.rst
102102
inputs_KineticEnergy.rst
103103
inputs_Enstrophy.rst
104+
inputs_FieldNorms.rst
104105
inputs_Actuator.rst
105106
inputs_multiphase.rst
106107
inputs_ocean_waves.rst

‎docs/sphinx/user/inputs_FieldNorms.rst‎

Lines changed: 44 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,44 @@
1+
.. _inputs_fieldnorms:
2+
3+
Section: FieldNorms
4+
~~~~~~~~~~~~~~~~~~~
5+
6+
This section controls field norms post-processing. Field norms
7+
output the global norm of every plot variable to a text file, writing
8+
a new line for each output time.
9+
The prefix is the label set in ``incflo.post_processing``. For example
10+
``incflo.post_processing = fieldnorms``. Inputs listed here only address
11+
options specific to field norms; the inputs controlling the output timing
12+
(which are shared by other post-processing types) are listed in the
13+
[:ref:`post-processing section <inputs_post_processing>`].
14+
15+
.. input_param:: fieldnorms.type
16+
17+
**type:** String, mandatory
18+
19+
To use field norms output specify with keyword ``FieldNorms``
20+
21+
.. input_param:: fieldnorms.mask_redundant_grids
22+
23+
**type:** Boolean, optional, default = true
24+
25+
Setting this option to true ensures that the field norm calculation
26+
does not include cells that are covered by finer cells from a different
27+
mesh level. This also ensures that no volume is double counted. Setting
28+
it to false makes sure that every cell, covered or not, is included.
29+
30+
.. input_param:: fieldnorms.use_vector_magnitude
31+
32+
**type:** Boolean, optional, default = false
33+
34+
By default, each component of vector output fields is considered separately.
35+
Setting this option to true ensures that the field norm calculation
36+
considers the local vector magnitude instead of separating each vector component.
37+
38+
.. input_param:: fieldnorms.norm_type
39+
40+
**type:** String, optional, default = ``2``
41+
42+
By default, this post-processing routine calculates the volume-weighted
43+
L2 norm, corresponding to norm type 2. Other options include the L1 norm,
44+
designated by ``1``, and the L-infinity norm, designated by ``infinity``.

0 commit comments

Comments
 (0)