Changeset 1923

Timestamp:

Jul 8, 2015 3:21:58 PM (8 years ago)

Author:

vondreele

Message:

further twin work

Location:

trunk

Files:

• GSASIIstrIO.py (modified) (3 diffs)
• GSASIIstrMath.py (modified) (9 diffs)

trunk/GSASIIstrIO.py

@@ -2265 +2265 @@
                 sumTwFr = 0.
                 controlDict[pfx+'TwinLaw'] = []
-                controlDict[pfx+'TwinInv'] = []               
+                controlDict[pfx+'TwinInv'] = []
+                NTL = 0            
                 for it,twin in enumerate(Twins):
                     if 'bool' in str(type(twin[0])):
@@ -2271 +2272 @@
                         controlDict[pfx+'TwinLaw'].append(np.zeros((3,3)))
                     else:
+                        NTL += 1
                         controlDict[pfx+'TwinInv'].append(False)
                         controlDict[pfx+'TwinLaw'].append(twin[0])
@@ -2277 +2279 @@
                         sumTwFr += twin[1]
                     else:
-                        hapDict[pfx+'TwinFr:'+str(it)] = twin[1][0]
+                        hapDict[pfx+'TwinFr:0'] = twin[1][0]
+                        controlDict[pfx+'TwinNMN'] = twin[1][2]
                     if Twins[0][1][1]:
                         hapVary.append(pfx+'TwinFr:'+str(it))
+                controlDict[pfx+'NTL'] = NTL
                 #need to make constraint on TwinFr
                 controlDict[pfx+'TwinLaw'] = np.array(controlDict[pfx+'TwinLaw'])

trunk/GSASIIstrMath.py

@@ -745 +745 @@
     TwDict = refDict.get('TwDict',{})           
     if 'S' in calcControls[hfx+'histType']:
+        NTL = calcControls[phfx+'NTL']
+        NM = calcControls[phfx+'TwinNMN']+1
         TwinLaw = calcControls[phfx+'TwinLaw']
         TwinFr = np.array([parmDict[phfx+'TwinFr:'+str(i)] for i in range(len(TwinLaw))])
@@ -785 +787 @@
                 if iref in TwDict:
                     for i in TwDict[iref]:
-                        H[ir][i] = np.array(TwDict[iref][i])*TwinInv[i]
+                        for n in range(NTL):
+                            H[ir][i+n*NM] = np.inner(TwinLaw[n*NM],np.array(TwDict[iref][i])*TwinInv[i+n*NM])
             TwMask = np.any(H,axis=-1)
         SQ = 1./(2.*refl.T[4])**2               #array(blkSize)
@@ -847 +850 @@
     TwDict = refDict.get('TwDict',{})           
     if 'S' in calcControls[hfx+'histType']:
+        NTL = calcControls[phfx+'NTL']
+        NM = calcControls[phfx+'TwinNMN']+1
         TwinLaw = calcControls[phfx+'TwinLaw']
         TwinFr = np.array([parmDict[phfx+'TwinFr:'+str(i)] for i in range(len(TwinLaw))])
@@ -891 +896 @@
             if iref in TwDict:
                 for i in TwDict[iref]:
-                    H[i] = np.array(TwDict[iref][i])*TwinInv[i]
+                    for n in range(NTL):
+                        H[i+n*NM] = np.inner(TwinLaw[n*NM],np.array(TwDict[iref][i])*TwinInv[i+n*NM])
             TwMask = np.any(H,axis=-1)
         SQ = 1./(2.*refl[4])**2             # or (sin(theta)/lambda)**2
@@ -971 +977 @@
                 dFdui[iref] = [2.*TwMask[it]*SA[it]*(dfadui[0][it]+dfbdui[1][it])+2.*TwMask[it]*SB[it]*(dfbdui[0][it]+dfadui[1][it]) for it in range(nTwin)]
                 dFdua[iref] = [2.*TwMask[it]*SA[it]*(dfadua[it][0]+dfbdua[it][1])+2.*TwMask[it]*SB[it]*(dfbdua[it][0]+dfadua[it][1]) for it in range(nTwin)]
-                dFdtw[iref] = TwMask*np.sum(fas,axis=0)**2+TwMask*np.sum(fbs,axis=0)**2
+                dFdtw[iref] = np.sum(TwMask*fas,axis=0)**2+np.sum(TwMask*fbs,axis=0)**2
+                
             else:   #these are good for no twin single crystals
                 dFdfr[iref] = 2.*SA*(dfadfr[0]+dfbdfr[1])*Mdata/len(Uniq)+ \
@@ -2456 +2463 @@
     G,g = G2lat.A2Gmat(A)       #recip & real metric tensors
     refDict = Histogram['Data']
+    if parmDict[phfx+'Scale'] < 0.:
+        parmDict[phfx+'Scale'] = .001
     if im:
         dFdvDict = SStructureFactorDerv(refDict,im,G,hfx,pfx,SGData,SSGData,calcControls,parmDict)
@@ -2481 +2490 @@
                             depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*ref[11+im]
                     if phfx+'Scale' in varylist:
-                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9+im]*ref[11+im]  #OK
+                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[7+im]*ref[11+im]/parmDict[phfx+'Scale']  #OK
                     elif phfx+'Scale' in dependentVars:
-                        depDerivDict[phfx+'Scale'][iref] = w*ref[9+im]*ref[11+im]   #OK
+                        depDerivDict[phfx+'Scale'][iref] = w*ref[7+im]*ref[11+im]/parmDict[phfx+'Scale']   #OK
                     for item in ['Ep','Es','Eg']:
                         if phfx+item in varylist and phfx+item in dervDict:
@@ -2510 +2519 @@
                             depDerivDict[var][iref] = w*dFdvDict[var][iref]*parmDict[phfx+'Scale']*ref[11+im]
                     if phfx+'Scale' in varylist:
-                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[9+im]*ref[11+im]  #OK
+                        dMdvh[varylist.index(phfx+'Scale')][iref] = w*ref[7+im]*ref[11+im]/parmDict[phfx+'Scale']  #OK
                     elif phfx+'Scale' in dependentVars:
-                        depDerivDict[phfx+'Scale'][iref] = w*ref[9+im]*ref[11+im]   #OK                    
+                        depDerivDict[phfx+'Scale'][iref] = w*ref[7+im]*ref[11+im]/parmDict[phfx+'Scale']   #OK                    
                     for item in ['Ep','Es','Eg']:   #OK!
                         if phfx+item in varylist and phfx+item in dervDict:
@@ -2728 +2737 @@
             SGData = Phase['General']['SGData']
             im = 0
+            if parmDict[phfx+'Scale'] < 0.:
+                parmDict[phfx+'Scale'] = .001                
             if Phase['General']['Type'] in ['modulated','magnetic']:
                 SSGData = Phase['General']['SSGData']