Lattice.cpp

   1 #include <stdexcept>
   2
   3 #include "Lattice.h"
   4 #include "SimpleIni.h"
   5
   6 // FIXME: LatticeType should derive operator= from Array<Occupation, 2>
   7 LatticeLayer& LatticeLayer::operator=(Occupation const& o) {
   8         for (LatticeLayer::iterator i = this->begin(); i != this->end(); ++i)
   9                 this->operator()(i.position()) = o;
  10         return *this;
  11 };
  12
  13 ostream& operator<<(ostream& output, const LatticeLayer& layer)
  14 {
  15         for (int y = 0; y < layer.extent(secondDim); ++y) {
  16                 for (int x = 0; x < layer.extent(firstDim); ++x) {
  17                         switch (layer(x,y)) {
  18                         case empty:     cout << "."; break;
  19                         case CO:        cout << "C"; break;
  20                         case Pd:        cout << "P"; break;
  21                         default:        throw out_of_range("Invalid occupation");
  22                         }
  23                 }
  24                 output << endl;
  25         }
  26
  27         return output;
  28 }
  29
  30 Lattice::Lattice(int UnitCellSizeX, int UnitCellSizeY, int GranularityX, int GranularityY) :
  31         GranularityX(GranularityX), GranularityY(GranularityY),
  32         UnitCellSizeX(UnitCellSizeX), UnitCellSizeY(UnitCellSizeY)
  33 {
  34         TinyVector<int, 2> extent(UnitCellSizeX * GranularityX, UnitCellSizeY * GranularityY);
  35         adsorbates.resize(extent);
  36         surface.resize(extent);
  37 }
  38
  39 void
  40 Lattice::assessInteractions(Interactions &interactions)
  41 {
  42         for (LatticeLayer::iterator i = adsorbates.begin(); i != adsorbates.end(); i++) {
  43                 if (adsorbates(i.position()) == empty)
  44                         continue;
  45                 cout << "POS" << i.position() << endl;
  46                 for (Interactions::iterator interaction = interactions.begin();
  47                                 interaction != interactions.end(); interaction++) {
  48                         cout << "checking " << interaction->name << ": ";
  49                         for (Directions::iterator direction = interaction->directions.begin();
  50                                         direction != interaction->directions.end(); direction++) {
  51                                 for (Directions::iterator symmetry_operation = _symmetryOperations.directions.begin();
  52                                         symmetry_operation != _symmetryOperations.directions.end(); symmetry_operation++) {
  53
  54                                         int x = (i.position()[firstDim] + (int)((direction->x / UnitCellSizeX)*
  55                                                                 (symmetry_operation->x * GranularityX)) %
  56                                                                 (UnitCellSizeX * GranularityX) );
  57                                         int y = (i.position()[secondDim] + (int)((direction->y / UnitCellSizeY)*
  58                                                                 (symmetry_operation->y * GranularityY)) %
  59                                                                 (UnitCellSizeY * GranularityY));
  60                                         cout << "(" << x << "/" << y << ") ";
  61                                         if (adsorbates(x, y) != empty) {
  62                                                 interaction->multiplicity++;
  63                                                 cout << "HIT ";
  64                                                 }
  65                                 }
  66                         }
  67                         cout << endl;
  68                 }
  69         }
  70 }
  71
  72
  73 ostream& operator<<(ostream& output, const Lattice& lattice)
  74 {
  75 #define OUTPUT(l) \
  76         for (int y = 0; y < l.extent(secondDim); ++y) { \
  77                 if (y && y % lattice.GranularityY == 0) { \
  78                         for (int i = 0; i < lattice.UnitCellSizeX * \
  79                                         (lattice.GranularityX+1)-1; ++i) { \
  80                                 if (i && (i+1) % (lattice.GranularityX+1) == 0) output << "+"; \
  81                                 else output << "-"; \
  82                         } \
  83                         output << endl; \
  84                 } \
  85                 for (int x = 0; x < l.extent(firstDim); ++x) { \
  86                         if (x && x % lattice.GranularityX == 0) output << "|"; \
  87                         switch (l(x,y)) { \
  88                         case empty:     cout << "."; break; \
  89                         case CO:        cout << "C"; break; \
  90                         case Pd:        cout << "P"; break; \
  91                         default:        throw out_of_range("Invalid occupation"); \
  92                         } \
  93                 } \
  94                 output << endl; \
  95         }
  96
  97
  98         output << "lattice: UnitCellSizeX = " << lattice.UnitCellSizeX
  99                 << ", UnitCellSizeY = " << lattice.UnitCellSizeY
 100                 << ", GranularityX = " << lattice.GranularityX
 101                 << ", GranularityY = " << lattice.GranularityY << endl
 102                 << endl << "surface: " << endl;
 103         OUTPUT(lattice.surface);
 104         output << endl << "adsorbates: " << endl ;
 105         OUTPUT(lattice.adsorbates);
 106         output << endl;
 107
 108         return output;
 109 }