Changeset 1988
 Timestamp:
 Oct 6, 2015 2:53:02 PM (8 years ago)
 Location:
 trunk
 Files:

 4 edited
Legend:
 Unmodified
 Added
 Removed

trunk/GSASIImath.py
r1987 r1988 926 926 ################################################################################ 927 927 928 def Modulation(waveTypes,Uniq,FSSdata,XSSdata,USSdata,Mast): 928 def Modulation(waveTypes,H,FSSdata,XSSdata,USSdata,Mast): 929 ''' 930 H: array ops X hklt 931 Ax: array atoms X waves X xyz 932 Bx: array atoms X waves X xyz 933 ''' 929 934 930 935 nxs = np.newaxis 931 936 glTau,glWt = pwd.pygauleg(0.,1.,32) #get GaussLegendre intervals & weights 932 933 def expModInt(H,Af,Bf,Ax,Bx,Au,Bu):934 '''935 H: array ops X hklt936 Ax: array atoms X waves X xyz937 Bx: array atoms X waves X xyz938 S: array ops939 '''940 if 'Fourier' in waveTypes:941 nf = 0942 nx = 0943 XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32))944 FmodZ = np.zeros((Af.shape[0],Af.shape[1],32))945 if 'Crenel' in waveTypes:946 nC = np.where('Crenel' in waveTypes)947 nf = 1948 #FmodZ = 0 replace949 else:950 nx = 1951 if 'Sawtooth' in waveTypes:952 nS = np.where('Sawtooth' in waveTypes)953 #XmodZ = 0 replace954 if Af.shape[1]:955 tauF = np.arange(1.,Af.shape[1]+1nf)[:,nxs]*glTau #Fwaves x 32956 FmodA = Af[:,nf:,nxs]*np.sin(twopi*tauF[nxs,:,:]) #atoms X Fwaves X 32957 FmodB = Bf[:,nf:,nxs]*np.cos(twopi*tauF[nxs,:,:])958 Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves959 else:960 Fmod = 1.0961 if Ax.shape[1]:962 tauX = np.arange(1.,Ax.shape[1]+1nx)[:,nxs]*glTau #Xwaves x 32963 XmodA = Ax[:,nx:,:,nxs]*np.sin(twopi*tauX[nxs,:,nxs,:]) #atoms X waves X pos X 32964 XmodB = Bx[:,nx:,:,nxs]*np.cos(twopi*tauX[nxs,:,nxs,:]) #ditto965 Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1) #atoms X pos X 32; sum waves966 if Au.shape[1]:967 tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 32968 UmodA = Au[:,:,:,:,nxs]*np.sin(twopi*tauU[nxs,:,nxs,nxs,:]) #atoms x waves x 3x3 x 32969 UmodB = Bu[:,:,:,:,nxs]*np.cos(twopi*tauU[nxs,:,nxs,nxs,:]) #ditto970 Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves971 HbH = np.exp(np.sum(H[:,:,nxs,nxs,:3]*np.inner(H[:,:,:3],Umod),axis=1)) # refBlk x ops x atoms x 32 add Overhauser corr.?972 else:973 HbH = 1.0974 D = H[:,:,3][:,:,nxs]*glTau[nxs,nxs,:] #m*tau;refBlk x ops X 32975 HdotX = np.inner(H[:,:,:3],np.swapaxes(Xmod,1,2))+D[:,:,nxs,:] #refBlk X ops x atoms X 32976 sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=1) #refBlk X ops x atoms; sum for GL integration977 cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=1) #ditto978 # GSASIIpath.IPyBreak()979 return cosHA,sinHA #each refBlk x ops X atoms980 981 937 Ax = np.array(XSSdata[:3]).T #atoms x waves x sin pos mods 982 938 Bx = np.array(XSSdata[3:]).T #...cos pos mods … … 985 941 Au = Mast*np.array(G2lat.U6toUij(USSdata[:6])).T #atoms x waves x sin Uij mods 986 942 Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T #...cos Uij mods 987 GpA = np.array(expModInt(Uniq,Af,Bf,Ax,Bx,Au,Bu)) 988 return GpA # 2 x refBlk x SGops x atoms 989 990 def ModulationDerv(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast): 943 944 if 'Fourier' in waveTypes: 945 nf = 0 946 nx = 0 947 XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32)) 948 FmodZ = np.zeros((Af.shape[0],Af.shape[1],32)) 949 if 'Crenel' in waveTypes: 950 nC = np.where('Crenel' in waveTypes) 951 nf = 1 952 #FmodZ = 0 replace 953 else: 954 nx = 1 955 if 'Sawtooth' in waveTypes: 956 nS = np.where('Sawtooth' in waveTypes) 957 #XmodZ = 0 replace 958 if Af.shape[1]: 959 tauF = np.arange(1.,Af.shape[1]+1nf)[:,nxs]*glTau #Fwaves x 32 960 FmodA = Af[:,nf:,nxs]*np.sin(twopi*tauF[nxs,:,:]) #atoms X Fwaves X 32 961 FmodB = Bf[:,nf:,nxs]*np.cos(twopi*tauF[nxs,:,:]) 962 Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves 963 else: 964 Fmod = 1.0 965 if Ax.shape[1]: 966 tauX = np.arange(1.,Ax.shape[1]+1nx)[:,nxs]*glTau #Xwaves x 32 967 XmodA = Ax[:,nx:,:,nxs]*np.sin(twopi*tauX[nxs,:,nxs,:]) #atoms X waves X pos X 32 968 XmodB = Bx[:,nx:,:,nxs]*np.cos(twopi*tauX[nxs,:,nxs,:]) #ditto 969 Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1) #atoms X pos X 32; sum waves 970 if Au.shape[1]: 971 tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 32 972 UmodA = Au[:,:,:,:,nxs]*np.sin(twopi*tauU[nxs,:,nxs,nxs,:]) #atoms x waves x 3x3 x 32 973 UmodB = Bu[:,:,:,:,nxs]*np.cos(twopi*tauU[nxs,:,nxs,nxs,:]) #ditto 974 Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves 975 HbH = np.exp(np.sum(H[:,:,nxs,nxs,:3]*np.inner(H[:,:,:3],Umod),axis=1)) # refBlk x ops x atoms x 32 add Overhauser corr.? 976 else: 977 HbH = 1.0 978 D = H[:,:,3][:,:,nxs]*glTau[nxs,nxs,:] #m*tau;refBlk x ops X 32 979 HdotX = np.inner(H[:,:,:3],np.swapaxes(Xmod,1,2))+D[:,:,nxs,:] #refBlk X ops x atoms X 32 980 sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=1) #refBlk X ops x atoms; sum for GL integration 981 cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=1) #ditto 982 # GSASIIpath.IPyBreak() 983 984 return np.array([cosHA,sinHA]) # 2 x refBlk x SGops x atoms 985 986 def ModulationDerv(waveTypes,H,FSSdata,XSSdata,USSdata,Mast): 987 ''' 988 H: array ops X hklt 989 Ax: array atoms X waves X xyz 990 Bx: array atoms X waves X xyz 991 ''' 991 992 992 993 nxs = np.newaxis 993 994 glTau,glWt = pwd.pygauleg(0.,1.,32) #get GaussLegendre intervals & weights 994 995 def expModInt(H,Af,Bf,Ax,Bx,Au,Bu):996 '''997 H: array ops X hklt998 Ax: array atoms X waves X xyz999 Bx: array atoms X waves X xyz1000 S: array ops1001 '''1002 if 'Fourier' in waveTypes:1003 nf = 01004 nx = 01005 XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32))1006 FmodZ = np.zeros((Af.shape[0],Af.shape[1],32))1007 if 'Crenel' in waveTypes:1008 nC = np.where('Crenel' in waveTypes)1009 nf = 11010 #FmodZ = 0 replace1011 else:1012 nx = 11013 if 'Sawtooth' in waveTypes:1014 nS = np.where('Sawtooth' in waveTypes)1015 #XmodZ = 0 replace1016 if Af.shape[1]:1017 tauF = np.arange(1.,Af.shape[1]+1nf)[:,nxs]*glTau #Fwaves x 321018 FmodA = Af[:,nf:,nxs]*np.sin(twopi*tauF[nxs,:,:]) #atoms X Fwaves X 321019 FmodB = Bf[:,nf:,nxs]*np.cos(twopi*tauF[nxs,:,:])1020 Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves1021 else:1022 Fmod = 1.01023 if Ax.shape[1]:1024 tauX = np.arange(1.,Ax.shape[1]+1nx)[:,nxs]*glTau #Xwaves x 321025 XmodA = Ax[:,nx:,:,nxs]*np.sin(twopi*tauX[nxs,:,nxs,:]) #atoms X waves X pos X 321026 XmodB = Bx[:,nx:,:,nxs]*np.cos(twopi*tauX[nxs,:,nxs,:]) #ditto1027 Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1) #atoms X pos X 32; sum waves1028 if Au.shape[1]:1029 tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 321030 UmodA = Au[:,:,:,:,nxs]*np.sin(twopi*tauU[nxs,:,nxs,nxs,:]) #atoms x waves x 3x3 x 321031 UmodB = Bu[:,:,:,:,nxs]*np.cos(twopi*tauU[nxs,:,nxs,nxs,:]) #ditto1032 Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves1033 HbH = np.exp(np.swapaxes(np.array([np.inner(h,np.inner(Umod,h)) for h in H[:,:3]]),1,2)).T #Overhauser corr.?1034 else:1035 HbH = 1.01036 D = H[:,3][:,nxs]*glTau[nxs,:] #m*tau; ops X 321037 HdotX = np.inner(np.swapaxes(Xmod,1,2),H[:,:3])+D.T[nxs,:,:] #atoms X 32 X ops1038 sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt[nxs,:,nxs],axis=1) #atoms X ops; sum for GL integration1039 cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt[nxs,:,nxs],axis=1) #ditto1040 # GSASIIpath.IPyBreak()1041 return cosHA.T,sinHA.T #ops X atoms1042 1043 995 Ax = np.array(XSSdata[:3]).T #atoms x waves x sin pos mods 1044 996 dGdAx = np.zeros_like(Ax) … … 1053 1005 Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T #...cos Uij mods 1054 1006 dGdBu = np.zeros_like(Bu) 1055 GpA = np.array(expModInt(SSUniq,Af,Bf,Ax,Bx,Au,Bu)) 1056 return GpA,dGdAf,dGdBf,dGdAx,dGdBx,dGdAu,dGdBu 1007 1008 if 'Fourier' in waveTypes: 1009 nf = 0 1010 nx = 0 1011 XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32)) 1012 FmodZ = np.zeros((Af.shape[0],Af.shape[1],32)) 1013 if 'Crenel' in waveTypes: 1014 nC = np.where('Crenel' in waveTypes) 1015 nf = 1 1016 #FmodZ = 0 replace 1017 else: 1018 nx = 1 1019 if 'Sawtooth' in waveTypes: 1020 nS = np.where('Sawtooth' in waveTypes) 1021 #XmodZ = 0 replace 1022 if Af.shape[1]: 1023 tauF = np.arange(1.,Af.shape[1]+1nf)[:,nxs]*glTau #Fwaves x 32 1024 StauF = np.ones_like(Af)[:,nf:,nxs]*np.sin(twopi*tauF)[nxs,:,:] #also dFmodA/dAf 1025 CtauF = np.ones_like(Bf)[:,nf:,nxs]*np.cos(twopi*tauF)[nxs,:,:] #also dFmodB/dBf 1026 FmodA = Af[:,nf:,nxs]*StauF #atoms X Fwaves X 32 1027 FmodB = Bf[:,nf:,nxs]*CtauF 1028 Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves 1029 else: 1030 Fmod = 1.0 1031 if Ax.shape[1]: 1032 tauX = np.arange(1.,Ax.shape[1]+1nx)[:,nxs]*glTau #Xwaves x 32 1033 StauX = np.ones_like(Ax)[:,nx:,:,nxs]*np.sin(twopi*tauX)[nxs,:,nxs,:] #also dXmodA/dAx 1034 CtauX = np.ones_like(Bx)[:,nx:,:,nxs]*np.cos(twopi*tauX)[nxs,:,nxs,:] #also dXmodB/dBx 1035 XmodA = Ax[:,nx:,:,nxs]*StauX #atoms X waves X pos X 32 1036 XmodB = Bx[:,nx:,:,nxs]*CtauX #ditto 1037 Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1) #atoms X pos X 32; sum waves 1038 if Au.shape[1]: 1039 tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 32 1040 StauU = np.ones_like(Au)[:,:,:,:,nxs]*np.sin(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmodA/dAu 1041 CtauU = np.ones_like(Bu)[:,:,:,:,nxs]*np.cos(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmodB/dBu 1042 UmodA = Au[:,:,:,:,nxs]*StauU #atoms x waves x 3x3 x 32 1043 UmodB = Bu[:,:,:,:,nxs]*CtauU #ditto 1044 Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves 1045 HbH = np.exp(np.sum(H[:,nxs,nxs,:3]*np.inner(H[:,:3],Umod),axis=1)) # ops x atoms x 32 add Overhauser corr.? 1046 else: 1047 HbH = 1.0 1048 D = H[:,3][:,nxs]*glTau[nxs,:] #m*tau; ops X 32 1049 HdotX = np.inner(H[:,:3],np.swapaxes(Xmod,1,2))+D[:,nxs,:] #ops x atoms X 32 1050 sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=1) #ops x atoms; sum for GL integration 1051 cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=1) #ditto 1052 # GSASIIpath.IPyBreak() 1053 return np.array([cosHA,sinHA]),dGdAf,dGdBf,dGdAx,dGdBx,dGdAu,dGdBu #ops X atoms 1057 1054 1058 1055 def posFourier(tau,psin,pcos,smul): 
trunk/GSASIIphsGUI.py
r1981 r1988 1872 1872 generalData = data['General'] 1873 1873 Type = generalData['Type'] 1874 if Type in ['nuclear','macromolecular']: 1875 choice = ['F  site fraction','X  coordinates','U  thermal parameters'] 1876 elif Type == 'magnetic': 1877 choice = ['F  site fraction','X  coordinates','U  thermal parameters'] 1874 choice = ['F  site fraction','X  coordinates','U  thermal parameters'] 1878 1875 dlg = wx.MultiChoiceDialog(G2frame,'Select','Refinement controls',choice) 1879 1876 if dlg.ShowModal() == wx.ID_OK: 
trunk/GSASIIrestrGUI.py
r1831 r1988 89 89 Types += [atom[ct] for atom in Atoms] 90 90 Coords += [atom[cx:cx+3] for atom in Atoms] 91 Ids += [atom[ 1] for atom in Atoms]91 Ids += [atom[cia+8] for atom in Atoms] 92 92 rama = G2data.ramachandranDist['All'] 93 93 ramaName = 'All' 
trunk/GSASIIstrMath.py
r1987 r1988 1290 1290 if nTwin > 1: 1291 1291 dfadx = np.array([np.sum(twopi*Uniq[it,:,:3]*np.swapaxes(fax,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)]) 1292 dfadua = np.array([np.sum(Hij[it]*np.swapaxes(fa ,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)])1292 dfadua = np.array([np.sum(Hij[it]*np.swapaxes(fag,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)]) 1293 1293 # array(nTwin,2,nAtom,3) & array(nTwin,2,nAtom,6) 1294 1294 else: 1295 1295 dfadx = np.sum(twopi*Uniq[:,:3]*np.swapaxes(fax,2,1)[:,:,:,np.newaxis],axis=2) 1296 dfadua = np.sum(Hij*np.swapaxes(fa ,2,1)[:,:,:,np.newaxis],axis=2)1296 dfadua = np.sum(Hij*np.swapaxes(fag,2,1)[:,:,:,np.newaxis],axis=2) 1297 1297 # array(2,nAtom,3) & array(2,nAtom,6) 1298 1298 #NB: the above have been checked against PA(1:10,1:2) in strfctr.for for al2O3! … … 1303 1303 if len(TwinLaw) > 1: 1304 1304 dfbdx = np.array([np.sum(twopi*Uniq[it]*np.swapaxes(fbx,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)]) 1305 dfbdua = np.array([np.sum(Hij[it]*np.swapaxes(fb ,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)])1305 dfbdua = np.array([np.sum(Hij[it]*np.swapaxes(fbg,2,1)[:,it,:,:,np.newaxis],axis=2) for it in range(nTwin)]) 1306 1306 else: 1307 1307 dfadfl = np.sum(FPP*Tcorr*sinp) 1308 1308 dfbdfl = np.sum(FPP*Tcorr*cosp) 1309 1309 dfbdx = np.sum(twopi*Uniq*np.swapaxes(fbx,2,1)[:,:,:,np.newaxis],axis=2) 1310 dfbdua = np.sum(Hij*np.swapaxes(fb ,2,1)[:,:,:,np.newaxis],axis=2)1310 dfbdua = np.sum(Hij*np.swapaxes(fbg,2,1)[:,:,:,np.newaxis],axis=2) 1311 1311 else: 1312 1312 dfbdfr = np.zeros_like(dfadfr)
Note: See TracChangeset
for help on using the changeset viewer.