]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add missing methods FEValuesViews::SymmetricTensor< 2, dim, spacedim >::value and...
authormcbride <mcbride@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 9 Sep 2010 08:51:27 +0000 (08:51 +0000)
committermcbride <mcbride@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 9 Sep 2010 08:51:27 +0000 (08:51 +0000)
git-svn-id: https://svn.dealii.org/trunk@21901 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 20af1d29696748f1f22c4b6e70bc25edcd310575..8ec03b79e47a4b9d093f7e5b195418d0d9631935 100644 (file)
@@ -4110,6 +4110,127 @@ namespace FEValuesViews
        return symmetrize(return_value);
       }
   }
+  
+  
+  
+  template <int dim, int spacedim>
+          inline
+          typename SymmetricTensor<2, dim, spacedim>::value_type
+          SymmetricTensor<2, dim, spacedim>::value (const unsigned int shape_function,
+                                                    const unsigned int q_point) const
+  {
+      typedef FEValuesBase<dim,spacedim> FVB;
+      Assert (shape_function < fe_values.fe->dofs_per_cell,
+              ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
+      Assert (fe_values.update_flags & update_values,
+              typename FVB::ExcAccessToUninitializedField());
+
+                                // similar to the vector case where
+                                // we have more then one index and we need
+                                // to convert between unrolled and component
+                                // indexing for tensors
+
+      const int snc = shape_function_data[shape_function].single_nonzero_component;
+
+      if (snc == -2)
+      {
+                                    // shape function is zero for the
+                                    // selected components
+          return value_type();
+
+      } else if (snc != -1)
+      {
+          value_type return_value;
+          const unsigned int comp =
+                  shape_function_data[shape_function].single_nonzero_component_index;
+          return_value[value_type::unrolled_to_component_indices(comp)]
+                  = fe_values.shape_values(snc,q_point);
+          return return_value;
+      }
+      else
+      {
+          value_type return_value;
+          for (unsigned int d = 0; d < value_type::n_independent_components; ++d)
+              if (shape_function_data[shape_function].is_nonzero_shape_function_component[d])
+                  return_value[value_type::unrolled_to_component_indices(d)]
+                          = fe_values.shape_values(shape_function_data[shape_function].row_index[d],q_point);
+          return return_value;
+      }
+  }
+
+  template <int dim, int spacedim>
+          inline
+          typename SymmetricTensor<2, dim, spacedim>::divergence_type
+          SymmetricTensor<2, dim, spacedim>::divergence(const unsigned int shape_function,
+             const unsigned int q_point) const
+  {
+      typedef FEValuesBase<dim,spacedim> FVB;
+      Assert (shape_function < fe_values.fe->dofs_per_cell,
+              ExcIndexRange (shape_function, 0, fe_values.fe->dofs_per_cell));
+      Assert (fe_values.update_flags & update_gradients,
+              typename FVB::ExcAccessToUninitializedField());
+
+      const int snc = shape_function_data[shape_function].single_nonzero_component;
+
+      if (snc == -2)
+      {
+                            // shape function is zero for the
+                            // selected components
+          return divergence_type();
+      } else if (snc != -1) {
+                            // have a single non-zero component when the
+                            // symmetric tensor is repsresented in unrolled form.
+                            // this implies we potentially have two non-zero
+                            // components when represented in component form!
+                            // we will only have one non-zero entry if the non-zero
+                            // component lies on the diagonal of the tensor.
+                            //
+                            // the divergence of a second-order tensor
+                            // is a first order tensor.
+                            //
+                            // assume the second-order tensor is A with componets A_{ij}.
+                            // then A_{ij} = A_{ji} and there is only one (if diagonal)
+                            // or two non-zero entries in the tensorial representation.
+                            // define the divergence as:
+                            // b_i := \dfrac{\partial A_{ij}}{\partial x_j}.
+                            //
+                            // Now, knowing the row ii and collumn jj of the non-zero entry
+                            // we compute the divergence as
+                            // b_ii = \dfrac{\partial A_{ij}}{\partial x_jj}  (no sum)
+                            // and if ii =! jj (not on a diagonal)
+                            // b_jj = \dfrac{\partial A_{ij}}{\partial x_ii}  (no sum)
+
+          divergence_type return_value;
+
+                            // non-zero index in unrolled format
+          const unsigned int comp =
+            shape_function_data[shape_function].single_nonzero_component_index;
+
+          const unsigned int ii = value_type::unrolled_to_component_indices(comp)[0];
+          const unsigned int jj = value_type::unrolled_to_component_indices(comp)[1];
+
+                            // value of the non-zero tensor component
+          const double A_ij = fe_values.shape_values(snc,q_point);
+
+                            // the gradient of the non-zero shape function
+          const Tensor<1, spacedim> phi_grad = fe_values.shape_gradients[snc][q_point];
+
+          return_value[ii] = A_ij * phi_grad[jj];
+
+                            // if we are not on a diagonal
+          if (ii != jj)
+              return_value[jj] = A_ij * phi_grad[ii];
+
+          return return_value;
+
+      } else
+      {
+        Assert (false, ExcNotImplemented());
+        divergence_type return_value;
+        return return_value;
+      }
+  }
+
 }
 
 

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.