*/
static
unsigned int
- unrolled_index (const TableIndices<rank> &indices);
+ component_to_unrolled_index (const TableIndices<rank> &indices);
+
+ /**
+ * The opposite of the previous
+ * function: given an index $i$
+ * in the unrolled form of the
+ * tensor, return what set of
+ * indices $(k,l)$ (for rank-2
+ * tensors) or $(k,l,m,n)$ (for
+ * rank-4 tensors) corresponds to
+ * it.
+ */
+ static
+ TableIndices<rank>
+ unrolled_to_component_indices (const unsigned int i);
/**
* Reset all values to zero.
template <>
inline
unsigned int
-SymmetricTensor<2,1>::unrolled_index (const TableIndices<2> &indices)
+SymmetricTensor<2,1>::component_to_unrolled_index (const TableIndices<2> &indices)
{
const unsigned int dim = 1;
Assert (indices[0] < dim, ExcIndexRange(indices[0], 0, dim));
template <>
inline
unsigned int
-SymmetricTensor<2,2>::unrolled_index (const TableIndices<2> &indices)
+SymmetricTensor<2,2>::component_to_unrolled_index (const TableIndices<2> &indices)
{
const unsigned int dim = 2;
Assert (indices[0] < dim, ExcIndexRange(indices[0], 0, dim));
template <>
inline
unsigned int
-SymmetricTensor<2,3>::unrolled_index (const TableIndices<2> &indices)
+SymmetricTensor<2,3>::component_to_unrolled_index (const TableIndices<2> &indices)
{
const unsigned int dim = 3;
Assert (indices[0] < dim, ExcIndexRange(indices[0], 0, dim));
+template <>
+inline
+TableIndices<2>
+SymmetricTensor<2,1>::unrolled_to_component_indices (const unsigned int i)
+{
+ Assert (i < n_independent_components, ExcIndexRange(i, 0, n_independent_components));
+
+ return TableIndices<2>(0,0);
+}
+
+
+
+template <>
+inline
+TableIndices<2>
+SymmetricTensor<2,2>::unrolled_to_component_indices (const unsigned int i)
+{
+ Assert (i < n_independent_components, ExcIndexRange(i, 0, n_independent_components));
+
+ static const TableIndices<2> table[n_independent_components] =
+ { TableIndices<2> (0,0),
+ TableIndices<2> (1,1),
+ TableIndices<2> (0,1) };
+
+ return table[i];
+}
+
+
+
+template <>
+inline
+TableIndices<2>
+SymmetricTensor<2,3>::unrolled_to_component_indices (const unsigned int i)
+{
+ Assert (i < n_independent_components, ExcIndexRange(i, 0, n_independent_components));
+
+ static const TableIndices<2> table[n_independent_components] =
+ { TableIndices<2> (0,0),
+ TableIndices<2> (1,1),
+ TableIndices<2> (2,2),
+ TableIndices<2> (0,1),
+ TableIndices<2> (0,2),
+ TableIndices<2> (1,2) };
+
+ return table[i];
+}
+
+
+
#endif // DOXYGEN
--- /dev/null
+//---------------------------- symmetric_tensor_24.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2006, 2008, 2009 by the deal.II authors
+//
+// This file is subject to QPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//---------------------------- symmetric_tensor_24.cc ---------------------------
+
+// check SymmetricTensor<2,dim>::component_to_unrolled_index and the
+// other way round
+
+#include "../tests.h"
+#include <base/symmetric_tensor.h>
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <fstream>
+#include <iomanip>
+
+
+template <int dim>
+void check ()
+{
+ typedef SymmetricTensor<2,dim> S;
+ for (unsigned int i=0; i<S::n_independent_components; ++i)
+ {
+ deallog << i << " -- "
+ << S::unrolled_to_component_indices (i)
+ << std::endl;
+ Assert (S::component_to_unrolled_index
+ (S::unrolled_to_component_indices (i))
+ ==
+ i,
+ ExcInternalError());
+ }
+}
+
+
+int main ()
+{
+ std::ofstream logfile("symmetric_tensor_24/output");
+ deallog << std::setprecision(3);
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ check<1> ();
+ check<2> ();
+ check<3> ();
+
+ deallog << "OK" << std::endl;
+}