Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #include "PiecewiseBilinear.h"
11 : #include "ColumnMajorMatrix.h"
12 : #include "BilinearInterpolation.h"
13 : #include "MooseUtils.h"
14 :
15 : #include <fstream>
16 :
17 : registerMooseObject("MooseApp", PiecewiseBilinear);
18 :
19 : InputParameters
20 14643 : PiecewiseBilinear::validParams()
21 : {
22 14643 : InputParameters params = Function::validParams();
23 58572 : params.addParam<FileName>(
24 : "data_file", "", "File holding csv data for use with PiecewiseBilinear");
25 58572 : params.addParam<std::vector<Real>>("x", "The x abscissa values");
26 58572 : params.addParam<std::vector<Real>>("y", "The y abscissa values");
27 58572 : params.addParam<std::vector<Real>>("z", "The ordinate values");
28 58572 : params.addParam<int>("axis", -1, "The axis used (0, 1, or 2 for x, y, or z).");
29 43929 : params.addParam<int>(
30 29286 : "xaxis", -1, "The coordinate used for x-axis data (0, 1, or 2 for x, y, or z).");
31 43929 : params.addParam<int>(
32 29286 : "yaxis", -1, "The coordinate used for y-axis data (0, 1, or 2 for x, y, or z).");
33 43929 : params.addParam<Real>(
34 29286 : "scale_factor", 1.0, "Scale factor to be applied to the axis, yaxis, or xaxis values");
35 43929 : params.addParam<bool>("radial",
36 29286 : false,
37 : "Set to true if you want to interpolate along a radius "
38 : "rather that along a specific axis, and note that you "
39 : "have to define xaxis and yaxis in the input file");
40 14643 : params.addClassDescription("Interpolates values from a csv file");
41 14643 : return params;
42 0 : }
43 :
44 180 : PiecewiseBilinear::PiecewiseBilinear(const InputParameters & parameters)
45 : : Function(parameters),
46 180 : _xaxis(getParam<int>("xaxis")),
47 360 : _yaxis(getParam<int>("yaxis")),
48 360 : _data_file_name(getParam<FileName>("data_file")),
49 360 : _axis(getParam<int>("axis")),
50 180 : _axisValid(_axis > -1 && _axis < 3),
51 180 : _yaxisValid(_yaxis > -1 && _yaxis < 3),
52 180 : _xaxisValid(_xaxis > -1 && _xaxis < 3),
53 360 : _scale_factor(getParam<Real>("scale_factor")),
54 720 : _radial(getParam<bool>("radial"))
55 : {
56 :
57 180 : if (!_axisValid && !_yaxisValid && !_xaxisValid)
58 0 : mooseError("In PiecewiseBilinear ",
59 0 : _name,
60 : ": None of axis, yaxis, or xaxis properly defined. Allowable range is 0-2");
61 :
62 180 : if (_axisValid && (_yaxisValid || _xaxisValid))
63 0 : mooseError("In PiecewiseBilinear ", _name, ": Cannot define axis with either yaxis or xaxis");
64 :
65 180 : if (_radial && (!_yaxisValid || !_xaxisValid))
66 0 : mooseError(
67 0 : "In PiecewiseBilinear ", _name, ": yaxis and xaxis must be defined when radial = true");
68 :
69 180 : std::vector<Real> x;
70 180 : std::vector<Real> y;
71 180 : ColumnMajorMatrix z;
72 180 : std::vector<Real> z_vec;
73 :
74 180 : if (!_data_file_name.empty())
75 : {
76 604 : if (parameters.isParamValid("x") || parameters.isParamValid("y") ||
77 256 : parameters.isParamValid("z"))
78 4 : mooseError("In PiecewiseBilinear: Cannot specify 'data_file' and 'x', 'y', or 'z' together.");
79 : else
80 84 : parse(_data_file_name, x, y, z, name());
81 : }
82 :
83 632 : else if (!(parameters.isParamValid("x") && parameters.isParamValid("y") &&
84 268 : parameters.isParamValid("z")))
85 4 : mooseError("In PiecewiseBilinear: 'x' and 'y' and 'z' must be specified if any one is "
86 : "specified or no 'data_file' is specified.");
87 :
88 : else
89 : {
90 176 : x = getParam<std::vector<Real>>("x");
91 176 : y = getParam<std::vector<Real>>("y");
92 176 : z_vec = getParam<std::vector<Real>>("z");
93 :
94 : // check that size of z = (size of x)*(size of y)
95 88 : if (z_vec.size() != x.size() * y.size())
96 4 : mooseError("In PiecewiseBilinear: Size of z should be the size of x times the size of y.");
97 :
98 : // reshape and populate z matrix
99 84 : z.reshape(y.size(), x.size());
100 84 : int idx = 0;
101 336 : for (unsigned int i = 0; i < y.size(); i++)
102 1008 : for (unsigned int j = 0; j < x.size(); j++)
103 : {
104 756 : z(i, j) = z_vec[idx];
105 756 : idx += 1;
106 : }
107 : }
108 :
109 168 : _bilinear_interp = std::make_unique<BilinearInterpolation>(x, y, z);
110 168 : }
111 :
112 168 : PiecewiseBilinear::~PiecewiseBilinear() {}
113 :
114 : Real
115 3648 : PiecewiseBilinear::value(Real t, const Point & p) const
116 : {
117 3648 : return valueInternal(t, p);
118 : }
119 :
120 : ADReal
121 0 : PiecewiseBilinear::value(const ADReal & t, const ADPoint & p) const
122 : {
123 0 : return valueInternal(t, p);
124 : }
125 :
126 : template <typename T, typename P>
127 : T
128 3648 : PiecewiseBilinear::valueInternal(T t, const P & p) const
129 : {
130 3648 : T retVal = 0.0;
131 3648 : if (_yaxisValid && _xaxisValid && _radial)
132 : {
133 608 : const auto rx = p(_xaxis) * p(_xaxis);
134 608 : const auto ry = p(_yaxis) * p(_yaxis);
135 : using std::sqrt;
136 608 : const auto r = sqrt(rx + ry);
137 608 : retVal = _bilinear_interp->sample(r, t);
138 608 : }
139 3040 : else if (_axisValid)
140 1216 : retVal = _bilinear_interp->sample(p(_axis), t);
141 1824 : else if (_yaxisValid && !_radial)
142 : {
143 1216 : if (_xaxisValid)
144 608 : retVal = _bilinear_interp->sample(p(_xaxis), p(_yaxis));
145 : else
146 608 : retVal = _bilinear_interp->sample(t, p(_yaxis));
147 : }
148 : else
149 608 : retVal = _bilinear_interp->sample(p(_xaxis), t);
150 :
151 3648 : return retVal * _scale_factor;
152 0 : }
153 :
154 : void
155 84 : PiecewiseBilinear::parse(const std::string & data_file_name,
156 : std::vector<Real> & x,
157 : std::vector<Real> & y,
158 : ColumnMajorMatrix & z,
159 : const std::string & object_name)
160 : {
161 84 : std::ifstream file(data_file_name.c_str());
162 84 : if (!file.good())
163 0 : ::mooseError(object_name, " : Error opening file '", data_file_name, "'.");
164 :
165 84 : std::size_t num_lines = 0;
166 84 : std::size_t num_cols = libMesh::invalid_uint;
167 84 : std::vector<Real> data;
168 :
169 84 : std::string line;
170 84 : std::vector<Real> line_data;
171 420 : while (std::getline(file, line))
172 : {
173 336 : num_lines++;
174 672 : if (!MooseUtils::tokenizeAndConvert<double>(line, line_data, ", "))
175 0 : ::mooseError(object_name, " : Error parsing file '", data_file_name, "' on line ", num_lines);
176 :
177 336 : data.insert(data.end(), line_data.begin(), line_data.end());
178 :
179 336 : if (num_cols == libMesh::invalid_uint)
180 84 : num_cols = line_data.size();
181 252 : else if (line_data.size() != num_cols + 1)
182 0 : ::mooseError(object_name,
183 : " : Read ",
184 0 : line_data.size(),
185 : " columns of data but expected ",
186 0 : num_cols + 1,
187 : " columns in line ",
188 : num_lines);
189 : }
190 :
191 84 : x.resize(num_cols);
192 84 : y.resize(num_lines - 1);
193 84 : z.reshape(num_lines - 1, num_cols);
194 84 : std::size_t offset = 0;
195 :
196 : // Extract the first line's data (the x axis data)
197 336 : for (unsigned int j = 0; j < num_cols; ++j)
198 252 : x[j] = data[offset++];
199 :
200 336 : for (unsigned int i = 0; i < num_lines - 1; ++i)
201 : {
202 : // Extract the y axis entry for this line
203 252 : y[i] = data[offset++];
204 :
205 : // Extract the function values for this row in the matrix
206 1008 : for (unsigned int j = 0; j < num_cols; ++j)
207 756 : z(i, j) = data[offset++];
208 : }
209 :
210 84 : if (data.size() != offset)
211 0 : ::mooseError(object_name, " : Inconsistency in data read from '", data_file_name, "'.");
212 84 : }
|