Changeset 1990 for trunk/GSASIImath.py


Ignore:
Timestamp:
Oct 7, 2015 1:06:38 PM (7 years ago)
Author:
vondreele
Message:

Structure factor derivs for SS continued...

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImath.py

    r1988 r1990  
    928928def Modulation(waveTypes,H,FSSdata,XSSdata,USSdata,Mast):
    929929    '''
    930     H: array ops X hklt
    931     Ax: array atoms X waves X xyz
    932     Bx: array atoms X waves X xyz
     930    H: array nRefBlk x ops X hklt
     931    FSSdata: array 2 x atoms x waves    (sin,cos terms)
     932    XSSdata: array 2x3 x atoms X waves (sin,cos terms)
     933    USSdata: array 2x6 x atoms X waves (sin,cos terms)
     934    Mast: array orthogonalization matrix for Uij
    933935    '''
    934936   
     
    987989    '''
    988990    H: array ops X hklt
    989     Ax: array atoms X waves X xyz
    990     Bx: array atoms X waves X xyz
     991    FSSdata: array 2 x atoms x waves    (sin,cos terms)
     992    XSSdata: array 2x3 x atoms X waves (sin,cos terms)
     993    USSdata: array 2x6 x atoms X waves (sin,cos terms)
     994    Mast: array orthogonalization matrix for Uij
    991995    '''
    992996   
     
    10201024            nS = np.where('Sawtooth' in waveTypes)
    10211025            #XmodZ = 0   replace
     1026    if Ax.shape[1]:
     1027        tauX = np.arange(1.,Ax.shape[1]+1-nx)[:,nxs]*glTau  #Xwaves x 32
     1028        StauX = np.ones_like(Ax)[:,nx:,:,nxs]*np.sin(twopi*tauX)[nxs,:,nxs,:]   #also dXmodA/dAx
     1029        CtauX = np.ones_like(Bx)[:,nx:,:,nxs]*np.cos(twopi*tauX)[nxs,:,nxs,:]   #also dXmodB/dBx
     1030        XmodA = Ax[:,nx:,:,nxs]*StauX #atoms X waves X pos X 32
     1031        XmodB = Bx[:,nx:,:,nxs]*CtauX #ditto
     1032        Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1)                #atoms X pos X 32; sum waves
     1033        D = H[:,3][:,nxs]*glTau[nxs,:]              #m*tau; ops X 32
     1034        HdotX = np.inner(H[:,:3],np.swapaxes(Xmod,1,2))+D[:,nxs,:]     #ops x atoms X 32
    10221035    if Af.shape[1]:
    10231036        tauF = np.arange(1.,Af.shape[1]+1-nf)[:,nxs]*glTau  #Fwaves x 32
     
    10281041        Fmod = np.sum(FmodA+FmodB+FmodC,axis=1)             #atoms X 32; sum waves
    10291042    else:
    1030         Fmod = 1.0           
    1031     if Ax.shape[1]:
    1032         tauX = np.arange(1.,Ax.shape[1]+1-nx)[:,nxs]*glTau  #Xwaves x 32
    1033         StauX = np.ones_like(Ax)[:,nx:,:,nxs]*np.sin(twopi*tauX)[nxs,:,nxs,:]   #also dXmodA/dAx
    1034         CtauX = np.ones_like(Bx)[:,nx:,:,nxs]*np.cos(twopi*tauX)[nxs,:,nxs,:]   #also dXmodB/dBx
    1035         XmodA = Ax[:,nx:,:,nxs]*StauX #atoms X waves X pos X 32
    1036         XmodB = Bx[:,nx:,:,nxs]*CtauX #ditto
    1037         Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1)                #atoms X pos X 32; sum waves
     1043        Fmod = np.ones_like(HdotX)           
    10381044    if Au.shape[1]:
    10391045        tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau     #Uwaves x 32
     
    10451051        HbH = np.exp(-np.sum(H[:,nxs,nxs,:3]*np.inner(H[:,:3],Umod),axis=-1)) # ops x atoms x 32 add Overhauser corr.?
    10461052    else:
    1047         HbH = 1.0
    1048     D = H[:,3][:,nxs]*glTau[nxs,:]              #m*tau; ops X 32
    1049     HdotX = np.inner(H[:,:3],np.swapaxes(Xmod,1,2))+D[:,nxs,:]     #ops x atoms X 32
    1050     sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=-1)       #ops x atoms; sum for G-L integration
     1053        HbH = np.ones_like(HdotX)
     1054    HdotXA = H[:,nxs,nxs,nxs,:3]*np.swapaxes(XmodA,-1,-2)[nxs,:,:,:,:]  #ops x atoms x waves x 32 x xyz
     1055    HdotXB = H[:,nxs,nxs,nxs,:3]*np.swapaxes(XmodB,-1,-2)[nxs,:,:,:,:]
     1056    dHdXA = H[:,nxs,nxs,nxs,:3]*np.swapaxes(StauX,-1,-2)[nxs,:,:,:,:]    #ops x atoms x waves x 32 x xyz
     1057    dHdXB = H[:,nxs,nxs,nxs,:3]*np.swapaxes(CtauX,-1,-2)[nxs,:,:,:,:]              #ditto
     1058    sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=-1)       #ops x atoms x waves x xyz; sum for G-L integration
    10511059    cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=-1)       #ditto
    1052 #   GSASIIpath.IPyBreak()                 
     1060    dGdAx = np.sum((Fmod*HbH)[:,:,nxs,:,nxs]*(twopi*dHdXA*np.sin(twopi*HdotXA)+np.cos(twopi*HdotXA))*glWt[nxs,nxs,nxs,:,nxs],axis=-2)
     1061# ops x atoms x waves x xyz
     1062    dGdBx = np.sum((Fmod*HbH)[:,:,nxs,:,nxs]*(twopi*dHdXB*np.cos(twopi*HdotXB)-np.sin(twopi*HdotXB))*glWt[nxs,nxs,nxs,:,nxs],axis=-2)
     1063#    GSASIIpath.IPyBreak()                 
    10531064    return np.array([cosHA,sinHA]),dGdAf,dGdBf,dGdAx,dGdBx,dGdAu,dGdBu      #ops X atoms
    10541065   
Note: See TracChangeset for help on using the changeset viewer.