F77GFMBoundary.h

Go to the documentation of this file.

// -*- C++ -*-

// Copyright (C) 2003-2007 California Institute of Technology
// Ralf Deiterding, ralf@amroc.net

#ifndef AMROC_F77_GFMBOUNDARY_H
#define AMROC_F77_GFMBOUNDARY_H

#include "GFMBoundary.h"

template <class VectorType, int dim>
class F77GFMBoundary : public GFMBoundary<VectorType,dim> {
  typedef typename VectorType::InternalDataType DataType;
  typedef GFMBoundary<VectorType,dim> base; 

public:
  typedef typename base::vec_grid_fct_type vec_grid_fct_type;  
  typedef typename base::vec_grid_data_type vec_grid_data_type;
  typedef typename base::grid_fct_type grid_fct_type;  
  typedef typename base::grid_data_type grid_data_type;
  typedef typename base::point_type point_type;

  typedef generic_fortran_func generic_func_type;

  typedef void (*bnd_1_func_type) ( FI(1,VectorType), const INTEGER &, const INTEGER &,
                                    const INTEGER *, const DOUBLE *, const DOUBLE *,
                                    const DOUBLE *, const DOUBLE *, const INTEGER &, 
                                    const DOUBLE *, const DOUBLE *, const DOUBLE &);

  typedef void (*bnd_2_func_type) ( FI(2,VectorType), const INTEGER &, const INTEGER &,
                                    const INTEGER *, const DOUBLE *, const DOUBLE *,
                                    const DOUBLE *, const DOUBLE *, const INTEGER &, 
                                    const DOUBLE *, const DOUBLE *, const DOUBLE &);

  typedef void (*bnd_3_func_type) ( FI(3,VectorType), const INTEGER &, const INTEGER &,
                                    const INTEGER *, const DOUBLE *, const DOUBLE *,
                                    const DOUBLE *, const DOUBLE *, const INTEGER &, 
                                    const DOUBLE *, const DOUBLE *, const DOUBLE &);

  typedef void (*trans_1_func_type) ( const INTEGER& mx, 
                                      const INTEGER& meqn, VectorType* q, VectorType* qt);

  typedef void (*trans_2_func_type) ( const INTEGER& mx, const INTEGER& my,
                                      const INTEGER& meqn, VectorType* q, VectorType* qt);

  typedef void (*trans_3_func_type) ( const INTEGER& mx, const INTEGER& my, const INTEGER& mz,
                                      const INTEGER& meqn, VectorType* q, VectorType* qt);

  typedef void (*bndaux_func_type) ( const INTEGER &, const INTEGER &, const DOUBLE *, 
                                     const DOUBLE *, const DOUBLE *, const DOUBLE *, 
                                     const INTEGER &, const DOUBLE *, const DOUBLE *, 
                                     const DOUBLE &);

public:
  F77GFMBoundary(generic_func_type bnd, generic_func_type trs) : 
    base(), f_intbnd(bnd), f_inttrans(trs), f_intaux(0) {} 
  F77GFMBoundary(generic_func_type bnd, generic_func_type trs, generic_func_type aux) : 
    base(), f_intbnd(bnd), f_inttrans(trs), f_intaux(aux) {} 

  virtual ~F77GFMBoundary() {}
  
  virtual void SetGrid(vec_grid_data_type& gdu, grid_data_type& gdphi, const BBox& bb,
                       const int& Level, double t, const int& nc, 
                       const int* idx, const point_type* xc, 
                       DataType* distance, point_type* normal) {

    Coords ex = gdu.extents();
    VectorType* u = (VectorType*) 0;
    DataType* aux = (DataType*) 0;
    DCoords dx = base::GH().worldStep(gdu.stepsize());
    if (base::NAux()>0)
      aux = new DataType[nc*base::NAux()];

    if (base::NormalFactor()>0)
      if (f_inttrans) {
        START_WATCH
        vec_grid_data_type gdhelp(gdu.bbox());  
        if (dim == 1) 
          ((trans_1_func_type) f_inttrans)(AA(1,ex()),base::NEquations(),gdu.data(),gdhelp.data());
        else if (dim == 2) 
          ((trans_2_func_type) f_inttrans)(AA(2,ex()),base::NEquations(),gdu.data(),gdhelp.data());
        else if (dim == 3) 
          ((trans_3_func_type) f_inttrans)(AA(3,ex()),base::NEquations(),gdu.data(),gdhelp.data());
        END_WATCH(GFM_TRANSFORM)

        base::Extrapolation(gdhelp,gdphi,u,Level,nc,idx,xc,
                            distance,normal,base::NormalFactor());
        if (base::NAux()>0) {
          START_WATCH
            SetBndryAux(gdhelp,gdphi,bb,u,aux,Level,t,nc,idx,xc,distance,normal);
          END_WATCH(GFM_AUXILIARY_VALUES)
        }
      }
      else {
        base::Extrapolation(gdu,gdphi,u,Level,nc,idx,xc,
                            distance,normal,base::NormalFactor());
        if (base::NAux()>0) {
          START_WATCH
            SetBndryAux(gdu,gdphi,bb,u,aux,Level,t,nc,idx,xc,distance,normal);
          END_WATCH(GFM_AUXILIARY_VALUES)
        }
      }
    else if (base::NAux()>0) {
      START_WATCH
        SetBndryAux(gdu,gdphi,bb,u,aux,Level,t,nc,idx,xc,distance,normal);
      END_WATCH(GFM_AUXILIARY_VALUES)
    }
    
    if (f_intbnd) {
        START_WATCH
        if (dim == 1) 
          ((bnd_1_func_type) f_intbnd)(FORTRAN_ARGS(1,gdu),base::NEquations(),nc,
                                       idx,u[0].data(),xc[0].data(),distance,
                                       normal[0].data(),base::NAux(),aux,dx(),t);
        else if (dim == 2) 
          ((bnd_2_func_type) f_intbnd)(FORTRAN_ARGS(2,gdu),base::NEquations(),nc,
                                       idx,u[0].data(),xc[0].data(),distance,
                                       normal[0].data(),base::NAux(),aux,dx(),t);
        else if (dim == 3) 
          ((bnd_3_func_type) f_intbnd)(FORTRAN_ARGS(3,gdu),base::NEquations(),nc,
                                       idx,u[0].data(),xc[0].data(),distance,
                                       normal[0].data(),base::NAux(),aux,dx(),t);
        END_WATCH(GFM_SETBNDRY)
    }    
    if (u) delete [] u;
    if (aux) delete [] aux;
  }

  virtual void SetBndryAux(vec_grid_data_type& gdu, grid_data_type& gdphi, const BBox& bb,
                           const VectorType* u, DataType* aux, const int& Level, 
                           double t, const int& nc, const int* idx, 
                           const point_type* xc, const DataType* distance, 
                           const point_type* normal) {
    if (f_intaux) {
      DCoords dx = base::GH().worldStep(gdu.stepsize());
      ((bndaux_func_type) f_intaux)(base::NEquations(),nc,u[0].data(),
                                    xc[0].data(),distance,normal[0].data(),
                                    base::NAux(),aux,dx(),t);
    }
    else 
      for (register int i=0; i<nc*base::NAux(); i++)
        aux[i] = static_cast<DataType>(0.);
  }

  inline void SetBndFunc(generic_func_type bnd) { f_intbnd = bnd; }
  generic_func_type GetBndFunc() const { return f_intbnd; }
  inline void SetTransFunc(generic_func_type trans) { f_inttrans = trans; }
  generic_func_type GetTransFunc() const { return f_inttrans; }
  inline void SetAuxFunc(generic_func_type aux) { f_intaux = aux; }
  generic_func_type GetAuxFunc() const { return f_intaux; }

protected:
  generic_func_type f_intbnd, f_inttrans, f_intaux;
};

#endif

---