Re: [pooma-dev] Revised Q: How to apply a Stencil to all leafs in an expression

[Richard Guenther]

As usually, I found a solution myself, but of course am going to share
it with you:

On Wed, 30 Oct 2002, Richard Guenther wrote:

> I.e. instead of writing
>
>   A(I) = (stencil(B)*stencil(C))(I);
>
> I'd like to write something like
>
>   A(I) = applyStencil(B*C)(I);
>
> for brevity (I really have huge expressions, they get completely
> unreadable otherwise).

Code along the following is needed to achieve this (where stencil
is called Dummy - still need to figure out how to template that thing
on a specific stencil, but at the moment I need this for polynomial
evaluation only):


struct ApplyDummy {};

template <class Mesh, class T, class EngineTag>
struct LeafFunctor< Field<Mesh, T, EngineTag>, ApplyDummy >
{
  typedef typename FieldStencilSimple<DummyStencil<Mesh::dimensions, T>,
    Field<Mesh, T, EngineTag> >::Type_t Type_t;
  inline static Type_t apply(const Field<Mesh, T, EngineTag> &f,
                             const ApplyDummy &)
  {
    return FieldStencilSimple<DummyStencil<Mesh::dimensions, T>,
    Field<Mesh, T, EngineTag> >::make(DummyStencil<Mesh::dimensions, T>(),
f);
  }
};

template <class Leaf>
struct LeafFunctor< Leaf, ApplyDummy >
{
  typedef Leaf Type_t;
  inline static Type_t apply(const Leaf &f,
                             const ApplyDummy &)
  {
    return f;
  }
};

template <class Mesh, class T, class Expr>
Field<Mesh, T, ExpressionTag< typename ForEach<Expr, ApplyDummy,
TreeCombine>::Type_t> >
applyDummy(const Field<Mesh, T, ExpressionTag<Expr> > &e)
{
  return Field<Mesh, T, ExpressionTag< typename ForEach<Expr, ApplyDummy,
TreeCombine>::Type_t> >
    (forEach(e.engine().expression(), ApplyDummy(), TreeCombine()));
}


So the idea is to just change the (Field) leaf engines to stencil
apply engines and stick them at the old place, for other leafs just
keep them, so f.i. applyStencil(A*B+0.5) still works.

Maybe this is of use to someone -- I use this to evaluate polynomial
expressions at certain points, f.e.

 peeval<0>((rhe_r + (gamma-1)*(rhe_r - 0.5*pow2(rhv_r)))
           * rhv_r/pow2(rh_r), 0.5)(I+dI)

where you can imagine the mess created by individual peval() calls at
every term of the expression.

Richard.

--
Richard Guenther <address@hidden>
WWW: http://www.tat.physik.uni-tuebingen.de/~rguenth/
The GLAME Project: http://www.glame.de/