Develop and Download Open Source Software

Browse Subversion Repository

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 33 - (show annotations) (download) (as text)
Wed Sep 7 04:38:51 2016 UTC (7 years, 6 months ago) by rtomiyasu
File MIME type: text/x-c++src
File size: 16419 byte(s)
The output format for base-centered monoclinic cells was corrected.
1 /*
2 * The MIT License
3
4 Conograph (powder auto-indexing program)
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 "../utility_func/chToDouble.hh"
28 #include "../utility_lattice_reduction/put_Buerger_reduced_lattice.hh"
29 #include "OutputInfo.hh"
30 #include "VCLatticeFigureOfMeritToCheckSymmetry.hh"
31
32 const string VCLatticeFigureOfMeritToCheckSymmetry::CS_LABEL[NUM_LS] =
33 { "01", "02", "03", "04", "05", "06", "07", "08", "09", "10", "11", "12", "13", "14" };
34
35 VCLatticeFigureOfMeritToCheckSymmetry::VCLatticeFigureOfMeritToCheckSymmetry(const Double& rhs)
36 : m_label(-1), m_latfom(rhs),
37 m_S_red( SymMat43_VCData( SymMat<VCData>(3), NRMat<Int4>(4,3) ) )
38 {
39 }
40
41
42 VCLatticeFigureOfMeritToCheckSymmetry::VCLatticeFigureOfMeritToCheckSymmetry(const BravaisType& brat,
43 const SymMat43_VCData& S,
44 const ePeakShiftFunctionType& type,
45 const Double& wave_length,
46 const vector<ZParawError>& peak_shift_param_rad)
47 : m_label(-1),
48 m_S_red( SymMat43_VCData( SymMat<VCData>(3), NRMat<Int4>(4,3) ) )
49 {
50 this->setLatticeConstants43(brat, S);
51 m_latfom.setPeakShiftParamRadian(type, wave_length, peak_shift_param_rad);
52 m_num_lattice_found = 0;
53 }
54
55
56 #ifdef DEBUG
57 static bool checkInitialLatticeParameters(
58 const BravaisType& brat,
59 const SymMat43_VCData& S_red)
60 {
61 const SymMat<Double> dbl_S_red( chToDouble(S_red.first) );
62
63 if( brat.enumLaueGroup() == Ci && brat.enumCentringType() == Prim )
64 {
65 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)
66 && fabs( dbl_S_red(0,1) ) * 1.9999 < dbl_S_red(1,1)
67 && fabs( dbl_S_red(0,2) ) * 1.9999 < dbl_S_red(2,2)
68 && fabs( dbl_S_red(1,2) ) * 1.9999 < dbl_S_red(2,2) );
69 }
70 else if( brat.enumLaueGroup() == C2h_Y && brat.enumCentringType() == Prim )
71 {
72 assert( 0.0 <= dbl_S_red(0,2) && dbl_S_red(2,2)*0.9999 < dbl_S_red(0,0)
73 && 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) );
74 }
75 else if( brat.enumLaueGroup() == C2h_Z && brat.enumCentringType() == Prim )
76 {
77 assert( 0.0 <= dbl_S_red(0,1) && dbl_S_red(1,1)*0.9999 < dbl_S_red(0,0)
78 && 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) );
79 }
80 else if( brat.enumLaueGroup() == C2h_X && brat.enumCentringType() == Prim )
81 {
82 assert( 0.0 <= dbl_S_red(1,2) && dbl_S_red(2,2)*0.9999 < dbl_S_red(1,1)
83 && 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) );
84 }
85 else if( brat.enumLaueGroup() == C2h_Y && brat.enumCentringType() == BaseZ )
86 {
87 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) );
88 }
89 else if( brat.enumLaueGroup() == C2h_Z && brat.enumCentringType() == BaseX )
90 {
91 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) );
92 }
93 else if( brat.enumLaueGroup() == C2h_X && brat.enumCentringType() == BaseY )
94 {
95 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) );
96 }
97 else if( brat.enumLaueGroup() == D2h
98 && brat.enumCentringType() != BaseX
99 && brat.enumCentringType() != BaseY
100 && brat.enumCentringType() != BaseZ )
101 {
102 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) );
103 }
104
105 const SymMat<VCData> S_super = transform_sym_matrix(S_red.second, S_red.first);
106 assert( S_super(0,1) <= VCData()
107 && S_super(0,2) <= VCData()
108 && S_super(0,3) <= VCData()
109 && S_super(1,2) <= VCData()
110 && S_super(1,3) <= VCData()
111 && S_super(2,3) <= VCData() );
112
113 SymMat<VCData> S_red_cp = S_red.first;
114 cal_average_crystal_system(brat.enumLaueGroup(), S_red_cp);
115 assert( S_red.first(0,0).putVecCoef() == S_red_cp(0,0).putVecCoef() );
116 assert( S_red.first(1,1).putVecCoef() == S_red_cp(1,1).putVecCoef() );
117 assert( S_red.first(2,2).putVecCoef() == S_red_cp(2,2).putVecCoef() );
118 assert( S_red.first(0,1).putVecCoef() == S_red_cp(0,1).putVecCoef() );
119 assert( S_red.first(0,2).putVecCoef() == S_red_cp(0,2).putVecCoef() );
120 assert( S_red.first(1,2).putVecCoef() == S_red_cp(1,2).putVecCoef() );
121
122 return true;
123 }
124 #endif
125
126
127 void VCLatticeFigureOfMeritToCheckSymmetry::setLatticeConstants43(const BravaisType& brat, const SymMat43_VCData& S)
128 {
129 m_S_red = S;
130 assert( checkInitialLatticeParameters(brat, m_S_red) );
131
132 m_latfom.setLatticeConstants43(brat, SymMat43_Double(chToDouble(m_S_red.first), m_S_red.second));
133 m_num_lattice_found = 0;
134 }
135
136
137
138
139 bool VCLatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsMonoclinic(const ePointGroup& epg_new,
140 const Double& cv2,
141 map< SymMat<VCData>, NRMat<Int4> >& ans) const
142 {
143 ans.clear();
144
145 SymMat<VCData> ans0 = m_S_red.first;
146 cal_average_crystal_system(C2h_X, ans0);
147
148 SymMat<VCData> S_red(3);
149 NRMat<Int4> trans_mat2;
150 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
151 {
152 if( epg_new == C2h_X )
153 {
154 S_red = ans0;
155 trans_mat2 = m_S_red.second;
156 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
157 }
158 else if( epg_new == C2h_Y )
159 {
160 S_red = transform_sym_matrix(put_matrix_YXZ(), ans0);
161 trans_mat2 = mprod(m_S_red.second, put_matrix_YXZ());
162 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
163 }
164 else // if( epg_new == C2h_Z )
165 {
166 S_red = transform_sym_matrix(put_matrix_YZX(), ans0);
167 trans_mat2 = mprod(m_S_red.second, put_matrix_ZXY());
168 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
169 }
170 ans.insert( SymMat43_VCData( S_red, trans_mat2) );
171 }
172
173 ans0 = m_S_red.first;
174 cal_average_crystal_system(C2h_Y, ans0);
175 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
176 {
177 if( epg_new == C2h_X )
178 {
179 S_red = transform_sym_matrix(put_matrix_YXZ(), ans0);
180 trans_mat2 = mprod(m_S_red.second, put_matrix_YXZ());
181 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
182 }
183 else if( epg_new == C2h_Y )
184 {
185 S_red = ans0;
186 trans_mat2 = m_S_red.second;
187 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
188 }
189 else // if( epg_new == C2h_Z )
190 {
191 S_red = transform_sym_matrix(put_matrix_XZY(), ans0);
192 trans_mat2 = mprod(m_S_red.second, put_matrix_XZY());
193 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
194 }
195 ans.insert( SymMat43_VCData( S_red, trans_mat2) );
196 }
197
198 ans0 = m_S_red.first;
199 cal_average_crystal_system(C2h_Z, ans0);
200 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
201 {
202 if( epg_new == C2h_X )
203 {
204 S_red = transform_sym_matrix(put_matrix_ZXY(), ans0);
205 trans_mat2 = mprod(m_S_red.second, put_matrix_YZX());
206 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
207 }
208 else if( epg_new == C2h_Y )
209 {
210 S_red = transform_sym_matrix(put_matrix_XZY(), ans0);
211 trans_mat2 = mprod(m_S_red.second, put_matrix_XZY());
212 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
213 }
214 else // if( epg_new == C2h_Z )
215 {
216 S_red = ans0;
217 trans_mat2 = m_S_red.second;
218 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
219 }
220 ans.insert( SymMat43_VCData( S_red, trans_mat2) );
221 }
222
223 return !( ans.empty() );
224 }
225
226
227 bool VCLatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsOrthorhombic(const Double& cv2,
228 map< SymMat<VCData>, NRMat<Int4> >& ans) const
229 {
230 ans.clear();
231
232 const BravaisType& brat = m_latfom.putBravaisType();
233
234 SymMat<VCData> ans0 = m_S_red.first;
235 cal_average_crystal_system(D2h, ans0);
236 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
237 {
238 if( brat.enumCentringType() == BaseX )
239 {
240 if( ans0(1,1) < ans0(2,2) )
241 {
242 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_ZYX(), ans0), mprod( m_S_red.second, put_matrix_ZYX() ) ) );
243 }
244 else
245 {
246 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_YZX(), ans0), mprod( m_S_red.second, put_matrix_ZXY() ) ) );
247 }
248 }
249 else if( brat.enumCentringType() == BaseY )
250 {
251 if( ans0(0,0) < ans0(2,2) )
252 {
253 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_ZXY(), ans0), mprod( m_S_red.second, put_matrix_YZX() ) ) );
254 }
255 else
256 {
257 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_XZY(), ans0), mprod( m_S_red.second, put_matrix_XZY() ) ) );
258 }
259 }
260 else if( brat.enumCentringType() == BaseZ )
261 {
262 if( ans0(0,0) < ans0(1,1) )
263 {
264 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_YXZ(), ans0), mprod( m_S_red.second, put_matrix_YXZ() ) ) );
265 }
266 else
267 {
268 ans.insert( SymMat43_VCData( transform_sym_matrix(put_matrix_XYZ(), ans0), m_S_red.second ) );
269 }
270 }
271 else
272 {
273 NRMat<Int4> trans_mat = m_S_red.second;
274 putBuergerReducedOrthorhombic(brat.enumCentringType(), ans0, trans_mat);
275 ans.insert( SymMat43_VCData(ans0, trans_mat ) );
276 }
277 return true;
278 }
279 return false;
280 }
281
282
283 bool VCLatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsTetragonal(const Double& cv2,
284 map< SymMat<VCData>, NRMat<Int4> >& ans) const
285 {
286 ans.clear();
287
288 SymMat<VCData> ans0 = m_S_red.first;
289 cal_average_crystal_system(D4h_X, ans0);
290 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
291 {
292 ans.insert( SymMat43_VCData(
293 transform_sym_matrix(put_matrix_YZX(), ans0), mprod( m_S_red.second, put_matrix_ZXY() ) ) );
294 }
295
296 ans0 = m_S_red.first;
297 cal_average_crystal_system(D4h_Y, ans0);
298 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
299 {
300 ans.insert( SymMat43_VCData(
301 transform_sym_matrix(put_matrix_XZY(), ans0), mprod( m_S_red.second, put_matrix_XZY() ) ) );
302 }
303
304 ans0 = m_S_red.first;
305 cal_average_crystal_system(D4h_Z, ans0);
306 if( check_equiv_m(ans0, m_S_red.first, cv2 ) )
307 {
308 ans.insert( SymMat43_VCData(ans0, m_S_red.second ) );
309 }
310
311 return !( ans.empty() );
312 }
313
314
315
316
317 bool VCLatticeFigureOfMeritToCheckSymmetry::checkIfLatticeIsHexagonal(const ePointGroup& epg_new, const Double& cv2,
318 map< SymMat<VCData>, NRMat<Int4> >& ans) const
319 {
320 ans.clear();
321 const BravaisType& brat = m_latfom.putBravaisType();
322
323 SymMat43_VCData ans2(SymMat<VCData>(3), NRMat<Int4>(3,3));
324
325 if( brat.enumLaueGroup() == C2h_X )
326 {
327 ans2.first = transform_sym_matrix(put_matrix_YZX(), m_S_red.first);
328 ans2.second = mprod( m_S_red.second, put_matrix_ZXY() );
329 }
330 else if( brat.enumLaueGroup() == C2h_Y )
331 {
332 ans2.first = transform_sym_matrix(put_matrix_ZXY(), m_S_red.first);
333 ans2.second = mprod( m_S_red.second, put_matrix_YZX() );
334 }
335 else // if( brat.enumLaueGroup() == C2h_Z )
336 {
337 ans2.first = transform_sym_matrix(put_matrix_XYZ(), m_S_red.first);
338 ans2.second = m_S_red.second;
339 }
340
341 if( ans2.first(0,1) < VCData() )
342 {
343 ans2.first(0,1) *= -1;
344 ans2.second[0][0] *= -1;
345 ans2.second[1][0] *= -1;
346 ans2.second[2][0] *= -1;
347 }
348
349 SymMat<VCData> ans0 = ans2.first;
350 cal_average_crystal_system(epg_new, ans2.first);
351 if( check_equiv_m(ans2.first, ans0, cv2 ) )
352 {
353 ans.insert( ans2 );
354 return true;
355 }
356 else return false;
357 }
358
359
360 bool VCLatticeFigureOfMeritToCheckSymmetry::checkLatticeSymmetry(const ePointGroup& epg_new, const Double& cv2,
361 map< SymMat<VCData>, NRMat<Int4> >& ans) const
362 {
363 ans.clear();
364 const BravaisType& brat = m_latfom.putBravaisType();
365 if( epg_new == Ci || epg_new == brat.enumLaueGroup() )
366 {
367 ans.insert( m_S_red );
368 return true;
369 }
370
371 if( epg_new == C2h_X || epg_new == C2h_Y || epg_new == C2h_Z )
372 {
373 assert( brat.enumLaueGroup() == Ci );
374 assert( brat.enumCentringType() == Prim );
375
376 return checkIfLatticeIsMonoclinic(epg_new, cv2, ans);
377 }
378 else if( epg_new == D4h_Z )
379 {
380 assert( brat.enumLaueGroup() == D2h );
381 assert( brat.enumCentringType() == Prim
382 || brat.enumCentringType() == Inner );
383
384 return checkIfLatticeIsTetragonal(cv2, ans);
385 }
386 else if( epg_new == D2h )
387 {
388 assert( brat.enumLaueGroup() != Ci || brat.enumCentringType() == Prim );
389 assert( brat.enumLaueGroup() != C2h_Z || brat.enumCentringType() == BaseX );
390 assert( brat.enumLaueGroup() != C2h_X || brat.enumCentringType() == BaseY );
391 assert( brat.enumLaueGroup() != C2h_Y || brat.enumCentringType() == BaseZ );
392 assert( brat.enumCentringType() != Rhom_hex );
393
394 return checkIfLatticeIsOrthorhombic(cv2, ans);
395 }
396 else if( epg_new == D6h )
397 {
398 assert( brat.enumCentringType() == Prim );
399 assert( brat.enumLaueGroup() == C2h_X
400 || brat.enumLaueGroup() == C2h_Y
401 || brat.enumLaueGroup() == C2h_Z );
402 return checkIfLatticeIsHexagonal(epg_new, cv2, ans);
403 }
404 else
405 {
406 assert( epg_new == Oh );
407 assert( brat.enumCentringType() == Prim
408 || brat.enumCentringType() == Inner
409 || brat.enumCentringType() == Face );
410
411 SymMat43_VCData ans2 = m_S_red;
412 cal_average_crystal_system(epg_new, ans2.first);
413 if( check_equiv_m(ans2.first, m_S_red.first, cv2 ) )
414 {
415 ans.insert( ans2 );
416 return true;
417 }
418 }
419 return !(ans.empty());
420 }
421
422
423 void VCLatticeFigureOfMeritToCheckSymmetry::putLatticesOfHigherSymmetry(
424 const ePointGroup& epg, const Double& cv2,
425 vector<VCLatticeFigureOfMeritToCheckSymmetry>& lattice_result) const
426 {
427 lattice_result.clear();
428 map< SymMat<VCData>, NRMat<Int4> > S_red_tray;
429 if( !this->checkLatticeSymmetry(epg, cv2, S_red_tray) ) return;
430
431 const BravaisType& ebrat_original = this->putLatticeFigureOfMerit().putBravaisType();
432 const eCentringType eblat = (ebrat_original.enumBravaisType()==Monoclinic_B?
433 (epg==D31d_rho?Prim:(epg==D3d_1_hex?Rhom_hex:BaseZ)):ebrat_original.enumCentringType());
434
435 const NRMat<Int4> matrix_min_to_sell = this->putInitialForm().second;
436
437 SymMat<Double> S_super(4);
438 NRMat<Int4> trans_mat(4,3);
439
440 for(map< SymMat<VCData>, NRMat<Int4> >::const_iterator it=S_red_tray.begin(); it!=S_red_tray.end(); it++)
441 {
442 // S_super = it->second * it->first * Transpose(it->second) is close to
443 // Delone-reduced form of the original lattice.
444 S_super = transform_sym_matrix(it->second, chToDouble(it->first) );
445
446 trans_mat = identity_matrix<Int4>(4);
447
448 // S_super = trans_mat * it->second * it->first * Transpose(trans_mat * it->second).
449 put_Selling_reduced_dim_less_than_4(S_super, trans_mat);
450 moveSmallerDiagonalLeftUpper(S_super, trans_mat);
451
452 lattice_result.push_back( VCLatticeFigureOfMeritToCheckSymmetry( BravaisType( pair<eCentringType, ePointGroup>(eblat, epg) ),
453 SymMat43_VCData(it->first, mprod(trans_mat, it->second) ),
454 this->putLatticeFigureOfMerit().putPeakShiftFunctionType(),
455 this->putLatticeFigureOfMerit().putWaveLength(),
456 this->putLatticeFigureOfMerit().putPeakShiftParamRadian() ) );
457 }
458 }
459
460
461 void VCLatticeFigureOfMeritToCheckSymmetry::printLatticeInformation(
462 const vector<VCLatticeFigureOfMeritToCheckSymmetry> lattice_result[],
463 const OutputInfo outinfo[],
464 const eABCaxis& abc_axis,
465 const eRHaxis& rh_axis,
466 const Double& resol,
467 const Int4& label_start0,
468 ostream* os) const
469 {
470 m_latfom.printLatticeInformation(abc_axis, rh_axis, resol, label_start0, os);
471
472 Int4 label_start = label_start0;
473 os->width(label_start);
474 *os << "" << "<NumberOfLatticesInNeighborhood>";
475 os->width(14);
476 *os << this->putNumberOfLatticesInNeighborhood();
477 *os << " </NumberOfLatticesInNeighborhood>\n\n";
478 }

Back to OSDN">Back to OSDN
ViewVC Help
Powered by ViewVC 1.1.26