PDEOperatorExpression.hpp

Go to the documentation of this file.

//  This file is part of ff3d - http://www.freefem.org/ff3d
//  Copyright (C) 2001, 2002, 2003 Stéphane Del Pino

//  This program is free software; you can redistribute it and/or modify
//  it under the terms of the GNU General Public License as published by
//  the Free Software Foundation; either version 2, or (at your option)
//  any later version.

//  This program is distributed in the hope that it will be useful,
//  but WITHOUT ANY WARRANTY; without even the implied warranty of
//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU General Public License for more details.

//  You should have received a copy of the GNU General Public License
//  along with this program; if not, write to the Free Software Foundation,
//  Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.  

//  $Id: PDEOperatorExpression.hpp,v 1.10 2007/05/20 23:02:48 delpinux Exp $

#ifndef PDE_OPERATOR_EXPRESSION_HPP
#define PDE_OPERATOR_EXPRESSION_HPP

#include <TinyVector.hpp>

#include <Expression.hpp>

#include <Variable.hpp>
#include <FunctionExpression.hpp>
#include <FunctionExpressionBinaryOperation.hpp>

#include <EmbededFunctions.hpp>

#include <ErrorHandler.hpp>

#include <limits>

class PDEOperator;

class PDEOperatorExpression
  : public Expression
{
public:
  enum PDEOperatorType {
    scalar,
    vectorial
  };

protected:
  PDEOperatorType __pdeOperatorType;

  virtual void __execute() = 0;

public:
  virtual bool __checkBoundaryExpression() const = 0;

  void execute();

  const PDEOperatorType& pdeOperatorType() const
  {
    return __pdeOperatorType;
  }

  PDEOperatorExpression(const PDEOperatorExpression& e);

  PDEOperatorExpression(PDEOperatorType t);

  virtual ~PDEOperatorExpression();
};

class PDEScalarOperatorExpression
  : public PDEOperatorExpression
{
public:
  enum PDEScalarOperatorType {
    divMuGrad,
    secondOrder,
    firstOrder,
    zerothOrder
  };

protected:
  const std::string __unknownName;

  ReferenceCounting<PDEOperator> __pdeOperator;

private:
  PDEScalarOperatorType __pdeScalarOperatorType;

public:
  const std::string& unknownName() const
  {
    return __unknownName;
  }

  ReferenceCounting<PDEOperator> pdeOperator();

  PDEScalarOperatorExpression(const PDEScalarOperatorExpression& e);

  PDEScalarOperatorExpression(const std::string& unknownName,
                              PDEScalarOperatorType t);

  virtual ~PDEScalarOperatorExpression();
};

class PDEVectorialOperatorExpression
  : public PDEOperatorExpression
{
public:
  enum PDEScalarOperatorType {
    firstOrder,
    secondOrder
  };

public:
  typedef std::map<std::string,
                   ReferenceCounting<PDEOperator> > PDEOperatorList;

protected:
  PDEOperatorList __pdeOperatorList;

private:
  PDEScalarOperatorType __pdeOperatorType;

public:

  PDEVectorialOperatorExpression::PDEOperatorList& pdeOperatorList()
  {
    return __pdeOperatorList;
  }

  PDEVectorialOperatorExpression(const PDEVectorialOperatorExpression::
                                 PDEScalarOperatorType& t);

  PDEVectorialOperatorExpression(const PDEVectorialOperatorExpression& e);

  virtual ~PDEVectorialOperatorExpression();
};

class PDEScalarOperatorExpressionDivMuGrad
  : public PDEScalarOperatorExpression
{
private:
  ReferenceCounting<FunctionExpression> __mu;

  std::ostream& put(std::ostream& os) const
  {
    os << "div(" << (*__mu) << "*grad(" << __unknownName << "))";
    return os;
  }

  void __execute();
public:

  bool __checkBoundaryExpression() const;

  PDEScalarOperatorExpressionDivMuGrad(const std::string& unknownName,
                                       ReferenceCounting<FunctionExpression> mu)
    : PDEScalarOperatorExpression(unknownName, PDEScalarOperatorExpression::divMuGrad),
      __mu(mu)
  {
    ;
  }

  PDEScalarOperatorExpressionDivMuGrad(const PDEScalarOperatorExpressionDivMuGrad& re)
    : PDEScalarOperatorExpression(re),
      __mu(re.__mu)
  {
    ;
  }

  ~PDEScalarOperatorExpressionDivMuGrad()
  {
    ;
  }
};


class PDEScalarOperatorExpressionOrderZero
  : public PDEScalarOperatorExpression
{
private:
  ReferenceCounting<FunctionExpression> __Alpha;

  std::ostream& put(std::ostream& os) const
  {
    os << (*__Alpha) << '*' << __unknownName;
    return os;
  }

  void __execute();

public:
  bool __checkBoundaryExpression() const;

  PDEScalarOperatorExpressionOrderZero(const std::string& unknownName,
                                       ReferenceCounting<FunctionExpression> alpha)
    : PDEScalarOperatorExpression(unknownName, PDEScalarOperatorExpression::zerothOrder),
      __Alpha(alpha)
  {
    ;
  }

  PDEScalarOperatorExpressionOrderZero(const PDEScalarOperatorExpressionOrderZero& re)
    : PDEScalarOperatorExpression(re),
      __Alpha(re.__Alpha)
  {
    ;
  }

  ~PDEScalarOperatorExpressionOrderZero()
  {
    ;
  }
};


class PDEScalarOperatorExpressionOrderOne
  : public PDEScalarOperatorExpression
{
protected:
  ReferenceCounting<FunctionExpression> __nu[3];

  std::ostream& put(std::ostream& os) const
  {
    const char d[3] = { 'x', 'y', 'z'};
    bool plus = false;
    for (size_t i=0; i<3; ++i)
      if (__nu[i] != 0) {
        if (plus)
          os << '+';
        os << *(__nu[i]) << "*d" << d[i] << '(' << __unknownName << ')';
        plus = true;
      }
    return os;
  }

  void __execute();

public:
  bool __checkBoundaryExpression() const;

  void operator*=(ReferenceCounting<FunctionExpression> f)
  {
    for (size_t i=0; i<3; ++i) {
      if (__nu[i] != 0)
        __nu[i] = new FunctionExpressionBinaryOperation(BinaryOperation::product,f,__nu[i]);
    }
  }

  virtual void preexec()
  {
    for(size_t i=0; i<3; ++i) {
      if (__nu[i] != 0)
        (*__nu[i]).execute();
    }
  }

  ReferenceCounting<FunctionExpression> nu(const size_t i) const
  {
    ASSERT (i<3);
    return __nu[i];
  }

  const std::string& unknownName() const
  {
    return __unknownName;
  }

  PDEScalarOperatorExpressionOrderOne(const std::string& operatorName,
                                      const std::string&unknownName,
                                      ReferenceCounting<FunctionExpression> alpha)
    : PDEScalarOperatorExpression(unknownName, PDEScalarOperatorExpression::firstOrder)
  {
    __nu[0] = 0; __nu[1] = 0; __nu[2] = 0;

    if (operatorName == "dx") {
      __nu[0] = alpha;
      return;
    }
    if (operatorName == "dy") {
      __nu[1] = alpha;
      return;
    }

    if (operatorName == "dz") {
      __nu[2] = alpha;
      return;
    }
  }

  PDEScalarOperatorExpressionOrderOne(const std::string& unknownName)
    : PDEScalarOperatorExpression(unknownName, PDEScalarOperatorExpression::firstOrder)
  {
    for (size_t i=0; i<3; ++i)
      __nu[i] = 0;
  }

  PDEScalarOperatorExpressionOrderOne(const PDEScalarOperatorExpressionOrderOne& re)
    : PDEScalarOperatorExpression(re)
  {
    for (size_t i=0; i<3; ++i)
      __nu[i] = re.__nu[i];
  }

  ~PDEScalarOperatorExpressionOrderOne()
  {
    ;
  }
};

class PDEVectorialOperatorExpressionOrderTwo;
class PDEVectorialOperatorExpressionOrderOne
  : public PDEVectorialOperatorExpression
{
  friend class PDEVectorialOperatorExpressionOrderTwo;
private:
  typedef
  std::multimap<std::string,
                ReferenceCounting<PDEScalarOperatorExpressionOrderOne> >
  FirstOrderSum;

  FirstOrderSum __sum;
  FirstOrderSum __diff;

protected:
  std::ostream& put(std::ostream& os) const
  {
    bool plus = false;
    for (FirstOrderSum::const_iterator i = __sum.begin();
         i != __sum.end(); ++i) {
      if (plus)
        os << '+';
      os << (*(*i).second);
      plus = true;
    }

    for (FirstOrderSum::const_iterator i = __diff.begin();
         i != __diff.end(); ++i) {
      os << '-' << (*(*i).second);
    }
    return os;
  }

  void __execute()
  {
    throw ErrorHandler(__FILE__,__LINE__,
                       "not implemented",
                       ErrorHandler::unexpected);
  }

public:
  bool __checkBoundaryExpression() const;

  void operator*=(ReferenceCounting<FunctionExpression> f)
  {
    for (FirstOrderSum::iterator i = __sum.begin();
         i != __sum.end(); ++i) {
      (*(*i).second) *= f;
    }

    for (FirstOrderSum::iterator i = __diff.begin();
         i != __diff.end(); ++i) {
      (*(*i).second) *= f;
    }
  }

  void operator+=(ReferenceCounting<PDEScalarOperatorExpressionOrderOne> pdeOp)
  {
    __sum.insert(std::make_pair(pdeOp->unknownName(),pdeOp));
  }

  void operator+=(ReferenceCounting<PDEVectorialOperatorExpressionOrderOne> pdeVectOp)
  {
    for (FirstOrderSum::iterator i = (*pdeVectOp).__sum.begin();
         i != (*pdeVectOp).__sum.end(); ++i) {
      __sum.insert(*i);
    }

    for (FirstOrderSum::iterator i = (*pdeVectOp).__diff.begin();
         i != (*pdeVectOp).__diff.end(); ++i) {
      __diff.insert(*i);
    }
  }

  void operator-=(ReferenceCounting<PDEVectorialOperatorExpressionOrderOne> pdeVectOp)
  {
    for (FirstOrderSum::iterator i = (*pdeVectOp).__sum.begin();
         i != (*pdeVectOp).__sum.end(); ++i) {
      __diff.insert(*i);
    }

    for (FirstOrderSum::iterator i = (*pdeVectOp).__diff.begin();
         i != (*pdeVectOp).__diff.end(); ++i) {
      __sum.insert(*i);
    }
  }

  void preexec()
  {
    for(FirstOrderSum::iterator i = __sum.begin();
        i != __sum.end(); ++i) {
      (*(*i).second).preexec();
    }
    for(FirstOrderSum::iterator i = __diff.begin();
        i != __diff.end(); ++i) {
      (*(*i).second).preexec();
    }
  }

  PDEVectorialOperatorExpressionOrderOne()
    : PDEVectorialOperatorExpression(PDEVectorialOperatorExpression::firstOrder)
  {
    ;
  }

  ~PDEVectorialOperatorExpressionOrderOne()
  {
    ;
  }
};

class PDEVectorialOperatorExpressionOrderTwo
  : public PDEVectorialOperatorExpression
{
private:
  typedef std::multimap<int,
                        ReferenceCounting<PDEVectorialOperatorExpressionOrderOne> > FirstOperatorList;
  FirstOperatorList __firstOrderOp;

  std::ostream& put(std::ostream& os) const
  {
    const char x[3] = {'x', 'y', 'z'};
    for (FirstOperatorList::const_iterator i = __firstOrderOp.begin();
         i != __firstOrderOp.end(); ++i) {
      os << 'd' << x[(*i).first] << '(' << (*(*i).second) << ')';
    }
    return os;
  }

  void __execute();

public:
  bool __checkBoundaryExpression() const;

  void add(const std::string& operatorName,
           ReferenceCounting<PDEVectorialOperatorExpressionOrderOne> pdeOp)
  {
    size_t operatorNumber = std::numeric_limits<size_t>::max();

    if (operatorName == "dx") {
      operatorNumber = 0;
    }
    if (operatorName == "dy") {
      operatorNumber = 1;
    }
    if (operatorName == "dz") {
      operatorNumber = 2;
    }

    if (operatorNumber > 2) {
      throw ErrorHandler(__FILE__,__LINE__,
                         "unknown operator",
                         ErrorHandler::unexpected);
    }

    __firstOrderOp.insert(std::pair<int,
                                    ReferenceCounting<PDEVectorialOperatorExpressionOrderOne> >(operatorNumber,
                                                                                                pdeOp));
  }

  PDEVectorialOperatorExpressionOrderTwo();

  ~PDEVectorialOperatorExpressionOrderTwo();
};

#endif // PDE_OPERATOR_EXPRESSION_HPP

---