Changeset 1993
- Timestamp:
- Oct 8, 2015 4:09:47 PM (7 years ago)
- Location:
- trunk
- Files:
-
- 4 edited
-
GSASIImath.py (modified) (2 diffs)
-
GSASIIphsGUI.py (modified) (2 diffs)
-
GSASIIplot.py (modified) (4 diffs)
-
GSASIIstrMath.py (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/GSASIImath.py
r1992 r1993 994 994 Mast: array orthogonalization matrix for Uij 995 995 ''' 996 997 nxs = np.newaxis 998 glTau,glWt = pwd.pygauleg(0.,1.,32) #get Gauss-Legendre intervals & weights 999 Ax = np.array(XSSdata[:3]).T #atoms x waves x sin pos mods 1000 dGdAx = np.zeros_like(Ax) 1001 Bx = np.array(XSSdata[3:]).T #...cos pos mods 1002 dGdBx = np.zeros_like(Bx) 1003 Af = np.array(FSSdata[0]).T #sin frac mods x waves x atoms 1004 dGdAf = np.zeros_like(Af) 1005 Bf = np.array(FSSdata[1]).T #cos frac mods... 1006 dGdBf = np.zeros_like(Bf) 1007 Au = Mast*np.array(G2lat.U6toUij(USSdata[:6])).T #atoms x waves x sin Uij mods 1008 dGdAu = np.zeros_like(Au) 1009 Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T #...cos Uij mods 1010 dGdBu = np.zeros_like(Bu) 1011 1012 if 'Fourier' in waveTypes: 1013 nf = 0 1014 nx = 0 1015 XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32)) 1016 FmodZ = np.zeros((Af.shape[0],Af.shape[1],32)) 1017 if 'Crenel' in waveTypes: 1018 nC = np.where('Crenel' in waveTypes) 1019 nf = 1 1020 #FmodZ = 0 replace 1021 else: 1022 nx = 1 1023 if 'Sawtooth' in waveTypes: 1024 nS = np.where('Sawtooth' in waveTypes) 1025 #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 1035 if Af.shape[1]: 1036 tauF = np.arange(1.,Af.shape[1]+1-nf)[:,nxs]*glTau #Fwaves x 32 1037 StauF = np.ones_like(Af)[:,nf:,nxs]*np.sin(twopi*tauF)[nxs,:,:] #also dFmodA/dAf 1038 CtauF = np.ones_like(Bf)[:,nf:,nxs]*np.cos(twopi*tauF)[nxs,:,:] #also dFmodB/dBf 1039 FmodA = Af[:,nf:,nxs]*StauF #atoms X Fwaves X 32 1040 FmodB = Bf[:,nf:,nxs]*CtauF 1041 Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves 1042 else: 1043 Fmod = np.ones_like(HdotX) 1044 if Au.shape[1]: 1045 tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 32 1046 StauU = np.ones_like(Au)[:,:,:,:,nxs]*np.sin(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmodA/dAu 1047 CtauU = np.ones_like(Bu)[:,:,:,:,nxs]*np.cos(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmodB/dBu 1048 UmodA = Au[:,:,:,:,nxs]*StauU #atoms x waves x 3x3 x 32 1049 UmodB = Bu[:,:,:,:,nxs]*CtauU #ditto 1050 Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves 1051 HbH = np.exp(-np.sum(H[:,nxs,nxs,:3]*np.inner(H[:,:3],Umod),axis=-1)) # ops x atoms x 32 add Overhauser corr.? 1052 else: 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 1059 cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=-1) #ditto 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) 996 import scipy.misc as sm 997 998 cosHA,sinHA = Modulation(waveTypes,np.array([H,]),FSSdata,XSSdata,USSdata,Mast) 999 deltx = 0.00001 1000 deltf = 0.0001 1001 deltu = 0.000001 1002 Mf = FSSdata.T.shape 1003 dGdMfC = np.zeros(Mf) 1004 dGdMfS = np.zeros(Mf) 1005 Mx = XSSdata.T.shape 1006 Mx = [H.shape[0],]+list(Mx) #atoms x waves x sin+cos pos mods 1007 dGdMxC = np.zeros(Mx) 1008 dGdMxS = np.zeros(Mx) 1009 for i in range(Mx[1]): #atoms 1010 for j in range(Mx[2]): #waves 1011 for k in range(Mx[3]): #coeff 1012 XSSdata[k,j,i] += deltx 1013 Cap,Sap = np.squeeze(Modulation(waveTypes,np.array([H,]),FSSdata,XSSdata,USSdata,Mast)) 1014 XSSdata[k,j,i] -= 2*deltx 1015 Cam,Sam = np.squeeze(Modulation(waveTypes,np.array([H,]),FSSdata,XSSdata,USSdata,Mast)) 1016 XSSdata[k,j,i] += deltx 1017 dGdMxC[:,:,j,k] += (Cap-Cam)/(2*deltx) 1018 dGdMxS[:,:,j,k] += (Sap-Sam)/(2*deltx) 1019 Mu = USSdata.T.shape 1020 dGdMuC = np.zeros(Mu) 1021 dGdMuS = np.zeros(Mu) 1063 1022 # GSASIIpath.IPyBreak() 1064 return np.array([cosHA,sinHA]),np.concatenate((dGdAf,dGdBf),-1), \ 1065 np.concatenate((dGdAx,dGdBx),-1),np.concatenate((dGdAu,dGdBu),-1) #ops X atoms 1023 1024 1025 # nxs = np.newaxis 1026 # glTau,glWt = pwd.pygauleg(0.,1.,32) #get Gauss-Legendre intervals & weights 1027 # Bx = np.array(XSSdata[3:]).T #...cos pos mods 1028 # dGdBx = np.zeros_like(Bx) 1029 # Af = np.array(FSSdata[0]).T #sin frac mods x waves x atoms 1030 # dGdAf = np.zeros_like(Af) 1031 # Bf = np.array(FSSdata[1]).T #cos frac mods... 1032 # dGdBf = np.zeros_like(Bf) 1033 # Au = Mast*np.array(G2lat.U6toUij(USSdata[:6])).T #atoms x waves x sin Uij mods 1034 # dGdAu = np.zeros_like(Au) 1035 # Bu = Mast*np.array(G2lat.U6toUij(USSdata[6:])).T #...cos Uij mods 1036 # dGdBu = np.zeros_like(Bu) 1037 # 1038 # GSASIIpath.IPyBreak() 1039 # if 'Fourier' in waveTypes: 1040 # nf = 0 1041 # nx = 0 1042 # XmodZ = np.zeros((Ax.shape[0],Ax.shape[1],3,32)) 1043 # FmodZ = np.zeros((Af.shape[0],Af.shape[1],32)) 1044 # if 'Crenel' in waveTypes: 1045 # nC = np.where('Crenel' in waveTypes) 1046 # nf = 1 1047 # #FmodZ = 0 replace 1048 # else: 1049 # nx = 1 1050 # if 'Sawtooth' in waveTypes: 1051 # nS = np.where('Sawtooth' in waveTypes) 1052 # #XmodZ = 0 replace 1053 # if Ax.shape[1]: 1054 # tauX = np.arange(1.,Ax.shape[1]+1-nx)[:,nxs]*glTau #Xwaves x 32 1055 # StauX = np.ones_like(Ax)[:,nx:,:,nxs]*np.sin(twopi*tauX)[nxs,:,nxs,:] #also dXmod/dAx 1056 # CtauX = np.ones_like(Bx)[:,nx:,:,nxs]*np.cos(twopi*tauX)[nxs,:,nxs,:] #also dXmod/dBx 1057 # XmodA = Ax[:,nx:,:,nxs]*StauX #atoms X waves X pos X 32 1058 # XmodB = Bx[:,nx:,:,nxs]*CtauX #ditto 1059 # Xmod = np.sum(XmodA+XmodB+XmodZ,axis=1) #atoms X pos X 32; sum waves 1060 # D = H[:,3][:,nxs]*glTau[nxs,:] #m*tau; ops X 32 1061 # HdotX = np.inner(H[:,:3],np.swapaxes(Xmod,1,2))+D[:,nxs,:] #ops x atoms X 32 1062 # if Af.shape[1]: 1063 # tauF = np.arange(1.,Af.shape[1]+1-nf)[:,nxs]*glTau #Fwaves x 32 1064 # StauF = np.ones_like(Af)[:,nf:,nxs]*np.sin(twopi*tauF)[nxs,:,:] #also dFmod/dAf 1065 # CtauF = np.ones_like(Bf)[:,nf:,nxs]*np.cos(twopi*tauF)[nxs,:,:] #also dFmod/dBf 1066 # FmodA = Af[:,nf:,nxs]*StauF #atoms X Fwaves X 32 1067 # FmodB = Bf[:,nf:,nxs]*CtauF 1068 # Fmod = np.sum(FmodA+FmodB+FmodC,axis=1) #atoms X 32; sum waves 1069 # else: 1070 # Fmod = np.ones_like(HdotX) 1071 # if Au.shape[1]: 1072 # tauU = np.arange(1.,Au.shape[1]+1)[:,nxs]*glTau #Uwaves x 32 1073 # StauU = np.ones_like(Au)[:,:,:,:,nxs]*np.sin(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmod/dAu 1074 # CtauU = np.ones_like(Bu)[:,:,:,:,nxs]*np.cos(twopi*tauU)[nxs,:,nxs,nxs,:] #also dUmod/dBu 1075 # UmodA = Au[:,:,:,:,nxs]*StauU #atoms x waves x 3x3 x 32 1076 # UmodB = Bu[:,:,:,:,nxs]*CtauU #ditto 1077 # Umod = np.swapaxes(np.sum(UmodA+UmodB,axis=1),1,3) #atoms x 3x3 x 32; sum waves 1078 # HbH = np.exp(-np.sum(H[:,nxs,nxs,:3]*np.inner(H[:,:3],Umod),axis=-1)) # ops x atoms x 32 add Overhauser corr.? 1079 # else: 1080 # HbH = np.ones_like(HdotX) 1081 # HdotXA = H[:,nxs,nxs,nxs,:3]*np.swapaxes(XmodA,-1,-2)[nxs,:,:,:,:] #ops x atoms x waves x 32 x xyz 1082 # HdotXB = H[:,nxs,nxs,nxs,:3]*np.swapaxes(XmodB,-1,-2)[nxs,:,:,:,:] 1083 # dHdXA = H[:,nxs,nxs,nxs,:3]*np.swapaxes(StauX,-1,-2)[nxs,:,:,:,:] #ops x atoms x waves x 32 x xyz 1084 # dHdXB = H[:,nxs,nxs,nxs,:3]*np.swapaxes(CtauX,-1,-2)[nxs,:,:,:,:] #ditto 1085 # sinHA = np.sum(Fmod*HbH*np.sin(twopi*HdotX)*glWt,axis=-1) #ops x atoms x waves x xyz; sum for G-L integration 1086 # cosHA = np.sum(Fmod*HbH*np.cos(twopi*HdotX)*glWt,axis=-1) #ditto 1087 # dGdAx = np.sum((Fmod*HbH)[:,:,nxs,:,nxs]*(twopi*dHdXA*np.sin(twopi*HdotXA)+np.cos(twopi*HdotXA))*glWt[nxs,nxs,nxs,:,nxs],axis=-2) 1088 ## ops x atoms x waves x xyz - imaginary part 1089 # dGdBx = np.sum((Fmod*HbH)[:,:,nxs,:,nxs]*(twopi*dHdXB*np.cos(twopi*HdotXB)-np.sin(twopi*HdotXB))*glWt[nxs,nxs,nxs,:,nxs],axis=-2) 1090 ## ops x atoms x waves x xyz - real part 1091 return np.array([cosHA,sinHA]),[dGdMfC,dGdMfS],[dGdMxC,dGdMxS],[dGdMuC,dGdMuS] 1066 1092 1067 1093 def posFourier(tau,psin,pcos,smul): … … 2296 2322 h,k,l,m = -hkl+Hmax 2297 2323 Fhkl[h,k,l,m] = dF*phasem 2298 SSrho = fft.fftn(fft.fftshift(Fhkl))/ cell[6] #4D map2324 SSrho = fft.fftn(fft.fftshift(Fhkl))/(10.*cell[6]) #4D map 2299 2325 rho = fft.fftn(fft.fftshift(Fhkl[:,:,:,maxM+1]))/cell[6] #3D map 2300 2326 map4DData['rho'] = np.real(SSrho) -
trunk/GSASIIphsGUI.py
r1988 r1993 1378 1378 atomData[r][c] = atomData[r][c].replace(rbAtmDict.get(atomData[r][ci+8],''),'') 1379 1379 if 'Atoms' in data['Drawing']: 1380 ci = colLabels.index('I/A') 1380 1381 DrawAtomsReplaceByID(data['Drawing'],ci+8,atomData[r],ID) 1381 1382 wx.CallAfter(Paint) … … 2483 2484 Map['MapType'] = mapType.GetValue() 2484 2485 2485 Map['Resolution'] = 0. 252486 Map['Resolution'] = 0.5 2486 2487 refsList = data['Histograms'].keys() 2487 2488 mapSizer = wx.BoxSizer(wx.HORIZONTAL) -
trunk/GSASIIplot.py
r1991 r1993 4475 4475 G2frame.tau = 0. 4476 4476 elif key in ['+','=']: 4477 G2frame.tau += 0. 054477 G2frame.tau += 0.1 4478 4478 elif key == '-': 4479 G2frame.tau -= 0. 054479 G2frame.tau -= 0.1 4480 4480 G2frame.tau %= 1. #force 0-1 range 4481 4481 G2frame.G2plotNB.status.SetStatusText('Modulation tau = %.2f'%(G2frame.tau),1) … … 5056 5056 # end of useful debug 5057 5057 mapData = generalData['Map'] 5058 D4mapData = generalData.get('4DmapData',{}) 5058 5059 pageName = '' 5059 5060 page = getSelection() … … 5061 5062 pageName = G2frame.dataDisplay.GetPageText(page) 5062 5063 rhoXYZ = [] 5063 if len(mapData['rho']): 5064 if len(D4mapData.get('rho',[])): #preferentially select 4D map if there 5065 rho = D4mapData['rho'][:,:,:,int(G2frame.tau*10)] #pick current tau 3D slice 5066 elif len(mapData['rho']): #ordinary 3D map 5067 rho = mapData['rho'] 5068 if len(rho): 5064 5069 VP = drawingData['viewPoint'][0]-np.array([.5,.5,.5]) 5065 5070 contLevel = drawingData['contourLevel']*mapData['rhoMax'] 5066 5071 if 'delt-F' in mapData['MapType'] or 'N' in mapData.get('Type',''): 5067 rho = ma.array( mapData['rho'],mask=(np.abs(mapData['rho'])<contLevel))5072 rho = ma.array(rho,mask=(np.abs(rho)<contLevel)) 5068 5073 else: 5069 rho = ma.array( mapData['rho'],mask=(mapData['rho']<contLevel))5074 rho = ma.array(rho,mask=(rho<contLevel)) 5070 5075 steps = 1./np.array(rho.shape) 5071 5076 incre = np.where(VP>=0,VP%steps,VP%steps-steps) … … 5077 5082 rcube = 2000.*Vol/(ForthirdPI*Nc) 5078 5083 rmax = math.exp(math.log(rcube)/3.)**2 5079 radius = min(drawingData ['mapSize']**2,rmax)5084 radius = min(drawingData.get('mapSize',10.)**2,rmax) 5080 5085 view = drawingData['viewPoint'][0] 5081 5086 Rok = np.sum(np.inner(Amat,rhoXYZ-view).T**2,axis=1)>radius -
trunk/GSASIIstrMath.py
r1992 r1993 1134 1134 dFdfl = np.zeros((nRef,nTwin)) 1135 1135 dFdtw = np.zeros((nRef,nTwin)) 1136 dFdGf = np.zeros((nRef,nTwin,mSize,FSSdata.shape[1],2)) 1137 dFdGx = np.zeros((nRef,nTwin,mSize,XSSdata.shape[1],6)) 1138 dFdGu = np.zeros((nRef,nTwin,mSize,USSdata.shape[1],12)) 1136 dFdGfA = np.zeros((nRef,nTwin,mSize,FSSdata.shape[1])) 1137 dFdGfB = np.zeros((nRef,nTwin,mSize,FSSdata.shape[1])) 1138 dFdGxA = np.zeros((nRef,nTwin,mSize,XSSdata.shape[1],3)) 1139 dFdGxB = np.zeros((nRef,nTwin,mSize,XSSdata.shape[1],3)) 1140 dFdGuA = np.zeros((nRef,nTwin,mSize,USSdata.shape[1],6)) 1141 dFdGuB = np.zeros((nRef,nTwin,mSize,USSdata.shape[1],6)) 1139 1142 else: 1140 1143 dFdfr = np.zeros((nRef,mSize)) … … 1190 1193 fotp = FPP*Tcorr #ops x atoms 1191 1194 GfpuA,dGdf,dGdx,dGdu = G2mth.ModulationDerv(waveTypes,Uniq,FSSdata,XSSdata,USSdata,Mast) 1195 # derivs are: ops x atoms x waves x 1,3,or 6 parms as [real,imag] parts 1192 1196 fa = np.array([((FF+FP).T-Bab).T*cosp*Tcorr,-Flack*FPP*sinp*Tcorr]) # array(2,nTwin,nEqv,nAtoms) 1193 1197 fb = np.array([((FF+FP).T-Bab).T*sinp*Tcorr,Flack*FPP*cosp*Tcorr]) … … 1208 1212 dfadui = np.sum(-SQfactor*fag,axis=1) 1209 1213 dfbdui = np.sum(-SQfactor*fbg,axis=1) 1214 # GSASIIpath.IPyBreak() 1210 1215 if nTwin > 1: 1211 dfadx = np.array([np.sum(twopi*Uniq[it,:,:3]*np.swapaxes(fax,-2,-1)[:,it,:,:,np.newaxis],axis=-2) for it in range(nTwin)]) 1212 dfbdx = np.array([np.sum(twopi*Uniq[it,:,:3]*np.swapaxes(fbx,-2,-1)[:,it,:,:,np.newaxis],axis=-2) for it in range(nTwin)]) 1213 dfadua = np.array([np.sum(-Hij[it]*np.swapaxes(fag,-2,-1)[:,it,:,:,np.newaxis],axis=-2) for it in range(nTwin)]) 1214 dfbdua = np.array([np.sum(-Hij[it]*np.swapaxes(fbg,-2,-1)[:,it,:,:,np.newaxis],axis=-2) for it in range(nTwin)]) 1215 # array(nTwin,2,nAtom,3) & array(nTwin,2,nAtom,6) 1216 dfadx = np.array([np.sum(twopi*Uniq[it,:,:3]*np.swapaxes(fax,-2,-1)[:,it,:,:,nxs],axis=-2) for it in range(nTwin)]) 1217 dfbdx = np.array([np.sum(twopi*Uniq[it,:,:3]*np.swapaxes(fbx,-2,-1)[:,it,:,:,nxs],axis=-2) for it in range(nTwin)]) 1218 dfadua = np.array([np.sum(-Hij[it]*np.swapaxes(fag,-2,-1)[:,it,:,:,nxs],axis=-2) for it in range(nTwin)]) 1219 dfbdua = np.array([np.sum(-Hij[it]*np.swapaxes(fbg,-2,-1)[:,it,:,:,nxs],axis=-2) for it in range(nTwin)]) 1220 # array(nTwin,2,nAtom,3) & array(2,nTwin,nAtom,6) 1221 dfadGx = np.sum(fa[:,it,:,:,nxs,nxs]*dGdx[0][nxs,nxs,:,:,:,:]-fb[:,it,:,:,nxs,nxs]*dGdx[1][nxs,nxs,:,:,:,:],axis=1) 1222 dfbdGx = np.sum(fb[:,it,:,:,nxs,nxs]*dGdx[0][nxs,nxs,:,:,:,:]+fa[:,it,:,:,nxs,nxs]*dGdx[1][nxs,nxs,:,:,:,:],axis=1) 1223 # array (2,nAtom,wave,3) 1216 1224 else: 1217 dfadx = np.sum(twopi*Uniq[:,:3]*np.swapaxes(fax,-2,-1)[:,:,:,n p.newaxis],axis=-2)1218 dfbdx = np.sum(twopi*Uniq[:,:3]*np.swapaxes(fbx,-2,-1)[:,:,:,n p.newaxis],axis=-2)1219 dfadua = np.sum(-Hij*np.swapaxes(fag,-2,-1)[:,:,:,n p.newaxis],axis=-2)1220 dfbdua = np.sum(-Hij*np.swapaxes(fbg,-2,-1)[:,:,:,n p.newaxis],axis=-2)1225 dfadx = np.sum(twopi*Uniq[:,:3]*np.swapaxes(fax,-2,-1)[:,:,:,nxs],axis=-2) 1226 dfbdx = np.sum(twopi*Uniq[:,:3]*np.swapaxes(fbx,-2,-1)[:,:,:,nxs],axis=-2) 1227 dfadua = np.sum(-Hij*np.swapaxes(fag,-2,-1)[:,:,:,nxs],axis=-2) 1228 dfbdua = np.sum(-Hij*np.swapaxes(fbg,-2,-1)[:,:,:,nxs],axis=-2) 1221 1229 # array(2,nAtom,3) & array(2,nAtom,6) 1230 dfadGx = np.sum(fa[:,:,:,nxs,nxs]*dGdx[0][nxs,:,:,:,:]-fb[:,:,:,nxs,nxs]*dGdx[1][nxs,:,:,:,:],axis=1) 1231 dfbdGx = np.sum(fb[:,:,:,nxs,nxs]*dGdx[0][nxs,:,:,:,:]+fa[:,:,:,nxs,nxs]*dGdx[1][nxs,:,:,:,:],axis=1) 1232 # array (2,nAtom,wave,3) 1222 1233 #NB: the above have been checked against PA(1:10,1:2) in strfctr.for for al2O3! 1223 1234 # GSASIIpath.IPyBreak() … … 1255 1266 dFdui[iref] = 2.*SA*(dfadui[0]+dfbdui[1])+2.*SB*(dfbdui[0]+dfadui[1]) #array(nRef,nAtom) 1256 1267 dFdua[iref] = 2.*SA*(dfadua[0]+dfbdua[1])+2.*SB*(dfbdua[0]+dfadua[1]) #array(nRef,nAtom,6) 1257 dFdfl[iref] = -SA*dfadfl-SB*dfbdfl #array(nRef,) 1268 dFdfl[iref] = -SA*dfadfl-SB*dfbdfl #array(nRef,) 1269 dFdGx[iref] = 2.*SA*(dfadGx[0]+dfbdGx[1])+2.*SB*(dfadGx[1]+dfbdGx[0]) #array(nRef,natom,nwave,6) 1258 1270 dFdbab[iref] = 2.*fas[0]*np.array([np.sum(dfadba*dBabdA),np.sum(-dfadba*parmDict[phfx+'BabA']*SQfactor*dBabdA)]).T+ \ 1259 1271 2.*fbs[0]*np.array([np.sum(dfbdba*dBabdA),np.sum(-dfbdba*parmDict[phfx+'BabA']*SQfactor*dBabdA)]).T 1260 1272 #loop over atoms - each dict entry is list of derivatives for all the reflections 1261 for i in range(len(Mdata)): 1273 for i in range(len(Mdata)): #loop over atoms 1262 1274 dFdvDict[pfx+'Afrac:'+str(i)] = dFdfr.T[i] 1263 1275 dFdvDict[pfx+'dAx:'+str(i)] = dFdx.T[0][i] … … 1271 1283 dFdvDict[pfx+'AU13:'+str(i)] = .5*dFdua.T[4][i] 1272 1284 dFdvDict[pfx+'AU23:'+str(i)] = .5*dFdua.T[5][i] 1273 #need dFdvDict[pfx+'Xsin:'+str[i]:str(m)], etc for modulations... 1285 for j in range(XSSdata.shape[1]): 1286 dFdvDict[pfx+'Xsin:'+str(i)+':'+str(j)] = dFdGx.T[0][j][i] 1287 dFdvDict[pfx+'Ysin:'+str(i)+':'+str(j)] = dFdGx.T[1][j][i] 1288 dFdvDict[pfx+'Zsin:'+str(i)+':'+str(j)] = dFdGx.T[2][j][i] 1289 dFdvDict[pfx+'Xcos:'+str(i)+':'+str(j)] = dFdGx.T[3][j][i] 1290 dFdvDict[pfx+'Ycos:'+str(i)+':'+str(j)] = dFdGx.T[4][j][i] 1291 dFdvDict[pfx+'Zcos:'+str(i)+':'+str(j)] = dFdGx.T[5][j][i] 1274 1292 dFdvDict[phfx+'BabA'] = dFdbab.T[0] 1275 1293 dFdvDict[phfx+'BabU'] = dFdbab.T[1]
Note: See TracChangeset
for help on using the changeset viewer.
