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 1876 by gezelter, Fri May 17 17:10:11 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;
55  
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 – }
191 –
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    
# Line 364 | Line 176 | int main(int argc, char *argv []) {
176  
177    createMdFile(inputFileName, outputFileName, Points.size());
178    
179 <  if (oldInfo != NULL)
368 <    delete oldInfo;
179 >  delete oldInfo;
180    
181    SimCreator newCreator;
182    SimInfo* NewInfo = newCreator.createSim(outputFileName, false);
# Line 383 | Line 194 | int main(int argc, char *argv []) {
194    Vector3d boxMax;
195    Vector3d boxMin;
196  
197 <  for (int n = 0; n < Points.size(); n++) {
197 >  for (unsigned int n = 0; n < Points.size(); n++) {
198      mol = NewInfo->getMoleculeByGlobalIndex(l);
199      Vector3d location = Points[n] * latticeConstant;
200      Vector3d orientation = Vector3d(0, 0, 1.0);
201 <
201 >    
202      if (n == 0) {
203        boxMax = location;
204        boxMin = location;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines