ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/applications/nanoRodBuilder/GeometryBuilder.cpp
Revision: 963
Committed: Wed May 17 21:51:42 2006 UTC (18 years, 11 months ago) by tim
File size: 27832 byte(s)
Log Message: 
Adding single precision capabilities to c++ side

File Contents

# User Rev Content
1 chuckv 467 /*
2     * GeometryBuilder.cpp
3     * nanorodBuilder
4     *
5     * Created by Charles Vardeman II on 4/4/05.
6     * Copyright 2005 University of Notre Dame. All rights reserved.
7     *
8     */
9    
10     /*
11     * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
12     *
13     * The University of Notre Dame grants you ("Licensee") a
14     * non-exclusive, royalty free, license to use, modify and
15     * redistribute this software in source and binary code form, provided
16     * that the following conditions are met:
17     *
18     * 1. Acknowledgement of the program authors must be made in any
19     * publication of scientific results based in part on use of the
20     * program. An acceptable form of acknowledgement is citation of
21     * the article in which the program was described (Matthew
22 gezelter 507 * A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
23     * J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
24 chuckv 467 * Parallel Simulation Engine for Molecular Dynamics,"
25 gezelter 507 * J. Comput. Chem. 26, pp. 252-271 (2005))
26 chuckv 467 *
27     * 2. Redistributions of source code must retain the above copyright
28     * notice, this list of conditions and the following disclaimer.
29     *
30     * 3. Redistributions in binary form must reproduce the above copyright
31     * notice, this list of conditions and the following disclaimer in the
32     * documentation and/or other materials provided with the
33     * distribution.
34     *
35     * This software is provided "AS IS," without a warranty of any
36     * kind. All express or implied conditions, representations and
37     * warranties, including any implied warranty of merchantability,
38     * fitness for a particular purpose or non-infringement, are hereby
39     * excluded. The University of Notre Dame and its licensors shall not
40     * be liable for any damages suffered by licensee as a result of
41     * using, modifying or distributing the software or its
42     * derivatives. In no event will the University of Notre Dame or its
43     * licensors be liable for any lost revenue, profit or data, or for
44     * direct, indirect, special, consequential, incidental or punitive
45     * damages, however caused and regardless of the theory of liability,
46     * arising out of the use of or inability to use software, even if the
47     * University of Notre Dame has been advised of the possibility of
48     * such damages.
49     */
50    
51     #include "GeometryBuilder.hpp"
52    
53     #include <CGAL/Simple_cartesian.h>
54     #include <CGAL/Polyhedron_incremental_builder_3.h>
55     #include <CGAL/Polyhedron_3.h>
56     #include <CGAL/Homogeneous.h>
57     #include <CGAL/Polyhedron_traits_with_normals_3.h>
58     #include <CGAL/Polyhedron_3.h>
59     #include <CGAL/Aff_transformation_3.h>
60     #include <CGAL/aff_transformation_tags.h>
61     #include <iostream>
62     #include <algorithm>
63 chuckv 518 #include <fstream>
64     #include <math.h>
65 chuckv 467
66 chuckv 518 using namespace std;
67 chuckv 467 using namespace oopse;
68    
69     //typedef CGAL::Homogeneous<int> Kernel;
70 tim 963 typedef CGAL::Simple_cartesian<RealType> Kernel;
71 chuckv 467 //typedef CGAL::Polyhedron_3<Kernel> Polyhedron;
72    
73     typedef Kernel::Point_3 Point_3;
74     typedef Kernel::Vector_3 Vector_3;
75     typedef CGAL::Polyhedron_traits_with_normals_3<Kernel> Traits;
76     //typedef CGAL::Polyhedron_3<Kernel> Polyhedron;
77     typedef CGAL::Polyhedron_3<Traits> Polyhedron;
78     typedef Polyhedron::HalfedgeDS HalfedgeDS;
79     typedef Polyhedron::Facet_iterator Facet_iterator;
80     typedef Polyhedron::Halfedge_around_facet_circulator Halfedge_facet_circulator;
81     typedef Polyhedron::Facet_iterator Facet_iterator;
82     typedef Polyhedron::Plane_iterator Plane_iterator;
83     typedef Polyhedron::Vertex_handle Vertex_handle;
84    
85     Polyhedron nanoRodPolyhedron;
86 chuckv 792 Polyhedron nanoRodTwinnedPolyhedron1;
87     Polyhedron nanoRodTwinnedPolyhedron2;
88     Polyhedron nanoRodTwinnedPolyhedron3;
89     Polyhedron nanoRodTwinnedPolyhedron4;
90     Polyhedron nanoRodTwinnedPolyhedron5;
91 chuckv 467
92    
93    
94    
95     //typedef CGAL::Scaling Scaling;
96     //typedef Aff_transformation_3<Kernel> A;( const Scaling,
97     // Kernel::RT s=RT(20),
98     // Kernel::RT hw = RT(1));
99    
100 chuckv 542
101    
102    
103    
104 chuckv 467 // A modifier creating a triangle with the incremental builder.
105     template <class HDS>
106 chuckv 521 class buildSingleCrystal : public CGAL::Modifier_base<HDS> {
107 gezelter 507 public:
108 chuckv 467 Vertex_handle end1;
109     Vertex_handle neight1;
110 gezelter 507 Vertex_handle end2;
111 chuckv 467 Vertex_handle neight2;
112 chuckv 518 Vertex_handle neight3;
113 chuckv 467
114 chuckv 521 buildSingleCrystal() {}
115 chuckv 467 void operator()( HDS& hds) {
116     // Postcondition: `hds' is a valid polyhedral surface.
117     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
118     B.begin_surface( 12, 15, 6);
119     typedef typename HDS::Vertex Vertex;
120     typedef typename Vertex::Point Point;
121    
122    
123    
124    
125    
126     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
127     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
128     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
129     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
130 chuckv 518 neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
131 chuckv 467 neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
132     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
133     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
134     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
135     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
136     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
137     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
138    
139     B.begin_facet();
140     B.add_vertex_to_facet( 7);
141     B.add_vertex_to_facet( 9);
142     B.add_vertex_to_facet( 11);
143     B.add_vertex_to_facet( 10);
144     B.end_facet();
145    
146     B.begin_facet();
147     B.add_vertex_to_facet( 8);
148     B.add_vertex_to_facet( 10);
149     B.add_vertex_to_facet( 11);
150     B.end_facet();
151    
152     B.begin_facet();
153     B.add_vertex_to_facet( 3);
154     B.add_vertex_to_facet( 9);
155     B.add_vertex_to_facet( 7);
156     B.end_facet();
157    
158     B.begin_facet();
159     B.add_vertex_to_facet( 9);
160     B.add_vertex_to_facet( 5);
161     B.add_vertex_to_facet( 6);
162     B.add_vertex_to_facet( 11);
163     B.end_facet();
164    
165     B.begin_facet();
166     B.add_vertex_to_facet( 8);
167     B.add_vertex_to_facet( 11);
168     B.add_vertex_to_facet( 6);
169     B.end_facet();
170    
171     B.begin_facet();
172     B.add_vertex_to_facet( 3);
173     B.add_vertex_to_facet( 5);
174     B.add_vertex_to_facet( 9);
175     B.end_facet();
176    
177     B.begin_facet();
178     B.add_vertex_to_facet( 5);
179     B.add_vertex_to_facet( 0);
180     B.add_vertex_to_facet( 2);
181     B.add_vertex_to_facet( 6);
182     B.end_facet();
183    
184     B.begin_facet();
185     B.add_vertex_to_facet( 8);
186     B.add_vertex_to_facet( 6);
187     B.add_vertex_to_facet( 2);
188     B.end_facet();
189    
190     B.begin_facet();
191     B.add_vertex_to_facet( 3);
192     B.add_vertex_to_facet( 0);
193     B.add_vertex_to_facet( 5);
194     B.end_facet();
195    
196     B.begin_facet();
197     B.add_vertex_to_facet( 0);
198     B.add_vertex_to_facet( 1);
199     B.add_vertex_to_facet( 4);
200     B.add_vertex_to_facet( 2);
201     B.end_facet();
202    
203     B.begin_facet();
204     B.add_vertex_to_facet( 8);
205     B.add_vertex_to_facet( 2);
206     B.add_vertex_to_facet( 4);
207     B.end_facet();
208    
209     B.begin_facet();
210     B.add_vertex_to_facet( 3);
211     B.add_vertex_to_facet( 1);
212     B.add_vertex_to_facet( 0);
213     B.end_facet();
214    
215     B.begin_facet();
216     B.add_vertex_to_facet( 1);
217     B.add_vertex_to_facet( 7);
218     B.add_vertex_to_facet( 10);
219     B.add_vertex_to_facet( 4);
220     B.end_facet();
221    
222     B.begin_facet();
223     B.add_vertex_to_facet( 8);
224     B.add_vertex_to_facet( 4);
225     B.add_vertex_to_facet( 10);
226     B.end_facet();
227    
228     B.begin_facet();
229     B.add_vertex_to_facet( 3);
230     B.add_vertex_to_facet( 7);
231     B.add_vertex_to_facet( 1);
232     B.end_facet();
233    
234    
235     B.end_surface();
236     }
237     };
238    
239    
240    
241 chuckv 792 /* First Crystal*/
242    
243     template <class HDS>
244     class buildtwinned1 : public CGAL::Modifier_base<HDS> {
245     public:
246     Vertex_handle end1;
247     Vertex_handle neight1;
248     Vertex_handle end2;
249     Vertex_handle neight2;
250     Vertex_handle neight3;
251    
252     buildtwinned1() {}
253     void operator()( HDS& hds) {
254     // Postcondition: `hds' is a valid polyhedral surface.
255     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
256     B.begin_surface( 12, 5, 6);
257     typedef typename HDS::Vertex Vertex;
258     typedef typename Vertex::Point Point;
259    
260    
261    
262    
263    
264     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
265     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
266     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
267     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
268     neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
269     neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
270     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
271     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
272     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
273     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
274     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
275     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
276    
277     B.begin_facet();
278     B.add_vertex_to_facet( 8);
279     B.add_vertex_to_facet( 4);
280     B.add_vertex_to_facet( 1);
281     B.add_vertex_to_facet( 3);
282     B.end_facet();
283    
284     B.begin_facet();
285     B.add_vertex_to_facet( 3);
286     B.add_vertex_to_facet( 0);
287     B.add_vertex_to_facet( 2);
288     B.add_vertex_to_facet( 8);
289     B.end_facet();
290    
291     B.begin_facet();
292     B.add_vertex_to_facet( 0);
293     B.add_vertex_to_facet( 1);
294     B.add_vertex_to_facet( 4);
295     B.add_vertex_to_facet( 2);
296     B.end_facet();
297    
298    
299     B.begin_facet();
300     B.add_vertex_to_facet( 3);
301     B.add_vertex_to_facet( 1);
302     B.add_vertex_to_facet( 0);
303     B.end_facet();
304    
305     B.begin_facet();
306     B.add_vertex_to_facet( 4);
307     B.add_vertex_to_facet( 8);
308     B.add_vertex_to_facet( 2);
309     B.end_facet();
310    
311     B.end_surface();
312     }
313     };
314    
315    
316    
317    
318     /*second crystal*/
319     template <class HDS>
320     class buildtwinned2 : public CGAL::Modifier_base<HDS> {
321     public:
322     Vertex_handle end1;
323     Vertex_handle neight1;
324     Vertex_handle end2;
325     Vertex_handle neight2;
326     Vertex_handle neight3;
327    
328     buildtwinned2() {}
329     void operator()( HDS& hds) {
330     // Postcondition: `hds' is a valid polyhedral surface.
331     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
332     B.begin_surface( 12, 5, 6);
333     typedef typename HDS::Vertex Vertex;
334     typedef typename Vertex::Point Point;
335    
336    
337    
338    
339    
340     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
341     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
342     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
343     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
344     neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
345     neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
346     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
347     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
348     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
349     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
350     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
351     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
352    
353     B.begin_facet();
354     B.add_vertex_to_facet( 8);
355     B.add_vertex_to_facet( 2);
356     B.add_vertex_to_facet( 0);
357     B.add_vertex_to_facet( 3);
358     B.end_facet();
359    
360     B.begin_facet();
361     B.add_vertex_to_facet( 3);
362     B.add_vertex_to_facet( 5);
363     B.add_vertex_to_facet( 6);
364     B.add_vertex_to_facet( 8);
365     B.end_facet();
366    
367     B.begin_facet();
368     B.add_vertex_to_facet( 5);
369     B.add_vertex_to_facet( 0);
370     B.add_vertex_to_facet( 2);
371     B.add_vertex_to_facet( 6);
372     B.end_facet();
373    
374    
375     B.begin_facet();
376     B.add_vertex_to_facet( 3);
377     B.add_vertex_to_facet( 0);
378     B.add_vertex_to_facet( 5);
379     B.end_facet();
380    
381     B.begin_facet();
382     B.add_vertex_to_facet( 8);
383     B.add_vertex_to_facet( 6);
384     B.add_vertex_to_facet( 2);
385     B.end_facet();
386    
387     B.end_surface();
388     }
389     };
390    
391     /* Third Crystal*/
392    
393     template <class HDS>
394     class buildtwinned3 : public CGAL::Modifier_base<HDS> {
395     public:
396     Vertex_handle end1;
397     Vertex_handle neight1;
398     Vertex_handle end2;
399     Vertex_handle neight2;
400     Vertex_handle neight3;
401    
402     buildtwinned3() {}
403     void operator()( HDS& hds) {
404     // Postcondition: `hds' is a valid polyhedral surface.
405     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
406     B.begin_surface( 12, 5, 6);
407     typedef typename HDS::Vertex Vertex;
408     typedef typename Vertex::Point Point;
409    
410    
411    
412    
413    
414     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
415     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
416     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
417     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
418     neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
419     neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
420     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
421     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
422     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
423     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
424     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
425     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
426    
427     B.begin_facet();
428     B.add_vertex_to_facet( 8);
429     B.add_vertex_to_facet( 6);
430     B.add_vertex_to_facet( 5);
431     B.add_vertex_to_facet( 3);
432     B.end_facet();
433    
434     B.begin_facet();
435     B.add_vertex_to_facet( 3);
436     B.add_vertex_to_facet( 9);
437     B.add_vertex_to_facet( 11);
438     B.add_vertex_to_facet( 8);
439     B.end_facet();
440    
441     B.begin_facet();
442     B.add_vertex_to_facet( 9);
443     B.add_vertex_to_facet( 5);
444     B.add_vertex_to_facet( 6);
445     B.add_vertex_to_facet( 11);
446     B.end_facet();
447    
448    
449     B.begin_facet();
450     B.add_vertex_to_facet( 3);
451     B.add_vertex_to_facet( 5);
452     B.add_vertex_to_facet( 9);
453     B.end_facet();
454    
455     B.begin_facet();
456     B.add_vertex_to_facet( 8);
457     B.add_vertex_to_facet( 11);
458     B.add_vertex_to_facet( 6);
459     B.end_facet();
460    
461     B.end_surface();
462     }
463     };
464    
465     /*Fourth Segment*/
466     template <class HDS>
467     class buildtwinned4 : public CGAL::Modifier_base<HDS> {
468     public:
469     Vertex_handle end1;
470     Vertex_handle neight1;
471     Vertex_handle end2;
472     Vertex_handle neight2;
473     Vertex_handle neight3;
474    
475     buildtwinned4() {}
476     void operator()( HDS& hds) {
477     // Postcondition: `hds' is a valid polyhedral surface.
478     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
479     B.begin_surface( 12, 5, 6);
480     typedef typename HDS::Vertex Vertex;
481     typedef typename Vertex::Point Point;
482    
483    
484    
485    
486    
487     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
488     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
489     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
490     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
491     neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
492     neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
493     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
494     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
495     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
496     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
497     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
498     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
499    
500     B.begin_facet();
501     B.add_vertex_to_facet( 8);
502     B.add_vertex_to_facet( 11);
503     B.add_vertex_to_facet( 9);
504     B.add_vertex_to_facet( 3);
505     B.end_facet();
506    
507     B.begin_facet();
508     B.add_vertex_to_facet( 3);
509     B.add_vertex_to_facet( 7);
510     B.add_vertex_to_facet( 10);
511     B.add_vertex_to_facet( 8);
512     B.end_facet();
513    
514     B.begin_facet();
515     B.add_vertex_to_facet( 7);
516     B.add_vertex_to_facet( 9);
517     B.add_vertex_to_facet( 11);
518     B.add_vertex_to_facet( 10);
519     B.end_facet();
520    
521    
522     B.begin_facet();
523     B.add_vertex_to_facet( 3);
524     B.add_vertex_to_facet( 9);
525     B.add_vertex_to_facet( 7);
526     B.end_facet();
527    
528     B.begin_facet();
529     B.add_vertex_to_facet( 8);
530     B.add_vertex_to_facet( 10);
531     B.add_vertex_to_facet( 11);
532     B.end_facet();
533    
534     B.end_surface();
535     }
536     };
537    
538    
539     /*Fifth Segment*/
540     template <class HDS>
541     class buildtwinned5 : public CGAL::Modifier_base<HDS> {
542     public:
543     Vertex_handle end1;
544     Vertex_handle neight1;
545     Vertex_handle end2;
546     Vertex_handle neight2;
547     Vertex_handle neight3;
548    
549     buildtwinned5() {}
550     void operator()( HDS& hds) {
551     // Postcondition: `hds' is a valid polyhedral surface.
552     CGAL::Polyhedron_incremental_builder_3<HDS> B( hds, true);
553     B.begin_surface( 12, 5, 6);
554     typedef typename HDS::Vertex Vertex;
555     typedef typename Vertex::Point Point;
556    
557    
558    
559    
560    
561     B.add_vertex( Point(-0.7887222926324, 0.4874571845315, -0.2562714077342));
562     B.add_vertex( Point(-0.4874571845316, 0.4874571845315, 0.6709272557930));
563     B.add_vertex( Point(-0.7887222926324, -0.4874571845315, -0.2562714077342)); //End vertex
564     end1 = B.add_vertex( Point( 0.0000000000000, 1.0000000000000, 0.0000000000000));
565     neight3 = B.add_vertex( Point(-0.4874571845315, -0.4874571845316, 0.6709272557930));
566     neight1 = B.add_vertex( Point(-0.0000000000000, 0.4874571845316, -0.8293116961175));
567     B.add_vertex( Point( 0.0000000000000, -0.4874571845316, -0.8293116961175));
568     B.add_vertex( Point( 0.4874571845315, 0.4874571845316, 0.6709272557930));
569     end2 = B.add_vertex( Point(-0.0000000000000, -1.0000000000000, 0.0000000000000)); //End Vertex
570     B.add_vertex( Point( 0.7887222926324, 0.4874571845315, -0.2562714077342));
571     neight2 = B.add_vertex( Point( 0.4874571845316, -0.4874571845315, 0.6709272557930));
572     B.add_vertex( Point( 0.7887222926324, -0.4874571845315, -0.2562714077342));
573    
574     B.begin_facet();
575     B.add_vertex_to_facet( 8);
576     B.add_vertex_to_facet( 4);
577     B.add_vertex_to_facet( 1);
578     B.add_vertex_to_facet( 3);
579     B.end_facet();
580    
581     B.begin_facet();
582     B.add_vertex_to_facet( 3);
583     B.add_vertex_to_facet( 7);
584     B.add_vertex_to_facet( 10);
585     B.add_vertex_to_facet( 8);
586     B.end_facet();
587    
588     B.begin_facet();
589     B.add_vertex_to_facet( 7);
590     B.add_vertex_to_facet( 1);
591     B.add_vertex_to_facet( 4);
592     B.add_vertex_to_facet( 10);
593     B.end_facet();
594    
595    
596     B.begin_facet();
597     B.add_vertex_to_facet( 1);
598     B.add_vertex_to_facet( 7);
599     B.add_vertex_to_facet( 3);
600     B.end_facet();
601    
602     B.begin_facet();
603     B.add_vertex_to_facet( 10);
604     B.add_vertex_to_facet( 4);
605     B.add_vertex_to_facet( 8);
606     B.end_facet();
607    
608     B.end_surface();
609     }
610     };
611    
612    
613    
614    
615    
616    
617 chuckv 467 struct Normal_vector {
618     template <class Facet>
619     typename Facet::Plane_3 operator()( Facet& f) {
620     typename Facet::Halfedge_handle h = f.halfedge();
621     // Facet::Plane_3 is the normal vector type. We assume the
622     // CGAL Kernel here and use its global functions.
623     return CGAL::cross_product(
624     h->next()->vertex()->point() - h->vertex()->point(),
625     h->next()->next()->vertex()->point() - h->next()->vertex()->point());
626     }
627     };
628    
629    
630 tim 963 bool GeometryBuilder::isInsidePolyhedron(RealType x, RealType y, RealType z) {
631 chuckv 467
632 gezelter 507 Point_3 point(x,y,z);
633 chuckv 467 Plane_iterator i;
634     Facet_iterator j;
635     for ( i =nanoRodPolyhedron.planes_begin(), j = nanoRodPolyhedron.facets_begin(); i != nanoRodPolyhedron.planes_end() && j !=nanoRodPolyhedron.facets_end(); ++i, ++j) {
636     Halfedge_facet_circulator k = j->facet_begin();
637    
638     Vector_3 newVector = point - k->vertex()->point();
639     Vector_3 normal = *i;
640 tim 963 RealType dot_product = newVector.x() * normal.x() + newVector.y() * normal.y() + newVector.z() * normal.z();
641 chuckv 467
642     if (dot_product < 0) {
643     return false;
644     }
645     }
646    
647     return true;
648     }
649    
650    
651 tim 963 GeometryBuilder::GeometryBuilder(RealType length,RealType width) {
652 gezelter 507 // Create the geometry for nanorod
653 chuckv 792 buildSingleCrystal<HalfedgeDS> singleCrystalNanorod;
654 chuckv 467
655 chuckv 792 nanoRodPolyhedron.delegate( singleCrystalNanorod);
656 chuckv 467
657 tim 963 RealType y1 = singleCrystalNanorod.end1->point().y() - singleCrystalNanorod.neight1->point().y();
658     RealType y2 = singleCrystalNanorod.end2->point().y() - singleCrystalNanorod.neight2->point().y();
659 chuckv 518
660 tim 963 RealType endDist = sqrt(pow(singleCrystalNanorod.neight2->point().x() - singleCrystalNanorod.neight3->point().x(),2)+
661 chuckv 792 pow(singleCrystalNanorod.neight2->point().y() - singleCrystalNanorod.neight3->point().y(),2)+
662     pow(singleCrystalNanorod.neight2->point().z() - singleCrystalNanorod.neight3->point().z(),2));
663 chuckv 518
664 tim 963 RealType endRatio1 = y1/endDist;
665     RealType endRatio2 = y2/endDist;
666 chuckv 518
667     std::cout << "End dist is " << endDist <<" ratio " << endRatio1 << std::endl;
668    
669 gezelter 507 CGAL::Aff_transformation_3<Kernel> aff_tranformation( width,
670     0.0,
671     0.0,
672     0.0,
673     0.0,
674     length,
675     0.0,
676     0.0,
677     0.0,
678     0.0,
679     width,
680     0.0);
681     std::transform( nanoRodPolyhedron.points_begin(), nanoRodPolyhedron.points_end(), nanoRodPolyhedron.points_begin(), aff_tranformation);
682 chuckv 467
683 chuckv 518
684 tim 963 RealType endDist2 = sqrt(pow(singleCrystalNanorod.neight2->point().x() -singleCrystalNanorod.neight3->point().x(),2)+
685 chuckv 792 pow(singleCrystalNanorod.neight2->point().y() - singleCrystalNanorod.neight3->point().y(),2)+
686     pow(singleCrystalNanorod.neight2->point().z() - singleCrystalNanorod.neight3->point().z(),2));
687 chuckv 518
688 chuckv 792 Point_3 point1(singleCrystalNanorod.end1->point().x(), endDist2*endRatio1 + singleCrystalNanorod.neight1->point().y(), singleCrystalNanorod.end1->point().z());
689     Point_3 point2(singleCrystalNanorod.end2->point().x(), endDist2*endRatio2 + singleCrystalNanorod.neight2->point().y(), singleCrystalNanorod.end2->point().z());
690     singleCrystalNanorod.end1->point() = point1;
691     singleCrystalNanorod.end2->point() = point2;
692 chuckv 467
693 gezelter 507 // Construct normal vector for each face.
694     std::transform( nanoRodPolyhedron.facets_begin(), nanoRodPolyhedron.facets_end(), nanoRodPolyhedron.planes_begin(),
695     Normal_vector());
696 chuckv 518 }
697 chuckv 792
698    
699    
700 tim 963 GeometryBuilder::GeometryBuilder(RealType length,RealType width, bool twinned) {
701 chuckv 792 // Create the geometry for nanorod
702    
703     buildtwinned1<HalfedgeDS> crystal1;
704     buildtwinned2<HalfedgeDS> crystal2;
705     buildtwinned3<HalfedgeDS> crystal3;
706     buildtwinned4<HalfedgeDS> crystal4;
707     buildtwinned5<HalfedgeDS> crystal5;
708    
709    
710     nanoRodTwinnedPolyhedron1.delegate( crystal1);
711     nanoRodTwinnedPolyhedron2.delegate( crystal2);
712     nanoRodTwinnedPolyhedron3.delegate( crystal3);
713     nanoRodTwinnedPolyhedron4.delegate( crystal4);
714     nanoRodTwinnedPolyhedron5.delegate( crystal5);
715    
716    
717    
718    
719    
720    
721 tim 963 RealType y1 = crystal1.end1->point().y() - crystal1.neight1->point().y();
722     RealType y2 = crystal1.end2->point().y() - crystal1.neight2->point().y();
723 chuckv 792
724 tim 963 RealType endDist = sqrt(pow(crystal1.neight2->point().x() - crystal1.neight3->point().x(),2)+
725 chuckv 792 pow(crystal1.neight2->point().y() - crystal1.neight3->point().y(),2)+
726     pow(crystal1.neight2->point().z() - crystal1.neight3->point().z(),2));
727    
728 tim 963 RealType endRatio1 = y1/endDist;
729     RealType endRatio2 = y2/endDist;
730 chuckv 792
731     std::cout << "End dist is " << endDist <<" ratio " << endRatio1 << std::endl;
732    
733     CGAL::Aff_transformation_3<Kernel> aff_tranformation( width,
734     0.0,
735     0.0,
736     0.0,
737     0.0,
738     length,
739     0.0,
740     0.0,
741     0.0,
742     0.0,
743     width,
744     0.0);
745     std::transform(nanoRodTwinnedPolyhedron1.points_begin(), nanoRodTwinnedPolyhedron1.points_end(), nanoRodTwinnedPolyhedron1.points_begin(), aff_tranformation);
746    
747    
748 tim 963 RealType endDist2 = sqrt(pow(crystal1.neight2->point().x() - crystal1.neight3->point().x(),2)+
749 chuckv 792 pow(crystal1.neight2->point().y() - crystal1.neight3->point().y(),2)+
750     pow(crystal1.neight2->point().z() - crystal1.neight3->point().z(),2));
751    
752     Point_3 point1(crystal1.end1->point().x(), endDist2*endRatio1 + crystal1.neight1->point().y(), crystal1.end1->point().z());
753     Point_3 point2(crystal1.end2->point().x(), endDist2*endRatio2 + crystal1.neight2->point().y(), crystal1.end2->point().z());
754     crystal1.end1->point() = point1;
755     crystal1.end2->point() = point2;
756    
757     // Construct normal vector for each face.
758     std::transform( nanoRodTwinnedPolyhedron1.facets_begin(), nanoRodTwinnedPolyhedron1.facets_end(), nanoRodTwinnedPolyhedron1.planes_begin(),
759     Normal_vector());
760     }
761    
762    
763    
764    
765    
766    
767    
768    
769    
770 chuckv 518 void GeometryBuilder::dumpGeometry(const std::string& geomFileName){
771    
772     std::ofstream newGeomFile;
773    
774     //create new .md file based on old .md file
775     newGeomFile.open(geomFileName.c_str());
776    
777     // Write polyhedron in Object File Format (OFF).
778     CGAL::set_ascii_mode( std::cout);
779     newGeomFile << "OFF" << std::endl << nanoRodPolyhedron.size_of_vertices() << ' '
780     << nanoRodPolyhedron.size_of_facets() << " 0" << std::endl;
781     std::copy( nanoRodPolyhedron.points_begin(), nanoRodPolyhedron.points_end(),
782     std::ostream_iterator<Point_3>( newGeomFile, "\n"));
783     for ( Facet_iterator i = nanoRodPolyhedron.facets_begin(); i != nanoRodPolyhedron.facets_end(); ++i) {
784     Halfedge_facet_circulator j = i->facet_begin();
785     // Facets in polyhedral surfaces are at least triangles.
786     CGAL_assertion( CGAL::circulator_size(j) >= 3);
787     newGeomFile << CGAL::circulator_size(j) << ' ';
788     do {
789     newGeomFile << ' ' << std::distance(nanoRodPolyhedron.vertices_begin(), j->vertex());
790     } while ( ++j != i->facet_begin());
791     newGeomFile << std::endl;
792     }
793    
794     newGeomFile.close();
795    
796    
797 chuckv 467
798    
799 chuckv 518 }
800 chuckv 467
801 chuckv 518