From 5c470e5e39c3b7acf0f1339056721453db3ee853 Mon Sep 17 00:00:00 2001 From: guido Date: Mon, 9 Aug 2004 21:34:29 +0000 Subject: [PATCH] new get_function_grads git-svn-id: https://svn.dealii.org/trunk@9544 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_values.h | 279 ++++++++++-------- deal.II/deal.II/source/fe/fe_values.cc | 102 +++++++ .../deal.II/source/fe/fe_values.instance.h | 8 + 3 files changed, 265 insertions(+), 124 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index f93a216e03..bdb6c40c67 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -121,7 +121,7 @@ class FEValuesData * derivatives. */ typedef std::vector > > GradGradVector; - + /** * Store the values of the shape * functions at the quadrature @@ -366,6 +366,9 @@ class FEValuesBase : protected FEValuesData * Destructor. */ ~FEValuesBase (); + /// @name ShapeAccess Access to shape function values + //@{ + /** * Value of a shape function at a * quadrature point on the cell, @@ -433,6 +436,138 @@ class FEValuesBase : protected FEValuesData const unsigned int point_no, const unsigned int component) const; + /** + * Compute the gradient of the + * @p ith shape function at the + * @p j quadrature point with + * respect to real cell + * coordinates. If you want to + * get the derivative in one of + * the coordinate directions, use + * the appropriate function of + * the Tensor class to + * extract one component. Since + * only a reference to the + * gradient's value is returned, + * there should be no major + * performance drawback. + * + * If the shape function is + * vector-valued, then this + * returns the only non-zero + * component. If the shape + * function has more than one + * non-zero component (i.e. it is + * not primitive), then throw an + * exception of type + * ExcShapeFunctionNotPrimitive. In + * that case, use the + * shape_grad_component() + * function. + */ + const Tensor<1,dim> & + shape_grad (const unsigned int function, + const unsigned int quadrature_point) const; + + /** + * Return one vector component of + * the gradient of a shape function + * at a quadrature point. If the + * finite element is scalar, then + * only component zero is allowed + * and the return value equals + * that of the shape_grad() + * function. If the finite + * element is vector valued but + * all shape functions are + * primitive (i.e. they are + * non-zero in only one + * component), then the value + * returned by shape_grad() + * equals that of this function + * for exactly one + * component. This function is + * therefore only of greater + * interest if the shape function + * is not primitive, but then it + * is necessary since the other + * function cannot be used. + */ + Tensor<1,dim> + shape_grad_component (const unsigned int function_no, + const unsigned int point_no, + const unsigned int component) const; + + /** + * Second derivatives of + * the @p function_noth shape function at + * the @p point_noth quadrature point + * with respect to real cell + * coordinates. If you want to + * get the derivatives in one of + * the coordinate directions, use + * the appropriate function of + * the @p Tensor class to + * extract one component. Since + * only a reference to the + * derivative values is returned, + * there should be no major + * performance drawback. + * + * If the shape function is + * vector-valued, then this + * returns the only non-zero + * component. If the shape + * function has more than one + * non-zero component (i.e. it is + * not primitive), then throw an + * exception of type + * @p ExcShapeFunctionNotPrimitive. In + * that case, use the + * shape_grad_grad_component() + * function. + */ + const Tensor<2,dim> & + shape_2nd_derivative (const unsigned int function_no, + const unsigned int point_no) const; + + + /** + * Return one vector component of + * the gradient of a shape + * function at a quadrature + * point. If the finite element + * is scalar, then only component + * zero is allowed and the return + * value equals that of the + * @p shape_2nd_derivative + * function. If the finite + * element is vector valued but + * all shape functions are + * primitive (i.e. they are + * non-zero in only one + * component), then the value + * returned by + * @p shape_2nd_derivative + * equals that of this function + * for exactly one + * component. This function is + * therefore only of greater + * interest if the shape function + * is not primitive, but then it + * is necessary since the other + * function cannot be used. + */ + Tensor<2,dim> + shape_2nd_derivative_component (const unsigned int function_no, + const unsigned int point_no, + const unsigned int component) const; + + + //@} + /// @name FunctionAccess Access to values of global finite element functions + //@{ + /** * Returns the values of the * finite element function @@ -570,68 +705,6 @@ class FEValuesBase : protected FEValuesData const std::vector& indices, std::vector >& values) const; - /** - * Compute the gradient of the - * @p ith shape function at the - * @p j quadrature point with - * respect to real cell - * coordinates. If you want to - * get the derivative in one of - * the coordinate directions, use - * the appropriate function of - * the Tensor class to - * extract one component. Since - * only a reference to the - * gradient's value is returned, - * there should be no major - * performance drawback. - * - * If the shape function is - * vector-valued, then this - * returns the only non-zero - * component. If the shape - * function has more than one - * non-zero component (i.e. it is - * not primitive), then throw an - * exception of type - * ExcShapeFunctionNotPrimitive. In - * that case, use the - * shape_grad_component() - * function. - */ - const Tensor<1,dim> & - shape_grad (const unsigned int function, - const unsigned int quadrature_point) const; - - /** - * Return one vector component of - * the gradient of a shape function - * at a quadrature point. If the - * finite element is scalar, then - * only component zero is allowed - * and the return value equals - * that of the shape_grad() - * function. If the finite - * element is vector valued but - * all shape functions are - * primitive (i.e. they are - * non-zero in only one - * component), then the value - * returned by shape_grad() - * equals that of this function - * for exactly one - * component. This function is - * therefore only of greater - * interest if the shape function - * is not primitive, but then it - * is necessary since the other - * function cannot be used. - */ - Tensor<1,dim> - shape_grad_component (const unsigned int function_no, - const unsigned int point_no, - const unsigned int component) const; - /** * Compute the gradients of the finite * element function characterized @@ -712,71 +785,28 @@ class FEValuesBase : protected FEValuesData void get_function_grads (const InputVector &fe_function, std::vector > > &gradients) const; - /** - * Second derivatives of - * the @p function_noth shape function at - * the @p point_noth quadrature point - * with respect to real cell - * coordinates. If you want to - * get the derivatives in one of - * the coordinate directions, use - * the appropriate function of - * the @p Tensor class to - * extract one component. Since - * only a reference to the - * derivative values is returned, - * there should be no major - * performance drawback. - * - * If the shape function is - * vector-valued, then this - * returns the only non-zero - * component. If the shape - * function has more than one - * non-zero component (i.e. it is - * not primitive), then throw an - * exception of type - * @p ExcShapeFunctionNotPrimitive. In - * that case, use the - * shape_grad_grad_component() - * function. + /** + * Function gradient access with + * more flexibility. see + * get_function_values() with + * corresponding arguments. */ - const Tensor<2,dim> & - shape_2nd_derivative (const unsigned int function_no, - const unsigned int point_no) const; - + template + void get_function_grads (const InputVector& fe_function, + const std::vector& indices, + std::vector >& gradients) const; /** - * Return one vector component of - * the gradient of a shape - * function at a quadrature - * point. If the finite element - * is scalar, then only component - * zero is allowed and the return - * value equals that of the - * @p shape_2nd_derivative - * function. If the finite - * element is vector valued but - * all shape functions are - * primitive (i.e. they are - * non-zero in only one - * component), then the value - * returned by - * @p shape_2nd_derivative - * equals that of this function - * for exactly one - * component. This function is - * therefore only of greater - * interest if the shape function - * is not primitive, but then it - * is necessary since the other - * function cannot be used. + * Function gradient access with + * more flexibility. see + * get_function_values() with + * corresponding arguments. */ - Tensor<2,dim> - shape_2nd_derivative_component (const unsigned int function_no, - const unsigned int point_no, - const unsigned int component) const; - + template + void get_function_grads (const InputVector& fe_function, + const std::vector& indices, + std::vector > >& gradients) const; + /** * Compute the tensor of second * derivatives of the finite @@ -855,6 +885,7 @@ class FEValuesBase : protected FEValuesData void get_function_2nd_derivatives (const InputVector &fe_function, std::vector > > &second_derivatives) const; + //@} /** * Position of the @p ith diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 2c6d8f27c7..9dff4e3b08 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -598,6 +598,43 @@ get_function_grads (const InputVector &fe_function, } +template +template +void FEValuesBase::get_function_grads ( + const InputVector& fe_function, + const std::vector& indices, + std::vector > &values) const +{ + Assert (this->update_flags & update_gradients, ExcAccessToUninitializedField()); + // This function fills a single + // component only + Assert (fe->n_components() == 1, + ExcWrongNoOfComponents()); + // One index for each dof + Assert (indices.size() == dofs_per_cell, + ExcDimensionMismatch(indices.size(), dofs_per_cell)); + // This vector has one entry for + // each quadrature point + Assert (values.size() == n_quadrature_points, + ExcWrongVectorSize(values.size(), n_quadrature_points)); + + // initialize with zero + std::fill_n (values.begin(), n_quadrature_points, Tensor<1,dim>()); + + // add up contributions of trial + // functions. note that here we + // deal with scalar finite + // elements, so no need to check + // for non-primitivity of shape + // functions + for (unsigned int point=0; pointshape_grad(shape_func, point)); +} + + + template template @@ -652,6 +689,71 @@ get_function_grads (const InputVector &fe_function, +template +template +void FEValuesBase::get_function_grads ( + const InputVector& fe_function, + const std::vector& indices, + std::vector > >& values) const +{ + // One value per quadrature point + Assert (n_quadrature_points == values.size(), + ExcWrongVectorSize(values.size(), n_quadrature_points)); + + const unsigned int n_components = fe->n_components(); + + // Size of indices must be a + // multiple of dofs_per_cell such + // that an integer number of + // function values is generated in + // each point. + Assert (indices.size() % dofs_per_cell == 0, + ExcNotMultiple(indices.size(), dofs_per_cell)); + + // The number of components of the + // result may be a multiple of the + // number of components of the + // finite element + const unsigned int result_components = indices.size() / dofs_per_cell; + + for (unsigned i=0;iupdate_flags & update_values, ExcAccessToUninitializedField()); + + // initialize with zero + for (unsigned i=0;i()); + + // add up contributions of trial + // functions. now check whether the + // shape function is primitive or + // not. if it is, then set its only + // non-zero component, otherwise + // loop over components + for (unsigned int mc = 0; mc < component_multiple; ++mc) + for (unsigned int point=0; pointis_primitive(shape_func)) + values[point][fe->system_to_component_index(shape_func).first + +mc * n_components] + += fe_function(indices[shape_func+mc*dofs_per_cell]) + * shape_grad(shape_func, point); + else + for (unsigned int c=0; c template void diff --git a/deal.II/deal.II/source/fe/fe_values.instance.h b/deal.II/deal.II/source/fe/fe_values.instance.h index 89158219ee..d76e138945 100644 --- a/deal.II/deal.II/source/fe/fe_values.instance.h +++ b/deal.II/deal.II/source/fe/fe_values.instance.h @@ -48,10 +48,18 @@ void FEValuesBase::get_function_values template void FEValuesBase::get_function_grads (const IN&, std::vector > &) const; +template +void FEValuesBase::get_function_grads +(const IN&, const std::vector&, + std::vector > &) const; template void FEValuesBase::get_function_grads (const IN&, std::vector > > &) const; +template +void FEValuesBase::get_function_grads +(const IN&, const std::vector&, + std::vector > > &) const; template void FEValuesBase::get_function_2nd_derivatives -- 2.39.5