#include <base/tensor_base.h>
-
/**
* Provide a general tensor class with an arbitrary rank, i.e. with
* an arbitrary number of indices. The Tensor class provides an
*/
Tensor<rank_,dim> operator - (const Tensor<rank_,dim> &) const;
+ /**
+ * Fill a vector with all tensor elements.
+ *
+ * Thsi function unrolls all
+ * tensor entries into a single,
+ * linearly numbered vector. As
+ * usual in C++, the rightmost
+ * index marches fastest.
+ */
+ void unroll(vector<double> & result) const;
+
+
/**
* Reset all values to zero.
*/
void clear ();
+ /**
+ * Exception
+ */
+ DeclException2(ExcWrongVectorSize, unsigned, int, << "Tensor has " << arg1
+ << " entries, but vector has size " << arg2);
/**
* Exception
* subelements.
*/
Tensor<rank_-1,dim> subtensor[dim];
+
+ /**
+ * Help function for unroll.
+ */
+ void unroll_recursion(vector<double> & result, unsigned& start_index) const;
+
+ template<>
+ friend void Tensor<rank_+1,dim>::unroll_recursion(vector<double> &, unsigned&) const;
};
#include <base/exceptions.h>
#include <iostream>
+#include <vector>
* Reset all values to zero.
*/
void clear ();
-
+
+ /**
+ * Fill a vector with all tensor elements.
+ *
+ * Thsi function unrolls all
+ * tensor entries into a single,
+ * linearly numbered vector. As
+ * usual in C++, the rightmost
+ * index marches fastest.
+ */
+ void unroll(vector<double> & result) const;
+
protected:
/**
* Stores the values in a simple array.
*/
double values[dim];
+
+ /**
+ * Help function for unroll.
+ */
+ void unroll_recursion(vector<double> & result, unsigned& start_index) const;
+
+ template<>
+ friend void Tensor<2,dim>::unroll_recursion(vector<double> &, unsigned&) const;
};
template<int rank_, int dim>
class Tensor;
+/**
+ * Base class for multi-valued functions.
+ *
+ *
+ */
+template <int dim>
+class VectorFunction :
+ public FunctionTime
+{
+ public:
+ /**
+ * Constructor. May take an initial vakue
+ * for the time variable, which defaults
+ * to zero.
+ */
+ VectorFunction (unsigned n_components, const double initial_time = 0.0);
+
+ /**
+ * Virtual destructor; absolutely
+ * necessary in this case.
+ */
+ virtual ~VectorFunction ();
+
+// /**
+// * Return the value of the function
+// * at the given point.
+// */
+// virtual double operator () (const Point<dim> &p, unsigned component) const;
+
+ /**
+ * Set #values# to the point values
+ * of the function at points #p#.
+ * It is assumed that #values#
+ * already has the right size, i.e.
+ * the same size as the #n_components#
+ * array.
+ */
+ virtual void value (const Point<dim> &p, vector<double> &values) const;
+
+ /**
+ * Set #values# to the point values
+ * of the function at the #points#.
+ * It is assumed that #values#
+ * already has the right size, i.e.
+ * the same size as the #points#
+ * array.
+ */
+ virtual void value_list (const vector<Point<dim> > &points,
+ vector<vector<double> > &values) const;
+
+ /**
+ * Set #gradients# to the gradients of
+ * the function at the #points#.
+ * It is assumed that #values#
+ * already has the right size, i.e.
+ * the same size as the #points# array.
+ */
+ virtual void gradient_list (const vector<Point<dim> > &points,
+ vector<vector<Tensor<1,dim> > > &gradients) const;
+
+ /**
+ * Number of vector components.
+ */
+ const unsigned n_components;
+
+ /**
+ * Exception
+ */
+ DeclException0 (ExcPureFunctionCalled);
+ /**
+ * Exception
+ */
+ DeclException2 (ExcVectorHasWrongSize,
+ int, int,
+ << "The vector has size " << arg1 << " but should have "
+ << arg2 << " elements.");
+
+};
+
/**
* This class is a model for a tensor valued function.
* It returns the value
*/
template <int rank_, int dim>
class TensorFunction :
- public FunctionTime
+ public VectorFunction<dim>
{
public:
/**
*/
virtual void gradient_list (const vector<Point<dim> > &points,
vector<Tensor<rank_+1,dim> > &gradients) const;
-
+
+ /**
+ * See #VectorFunction#.
+ */
+ virtual void value_list (const vector<Point<dim> > &points,
+ vector<vector<double> > &values) const;
+
+ /**
+ * See #VectorFunction#.
+ */
+ virtual void gradient_list (const vector<Point<dim> > &points,
+ vector<vector<Tensor<1,dim> > > &gradients) const;
+
/**
* Exception
*/
--- /dev/null
+// $Id$
+
+#include <base/tensor.h>
+#include <cmath>
+
+template <int rank_, int dim> void
+Tensor<rank_, dim>::unroll( vector<double>& result) const
+{
+ Assert(result.size()==pow(dim,rank_),
+ ExcWrongVectorSize(pow(dim,rank_), result.size()));
+
+ unsigned index = 0;
+ unroll_recursion(result,index);
+}
+
+
+template <int rank_, int dim> void
+Tensor<rank_, dim>::unroll_recursion( vector<double>& result, unsigned& index) const
+{
+ for (unsigned i=0; i<dim; ++i)
+ {
+ operator[](i).unroll_recursion(result, index);
+ }
+
+}
+
+
+template<int dim> void
+Tensor<1,dim>::unroll_recursion( vector<double>& result, unsigned& index) const
+{
+ for (unsigned i=0; i<dim; ++i)
+ {
+ cerr << "[" << index << ',' << operator[](i) << ']' << endl;
+ result[index++] = operator[](i);
+ }
+
+}
+
+
+template class Tensor<1, 1>;
+template class Tensor<1, 2>;
+template class Tensor<1, 3>;
+template class Tensor<2, 1>;
+template class Tensor<2, 2>;
+template class Tensor<2, 3>;
+template class Tensor<3, 1>;
+template class Tensor<3, 2>;
+template class Tensor<3, 3>;
+template class Tensor<4, 1>;
+template class Tensor<4, 2>;
+template class Tensor<4, 3>;
#include <base/tensorfunction.h>
#include <vector>
#include <base/tensor.h>
+#include <cmath>
+template <int dim>
+VectorFunction<dim>::VectorFunction(unsigned n_components, const double initial_time)
+ :
+ FunctionTime(initial_time),
+ n_components(n_components)
+{}
+
+
+template <int dim>
+VectorFunction<dim>::~VectorFunction()
+{}
+
+/*
+template <int dim> double
+VectorFunction<dim>::operator () (const Point<dim> &, unsigned) const
+
+{
+ Assert (false, ExcPureFunctionCalled());
+ return 0.;
+}
+*/
+
+template <int dim> void
+VectorFunction<dim>::value (const Point<dim> &, vector<double> &) const
+{
+ Assert (false, ExcPureFunctionCalled());
+}
+
+
+template <int dim> void
+VectorFunction<dim>::value_list (const vector<Point<dim> > &,
+ vector<vector<double> > &) const
+{
+ Assert (false, ExcPureFunctionCalled());
+}
+
+
+template <int dim> void
+VectorFunction<dim>::gradient_list (const vector<Point<dim> > &,
+ vector<vector<Tensor<1,dim> > > &) const
+{
+ Assert (false, ExcPureFunctionCalled());
+}
+
+
+//////////////////////////////////////////////////////////////////////
+// TensorFunction
+//////////////////////////////////////////////////////////////////////
template <int rank_, int dim>
-TensorFunction<rank_, dim>::TensorFunction (const double initial_time) :
- FunctionTime(initial_time)
+TensorFunction<rank_, dim>::TensorFunction (const double initial_time)
+ :
+ VectorFunction<dim>(pow(dim,rank_), initial_time)
{};
gradients[i] = gradient(points[i]);
};
+/*
+template <int rank_, int dim> void
+TensorFunction<rank_, dim>::value (const Point<dim> &p, vector<double> &erg) const
+{
+ Tensor<rank_,dim> h = operator()(p);
+ h.unroll(erg);
+}
+*/
+
+template <int rank_, int dim> void
+TensorFunction<rank_, dim>::value_list (const vector<Point<dim> > & points,
+ vector<vector<double> > & values) const
+{
+ Assert (values.size() == points.size(),
+ ExcVectorHasWrongSize(values.size(), points.size()));
+
+ for (unsigned int i=0; i<points.size(); ++i)
+ operator() (points[i]).unroll(values[i]);
+
+}
+
+
+template <int rank_, int dim> void
+TensorFunction<rank_, dim>::gradient_list (const vector<Point<dim> > &,
+ vector<vector<Tensor<1,dim> > > &) const
+{
+ Assert (false, ExcPureFunctionCalled());
+}
+
+
+
template class TensorFunction<1,1>;
template class TensorFunction<2,1>;
template class TensorFunction<3,1>;