Ignore:
Timestamp:
Feb 6, 2019 2:11:30 PM (3 years ago)
Author:
vondreele
Message:

incommensurate mag str fctr. closer but still not right

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIstrMath.py

    r3808 r3812  
    11081108        fam = Q*TMcorr[nxs,:,nxs,:]*cosm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
    11091109        fbm = Q*TMcorr[nxs,:,nxs,:]*sinm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
    1110         fams = np.sum(np.sum(fam,axis=-1),axis=-1)                          #Mxyz,Nref
    1111         fbms = np.sum(np.sum(fbm,axis=-1),axis=-1)                          #ditto
    1112         refl.T[9] = np.sum(fams**2,axis=0)+np.sum(fbms**2,axis=0)
     1110        fams = np.sum(np.sum(fam,axis=-1),axis=-1)     #Mxyz,Nref  Sum(sum(fam,atoms),ops)
     1111        fbms = np.sum(np.sum(fbm,axis=-1),axis=-1)     #ditto
     1112        refl.T[9] = np.sum(fams**2,axis=0)+np.sum(fbms**2,axis=0)   #Sum(fams**2,Mxyz) Re + Im
    11131113        refl.T[7] = np.copy(refl.T[9])               
    11141114        refl.T[10] = atan2d(fbms[0],fams[0]) #- what is phase for mag refl?
     
    15111511    modQ = np.array([parmDict[pfx+'mV0'],parmDict[pfx+'mV1'],parmDict[pfx+'mV2']])
    15121512
    1513     if parmDict[pfx+'isMag']:       #TODO: fix the math
     1513    if parmDict[pfx+'isMag']:       #TODO: fix the math       
     1514#1st try at this
    15141515        GSdata = Gdata[nxs,nxs,:,:]+Mmod    #ReIm,ntau,Mxyz,Natm
    15151516        GSdata = np.inner(np.swapaxes(GSdata,2,3),np.swapaxes(SGMT,1,2)).T  #apply sym. ops.--> Mxyz,Nops,Natm,Ntau,ReIm
     
    15181519        GSdata = np.hstack([GSdata for icen in range(Ncen)])        #dup over cell centering
    15191520        GSdata = SGData['MagMom'][nxs,:,nxs,nxs,nxs]*GSdata   #flip vectors according to spin flip * det(opM)
    1520         SMag = np.sqrt(np.sum((np.inner(GSdata.T,Ginv)*GSdata.T),axis=-1)).T    #Nops,Natm,Ntau,ReIm
    15211521        Kdata = np.inner(GSdata.T,uAmat).T     #Cartesian unit vectors
    1522         Kmean = np.mean(np.sqrt(np.sum(np.sum(Kdata**2,axis=0),axis=-1)),axis=0)    #normalization --> unit vectors
    1523         Kdata /= Kmean.T      #mxyz,nops,natm,ntau,ReIm
     1522        SMag = np.sqrt(np.sum(np.sum(Kdata**2,axis=0),axis=-1))
     1523        Kmean = np.mean(SMag,axis=0)    #normalization --> unit vectors
     1524        Kdata /= Kmean[nxs,nxs,:,:,nxs]      #mxyz,nops,natm,ntau,ReIm
    15241525
    15251526    FF = np.zeros(len(Tdata))
     
    15851586        Tuij = np.where(HbH<1.,np.exp(HbH),1.0)
    15861587        Tcorr = np.reshape(Tiso,Tuij.shape)*Tuij*Mdata*Fdata/Uniq.shape[1]  #refBlk x ops x atoms
    1587         GfpuA = G2mth.Modulation(Uniq,UniqP,nWaves,Fmod,Xmod,Umod,glTau,glWt) #2 x refBlk x sym X atoms
    1588         if 'T' in calcControls[hfx+'histType']:
    1589             fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-np.reshape(Flack*FPP,sinp.shape)*sinp*Tcorr])
    1590             fb = np.array([np.reshape(Flack*FPP,cosp.shape)*cosp*Tcorr,np.reshape(((FF+FP).T-Bab).T,sinp.shape)*sinp*Tcorr])
    1591         else:
    1592             fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-Flack*FPP*sinp*Tcorr])
    1593             fb = np.array([Flack*FPP*cosp*Tcorr,np.reshape(((FF+FP).T-Bab).T,sinp.shape)*sinp*Tcorr])
    1594         fams = 0.
    1595         fbms = 0.
    15961588
    15971589        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:       #TODO: mag math here??
    15981590            MF = refDict['FF']['MF'][iBeg:iFin].T[Tindx].T   #Nref,Natm
    15991591            TMcorr = 0.539*(np.reshape(Tiso,Tuij.shape)*Tuij)[:,0,:]*Fdata*Mdata*MF/(2*Nops)     #Nref,Natm
     1592                     
    16001593            D = twopi*H.T[:,3:]*glTau[nxs,:]
    16011594            mphase = phase[:,:,nxs,:]+D[:,nxs,:,nxs]
     
    16101603            Q = HM[:,:,nxs,nxs,nxs,nxs]*eDotK[nxs,:,:,:,:,:]-Kdata[:,nxs,:,:,:,:] #Mxyz,Nref,Nop,Natm,Ntau,ReIm
    16111604
    1612 #flawed but "close"           
    1613             fam = (Q*TMcorr[nxs,:,nxs,:,nxs,nxs]*cosm[nxs,:,:,:,:,nxs]*SMag[nxs,nxs,:,:,:,:]).T   #Mxyz,Nref,Nop,Natm,Ntau,ReIm
    1614             fbm = (Q*TMcorr[nxs,:,nxs,:,nxs,nxs]*sinm[nxs,:,:,:,:,nxs]*SMag[nxs,nxs,:,:,:,:]).T  #--> ReIm,Ntau,Natm,Nop,Nref,Mxyzditto
     1605            fam = (Q*TMcorr[nxs,:,nxs,:,nxs,nxs]*cosm[nxs,:,:,:,:,nxs]*SMag[nxs,nxs,:,:,:,nxs])   #Mxyz,Nref,Nop,Natm,Ntau,ReIm
     1606            fbm = (Q*TMcorr[nxs,:,nxs,:,nxs,nxs]*sinm[nxs,:,:,:,:,nxs]*SMag[nxs,nxs,:,:,:,nxs])
    16151607           
    1616             famg = (fam[0]*fbm[0]-fam[1]*fbm[1]).T          #--> Mxyz,Nref,Nop,Natm,Ntau
    1617             fbmg = (fam[0]*fbm[1]+fam[1]*fbm[1]).T
     1608            fams = np.sum(np.sum(np.sum(fam,axis=-1),axis=2),axis=2)      #xyz,Nref,ntau; sum ops & atoms
     1609            fbms = np.sum(np.sum(np.sum(fbm,axis=-1),axis=2),axis=2)      #ditto
    16181610           
    1619             fams = np.sum(np.sum(famg,axis=2),axis=2)      #xyz,Nref,ntau; sum ops & atoms
    1620             fbms = np.sum(np.sum(fbmg,axis=2),axis=2)      #ditto
     1611            fas = np.sum(fams*glWt[nxs,nxs,:],axis=-1)
     1612            fbs = np.sum(fbms*glWt[nxs,nxs,:],axis=-1)
     1613
     1614        else:
     1615            GfpuA = G2mth.Modulation(Uniq,UniqP,nWaves,Fmod,Xmod,Umod,glTau,glWt) #2 x refBlk x sym X atoms
     1616            if 'T' in calcControls[hfx+'histType']:
     1617                fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-np.reshape(Flack*FPP,sinp.shape)*sinp*Tcorr])
     1618                fb = np.array([np.reshape(Flack*FPP,cosp.shape)*cosp*Tcorr,np.reshape(((FF+FP).T-Bab).T,sinp.shape)*sinp*Tcorr])
     1619            else:
     1620                fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-Flack*FPP*sinp*Tcorr])
     1621                fb = np.array([Flack*FPP*cosp*Tcorr,np.reshape(((FF+FP).T-Bab).T,sinp.shape)*sinp*Tcorr])
     1622            fag = fa*GfpuA[0]-fb*GfpuA[1]   #real; 2 x refBlk x sym x atoms
     1623            fbg = fb*GfpuA[0]+fa*GfpuA[1]
     1624            fas = np.sum(np.sum(fag,axis=-1),axis=-1)   #2 x refBlk; sum sym & atoms
     1625            fbs = np.sum(np.sum(fbg,axis=-1),axis=-1)
    16211626           
    1622             fams = np.sum(np.sum(fams**2,axis=0)*glWt,axis=-1)
    1623             fbms = np.sum(np.sum(fbms**2,axis=0)*glWt,axis=-1)
    1624 
    1625 #why doesn't this work?
    1626            
    1627 #            fam = Q*TMcorr[nxs,:,nxs,nxs,:]*cosm[nxs,:,:,nxs,:]*SMag[nxs,nxs,:,:,:]    #RI,Mxyz,Nref,Nop,Ntau,Natm
    1628 #            fbm = Q*TMcorr[nxs,:,nxs,nxs,:]*sinm[nxs,:,:,nxs,:]*SMag[nxs,nxs,:,:,:]    #ditto
    1629 #           
    1630 #            famg = fa[:,nxs,:,:,nxs,:]*fam[nxs,:,:,:,:]-fb[:,nxs,:,:,nxs,:]*fbm[nxs,:,:,:,:]
    1631 #            fbmg = fb[:,nxs,:,:,nxs,:]*fam[nxs,:,:,:,:]+fa[:,nxs,:,:,nxs,:]*fbm[nxs,:,:,:,:]
    1632 #           
    1633 #            fams = np.sum(np.sum(famg,axis=3),axis=-1)     #RI,xyz,Nref,ntau
    1634 #            fbms = np.sum(np.sum(fbmg,axis=3),axis=-1)     #ditto
    1635 #           
    1636 #            fams = np.sum(np.sum(np.sum(fams**2,axis=0),axis=0)*glWt[nxs,nxs,nxs,:],axis=-1)
    1637 #            fbms = np.sum(np.sum(np.sum(fbms**2,axis=0),axis=0)*glWt[nxs,nxs,nxs,:],axis=-1)
    1638 
    1639 #or this
    1640 
    1641 #            fam = TMcorr[:,nxs,:]*cosm    #ditto
    1642 #            fbm = TMcorr[:,nxs,:]*sinm    #ditto
    1643 #
    1644 #            cosQ = np.sum(np.cos(twopi*Q)*SMag[nxs,nxs,:,:,:]*glWt[nxs,nxs,nxs,:,nxs],axis=3)           
    1645 #            sinQ = np.sum(np.sin(twopi*Q)*SMag[nxs,nxs,:,:,:]*glWt[nxs,nxs,nxs,:,nxs],axis=3)
    1646 #           
    1647 #            fagm = fam[nxs,:,:,:]*cosQ-fbm[nxs,:,:,:]*sinQ
    1648 #            fbgm = fbm[nxs,:,:,:]*cosQ+fam[nxs,:,:,:]*sinQ
    1649 #           
    1650 #            fams = np.sum(np.sum(fagm,axis=-1),axis=-1)                          #xyz,Nref,ntau
    1651 #            fbms = np.sum(np.sum(fbgm,axis=-1),axis=-1)                          #ditto
    1652 #           
    1653 #            fams = np.sum(fams**2,axis=0)
    1654 #            fbms = np.sum(fbms**2,axis=0)
    1655 
    1656         fag = fa*GfpuA[0]-fb*GfpuA[1]   #real; 2 x refBlk x sym x atoms
    1657         fbg = fb*GfpuA[0]+fa*GfpuA[1]
    1658         fas = np.sum(np.sum(fag,axis=-1),axis=-1)   #2 x refBlk; sum sym & atoms
    1659         fbs = np.sum(np.sum(fbg,axis=-1),axis=-1)
    1660         refl.T[10] = np.sum(fas,axis=0)**2+np.sum(fbs,axis=0)**2+fams+fbms    #square of sums
     1627        refl.T[10] = np.sum(fas,axis=0)**2+np.sum(fbs,axis=0)**2    #square of sums
    16611628        refl.T[11] = atan2d(fbs[0],fas[0])  #ignore f' & f"
    16621629        if 'P' not in calcControls[hfx+'histType']:
Note: See TracChangeset for help on using the changeset viewer.