Changeset 1977

Timestamp:

Sep 23, 2015 3:28:12 PM (7 years ago)

Author:

vondreele

Message:

More modulation structure factor stuff, begin SSF derivatives

Location:

trunk

Files:

• GSASIImath.py (modified) (4 diffs)
• GSASIIstrMath.py (modified) (4 diffs)

trunk/GSASIImath.py

@@ -34 +34 @@
 import numpy.fft as fft
 import scipy.optimize as so
-import pypowder as pyd
+import pypowder as pwd
 
 sind = lambda x: np.sin(x*np.pi/180.)
@@ -924 +924 @@
 ################################################################################
     
-def Modulation(waveTypes,SSUniq,SGT,FSSdata,XSSdata,USSdata,Mast):
-    import pypowder as pwd
+def Modulation(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast):
     
     nxs = np.newaxis
     glTau,glWt = pwd.pygauleg(0.,1.,32)         #get Gauss-Legendre intervals & weights
     
-    def expModInt(H,Af,Bf,Ax,Bx,Au,Bu,S):
+    def expModInt(H,Af,Bf,Ax,Bx,Au,Bu):
         '''
         H: array ops X hklt
@@ -971 +970 @@
         else:
             HbH = 1.0
-        D = H[:,3][:,nxs]*(glTau[nxs,:]+S[:,nxs])              #m*tau; ops X 32
+        D = H[:,3][:,nxs]*glTau[nxs,:]              #m*tau; ops X 32
         HdotX = np.inner(np.swapaxes(Xmod,1,2),H[:,:3])+D.T[nxs,:,:]     #atoms X 32 X ops
         sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt[nxs,:,nxs],axis=1)       #atoms X ops; sum for G-L integration
         cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt[nxs,:,nxs],axis=1)       #ditto
-        GSASIIpath.IPyBreak()                  
+#        GSASIIpath.IPyBreak()                  
         return cosHA.T,sinHA.T      #ops X atoms
 
@@ -984 +983 @@
     Au = Mast*np.array(G2lat.U6toUij(USSdata[:6])).T   #atoms x waves x sin Uij mods
     Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T   #...cos Uij mods
-    GpA = np.array(expModInt(SSUniq,Af,Bf,Ax,Bx,Au,Bu,SGT))
+    GpA = np.array(expModInt(SSUniq,Af,Bf,Ax,Bx,Au,Bu))
     return GpA             # SGops x atoms
     
-def ModulationDerv(waveTypes,SSUniq,SSPhi,FSSdata,XSSdata,USSdata,SStauM,Mast):
-    import scipy.special as sp
-    Smult,TauT = SStauM             # both atoms x SGops
-    nh = np.zeros(1)
-    if XSSdata.ndim > 2:
-        nh = np.arange(XSSdata.shape[1])        
-    M = np.where(m>0,m+nh[:,np.newaxis],m-nh[:,np.newaxis])
-    A = np.array([[a,b] for a,b in zip(XSSdata[:3],XSSdata[3:])])
-    HdotA = twopi*(np.inner(A.T,SSUniq.T[:3].T)+SSPhi+TauT[:,np.newaxis,:])
-    Gpm = sp.jn(M-1,HdotA)
-    Gpp = sp.jn(M+1,HdotA)
-    if Gpm.ndim > 3: #sum over multiple harmonics
-        Gpm = np.sum(Gpm,axis=0)
-        Gpp = np.sum(Gpp,axis=0)
-    dGpdk = 0.5*(Gpm+Gpp)
-    return np.real(dGpdk),np.imag(dGpdk)
+def ModulationDerv(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast):
+    
+    nxs = np.newaxis
+    glTau,glWt = pwd.pygauleg(0.,1.,32)         #get Gauss-Legendre intervals & weights
+    
+    Ax = np.array(XSSdata[:3]).T   #atoms x waves x sin pos mods
+    Bx = np.array(XSSdata[3:]).T   #...cos pos mods
+    Af = np.array(FSSdata[0]).T    #sin frac mods x waves x atoms
+    Bf = np.array(FSSdata[1]).T    #cos frac mods...
+    Au = Mast*np.array(G2lat.U6toUij(USSdata[:6])).T   #atoms x waves x sin Uij mods
+    Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T   #...cos Uij mods
+
+    return dGpdk
     
 def posFourier(tau,psin,pcos,smul):

trunk/GSASIIstrMath.py

@@ -983 +983 @@
     waveTypes,FSSdata,XSSdata,USSdata,MSSdata = GetAtomSSFXU(pfx,calcControls,parmDict)
     SStauM = list(GetSSTauM(SGData['SGOps'],SSGData['SSGOps'],pfx,calcControls,Xdata))
-    SST = SSGT[:,3]
     if SGData['SGInv']:
         SStauM[0] = np.hstack((SStauM[0],SStauM[0]))
         SStauM[1] = np.hstack((SStauM[1],SStauM[1]))
-        SST = np.hstack((SST,-SST))
     modQ = np.array([parmDict[pfx+'mV0'],parmDict[pfx+'mV1'],parmDict[pfx+'mV2']])
     FF = np.zeros(len(Tdata))
@@ -1036 +1034 @@
             SSPhi = np.hstack((SSPhi,-SSPhi))
 #        GSASIIpath.IPyBreak()
-        GfpuA = G2mth.Modulation(waveTypes,SSUniq,SST,FSSdata,XSSdata,USSdata,Mast)
-        phase = twopi*(np.inner(Uniq,(dXdata.T+Xdata.T))+Phi[:,np.newaxis])
+        GfpuA = G2mth.Modulation(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast)
+        phase = twopi*(np.inner(Uniq,(dXdata.T+Xdata.T))+SSPhi[:,np.newaxis])
         sinp = np.sin(phase)
         cosp = np.cos(phase)
@@ -1103 +1101 @@
         Phi = np.inner(H[:3],SGT)
         SSPhi = np.inner(H,SSGT)
-        GfpuA,GfpuB = G2mth.Modulation(waveTypes,SSUniq,SSPhi,FSSdata,XSSdata,USSdata,SStauM,Mast)
-        dGAdk,dGBdk = G2mth.ModulationDerv(waveTypes,SSUniq,SSPhi,FSSdata,XSSdata,USSdata,SStauM,Mast)
+        GfpuA = G2mth.Modulation(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast)
+        dGAdk = G2mth.ModulationDerv(waveTypes,SSUniq,FSSdata,XSSdata,USSdata,Mast)
         #need ModulationDerv here dGAdXsin, etc  
-        phase = twopi*(np.inner((dXdata.T+Xdata.T),Uniq)+Phi[np.newaxis,:])
+        phase = twopi*(np.inner((dXdata.T+Xdata.T),Uniq)+SSPhi[np.newaxis,:])
         sinp = np.sin(phase)
         cosp = np.cos(phase)
@@ -1119 +1117 @@
         fot = (FF+FP-Bab)*occ*Tcorr
         fotp = FPP*occ*Tcorr
-        fa = np.array([fot[:,np.newaxis]*cosp,fotp[:,np.newaxis]*cosp])       #non positions
-        fb = np.array([fot[:,np.newaxis]*sinp,-fotp[:,np.newaxis]*sinp])
-        fa = fa*GfpuA[:,:,np.newaxis]-fb*GfpuB[:,:,np.newaxis]
-        fb = fb*GfpuA[:,:,np.newaxis]+fa*GfpuB[:,:,np.newaxis]
+        fa *= GfpuA
+        fb *= GfpuA       
         
         fas = np.sum(np.sum(fa,axis=1),axis=1)