---

SymbolicUtilities.h

Go to the documentation of this file.

#ifndef SYMBOLICUTILITIES_H
#define SYMBOLICUTILITIES_H

#include "SundanceDefs.h"

#include "Expr.h"
#include "Derivative.h"


namespace Sundance
{

  using namespace TSF;

  template<int dim> Expr nabla();

  template<> inline Expr nabla<1>()
    {
      static Expr dx = new Derivative(0);
      return dx;
    }

  template<> inline Expr nabla<2>()
    {
      static Expr dx = new Derivative(0);
      static Expr dy = new Derivative(1);
      return List(dx, dy);
    }

  template<> inline Expr nabla<3>()
    {
      static Expr dx = new Derivative(0);
      static Expr dy = new Derivative(1);
      static Expr dz = new Derivative(2);
      return List(dx, dy, dz);
    }

  template<int dim> Expr gradient(const Expr& f);

  template<> inline Expr gradient<1>(const Expr& f)
    {
      return nabla<1>()*f;
    }

  template<> inline Expr gradient<2>(const Expr& f)
    {
      return nabla<2>()*f;
    }

  template<> inline Expr gradient<3>(const Expr& f)
    {
      return nabla<2>()*f;
    }

  template<int dim> Expr divergence(const Expr& f);

  template<> inline Expr divergence<1>(const Expr& f)
    {
      if (f.length() != 1)
        {
          TSFError::raise("1D divergence operator called on a vector "
                          "of dimension " + TSF::toString(f.length())
                          + ": " + f.toString());
        }
      return nabla<1>()*f;
    }

  template<> inline Expr divergence<2>(const Expr& f)
    {
      if (f.length() != 2)
        {
          TSFError::raise("2D divergence operator called on a vector "
                          "of dimension " + TSF::toString(f.length())
                          + ": " + f.toString());
        }
      return nabla<2>() * f;
    }

  template<>  inline Expr divergence<3>(const Expr& f)
    {
      if (f.length() != 3)
        {
          TSFError::raise("3D divergence operator called on a vector "
                          "of dimension " + TSF::toString(f.length())
                          + ": " + f.toString());
        }
      return nabla<3>() * f;
    }

  inline Expr cross(const Expr& a, const Expr& b)
    {
      if (a.length() != b.length())
        {
          TSFError::raise("Unequal vector lengths in cross product: "
                          "first arg= " + a.toString()
                          + " second arg= " + b.toString());
        }
      if (a.length()==2)
        {
          return a[0]*b[1] - a[1]*b[0];
        }
      else if (a.length() == 3)
        {
          return List(a[1]*b[2] - a[2]*b[1],
                      -a[0]*b[2] + a[2]*b[0],
                      a[0]*b[1] - a[1]*b[0]);
        }
      else
        {
          TSFError::raise("vector lengths in cross product must be 2 or 3: "
                          "first arg= " + a.toString()
                          + " second arg= " + b.toString());
          return 0.0; // -Wall
        }
    }

  template<int dim> Expr curl(const Expr& f);

  template<> inline Expr curl<2>(const Expr& f)
    {
      return cross(nabla<2>(), f);
    }

  template<> inline Expr curl<3>(const Expr& f)
    {
      return cross(nabla<3>(), f);
    }

  template<int dim> Expr grad(const Expr& f);

  template<> inline Expr grad<1>(const Expr& f) {return gradient<1>(f);}
  template<> inline Expr grad<2>(const Expr& f) {return gradient<2>(f);}
  template<> inline Expr grad<3>(const Expr& f) {return gradient<3>(f);}

  template<int dim> Expr div(const Expr& f);

  template<> inline Expr div<1>(const Expr& f) {return divergence<1>(f);}
  template<> inline Expr div<2>(const Expr& f) {return divergence<2>(f);}
  template<> inline Expr div<3>(const Expr& f) {return divergence<3>(f);}
}


#endif

---

Contact:

Kevin Long (krlong@ca.sandia.gov)


Documentation generated by