ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/applications/nanoparticleBuilder/icosahedralBuilder.cpp
(Generate patch)

Comparing branches/development/src/applications/nanoparticleBuilder/icosahedralBuilder.cpp (file contents):
Revision 1861 by gezelter, Tue Apr 9 19:45:54 2013 UTC vs.
Revision 1862 by gezelter, Fri Apr 12 20:45:04 2013 UTC

# Line 47 | Line 47
47   #include "brains/SimInfo.hpp"
48   #include "brains/SimCreator.hpp"
49   #include "io/DumpWriter.hpp"
50 < #include "math/SquareMatrix3.hpp"
50 > #include "clusters/Icosahedron.hpp"
51 > #include "clusters/Decahedron.hpp"
52  
53   using namespace OpenMD;
54   using namespace std;
54
55 //
56 // Create Mackay icosaheron structure.
57 //
58 // Heavily modified from a code created by: Yanting Wang    07/21/2003
59 //
60
61 vector<Vector3d> Points;
62 vector<pair<int, int> > Edges;
63 vector<tuple3<int, int, int> > Facets;
64 vector<Vector3d> Basis; // Basis vectors of the edges
65
66 //
67 // function np
68 //
69 // Calculate number of particles on the nth layer.
70 //
71
72 int np( int n ) {
73  if( n<0 ) return -1;
74  else if( n==0 ) return 1;
75  else if( n==1 ) return 12;
76  else if( n==2 ) return 42;
77  else {
78    int count = 0;    
79    count += 12;   // edge particles
80    count += (n-1)*30;   // side particles
81    for( int i = 1; i <= n-2; i++ ) count += i*20;   // body particles
82    return count;
83  }
84 }
85
86 //
87 // function init
88 //
89 // Initialize some constants.
90 //
91
92 void init() {
93  
94  Basis.clear();
95  Edges.clear();
96  Facets.clear();
97  Points.clear();
98  
99  //
100  // Initialize Basis vectors.
101  //
102  const RealType HGR = ( sqrt(5.0) + 1.0 ) / 4.0;   // half of the golden ratio
103  
104  Basis.push_back( Vector3d(  HGR,  0.0,  0.5 ));
105  Basis.push_back( Vector3d(  HGR,  0.0, -0.5 ));
106  Basis.push_back( Vector3d(  0.5,  HGR,  0.0 ));
107  Basis.push_back( Vector3d( -0.5,  HGR,  0.0 ));
108  Basis.push_back( Vector3d(  0.0,  0.5,  HGR ));
109  Basis.push_back( Vector3d(  0.0, -0.5,  HGR ));
110  Basis.push_back( Vector3d(  0.5, -HGR,  0.0 ));
111  Basis.push_back( Vector3d(  0.0,  0.5, -HGR ));
112  Basis.push_back( Vector3d( -HGR,  0.0,  0.5 ));
113  Basis.push_back( Vector3d(  0.0, -0.5, -HGR ));
114  Basis.push_back( Vector3d( -HGR,  0.0, -0.5 ));
115  Basis.push_back( Vector3d( -0.5, -HGR,  0.0 ));
116  
117  //
118  // Initialize 30 edges
119  //
120  
121  Edges.push_back(std::make_pair(0, 1));
122  Edges.push_back(std::make_pair(0, 2));
123  Edges.push_back(std::make_pair(0, 4));
124  Edges.push_back(std::make_pair(0, 5));
125  Edges.push_back(std::make_pair(0, 6));
126  
127  Edges.push_back(std::make_pair(10, 3));
128  Edges.push_back(std::make_pair(10, 7));
129  Edges.push_back(std::make_pair(10, 8));
130  Edges.push_back(std::make_pair(10, 9));
131  Edges.push_back(std::make_pair(10, 11));
132  
133  Edges.push_back(std::make_pair(1, 2));
134  Edges.push_back(std::make_pair(1, 6));
135  Edges.push_back(std::make_pair(1, 7));
136  Edges.push_back(std::make_pair(1, 9));
137  
138  Edges.push_back(std::make_pair(8, 3));
139  Edges.push_back(std::make_pair(8, 4));
140  Edges.push_back(std::make_pair(8, 5));
141  Edges.push_back(std::make_pair(8, 11));
142  
143  Edges.push_back(std::make_pair(2, 3));
144  Edges.push_back(std::make_pair(2, 4));
145  Edges.push_back(std::make_pair(2, 7));
146  
147  Edges.push_back(std::make_pair(11, 5));
148  Edges.push_back(std::make_pair(11, 6));
149  Edges.push_back(std::make_pair(11, 9));
150  
151  Edges.push_back(std::make_pair(6, 5));
152  Edges.push_back(std::make_pair(6, 9));
153  
154  Edges.push_back(std::make_pair(3, 4));
155  Edges.push_back(std::make_pair(3, 7));
156  
157  Edges.push_back(std::make_pair(7, 9));
158  
159  Edges.push_back(std::make_pair(5, 4));
160  
161  //
162  // Initialize 20 facets
163  //
164  
165  Facets.push_back(make_tuple3(0, 1, 2));
166  Facets.push_back(make_tuple3(0, 2, 4));
167  Facets.push_back(make_tuple3(0, 4, 5));
168  Facets.push_back(make_tuple3(0, 5, 6));
169  Facets.push_back(make_tuple3(0, 1, 6));
170  
171  Facets.push_back(make_tuple3(10, 3, 7));
172  Facets.push_back(make_tuple3(10, 3, 8));
173  Facets.push_back(make_tuple3(10, 8, 11));
174  Facets.push_back(make_tuple3(10, 9, 11));
175  Facets.push_back(make_tuple3(10, 7, 9));
176                                    
177  Facets.push_back(make_tuple3(1, 2, 7));
178  Facets.push_back(make_tuple3(1, 7, 9));
179  Facets.push_back(make_tuple3(1, 6, 9));
180  
181  Facets.push_back(make_tuple3(8, 5, 11));
182  Facets.push_back(make_tuple3(8, 4, 5));
183  Facets.push_back(make_tuple3(8, 3, 4));
184  
185  Facets.push_back(make_tuple3(2, 3, 7));
186  Facets.push_back(make_tuple3(2, 3, 4));
187  
188  Facets.push_back(make_tuple3(11, 5, 6));
189  Facets.push_back(make_tuple3(11, 6, 9));
190 }
55  
192 //
193 // function ih
194 //
195 // Create nth layer particles.
196 // The distance between nearest neighbors has the unit length of 1.
197
198 vector<Vector3d> ih( int n ) {
199
200  if( n < 0 ) return Points;
201  
202  if( n==0 ) {
203    // center particle only
204    Points.push_back(Vector3d( 0.0, 0.0, 0.0 ));
205    return Points;
206  }
207  
208  //
209  // Generate edge particles
210  //
211  for( vector<Vector3d>::iterator i = Basis.begin(); i != Basis.end(); ++i ) {
212    
213    Points.push_back( (*i) * RealType(n) );
214  }
215  
216  //
217  // Generate side particles
218  //
219  
220  if( n<2 ) return Points;
221  
222  for( vector<pair<int,int> >::iterator i=Edges.begin();
223       i != Edges.end(); ++i ) {
224    
225    Vector3d e1 = Basis[ (*i).first  ] * RealType(n);
226    Vector3d e2 = Basis[ (*i).second ] * RealType(n);
227    
228    for( int j = 1; j <= n-1; j++ ) {
229      Points.push_back( e1 + (e2-e1) * RealType(j) / RealType(n));
230    }      
231  }
232  
233  //
234  // Generate body particles
235  //
236  
237  if( n<3 ) return Points;
238  
239  for( vector<tuple3<int,int,int> >::iterator i = Facets.begin();
240       i != Facets.end(); ++i) {
241    
242    Vector3d e1 = Basis[ (*i).first  ] * RealType(n);
243    Vector3d e2 = Basis[ (*i).second ] * RealType(n);
244    Vector3d e3 = Basis[ (*i).third  ] * RealType(n);
245    
246    for( int j=1; j<=n-2; j++ ) {
247      
248      Vector3d v1 = e1 + (e2-e1) * RealType(j+1) / RealType(n);
249      Vector3d v2 = e1 + (e3-e1) * RealType(j+1) / RealType(n);
250      
251      for( int k=1; k<=j; k++ ) {
252        Points.push_back(v1 + (v2-v1) * RealType(k) / RealType(j+1));
253      }
254    }
255  }
256  return Points;
257 }
258
259
56   void createMdFile(const std::string&oldMdFileName,
57                    const std::string&newMdFileName,
58                    int nMol) {
# Line 297 | Line 93 | int main(int argc, char *argv []) {
93    MoLocator* locator;
94    RealType latticeConstant;
95    int nShells;
300
301  DumpWriter *writer;
96  
97 <  init();
97 >  DumpWriter *writer;
98    
99    // Parse Command Line Arguments
100    if (cmdline_parser(argc, argv, &args_info) != 0)
# Line 336 | Line 130 | int main(int argc, char *argv []) {
130  
131    if (args_info.latticeConstant_given) {    
132      latticeConstant = args_info.latticeConstant_arg;
133 <  } else  {
340 <    
341 <    int count=0;
342 <    for( int i = 0; i <= nShells; i++ ) count += np( i );
133 >  } else  {  
134      sprintf(painCave.errMsg, "icosahedralBuilder:  No lattice constant\n"
135 <            "\tgiven.  Total number of atoms in a Mackay Icosahedron with\n"
345 <            "\t%d shells is %d.", nShells, count);
135 >            "\tgiven.");
136      painCave.isFatal = 1;
347    painCave.severity = OPENMD_INFO;
137      cmdline_parser_print_help();
138      simError();
139    }
# Line 352 | Line 141 | int main(int argc, char *argv []) {
141    /* parse md file and set up the system */
142    SimCreator oldCreator;
143    SimInfo* oldInfo = oldCreator.createSim(inputFileName, false);
355
144    Globals* simParams = oldInfo->getSimParams();
145  
146 <  //generate the coordinates
147 <  for( int i = 0; i <= nShells; i++ ) ih( i );
146 >  vector<Vector3d> Points;
147 >  if (args_info.ico_given){
148 >    Icosahedron* ico = new Icosahedron();
149 >    Points = ico->getPoints(nShells);
150 >  } else if (args_info.deca_given) {
151 >    RegularDecahedron* deca = new RegularDecahedron(nShells);
152 >    Points = deca->getPoints();
153 >  } else if (args_info.ino_given) {
154 >    int columnAtoms = args_info.columnAtoms_arg;    
155 >    InoDecahedron* ino = new InoDecahedron(columnAtoms, nShells);
156 >    Points = ino->getPoints();
157 >  } else if (args_info.marks_given) {
158 >    int columnAtoms = args_info.columnAtoms_arg;    
159 >    int twinAtoms = args_info.twinAtoms_arg;
160 >    Decahedron* marks = new Decahedron(columnAtoms, nShells, twinAtoms);
161 >    Points = marks->getPoints();
162 >  }
163 >  else if (args_info.stone_given) {
164 >    int columnAtoms = args_info.columnAtoms_arg;    
165 >    int twinAtoms = args_info.twinAtoms_arg;
166 >    int truncatedPlanes = args_info.truncatedPlanes_arg;
167 >    CurlingStoneDecahedron* csd = new CurlingStoneDecahedron(columnAtoms, nShells, twinAtoms, truncatedPlanes);
168 >    Points = csd->getPoints();
169 >  }
170 >  
171 >    
172  
173    outputFileName = args_info.output_arg;
174    

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines