23 #ifdef DEAL_II_WITH_GSL
24 # include <gsl/gsl_sf_legendre.h>
35 <<
"x[" << arg1 <<
"] = " << arg2 <<
" is not in [0,1]");
45 #ifdef DEAL_II_WITH_GSL
47 for (
unsigned int d = 0;
d < dim;
d++)
49 const double x = 2.0 * (x_q[
d] - 0.5);
50 Assert((x_q[
d] <= 1.0) && (x_q[
d] >= 0.), ExcLegendre(
d, x_q[
d]));
51 const int ind = indices[
d];
52 res *=
std::sqrt(2.0) * gsl_sf_legendre_Pl(ind, x);
61 ExcMessage(
"deal.II has to be configured with GSL "
62 "in order to use Legendre transformation."));
77 for (
unsigned int d = 0;
d < dim;
d++)
78 res *= (0.5 + indices[
d]);
85 template <
int dim,
int spacedim>
90 const unsigned int dof)
93 for (
unsigned int q = 0; q < quadrature.
size(); ++q)
99 return sum * multiplier(indices);
108 template <
int spacedim>
111 const std::vector<unsigned int> & n_coefficients_per_direction,
114 const unsigned int fe,
119 if (legendre_transform_matrices[fe].m() == 0)
121 legendre_transform_matrices[fe].reinit(n_coefficients_per_direction[fe],
122 fe_collection[fe].dofs_per_cell);
124 for (
unsigned int k = 0; k < n_coefficients_per_direction[fe]; ++k)
125 for (
unsigned int j = 0; j < fe_collection[fe].dofs_per_cell; ++j)
126 legendre_transform_matrices[fe](k, j) = integrate(
131 template <
int spacedim>
134 const std::vector<unsigned int> & n_coefficients_per_direction,
137 const unsigned int fe,
142 if (legendre_transform_matrices[fe].m() == 0)
144 legendre_transform_matrices[fe].reinit(
145 Utilities::fixed_power<2>(n_coefficients_per_direction[fe]),
146 fe_collection[fe].dofs_per_cell);
149 for (
unsigned int k1 = 0; k1 < n_coefficients_per_direction[fe]; ++k1)
150 for (
unsigned int k2 = 0; k2 < n_coefficients_per_direction[fe];
152 for (
unsigned int j = 0; j < fe_collection[fe].dofs_per_cell; ++j)
153 legendre_transform_matrices[fe](k, j) =
154 integrate(fe_collection[fe],
161 template <
int spacedim>
164 const std::vector<unsigned int> & n_coefficients_per_direction,
167 const unsigned int fe,
172 if (legendre_transform_matrices[fe].m() == 0)
174 legendre_transform_matrices[fe].reinit(
175 Utilities::fixed_power<3>(n_coefficients_per_direction[fe]),
176 fe_collection[fe].dofs_per_cell);
179 for (
unsigned int k1 = 0; k1 < n_coefficients_per_direction[fe]; ++k1)
180 for (
unsigned int k2 = 0; k2 < n_coefficients_per_direction[fe]; ++k2)
181 for (
unsigned int k3 = 0; k3 < n_coefficients_per_direction[fe];
183 for (
unsigned int j = 0; j < fe_collection[fe].dofs_per_cell; ++j)
184 legendre_transform_matrices[fe](k, j) =
185 integrate(fe_collection[fe],
197 template <
int dim,
int spacedim>
199 const std::vector<unsigned int> & n_coefficients_per_direction,
202 : n_coefficients_per_direction(n_coefficients_per_direction)
203 , fe_collection(&fe_collection)
204 , q_collection(q_collection)
205 , legendre_transform_matrices(fe_collection.size())
209 ExcMessage(
"All parameters are supposed to have the same size."));
212 const unsigned int max_n_coefficients_per_direction =
216 Utilities::fixed_power<dim>(max_n_coefficients_per_direction));
221 template <
int dim,
int spacedim>
223 const unsigned int n_coefficients_per_direction,
227 std::vector<unsigned
int>(fe_collection.size(),
228 n_coefficients_per_direction),
235 template <
int dim,
int spacedim>
249 template <
int dim,
int spacedim>
254 for (
unsigned int fe = 0; fe < fe_collection->
size(); ++fe)
256 ensure_existence(n_coefficients_per_direction,
260 legendre_transform_matrices);
268 template <
int dim,
int spacedim>
271 const unsigned int index)
const
273 return n_coefficients_per_direction[index];
278 template <
int dim,
int spacedim>
279 template <
typename Number>
282 const ::Vector<Number> &local_dof_values,
283 const unsigned int cell_active_fe_index,
286 for (
unsigned int d = 0;
d < dim; ++
d)
288 n_coefficients_per_direction[cell_active_fe_index]);
290 ensure_existence(n_coefficients_per_direction,
293 cell_active_fe_index,
294 legendre_transform_matrices);
297 legendre_transform_matrices[cell_active_fe_index];
299 unrolled_coefficients.resize(Utilities::fixed_power<dim>(
300 n_coefficients_per_direction[cell_active_fe_index]));
301 std::fill(unrolled_coefficients.begin(),
302 unrolled_coefficients.end(),
310 for (
unsigned int i = 0; i < unrolled_coefficients.size(); i++)
311 for (
unsigned int j = 0; j < local_dof_values.size(); j++)
312 unrolled_coefficients[i] +=
matrix[i][j] * local_dof_values[j];
314 legendre_coefficients.
fill(unrolled_coefficients.begin());
320 #include "fe_series_legendre.inst"