]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Change the data type of one argument. Add get_function_grads for vector valued finite...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 Mar 1999 18:46:56 +0000 (18:46 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 Mar 1999 18:46:56 +0000 (18:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@1059 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 301451e1d04c25a81d72f47966a91a3b0a037722..72a0f6393b008bc92d59637a4fa1d9b4846d4dbf 100644 (file)
@@ -296,16 +296,16 @@ class FEValuesBase {
     void get_function_values (const Vector<double> &fe_function,
                              vector<double>       &values) const;
 
-                                  /**
-                                   * Access to vector valued finite element functions.
-                                   *
-                                   * This function does the same as
-                                   * the other #get_function_values#,
-                                   * but applied to multi-component
-                                   * elements.
-                                   */
+                                    /**
+                                     * Access to vector valued finite element functions.
+                                     *
+                                     * This function does the same as
+                                     * the other #get_function_values#,
+                                     * but applied to multi-component
+                                     * elements.
+                                     */
     void get_function_values (const Vector<double>    &fe_function,
-                             vector<Vector<double> > &values) const;
+                             vector<vector<double> > &values) const;
 
                                     /**
                                      * Return the gradient of the #i#th shape
@@ -345,6 +345,24 @@ class FEValuesBase {
     void get_function_grads (const Vector<double>   &fe_function,
                             vector<Tensor<1,dim> > &gradients) const;
 
+                                    /**
+                                     * Return the gradients of the finite
+                                     * element function characterized
+                                     * by #fe_function# restricted to
+                                     * #cell# at the quadrature points.
+                                     *
+                                     * The function assumes that the
+                                     * #gradients# object already has the
+                                     * right size.
+                                     *
+                                     * This function does the same as
+                                     * the other #get_function_values#,
+                                     * but applied to multi-component
+                                     * elements.
+                                     */
+    void get_function_grads (const Vector<double>   &fe_function,
+                            vector<vector<Tensor<1,dim> > > &gradients) const;
+
                                     /**
                                      * Return the 2nd derivatives of
                                      * the #i#th shape function at
index 158dd8043d732f646bfe8f0898e128619c49131c..41f4dec03919ef31e9b62213756fc1abaa24a200 100644 (file)
@@ -90,7 +90,7 @@ void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
 
 template <int dim>
 void FEValuesBase<dim>::get_function_values (const Vector<double>     &fe_function,
-                                            vector< Vector<double> > &values) const
+                                            vector< vector<double> > &values) const
 {
   Assert (n_quadrature_points == values.size(),
          ExcWrongNoOfComponents());
@@ -113,7 +113,7 @@ void FEValuesBase<dim>::get_function_values (const Vector<double>     &fe_functi
                                   // functions
   for (unsigned int point=0; point<n_quadrature_points; ++point)
     for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
-      values[point](fe->system_to_component_index(shape_func).first)
+      values[point][fe->system_to_component_index(shape_func).first]
        += (dof_values(shape_func) * shape_values[selected_dataset](shape_func, point));
 };
 
@@ -160,6 +160,38 @@ void FEValuesBase<dim>::get_function_grads (const Vector<double>   &fe_function,
       };
 };
 
+template <int dim>
+void FEValuesBase<dim>::get_function_grads (const Vector<double>            &fe_function,
+                                           vector<vector<Tensor<1,dim> > > &gradients) const
+{
+  Assert (n_quadrature_points == gradients.size(),
+         ExcWrongNoOfComponents());
+  Assert (selected_dataset<shape_values.size(),
+         ExcInvalidIndex (selected_dataset, shape_values.size()));
+  for (unsigned i=0;i<gradients.size();++i)
+    Assert (gradients[i].size() == fe->n_components,
+           ExcWrongVectorSize(gradients.size(), n_quadrature_points));
+
+                                  // get function values of dofs
+                                  // on this cell
+  Vector<double> dof_values (total_dofs);
+  present_cell->get_dof_values (fe_function, dof_values);
+
+                                  // initialize with zero
+  for (unsigned i=0;i<gradients.size();++i)
+    fill_n (gradients[i].begin(), gradients[i].size(), Tensor<1,dim>());
+
+                                  // add up contributions of trial
+                                  // functions
+  for (unsigned int point=0; point<n_quadrature_points; ++point)
+    for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
+      {
+       Tensor<1,dim> tmp(shape_gradients[shape_func][point]);
+       tmp *= dof_values(shape_func);
+       gradients[point][fe->system_to_component_index(shape_func).first]
+         += tmp;
+      };
+};
 
 
 template <int dim>

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.