Changeset 3877


Ignore:
Timestamp:
Apr 5, 2019 3:36:26 PM (3 years ago)
Author:
vondreele
Message:

micsf again

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImath.py

    r3867 r3877  
    13871387    psin = np.sin(twopi*phase)
    13881388    pcos = np.cos(twopi*phase)
    1389     Mmod = np.sum(Am[nxs,nxs,:,:,:]*pcos[:,:,:,nxs,nxs]+Bm[nxs,nxs,:,:,:]*psin[:,:,:,nxs,nxs],axis=3)
     1389    MmodA = np.sum(Am[nxs,nxs,:,:,:]*pcos[:,:,:,nxs,nxs],axis=3)
     1390    MmodB = np.sum(Bm[nxs,nxs,:,:,:]*psin[:,:,:,nxs,nxs],axis=3)
    13901391    if SGData['SGGray']:
    1391         Mmod = -np.sum(SGMT[nxs,:,nxs,:,:]*Mmod[:,:,:,nxs,:],axis=-1)*detSM[nxs,:,nxs,nxs]
     1392        MmodA = -np.sum(SGMT[nxs,:,nxs,:,:]*MmodA[:,:,:,nxs,:],axis=-1)*detSM[nxs,:,nxs,nxs]
     1393        MmodB = -np.sum(SGMT[nxs,:,nxs,:,:]*MmodB[:,:,:,nxs,:],axis=-1)*detSM[nxs,:,nxs,nxs]
    13921394    else:
    1393         Mmod = np.sum(SGMT[nxs,:,nxs,:,:]*Mmod[:,:,:,nxs,:],axis=-1)*SGData['MagMom'][nxs,:,nxs,nxs]
    1394     return Mmod    #Ntau,Nops,Natm,,Mxyz
     1395        MmodA = np.sum(SGMT[nxs,:,nxs,:,:]*MmodA[:,:,:,nxs,:],axis=-1)*SGData['MagMom'][nxs,:,nxs,nxs]
     1396        MmodB = np.sum(SGMT[nxs,:,nxs,:,:]*MmodB[:,:,:,nxs,:],axis=-1)*SGData['MagMom'][nxs,:,nxs,nxs]
     1397    Mmod = MmodA+MmodB
     1398    return Mmod,MmodA,MmodB    #Ntau,Nops,Natm,,Mxyz; sum,Re & Im parts
    13951399       
    13961400def Modulation(H,HP,nWaves,Fmod,Xmod,Umod,glTau,glWt):
     
    14231427    return np.array([cosHA,sinHA])             # 2 x refBlk x SGops x atoms
    14241428   
     1429def MagModulation(H,HP,nWaves,Fmod,Xmod,Umod,Mmod,glTau,glWt):
     1430    '''
     1431    H: array nRefBlk x ops X hklt
     1432    HP: array nRefBlk x ops X hklt proj to hkl
     1433    nWaves: list number of waves for frac, pos, uij & mag
     1434    Fmod: array 2 x atoms x waves    (sin,cos terms)
     1435    Xmod: array atoms X 3 X ngl
     1436    Umod: array atoms x 3x3 x ngl
     1437    glTau,glWt: arrays Gauss-Lorentzian pos & wts
     1438    '''
     1439   
     1440    if nWaves[2]:       #uij (adp) waves
     1441        if len(HP.shape) > 2:
     1442            HbH = np.exp(-np.sum(HP[:,:,nxs,nxs,:]*np.inner(HP,Umod),axis=-1)) # refBlk x ops x atoms x ngl add Overhauser corr.?
     1443        else:
     1444            HbH = np.exp(-np.sum(HP[:,nxs,nxs,:]*np.inner(HP,Umod),axis=-1)) # refBlk x ops x atoms x ngl add Overhauser corr.?
     1445    else:
     1446        HbH = 1.0
     1447    HdotX = np.inner(HP,Xmod)                   #refBlk x ops x atoms X ngl
     1448    if len(H.shape) > 2:
     1449        D = H[:,:,3:]*glTau[nxs,nxs,:]              #m*e*tau; refBlk x ops X ngl
     1450        HdotXD = twopi*(HdotX+D[:,:,nxs,:])
     1451    else:
     1452        D = H[:,3:]*glTau[nxs,:]              #m*e*tau; refBlk x ops X ngl
     1453        HdotXD = twopi*(HdotX+D[:,nxs,:])
     1454    M = np.swapaxes(Mmod,1,2)
     1455    cosHA = np.sum(M[nxs,nxs,:,:,:]*(Fmod*HbH*np.cos(HdotXD)[:,:,:,nxs,:]*glWt),axis=-1)       #real part; refBlk X ops x atoms; sum for G-L integration
     1456    sinHA = np.sum(M[nxs,nxs,:,:,:]*(Fmod*HbH*np.sin(HdotXD)[:,:,:,nxs,:]*glWt),axis=-1)       #imag part; ditto
     1457    return np.array([cosHA,sinHA])             # 2 x refBlk x SGops x atoms
     1458
    14251459def ModulationTw(H,HP,nWaves,Fmod,Xmod,Umod,glTau,glWt):
    14261460    '''
  • trunk/GSASIIstrMath.py

    r3870 r3877  
    15091509    if parmDict[pfx+'isMag']:       #This part correct for making modulated mag moments on equiv atoms
    15101510       
    1511         mXYZ = np.array([[xyz[0] for xyz in list(G2spc.GenAtom(xyz,SGData,All=True,Move=True))] for xyz in (Xdata+dXdata).T])%1. #Natn,Nop,xyz
    1512         Tmag = G2mth.MagMod(ngl,mXYZ,modQ,MSSdata,SGData,SSGData)   #Ntau,Nops,Natm,Mxyz
     1511#        mXYZ = np.array([[xyz[0] for xyz in list(G2spc.GenAtom(xyz,SGData,All=True,Move=True))] for xyz in (Xdata+dXdata).T])%1. #Natn,Nop,xyz
     1512#        Tmag,TmagA,TmagB = G2mth.MagMod(ngl,mXYZ,modQ,MSSdata,SGData,SSGData)   #Ntau,Nops,Natm,Mxyz-Tmag matches drawing moments
    15131513       
    15141514        if not SGData['SGGray']:
     
    15191519            GSdata = SGData['MagMom'][nxs,:,nxs]*GSdata   #flip vectors according to spin flip * det(opM)
    15201520            GSdata = np.swapaxes(GSdata,0,1)    #Nop,Natm,Mxyz
    1521             Tmag += Gdata.T[nxs,nxs,:,:]
     1521#            Tmag += Gdata.T[nxs,nxs,:,:]
    15221522           
    1523         TmagC = np.inner(Tmag,uAmat.T)   #make cartesian; Ntau,Nops,Natm,,Mxyz
    1524         Smag = np.sqrt(np.sum(TmagC**2,axis=-1))
    1525         Kmag = TmagC/Smag[:,:,:,nxs]
     1523#        TmagC = np.inner(Tmag,uAmat.T)   #make cartesian; Ntau,Nops,Natm,,Mxyz
     1524#        Smag = np.sqrt(np.sum(TmagC**2,axis=-1))
     1525#        Kmag = TmagC/Smag[:,:,:,nxs]
    15261526
    15271527    FF = np.zeros(len(Tdata))
     
    15551555    while iBeg < nRef:
    15561556        iFin = min(iBeg+blkSize,nRef)
    1557         mRef = iFin-iBeg
    15581557        refl = refDict['RefList'][iBeg:iFin]    #array(blkSize,nItems)
    15591558        H = refl.T[:4]                          #array(blkSize,4)
     
    15911590
    15921591        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:       #TODO: mag math here??
    1593             phasem = twopi*np.inner(H.T[:,:3],mXYZ)
    1594             phasem = np.swapaxes(phasem,1,2)
    1595             cosm = np.cos(phasem)
    1596             sinm = np.sin(phasem)
     1592            GmfpuA = G2mth.MagModulation(Uniq,UniqP,nWaves,Fmod,Xmod,Umod,Mmod,glTau,glWt) #2 x refBlk x sym X atoms
     1593           
     1594           
     1595           
     1596#            phasem = twopi*np.inner(H.T[:,:3],mXYZ)
     1597#            phasem = np.swapaxes(phasem,1,2)
     1598#            cosm = np.cos(phasem)
     1599#            sinm = np.sin(phasem)
    15971600            MF = refDict['FF']['MF'][iBeg:iFin].T[Tindx].T   #Nref,Natm
    15981601            TMcorr = 0.539*(np.reshape(Tiso,Tuij.shape)*Tuij)[:,0,:]*Fdata*Mdata*MF/(2*Nops)     #Nref,Natm
    1599                      
     1602#                     
    16001603#            HM = np.inner(Bmat,HP.T)                            #put into cartesian space
    1601 #            HM = HM/np.sqrt(np.sum(HM**2,axis=0))               #& normalize
     1604#            eM = HM/np.sqrt(np.sum(HM**2,axis=0))               #& normalize
     1605#            cosq = np.inner(eM.T,Kmag)
    16021606#
    1603 #            fam = TMcorr[:,nxs,nxs,:,nxs]*Tmag[nxs,:,:,:,:]*cosm[:,nxs,:,:,nxs]    #Nref,Ntau,Nops,Natm,Mxyz
    1604 #            fbm = TMcorr[:,nxs,nxs,:,nxs]*Tmag[nxs,:,:,:,:]*sinm[:,nxs,:,:,nxs]   
    1605 #                       
    1606 #            famq = np.sum(np.sum(fam,axis=-2),axis=-2)        #Nref,Mxyz; sum ops & atoms
    1607 #            fbmq = np.sum(np.sum(fbm,axis=-2),axis=-2)
     1607#            fam = TMcorr[:,nxs,nxs,:,nxs]*(TmagA[nxs,:,:,:,:]*cosm[:,nxs,:,:,nxs]-TmagB[nxs,:,:,:,:]*sinm[:,nxs,:,:,nxs])/2.    #Nref,Ntau,Nops,Natm,Mxyz
     1608#            fbm = TMcorr[:,nxs,nxs,:,nxs]*(TmagA[nxs,:,:,:,:]*sinm[:,nxs,:,:,nxs]+TmagB[nxs,:,:,:,:]*cosm[:,nxs,:,:,nxs])/2.
    16081609#           
    1609 #            fas = np.sum(famq,axis=-1)**2-np.sum(HM.T[:,nxs,:]*famq,axis=-1)**2   #mag intensity calc F^2-(e.F)^2
    1610 #            fbs = np.sum(fbmq,axis=-1)**2-np.sum(HM.T[:,nxs,:]*fbmq,axis=-1)**2
     1610#           
     1611#            famq = np.sum(np.sum(fam,axis=2),axis=2)        #Nref,Mxyz; sum ops & atoms
     1612#            fbmq = np.sum(np.sum(fbm,axis=2),axis=2)           
     1613#           
     1614##            fas = np.sum(famq,axis=-1)**2-np.sum(HM.T[:,nxs,:]*famq,axis=-1)**2   #mag intensity calc F^2-(e.F)^2
     1615##            fbs = np.sum(fbmq,axis=-1)**2-np.sum(HM.T[:,nxs,:]*fbmq,axis=-1)**2
     1616#            fas = famq**2-(cosq*famq)**2
     1617#            fbs = fbmq**2-(cosq*fbmq)**2                       
     1618#
     1619#            refl.T[10] = np.sum(fas/ngl,axis=-1)+np.sum(fbs/ngl,axis=-1)
     1620#            refl.T[11] = atan2d(fbs[:,0],fas[:,0])
    16111621#for modulated moments --> m != 0 reflections
    16121622#            M = np.array(np.abs(H[3]),dtype=np.int)-1
     
    16261636#            sinm = np.swapaxes(np.sin(mphase),1,2)    #--> Nref,Nop,Natm,Ntau
    16271637#            cosm = np.swapaxes(np.cos(mphase),1,2)                               #ditto
    1628            
    1629             HM = np.inner(Bmat,HP.T)                             #put into cartesian space
    1630             HM = HM/np.sqrt(np.sum(HM**2,axis=0))               #normalize
    1631             eDotK = np.sum(HM.T[:,nxs,nxs,nxs,:]*Kmag[nxs,:,:,:,:],axis=-1)
    1632             Q = HM.T[:,nxs,nxs,nxs,:]*eDotK[:,:,:,:,nxs]-Kmag[nxs,:,:,:,:] #Nref,Ntau,Nop,Natm,Mxyz
    1633 
     1638#           
     1639#            HM = np.inner(Bmat,HP.T)                             #put into cartesian space
     1640#            HM = HM/np.sqrt(np.sum(HM**2,axis=0))               #normalize
     1641#            eDotK = np.sum(HM.T[:,nxs,nxs,nxs,:]*Kmag[nxs,:,:,:,:],axis=-1)
     1642#            Q = HM.T[:,nxs,nxs,nxs,:]*eDotK[:,:,:,:,nxs]-Kmag[nxs,:,:,:,:] #Nref,Ntau,Nop,Natm,Mxyz
     1643#
    16341644            fam = (Q*TMcorr[:,nxs,nxs,:,nxs]*cosm[:,nxs,:,:,nxs]*Smag[nxs,:,:,:,nxs])   #Nref,Ntau,Nop,Natm,Mxyz
    16351645            fbm = (Q*TMcorr[:,nxs,nxs,:,nxs]*sinm[:,nxs,:,:,nxs]*Smag[nxs,:,:,:,nxs])
    1636            
    1637             fams = np.sqrt(np.sum(fam**2,axis=-1))
    1638             fbms = np.sqrt(np.sum(fbm**2,axis=-1))
    1639                        
    1640             fas = np.sum(np.sum(fams,axis=2),axis=2)      #Nref,ntau,Mxyz; sum ops & atoms
    1641             fbs = np.sum(np.sum(fbms,axis=2),axis=2)      #ditto
    1642            
    1643             refl.T[10] = np.sum(fas/ngl,axis=-1)**2+np.sum(fbs/ngl,axis=-1)**2
    1644             refl.T[11] = atan2d(fbs[:,0],fas[:,0])
     1646#           
     1647#            fams = np.sqrt(np.sum(fam**2,axis=-1))
     1648#            fbms = np.sqrt(np.sum(fbm**2,axis=-1))
     1649#                       
     1650#            fas = np.sum(np.sum(fams,axis=2),axis=2)      #Nref,ntau,Mxyz; sum ops & atoms
     1651#            fbs = np.sum(np.sum(fbms,axis=2),axis=2)      #ditto
     1652#           
     1653#            refl.T[10] = np.sum(fas/ngl,axis=-1)**2+np.sum(fbs/ngl,axis=-1)**2
     1654#            refl.T[11] = atan2d(fbs[:,0],fas[:,0])
    16451655           
    16461656        else:
Note: See TracChangeset for help on using the changeset viewer.