Develop and Download Open Source Software

Powder auto-indexing software Conograph

Browse Subversion Repository

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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 33 - (hide 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 rtomiyasu 25 /*
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
157 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
163 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
169 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
182 rtomiyasu 25 }
183     else if( epg_new == C2h_Y )
184     {
185     S_red = ans0;
186     trans_mat2 = m_S_red.second;
187 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
188 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
194 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(1, 2, S_red, trans_mat2);
207 rtomiyasu 25 }
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 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 2, S_red, trans_mat2);
213 rtomiyasu 25 }
214     else // if( epg_new == C2h_Z )
215     {
216     S_red = ans0;
217     trans_mat2 = m_S_red.second;
218 rtomiyasu 33 putBuergerReducedMonoclinicP<VCData>(0, 1, S_red, trans_mat2);
219 rtomiyasu 25 }
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  
About OSDN
Find Software
Develop Software
Help