virtual Tensor<1,dim> gradient (const Point<dim> &p,
const unsigned int component = 0) const;
};
+
+
+
+/**
+ * A class that represents a function object for a monomial. Monomials are
+ * polynomials with only a single term, i.e. in 1-d they have the form
+ * $x^\alpha$, in 2-d the form $x_1^{\alpha_1}x_2^{\alpha_2}$, and in 3-d
+ * $x_1^{\alpha_1}x_2^{\alpha_2}x_3^{\alpha_3}$. Monomials are therefore
+ * described by a $dim$-tuple of exponents. Consequently, the class's
+ * constructor takes a Tensor<1,dim> to describe the set of exponents. Most of
+ * the time these exponents will of course be integers, but real exponents are
+ * of course equally valid.
+ *
+ * @author Wolfgang Bangerth, 2006
+ */
+ template <int dim>
+ class Monomial : public Function<dim>
+ {
+ public:
+ /**
+ * Constructor. The arguments are
+ * described in the general description
+ * of the class.
+ */
+ Monomial (const Tensor<1,dim> &exponents);
+
+ /**
+ * Function value at one point.
+ */
+ virtual double value (const Point<dim> &p,
+ const unsigned int component = 0) const;
+
+ /**
+ * Function values at multiple points.
+ */
+ virtual void value_list (const std::vector<Point<dim> > &points,
+ std::vector<double> &values,
+ const unsigned int component = 0) const;
+
+ /**
+ * Function gradient at one point.
+ */
+ virtual Tensor<1,dim> gradient (const Point<dim> &p,
+ const unsigned int component = 0) const;
+
+ private:
+ /**
+ * The set of exponents.
+ */
+ const Tensor<1,dim> exponents;
+ };
+
}
DEAL_II_NAMESPACE_CLOSE
return sum;
}
+
+
+
+
+/* ---------------------- Monomial ----------------------- */
+
+
+
+ template <int dim>
+ Monomial<dim>::
+ Monomial (const Tensor<1,dim> &exponents)
+ :
+ exponents (exponents)
+ {}
+
+
+
+ template <int dim>
+ double
+ Monomial<dim>::value (const Point<dim> &p,
+ const unsigned int component) const
+ {
+ Assert (component==0, ExcIndexRange(component,0,1)) ;
+
+ double prod = 1;
+ for (unsigned int s=0; s<dim; ++s)
+ prod *= std::pow(p[s], exponents[s]);
+
+ return prod;
+ }
+
+
+
+ template <int dim>
+ Tensor<1,dim>
+ Monomial<dim>::gradient (const Point<dim> &p,
+ const unsigned int component) const
+ {
+ Assert (component==0, ExcIndexRange(component,0,1)) ;
+
+ Tensor<1,dim> r;
+ for (unsigned int d=0; d<dim; ++d)
+ {
+ double prod = 1;
+ for (unsigned int s=0; s<dim; ++s)
+ prod *= (s==d
+ ?
+ exponents[s] * std::pow(p[s], exponents[s]-1)
+ :
+ std::pow(p[s], exponents[s]));
+
+ r[d] = prod;
+ }
+
+ return r;
+ }
+
+
+
+ template<int dim>
+ void
+ Monomial<dim>::value_list (const std::vector<Point<dim> > &points,
+ std::vector<double> &values,
+ const unsigned int component) const
+ {
+ Assert (values.size() == points.size(),
+ ExcDimensionMismatch(values.size(), points.size()));
+
+ for (unsigned int i=0; i<points.size(); ++i)
+ values[i] = Monomial<dim>::value (points[i], component);
+ }
+
// explicit instantiations
template class FourierSineSum<3>;
template class SlitSingularityFunction<2>;
template class SlitSingularityFunction<3>;
+ template class Monomial<1>;
+ template class Monomial<2>;
+ template class Monomial<3>;
}
DEAL_II_NAMESPACE_CLOSE
<h3>base</h3>
<ol>
+ <li> <p> New: The new <code>Functions::Monomial</code> class implements
+ monomials as a function object.
+ <br>
+ (WB 2006/12/15)
+ </p>
+
<li> <p> Fixed: If no substring is found with a width smaller than the given
threshold, the <code>Utilities::break_text_into_lines</code> function now
returns the smallest possible substring (larger than the threshold).
- <br> (Tobias Leicht 2006/12/13)
+ <br>
+ (Tobias Leicht 2006/12/13)
</p>
<li> <p> Changed: We used to set preprocessor variables