16 #ifndef dealii_fe_series_H 17 #define dealii_fe_series_H 21 #include <deal.II/base/config.h> 23 #include <deal.II/base/exceptions.h> 24 #include <deal.II/base/subscriptor.h> 25 #include <deal.II/base/table.h> 26 #include <deal.II/base/table_indices.h> 27 #include <deal.II/base/tensor.h> 29 #include <deal.II/hp/fe_collection.h> 30 #include <deal.II/hp/q_collection.h> 32 #include <deal.II/lac/full_matrix.h> 33 #include <deal.II/lac/vector.h> 35 #include <deal.II/numerics/vector_tools.h> 42 DEAL_II_NAMESPACE_OPEN
101 Fourier(
const unsigned int size_in_each_direction,
111 calculate(const ::Vector<double> & local_dof_values,
112 const unsigned int cell_active_fe_index,
113 Table<dim, std::complex<double>> &fourier_coefficients);
203 Legendre(
const unsigned int size_in_each_direction,
213 calculate(const ::Vector<double> &local_dof_values,
214 const unsigned int cell_active_fe_index,
221 const unsigned int N;
266 template <
int dim,
typename T>
267 std::pair<std::vector<unsigned int>, std::vector<double>>
269 const std::function<std::pair<bool, unsigned int>(
280 std::pair<double, double>
293 namespace FESeriesImplementation
295 template <
int dim,
typename T>
299 const std::function<std::pair<bool, unsigned int>(
301 std::map<
unsigned int, std::vector<T>> &pred_to_values)
303 const std::pair<bool, unsigned int> pred_pair = predicate(ind);
305 if (pred_pair.first ==
false)
308 const unsigned int &pred_value = pred_pair.second;
309 const T & coeff_value = coefficients(ind);
312 pred_to_values[pred_value].push_back(coeff_value);
315 template <
typename T>
318 const std::function<std::pair<bool, unsigned int>(
320 std::map<
unsigned int, std::vector<T>> &pred_to_values)
322 for (
unsigned int i = 0; i < coefficients.
size(0); i++)
325 fill_map_index(coefficients, ind, predicate, pred_to_values);
329 template <
typename T>
332 const std::function<std::pair<bool, unsigned int>(
334 std::map<
unsigned int, std::vector<T>> &pred_to_values)
336 for (
unsigned int i = 0; i < coefficients.
size(0); i++)
337 for (
unsigned int j = 0; j < coefficients.
size(1); j++)
340 fill_map_index(coefficients, ind, predicate, pred_to_values);
344 template <
typename T>
347 const std::function<std::pair<bool, unsigned int>(
349 std::map<
unsigned int, std::vector<T>> &pred_to_values)
351 for (
unsigned int i = 0; i < coefficients.
size(0); i++)
352 for (
unsigned int j = 0; j < coefficients.
size(1); j++)
353 for (
unsigned int k = 0; k < coefficients.
size(2); k++)
356 fill_map_index(coefficients, ind, predicate, pred_to_values);
361 template <
typename T>
363 complex_mean_value(
const T &value)
368 template <
typename T>
370 complex_mean_value(
const std::complex<T> &value)
374 "FESeries::process_coefficients() can not be used with" 375 "complex-valued coefficients and VectorTools::mean norm."));
376 return std::abs(value);
382 template <
int dim,
typename T>
383 std::pair<std::vector<unsigned int>, std::vector<double>>
390 std::vector<unsigned int> predicate_values;
391 std::vector<double> norm_values;
396 std::map<unsigned int, std::vector<T>> pred_to_values;
397 internal::FESeriesImplementation::fill_map(coefficients,
402 for (
typename std::map<
unsigned int, std::vector<T>>::const_iterator it =
403 pred_to_values.begin();
404 it != pred_to_values.end();
407 predicate_values.push_back(it->first);
408 Vector<T> values(it->second.begin(), it->second.end());
414 norm_values.push_back(values.l2_norm());
419 norm_values.push_back(values.l1_norm());
424 norm_values.push_back(values.linfty_norm());
429 norm_values.push_back(
430 internal::FESeriesImplementation::complex_mean_value(
431 values.mean_value()));
440 return std::make_pair(predicate_values, norm_values);
446 DEAL_II_NAMESPACE_CLOSE
448 #endif // dealii_fe_series_H
Table< dim, Tensor< 1, dim > > k_vectors
void calculate(const ::Vector< double > &local_dof_values, const unsigned int cell_active_fe_index, Table< dim, std::complex< double >> &fourier_coefficients)
SmartPointer< const hp::FECollection< dim > > fe_collection
std::pair< double, double > linear_regression(const std::vector< double > &x, const std::vector< double > &y)
SmartPointer< const hp::QCollection< dim > > q_collection
#define AssertThrow(cond, exc)
Fourier(const unsigned int size_in_each_direction, const hp::FECollection< dim > &fe_collection, const hp::QCollection< dim > &q_collection)
std::pair< std::vector< unsigned int >, std::vector< double > > process_coefficients(const Table< dim, T > &coefficients, const std::function< std::pair< bool, unsigned int >(const TableIndices< dim > &)> &predicate, const VectorTools::NormType norm)
std::vector< std::complex< double > > unrolled_coefficients
SmartPointer< const hp::QCollection< dim > > q_collection
static::ExceptionBase & ExcMessage(std::string arg1)
size_type size(const unsigned int i) const
SmartPointer< const hp::FECollection< dim > > fe_collection
void ensure_existence(const unsigned int fe_index)
std::vector< double > unrolled_coefficients
static::ExceptionBase & ExcNotImplemented()
std::vector< FullMatrix< std::complex< double > > > fourier_transform_matrices
std::vector< FullMatrix< double > > legendre_transform_matrices