]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new get_function_grads
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Aug 2004 21:34:29 +0000 (21:34 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Aug 2004 21:34:29 +0000 (21:34 +0000)
git-svn-id: https://svn.dealii.org/trunk@9544 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_values.h
deal.II/deal.II/source/fe/fe_values.cc
deal.II/deal.II/source/fe/fe_values.instance.h

index f93a216e03251f7fa45856b913409dc3639db771..bdb6c40c671aed8ae013b886a47bdeb00811eefc 100644 (file)
@@ -121,7 +121,7 @@ class FEValuesData
                                      * derivatives.
                                      */
     typedef std::vector<std::vector<Tensor<2,dim> > > GradGradVector;
-    
+
                                     /**
                                      * Store the values of the shape
                                      * functions at the quadrature
@@ -366,6 +366,9 @@ class FEValuesBase : protected FEValuesData<dim>
                                      * 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<dim>
                                  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<dim>
                              const std::vector<unsigned int>& indices,
                              std::vector<Vector<number> >& 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<dim>
     void get_function_grads (const InputVector               &fe_function,
                             std::vector<std::vector<Tensor<1,dim> > > &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 <class InputVector>
+    void get_function_grads (const InputVector& fe_function,
+                            const std::vector<unsigned int>& indices,
+                            std::vector<Tensor<1,dim> >& 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 <class InputVector>
+    void get_function_grads (const InputVector& fe_function,
+                            const std::vector<unsigned int>& indices,
+                            std::vector<std::vector<Tensor<1,dim> > >& gradients) const;
+
                                     /**
                                      * Compute the tensor of second
                                      * derivatives of the finite
@@ -855,6 +885,7 @@ class FEValuesBase : protected FEValuesData<dim>
     void
     get_function_2nd_derivatives (const InputVector      &fe_function,
                                   std::vector<std::vector<Tensor<2,dim> > > &second_derivatives) const;
+                                    //@}
     
                                     /**
                                      * Position of the @p ith
index 2c6d8f27c7f82a011fa232ea4831006d2b9b1669..9dff4e3b08580b90f6a43057d491af21151ead5a 100644 (file)
@@ -598,6 +598,43 @@ get_function_grads (const InputVector           &fe_function,
 }
 
 
+template <int dim>
+template <class InputVector>
+void FEValuesBase<dim>::get_function_grads (
+  const InputVector& fe_function,
+  const std::vector<unsigned int>& indices,
+  std::vector<Tensor<1,dim> > &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; point<n_quadrature_points; ++point)
+    for (unsigned int shape_func=0; shape_func<dofs_per_cell; ++shape_func)
+      values[point] += (fe_function(indices[shape_func]) *
+                       this->shape_grad(shape_func, point));
+}
+
+
+
 
 template <int dim>
 template <class InputVector>
@@ -652,6 +689,71 @@ get_function_grads (const InputVector                         &fe_function,
 
 
 
+template <int dim>
+template <class InputVector>
+void FEValuesBase<dim>::get_function_grads (
+  const InputVector& fe_function,
+  const std::vector<unsigned int>& indices,
+  std::vector<std::vector<Tensor<1,dim> > >& 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;i<values.size();++i)
+    Assert (values[i].size() == result_components,
+           ExcDimensionMismatch(values[i].size(), result_components));
+
+                                  // If the result has more
+                                  // components than the finite
+                                  // element, we need this number for
+                                  // loops filling all components
+  const unsigned int component_multiple = result_components / n_components;
+  
+  Assert (this->update_flags & update_values, ExcAccessToUninitializedField());
+    
+                                  // initialize with zero
+  for (unsigned i=0;i<values.size();++i)
+    std::fill_n (values[i].begin(), values[i].size(), Tensor<1,dim>());
+
+                                  // 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; point<n_quadrature_points; ++point)
+      for (unsigned int shape_func=0; shape_func<dofs_per_cell; ++shape_func)
+       if (fe->is_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<n_components; ++c)
+           values[point][c] += (fe_function(indices[shape_func]) *
+                                shape_grad_component(shape_func, point, c));
+}
+
+
+
 template <int dim>
 template <class InputVector>
 void
index 89158219ee1e69d2ac0c21901822b27e663bef31..d76e13894585cd90f69027085c066bf9de89553c 100644 (file)
@@ -48,10 +48,18 @@ void FEValuesBase<deal_II_dimension>::get_function_values<IN>
 template
 void FEValuesBase<deal_II_dimension>::get_function_grads<IN>
 (const IN&, std::vector<Tensor<1,deal_II_dimension> > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads<IN>
+(const IN&, const std::vector<unsigned int>&,
+ std::vector<Tensor<1,deal_II_dimension> > &) const;
 
 template
 void FEValuesBase<deal_II_dimension>::get_function_grads<IN>
 (const IN&, std::vector<std::vector<Tensor<1,deal_II_dimension> > > &) const;
+template
+void FEValuesBase<deal_II_dimension>::get_function_grads<IN>
+(const IN&, const std::vector<unsigned int>&,
+ std::vector<std::vector<Tensor<1,deal_II_dimension> > > &) const;
 
 template
 void FEValuesBase<deal_II_dimension>::get_function_2nd_derivatives<IN>

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.