From: guido Date: Thu, 12 May 2005 07:13:50 +0000 (+0000) Subject: new fe_values::get_function* for nonprimitive elements X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=51d8918a6a71770fc7f3851954ff7d6900506c5c;p=dealii-svn.git new fe_values::get_function* for nonprimitive elements git-svn-id: https://svn.dealii.org/trunk@10665 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index c9fbc2505b..d40864b352 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -706,6 +706,69 @@ class FEValuesBase : protected FEValuesData const VectorSlice >& indices, std::vector >& values) const; + + /** + * Generate vector function + * values from an arbitrary + * vector. + * + * This function offers the + * possibility to extract + * function values in quadrature + * points from vectors not + * corresponding to a whole + * discretization. + * + * The length of the vector + * indices may even be a + * multiple of the number of dofs + * per cell. Then, the vectors in + * value should allow + * for the same multiple of the + * components of the finite + * element. + * + * Depending on the last + * argument, the outer vector of + * values has either the + * length of the quadrature rule + * (quadrature_points_fastest + * == false) or the length + * of components to be filled + * quadrature_points_fastest + * == true. If p is + * the xurrent quadrature point + * number and i is the + * vector component of the + * solution desired, the access + * to values is + * values[p][i] if + * quadrature_points_fastest + * == false, and + * values[i][p] + * otherwise. + * + * You may want to use this + * function, if you want to + * access just a single block + * from a BlockVector, if you + * have a multi-level vector or + * if you already have a local + * representation of your finite + * element data. + * + * Since this function allows for + * fairly general combinations of + * argument sizes, be aware that + * the checks on the arguments + * may not detect errors. + */ + template + void get_function_values (const InputVector& fe_function, + const VectorSlice >& indices, + std::vector >& values, + bool quadrature_points_fastest) const; + /** * Compute the gradients of the finite * element function characterized @@ -806,7 +869,8 @@ class FEValuesBase : protected FEValuesData template void get_function_grads (const InputVector& fe_function, const VectorSlice >& indices, - std::vector > >& gradients) const; + std::vector > >& gradients, + bool quadrature_points_fastest = false) const; /** * Compute the tensor of second @@ -885,7 +949,8 @@ class FEValuesBase : protected FEValuesData template void get_function_2nd_derivatives (const InputVector &fe_function, - std::vector > > &second_derivatives) const; + std::vector > > &second_derivatives, + bool quadrature_points_fastest = false) const; //@} /** diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 3e9f1a6c98..477e360bc3 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -517,7 +517,7 @@ void FEValuesBase::get_function_values ( // result may be a multiple of the // number of components of the // finite element - const unsigned int result_components = indices.size() / dofs_per_cell; + const unsigned int result_components = indices.size() * n_components / dofs_per_cell; for (unsigned i=0;i::get_function_values ( * shape_value(shape_func, point)); else for (unsigned int c=0; c +template +void FEValuesBase::get_function_values ( + const InputVector& fe_function, + const VectorSlice >& indices, + std::vector >& values, + bool quadrature_points_fastest) const +{ + 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() * n_components / dofs_per_cell; + + // Check if the value argument is + // initialized to the correct sizes + if (quadrature_points_fastest) + { + Assert (values.size() == result_components, + ExcDimensionMismatch(values.size(), result_components)); + for (unsigned i=0;iupdate_flags & update_values, ExcAccessToUninitializedField()); + + // initialize with zero + for (unsigned i=0;iis_primitive(shape_func)) + values[fe->system_to_component_index(shape_func).first + +mc * n_components][point] + += (fe_function(indices[shape_func+mc*dofs_per_cell]) + * shape_value(shape_func, point)); + else + for (unsigned int c=0; cis_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_value(shape_func, point)); + else + for (unsigned int c=0; c tmp = this->shape_grad(shape_func,point); tmp *= dof_values(shape_func); - gradients[point][fe->system_to_component_index(shape_func).first] + gradients[fe->system_to_component_index(shape_func).first][point] += tmp; } else @@ -683,8 +777,8 @@ get_function_grads (const InputVector &fe_function, { Tensor<1,dim> tmp = this->shape_grad_component(shape_func,point,c); tmp *= dof_values(shape_func); - gradients[point][c] += tmp; - }; + gradients[c][point] += tmp; + } } @@ -694,12 +788,9 @@ template void FEValuesBase::get_function_grads ( const InputVector& fe_function, const VectorSlice >& indices, - std::vector > >& values) const + std::vector > >& values, + bool quadrature_points_fastest) const { - // One value per quadrature point - Assert (n_quadrature_points == values.size(), - ExcDimensionMismatch(values.size(), n_quadrature_points)); - const unsigned int n_components = fe->n_components(); // Size of indices must be a @@ -714,11 +805,26 @@ void FEValuesBase::get_function_grads ( // result may be a multiple of the // number of components of the // finite element - const unsigned int result_components = indices.size() / dofs_per_cell; + const unsigned int result_components = indices.size() * n_components / dofs_per_cell; - for (unsigned i=0;i::get_function_grads ( // 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; cis_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; cis_primitive(shape_func)) + values[fe->system_to_component_index(shape_func).first + +mc * n_components][point] + += fe_function(indices[shape_func+mc*dofs_per_cell]) + * shape_grad(shape_func, point); + else + for (unsigned int c=0; c void FEValuesBase:: get_function_2nd_derivatives (const InputVector &fe_function, - std::vector > > &second_derivs) const + std::vector > > &second_derivs, + bool quadrature_points_fastest) const { Assert (n_quadrature_points == second_derivs.size(), ExcDimensionMismatch(second_derivs.size(), n_quadrature_points)); @@ -824,22 +945,40 @@ get_function_2nd_derivatives (const InputVector &fe_func // add up contributions of trial // functions - for (unsigned int point=0; pointis_primitive(shape_func)) - { - Tensor<2,dim> tmp(shape_2nd_derivative(shape_func,point)); - tmp *= dof_values(shape_func); - second_derivs[point][fe->system_to_component_index(shape_func).first] - += tmp; - } - else - for (unsigned int c=0; cis_primitive(shape_func)) + { + Tensor<2,dim> tmp(shape_2nd_derivative(shape_func,point)); + tmp *= dof_values(shape_func); + second_derivs[point][fe->system_to_component_index(shape_func).first] + += tmp; + } + else + for (unsigned int c=0; c tmp = this->shape_2nd_derivative_component(shape_func,point,c); + tmp *= dof_values(shape_func); + second_derivs[point][c] += tmp; + } + else + for (unsigned int point=0; pointis_primitive(shape_func)) { - Tensor<2,dim> tmp = this->shape_2nd_derivative_component(shape_func,point,c); + Tensor<2,dim> tmp(shape_2nd_derivative(shape_func,point)); tmp *= dof_values(shape_func); - second_derivs[point][c] += tmp; - }; + second_derivs[fe->system_to_component_index(shape_func).first][point] + += tmp; + } + else + for (unsigned int c=0; c tmp = this->shape_2nd_derivative_component(shape_func,point,c); + tmp *= dof_values(shape_func); + second_derivs[c][point] += tmp; + } } 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 ac9e17c272..f02bebe622 100644 --- a/deal.II/deal.II/source/fe/fe_values.instance.h +++ b/deal.II/deal.II/source/fe/fe_values.instance.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -45,6 +45,15 @@ void FEValuesBase::get_function_values (const IN&, const VectorSlice >&, std::vector > &) const; +template +void FEValuesBase::get_function_values +(const IN&, const VectorSlice >&, + std::vector > &, bool) const; +template +void FEValuesBase::get_function_values +(const IN&, const VectorSlice >&, + std::vector > &, bool) const; + template void FEValuesBase::get_function_grads (const IN&, std::vector > &) const; @@ -59,7 +68,7 @@ void FEValuesBase::get_function_grads template void FEValuesBase::get_function_grads (const IN&, const VectorSlice >&, - std::vector > > &) const; + std::vector > > &, bool) const; template void FEValuesBase::get_function_2nd_derivatives @@ -67,4 +76,4 @@ void FEValuesBase::get_function_2nd_derivatives template void FEValuesBase::get_function_2nd_derivatives -(const IN&, std::vector > > &) const; +(const IN&, std::vector > > &, bool) const;