Changeset 3777

Timestamp:

Jan 9, 2019 10:51:24 AM (5 years ago)

Author:

vondreele

Message:

fix wx.SAVE --> wx.FD_SAVE bug
work on incommensurate magnetic structure factor calcs. close but no cigar yet

Location:

trunk

Files:

• GSASIIdataGUI.py (modified) (1 diff)
• GSASIImath.py (modified) (3 diffs)
• GSASIIstrMath.py (modified) (12 diffs)

trunk/GSASIIdataGUI.py

@@ -3490 +3490 @@
                         dlg = wx.FileDialog(self, 'Choose sum image filename', pth,outname.split('IMG ')[1], 
                             'G2img files (*.G2img)|*.G2img', 
-                            wx.SAVE|wx.FD_OVERWRITE_PROMPT)
+                            wx.FD_SAVE|wx.FD_OVERWRITE_PROMPT)
                         if dlg.ShowModal() == wx.ID_OK:
                             newimagefile = dlg.GetPath()

trunk/GSASIImath.py

@@ -1346 +1346 @@
         Umod = 1.0
     if nWaves[3]:
-        tauM = np.arange(1.,nWaves[1]+1-nx)[:,nxs]*glTau  #Mwaves x ngl
+        tauM = np.arange(1.,nWaves[3]+1-nx)[:,nxs]*glTau  #Mwaves x ngl
         MmodA = Am[:,:,:,nxs]*np.sin(twopi*tauM[nxs,:,nxs,:]) #atoms X waves X 3 X ngl
         MmodB = Bm[:,:,:,nxs]*np.cos(twopi*tauM[nxs,:,nxs,:]) #ditto
@@ -1359 +1359 @@
     H: array nRefBlk x ops X hklt
     HP: array nRefBlk x ops X hklt proj to hkl
+    nWaves: list number of waves for frac, pos, uij & mag
     Fmod: array 2 x atoms x waves    (sin,cos terms)
     Xmod: array atoms X 3 X ngl
@@ -1366 +1367 @@
     '''
     
-    if nWaves[2]:
+    if nWaves[2]:       #uij (adp) waves
         if len(HP.shape) > 2:
             HbH = np.exp(-np.sum(HP[:,:,nxs,nxs,:]*np.inner(HP,Umod),axis=-1)) # refBlk x ops x atoms x ngl add Overhauser corr.?

trunk/GSASIIstrMath.py

@@ -1412 +1412 @@
     '''
     phfx = pfx.split(':')[0]+hfx
+    g = nl.inv(G)
     ast = np.sqrt(np.diag(G))
+    ainv = np.sqrt(np.diag(g))
+    GS = G/np.outer(ast,ast)
+    Ginv = g/np.outer(ainv,ainv)
+    VGi = np.sqrt(nl.det(Ginv))
+    uAmat = G2lat.Gmat2AB(GS)[0]
     Mast = twopisq*np.multiply.outer(ast,ast)    
     SGInv = SGData['SGInv']
@@ -1431 +1437 @@
     if not Xdata.size:          #no atoms in phase!
         return
-
-    if parmDict[pfx+'isMag']:       #TODO: fix the math
-        Mag = np.sqrt(np.sum(Gdata**2,axis=0))      #magnitude of moments for uniq atoms
-        Gdata = np.where(Mag>0.,Gdata/Mag,0.)       #normalze mag. moments
-        Gdata = np.inner(Gdata.T,SGMT).T            #apply sym. ops.
-        if SGData['SGInv'] and not SGData['SGFixed']:
-            Gdata = np.hstack((Gdata,-Gdata))       #inversion if any
-        Gdata = np.hstack([Gdata for icen in range(Ncen)])        #dup over cell centering
-        Gdata = SGData['MagMom'][nxs,:,nxs]*Gdata   #flip vectors according to spin flip * det(opM)
-        Mag = np.tile(Mag[:,nxs],len(SGMT)*Ncen).T
-        if SGData['SGInv'] and not SGData['SGFixed']:
-            Mag = np.repeat(Mag,2,axis=0)                  #Mag same shape as Gdata
-
-
     waveTypes,FSSdata,XSSdata,USSdata,MSSdata = GetAtomSSFXU(pfx,calcControls,parmDict)
     ngl,nWaves,Fmod,Xmod,Umod,Mmod,glTau,glWt = G2mth.makeWaves(waveTypes,FSSdata,XSSdata,USSdata,MSSdata,Mast)
     modQ = np.array([parmDict[pfx+'mV0'],parmDict[pfx+'mV1'],parmDict[pfx+'mV2']])
+
+    if parmDict[pfx+'isMag']:       #TODO: fix the math
+        GSdata = Gdata[:,nxs,:]+Mmod.T 
+        GSdata = np.inner(GSdata.T,SGMT).T            #apply sym. ops.
+        if SGData['SGInv'] and not SGData['SGFixed']:
+            GSdata = np.hstack((GSdata,-GSdata))       #inversion if any
+        GSdata = np.hstack([GSdata for icen in range(Ncen)])        #dup over cell centering
+        GSdata = SGData['MagMom'][nxs,:,nxs,nxs]*GSdata   #flip vectors according to spin flip * det(opM)
+        SMag = np.sqrt(np.sum((np.inner(GSdata.T,Ginv)*GSdata.T),axis=-1)).T
+        Kdata = np.inner(GSdata.T,uAmat).T*VGi/SMag[nxs,:,:,:]     #Cartesian unit vectors = 0.9626 for hexagonal???
+
     FF = np.zeros(len(Tdata))
     if 'NC' in calcControls[hfx+'histType']:
@@ -1456 +1459 @@
     Uij = np.array(G2lat.U6toUij(Uijdata)).T
     bij = Mast*Uij
-    blkSize = 32       #no. of reflections in a block
+    blkSize = 48       #no. of reflections in a block
     nRef = refDict['RefList'].shape[0]
     SQ = 1./(2.*refDict['RefList'].T[5])**2
@@ -1511 +1514 @@
         Tcorr = np.reshape(Tiso,Tuij.shape)*Tuij*Mdata*Fdata/Uniq.shape[1]  #refBlk x ops x atoms
         GfpuA = G2mth.Modulation(Uniq,UniqP,nWaves,Fmod,Xmod,Umod,Mmod,glTau,glWt) #2 x refBlk x sym X atoms
-        fams = np.zeros(32)
-        fbms = np.zeros(32)
+        fams = 0.
+        fbms = 0.
 
         if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:       #TODO: mag math here??
@@ -1525 +1528 @@
             sinm = np.sin(mphase)                               #ditto - match magstrfc.for
             cosm = np.cos(mphase)                               #ditto
+            
             HM = np.inner(Bmat,HP.T)                             #put into cartesian space
             HM = HM/np.sqrt(np.sum(HM**2,axis=0))               #Gdata = MAGS & HM = UVEC in magstrfc.for both OK
-            eDotK = np.sum(HM[:,:,nxs,nxs]*Gdata[:,nxs,:,:],axis=0)
-            Q = HM[:,:,nxs,nxs]*eDotK[nxs,:,:,:]-Gdata[:,nxs,:,:] #xyz,Nref,Nop,Natm = BPM in magstrfc.for OK
-            fam = Q*TMcorr[nxs,:,nxs,:]*cosm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
-            fbm = Q*TMcorr[nxs,:,nxs,:]*sinm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
-            fagm = fam*GfpuA[0]-fbm*GfpuA[1]   #real; 2 x refBlk x sym x atoms
-            fbgm = fbm*GfpuA[0]+fam*GfpuA[1]
-            fams = np.sum(np.sum(fagm,axis=-1),axis=-1)                          #xyz,Nref
-            fbms = np.sum(np.sum(fbgm,axis=-1),axis=-1)                          #ditto
-
+            eDotK = np.sum(HM[:,:,nxs,nxs,nxs]*Kdata[:,nxs,:,:,:],axis=0)
+            Q = HM[:,:,nxs,nxs,nxs]*eDotK[nxs,:,:,:,:]-Kdata[:,nxs,:,:,:] #Mxyz,Nref,Nop,Ntau,Natm
+            
+            fam = Q*TMcorr[nxs,:,nxs,nxs,:]*cosm[nxs,:,:,nxs,:]*SMag[nxs,nxs,:,:,:]    #ditto
+            fbm = Q*TMcorr[nxs,:,nxs,nxs,:]*sinm[nxs,:,:,nxs,:]*SMag[nxs,nxs,:,:,:]    #ditto
+            
+            fams = np.sum(np.sum(fam,axis=2),axis=-1)                          #xyz,Nref,ntau
+            fbms = np.sum(np.sum(fbm,axis=2),axis=-1)                          #ditto
+            
+            fams = np.sum(np.sum(fams**2,axis=0)*glWt[nxs,:],axis=-1)
+            fbms = np.sum(np.sum(fbms**2,axis=0)*glWt[nxs,:],axis=-1)
 
         if 'T' in calcControls[hfx+'histType']:
@@ -1547 +1553 @@
         fas = np.sum(np.sum(fag,axis=-1),axis=-1)   #2 x refBlk; sum sym & atoms
         fbs = np.sum(np.sum(fbg,axis=-1),axis=-1)
-        refl.T[10] = np.sum(fas,axis=0)**2+np.sum(fbs,axis=0)**2+np.sum(fams,axis=0)**2+np.sum(fbms,axis=0)**2    #square of sums
+        refl.T[10] = np.sum(fas,axis=0)**2+np.sum(fbs,axis=0)**2+fams+fbms    #square of sums
         refl.T[11] = atan2d(fbs[0],fas[0])  #ignore f' & f"
         if 'P' not in calcControls[hfx+'histType']:
@@ -1676 +1682 @@
         Tuij = np.where(HbH<1.,np.exp(HbH),1.0)
         Tcorr = np.reshape(Tiso,Tuij.shape)*Tuij*Mdata*Fdata/Uniq.shape[1]  #refBlk x ops x atoms
-#        GSASIIpath.IPyBreak()
         if 'T' in calcControls[hfx+'histType']:
             fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-np.reshape(Flack*FPP,sinp.shape)*sinp*Tcorr])
@@ -1758 +1763 @@
     dFdGf = np.zeros((nRef,mSize,FSSdata.shape[1],2))
     dFdGx = np.zeros((nRef,mSize,XSSdata.shape[1],6))
-#    dFdGz = np.zeros((nRef,mSize,5))
     dFdGu = np.zeros((nRef,mSize,USSdata.shape[1],12))
     Flack = 1.0
@@ -1910 +1914 @@
             dFdvDict[pfx+'U23cos:'+str(i)+':'+str(j)] = 2.*dFdGu.T[11][j][i]
             
-#        GSASIIpath.IPyBreak()
     dFdvDict[phfx+'Flack'] = 4.*dFdfl.T
     dFdvDict[phfx+'BabA'] = dFdbab.T[0]
@@ -1942 +1945 @@
         if parm == '0':
             continue
-        if parm.split(':')[2] in ['Tmin','Tmax','Xmax','Ymax','Zmax','Fzero','Fwid',]:
+        if parm.split(':')[2] in ['Tmin','Tmax','Xmax','Ymax','Zmax','Fzero','Fwid',
+            'MXsin','MXcos','MYsin','MYcos','MZsin','MZcos','AMx','AMy','AMz',]:
             parmDict[parm] += dM
             prefList = SStructureFactor(trefDict,G,hfx,pfx,SGData,SSGData,calcControls,parmDict)
@@ -3326 +3330 @@
                     aname = name.split(pfx)[1][:2]
                     if aname not in ['Af','dA','AU','RB','AM','Xs','Xc','Ys','Yc','Zs','Zc',    \
-                        'Tm','Xm','Ym','Zm','U1','U2','U3']: continue # skip anything not an atom or rigid body param
+                        'Tm','Xm','Ym','Zm','U1','U2','U3','MX','MY','MZ']: continue # skip anything not an atom or rigid body param
                 except IndexError:
                     continue
@@ -3336 +3340 @@
                     aname = name.split(pfx)[1][:2]
                     if aname not in ['Af','dA','AU','RB','AM','Xs','Xc','Ys','Yc','Zs','Zc',    \
-                        'Tm','Xm','Ym','Zm','U1','U2','U3']: continue # skip anything not an atom or rigid body param
+                        'Tm','Xm','Ym','Zm','U1','U2','U3','MX','MY','MZ']: continue # skip anything not an atom or rigid body param
                 except IndexError:
                     continue