Develop and Download Open Source Software

Powder auto-indexing software Conograph

Browse Subversion Repository

/[conograph]/Conograph/trunk/src/lattice_symmetry/LatticeFigureOfMeritToCheckSymmetry.cc

Contents of /Conograph/trunk/src/lattice_symmetry/LatticeFigureOfMeritToCheckSymmetry.cc

Parent Directory Parent Directory | Revision Log Revision Log

Revision 25 - (show annotations) (download) (as text)
Mon Jul 7 02:35:51 2014 UTC (9 years, 8 months ago) by rtomiyasu
File MIME type: text/x-c++src
File size: 15205 byte(s) 
Source codes of version 0.9.99
1 /*
2 * The MIT License
3
4 BLDConograph (Bravais lattice determination module in Conograph)
5
6 Copyright (c) <2012> <Ryoko Oishi-Tomiyasu, KEK>
7
8 Permission is hereby granted, free of charge, to any person obtaining a copy
9 of this software and associated documentation files (the "Software"), to deal
10 in the Software without restriction, including without limitation the rights
11 to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
12 copies of the Software, and to permit persons to whom the Software is
13 furnished to do so, subject to the following conditions:
14
15 The above copyright notice and this permission notice shall be included in
16 all copies or substantial portions of the Software.
17
18 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
19 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
20 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
21 AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
22 LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
23 OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
24 THE SOFTWARE.
25 *
26 */
27 #include <limits>
28 #include "../utility_lattice_reduction/put_Buerger_reduced_lattice.hh"
29 #include "check_equiv.hh"
30 #include "LatticeFigureOfMeritToCheckSymmetry.hh"
31
32 LatticeFigureOfMeritToCheckSymmetry::LatticeFigureOfMeritToCheckSymmetry()
33 : m_S_red( SymMat43_Double( SymMat<Double>(3), NRMat<Int4>(4,3) ) )
34 {
35 }
36
37
38 LatticeFigureOfMeritToCheckSymmetry::LatticeFigureOfMeritToCheckSymmetry(const Double& rhs)
39 : m_latfom(rhs),
40 m_S_red( SymMat43_Double( SymMat<Double>(3), NRMat<Int4>(4,3) ) )
41 {
42 }
43
44
45 LatticeFigureOfMeritToCheckSymmetry::LatticeFigureOfMeritToCheckSymmetry(const BravaisType& brat,
46 const SymMat43_Double& S)
47 : m_S_red( SymMat43_Double( SymMat<Double>(3), NRMat<Int4>(4,3) ) )
48 {
49 setLatticeConstants43(brat, S);
50 }
51
52
53
54
55 #ifdef DEBUG
56 static bool checkInitialLatticeParameters(
57 const BravaisType& brat,
58 const SymMat43_Double& S_red)
59 {
60 const SymMat<Double> dbl_S_red( S_red.first );
61
62 if( brat.enumLaueGroup() == Ci && brat.enumCentringType() == Prim )
63 {
64 assert( dbl_S_red(2,2)*0.9999 < dbl_S_red(1,1) && dbl_S_red(1,1)*0.9999 < dbl_S_red(0,0)
65 && fabs( dbl_S_red(0,1) ) * 1.9999 < dbl_S_red(1,1)
66 && fabs( dbl_S_red(0,2) ) * 1.9999 < dbl_S_red(2,2)
67 && fabs( dbl_S_red(1,2) ) * 1.9999 < dbl_S_red(2,2) );
68 }
69 else if( brat.enumLaueGroup() == C2h_Y && brat.enumCentringType() == Prim )
70 {
71 assert( 0.0 <= dbl_S_red(0,2) && dbl_S_red(2,2)*0.9999 < dbl_S_red(0,0)
72 && fabs( dbl_S_red(0,2) ) * 1.9999 < dbl_S_red(2,2) && fabs( dbl_S_red(0,2) ) * 1.9999 < dbl_S_red(0,0) );
73 }
74 else if( brat.enumLaueGroup() == C2h_Z && brat.enumCentringType() == Prim )
75 {
76 assert( 0.0 <= dbl_S_red(0,1) && dbl_S_red(1,1)*0.9999 < dbl_S_red(0,0)
77 && fabs( dbl_S_red(0,1) ) * 1.9999 < dbl_S_red(0,0) && fabs( dbl_S_red(0,1) ) * 1.9999 < dbl_S_red(1,1) );
78 }
79 else if( brat.enumLaueGroup() == C2h_X && brat.enumCentringType() == Prim )
80 {
81 assert( 0.0 <= dbl_S_red(1,2) && dbl_S_red(2,2)*0.9999 < dbl_S_red(1,1)
82 && fabs( dbl_S_red(1,2) ) * 1.9999 < dbl_S_red(1,1) && fabs( dbl_S_red(1,2) ) * 1.9999 < dbl_S_red(2,2) );
83 }
84 else if( brat.enumLaueGroup() == C2h_Y && brat.enumCentringType() == BaseZ )
85 {
86 assert( 0.0 <= dbl_S_red(0,2) && fabs( dbl_S_red(0,2) ) * 1.9999 < dbl_S_red(2,2) && fabs( dbl_S_red(0,2) ) * 0.9999 < dbl_S_red(0,0) );
87 }
88 else if( brat.enumLaueGroup() == C2h_Z && brat.enumCentringType() == BaseX )
89 {
90 assert( 0.0 <= dbl_S_red(0,1) && fabs( dbl_S_red(0,1) ) * 1.9999 < dbl_S_red(0,0) && fabs( dbl_S_red(0,1) ) * 0.9999 < dbl_S_red(1,1) );
91 }
92 else if( brat.enumLaueGroup() == C2h_X && brat.enumCentringType() == BaseY )
93 {
94 assert( 0.0 <= dbl_S_red(1,2) && fabs( dbl_S_red(1,2) ) * 1.9999 < dbl_S_red(1,1) && fabs( dbl_S_red(1,2) ) * 0.9999 < dbl_S_red(2,2) );
95 }
96 else if( brat.enumLaueGroup() == D2h
97 && brat.enumCentringType() != BaseX
98 && brat.enumCentringType() != BaseY
99 && brat.enumCentringType() != BaseZ )
100 {
101 assert( dbl_S_red(2,2)*0.9999 < dbl_S_red(1,1) && dbl_S_red(1,1)*0.9999 < dbl_S_red(0,0) );
102 }
103
104 const SymMat<Double> S_super = transform_sym_matrix(S_red.second, S_red.first);
105 assert( S_super(0,1) <= 0.0
106 && S_super(0,2) <= 0.0
107 && S_super(0,3) <= 0.0
108 && S_super(1,2) <= 0.0
109 && S_super(1,3) <= 0.0
110 && S_super(2,3) <= 0.0 );
111
112 return true;
113 }
114 #endif
115
116
117 void LatticeFigureOfMeritToCheckSymmetry::setLatticeConstants43(const BravaisType& brat, const SymMat43_Double& S)
118 {
119 m_S_red = S;
120 assert( checkInitialLatticeParameters(brat, m_S_red) );
121
122 m_latfom.setLatticeConstants43(brat, S);
123 }
124
125
126
127 static bool operator<(const SymMat<Double>& lhs, const SymMat<Double>& rhs)
128 {
129 static const Double EPS_1 = 1.0+sqrt( numeric_limits<double>::epsilon() );
130 assert( lhs.size() == 3 );
131 assert( rhs.size() == 3 );
132
133 static const Int4 ISIZE = 3;
134 for(Int4 i=0; i<ISIZE; i++)
135 {
136 if( lhs(i,i)*EPS_1 < rhs(i,i) ) return true;
137 if( rhs(i,i)*EPS_1 < lhs(i,i) ) return false;
138
139 for(Int4 j=0; j<i; j++)
140 {
141 const Double lhs_ij = lhs(i,i)+lhs(j,j)+lhs(i,j)*2.0;
142 const Double rhs_ij = rhs(i,i)+rhs(j,j)+rhs(i,j)*2.0;
143
144 if( lhs_ij*EPS_1 < rhs_ij ) return true;
145 if( rhs_ij*EPS_1 < lhs_ij ) return false;
146 }
147 }
148
149 return false;
150 }
151
152
153
154 bool LatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsMonoclinic(const ePointGroup& epg_new,
155 const Double& resol,
156 map< SymMat<Double>, NRMat<Int4> >& ans) const
157 {
158 ans.clear();
159
160 SymMat<Double> ans0 = m_S_red.first;
161 cal_average_crystal_system(C2h_X, ans0);
162
163 SymMat<Double> S_red(3);
164 NRMat<Int4> trans_mat2;
165 if( check_equiv_m(ans0, m_S_red.first, resol ) )
166 {
167 if( epg_new == C2h_X )
168 {
169 S_red = ans0;
170 trans_mat2 = m_S_red.second;
171 putBuergerReducedMonoclinicP(1, 2, S_red, trans_mat2);
172 }
173 else if( epg_new == C2h_Y )
174 {
175 S_red = transform_sym_matrix(put_matrix_YXZ(), ans0);
176 trans_mat2 = mprod(m_S_red.second, put_matrix_YXZ());
177 putBuergerReducedMonoclinicP(0, 2, S_red, trans_mat2);
178 }
179 else // if( epg_new == C2h_Z )
180 {
181 S_red = transform_sym_matrix(put_matrix_YZX(), ans0);
182 trans_mat2 = mprod(m_S_red.second, put_matrix_ZXY());
183 putBuergerReducedMonoclinicP(0, 1, S_red, trans_mat2);
184 }
185 ans.insert( SymMat43_Double( S_red, trans_mat2) );
186 }
187
188 ans0 = m_S_red.first;
189 cal_average_crystal_system(C2h_Y, ans0);
190 if( check_equiv_m(ans0, m_S_red.first, resol ) )
191 {
192 if( epg_new == C2h_X )
193 {
194 S_red = transform_sym_matrix(put_matrix_YXZ(), ans0);
195 trans_mat2 = mprod(m_S_red.second, put_matrix_YXZ());
196 putBuergerReducedMonoclinicP(1, 2, S_red, trans_mat2);
197 }
198 else if( epg_new == C2h_Y )
199 {
200 S_red = ans0;
201 trans_mat2 = m_S_red.second;
202 putBuergerReducedMonoclinicP(0, 2, S_red, trans_mat2);
203 }
204 else // if( epg_new == C2h_Z )
205 {
206 S_red = transform_sym_matrix(put_matrix_XZY(), ans0);
207 trans_mat2 = mprod(m_S_red.second, put_matrix_XZY());
208 putBuergerReducedMonoclinicP(0, 1, S_red, trans_mat2);
209 }
210 ans.insert( SymMat43_Double( S_red, trans_mat2) );
211 }
212
213 ans0 = m_S_red.first;
214 cal_average_crystal_system(C2h_Z, ans0);
215 if( check_equiv_m(ans0, m_S_red.first, resol ) )
216 {
217 if( epg_new == C2h_X )
218 {
219 S_red = transform_sym_matrix(put_matrix_ZXY(), ans0);
220 trans_mat2 = mprod(m_S_red.second, put_matrix_YZX());
221 putBuergerReducedMonoclinicP(1, 2, S_red, trans_mat2);
222 }
223 else if( epg_new == C2h_Y )
224 {
225 S_red = transform_sym_matrix(put_matrix_XZY(), ans0);
226 trans_mat2 = mprod(m_S_red.second, put_matrix_XZY());
227 putBuergerReducedMonoclinicP(0, 2, S_red, trans_mat2);
228 }
229 else // if( epg_new == C2h_Z )
230 {
231 S_red = ans0;
232 trans_mat2 = m_S_red.second;
233 putBuergerReducedMonoclinicP(0, 1, S_red, trans_mat2);
234 }
235 ans.insert( SymMat43_Double( S_red, trans_mat2) );
236 }
237
238 return !( ans.empty() );
239 }
240
241
242 bool LatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsOrthorhombic(const Double& resol,
243 map< SymMat<Double>, NRMat<Int4> >& ans) const
244 {
245 ans.clear();
246
247 const BravaisType& brat = m_latfom.putBravaisType();
248
249 SymMat<Double> ans0 = m_S_red.first;
250 cal_average_crystal_system(D2h, ans0);
251 if( check_equiv_m(ans0, m_S_red.first, resol ) )
252 {
253 if( brat.enumCentringType() == BaseX )
254 {
255 if( ans0(1,1) < ans0(2,2) )
256 {
257 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_ZYX(), ans0), mprod( m_S_red.second, put_matrix_ZYX() ) ) );
258 }
259 else
260 {
261 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_YZX(), ans0), mprod( m_S_red.second, put_matrix_ZXY() ) ) );
262 }
263 }
264 else if( brat.enumCentringType() == BaseY )
265 {
266 if( ans0(0,0) < ans0(2,2) )
267 {
268 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_ZXY(), ans0), mprod( m_S_red.second, put_matrix_YZX() ) ) );
269 }
270 else
271 {
272 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_XZY(), ans0), mprod( m_S_red.second, put_matrix_XZY() ) ) );
273 }
274 }
275 else if( brat.enumCentringType() == BaseZ )
276 {
277 if( ans0(0,0) < ans0(1,1) )
278 {
279 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_YXZ(), ans0), mprod( m_S_red.second, put_matrix_YXZ() ) ) );
280 }
281 else
282 {
283 ans.insert( SymMat43_Double( transform_sym_matrix(put_matrix_XYZ(), ans0), m_S_red.second ) );
284 }
285 }
286 else
287 {
288 NRMat<Int4> trans_mat = m_S_red.second;
289 putBuergerReducedOrthorhombic(brat.enumCentringType(), ans0, trans_mat);
290 ans.insert( SymMat43_Double(ans0, trans_mat ) );
291 }
292 return true;
293 }
294 return false;
295 }
296
297
298 bool LatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsTetragonal(const Double& resol,
299 map< SymMat<Double>, NRMat<Int4> >& ans) const
300 {
301 ans.clear();
302
303 SymMat<Double> ans0 = m_S_red.first;
304 cal_average_crystal_system(D4h_X, ans0);
305 if( check_equiv_m(ans0, m_S_red.first, resol ) )
306 {
307 ans.insert( SymMat43_Double(
308 transform_sym_matrix(put_matrix_YZX(), ans0), mprod( m_S_red.second, put_matrix_ZXY() ) ) );
309 }
310
311 ans0 = m_S_red.first;
312 cal_average_crystal_system(D4h_Y, ans0);
313 if( check_equiv_m(ans0, m_S_red.first, resol ) )
314 {
315 ans.insert( SymMat43_Double(
316 transform_sym_matrix(put_matrix_XZY(), ans0), mprod( m_S_red.second, put_matrix_XZY() ) ) );
317 }
318
319 ans0 = m_S_red.first;
320 cal_average_crystal_system(D4h_Z, ans0);
321 if( check_equiv_m(ans0, m_S_red.first, resol ) )
322 {
323 ans.insert( SymMat43_Double(ans0, m_S_red.second ) );
324 }
325
326 return !( ans.empty() );
327 }
328
329
330
331
332 bool LatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsHexagonal(const ePointGroup& epg_new, const Double& resol,
333 map< SymMat<Double>, NRMat<Int4> >& ans) const
334 {
335 ans.clear();
336 const BravaisType& brat = m_latfom.putBravaisType();
337
338 SymMat43_Double ans2(SymMat<Double>(3), NRMat<Int4>(3,3));
339
340 if( brat.enumLaueGroup() == C2h_X )
341 {
342 ans2.first = transform_sym_matrix(put_matrix_YZX(), m_S_red.first);
343 ans2.second = mprod( m_S_red.second, put_matrix_ZXY() );
344 }
345 else if( brat.enumLaueGroup() == C2h_Y )
346 {
347 ans2.first = transform_sym_matrix(put_matrix_ZXY(), m_S_red.first);
348 ans2.second = mprod( m_S_red.second, put_matrix_YZX() );
349 }
350 else // if( brat.enumLaueGroup() == C2h_Z )
351 {
352 ans2.first = transform_sym_matrix(put_matrix_XYZ(), m_S_red.first);
353 ans2.second = m_S_red.second;
354 }
355
356 if( ans2.first(0,1) < 0.0 )
357 {
358 ans2.first(0,1) *= -1;
359 ans2.second[0][0] *= -1;
360 ans2.second[1][0] *= -1;
361 ans2.second[2][0] *= -1;
362 }
363
364 SymMat<Double> ans0 = ans2.first;
365 cal_average_crystal_system(epg_new, ans2.first);
366 if( check_equiv_m(ans2.first, ans0, resol ) )
367 {
368 ans.insert( ans2 );
369 return true;
370 }
371 else return false;
372 }
373
374
375 bool LatticeFigureOfMeritToCheckSymmetry::checkLatticeSymmetry(const ePointGroup& epg_new, const Double& resol,
376 map< SymMat<Double>, NRMat<Int4> >& ans) const
377 {
378 ans.clear();
379 const BravaisType& brat = m_latfom.putBravaisType();
380 if( epg_new == Ci || epg_new == brat.enumLaueGroup() )
381 {
382 ans.insert( m_S_red );
383 return true;
384 }
385
386 if( epg_new == C2h_X || epg_new == C2h_Y || epg_new == C2h_Z )
387 {
388 assert( brat.enumLaueGroup() == Ci );
389 assert( brat.enumCentringType() == Prim );
390
391 return checkIfLatticeIsMonoclinic(epg_new, resol, ans);
392 }
393 else if( epg_new == D4h_Z )
394 {
395 assert( brat.enumLaueGroup() == D2h );
396 assert( brat.enumCentringType() == Prim
397 || brat.enumCentringType() == Inner );
398
399 return checkIfLatticeIsTetragonal(resol, ans);
400 }
401 else if( epg_new == D2h )
402 {
403 assert( brat.enumLaueGroup() != Ci || brat.enumCentringType() == Prim );
404 assert( brat.enumLaueGroup() != C2h_Z || brat.enumCentringType() == BaseX );
405 assert( brat.enumLaueGroup() != C2h_X || brat.enumCentringType() == BaseY );
406 assert( brat.enumLaueGroup() != C2h_Y || brat.enumCentringType() == BaseZ );
407 assert( brat.enumCentringType() != Rhom_hex );
408
409 return checkIfLatticeIsOrthorhombic(resol, ans);
410 }
411 else if( epg_new == D6h )
412 {
413 assert( brat.enumCentringType() == Prim );
414 assert( brat.enumLaueGroup() == C2h_X
415 || brat.enumLaueGroup() == C2h_Y
416 || brat.enumLaueGroup() == C2h_Z );
417 return checkIfLatticeIsHexagonal(epg_new, resol, ans);
418 }
419 else
420 {
421 assert( epg_new == Oh );
422 assert( brat.enumCentringType() == Prim
423 || brat.enumCentringType() == Inner
424 || brat.enumCentringType() == Face );
425
426 SymMat43_Double ans2 = m_S_red;
427 cal_average_crystal_system(epg_new, ans2.first);
428 if( check_equiv_m(ans2.first, m_S_red.first, resol ) )
429 {
430 ans.insert( ans2 );
431 return true;
432 }
433 }
434 return !(ans.empty());
435 }
436
437
438 void LatticeFigureOfMeritToCheckSymmetry::putLatticesOfHigherSymmetry(
439 const ePointGroup& epg, const Double& resol,
440 vector<LatticeFigureOfMeritToCheckSymmetry>& lattice_result) const
441 {
442 lattice_result.clear();
443 map< SymMat<Double>, NRMat<Int4> > S_red_tray;
444 if( !this->checkLatticeSymmetry(epg, resol, S_red_tray) ) return;
445
446 const BravaisType& ebrat_original = this->putLatticeFigureOfMerit().putBravaisType();
447 const eCentringType eblat = (ebrat_original.enumBravaisType()==Monoclinic_B?
448 (epg==D31d_rho?Prim:(epg==D3d_1_hex?Rhom_hex:BaseZ)):ebrat_original.enumCentringType());
449
450 const NRMat<Int4> matrix_min_to_sell = this->putInitialForm().second;
451
452 SymMat<Double> S_super(4);
453 NRMat<Int4> trans_mat(4,3);
454
455 for(map< SymMat<Double>, NRMat<Int4> >::const_iterator it=S_red_tray.begin(); it!=S_red_tray.end(); it++)
456 {
457 // S_super = it->second * it->first * Transpose(it->second) is close to
458 // Delone-reduced form of the original lattice.
459 S_super = transform_sym_matrix(it->second, it->first );
460
461 trans_mat = identity_matrix<Int4>(4);
462
463 // S_super = trans_mat * it->second * it->first * Transpose(trans_mat * it->second).
464 put_Selling_reduced_dim_less_than_4(S_super, trans_mat);
465 moveSmallerDiagonalLeftUpper(S_super, trans_mat);
466
467 lattice_result.push_back( LatticeFigureOfMeritToCheckSymmetry( BravaisType( pair<eCentringType, ePointGroup>(eblat, epg) ),
468 SymMat43_Double(it->first, mprod(trans_mat, it->second) ) ) );
469 }
470 }
Back to OSDN">Back to OSDN
 ViewVC Help
Powered by ViewVC 1.1.26  
About OSDN
Find Software
Develop Software
Help