Changeset 3784


Ignore:
Timestamp:
Jan 16, 2019 7:03:59 PM (3 years ago)
Author:
toby
Message:

generate lattice constraints for mag. xform

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIconstrGUI.py

    r3760 r3784  
    913913                                eqString[-1] += ' - '
    914914                                m = abs(m)
    915                         eqString[-1] += '%.3f*%s '%(m,var)
     915                        eqString[-1] += '{:g}*{:} '.format(m,var)
    916916                        varMean = G2obj.fmtVarDescr(var)
    917917                        helptext += "\n" + var + " ("+ varMean + ")"
     
    937937                        if len(eqString[-1]) > maxlen:
    938938                            eqString.append(' ')
     939                        m = term[0]
    939940                        if eqString[-1] != '':
    940941                            if term[0] > 0:
     
    942943                            else:
    943944                                eqString[-1] += ' - '
    944                         eqString[-1] += '%.3f*%s '%(abs(term[0]),var)
     945                                m = -term[0]
     946                        eqString[-1] += '{:g}*{:} '.format(m,var)
    945947                        varMean = G2obj.fmtVarDescr(var)
    946948                        helptext += "\n" + var + " ("+ varMean + ")"
     
    961963                            eqString[-1] += ' = '
    962964                        if normval/term[0] == 1:
    963                             eqString[-1] += '%s'% var
     965                            eqString[-1] += '{:}'.format(var)
    964966                        else:
    965                             eqString[-1] += '%.3f*%s'%(normval/term[0],var)
     967                            eqString[-1] += '{:g}*{:} '.format(normval/term[0],var)
    966968                        varMean = G2obj.fmtVarDescr(var)
    967969                        helptext += "\n" + var + " ("+ varMean + ")"
    968970                    if normval/indepterm[0] == 1:
    969                         eqString[-1] += ' = %s'% str(indepterm[1])
     971                        eqString[-1] += ' = {:} '.format(indepterm[1])
    970972                    else:
    971                         eqString[-1] += ' = %.3f*%s'%(normval/indepterm[0],str(indepterm[1]))
     973                        eqString[-1] += ' = {:g}*{:} '.format(normval/indepterm[0],str(indepterm[1]))
    972974                    typeString = 'EQUIV'
    973975                else:
     
    10141016        del data[name][Id]
    10151017        allcons = FindAllCons(data)     #should I call CheckChangedConstraint() instead?
    1016         if not len(allcons): return
    1017         CheckConstraints(allcons)
     1018        if len(allcons):
     1019            CheckConstraints(allcons)
    10181020        wx.CallAfter(OnPageChanged,None)
    1019        
     1021
    10201022    def OnConstEdit(event):
    10211023        '''Called to edit an individual contraint in response to a
     
    13321334    item = G2gd.GetGPXtreeItemId(G2frame,G2frame.root,'Constraints')
    13331335    if not item:
     1336        print('Error: no constraints in Data Tree')
    13341337        return
    13351338    constraints = G2frame.GPXtree.GetItemPyData(item)
     
    13371340    varyList = []
    13381341    xnames = ['dAx','dAy','dAz']
    1339 #    invTrans = nl.inv(Trans)
    1340 #    Us = ['AU11','AU22','AU33','AU12','AU13','AU23']
    1341 #    Uids = [[0,0,'AU11'],[1,1,'AU22'],[2,2,'AU33'],[0,1,'AU12'],[0,2,'AU13'],[1,2,'AU23']]
     1342    # constraints on matching atom params between phases
    13421343    for ia,code in enumerate(atCodes):
    13431344        atom = nAtoms[ia]
     
    13771378            DepCons = [1.0,G2obj.G2VarObj('%d::%s:%s'%(opId,name,iat))]
    13781379            constraints['Phase'].append([DepCons,IndpCon,None,None,'e'])
    1379        
     1380           
     1381        # unfinished Anisotropic constraint generation
     1382#        Uids = [[0,0,'AU11'],[1,1,'AU22'],[2,2,'AU33'],[0,1,'AU12'],[0,2,'AU13'],[1,2,'AU23']]
    13801383#        DepConsDict = dict(zip(Us,[[],[],[],[],[],[]]))
    13811384#        for iu,Uid in enumerate(Uids):
     
    14091412        #how do I do Uij's for most Trans?
    14101413       
    1411 #unfortunately, this doesn't always work!
    1412 #    As = ['A0','A1','A2','A3','A4','A5']
    1413 #    Aids = [[0,0,'A0'],[1,1,'A1'],[2,2,'A2'],[0,1,'A3'],[0,2,'A4'],[1,2,'A5']]
    1414 #    DepConsDict = dict(zip(As,[[],[],[],[],[],[]]))
    1415 #    T = Trans
    1416 ##Symbolic code:
    1417 #    '''
     1414    T = nl.inv(Trans)
     1415    conMat = [
     1416        [T[0,0]**2,T[0,1]**2,T[0,2]**2,T[0,0]*T[0,1],T[0,0]*T[0,2],T[0,1]*T[0,2]],
     1417        [T[1,0]**2,T[1,1]**2,T[1,2]**2,T[1,0]*T[1,1],T[1,0]*T[1,2],T[1,1]*T[1,2]],
     1418        [T[2,0]**2,T[2,1]**2,T[2,2]**2,T[2,0]*T[2,1],T[2,0]*T[2,2],T[2,1]*T[2,2]],
     1419        [2.*T[0,0]*T[1,0],2.*T[0,1]*T[1,1],2.*T[0,2]*T[1,2],T[0,0]*T[1,1]+T[0,1]*T[1,0],T[0,0]*T[1,2]+T[0,2]*T[1,0],T[0,1]*T[1,2]+T[0,2]*T[1,1]],
     1420        [2.*T[0,0]*T[2,0],2.*T[0,1]*T[2,1],2.*T[0,2]*T[2,2],T[0,0]*T[2,1]+T[0,1]*T[2,0],T[0,0]*T[2,2]+T[0,2]*T[2,0],T[0,1]*T[2,2]+T[0,2]*T[2,1]],
     1421        [2.*T[1,0]*T[2,0],2.*T[1,1]*T[2,1],2.*T[1,2]*T[2,2],T[1,0]*T[2,1]+T[1,1]*T[2,0],T[1,0]*T[2,2]+T[1,2]*T[2,0],T[1,1]*T[2,2]+T[1,2]*T[2,1]]]
     1422    # Gnew = conMat * A:
    14181423#         T00**2*a0  T01**2*a1 T02**2*a2 T00*T01*a3    T00*T02*a4    T01*T02*a5
    14191424#         T10**2*a0  T11**2*a1 T12**2*a2 T10*T11*a3    T10*T12*a4    T11*T12*a5
     
    14221427#         2*T00*T20*a0      2*T01*T21*a1     2*T02*T22*a2     (T00*T21 + T01*T20)*a3      (T00*T22 + T02*T20)*a4      (T01*T22 + T02*T21)*a5
    14231428#         2*T10*T20*a0      2*T11*T21*a1     2*T12*T22*a2     (T10*T21 + T11*T20)*a3      (T10*T22 + T12*T20)*a4      (T11*T22 + T12*T21)*a5
    1424 #    '''
    1425 #    conMat = [
    1426 #        [T[0,0]**2,T[0,1]**2,T[0,2]**2,T[0,0]*T[0,1],T[0,0]*T[0,2],T[0,1]*T[0,2]],
    1427 #        [T[1,0]**2,T[1,1]**2,T[1,2]**2,T[1,0]*T[1,1],T[1,0]*T[1,2],T[1,1]*T[1,2]],
    1428 #        [T[2,0]**2,T[2,1]**2,T[2,2]**2,T[2,0]*T[2,1],T[2,0]*T[2,2],T[2,1]*T[2,2]],
    1429 #        [2.*T[0,0]*T[1,0],2.*T[0,1]*T[1,1],2.*T[0,2]*T[1,2],T[0,0]*T[1,1]+T[0,1]*T[1,0],T[0,0]*T[1,2]+T[0,2]*T[1,0],T[0,1]*T[1,2]+T[0,2]*T[1,1]],
    1430 #        [2.*T[0,0]*T[2,0],2.*T[0,1]*T[2,1],2.*T[0,2]*T[2,2],T[0,0]*T[2,1]+T[0,1]*T[2,0],T[0,0]*T[2,2]+T[0,2]*T[2,0],T[0,1]*T[2,2]+T[0,2]*T[2,1]],
    1431 #        [2.*T[1,0]*T[2,0],2.*T[1,1]*T[2,1],2.*T[1,2]*T[2,2],T[1,0]*T[2,1]+T[1,1]*T[2,0],T[1,0]*T[2,2]+T[1,2]*T[2,0],T[1,1]*T[2,2]+T[1,2]*T[2,1]]]
    1432 #   
    1433 #    for iA,Aid in enumerate(Aids):
    1434 #        if abs(nAcof[iA]) > 1.e-8:
    1435 #            for ia in [0,1,2,3,4,5]:
    1436 #                cA = conMat[ia][iA]
    1437 #                if abs(cA) > 1.e-8:
    1438 #                    parm = SetUniqAj(npId,ia,nSGData)
    1439 #                    DepConsDict[Aid[2]].append([cA,G2obj.G2VarObj(parm)])
    1440 #    conStrings = []
    1441 #    for iA,Asi in enumerate(As):
    1442 #        parm = SetUniqAj(opId,iA,oSGData)
    1443 #        parmDict[parm] = oAcof[iA]
    1444 #        varyList.append(parm)
    1445 #        IndpCon = [1.0,G2obj.G2VarObj(parm)]
    1446 #        conStr = str([IndpCon,DepConsDict[Asi]])
    1447 #        if conStr in conStrings:
    1448 #            continue
    1449 #        conStrings.append(conStr)
    1450 #        if len(DepConsDict[Asi]) == 1:
    1451 #            if DepConsDict[Asi][0]:
    1452 #                constraints['Phase'].append([IndpCon,DepConsDict[Asi][0],None,None,'e'])
    1453 #        elif len(DepConsDict[Asi]) > 1:       
    1454 #            for Dep in DepConsDict[Asi]:
    1455 #                Dep[0] *= -1
    1456 #            constraints['Phase'].append([IndpCon]+DepConsDict[Asi]+[0.0,None,'c'])
    1457            
     1429    # Generated as symbolic code using:
     1430    # import sym
     1431    # A0, A1, A2, A3, A4, A5 = sym.symbols('A0, A1, A2, A3, A4, A5')
     1432    # G = sym.Matrix([ [A0,    A3/2.,  A4/2.], [A3/2.,  A1,    A5/2.], [A4/2., A5/2.,    A2]])
     1433    # transformation matrix
     1434    # T00, T10, T20, T01, T11, T21, T02, T12, T22 = sym.symbols('T00, T10, T20, T01, T11, T21, T02, T12, T22')
     1435    # Tr = sym.Matrix([ [T00, T10, T20], [T01, T11, T21], [T02, T12, T22],])
     1436    # Gnew = (Tr.T*G)*Tr
     1437   
     1438    for iAnew,Asi in enumerate(['A0','A1','A2','A3','A4','A5']): # loop through A[i] for new cell
     1439        Nparm = str(npId) + '::' + Asi
     1440        if Nparm != SetUniqAj(npId,iAnew,nSGData): continue # skip if already constrainted
     1441        multDict = {}
     1442        for iAorg in range(6):
     1443            cA = conMat[iAorg][iAnew] # coeff for A[i] in constraint matrix
     1444            if abs(cA) < 1.e-8: continue
     1445            parm = SetUniqAj(opId,iAorg,oSGData) # translate to unique A[i] in original cell
     1446            if not parm: continue
     1447            # sum coeff
     1448            if parm in multDict:
     1449                multDict[parm] += cA
     1450            else:
     1451                multDict[parm] = cA
     1452        # any non-zero multipliers?
     1453        maxMult = 0
     1454        for i in multDict:
     1455            maxMult = max(maxMult,abs(multDict[i]))
     1456        if maxMult <= 0: continue
     1457           
     1458        # create constraint (if needed) or equivalence
     1459        if len(multDict) == 1:
     1460            key = list(multDict.keys())[0]
     1461            constr = [[multDict[key],G2obj.G2VarObj(Nparm)],
     1462                          [1.0,G2obj.G2VarObj(key)],None,None,'e']
     1463        else:
     1464            constr = [[-1.0,G2obj.G2VarObj(Nparm)]]
     1465        for key in multDict:
     1466            constr += [[multDict[key],G2obj.G2VarObj(key)]]
     1467        constr += [0.0,None,'c']
     1468        constraints['Phase'] += [constr]
     1469   
     1470    # constraints on HAP Scale, etc.
    14581471    for hId,hist in enumerate(UseList):    #HAP - seems OK
    14591472        ohapkey = '%d:%d:'%(opId,hId)
Note: See TracChangeset for help on using the changeset viewer.