moab
moab::DataCoupler Class Reference

This class couples data between meshes. More...

#include <DataCoupler.hpp>

List of all members.

Public Types

enum  Method { CONSTANT, LINEAR_FE, QUADRATIC_FE, SPECTRAL }
enum  IntegType { VOLUME }

Public Member Functions

 DataCoupler (Interface *impl, Range &source_ents, int coupler_id, ParallelComm *pc=NULL, bool init_locator=true, int dim=-1)
virtual ~DataCoupler ()
ErrorCode locate_points (double *xyz, int num_points, const double rel_iter_tol=1.0e-10, const double abs_iter_tol=1.0e-10, const double inside_tol=1.0e-6)
ErrorCode locate_points (Range &ents, const double rel_iter_tol=1.0e-10, const double abs_iter_tol=1.0e-10, const double inside_tol=1.0e-6)
ErrorCode interpolate (int method, Tag tag, double *interp_vals=NULL, std::vector< int > *point_indices=NULL, bool normalize=true)
ErrorCode interpolate (int method, const std::string &tag_name, double *interp_vals=NULL, std::vector< int > *point_indices=NULL, bool normalize=true)
ErrorCode interpolate (int *methods, const std::string *tag_names, int *points_per_method, int num_methods, double *interp_vals=NULL, std::vector< int > *point_indices=NULL, bool normalize=true)
ErrorCode interpolate (int *methods, Tag *tag_names, int *points_per_method, int num_methods, double *interp_vals=NULL, std::vector< int > *point_indices=NULL, bool normalize=true)
SpatialLocatorspatial_locator ()
int my_id () const
const Rangetarget_ents () const
Rangetarget_ents ()
int get_dim () const

Private Attributes

InterfacembImpl
ParallelCommmyPcomm
SpatialLocatormyLocator
ErrormError
int myId
Range targetEnts
int myDim

Detailed Description

This class couples data between meshes.

The coupler interpolates solution data at a set of points. Data being interpolated resides on a "source" mesh, in a tag or in vertex coords. Applications calling this coupler send in entities, and receive back data interpolated at those points. Entities in the source mesh containing those points do not have to reside on the same processor.

To use, an application should:

  • instantiate this DataCoupler by calling the constructor collectively on all processors in the communicator
  • call locate_points, which locates the points to be interpolated and (optionally) caches the results in this class and SpatialLocator
  • call interpolate, which does the interpolation

Multiple interpolations (of multiple tags, or element-average vs. true interpolation) can be done after locating the points.

SpatialLocator is used for the spatial location portion of this work.

This class is a next-generation implementation of Coupler.

Examples:
DeformMeshRemap.cpp.

Definition at line 38 of file DataCoupler.hpp.


Member Enumeration Documentation

Enumerator:
VOLUME 

Definition at line 44 of file DataCoupler.hpp.

{VOLUME};
Enumerator:
CONSTANT 
LINEAR_FE 
QUADRATIC_FE 
SPECTRAL 

Definition at line 42 of file DataCoupler.hpp.


Constructor & Destructor Documentation

moab::DataCoupler::DataCoupler ( Interface impl,
Range source_ents,
int  coupler_id,
ParallelComm pc = NULL,
bool  init_locator = true,
int  dim = -1 
)

Definition at line 21 of file DataCoupler.cpp.

        : mbImpl(impl), myPcomm(pc), myId(coupler_id), myDim(dim)
{
  assert(NULL != mbImpl && (myPcomm || !source_ents.empty()));

    // now initialize the tree
  if (init_locator) {
    myLocator = new SpatialLocator(mbImpl, source_ents);
    myLocator->elem_eval(new ElemEvaluator(mbImpl));

      // initialize element evaluator with the default for the entity types in source_ents; 
      // can be replaced later by application if desired
    if (!source_ents.empty()) {
      Range::pair_iterator pit = source_ents.pair_begin();
      EntityType last_type = MBMAXTYPE;
      for (; pit != source_ents.pair_end(); pit++) {
        EntityType this_type = mbImpl->type_from_handle(pit->first);
        if (last_type == this_type) continue;
        ErrorCode rval = myLocator->elem_eval()->set_eval_set(pit->first);
        if (MB_SUCCESS != rval) throw(rval);
        last_type = this_type;
      }
    }
  }
  
  if (-1 == dim && !source_ents.empty()) 
    dim = mbImpl->dimension_from_handle(*source_ents.rbegin());

  ErrorCode rval = impl->query_interface(mError);
  if (MB_SUCCESS != rval) throw(rval);
}

Definition at line 60 of file DataCoupler.cpp.

{
  delete myLocator;
}

Member Function Documentation

int moab::DataCoupler::get_dim ( ) const [inline]

Definition at line 213 of file DataCoupler.hpp.

{return myDim;}
ErrorCode moab::DataCoupler::interpolate ( int  method,
Tag  tag,
double *  interp_vals = NULL,
std::vector< int > *  point_indices = NULL,
bool  normalize = true 
) [inline]

Definition at line 246 of file DataCoupler.hpp.

{
    // no point indices input, 
  int num_pts = (point_indices ? point_indices->size() : targetEnts.size());
  return interpolate(&method, &tag, &num_pts, 1,
                     interp_vals, point_indices, normalize);
}
ErrorCode moab::DataCoupler::interpolate ( int  method,
const std::string &  tag_name,
double *  interp_vals = NULL,
std::vector< int > *  point_indices = NULL,
bool  normalize = true 
)

Definition at line 91 of file DataCoupler.cpp.

{
    // tag name input, translate to tag handle and pass down the chain

    // not inlined because of call to set_last_error, class Error isn't in public interface
  Tag tag;
  ErrorCode result = mbImpl->tag_get_handle(interp_tag.c_str(), tag);
  if (MB_SUCCESS != result) {
    std::ostringstream str;
    str << "Failed to get handle for interpolation tag \"" << interp_tag << "\"";
    mError->set_last_error(str.str());
    return result;
  }
  return interpolate(method, tag, interp_vals, point_indices, normalize);
}
ErrorCode moab::DataCoupler::interpolate ( int *  methods,
const std::string *  tag_names,
int *  points_per_method,
int  num_methods,
double *  interp_vals = NULL,
std::vector< int > *  point_indices = NULL,
bool  normalize = true 
)
ErrorCode moab::DataCoupler::interpolate ( int *  methods,
Tag tag_names,
int *  points_per_method,
int  num_methods,
double *  interp_vals = NULL,
std::vector< int > *  point_indices = NULL,
bool  normalize = true 
)

Definition at line 111 of file DataCoupler.cpp.

{
    // lowest-level interpolate function, does actual interpolation using calls to ElemEvaluator

  ErrorCode result = MB_SUCCESS;

  unsigned int pts_total = 0;
  for (int i = 0; i < num_methods; i++) pts_total += (points_per_method ? points_per_method[i] : targetEnts.size());

  unsigned int num_indices = (point_indices ? point_indices->size() : targetEnts.size());
    // # points and indices should be identical
  if (pts_total != num_indices)
    return MB_FAILURE;

    // since each tuple contains one interpolated tag, if we're interpolating multiple tags, every tuple
    // needs to be able to store up to max tag size
  int max_tsize = -1;
  for (int i = 0; i < num_methods; i++) {
    int tmp_tsize;
    result = mbImpl->tag_get_length(tags[i], tmp_tsize);
    if (MB_SUCCESS != result) return MB_FAILURE;
    max_tsize = std::max(max_tsize, tmp_tsize);
  }

  if (myPcomm && myPcomm->size() > 1) {
      // build up the tuple list to distribute from my target points; assume that
      // all procs use same method/tag input
    TupleList TLob; // TLob structure: (pto_i, ridx_i, lidx_i, meth_tagidx_i)
    TLob.initialize(4, 0, 0, 0, num_indices);
    int tn = 0; // the tuple number we're currently on
    TLob.enableWriteAccess();
    for (int m = 0; m < num_methods; m++) {
      int num_points = (points_per_method ? points_per_method[m] : targetEnts.size());
      for (int j = 0; j < num_points; j++) {
        int idx = (point_indices ? (*point_indices)[j] : j); // the index in my targetEnts for this interpolation point

          // remote proc/idx from myLocator->parLocTable
        TLob.vi_wr[4*tn]   = myLocator->par_loc_table().vi_rd[2*idx]; // proc
        TLob.vi_wr[4*tn+1] = myLocator->par_loc_table().vi_rd[2*idx+1]; // remote idx
  
          // local entity index, tag/method index from my data
        TLob.vi_wr[4*tn+2] = idx;
        TLob.vi_wr[4*tn+3] = m;
        TLob.inc_n();
        tn++;
      }
    }

      // scatter/gather interpolation points
    myPcomm->proc_config().crystal_router()->gs_transfer(1, TLob, 0);

      // perform interpolation on local source mesh; put results into TLinterp
    TupleList TLinterp; // TLinterp structure: (pto_i, ridx_i, vals[max_tsize]_d)
    TLinterp.initialize(2, 0, 0, max_tsize, TLob.get_n());
    TLinterp.set_n(TLob.get_n()); // set the size right away
    TLinterp.enableWriteAccess();

    for (unsigned int i = 0; i < TLob.get_n(); i++) {
      int lidx = TLob.vi_rd[4*i+1]; // index into myLocator's local table
        // method and tag indexed with same index
//      /*Method*/ int method = methods[TLob.vi_rd[4*i+3]];
      Tag tag = tags[TLob.vi_rd[4*i+3]];
        // copy proc/remote index from TLob
      TLinterp.vi_wr[2*i] = TLob.vi_rd[4*i];
      TLinterp.vi_wr[2*i+1] = TLob.vi_rd[4*i+2]; 

        // set up the evaluator for the tag and entity, then interpolate, putting result in TLinterp
      myLocator->elem_eval()->set_tag_handle(tag);
      myLocator->elem_eval()->set_ent_handle(myLocator->loc_table().vul_rd[lidx]);
      result = myLocator->elem_eval()->eval(myLocator->loc_table().vr_rd+3*lidx, TLinterp.vr_rd+i*max_tsize);
      if (MB_SUCCESS != result) return result;
    }

      // scatter/gather interpolation data
    myPcomm->proc_config().crystal_router()->gs_transfer(1, TLinterp, 0);

      // copy the interpolated field as a unit
    std::copy(TLinterp.vr_rd, TLinterp.vr_rd+TLinterp.get_n()*max_tsize, interp_vals);
  }
  else {
      // if called serially, interpolate directly from local locations/entities,
      // into either interp_vals or by setting data directly on tags on those entities
    std::vector<double> tmp_vals;
    std::vector<EntityHandle> tmp_ents;
    double *tmp_dbl = interp_vals;
    for (int i = 0; i < num_methods; i++) {
      int num_points = (points_per_method ? points_per_method[i] : targetEnts.size());

        // interpolated data is tsize long, which is either max size (if data passed back to caller in interp_vals)
        // or tag size (if data will be set on entities, in which case it shouldn't have padding)
        // set sizes here, inside loop over methods, to reduce memory usage
      int tsize = max_tsize, tsize_bytes = 0;
      if (!interp_vals) {
        tmp_vals.resize(num_points*max_tsize);
        tmp_dbl = &tmp_vals[0];
        tmp_ents.resize(num_points);
        result = mbImpl->tag_get_length(tags[i], tsize);
        if (MB_SUCCESS != result) return result;
        result = mbImpl->tag_get_bytes(tags[i], tsize_bytes);
        if (MB_SUCCESS != result) return result;
      }
      
      for (int j = 0; j < num_points; j++) {
        int lidx;
        if (point_indices) {
          lidx = (*point_indices)[j];
        }
        else {
          lidx = j;
        }

        myLocator->elem_eval()->set_tag_handle(tags[i]);
        myLocator->elem_eval()->set_ent_handle(myLocator->loc_table().vul_rd[lidx]);
        if (!interp_vals) tmp_ents[j] = targetEnts[lidx]; // could be performance-sensitive, thus the if test
        result = myLocator->elem_eval()->eval(myLocator->loc_table().vr_rd+3*lidx, tmp_dbl);
        tmp_dbl += tsize;
        if (MB_SUCCESS != result) return result;
      } // for j over points

      if (!interp_vals) {
          // set tags on tmp_ents; data is already w/o padding, due to tsize setting above
        result = mbImpl->tag_set_data(tags[i], &tmp_ents[0], tmp_ents.size(), &tmp_vals[0]);
        if (MB_SUCCESS != result) return result;
      }

    } // for i over methods
  } // if myPcomm
  
      // done
  return MB_SUCCESS;
}
ErrorCode moab::DataCoupler::locate_points ( double *  xyz,
int  num_points,
const double  rel_iter_tol = 1.0e-10,
const double  abs_iter_tol = 1.0e-10,
const double  inside_tol = 1.0e-6 
)

Definition at line 79 of file DataCoupler.cpp.

{
#ifdef USE_MPI
  if (myPcomm && myPcomm->size() > 1) 
    return myLocator->par_locate_points(myPcomm, xyz, num_points, rel_iter_tol, abs_iter_tol, inside_tol);
#endif

  return myLocator->locate_points(xyz, num_points, rel_iter_tol, abs_iter_tol, inside_tol);
}
ErrorCode moab::DataCoupler::locate_points ( Range ents,
const double  rel_iter_tol = 1.0e-10,
const double  abs_iter_tol = 1.0e-10,
const double  inside_tol = 1.0e-6 
)

Definition at line 65 of file DataCoupler.cpp.

{
  targetEnts = targ_ents;

#ifdef USE_MPI
  if (myPcomm && myPcomm->size() > 1) 
    return myLocator->par_locate_points(myPcomm, targ_ents, rel_iter_tol, abs_iter_tol, inside_tol);
#endif

  return myLocator->locate_points(targ_ents, rel_iter_tol, abs_iter_tol, inside_tol);
}
int moab::DataCoupler::my_id ( ) const [inline]

Definition at line 210 of file DataCoupler.hpp.

{return myId;}

Definition at line 209 of file DataCoupler.hpp.

{return myLocator;}
const Range& moab::DataCoupler::target_ents ( ) const [inline]

Definition at line 211 of file DataCoupler.hpp.

{return targetEnts;}

Definition at line 212 of file DataCoupler.hpp.

{return targetEnts;}

Member Data Documentation

Definition at line 219 of file DataCoupler.hpp.

Definition at line 231 of file DataCoupler.hpp.

int moab::DataCoupler::myDim [private]

Definition at line 242 of file DataCoupler.hpp.

int moab::DataCoupler::myId [private]

Definition at line 235 of file DataCoupler.hpp.

Definition at line 227 of file DataCoupler.hpp.

Definition at line 223 of file DataCoupler.hpp.

Definition at line 239 of file DataCoupler.hpp.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines