Changeset 2484


Ignore:
Timestamp:
Oct 4, 2016 3:11:45 PM (7 years ago)
Author:
vondreele
Message:

fix transform to properly move magnetic moments (if space group is not changed)
some work on magnetic derivatives

Location:
trunk
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/GSASIIlattice.py

    r2482 r2484  
    246246   
    247247    cx,ct,cs,cia = oldPhase['General']['AtomPtrs']
     248    cm = 0
     249    if oldPhase['General']['Type'] == 'magnetic':
     250        cm = cx+4
     251    oAmat,oBmat = cell2AB(oldPhase['General']['Cell'][1:7])
     252    nAmat,nBmat = cell2AB(newPhase['General']['Cell'][1:7])
    248253    SGData = newPhase['General']['SGData']
    249254    invTrans = nl.inv(Trans)
    … …  
    275280        atom[cs:cs+2] = G2spc.SytSym(atom[cx:cx+3],SGData)[:2]
    276281        atom[cia+8] = ran.randint(0,sys.maxint)
     282        if cm:
     283            mag = np.sqrt(np.sum(np.array(atom[cm:cm+3])**2))
     284            mom = np.inner(np.array(atom[cm:cm+3]),oBmat)
     285            mom = np.inner(mom,invTrans.T)
     286            mom = np.inner(mom,nAmat)
     287            mom /= np.sqrt(np.sum(mom**2))
     288            atom[cm:cm+3] = mom*mag
    277289    newPhase['Atoms'] = newAtoms
    278290    newPhase['Atoms'] = GetUnique(newPhase)
    … …  
    284296    atomData = []
    285297    SGData = Phase['General']['SGData']
     298    SpnFlp = SGData.get('SpnFlp',[])
     299    Amat,Bmat = cell2AB(Phase['General']['Cell'][1:7])
    286300    cx,ct,cs,cia = Phase['General']['AtomPtrs']
     301    cm = 0
     302    if Phase['General']['Type'] == 'magnetic':
     303        cm = cx+4
    287304    unit = np.zeros(3)
    288305    for atom in Atoms:
    … …  
    297314                atom[cx:cx+3] = item[0]
    298315                atom[cia+2:cia+8] = item[1]
     316                if cm:
     317                    Opr = abs(item[2])%100
     318                    M = SGData['SGOps'][Opr-1][0]
     319                    opNum = G2spc.GetOpNum(item[2],SGData)
     320                    mom = np.inner(np.array(atom[cm:cm+3]),Bmat)
     321                    atom[cm:cm+3] = np.inner(np.inner(mom,M),Amat)*nl.det(M)*SpnFlp[opNum-1]
    299322                atomData.append(atom[:cia+9])  #not SS stuff
    300323        else:
    … …  
    303326                if item[0][2] >= .95: item[0][2] -= 1.
    304327                atom[cx:cx+3] = item[0]
     328                if cm:
     329                    Opr = abs(item[1])%100
     330                    M = SGData['SGOps'][Opr-1][0]
     331                    opNum = G2spc.GetOpNum(item[1],SGData)
     332                    mom = np.inner(np.array(atom[cm:cm+3]),Bmat)
     333                    atom[cm:cm+3] = np.inner(np.inner(mom,M),Amat)*nl.det(M)*SpnFlp[opNum-1]
     334#                GSASIIpath.IPyBreak()
    305335                atomData.append(atom[:cia+9])  #not SS stuff
     336           
    306337    return atomData
    307338       

  • trunk/GSASIIplot.py

    r2480 r2484  
    52985298        xyz = np.array([x,y,z])
    52995299        glEnable(GL_COLOR_MATERIAL)
    5300         glLineWidth(1)
     5300        glLineWidth(2)
     5301        glEnable(GL_BLEND)
     5302        glBlendFunc(GL_SRC_ALPHA,GL_ONE_MINUS_SRC_ALPHA)
     5303        glEnable(GL_LINE_SMOOTH)
    53015304        glPushMatrix()
    53025305        glTranslate(x,y,z)
    … …  
    53105313        glPopMatrix()
    53115314        glColor4ubv([0,0,0,0])
     5315        glDisable(GL_LINE_SMOOTH)
     5316        glDisable(GL_BLEND)
    53125317        glDisable(GL_COLOR_MATERIAL)
    53135318       

  • trunk/GSASIIstrMath.py

    r2483 r2484  
    766766        FF = np.repeat(refDict['FF']['FF'][iBeg:iFin].T[Tindx].T,len(SGT)*len(TwinLaw),axis=0)
    767767        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
    768             MF = np.repeat(refDict['FF']['MF'][iBeg:iFin].T[Tindx].T,len(TwinLaw),axis=0)   #Nref,Natm
    769             TMcorr = 0.5*0.539*Tcorr[:,0,:]*MF*len(SGMT)/Mdata                                  #Nref,Natm
     768            MF = refDict['FF']['MF'][iBeg:iFin]   #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))
    … …  
    773773                mphase = phase
    774774            mphase = np.array([mphase+twopi*np.inner(cen,H)[:,nxs,nxs] for cen in SGData['SGCen']])
    775             mphase = np.concatenate(mphase,axis=1)              #Nref,Nop,Natm
     775            mphase = np.concatenate(mphase,axis=1)              #Nref,full Nop,Natm
    776776            sinm = np.sin(mphase)                               #ditto - match magstrfc.for
    777777            cosm = np.cos(mphase)                               #ditto
    … …  
    10971097        iFin = min(iBeg+blkSize,nRef)
    10981098        refl = refDict['RefList'][iBeg:iFin]    #array(blkSize,nItems)
    1099         H = refl.T[:3]
     1099        H = refl.T[:3].T
    11001100        SQ = 1./(2.*refl.T[4])**2             # or (sin(theta)/lambda)**2
    11011101        SQfactor = 8.0*SQ*np.pi**2
    … …  
    11121112        FF = np.repeat(refDict['FF']['FF'][iBeg:iFin].T[Tindx].T,len(SGT),axis=0)
    11131113#        GSASIIpath.IPyBreak()
    1114         Uniq = np.inner(H.T,SGMT)             # array(nSGOp,3)
    1115         Phi = np.inner(H.T,SGT)
     1114        Uniq = np.inner(H,SGMT)             # array(nSGOp,3)
     1115        Phi = np.inner(H,SGT)
    11161116        phase = twopi*(np.inner(Uniq,(dXdata+Xdata).T).T+Phi.T).T
    11171117        sinp = np.sin(phase)        #refBlk x nOps x nAtoms
    … …  
    11491149            fams = np.sum(np.sum(fam,axis=-1),axis=-1)                          #xyz,Nref
    11501150            fbms = np.sum(np.sum(fbm,axis=-1),axis=-1)                          #ditto
     1151            Hij = np.hstack((Hij,Hij))
    11511152#            GSASIIpath.IPyBreak()
    11521153        if 'T' in calcControls[hfx+'histType']:
    … …  
    11601161        fax = np.array([-fot*sinp,-fotp*cosp])   #positions array(2,refBlk,nEqv,nAtoms)
    11611162        fbx = np.array([fot*cosp,-fotp*sinp])
    1162         #sum below is over Uniq
    1163         dfadfr = np.sum(fa/occ,axis=-2)        #array(2,refBlk,nAtom) Fdata != 0 avoids /0. problem
    1164         dfadba = np.sum(-cosp*Tcorr,axis=-2)  #array(refBlk,nAtom)
    1165         dfadx = np.sum(twopi*Uniq[nxs,:,nxs,:,:]*np.swapaxes(fax,-2,-1)[:,:,:,:,nxs],axis=-2)
    1166         dfadui = np.sum(-SQfactor[nxs,:,nxs,nxs]*fa,axis=-2) #array(Ops,refBlk,nAtoms)
    1167         dfadua = np.sum(-Hij[nxs,:,nxs,:,:]*np.swapaxes(fa,-2,-1)[:,:,:,:,nxs],axis=-2)
     1163        #sum below is over Uniq
     1164        if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
     1165            dfadfr = np.sum(fam/occ,axis=-2)        #array(refBlk,nAtom) Fdata != 0 avoids /0. problem
     1166            dfadba = np.sum(-cosm*TMcorr[:,nxs,:],axis=-2)    #array(refBlk,nAtom)
     1167            dfadx = np.sum(twopi*Uniq[:,nxs,:,:]*np.swapaxes(fax[0],-2,-1)[:,:,:,nxs],axis=-2)
     1168            dfadui = np.sum(-SQfactor[:,nxs,nxs]*fam,axis=-2) #array(Ops,refBlk,nAtoms)
     1169            dfadua = np.sum(-Hij[nxs,:,:,nxs,:]*fam[:,:,:,:,nxs],axis=-3)
     1170        else:
     1171            dfadfr = np.sum(fa/occ,axis=-2)        #array(2,refBlk,nAtom) Fdata != 0 avoids /0. problem
     1172            dfadba = np.sum(-cosp*Tcorr,axis=-2)  #array(refBlk,nAtom)
     1173            dfadx = np.sum(twopi*Uniq[nxs,:,nxs,:,:]*np.swapaxes(fax,-2,-1)[:,:,:,:,nxs],axis=-2)
     1174            dfadui = np.sum(-SQfactor[nxs,:,nxs,nxs]*fa,axis=-2) #array(Ops,refBlk,nAtoms)
     1175            dfadua = np.sum(-Hij[nxs,:,nxs,:,:]*np.swapaxes(fa,-2,-1)[:,:,:,:,nxs],axis=-2)
    11681176        # array(2,refBlk,nAtom,3) & array(2,refBlk,nAtom,6)
    11691177        if not SGData['SGInv']:
    1170             dfbdfr = np.sum(fb/occ,axis=-2)        #Fdata != 0 avoids /0. problem
    1171             dfbdba = np.sum(-sinp*Tcorr,axis=-2)
    1172             dfadfl = np.sum(np.sum(-fotp*sinp,axis=-1),axis=-1)
    1173             dfbdfl = np.sum(np.sum(fotp*cosp,axis=-1),axis=-1)
    1174             dfbdx = np.sum(twopi*Uniq[nxs,:,nxs,:,:]*np.swapaxes(fbx,-2,-1)[:,:,:,:,nxs],axis=-2)           
    1175             dfbdui = np.sum(-SQfactor[nxs,:,nxs,nxs]*fb,axis=-2)
    1176             dfbdua = np.sum(-Hij[nxs,:,nxs,:,:]*np.swapaxes(fb,-2,-1)[:,:,:,:,nxs],axis=-2)
     1178            if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
     1179                dfbdfr = np.sum(fbm/occ,axis=-2)        #array(refBlk,nAtom) Fdata != 0 avoids /0. problem
     1180                dfbdba = np.sum(-sinm*TMcorr[:,nxs,:],axis=-2)    #array(refBlk,nAtom)
     1181                dfbdx = np.sum(twopi*Uniq[:,nxs,:,:]*np.swapaxes(fbx[0],-2,-1)[:,:,:,nxs],axis=-2)
     1182                dfbdui = np.sum(-SQfactor[:,nxs,nxs]*fbm,axis=-2) #array(Ops,refBlk,nAtoms)
     1183                dfbdua = np.sum(-Hij[nxs,:,:,nxs,:]*fbm[:,:,:,:,nxs],axis=-3)
     1184            else:
     1185                dfbdfr = np.sum(fb/occ,axis=-2)        #Fdata != 0 avoids /0. problem
     1186                dfbdba = np.sum(-sinp*Tcorr,axis=-2)
     1187                dfadfl = np.sum(np.sum(-fotp*sinp,axis=-1),axis=-1)
     1188                dfbdfl = np.sum(np.sum(fotp*cosp,axis=-1),axis=-1)
     1189                dfbdx = np.sum(twopi*Uniq[nxs,:,nxs,:,:]*np.swapaxes(fbx,-2,-1)[:,:,:,:,nxs],axis=-2)           
     1190                dfbdui = np.sum(-SQfactor[nxs,:,nxs,nxs]*fb,axis=-2)
     1191                dfbdua = np.sum(-Hij[nxs,:,nxs,:,:]*np.swapaxes(fb,-2,-1)[:,:,:,:,nxs],axis=-2)
    11771192        else:
    11781193            dfbdfr = np.zeros_like(dfadfr)
    … …  
    11871202        SB = fbs[0]+fbs[1]
    11881203        if 'P' in calcControls[hfx+'histType']: #checked perfect for centro & noncentro
    1189             dFdfr[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadfr+fbs[:,:,nxs]*dfbdfr,axis=0)*Mdata/len(SGMT)
    1190             dFdx[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadx+fbs[:,:,nxs,nxs]*dfbdx,axis=0)
    1191             dFdui[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadui+fbs[:,:,nxs]*dfbdui,axis=0)
    1192             dFdua[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadua+fbs[:,:,nxs,nxs]*dfbdua,axis=0)
     1204            if 'N' in calcControls[hfx+'histType'] and parmDict[pfx+'isMag']:
     1205                dFdfr[iBeg:iFin] = 2.*(fams[:,nxs]*dfadfr+fbms[:,nxs]*dfbdfr)*Mdata/len(SGMT)
     1206                dFdx[iBeg:iFin] = 2.*(fams[:,nxs,nxs]*dfadx+fbms[:,nxs,nxs]*dfbdx)
     1207                dFdui[iBeg:iFin] = 2.*(fams[:,nxs]*dfadui+fbms[:,nxs]*dfbdui)
     1208                dFdua[iBeg:iFin] = 2.*(fams[:,nxs,nxs]*dfadua+fbms[:,nxs,nxs]*dfbdua)
     1209            else:
     1210                dFdfr[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadfr+fbs[:,:,nxs]*dfbdfr,axis=0)*Mdata/len(SGMT)
     1211                dFdx[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadx+fbs[:,:,nxs,nxs]*dfbdx,axis=0)
     1212                dFdui[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs]*dfadui+fbs[:,:,nxs]*dfbdui,axis=0)
     1213                dFdua[iBeg:iFin] = 2.*np.sum(fas[:,:,nxs,nxs]*dfadua+fbs[:,:,nxs,nxs]*dfbdua,axis=0)
    11931214        else:
    11941215            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.