Changeset 2481


Ignore:
Timestamp:
Sep 30, 2016 10:33:12 AM (7 years ago)
Author:
vondreele
Message:

force Transform to delete nonmagnetic atoms if phase made magnetic & add 'mag' to new phase name
fix TOF cosine background function
magnetic structure refinement works (in numeric mode only)
magnetic structure factor calculation correct

Location:
trunk
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • trunk/GSASII.py

    r2479 r2481  
    36563656        rbVary,rbDict = G2stIO.GetRigidBodyModels(rigidbodyDict,Print=False)
    36573657        rbIds = rigidbodyDict.get('RBIds',{'Vector':[],'Residue':[]})
    3658         Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,maxSSwave = G2stIO.GetPhaseData(Phases,RestraintDict=None,rbIds=rbIds,Print=False)       
     3658        Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,MFtable,maxSSwave = G2stIO.GetPhaseData(Phases,RestraintDict=None,rbIds=rbIds,Print=False)       
    36593659        hapVary,hapDict,controlDict = G2stIO.GetHistogramPhaseData(Phases,Histograms,Print=False)
    36603660        histVary,histDict,controlDict = G2stIO.GetHistogramData(Histograms,Print=False)

  • trunk/GSASIIconstrGUI.py

    r2480 r2481  
    239239
    240240    # create a list of the phase variables
    241     Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,maxSSwave = G2stIO.GetPhaseData(Phases,rbIds=rbIds,Print=False)
     241    Natoms,atomIndx,phaseVary,phaseDict,pawleyLookup,FFtable,BLtable,MFtable,maxSSwave = G2stIO.GetPhaseData(Phases,rbIds=rbIds,Print=False)
    242242    phaseList = []
    243243    for item in phaseDict:

  • trunk/GSASIIlattice.py

    r2480 r2481  
    3939import GSASIImath as G2mth
    4040import GSASIIspc as G2spc
     41import GSASIIElem as G2elem
    4142GSASIIpath.SetVersionNumber("$Revision$")
    4243# trig functions in degrees
    … …  
    238239    :param Trans: array transformation matrix
    239240    :param Vec: array transformation vector
    240     :param ifMag: bool True if convert to magnetic phase
     241    :param ifMag: bool True if convert to magnetic phase;
     242        if True all nonmagnetic atoms will be removed
    241243    '''
    242244   
    … …  
    259261        newPhase['General']['AtomPtrs'] = [cx,ct,cs,cia]
    260262        magAtoms = []
     263        Landeg = 2.0
    261264        for atom in newAtoms:
    262             magAtoms.append(atom[:cx+4]+[0.,0.,0.]+atom[cx+4:])
     265            if len(G2elem.GetMFtable([atom[ct],],[Landeg,])):
     266                magAtoms.append(atom[:cx+4]+[0.,0.,0.]+atom[cx+4:])
    263267        newAtoms = magAtoms
    264268        newPhase['Draw Atoms'] = []

  • trunk/GSASIIphsGUI.py

    r2480 r2481  
    13831383            dlg.Destroy()
    13841384        phaseName = newPhase['General']['Name']
     1385        if ifMag:
     1386            phaseName += ' mag'
    13851387        newPhase = G2lat.TransformPhase(data,newPhase,Trans,Vec,ifMag)
    13861388        newPhase['General']['Map'] = mapDefault.copy()

  • trunk/GSASIIpwd.py

    r2473 r2481  
    622622            key = pfx+'Back;'+str(iBak)
    623623            if bakType == 'chebyschev':
    624                 ybi = parmDict[key]*(2.*(xdata-xdata[0])/dt-1.)**iBak
     624                ybi = parmDict[key]*(-1.+2.*(xdata-xdata[0])/dt)**iBak
    625625            elif bakType == 'cosine':
    626                 ybi = parmDict[key]*npcosd(xdata*iBak)
     626                ybi = parmDict[key]*npcosd(180.*xdata*iBak/xdata[-1])
    627627            yb += ybi
    628628        sumBk[0] = np.sum(yb)
    … …  
    744744        for iBak in range(nBak):   
    745745            if bakType == 'chebyschev':
    746                 dydb[iBak] = (2.*(xdata-xdata[0])/dt-1.)**iBak
     746                dydb[iBak] = (-1.+2.*(xdata-xdata[0])/dt)**iBak
    747747            elif bakType == 'cosine':
    748                 dydb[iBak] = npcosd(xdata*iBak)
     748                dydb[iBak] = npcosd(180.*xdata*iBak/xdata[-1])
    749749    elif bakType in ['Q^2 power series','Q^-2 power series']:
    750750        QT = 1.

  • trunk/GSASIIstrIO.py

    r2478 r2481  
    10991099        BLtables.update(BLtable)
    11001100        phaseDict[pfx+'isMag'] = False
     1101        SGData = General['SGData']
     1102        SGtext,SGtable = G2spc.SGPrint(SGData)
    11011103        if General['Type'] == 'magnetic':
    11021104            MFtable = G2el.GetMFtable(General['AtomTypes'],General['Lande g'])
    11031105            MFtables.update(MFtable)
    11041106            phaseDict[pfx+'isMag'] = True
     1107            SpnFlp = SGData['SpnFlp']
    11051108        Atoms = PhaseData[name]['Atoms']
    11061109        if Atoms and not General.get('doPawley'):
    … …  
    11081111            AtLookup = G2mth.FillAtomLookUp(Atoms,cia+8)
    11091112        PawleyRef = PhaseData[name].get('Pawley ref',[])
    1110         SGData = General['SGData']
    1111         SGtext,SGtable = G2spc.SGPrint(SGData)
    11121113        cell = General['Cell']
    11131114        A = G2lat.cell2A(cell[1:7])
    … …  
    12081209                                G2mv.StoreEquivalence(name,equiv[1:])
    12091210                if 'M' in at[ct+1]:
    1210                     pass    #magnetic moment vary here
     1211                    SytSym,Mul,Nop = G2spc.SytSym(at[cx:cx+3],SGData)
     1212                    if SpnFlp[Nop] > 0.:    #black use p
     1213                        mId,mCoef = G2spc.GetCSpqinel(SytSym)[0]
     1214                    else:                   #red use q
     1215                        mId,mCoef = G2spc.GetCSpqinel(SytSym)[1]
     1216                    names = [pfx+'AMx:'+str(i),pfx+'AMy:'+str(i),pfx+'AMz:'+str(i)]
     1217                    equivs = [[],[],[]]
     1218                    for j in range(3):
     1219                        if mId[j] > 0:
     1220                            phaseVary.append(names[j])
     1221                            equivs[mId[j]-1].append([names[j],mCoef[j]])
     1222                    for equiv in equivs:
     1223                        if len(equiv) > 1:
     1224                            name = equiv[0][0]
     1225                            coef = equiv[0][1]
     1226                            for eqv in equiv[1:]:
     1227                                eqv[1] /= coef
     1228                                G2mv.StoreEquivalence(name,(eqv,))
    12111229                if General.get('Modulated',False):
    12121230                    AtomSS = at[-1]['SS1']
    … …  
    19691987                names = {cx:pfx+'Ax:'+str(i),cx+1:pfx+'Ay:'+str(i),cx+2:pfx+'Az:'+str(i),cx+3:pfx+'Afrac:'+str(i),
    19701988                    cia+1:pfx+'AUiso:'+str(i),cia+2:pfx+'AU11:'+str(i),cia+3:pfx+'AU22:'+str(i),cia+4:pfx+'AU33:'+str(i),
    1971                     cia+5:pfx+'AU12:'+str(i),cia+6:pfx+'AU13:'+str(i),cia+7:pfx+'AU23:'+str(i)}
     1989                    cia+5:pfx+'AU12:'+str(i),cia+6:pfx+'AU13:'+str(i),cia+7:pfx+'AU23:'+str(i),
     1990                    cx+4:pfx+'AMx:'+str(i),cx+5:pfx+'AMy:'+str(i),cx+6:pfx+'AMz:'+str(i)}
    19721991                for ind in range(cx,cx+4):
    19731992                    at[ind] = parmDict[names[ind]]
    … …  
    19842003                else:
    19852004                    for ind in range(cia+2,cia+8):
     2005                        at[ind] = parmDict[names[ind]]
     2006                        if names[ind] in sigDict:
     2007                            atomsSig[str(i)+':'+str(ind)] = sigDict[names[ind]]
     2008                if General['Type'] == 'magnetic':
     2009                    for ind in range(cx+4,cx+7):
    19862010                        at[ind] = parmDict[names[ind]]
    19872011                        if names[ind] in sigDict:

  • trunk/GSASIIstrMain.py

    r2478 r2481  
    200200    print >>printFile,135*'-'
    201201    try:
     202        covData = {}
    202203        IfOK,Rvals,result,covMatrix,sig = RefineCore(Controls,Histograms,Phases,restraintDict,
    203204            rigidbodyDict,parmDict,varyList,calcControls,pawleyLookup,ifPrint,printFile,dlg)

  • trunk/GSASIIstrMath.py

    r2480 r2481  
    767767        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
    768768            MF = np.repeat(refDict['FF']['MF'][iBeg:iFin].T[Tindx].T,len(TwinLaw),axis=0)   #Nref,Natm
    769             TMcorr = 0.539*Tcorr[:,0,:]*MF                                                  #Nref,Natm
     769            TMcorr = 0.5*0.539*Tcorr[:,0,:]*MF*len(SGMT)/Mdata                                  #Nref,Natm
    770770            if SGData['SGInv']:
    771771                mphase = np.hstack((phase,-phase))
    … …  
    780780            eDotK = np.sum(HM[:,:,nxs,nxs]*Gdata[:,nxs,:,:],axis=0)
    781781            Q = -HM[:,:,nxs,nxs]*eDotK[nxs,:,:,:]+Gdata[:,nxs,:,:] #xyz,Nref,Nop,Natm = BPM in magstrfc.for OK
    782             fam = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*cosm[nxs,:,:,:]    #ditto
    783             fbm = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*sinm[nxs,:,:,:]    #ditto
     782            fam = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*cosm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
     783            fbm = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*sinm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
    784784            fams = np.sum(np.sum(fam,axis=-1),axis=-1)                          #xyz,Nref
    785785            fbms = np.sum(np.sum(fbm,axis=-1),axis=-1)                          #ditto
    786             GSASIIpath.IPyBreak()
     786#            GSASIIpath.IPyBreak()
    787787        if 'T' in calcControls[hfx+'histType']: #fa,fb are 2 X blkSize X nTwin X nOps x nAtoms
    788788            fa = np.array([np.reshape(((FF+FP).T-Bab).T,cosp.shape)*cosp*Tcorr,-np.reshape(Flack*FPP,sinp.shape)*sinp*Tcorr])
    … …  
    800800            refl.T[10] = atan2d(fbs[0],fas[0])  #ignore f' & f"
    801801            if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
    802                 refl.T[9] += np.sum(fams**2,axis=0)+np.sum(fbms**2,axis=0)
     802                refl.T[9] = np.sum(fams**2,axis=0)+np.sum(fbms**2,axis=0)
    803803        else:                                       #HKLF: F^2 = (A[0]+A[1])^2 + (B[0]+B[1])^2
    804804            if len(TwinLaw) > 1:
    … …  
    811811                refl.T[7] = np.copy(refl.T[9])               
    812812                refl.T[10] = atan2d(fbs[0],fas[0])  #ignore f' & f"
     813                if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
     814                    refl.T[9] = np.sum(fams**2,axis=0)+np.sum(fbms**2,axis=0)
    813815#        GSASIIpath.IPyBreak()
    814816#                refl.T[10] = atan2d(np.sum(fbs,axis=0),np.sum(fas,axis=0)) #include f' & f"
    … …  
    10501052        GetAtomFXU(pfx,calcControls,parmDict)
    10511053    mSize = len(Mdata)
     1054    if parmDict[pfx+'isMag']:
     1055        Mag = np.sqrt(np.sum(Gdata**2,axis=0))      #magnitude of moments for uniq atoms
     1056        Gdata = np.where(Mag>0.,Gdata/Mag,0.)       #normalze mag. moments
     1057        Gdata = np.inner(Bmat,Gdata.T)              #convert to crystal space
     1058        Gdata = np.inner(Gdata.T,SGMT).T            #apply sym. ops.
     1059        if SGData['SGInv']:
     1060            Gdata = np.hstack((Gdata,-Gdata))       #inversion if any
     1061        Gdata = np.repeat(Gdata,Ncen,axis=1)    #dup over cell centering
     1062        Gdata = SGData['MagMom'][nxs,:,nxs]*Gdata   #flip vectors according to spin flip
     1063        Gdata = np.inner(Amat,Gdata.T)              #convert back to cart. space MXYZ, Natoms, NOps*Inv*Ncen
     1064        Gdata = np.swapaxes(Gdata,1,2)              # put Natoms last
     1065#        GSASIIpath.IPyBreak()
     1066        Mag = np.tile(Mag[:,nxs],len(SGMT)*Ncen).T
     1067        if SGData['SGInv']:
     1068            Mag = np.repeat(Mag,2,axis=0)                  #Mag same length as Gdata
     1069        dFdMx = np.zeros((nRef,mSize,3))
    10521070    FF = np.zeros(len(Tdata))
    10531071    if 'NC' in calcControls[hfx+'histType']:
    … …  
    11091127        else:
    11101128            fotp = FPP*Tcorr     
     1129        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
     1130            MF = np.repeat(refDict['FF']['MF'][iBeg:iFin].T[Tindx].T,len(TwinLaw),axis=0)   #Nref,Natm
     1131            TMcorr = 0.5*0.539*Tcorr[:,0,:]*MF*len(SGMT)/Mdata                                  #Nref,Natm
     1132            if SGData['SGInv']:
     1133                mphase = np.hstack((phase,-phase))
     1134            else:
     1135                mphase = phase
     1136            mphase = np.array([mphase+twopi*np.inner(cen,H)[:,nxs,nxs] for cen in SGData['SGCen']])
     1137            mphase = np.concatenate(mphase,axis=1)              #Nref,Nop,Natm
     1138            sinm = np.sin(mphase)                               #ditto - match magstrfc.for
     1139            cosm = np.cos(mphase)                               #ditto
     1140            HM = np.inner(Bmat.T,H)                             #put into cartesian space
     1141            HM = HM/np.sqrt(np.sum(HM**2,axis=0))               #Gdata = MAGS & HM = UVEC in magstrfc.for both OK
     1142            eDotK = np.sum(HM[:,:,nxs,nxs]*Gdata[:,nxs,:,:],axis=0)
     1143            Q = -HM[:,:,nxs,nxs]*eDotK[nxs,:,:,:]+Gdata[:,nxs,:,:] #xyz,Nref,Nop,Natm = BPM in magstrfc.for OK
     1144            fam = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*cosm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
     1145            fbm = Q*TMcorr[nxs,:,nxs,:]*SGData['MagMom'][nxs,nxs,:,nxs]*sinm[nxs,:,:,:]*Mag[nxs,nxs,:,:]    #ditto
     1146            fams = np.sum(np.sum(fam,axis=-1),axis=-1)                          #xyz,Nref
     1147            fbms = np.sum(np.sum(fbm,axis=-1),axis=-1)                          #ditto
     1148#            GSASIIpath.IPyBreak()
    11111149        if 'T' in calcControls[hfx+'histType']:
    11121150            fa = np.array([fot*cosp,-np.reshape(Flack*FPP,sinp.shape)*sinp*Tcorr])
    … …  
    11471185        if 'P' in calcControls[hfx+'histType']: #checked perfect for centro & noncentro
    11481186            dFdfr[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadfr+fbs[:,:,nxs]*dfbdfr,axis=0)*Mdata/len(SGMT)
    1149 #            dFdfr[iBeg:iFin] = 2.*(fas[0,:,nxs]*dfadfr[0]+fas[1,:,nxs]*dfadfr[1])*Mdata/len(SGMT)+   \
    1150 #                2.*(fbs[0,:,nxs]*dfbdfr[0]+fbs[1,:,nxs]*dfbdfr[1])*Mdata/len(SGMT)
    11511187            dFdx[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadx+fbs[:,:,nxs,nxs]*dfbdx,axis=0)
    1152 #            dFdx[iBeg:iFin] = 2.*(fas[0,:,nxs,nxs]*dfadx[0]+fas[1,:,nxs,nxs]*dfadx[1])+  \
    1153 #                2.*(fbs[0,:,nxs,nxs]*dfbdx[0]+fbs[1,:,nxs,nxs]*dfbdx[1])
    11541188            dFdui[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadui+fbs[:,:,nxs]*dfbdui,axis=0)
    1155 #            dFdui[iBeg:iFin] = 2.*(fas[0,:,nxs]*dfadui[0]+fas[1,:,nxs]*dfadui[1])+   \
    1156 #                2.*(fbs[0,:,nxs]*dfbdui[0]+fbs[1,:,nxs]*dfbdui[1])
    11571189            dFdua[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadua+fbs[:,:,nxs,nxs]*dfbdua,axis=0)
    1158 #            dFdua[iBeg:iFin] = 2.*(fas[0,:,nxs,nxs]*dfadua[0]+fas[1,:,nxs,nxs]*dfadua[1])+   \
    1159 #                2.*(fbs[0,:,nxs,nxs]*dfbdua[0]+fbs[1,:,nxs,nxs]*dfbdua[1])
    11601190        else:
    11611191            dFdfr[iBeg:iFin] = (2.*SA[:,nxs]*(dfadfr[0]+dfadfr[1])+2.*SB[:,nxs]*(dfbdfr[0]+dfbdfr[1]))*Mdata/len(SGMT)
Note: See TracChangeset for help on using the changeset viewer.