source_sampling.h

Go to the documentation of this file.

 
#ifndef PYNE_6OR6BJURKJHHTOFWXO2VMQM5EY
#define PYNE_6OR6BJURKJHHTOFWXO2VMQM5EY

#include <assert.h>
#include <iostream>
#include <fstream>
#include <stdio.h>
#include <stdlib.h>
#include <vector>
#include <stdexcept> 
#include <sstream>
#include <string>

#include "moab/Range.hpp"
#include "moab/Core.hpp"
#include "measure.h"
#include "moab/CartVect.hpp"

#ifdef __cplusplus
extern "C" {
#endif

namespace pyne {

  void sampling_setup_(int* mode);
  void particle_birth_(double* rands,
                            double* x,
                            double* y,
                            double* z,
                            double* e,
                            double* w);
  std::vector<double> read_e_bounds(std::string e_bounds_file);

  struct edge_points {
    moab::CartVect o_point;
    moab::CartVect x_vec;
    moab::CartVect y_vec;
    moab::CartVect z_vec;
  };
  
  class AliasTable {
  public:
    AliasTable(std::vector<double> p);
    int sample_pdf(double rand1, double rand2);
    ~AliasTable(){};
    int n; 
    std::vector<double> prob; 
    std::vector<int> alias; 
  };
  
  enum Mode {USER, ANALOG, UNIFORM};
  
  class Sampler {
  public:
    Sampler(std::string filename, 
            std::string src_tag_name, 
            std::vector<double> e_bounds, 
            bool uniform);
    Sampler(std::string filename, 
            std::string src_tag_name, 
            std::vector<double> e_bounds, 
            std::string bias_tag_name);
    std::vector<double> particle_birth(std::vector<double> rands);
    ~Sampler() {
      delete mesh;
      delete at;
    };
  
  // member variables
  private:
    // problem parameters
    std::string filename; 
    std::string src_tag_name; 
    std::string bias_tag_name; 
    std::vector<double> e_bounds;  
    int num_e_groups; 
    int num_bias_groups; 
    Mode mode; 
    // mesh
    moab::Interface* mesh; 
    int num_ves; 
    moab::EntityType ve_type; 
    int verts_per_ve; 
    // sampling
    std::vector<edge_points> all_edge_points; 
    std::vector<double> biased_weights; 
    AliasTable* at; 
  
  // member functions
  private:
    // instantiation
    void setup();
    void mesh_geom_data(moab::Range ves, std::vector<double> &volumes);
    void mesh_tag_data(moab::Range ves, const std::vector<double> volumes);
    // select birth parameters
    moab::CartVect sample_xyz(int ve_idx, std::vector<double> rands);
    double sample_e(int e_idx, double rand);
    double sample_w(int pdf_idx);
    // helper functions
    void normalize_pdf(std::vector<double> & pdf);
    int num_groups(moab::Tag tag);
    std::vector<double> read_bias_pdf(moab::Range ves, std::vector<double> volumes, 
                                      std::vector<double> pdf);
  };
} //end namespace pyne

#ifdef __cplusplus
} // extern "C"
#endif

#endif // PYNE_6OR6BJURKJHHTOFWXO2VMQM5EY