source: trunk/GSASIIspc.py @ 946

Last change on this file since 946 was 946, checked in by toby, 10 years ago

update self-docs, start work on constraints object

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 54.2 KB
Line 
1# -*- coding: utf-8 -*-
2"""
3*GSASIIspc: Space group module*
4-------------------------------
5
6Space group interpretation routines. Note that space group information is
7stored in a :ref:`Space Group (SGData)<SGData_table>` object.
8
9"""
10########### SVN repository information ###################
11# $Date: 2013-06-11 15:08:43 +0000 (Tue, 11 Jun 2013) $
12# $Author: toby $
13# $Revision: 946 $
14# $URL: trunk/GSASIIspc.py $
15# $Id: GSASIIspc.py 946 2013-06-11 15:08:43Z toby $
16########### SVN repository information ###################
17import numpy as np
18import numpy.ma as ma
19import numpy.linalg as nl
20import math
21import sys
22import os.path as ospath
23
24import GSASIIpath
25GSASIIpath.SetVersionNumber("$Revision: 946 $")
26import pyspg
27
28def SpcGroup(SGSymbol):
29    """
30    Determines cell and symmetry information from a short H-M space group name
31
32    :param SGSymbol: space group symbol (string) with spaces between axial fields
33    :returns: (SGError,SGData)
34       * SGError = 0 for no errors; >0 for errors (see SGErrors below for details)
35       * SGData - is a dict (see :ref:`Space Group object<SGData_table>`) with entries:
36       
37             * 'SpGrp': space group symbol, slightly cleaned up
38             * 'Laue':  one of '-1', '2/m', 'mmm', '4/m', '4/mmm', '3R',
39               '3mR', '3', '3m1', '31m', '6/m', '6/mmm', 'm3', 'm3m'
40             * 'SGInv': boolean; True if centrosymmetric, False if not
41             * 'SGLatt': one of 'P', 'A', 'B', 'C', 'I', 'F', 'R'
42             * 'SGUniq': one of 'a', 'b', 'c' if monoclinic, '' otherwise
43             * 'SGCen': cell centering vectors [0,0,0] at least
44             * 'SGOps': symmetry operations as [M,T] so that M*x+T = x'
45             * 'SGSys': one of 'triclinic', 'monoclinic', 'orthorhombic',
46               'tetragonal', 'rhombohedral', 'trigonal', 'hexagonal', 'cubic'
47             * 'SGPolax': one of '', 'x', 'y', 'x y', 'z', 'x z', 'y z',
48               'xyz', '111' for arbitrary axes
49
50    """
51    LaueSym = ('-1','2/m','mmm','4/m','4/mmm','3R','3mR','3','3m1','31m','6/m','6/mmm','m3','m3m')
52    LattSym = ('P','A','B','C','I','F','R')
53    UniqSym = ('','','a','b','c','',)
54    SysSym = ('triclinic','monoclinic','orthorhombic','tetragonal','rhombohedral','trigonal','hexagonal','cubic')
55    SGData = {}
56    SGData['SpGrp'] = SGSymbol.strip().lower().capitalize()
57    SGInfo = pyspg.sgforpy(SGSymbol)
58    SGData['SGLaue'] = LaueSym[SGInfo[0]-1]
59    SGData['SGInv'] = bool(SGInfo[1])
60    SGData['SGLatt'] = LattSym[SGInfo[2]-1]
61    SGData['SGUniq'] = UniqSym[SGInfo[3]+1]
62    if SGData['SGLatt'] == 'P':
63        SGData['SGCen'] = np.array(([0,0,0],))
64    elif SGData['SGLatt'] == 'A':
65        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5]))
66    elif SGData['SGLatt'] == 'B':
67        SGData['SGCen'] = np.array(([0,0,0],[.5,0,.5]))
68    elif SGData['SGLatt'] == 'C':
69        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,0,]))
70    elif SGData['SGLatt'] == 'I':
71        SGData['SGCen'] = np.array(([0,0,0],[.5,.5,.5]))
72    elif SGData['SGLatt'] == 'F':
73        SGData['SGCen'] = np.array(([0,0,0],[0,.5,.5],[.5,0,.5],[.5,.5,0,]))
74    elif SGData['SGLatt'] == 'R':
75        SGData['SGCen'] = np.array(([0,0,0],[1./3.,2./3.,2./3.],[2./3.,1./3.,1./3.]))
76    SGData['SGOps'] = []
77    for i in range(SGInfo[5]):
78        Mat = np.array(SGInfo[6][i])
79        Trns = np.array(SGInfo[7][i])
80        SGData['SGOps'].append([Mat,Trns])
81    if SGData['SGLaue'] in '-1':
82        SGData['SGSys'] = SysSym[0]
83    elif SGData['SGLaue'] in '2/m':
84        SGData['SGSys'] = SysSym[1]
85    elif SGData['SGLaue'] in 'mmm':
86        SGData['SGSys'] = SysSym[2]
87    elif SGData['SGLaue'] in ['4/m','4/mmm']:
88        SGData['SGSys'] = SysSym[3]
89    elif SGData['SGLaue'] in ['3R','3mR']:
90        SGData['SGSys'] = SysSym[4]
91    elif SGData['SGLaue'] in ['3','3m1','31m']:
92        SGData['SGSys'] = SysSym[5]
93    elif SGData['SGLaue'] in ['6/m','6/mmm']:
94        SGData['SGSys'] = SysSym[6]
95    elif SGData['SGLaue'] in ['m3','m3m']:
96        SGData['SGSys'] = SysSym[7]
97    SGData['SGPolax'] = SGpolar(SGData)
98    return SGInfo[8],SGData
99
100def SGErrors(IErr):
101    '''
102    Interprets the error message code from SpcGroup. Used in SpaceGroup.
103   
104    :param IErr: see SGError in :func:`SpcGroup`
105    :returns:
106        ErrString - a string with the error message or "Unknown error"
107    '''
108
109    ErrString = [' ',
110        'Less than 2 operator fields were found',
111        'Illegal Lattice type, not P, A, B, C, I, F or R',
112        'Rhombohedral lattice requires a 3-axis',
113        'Minus sign does not preceed 1, 2, 3, 4 or 6',
114        'Either a 5-axis anywhere or a 3-axis in field not allowed',
115        ' ',
116        'I for COMPUTED GO TO out of range.',
117        'An a-glide mirror normal to A not allowed',
118        'A b-glide mirror normal to B not allowed',
119        'A c-glide mirror normal to C not allowed',
120        'D-glide in a primitive lattice not allowed',
121        'A 4-axis not allowed in the 2nd operator field',
122        'A 6-axis not allowed in the 2nd operator field',
123        'More than 24 matrices needed to define group',
124        ' ',
125        'Improper construction of a rotation operator',
126        'Mirror following a / not allowed',
127        'A translation conflict between operators',
128        'The 2bar operator is not allowed',
129        '3 fields are legal only in R & m3 cubic groups',
130        'Syntax error. Expected I -4 3 d at this point',
131        ' ',
132        'A or B centered tetragonal not allowed',
133        ' ','unknown error in sgroup',' ',' ',' ',
134        'Illegal character in the space group symbol',
135        ]
136    try:
137        return ErrString[IErr]
138    except:
139        return "Unknown error"
140
141def SGpolar(SGData):
142    '''
143    Determine identity of polar axes if any
144    '''
145    POL = ('','x','y','x y','z','x z','y z','xyz','111')
146    NP = [1,2,4]
147    NPZ = [0,1]
148    for M,T in SGData['SGOps']:
149        for i in range(3):
150            if M[i][i] <= 0.: NP[i] = 0
151        if M[0][2] > 0: NPZ[0] = 8
152        if M[1][2] > 0: NPZ[1] = 0
153    NPol = (NP[0]+NP[1]+NP[2]+NPZ[0]*NPZ[1])*(1-int(SGData['SGInv']))
154    return POL[NPol]
155   
156def SGPrint(SGData):
157    '''
158    Print the output of SpcGroup in a nicely formatted way. Used in SpaceGroup
159
160    :param SGData: from :func:`SpcGroup`
161    :returns:
162        SGText - list of strings with the space group details
163    '''
164    Mult = len(SGData['SGCen'])*len(SGData['SGOps'])*(int(SGData['SGInv'])+1)
165    SGText = []
166    SGText.append(' Space Group: '+SGData['SpGrp'])
167    CentStr = 'centrosymmetric'
168    if not SGData['SGInv']:
169        CentStr = 'non'+CentStr
170    if SGData['SGLatt'] in 'ABCIFR':
171        SGText.append(' The lattice is '+CentStr+' '+SGData['SGLatt']+'-centered '+SGData['SGSys'].lower())
172    else:
173        SGText.append(' The lattice is '+CentStr+' '+'primitive '+SGData['SGSys'].lower())       
174    SGText.append(' Multiplicity of a general site is '+str(Mult))
175    SGText.append(' The Laue symmetry is '+SGData['SGLaue'])
176    if SGData['SGUniq'] in ['a','b','c']:
177        SGText.append(' The unique monoclinic axis is '+SGData['SGUniq'])
178    if SGData['SGInv']:
179        SGText.append(' The inversion center is located at 0,0,0')
180    if SGData['SGPolax']:
181        SGText.append(' The location of the origin is arbitrary in '+SGData['SGPolax'])
182    SGText.append('\n'+' The equivalent positions are:')
183    if SGData['SGLatt'] != 'P':
184        SGText.append('\n ('+Latt2text(SGData['SGLatt'])+')+\n')
185    Ncol = 2
186    line = ' '
187    col = 0
188    for iop,[M,T] in enumerate(SGData['SGOps']):
189        OPtxt = MT2text(M,T)
190        Fld = '(%2i) '%(iop+1)+OPtxt+'\t'
191        line += Fld
192        if '/' not in Fld:
193            line += '\t'
194        col += 1
195        if col == Ncol:
196            SGText.append(line)       
197            line = ' '
198            col = 0
199    SGText.append(line)       
200    return SGText
201   
202def MT2text(M,T):
203    "From space group matrix/translation operator returns text version"
204    XYZ = ('-Z ','-Y ','-X ','X-Y','ERR','Y-X',' X ',' Y ',' Z ','+X ','+Y ','+Z ')
205    TRA = ('   ','ERR','1/6','1/4','1/3','ERR','1/2','ERR','2/3','3/4','5/6','ERR')
206    Fld = ''
207    for j in range(3):
208        IJ = int(round(2*M[j][0]+3*M[j][1]+4*M[j][2]+4))%12
209        IK = int(round(T[j]*12))%12
210        if IK > 0 and IJ > 4: IJ += 3
211        Fld += TRA[IK]+XYZ[IJ]
212        if j != 2: Fld += ','
213    return Fld
214   
215def Latt2text(Latt):
216    "From lattice type ('P',A', etc.) returns ';' delimited cell centering vectors"
217    lattTxt = {'A':'0,0,0; 0,1/2,1/2','B':'0,0,0; 1/2,0,1/2',
218        'C':'0,0,0; 1/2,1/2,0','I':'0,0,0; 1/2,1/2,1/2',
219        'F':'0,0,0; 0,1/2,1/2; 1/2,0,1/2; 1/2,1/2,0',
220        'R':'0,0,0; 1/3,2/3,2/3; 2/3,1/3,1/3','P':'0,0,0'}
221    return lattTxt[Latt]   
222       
223def SpaceGroup(SGSymbol):
224    '''
225    Print the output of SpcGroup in a nicely formatted way.
226
227    :param SGSymbol: space group symbol (string) with spaces between axial fields
228    :returns: nothing
229    '''
230    E,A = SpcGroup(SGSymbol)
231    if E > 0:
232        print SGErrors(E)
233        return
234    for l in SGPrint(A):
235        print l
236
237def MoveToUnitCell(xyz):
238    '''
239    Translates a set of coordinates so that all values are >=0 and < 1
240
241    :param xyz: a list or numpy array of fractional coordinates
242    :returns: XYZ - numpy array of new coordinates now 0 or greater and less than 1
243    '''
244    XYZ = np.zeros(3)
245    for i,x in enumerate(xyz):
246        XYZ[i] = (x-int(x))%1.0
247    return XYZ
248       
249def Opposite(XYZ,toler=0.0002):
250    '''
251    Gives opposite corner, edge or face of unit cell for position within tolerance.
252        Result may be just outside the cell within tolerance
253
254    :param XYZ: 0 >= np.array[x,y,z] > 1 as by MoveToUnitCell
255    :param toler: unit cell fraction tolerance making opposite
256    :returns:
257        XYZ: array of opposite positions; always contains XYZ
258    '''
259    perm3 = [[1,1,1],[0,1,1],[1,0,1],[1,1,0],[1,0,0],[0,1,0],[0,0,1],[0,0,0]]
260    TB = np.where(abs(XYZ-1)<toler,-1,0)+np.where(abs(XYZ)<toler,1,0)
261    perm = TB*perm3
262    cperm = ['%d%d%d'%(i,j,k) for i,j,k in perm]
263    D = dict(zip(cperm,perm))
264    new = []
265    for key in D:
266        new.append(np.array(D[key])+np.array(XYZ))
267    return new
268       
269def GenAtom(XYZ,SGData,All=False,Uij=[],Move=True):
270    '''
271    Generates the equivalent positions for a specified coordinate and space group
272
273    :param XYZ: an array, tuple or list containing 3 elements: x, y & z
274    :param SGData: from :func:`SpcGroup`
275    :param All: True return all equivalent positions including duplicates;
276      False return only unique positions
277    :param Uij: [U11,U22,U33,U12,U13,U23] or [] if no Uij
278    :param Move: True move generated atom positions to be inside cell
279      False do not move atoms       
280    :return: [[XYZEquiv],Idup,[UijEquiv]]
281
282      *  [XYZEquiv] is list of equivalent positions (XYZ is first entry)
283      *  Idup = [-][C]SS where SS is the symmetry operator number (1-24), C (if not 0,0,0)
284      * is centering operator number (1-4) and - is for inversion
285        Cell = unit cell translations needed to put new positions inside cell
286        [UijEquiv] - equivalent Uij; absent if no Uij given
287       
288    '''
289    XYZEquiv = []
290    UijEquiv = []
291    Idup = []
292    Cell = []
293    X = np.array(XYZ)
294    if Move:
295        X = MoveToUnitCell(X)
296    for ic,cen in enumerate(SGData['SGCen']):
297        C = np.array(cen)
298        for invers in range(int(SGData['SGInv']+1)):
299            for io,[M,T] in enumerate(SGData['SGOps']):
300                idup = ((io+1)+100*ic)*(1-2*invers)
301                XT = np.inner(M,X)+T
302                if len(Uij):
303                    U = Uij2U(Uij)
304                    U = np.inner(M,np.inner(U,M).T)
305                    newUij = U2Uij(U)
306                if invers:
307                    XT = -XT
308                XT += C
309                newX = MoveToUnitCell(XT)
310                cell = np.asarray(np.rint(newX-XT),dtype=np.int32)
311                if All:
312                    if np.allclose(newX,X,atol=0.0002):
313                        idup = False
314                else:
315                    if True in [np.allclose(newX,oldX,atol=0.0002) for oldX in XYZEquiv]:
316                        idup = False
317                if All or idup:
318                    XYZEquiv.append(newX)
319                    Idup.append(idup)
320                    Cell.append(cell)
321                    if len(Uij):
322                        UijEquiv.append(newUij)                   
323    if len(Uij):
324        return zip(XYZEquiv,UijEquiv,Idup,Cell)
325    else:
326        return zip(XYZEquiv,Idup,Cell)
327
328def GenHKLf(HKL,SGData):
329    '''
330    Uses old GSAS Fortran routine genhkl.for
331
332    :param HKL:  [h,k,l]
333    :param SGData: space group data obtained from SpcGroup
334    :returns: iabsnt,mulp,Uniq,phi
335
336     *   iabsnt = True is reflection is forbidden by symmetry
337     *   mulp = reflection multiplicity including Friedel pairs
338     *   Uniq = numpy array of equivalent hkl in descending order of h,k,l
339
340    '''
341    hklf = HKL+[0,]
342    Ops = SGData['SGOps']
343    OpM = np.array([op[0] for op in Ops])
344    OpT = np.array([op[1] for op in Ops])
345    Inv = SGData['SGInv']
346    Cen = np.array([cen for cen in SGData['SGCen']])
347   
348    Nuniq,Uniq,iabsnt,mulp = pyspg.genhklpy(hklf,len(Ops),OpM,OpT,SGData['SGInv'],len(Cen),Cen)
349    h,k,l,f = Uniq
350    Uniq=np.array(zip(h[:Nuniq],k[:Nuniq],l[:Nuniq]))
351    phi = f[:Nuniq]
352   
353    return iabsnt,mulp,Uniq,phi
354                                 
355def GetOprPtrName(key):
356    'Needs a doc string'
357    OprPtrName = {
358        '-6643':[   2,' 1bar ', 1],'6479' :[  10,'  2z  ', 2],'-6479':[   9,'  mz  ', 3],
359        '6481' :[   7,'  my  ', 4],'-6481':[   6,'  2y  ', 5],'6641' :[   4,'  mx  ', 6],
360        '-6641':[   3,'  2x  ', 7],'6591' :[  28,' m+-0 ', 8],'-6591':[  27,' 2+-0 ', 9],
361        '6531' :[  25,' m110 ',10],'-6531':[  24,' 2110 ',11],'6537' :[  61,'  4z  ',12],
362        '-6537':[  62,' -4z  ',13],'975'  :[  68,' 3+++1',14],'6456' :[ 114,'  3z1 ',15],
363        '-489' :[  73,' 3+-- ',16],'483'  :[  78,' 3-+- ',17],'-969' :[  83,' 3--+ ',18],
364        '819'  :[  22,' m+0- ',19],'-819' :[  21,' 2+0- ',20],'2431' :[  16,' m0+- ',21],
365        '-2431':[  15,' 20+- ',22],'-657' :[  19,' m101 ',23],'657'  :[  18,' 2101 ',24],
366        '1943' :[  48,' -4x  ',25],'-1943':[  47,'  4x  ',26],'-2429':[  13,' m011 ',27],
367        '2429' :[  12,' 2011 ',28],'639'  :[  55,' -4y  ',29],'-639' :[  54,'  4y  ',30],
368        '-6484':[ 146,' 2010 ', 4],'6484' :[ 139,' m010 ', 5],'-6668':[ 145,' 2100 ', 6],
369        '6668' :[ 138,' m100 ', 7],'-6454':[ 148,' 2120 ',18],'6454' :[ 141,' m120 ',19],
370        '-6638':[ 149,' 2210 ',20],'6638' :[ 142,' m210 ',21],              #search ends here
371        '2223' :[  68,' 3+++2',39],
372        '6538' :[ 106,'  6z1 ',40],'-2169':[  83,' 3--+2',41],'2151' :[  73,' 3+--2',42],
373        '2205' :[  79,'-3-+-2',43],'-2205':[  78,' 3-+-2',44],'489'  :[  74,'-3+--1',45],
374        '801'  :[  53,'  4y1 ',46],'1945' :[  47,'  4x3 ',47],'-6585':[  62,' -4z3 ',48],
375        '6585' :[  61,'  4z3 ',49],'6584' :[ 114,'  3z2 ',50],'6666' :[ 106,'  6z5 ',51],
376        '6643' :[   1,' Iden ',52],'-801' :[  55,' -4y1 ',53],'-1945':[  48,' -4x3 ',54],
377        '-6666':[ 105,' -6z5 ',55],'-6538':[ 105,' -6z1 ',56],'-2223':[  69,'-3+++2',57],
378        '-975' :[  69,'-3+++1',58],'-6456':[ 113,' -3z1 ',59],'-483' :[  79,'-3-+-1',60],
379        '969'  :[  84,'-3--+1',61],'-6584':[ 113,' -3z2 ',62],'2169' :[  84,'-3--+2',63],
380        '-2151':[  74,'-3+--2',64],'0':[0,' ????',0]
381        }
382    return OprPtrName[key]
383
384def GetKNsym(key):
385    'Needs a doc string'
386    KNsym = {
387        '0'         :'    1   ','1'         :'   -1   ','64'        :'  2(100)','32'        :'  m(100)',
388        '97'        :'2/m(100)','16'        :'  2(010)','8'         :'  m(010)','25'        :'2/m(010)',
389        '2'         :'  2(001)','4'         :'  m(001)','7'         :'2/m(001)','134217728' :'  2(011)',
390        '67108864'  :'  m(011)','201326593' :'2/m(011)','2097152'   :'  2(0+-)','1048576'   :'  m(0+-)',
391        '3145729'   :'2/m(0+-)','8388608'   :'  2(101)','4194304'   :'  m(101)','12582913'  :'2/m(101)',
392        '524288'    :'  2(+0-)','262144'    :'  m(+0-)','796433'    :'2/m(+0-)','1024'      :'  2(110)',
393        '512'       :'  m(110)','1537'      :'2/m(110)','256'       :'  2(+-0)','128'       :'  m(+-0)',
394        '385'       :'2/m(+-0)','76'        :'mm2(100)','52'        :'mm2(010)','42'        :'mm2(001)',
395        '135266336' :'mm2(011)','69206048'  :'mm2(0+-)','8650760'   :'mm2(101)','4718600'   :'mm2(+0-)',
396        '1156'      :'mm2(110)','772'       :'mm2(+-0)','82'        :'  222   ','136314944' :'222(100)',
397        '8912912'   :'222(010)','1282'      :'222(001)','127'       :'  mmm   ','204472417' :'mmm(100)',
398        '13369369'  :'mmm(010)','1927'      :'mmm(001)','33554496'  :'  4(100)','16777280'  :' -4(100)',
399        '50331745'  :'4/m(100)','169869394' :'422(100)','84934738'  :'-42m 100','101711948' :'4mm(100)',
400        '254804095' :'4/mmm100','536870928 ':'  4(010)','268435472' :' -4(010)','805306393' :'4/m (10)',
401        '545783890' :'422(010)','272891986' :'-42m 010','541327412' :'4mm(010)','818675839' :'4/mmm010',
402        '2050'      :'  4(001)','4098'      :' -4(001)','6151'      :'4/m(001)','3410'      :'422(001)',
403        '4818'      :'-42m 001','2730'      :'4mm(001)','8191'      :'4/mmm001','8192'      :'  3(111)',
404        '8193'      :' -3(111)','2629888'   :' 32(111)','1319040'   :' 3m(111)','3940737'   :'-3m(111)',
405        '32768'     :'  3(+--)','32769'     :' -3(+--)','10519552'  :' 32(+--)','5276160'   :' 3m(+--)',
406        '15762945'  :'-3m(+--)','65536'     :'  3(-+-)','65537'     :' -3(-+-)','134808576' :' 32(-+-)',
407        '67437056'  :' 3m(-+-)','202180097' :'-3m(-+-)','131072'    :'  3(--+)','131073'    :' -3(--+)',
408        '142737664' :' 32(--+)','71434368'  :' 3m(--+)','214040961' :'-3m(--+)','237650'    :'   23   ',
409        '237695'    :'   m3   ','715894098' :'   432  ','358068946' :'  -43m  ','1073725439':'   m3m  ',
410        '68157504'  :' mm2d100','4456464'   :' mm2d010','642'       :' mm2d001','153092172' :'-4m2 100',
411        '277348404' :'-4m2 010','5418'      :'-4m2 001','1075726335':'  6/mmm ','1074414420':'-6m2 100',
412        '1075070124':'-6m2 120','1075069650':'   6mm  ','1074414890':'   622  ','1073758215':'   6/m  ',
413        '1073758212':'   -6   ','1073758210':'    6   ','1073759865':'-3m(100)','1075724673':'-3m(120)',
414        '1073758800':' 3m(100)','1075069056':' 3m(120)','1073759272':' 32(100)','1074413824':' 32(120)',
415        '1073758209':'   -3   ','1073758208':'    3   ','1074135143':'mmm(100)','1075314719':'mmm(010)',
416        '1073743751':'mmm(110)','1074004034':' mm2z100','1074790418':' mm2z010','1073742466':' mm2z110',
417        '1074004004':'mm2(100)','1074790412':'mm2(010)','1073742980':'mm2(110)','1073872964':'mm2(120)',
418        '1074266132':'mm2(210)','1073742596':'mm2(+-0)','1073872930':'222(100)','1074266122':'222(010)',
419        '1073743106':'222(110)','1073741831':'2/m(001)','1073741921':'2/m(100)','1073741849':'2/m(010)',
420        '1073743361':'2/m(110)','1074135041':'2/m(120)','1075314689':'2/m(210)','1073742209':'2/m(+-0)',
421        '1073741828':' m(001) ','1073741888':' m(100) ','1073741840':' m(010) ','1073742336':' m(110) ',
422        '1074003968':' m(120) ','1074790400':' m(210) ','1073741952':' m(+-0) ','1073741826':' 2(001) ',
423        '1073741856':' 2(100) ','1073741832':' 2(010) ','1073742848':' 2(110) ','1073872896':' 2(120) ',
424        '1074266112':' 2(210) ','1073742080':' 2(+-0) ','1073741825':'   -1   '
425        }
426    return KNsym[key]       
427
428def GetNXUPQsym(siteSym):       
429    'Needs a doc string'
430    NXUPQsym = {
431        '    1   ':(28,29,28,28),'   -1   ':( 1,29,28, 0),'  2(100)':(12,18,12,25),'  m(100)':(25,18,12,25),
432        '2/m(100)':( 1,18, 0,-1),'  2(010)':(13,17,13,24),'  m(010)':(24,17,13,24),'2/m(010)':( 1,17, 0,-1),
433        '  2(001)':(14,16,14,23),'  m(001)':(23,16,14,23),'2/m(001)':( 1,16, 0,-1),'  2(011)':(10,23,10,22),
434        '  m(011)':(22,23,10,22),'2/m(011)':( 1,23, 0,-1),'  2(0+-)':(11,24,11,21),'  m(0+-)':(21,24,11,21),
435        '2/m(0+-)':( 1,24, 0,-1),'  2(101)':( 8,21, 8,20),'  m(101)':(20,21, 8,20),'2/m(101)':( 1,21, 0,-1),
436        '  2(+0-)':( 9,22, 9,19),'  m(+0-)':(19,22, 9,19),'2/m(+0-)':( 1,22, 0,-1),'  2(110)':( 6,19, 6,18),
437        '  m(110)':(18,19, 6,18),'2/m(110)':( 1,19, 0,-1),'  2(+-0)':( 7,20, 7,17),'  m(+-0)':(17,20, 7,17),
438        '2/m(+-0)':( 1,20, 0,-1),'mm2(100)':(12,10, 0,-1),'mm2(010)':(13,10, 0,-1),'mm2(001)':(14,10, 0,-1),
439        'mm2(011)':(10,13, 0,-1),'mm2(0+-)':(11,13, 0,-1),'mm2(101)':( 8,12, 0,-1),'mm2(+0-)':( 9,12, 0,-1),
440        'mm2(110)':( 6,11, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'  222   ':( 1,10, 0,-1),'222(100)':( 1,13, 0,-1),
441        '222(010)':( 1,12, 0,-1),'222(001)':( 1,11, 0,-1),'  mmm   ':( 1,10, 0,-1),'mmm(100)':( 1,13, 0,-1),
442        'mmm(010)':( 1,12, 0,-1),'mmm(001)':( 1,11, 0,-1),'  4(100)':(12, 4,12, 0),' -4(100)':( 1, 4,12, 0),
443        '4/m(100)':( 1, 4,12,-1),'422(100)':( 1, 4, 0,-1),'-42m 100':( 1, 4, 0,-1),'4mm(100)':(12, 4, 0,-1),
444        '4/mmm100':( 1, 4, 0,-1),'  4(010)':(13, 3,13, 0),' -4(010)':( 1, 3,13, 0),'4/m (10)':( 1, 3,13,-1),
445        '422(010)':( 1, 3, 0,-1),'-42m 010':( 1, 3, 0,-1),'4mm(010)':(13, 3, 0,-1),'4/mmm010':(1, 3, 0,-1,),
446        '  4(001)':(14, 2,14, 0),' -4(001)':( 1, 2,14, 0),'4/m(001)':( 1, 2,14,-1),'422(001)':( 1, 2, 0,-1),
447        '-42m 001':( 1, 2, 0,-1),'4mm(001)':(14, 2, 0,-1),'4/mmm001':( 1, 2, 0,-1),'  3(111)':( 2, 5, 2, 0),
448        ' -3(111)':( 1, 5, 2, 0),' 32(111)':( 1, 5, 0, 2),' 3m(111)':( 2, 5, 0, 2),'-3m(111)':( 1, 5, 0,-1),
449        '  3(+--)':( 5, 8, 5, 0),' -3(+--)':( 1, 8, 5, 0),' 32(+--)':( 1, 8, 0, 5),' 3m(+--)':( 5, 8, 0, 5),
450        '-3m(+--)':( 1, 8, 0,-1),'  3(-+-)':( 4, 7, 4, 0),' -3(-+-)':( 1, 7, 4, 0),' 32(-+-)':( 1, 7, 0, 4),
451        ' 3m(-+-)':( 4, 7, 0, 4),'-3m(-+-)':( 1, 7, 0,-1),'  3(--+)':( 3, 6, 3, 0),' -3(--+)':( 1, 6, 3, 0),
452        ' 32(--+)':( 1, 6, 0, 3),' 3m(--+)':( 3, 6, 0, 3),'-3m(--+)':( 1, 6, 0,-1),'   23   ':( 1, 1, 0, 0),
453        '   m3   ':( 1, 1, 0, 0),'   432  ':( 1, 1, 0, 0),'  -43m  ':( 1, 1, 0, 0),'   m3m  ':( 1, 1, 0, 0),
454        ' mm2d100':(12,13, 0,-1),' mm2d010':(13,12, 0,-1),' mm2d001':(14,11, 0,-1),'-4m2 100':( 1, 4, 0,-1),
455        '-4m2 010':( 1, 3, 0,-1),'-4m2 001':( 1, 2, 0,-1),'  6/mmm ':( 1, 9, 0,-1),'-6m2 100':( 1, 9, 0,-1),
456        '-6m2 120':( 1, 9, 0,-1),'   6mm  ':(14, 9, 0,-1),'   622  ':( 1, 9, 0,-1),'   6/m  ':( 1, 9,14,-1),
457        '   -6   ':( 1, 9,14, 0),'    6   ':(14, 9,14, 0),'-3m(100)':( 1, 9, 0,-1),'-3m(120)':( 1, 9, 0,-1),
458        ' 3m(100)':(14, 9, 0,14),' 3m(120)':(14, 9, 0,14),' 32(100)':( 1, 9, 0,14),' 32(120)':( 1, 9, 0,14),
459        '   -3   ':( 1, 9,14, 0),'    3   ':(14, 9,14, 0),'mmm(100)':( 1,14, 0,-1),'mmm(010)':( 1,15, 0,-1),
460        'mmm(110)':( 1,11, 0,-1),' mm2z100':(14,14, 0,-1),' mm2z010':(14,15, 0,-1),' mm2z110':(14,11, 0,-1),
461        'mm2(100)':(12,14, 0,-1),'mm2(010)':(13,15, 0,-1),'mm2(110)':( 6,11, 0,-1),'mm2(120)':(15,14, 0,-1),
462        'mm2(210)':(16,15, 0,-1),'mm2(+-0)':( 7,11, 0,-1),'222(100)':( 1,14, 0,-1),'222(010)':( 1,15, 0,-1),
463        '222(110)':( 1,11, 0,-1),'2/m(001)':( 1,16,14,-1),'2/m(100)':( 1,25,12,-1),'2/m(010)':( 1,28,13,-1),
464        '2/m(110)':( 1,19, 6,-1),'2/m(120)':( 1,27,15,-1),'2/m(210)':( 1,26,16,-1),'2/m(+-0)':( 1,20,17,-1),
465        ' m(001) ':(23,16,14,23),' m(100) ':(26,25,12,26),' m(010) ':(27,28,13,27),' m(110) ':(18,19, 6,18),
466        ' m(120) ':(24,27,15,24),' m(210) ':(25,26,16,25),' m(+-0) ':(17,20, 7,17),' 2(001) ':(14,16,14,23),
467        ' 2(100) ':(12,25,12,26),' 2(010) ':(13,28,13,27),' 2(110) ':( 6,19, 6,18),' 2(120) ':(15,27,15,24),
468        ' 2(210) ':(16,26,16,25),' 2(+-0) ':( 7,20, 7,17),'   -1   ':( 1,29,28, 0)
469        }
470    return NXUPQsym[siteSym]
471
472def GetCSxinel(siteSym): 
473    'Needs a doc string'
474    CSxinel = [[],                         # 0th empty - indices are Fortran style
475        [[0,0,0],[ 0.0, 0.0, 0.0]],      #  0  0  0
476        [[1,1,1],[ 1.0, 1.0, 1.0]],      #  X  X  X
477        [[1,1,1],[ 1.0, 1.0,-1.0]],      #  X  X -X
478        [[1,1,1],[ 1.0,-1.0, 1.0]],      #  X -X  X
479        [[1,1,1],[ 1.0,-1.0,-1.0]],      # -X  X  X
480        [[1,1,0],[ 1.0, 1.0, 0.0]],      #  X  X  0
481        [[1,1,0],[ 1.0,-1.0, 0.0]],      #  X -X  0
482        [[1,0,1],[ 1.0, 0.0, 1.0]],      #  X  0  X
483        [[1,0,1],[ 1.0, 0.0,-1.0]],      #  X  0 -X
484        [[0,1,1],[ 0.0, 1.0, 1.0]],      #  0  Y  Y
485        [[0,1,1],[ 0.0, 1.0,-1.0]],      #  0  Y -Y
486        [[1,0,0],[ 1.0, 0.0, 0.0]],      #  X  0  0
487        [[0,1,0],[ 0.0, 1.0, 0.0]],      #  0  Y  0
488        [[0,0,1],[ 0.0, 0.0, 1.0]],      #  0  0  Z
489        [[1,1,0],[ 1.0, 2.0, 0.0]],      #  X 2X  0
490        [[1,1,0],[ 2.0, 1.0, 0.0]],      # 2X  X  0
491        [[1,1,2],[ 1.0, 1.0, 1.0]],      #  X  X  Z
492        [[1,1,2],[ 1.0,-1.0, 1.0]],      #  X -X  Z
493        [[1,2,1],[ 1.0, 1.0, 1.0]],      #  X  Y  X
494        [[1,2,1],[ 1.0, 1.0,-1.0]],      #  X  Y -X
495        [[1,2,2],[ 1.0, 1.0, 1.0]],      #  X  Y  Y
496        [[1,2,2],[ 1.0, 1.0,-1.0]],      #  X  Y -Y
497        [[1,2,0],[ 1.0, 1.0, 0.0]],      #  X  Y  0
498        [[1,0,2],[ 1.0, 0.0, 1.0]],      #  X  0  Z
499        [[0,1,2],[ 0.0, 1.0, 1.0]],      #  0  Y  Z
500        [[1,1,2],[ 1.0, 2.0, 1.0]],      #  X 2X  Z
501        [[1,1,2],[ 2.0, 1.0, 1.0]],      # 2X  X  Z
502        [[1,2,3],[ 1.0, 1.0, 1.0]],      #  X  Y  Z
503        ]
504    indx = GetNXUPQsym(siteSym)
505    return CSxinel[indx[0]]
506   
507def GetCSuinel(siteSym):
508    "returns Uij terms, multipliers, GUI flags & Uiso2Uij multipliers"
509    CSuinel = [[],                                             # 0th empty - indices are Fortran style
510        [[1,1,1,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,0,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  0  0  0
511        [[1,1,2,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,0,1,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  0  0  0
512        [[1,2,1,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  0  0  0
513        [[1,2,2,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,0,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  0  0  0
514        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D  D  D
515        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0,-1.0,-1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D -D -D
516        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0,-1.0, 1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D -D  D
517        [[1,1,1,2,2,2],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,0,0,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  A  D  D -D
518        [[1,1,2,1,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,0,1,0,0,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  A  C A/2 0  0
519        [[1,2,3,0,0,0],[ 1.0, 1.0, 1.0, 0.0, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  0  0
520        [[1,1,2,3,0,0],[ 1.0, 1.0, 1.0, 1.0, 0.0, 0.0],[1,0,1,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  0  0
521        [[1,2,1,0,3,0],[ 1.0, 1.0, 1.0, 0.0, 1.0, 0.0],[1,1,0,0,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  0  E  0
522        [[1,2,2,0,0,3],[ 1.0, 1.0, 1.0, 0.0, 0.0, 1.0],[1,1,0,0,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  0  0  F
523        [[1,2,3,2,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.5,0.0]],    #  A  B  C B/2 0  0
524        [[1,2,3,1,0,0],[ 1.0, 1.0, 1.0, 0.5, 0.0, 0.0],[1,1,1,0,0,0],[1.0,1.0,1.0,0.0,0.5,0.0]],    #  A  B  C A/2 0  0
525        [[1,2,3,4,0,0],[ 1.0, 1.0, 1.0, 1.0, 0.0, 0.0],[1,1,1,1,0,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  D  0  0
526        [[1,2,3,0,4,0],[ 1.0, 1.0, 1.0, 0.0, 1.0, 0.0],[1,1,1,0,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  E  0
527        [[1,2,3,0,0,4],[ 1.0, 1.0, 1.0, 0.0, 0.0, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  0  0  F
528        [[1,1,2,3,4,4],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,0,1,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  E -E
529        [[1,1,2,3,4,4],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,0,1,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  A  C  D  E  E
530        [[1,2,1,3,4,3],[ 1.0, 1.0, 1.0, 1.0, 1.0,-1.0],[1,1,0,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  D  E -D
531        [[1,2,1,3,4,3],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,0,1,1,0],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  A  D  E  D
532        [[1,2,2,3,3,4],[ 1.0, 1.0, 1.0, 1.0,-1.0, 1.0],[1,1,0,1,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  D -D  F
533        [[1,2,2,3,3,4],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,0,1,0,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  B  D  D  F
534        [[1,2,3,2,4,4],[ 1.0, 1.0, 1.0, 0.5, 0.5, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C B/2 F/2 F
535        [[1,2,3,1,0,4],[ 1.0, 1.0, 1.0, 0.5, 0.0, 1.0],[1,1,1,0,0,1],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C A/2  0  F
536        [[1,2,3,2,4,0],[ 1.0, 1.0, 1.0, 0.5, 1.0, 0.0],[1,1,1,0,1,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C B/2  E  0
537        [[1,2,3,1,4,4],[ 1.0, 1.0, 1.0, 0.5, 1.0, 0.5],[1,1,1,0,1,0],[1.0,1.0,1.0,0.5,0.0,0.0]],    #  A  B  C A/2  E E/2
538        [[1,2,3,4,5,6],[ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0],[1,1,1,1,1,1],[1.0,1.0,1.0,0.0,0.0,0.0]],    #  A  B  C  D  E   F
539        ]
540    indx = GetNXUPQsym(siteSym)
541    return CSuinel[indx[1]]
542   
543def MustrainNames(SGData):
544    'Needs a doc string'
545    laue = SGData['SGLaue']
546    uniq = SGData['SGUniq']
547    if laue in ['m3','m3m']:
548        return ['S400','S220']
549    elif laue in ['6/m','6/mmm','3m1']:
550        return ['S400','S004','S202']
551    elif laue in ['31m','3']:
552        return ['S400','S004','S202','S211']
553    elif laue in ['3R','3mR']:
554        return ['S400','S220','S310','S211']
555    elif laue in ['4/m','4/mmm']:
556        return ['S400','S004','S220','S022']
557    elif laue in ['mmm']:
558        return ['S400','S040','S004','S220','S202','S022']
559    elif laue in ['2/m']:
560        SHKL = ['S400','S040','S004','S220','S202','S022']
561        if uniq == 'a':
562            SHKL += ['S013','S031','S211']
563        elif uniq == 'b':
564            SHKL += ['S301','S103','S121']
565        elif uniq == 'c':
566            SHKL += ['S130','S310','S112']
567        return SHKL
568    else:
569        SHKL = ['S400','S040','S004','S220','S202','S022']
570        SHKL += ['S310','S103','S031','S130','S301','S013']
571        SHKL += ['S211','S121','S112']
572        return SHKL
573
574def HStrainNames(SGData):
575    'Needs a doc string'
576    laue = SGData['SGLaue']
577    uniq = SGData['SGUniq']
578    if laue in ['m3','m3m']:
579        return ['D11','eA']         #add cubic strain term
580    elif laue in ['6/m','6/mmm','3m1','31m','3']:
581        return ['D11','D33']
582    elif laue in ['3R','3mR']:
583        return ['D11','D12']
584    elif laue in ['4/m','4/mmm']:
585        return ['D11','D33']
586    elif laue in ['mmm']:
587        return ['D11','D22','D33']
588    elif laue in ['2/m']:
589        Dij = ['D11','D22','D33']
590        if uniq == 'a':
591            Dij += ['D23']
592        elif uniq == 'b':
593            Dij += ['D13']
594        elif uniq == 'c':
595            Dij += ['D12']
596        return Dij
597    else:
598        Dij = ['D11','D22','D33','D12','D13','D23']
599        return Dij
600   
601def MustrainCoeff(HKL,SGData):
602    'Needs a doc string'
603    #NB: order of terms is the same as returned by MustrainNames
604    laue = SGData['SGLaue']
605    uniq = SGData['SGUniq']
606    h,k,l = HKL
607    Strm = []
608    if laue in ['m3','m3m']:
609        Strm.append(h**4+k**4+l**4)
610        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
611    elif laue in ['6/m','6/mmm','3m1']:
612        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
613        Strm.append(l**4)
614        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
615    elif laue in ['31m','3']:
616        Strm.append(h**4+k**4+2.0*k*h**3+2.0*h*k**3+3.0*(h*k)**2)
617        Strm.append(l**4)
618        Strm.append(3.0*((h*l)**2+(k*l)**2+h*k*l**2))
619        Strm.append(4.0*h*k*l*(h+k))
620    elif laue in ['3R','3mR']:
621        Strm.append(h**4+k**4+l**4)
622        Strm.append(3.0*((h*k)**2+(h*l)**2+(k*l)**2))
623        Strm.append(2.0*(h*l**3+l*k**3+k*h**3)+2.0*(l*h**3+k*l**3+l*k**3))
624        Strm.append(4.0*(k*l*h**2+h*l*k**2+h*k*l**2))
625    elif laue in ['4/m','4/mmm']:
626        Strm.append(h**4+k**4)
627        Strm.append(l**4)
628        Strm.append(3.0*(h*k)**2)
629        Strm.append(3.0*((h*l)**2+(k*l)**2))
630    elif laue in ['mmm']:
631        Strm.append(h**4)
632        Strm.append(k**4)
633        Strm.append(l**4)
634        Strm.append(3.0*(h*k)**2)
635        Strm.append(3.0*(h*l)**2)
636        Strm.append(3.0*(k*l)**2)
637    elif laue in ['2/m']:
638        Strm.append(h**4)
639        Strm.append(k**4)
640        Strm.append(l**4)
641        Strm.append(3.0*(h*k)**2)
642        Strm.append(3.0*(h*l)**2)
643        Strm.append(3.0*(k*l)**2)
644        if uniq == 'a':
645            Strm.append(2.0*k*l**3)
646            Strm.append(2.0*l*k**3)
647            Strm.append(4.0*k*l*h**2)
648        elif uniq == 'b':
649            Strm.append(2.0*l*h**3)
650            Strm.append(2.0*h*l**3)
651            Strm.append(4.0*h*l*k**2)
652        elif uniq == 'c':
653            Strm.append(2.0*h*k**3)
654            Strm.append(2.0*k*h**3)
655            Strm.append(4.0*h*k*l**2)
656    else:
657        Strm.append(h**4)
658        Strm.append(k**4)
659        Strm.append(l**4)
660        Strm.append(3.0*(h*k)**2)
661        Strm.append(3.0*(h*l)**2)
662        Strm.append(3.0*(k*l)**2)
663        Strm.append(2.0*k*h**3)
664        Strm.append(2.0*h*l**3)
665        Strm.append(2.0*l*k**3)
666        Strm.append(2.0*h*k**3)
667        Strm.append(2.0*l*h**3)
668        Strm.append(2.0*k*l**3)
669        Strm.append(4.0*k*l*h**2)
670        Strm.append(4.0*h*l*k**2)
671        Strm.append(4.0*k*h*l**2)
672    return Strm
673   
674def Muiso2Shkl(muiso,SGData,cell):
675    "this is to convert isotropic mustrain to generalized Shkls - doesn't work just now"
676    import GSASIIlattice as G2lat
677    from scipy.optimize import fmin
678    A = G2lat.cell2AB(cell)[0]
679    def minMus(Shkl,H,muiso,SGData,A):
680        U = np.inner(A.T,H)
681        S = np.array(MustrainCoeff(H.T,SGData))
682        sum = np.sqrt(np.sum(np.multiply(S,Shkl)))
683        return abs(muiso-sum*H)
684    laue = SGData['SGLaue']
685    if laue in ['m3','m3m']:
686        H = [[1,0,0],[1,1,0]]
687        S0 = [0.01,0.01]
688    elif laue in ['6/m','6/mmm','3m1']:
689        H = [[1,0,0],[0,0,1],[1,0,1]]
690        S0 = [0.01,0.01,0.01]
691    elif laue in ['31m','3']:
692        H = [[1,0,0],[0,0,1],[1,0,1],[1,1,1]]
693        S0 = [0.01,0.01,0.01,0.01]
694    elif laue in ['3R','3mR']:
695        H = [[1,0,0],[1,1,0],[1,0,1],[1,1,1]]
696        S0 = [0.01,0.01,0.01,0.01]
697    elif laue in ['4/m','4/mmm']:
698        H = [[1,0,0],[0,0,1],[1,1,0],[1,0,1]]
699        S0 = [0.01,0.01,0.01,0.01]
700    elif laue in ['mmm']:
701        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1]]
702        S0 = [0.01,0.01,0.01,0.01,0.01,0.01]
703    elif laue in ['2/m']:
704        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1]]
705        if uniq == 'a':
706            H.append([0,1,-1])
707            H.append([0,-2,1])
708        elif uniq == 'b':
709            H.append([1,0,-1])
710            H.append([-2,0,1])
711        elif uniq == 'c':
712            H.append([1,-1,0])
713            H.append([-2,1,0])
714        H.append([1,1,1])
715        S0 = [9*[0.01,]]
716    else:
717        H = [[1,0,0],[0,1,0],[0,0,1],[1,1,0],[1,0,1],[0,1,1],
718            [-1,1,0],[1,0,-1],[0,-1,1],[1,-2,0],[-2,0,1],[0,1,-2],
719            [1,-1,1],[-1, 1, 1],[1,-1,1]]
720        S0 = [15*[0.01,]]
721    H = np.array(H)
722    S0 = np.array(S0)
723    return fmin(minMus,S0,(H,muiso,SGData,A))
724       
725def SytSym(XYZ,SGData):
726    '''
727    Generates the number of equivalent positions and a site symmetry code for a specified coordinate and space group
728
729    :param XYZ: an array, tuple or list containing 3 elements: x, y & z
730    :param SGData: from SpcGroup
731    :Returns: a two element tuple:
732
733     * The 1st element is a code for the site symmetry (see GetKNsym)
734     * The 2nd element is the site multiplicity
735
736    '''
737    def PackRot(SGOps):
738        IRT = []
739        for ops in SGOps:
740            M = ops[0]
741            irt = 0
742            for j in range(2,-1,-1):
743                for k in range(2,-1,-1):
744                    irt *= 3
745                    irt += M[k][j]
746            IRT.append(int(irt))
747        return IRT
748       
749    SymName = ''
750    Mult = 1
751    Isym = 0
752    if SGData['SGLaue'] in ['3','3m1','31m','6/m','6/mmm']:
753        Isym = 1073741824
754    Jdup = 0
755    Xeqv = GenAtom(XYZ,SGData,True)
756    IRT = PackRot(SGData['SGOps'])
757    L = -1
758    for ic,cen in enumerate(SGData['SGCen']):
759        for invers in range(int(SGData['SGInv']+1)):
760            for io,ops in enumerate(SGData['SGOps']):
761                irtx = (1-2*invers)*IRT[io]
762                L += 1
763                if not Xeqv[L][1]:
764                    Jdup += 1
765                    jx = GetOprPtrName(str(irtx))
766                    if jx[2] < 39:
767                        Isym += 2**(jx[2]-1)
768    if Isym == 1073741824: Isym = 0
769    Mult = len(SGData['SGOps'])*len(SGData['SGCen'])*(int(SGData['SGInv'])+1)/Jdup
770         
771    return GetKNsym(str(Isym)),Mult
772   
773def ElemPosition(SGData):
774    ''' Under development.
775    Object here is to return a list of symmetry element types and locations suitable
776    for say drawing them.
777    So far I have the element type... getting all possible locations without lookup may be impossible!
778    '''
779    SymElements = []
780    Inv = SGData['SGInv']
781    Cen = SGData['SGCen']
782    eleSym = {-3:['','-1'],-2:['',-6],-1:['2','-4'],0:['3','-3'],1:['4','m'],2:['6',''],3:['1','']}
783    # get operators & expand if centrosymmetric
784    Ops = SGData['SGOps']
785    opM = np.array([op[0].T for op in Ops])
786    opT = np.array([op[1] for op in Ops])
787    if Inv:
788        opM = np.concatenate((opM,-opM))
789        opT = np.concatenate((opT,-opT))
790    opMT = zip(opM,opT)
791    for M,T in opMT[1:]:        #skip I
792        Dt = int(nl.det(M))
793        Tr = int(np.trace(M))
794        Dt = -(Dt-1)/2
795        Es = eleSym[Tr][Dt]
796        if Dt:              #rotation-inversion
797            I = np.eye(3)
798            if Tr == 1:     #mirrors/glides
799                if np.any(T):       #glide
800                    M2 = np.inner(M,M)
801                    MT = np.inner(M,T)+T
802                    print 'glide',Es,MT
803                    print M2
804                else:               #mirror
805                    print 'mirror',Es,T
806                    print I-M
807                X = [-1,-1,-1]
808            elif Tr == -3:  # pure inversion
809                X = np.inner(nl.inv(I-M),T)
810                print 'inversion',Es,X
811            else:           #other rotation-inversion
812                M2 = np.inner(M,M)
813                MT = np.inner(M,T)+T
814                print 'rot-inv',Es,MT
815                print M2
816                X = [-1,-1,-1]
817               
818           
819           
820        else:               #rotations
821            print 'rotation',Es
822            X = [-1,-1,-1]
823        #SymElements.append([Es,X])
824       
825    return #SymElements
826   
827def ApplyStringOps(A,SGData,X,Uij=[]):
828    'Needs a doc string'
829    SGOps = SGData['SGOps']
830    SGCen = SGData['SGCen']
831    Ax = A.split('+')
832    Ax[0] = int(Ax[0])
833    iC = 0
834    if Ax[0] < 0:
835        iC = 1
836    Ax[0] = abs(Ax[0])
837    nA = Ax[0]%100-1
838    cA = Ax[0]/100
839    Cen = SGCen[cA]
840    M,T = SGOps[nA]
841    if len(Ax)>1:
842        cellA = Ax[1].split(',')
843        cellA = np.array([int(a) for a in cellA])
844    else:
845        cellA = np.zeros(3)
846    newX = (1-2*iC)*(Cen+np.inner(M,X)+T)+cellA
847    if len(Uij):
848        U = Uij2U(Uij)
849        U = np.inner(M,np.inner(U,M).T)
850        newUij = U2Uij(U)
851        return [newX,newUij]
852    else:
853        return newX
854       
855def StringOpsProd(A,B,SGData):
856    """
857    Find A*B where A & B are in strings '-' + '100*c+n' + '+ijk'
858    where '-' indicates inversion, c(>0) is the cell centering operator,
859    n is operator number from SgOps and ijk are unit cell translations (each may be <0).
860    Should return resultant string - C. SGData - dictionary using entries:
861
862       *  'SGCen': cell centering vectors [0,0,0] at least
863       *  'SGOps': symmetry operations as [M,T] so that M*x+T = x'
864
865    """
866    SGOps = SGData['SGOps']
867    SGCen = SGData['SGCen']
868    #1st split out the cell translation part & work on the operator parts
869    Ax = A.split('+'); Bx = B.split('+')
870    Ax[0] = int(Ax[0]); Bx[0] = int(Bx[0])
871    iC = 0
872    if Ax[0]*Bx[0] < 0:
873        iC = 1
874    Ax[0] = abs(Ax[0]); Bx[0] = abs(Bx[0])
875    nA = Ax[0]%100-1;  nB = Bx[0]%100-1
876    cA = Ax[0]/100;  cB = Bx[0]/100
877    Cen = (SGCen[cA]+SGCen[cB])%1.0
878    cC = np.nonzero([np.allclose(C,Cen) for C in SGCen])[0][0]
879    Ma,Ta = SGOps[nA]; Mb,Tb = SGOps[nB]
880    Mc = np.inner(Ma,Mb.T)
881#    print Ma,Mb,Mc
882    Tc = (np.add(np.inner(Mb,Ta)+1.,Tb))%1.0
883#    print Ta,Tb,Tc
884#    print [np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps]
885    nC = np.nonzero([np.allclose(M,Mc)&np.allclose(T,Tc) for M,T in SGOps])[0][0]
886    #now the cell translation part
887    if len(Ax)>1:
888        cellA = Ax[1].split(',')
889        cellA = [int(a) for a in cellA]
890    else:
891        cellA = [0,0,0]
892    if len(Bx)>1:
893        cellB = Bx[1].split(',')
894        cellB = [int(b) for b in cellB]
895    else:
896        cellB = [0,0,0]
897    cellC = np.add(cellA,cellB)
898    C = str(((nC+1)+(100*cC))*(1-2*iC))+'+'+ \
899        str(int(cellC[0]))+','+str(int(cellC[1]))+','+str(int(cellC[2]))
900    return C
901           
902def U2Uij(U):
903    #returns the UIJ vector U11,U22,U33,U12,U13,U23 from tensor U
904    return [U[0][0],U[1][1],U[2][2],U[0][1],U[0][2],U[1][2]]
905   
906def Uij2U(Uij):
907    #returns the thermal motion tensor U from Uij as numpy array
908    return np.array([[Uij[0],Uij[3],Uij[4]],[Uij[3],Uij[1],Uij[5]],[Uij[4],Uij[5],Uij[2]]])
909
910def StandardizeSpcName(spcgroup):
911    '''Accept a spacegroup name where spaces may have not been used
912    in the names according to the GSAS convention (spaces between symmetry
913    for each axis) and return the space group name as used in GSAS
914    '''
915    rspc = spcgroup.replace(' ','').upper()
916    # deal with rhombohedral and hexagonal setting designations
917    rhomb = ''
918    if rspc[-1:] == 'R':
919        rspc = rspc[:-1]
920        rhomb = ' R'
921    if rspc[-1:] == 'H': # hexagonal is assumed and thus can be ignored
922        rspc = rspc[:-1]
923    # look for a match in the spacegroup lists
924    for i in spglist.values():
925        for spc in i:
926            if rspc == spc.replace(' ','').upper():
927                return spc + rhomb
928    # how about the post-2002 orthorhombic names?
929    for i,spc in sgequiv_2002_orthorhombic:
930        if rspc == i.replace(' ','').upper():
931            return spc
932    # not found
933    return ''
934
935   
936'''A dictionary of space groups as ordered and named in the pre-2002 International
937Tables Volume A, except that spaces are used following the GSAS convention to
938separate the different crystallographic directions.
939Note that the symmetry codes here will recognize many non-standard space group
940symbols with different settings. They are ordered by Laue group
941'''
942spglist = {
943    'P1' : ('P 1','P -1',), # 1-2
944    'P2/m': ('P 2','P 21','P m','P a','P c','P n',
945        'P 2/m','P 21/m','P 2/c','P 2/a','P 2/n','P 21/c','P 21/a','P 21/n',), #3-15
946    'C2/m':('C 2','C m','C c','C n',
947        'C 2/m','C 2/c','C 2/n',),
948    'Pmmm':('P 2 2 2',
949        'P 2 2 21','P 2 21 2','P 21 2 2',
950        'P 21 21 2','P 21 2 21','P 2 21 21',
951        'P 21 21 21',
952        'P m m 2','P m 2 m','P 2 m m',
953        'P m c 21','P c m 21','P 21 m a','P 21 a m','P b 21 m','P m 21 b',
954        'P c c 2','P 2 a a','P b 2 b',
955        'P m a 2','P b m 2','P 2 m b','P 2 c m','P c 2 m','P m 2 a',
956        'P c a 21','P b c 21','P 21 a b','P 21 c a','P c 21 b','P b 21 a',
957        'P n c 2','P c n 2','P 2 n a','P 2 a n','P b 2 n','P n 2 b',
958        'P m n 21','P n m 21','P 21 m n','P 21 n m','P n 21 m','P m 21 n',
959        'P b a 2','P 2 c b','P c 2 a',
960        'P n a 21','P b n 21','P 21 n b','P 21 c n','P c 21 n','P n 21 a',
961        'P n n 2','P 2 n n','P n 2 n',
962        'P m m m','P n n n',
963        'P c c m','P m a a','P b m b',
964        'P b a n','P n c b','P c n a',
965        'P m m a','P m m b','P b m m','P c m m','P m c m','P m a m',
966        'P n n a','P n n b','P b n n','P c n n','P n c n','P n a n',
967        'P m n a','P n m b','P b m n','P c n m','P n c m','P m a n',
968        'P c c a','P c c b','P b a a','P c a a','P b c b','P b a b',
969        'P b a m','P m c b','P c m a',
970        'P c c n','P n a a','P b n b',
971        'P b c m','P c a m','P m c a','P m a b','P b m a','P c m b',
972        'P n n m','P m n n','P n m n',
973        'P m m n','P n m m','P m n m',
974        'P b c n','P c a n','P n c a','P n a b','P b n a','P c n b',
975        'P b c a','P c a b',
976        'P n m a','P m n b','P b n m','P c m n','P m c n','P n a m',
977        ),
978    'Cmmm':('C 2 2 21','C 2 2 2','C m m 2','C m c 21','C c c 2','C m 2 m','C 2 m m',
979        'C m 2 a','C 2 m b','C 2 c m','C c 2 m','C 2 c m','C c 2 m', # check: C c 2 m & C c 2 m twice
980        'C m c a','C m m m','C c c m','C m m a','C c c a','C m c m',),
981    'Immm':('I 2 2 2','I 21 21 21','I m m m',
982        'I m m 2','I m 2 m','I 2 m m',
983        'I b a 2','I 2 c b','I c 2 a',
984        'I m a 2','I b m 2','I 2 m b','I 2 c m','I c 2 m','I m 2 a',
985        'I b a m','I m c b','I c m a',
986        'I b c a','I c a b',
987        'I m m a','I m m b','I b m m ','I c m m','I m c m','I m a m',),
988    'Fmmm':('F 2 2 2','F m m m', 'F d d d',
989        'F m m 2','F m 2 m','F 2 m m',
990        'F d d 2','F d 2 d','F 2 d d',),
991    'P4/mmm':('P 4','P 41','P 42','P 43','P -4','P 4/m','P 42/m','P 4/n','P 42/n',
992        'P 4 2 2','P 4 21 2','P 41 2 2','P 41 21 2','P 42 2 2',
993        'P 42 21 2','P 43 2 2','P 43 21 2','P 4 m m','P 4 b m','P 42 c m',
994        'P 42 n m','P 4 c c','P 4 n c','P 42 m c','P 42 b c','P -4 2 m',
995        'P -4 2 c','P -4 21 m','P -4 21 c','P -4 m 2','P -4 c 2','P -4 b 2',
996        'P -4 n 2','P 4/m m m','P 4/m c c','P 4/n b m','P 4/n n c','P 4/m b m',
997        'P 4/m n c','P 4/n m m','P 4/n c c','P 42/m m c','P 42/m c m',
998        'P 42/n b c','P 42/n n m','P 42/m b c','P 42/m n m','P 42/n m c',
999        'P 42/n c m',),
1000    'I4/mmm':('I 4','I 41','I -4','I 4/m','I 41/a','I 4 2 2','I 41 2 2','I 4 m m',
1001        'I 4 c m','I 41 m d','I 41 c d',
1002        'I -4 m 2','I -4 c 2','I -4 2 m','I -4 2 d','I 4/m m m','I 4/m c m',
1003        'I 41/a m d','I 41/a c d'),
1004    'R3-H':('R 3','R -3','R 3 2','R 3 m','R 3 c','R -3 m','R -3 c',),
1005    'P6/mmm': ('P 3','P 31','P 32','P -3','P 3 1 2','P 3 2 1','P 31 1 2',
1006        'P 31 2 1','P 32 1 2','P 32 2 1', 'P 3 m 1','P 3 1 m','P 3 c 1',
1007        'P 3 1 c','P -3 1 m','P -3 1 c','P -3 m 1','P -3 c 1','P 6','P 61',
1008        'P 65','P 62','P 64','P 63','P -6','P 6/m','P 63/m','P 6 2 2',
1009        'P 61 2 2','P 65 2 2','P 62 2 2','P 64 2 2','P 63 2 2','P 6 m m',
1010        'P 6 c c','P 63 c m','P 63 m c','P -6 m 2','P -6 c 2','P -6 2 m',
1011        'P -6 2 c','P 6/m m m','P 6/m c c','P 63/m c m','P 63/m m c',),
1012    'Pm3m': ('P 2 3','P 21 3','P m 3','P n 3','P a 3','P 4 3 2','P 42 3 2',
1013        'P 43 3 2','P 41 3 2','P -4 3 m','P -4 3 n','P m 3 m','P n 3 n',
1014        'P m 3 n','P n 3 m',),
1015    'Im3m':('I 2 3','I 21 3','I m -3','I a -3', 'I 4 3 2','I 41 3 2',
1016        'I -4 3 m', 'I -4 3 d','I m -3 m','I m 3 m','I a -3 d',),
1017    'Fm3m':('F 2 3','F m -3','F d -3','F 4 3 2','F 41 3 2','F -4 3 m',
1018        'F -4 3 c','F m -3 m','F m 3 m','F m -3 c','F d -3 m','F d -3 c',),
1019}
1020
1021#'A few non-standard space groups for test use'
1022nonstandard_sglist = ('P 21 1 1','P 1 21 1','P 1 1 21','R 3 r','R 3 2 h', 
1023                      'R -3 r', 'R 3 2 r','R 3 m h', 'R 3 m r',
1024                      'R 3 c r','R -3 c r','R -3 m r',),
1025
1026#A list of orthorhombic space groups that were renamed in the 2002 Volume A,
1027# along with the pre-2002 name. The e designates a double glide-plane'''
1028sgequiv_2002_orthorhombic= (('A e m 2', 'A b m 2',),
1029                            ('A e a 2', 'A b a 2',),
1030                            ('C m c e', 'C m c a',),
1031                            ('C m m e', 'C m m a',),
1032                            ('C c c e', 'C c c a'),)
1033#Use the space groups types in this order to list the symbols in the
1034#order they are listed in the International Tables, vol. A'''
1035symtypelist = ('triclinic', 'monoclinic', 'orthorhombic', 'tetragonal', 
1036               'trigonal', 'hexagonal', 'cubic')
1037
1038# self-test materials follow. Requires files in directory testinp
1039selftestlist = []
1040'''Defines a list of self-tests'''
1041selftestquiet = True
1042def _ReportTest():
1043    'Report name and doc string of current routine when ``selftestquiet`` is False'
1044    if not selftestquiet:
1045        import inspect
1046        caller = inspect.stack()[1][3]
1047        doc = eval(caller).__doc__
1048        if doc is not None:
1049            print('testing '+__file__+' with '+caller+' ('+doc+')')
1050        else:
1051            print('testing '+__file__()+" with "+caller)
1052def test0():
1053    '''self-test #0: exercise MoveToUnitCell'''
1054    _ReportTest()
1055    msg = "MoveToUnitCell failed"
1056    assert (MoveToUnitCell([1,2,3]) == [0,0,0]).all, msg
1057    assert (MoveToUnitCell([2,-1,-2]) == [0,0,0]).all, msg
1058    assert abs(MoveToUnitCell(np.array([-.1]))[0]-0.9) < 1e-6, msg
1059    assert abs(MoveToUnitCell(np.array([.1]))[0]-0.1) < 1e-6, msg
1060selftestlist.append(test0)
1061
1062def test1():
1063    '''self-test #1: SpcGroup and SGPrint against previous results'''
1064    _ReportTest()
1065    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
1066    if ospath.exists(testdir):
1067        if testdir not in sys.path: sys.path.insert(0,testdir)
1068    import spctestinp
1069    def CompareSpcGroup(spc, referr, refdict, reflist): 
1070        'Compare output from GSASIIspc.SpcGroup with results from a previous run'
1071        # if an error is reported, the dictionary can be ignored
1072        msg0 = "CompareSpcGroup failed on space group %s" % spc
1073        result = SpcGroup(spc)
1074        if result[0] == referr and referr > 0: return True
1075        keys = result[1].keys()
1076        #print result[1]['SpGrp']
1077        msg = msg0 + " in list lengths"
1078        assert len(keys) == len(refdict.keys()), msg
1079        for key in keys:
1080            if key == 'SGOps' or  key == 'SGCen':
1081                msg = msg0 + (" in key %s length" % key)
1082                assert len(refdict[key]) == len(result[1][key]), msg
1083                for i in range(len(refdict[key])):
1084                    msg = msg0 + (" in key %s level 0" % key)
1085                    assert np.allclose(result[1][key][i][0],refdict[key][i][0]), msg
1086                    msg = msg0 + (" in key %s level 1" % key)
1087                    assert np.allclose(result[1][key][i][1],refdict[key][i][1]), msg
1088            else:
1089                msg = msg0 + (" in key %s" % key)
1090                assert result[1][key] == refdict[key], msg
1091        msg = msg0 + (" in key %s reflist" % key)
1092        #for (l1,l2) in zip(reflist, SGPrint(result[1])):
1093        #    assert l2.replace('\t','').replace(' ','') == l1.replace(' ',''), 'SGPrint ' +msg
1094        assert reflist == SGPrint(result[1]), 'SGPrint ' +msg
1095    for spc in spctestinp.SGdat:
1096        CompareSpcGroup(spc, 0, spctestinp.SGdat[spc], spctestinp.SGlist[spc] )
1097selftestlist.append(test1)
1098
1099def test2():
1100    '''self-test #2: SpcGroup against cctbx (sgtbx) computations'''
1101    _ReportTest()
1102    testdir = ospath.join(ospath.split(ospath.abspath( __file__ ))[0],'testinp')
1103    if ospath.exists(testdir):
1104        if testdir not in sys.path: sys.path.insert(0,testdir)
1105    import sgtbxtestinp
1106    def CompareWcctbx(spcname, cctbx_in, debug=0):
1107        'Compare output from GSASIIspc.SpcGroup with results from cctbx.sgtbx'
1108        cctbx = cctbx_in[:] # make copy so we don't delete from the original
1109        spc = (SpcGroup(spcname))[1]
1110        if debug: print spc['SpGrp']
1111        if debug: print spc['SGCen']
1112        latticetype = spcname.strip().upper()[0]
1113        # lattice type of R implies Hexagonal centering", fix the rhombohedral settings
1114        if latticetype == "R" and len(spc['SGCen']) == 1: latticetype = 'P'
1115        assert latticetype == spc['SGLatt'], "Failed: %s does not match Lattice: %s" % (spcname, spc['SGLatt'])
1116        onebar = [1]
1117        if spc['SGInv']: onebar.append(-1)
1118        for (op,off) in spc['SGOps']:
1119            for inv in onebar:
1120                for cen in spc['SGCen']:
1121                    noff = off + cen
1122                    noff = MoveToUnitCell(noff)
1123                    mult = tuple((op*inv).ravel().tolist())
1124                    if debug: print "\n%s: %s + %s" % (spcname,mult,noff)
1125                    for refop in cctbx:
1126                        if debug: print refop
1127                        # check the transform
1128                        if refop[:9] != mult: continue
1129                        if debug: print "mult match"
1130                        # check the translation
1131                        reftrans = list(refop[-3:])
1132                        reftrans = MoveToUnitCell(reftrans)
1133                        if all(abs(noff - reftrans) < 1.e-5):
1134                            cctbx.remove(refop)
1135                            break
1136                    else:
1137                        assert False, "failed on %s:\n\t %s + %s" % (spcname,mult,noff)
1138    for key in sgtbxtestinp.sgtbx:
1139        CompareWcctbx(key, sgtbxtestinp.sgtbx[key])
1140selftestlist.append(test2)
1141
1142def test3(): 
1143    '''self-test #3: exercise SytSym (includes GetOprPtrName, GenAtom, GetKNsym)
1144     for selected space groups against info in IT Volume A '''
1145    _ReportTest()
1146    def ExerciseSiteSym (spc, crdlist):
1147        'compare site symmetries and multiplicities for a specified space group'
1148        msg = "failed on site sym test for %s" % spc
1149        (E,S) = SpcGroup(spc)
1150        assert not E, msg
1151        for t in crdlist:
1152            symb, m = SytSym(t[0],S)
1153            if symb.strip() != t[2].strip() or m != t[1]:
1154                print spc,t[0],m,symb
1155            assert m == t[1]
1156            #assert symb.strip() == t[2].strip()
1157
1158    ExerciseSiteSym('p 1',[
1159            ((0.13,0.22,0.31),1,'1'),
1160            ((0,0,0),1,'1'),
1161            ])
1162    ExerciseSiteSym('p -1',[
1163            ((0.13,0.22,0.31),2,'1'),
1164            ((0,0.5,0),1,'-1'),
1165            ])
1166    ExerciseSiteSym('C 2/c',[
1167            ((0.13,0.22,0.31),8,'1'),
1168            ((0.0,.31,0.25),4,'2(010)'),
1169            ((0.25,.25,0.5),4,'-1'),
1170            ((0,0.5,0),4,'-1'),
1171            ])
1172    ExerciseSiteSym('p 2 2 2',[
1173            ((0.13,0.22,0.31),4,'1'),
1174            ((0,0.5,.31),2,'2(001)'),
1175            ((0.5,.31,0.5),2,'2(010)'),
1176            ((.11,0,0),2,'2(100)'),
1177            ((0,0.5,0),1,'222'),
1178            ])
1179    ExerciseSiteSym('p 4/n',[
1180            ((0.13,0.22,0.31),8,'1'),
1181            ((0.25,0.75,.31),4,'2(001)'),
1182            ((0.5,0.5,0.5),4,'-1'),
1183            ((0,0.5,0),4,'-1'),
1184            ((0.25,0.25,.31),2,'4(001)'),
1185            ((0.25,.75,0.5),2,'-4(001)'),
1186            ((0.25,.75,0.0),2,'-4(001)'),
1187            ])
1188    ExerciseSiteSym('p 31 2 1',[
1189            ((0.13,0.22,0.31),6,'1'),
1190            ((0.13,0.0,0.833333333),3,'2(100)'),
1191            ((0.13,0.13,0.),3,'2(110)'),
1192            ])
1193    ExerciseSiteSym('R 3 c',[
1194            ((0.13,0.22,0.31),18,'1'),
1195            ((0.0,0.0,0.31),6,'3'),
1196            ])
1197    ExerciseSiteSym('R 3 c R',[
1198            ((0.13,0.22,0.31),6,'1'),
1199            ((0.31,0.31,0.31),2,'3(111)'),
1200            ])
1201    ExerciseSiteSym('P 63 m c',[
1202            ((0.13,0.22,0.31),12,'1'),
1203            ((0.11,0.22,0.31),6,'m(100)'),
1204            ((0.333333,0.6666667,0.31),2,'3m(100)'),
1205            ((0,0,0.31),2,'3m(100)'),
1206            ])
1207    ExerciseSiteSym('I a -3',[
1208            ((0.13,0.22,0.31),48,'1'),
1209            ((0.11,0,0.25),24,'2(100)'),
1210            ((0.11,0.11,0.11),16,'3(111)'),
1211            ((0,0,0),8,'-3(111)'),
1212            ])
1213selftestlist.append(test3)
1214
1215if __name__ == '__main__':
1216    # run self-tests
1217    selftestquiet = False
1218    for test in selftestlist:
1219        test()
1220    print "OK"
Note: See TracBrowser for help on using the repository browser.