From: Denis Davydov Date: Mon, 22 Feb 2016 12:41:11 +0000 (+0100) Subject: extend FEFieldFunction to complex-valued vectors X-Git-Tag: v8.5.0-rc1~1297^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0dc4d47;p=dealii.git extend FEFieldFunction to complex-valued vectors --- diff --git a/include/deal.II/numerics/fe_field_function.h b/include/deal.II/numerics/fe_field_function.h index 65e04efe3e..e6beccc6e5 100644 --- a/include/deal.II/numerics/fe_field_function.h +++ b/include/deal.II/numerics/fe_field_function.h @@ -162,7 +162,7 @@ namespace Functions template , typename VectorType=Vector > - class FEFieldFunction : public Function + class FEFieldFunction : public Function { public: /** @@ -200,7 +200,7 @@ namespace Functions * information. */ virtual void vector_value (const Point &p, - Vector &values) const; + Vector &values) const; /** * Return the value of the function at the given point. Unless there is @@ -219,8 +219,8 @@ namespace Functions * See the section in the general documentation of this class for more * information. */ - virtual double value (const Point< dim > &p, - const unsigned int component = 0) const; + virtual typename VectorType::value_type value (const Point< dim > &p, + const unsigned int component = 0) const; /** * Set @p values to the point values of the specified component of the @@ -238,7 +238,7 @@ namespace Functions * information. */ virtual void value_list (const std::vector > &points, - std::vector< double > &values, + std::vector &values, const unsigned int component = 0) const; @@ -258,7 +258,7 @@ namespace Functions * information. */ virtual void vector_value_list (const std::vector > &points, - std::vector< Vector > &values) const; + std::vector > &values) const; /** * Return the gradient of all components of the function at the given @@ -277,7 +277,7 @@ namespace Functions */ virtual void vector_gradient (const Point< dim > &p, - std::vector< Tensor< 1, dim > > &gradients) const; + std::vector< Tensor< 1, dim,typename VectorType::value_type > > &gradients) const; /** * Return the gradient of the specified component of the function at the @@ -294,8 +294,8 @@ namespace Functions * See the section in the general documentation of this class for more * information. */ - virtual Tensor<1,dim> gradient(const Point< dim > &p, - const unsigned int component = 0)const; + virtual Tensor<1,dim,typename VectorType::value_type> gradient(const Point< dim > &p, + const unsigned int component = 0)const; /** * Return the gradient of all components of the function at all the given @@ -313,7 +313,7 @@ namespace Functions virtual void vector_gradient_list (const std::vector< Point< dim > > &p, std::vector< - std::vector< Tensor< 1, dim > > > &gradients) const; + std::vector< Tensor< 1, dim,typename VectorType::value_type > > > &gradients) const; /** * Return the gradient of the specified component of the function at all @@ -330,7 +330,7 @@ namespace Functions */ virtual void gradient_list (const std::vector< Point< dim > > &p, - std::vector< Tensor< 1, dim > > &gradients, + std::vector< Tensor< 1, dim,typename VectorType::value_type > > &gradients, const unsigned int component=0) const; @@ -345,7 +345,7 @@ namespace Functions * See the section in the general documentation of this class for more * information. */ - virtual double + virtual typename VectorType::value_type laplacian (const Point &p, const unsigned int component = 0) const; @@ -363,7 +363,7 @@ namespace Functions */ virtual void vector_laplacian (const Point &p, - Vector &values) const; + Vector &values) const; /** * Compute the Laplacian of one component at a set of points. @@ -378,7 +378,7 @@ namespace Functions */ virtual void laplacian_list (const std::vector > &points, - std::vector &values, + std::vector &values, const unsigned int component = 0) const; /** @@ -394,7 +394,7 @@ namespace Functions */ virtual void vector_laplacian_list (const std::vector > &points, - std::vector > &values) const; + std::vector > &values) const; /** * Create quadrature rules. This function groups the points into blocks diff --git a/include/deal.II/numerics/fe_field_function.templates.h b/include/deal.II/numerics/fe_field_function.templates.h index df5bbd45d5..025d9cf59f 100644 --- a/include/deal.II/numerics/fe_field_function.templates.h +++ b/include/deal.II/numerics/fe_field_function.templates.h @@ -37,7 +37,7 @@ namespace Functions const VectorType &myv, const Mapping &mymapping) : - Function(mydh.get_fe().n_components()), + Function(mydh.get_fe().n_components()), dh(&mydh, "FEFieldFunction"), data_vector(myv), mapping(mymapping), @@ -60,7 +60,7 @@ namespace Functions template void FEFieldFunction::vector_value (const Point &p, - Vector &values) const + Vector &values) const { Assert (values.size() == n_components, ExcDimensionMismatch(values.size(), n_components)); @@ -97,11 +97,11 @@ namespace Functions template - double + typename VectorType::value_type FEFieldFunction::value (const Point &p, const unsigned int comp) const { - Vector values(n_components); + Vector values(n_components); vector_value(p, values); return values(comp); } @@ -112,7 +112,7 @@ namespace Functions void FEFieldFunction:: vector_gradient (const Point &p, - std::vector > &gradients) const + std::vector > &gradients) const { Assert (gradients.size() == n_components, ExcDimensionMismatch(gradients.size(), n_components)); @@ -140,25 +140,18 @@ namespace Functions FEValues fe_v(mapping, cell->get_fe(), quad, update_gradients); fe_v.reinit(cell); - - // FIXME: we need a temp argument because get_function_values wants to put - // its data into an object storing the correct scalar type, but this - // function wants to return everything in a vector - std::vector< std::vector > > vgrads - (1, std::vector >(n_components) ); - fe_v.get_function_gradients(data_vector, vgrads); - gradients = std::vector >(vgrads[0].begin(), vgrads[0].end()); + fe_v.get_function_gradients(data_vector, gradients); } template - Tensor<1,dim> + Tensor<1,dim,typename VectorType::value_type> FEFieldFunction:: gradient (const Point &p, const unsigned int comp) const { - std::vector > grads(n_components); + std::vector > grads(n_components); vector_gradient(p, grads); return grads[comp]; } @@ -169,7 +162,7 @@ namespace Functions void FEFieldFunction:: vector_laplacian (const Point &p, - Vector &values) const + Vector &values) const { Assert (values.size() == n_components, ExcDimensionMismatch(values.size(), n_components)); @@ -206,11 +199,11 @@ namespace Functions template - double FEFieldFunction::laplacian + typename VectorType::value_type FEFieldFunction::laplacian (const Point &p, const unsigned int comp) const { - Vector lap(n_components); + Vector lap(n_components); vector_laplacian(p, lap); return lap[comp]; } @@ -223,7 +216,7 @@ namespace Functions void FEFieldFunction:: vector_value_list (const std::vector > &points, - std::vector< Vector > &values) const + std::vector< Vector > &values) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); @@ -267,12 +260,12 @@ namespace Functions void FEFieldFunction:: value_list (const std::vector > &points, - std::vector< double > &values, + std::vector< typename VectorType::value_type > &values, const unsigned int component) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); - std::vector< Vector > vvalues(points.size(), Vector(n_components)); + std::vector< Vector > vvalues(points.size(), Vector(n_components)); vector_value_list(points, vvalues); for (unsigned int q=0; q:: vector_gradient_list (const std::vector > &points, - std::vector > > &values) const + std::vector > > &values) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); @@ -332,13 +325,13 @@ namespace Functions void FEFieldFunction:: gradient_list (const std::vector > &points, - std::vector< Tensor<1,dim> > &values, + std::vector< Tensor<1,dim, typename VectorType::value_type> > &values, const unsigned int component) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); - std::vector< std::vector > > - vvalues(points.size(), std::vector >(n_components)); + std::vector< std::vector > > + vvalues(points.size(), std::vector >(n_components)); vector_gradient_list(points, vvalues); for (unsigned int q=0; q:: vector_laplacian_list (const std::vector > &points, - std::vector< Vector > &values) const + std::vector< Vector > &values) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); @@ -391,12 +384,12 @@ namespace Functions void FEFieldFunction:: laplacian_list (const std::vector > &points, - std::vector &values, + std::vector &values, const unsigned int component) const { Assert(points.size() == values.size(), ExcDimensionMismatch(points.size(), values.size())); - std::vector< Vector > vvalues(points.size(), Vector(n_components)); + std::vector< Vector > vvalues(points.size(), Vector(n_components)); vector_laplacian_list(points, vvalues); for (unsigned int q=0; q