material.h

#ifndef PYNE_MR34UE5INRGMZK2QYRDWICFHVM
#define PYNE_MR34UE5INRGMZK2QYRDWICFHVM

#include <iostream>
#include <fstream>
#include <string>
#include <map>
#include <set>
#include <stdio.h>
#include <stdlib.h>
#include <sstream>  // std::ostringstream

#if !defined(JSON_IS_AMALGAMATION)
  #define JSON_IS_AMALGAMATION
#endif

#ifndef PYNE_IS_AMALGAMATED
#include "json-forwards.h"
#include "json.h"
#include "h5wrap.h"
#include "utils.h"
#include "nucname.h"
#include "data.h"
#include "decay.h"
#endif

namespace pyne
{
  // Set Type Definitions
  typedef std::map<int, double> comp_map; 
  typedef comp_map::iterator comp_iter;   

  #ifdef PYNE_IS_AMALGAMATED
  namespace decayers {
    extern comp_map decay(comp_map, double);
  }  // namespace decayers
  #endif


  // These 37 strings are predefined FLUKA materials. 
  // Materials not on this list requires a MATERIAL card. 
  static std::string fluka_mat_strings[] = {
   "BLCKHOLE", "VACUUM",   "HYDROGEN", "HELIUM",   "BERYLLIU", "CARBON", 
   "NITROGEN", "OXYGEN",   "MAGNESIU", "ALUMINUM", "IRON",     "COPPER", 
   "SILVER",   "SILICON",  "GOLD",     "MERCURY",  "LEAD",     "TANTALUM", 
   "SODIUM",   "ARGON",    "CALCIUM",  "TIN",      "TUNGSTEN", "TITANIUM", 
   "NICKEL",   "WATER",    "POLYSTYR", "PLASCINT", "PMMA",     "BONECOMP", 
   "BONECORT", "MUSCLESK", "MUSCLEST", "ADTISSUE", "KAPTON", "POLYETHY", "AIR"
  };

  static int FLUKA_MAT_NUM = 37;

  class Material
  {
  protected:

    double get_comp_sum ();

  public:

    // Material Constructors
    Material ();  
    Material(comp_map cm, double m=-1.0, double d=-1.0, double apm=-1.0,
             Json::Value attributes=Json::Value(Json::objectValue));
    Material(char * filename, double m=-1.0, double d=-1.0, double apm=-1.0,
             Json::Value attributes=Json::Value(Json::objectValue));
    Material(std::string filename, double m=-1.0, double d=-1.0, double apm=-1.0,
             Json::Value attributes=Json::Value(Json::objectValue));
    ~Material (); 

    void norm_comp ();

    // Persistence functions.

    void _load_comp_protocol0(hid_t db, std::string datapath, int row);

    void _load_comp_protocol1(hid_t db, std::string datapath, int row);

    void from_hdf5(char * filename, char * datapath, int row=-1, int protocol=1);

    void from_hdf5(std::string filename, std::string datapath="/material",
                                                          int row=-1, int protocol=1);

    void write_hdf5(char * filename, char * datapath, char * nucpath, float row=-0.0,
                                                                    int chunksize=100);
    void write_hdf5(std::string filename, std::string datapath="/material",
                    std::string nucpath="/nucid", float row=-0.0, int chunksize=100);

    std::string mcnp(std::string frac_type = "mass");
    std::string fluka(int id, std::string frac_type = "mass");
    bool not_fluka_builtin(std::string fluka_name);
    std::string fluka_material_str(int id);
    std::string fluka_material_component(int fid, int nucid, 
                                         std::string fluka_name);
    std::string fluka_material_line(int znum, double atomic_mass, 
                              int fid, std::string fluka_name);
    std::string fluka_format_field(float field);
    std::string fluka_compound_str(int id, std::string frac_type = "mass");

    void from_text(char * filename);
    void from_text(std::string filename);

    void write_text(char * filename);
    void write_text(std::string filename);

    void load_json(Json::Value);
    Json::Value dump_json();
    void from_json(char * filename);
    void from_json(std::string filname);
    void write_json(char * filename);
    void write_json(std::string filename);

    // Fundemental mass stream data
    comp_map comp;
    double mass;  
    double density; 
    double atoms_per_molecule; 
    Json::Value metadata;

    // Material function definitions
    void normalize ();  
    comp_map mult_by_mass();
    double molecular_mass(double apm=-1.0);
    comp_map activity();
    comp_map decay_heat();
    comp_map dose_per_g(std::string dose_type, int source=0);
    Material expand_elements();
    // Returns a copy of the current material where all the isotopes of the elements
    // are added up, atomic-fraction-wise, unless they are in the exception set
    Material collapse_elements(std::set<int> exception_znum);
    // Wrapped version to facilitate calling from python
    Material collapse_elements(int **int_ptr_arry);
    // void print_material( pyne::Material test_mat);
    double mass_density(double num_dens=-1.0, double apm=-1.0);
    double number_density(double mass_dens=-1.0, double apm=-1.0);

    // Sub-Stream Computation
    Material sub_mat(std::set<int> nucset);
    Material sub_mat(std::set<std::string> nucset);

    Material set_mat(std::set<int> nucset, double value);
    Material set_mat(std::set<std::string> nucset, double value);

    Material del_mat(std::set<int> nucset);
    Material del_mat(std::set<std::string> nucset);

    Material sub_range(int lower=0, int upper=10000000);
    Material set_range(int lower=0, int upper=10000000, double value=0.0);
    Material del_range(int lower=0, int upper=10000000);

    Material sub_elem(int element);
    Material sub_lan();
    Material sub_act();
    Material sub_tru();
    Material sub_ma();
    Material sub_fp();

    // Atom fraction functions
    std::map<int, double> to_atom_frac();
    void from_atom_frac(std::map<int, double> atom_fracs);
    
    std::map<int, double> to_atom_dens();

    // Radioactive Material functions
    std::vector<std::pair<double, double> > gammas();
    std::vector<std::pair<double, double> > xrays();
    std::vector<std::pair<double, double> > photons(bool norm);
    std::vector<std::pair<double, double> > normalize_radioactivity(
      std::vector<std::pair<double, double> > unnormed);

    Material decay(double t);

    // Overloaded Operators
    Material operator+ (double);
    Material operator+ (Material);
    Material operator* (double);
    Material operator/ (double);
  };

  std::ostream& operator<< (std::ostream& os, Material mat);

  typedef struct material_data {
    double mass;  
    double density; 
    double atoms_per_mol; 
    double comp[1]; 
  } material_data;

  class MaterialProtocolError: public std::exception
  {
    virtual const char* what() const throw()
    {
      return "Invalid loading protocol number; please use 0 or 1.";
    }
  };

// End pyne namespace
}

#endif  // PYNE_MR34UE5INRGMZK2QYRDWICFHVM