# Changeset 961

Ignore:
Timestamp:
Jun 20, 2013 11:57:22 AM (10 years ago)
Message:

fix rb constraints - useCount
fix Uij2Ueqv
mods to anneal algorithms
allow save of MC/SA results from one set of runs to another
a callafter in UpdatePDFGrid - may need more
fixes for when RBs defined but not used

Location:
trunk
Files:
7 edited

Unmodified
Removed
• ## trunk/GSASIIconstrGUI.py

 r934 del data['Vector'][rbId] data['RBIds']['Vector'].remove(rbId) rbData['useCount'] -= 1 wx.CallAfter(UpdateVectorRB)
• ## trunk/GSASIIlattice.py

 r959 U = np.inner(Amat,np.inner(U,Amat).T) E,R = nl.eigh(U) return (E[0]+E[2]+E[5])/3.      #lower triangle? return np.sum(E)/3. def CosSinAngle(U,V,G):
• ## trunk/GSASIImath.py

 r960 #  A schedule due to Lester Ingber #  A schedule due to Lester Ingber modified to use bounds - OK class fast_sa(base_schedule): def init(self, **options): u = squeeze(random.uniform(0.0, 1.0, size=self.dims)) T = self.T y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)+1.0 xc = y*(self.upper - self.lower)/2.0+self.lower return xc xc = (sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)+1.0)/2.0 xnew = xc*(self.upper - self.lower)+self.lower return xnew #        y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0) #        xc = y*(self.upper - self.lower) return class cauchy_sa(base_schedule): #    def update_guess(self, x0): #        x0 = asarray(x0) class cauchy_sa(base_schedule):     #modified to use bounds - not good def update_guess(self, x0): x0 = asarray(x0) numbers = squeeze(random.uniform(-pi/4, pi/4, size=self.dims)) xc = (1.+(self.learn_rate * self.T * tan(numbers))%1.) xnew = xc*(self.upper - self.lower)+self.lower return xnew #        numbers = squeeze(random.uniform(-pi/2, pi/2, size=self.dims)) #        xc = self.learn_rate * self.T * tan(numbers) return class log_sa(base_schedule): class log_sa(base_schedule):        #OK def init(self,**options): if abs(af[-1]-best_state.cost) > feps*10: retval = 5 print "Warning: Cooled to %f at %s but this is not" \ % (squeeze(last_state.cost), str(squeeze(last_state.x))) \ + " the smallest point found." #                print "Warning: Cooled to %f at %s but this is not" \ #                      % (squeeze(last_state.cost), str(squeeze(last_state.x))) \ #                      + " the smallest point found." break if (Tf is not None) and (schedule.T < Tf): ''' gamFW = lambda s,g: math.exp(math.log(s**5+2.69269*s**4*g+2.42843*s**3*g**2+4.47163*s**2*g**3+0.07842*s*g**4+g**5)/5.) twopi = 2.0*np.pi global tsum tsum = 0. def getMDparms(item,pfx,parmDict,varyList): def GetAtomM(Xdata,SGData): Mdata = [] for xyz in Xdata.T: for xyz in Xdata: Mdata.append(float(len(G2spc.GenAtom(xyz,SGData)))) return np.array(Mdata) if parmDict[pfx+'Type'] in ['Vector','Residue']: if parmDict[pfx+'Type'] == 'Vector': RBId = parmDict[pfx+'RBId'] RBRes = RBdata['Vector'][RBId] RBRes = RBdata['Vector'][parmDict[pfx+'RBId']] aTypes = RBRes['rbTypes'] vecs = RBRes['rbVect'] Cart += vec*mag elif parmDict[pfx+'Type'] == 'Residue': RBId = parmDict[pfx+'RBId'] RBRes = RBdata['Residue'][RBId] RBRes = RBdata['Residue'][parmDict[pfx+'RBId']] aTypes = RBRes['rbTypes'] Cart = np.array(RBRes['rbXYZ']) if parm in parmDict: keys[key][atNo] = parmDict[parm] iatm += 1 else: continue        #skips March Dollase return Tdata,Xdata return Tdata,Xdata.T def calcMDcorr(MDval,MDaxis,Uniq,G): #            ParmDict: #        puts result F^2 in each ref[8] in refList #        ''' twopi = 2.0*np.pi #        ''' global tsum parmDict.update(dict(zip(varyList,values))) Tdata,Xdata = GetAtomTX(RBdata,parmDict) Srefs = np.zeros(len(refList)) sumFcsq = 0. t0 = time.time() for refl in refList: fbs = 0. H = refl[:3] for i,El in enumerate(Tdata): for i,El in enumerate(Tdata):       #NB: generator here is slower! FF[i] = refl[7][El] Uniq = refl[8] phi = refl[9] phase = twopi*(np.inner(Uniq,(Xdata.T))+phi[:,np.newaxis]) sinp = np.sin(phase) cosp = np.cos(phase) occ = Mdata/len(Uniq) fa = np.array(FF*occ*cosp) fas = np.sum(fa) if not ifInv: fb = np.array(FF*occ*sinp) fbs = np.sum(fb) fcsq = (fas**2+fbs**2)*refl[3]      #*calcMDcorr(MDval,MDaxis,Uniq,Gmat) sumFcsq += fcsq refl[5] = fcsq FF *= Mdata/len(refl[8])            #FF*occ phase = np.inner(refl[8],Xdata)     #hx+ky+lz phase += refl[9][:,np.newaxis]      #hx+ky+lz+phi cosp = np.cos(twopi*phase) fas = np.sum(FF*cosp) if ifInv: fbs = 0. else: sinp = np.sin(twopi*phase) fbs = np.sum(FF*sinp) refl[5] = (fas**2+fbs**2)*refl[3]  #*calcMDcorr(MDval,MDaxis,Uniq,Gmat) sumFcsq += refl[5] tsum += (time.time()-t0) scale = (parmDict['sumFosq']/sumFcsq) for iref,refl in enumerate(refList): quench=MCSA['fast parms'][0], m=MCSA['fast parms'][1], n=MCSA['fast parms'][2], lower=lower, upper=upper, slope=MCSA['log slope'],dlg=pgbar) return [False,results[1],results[2],results[0],varyList] Result = [False,False,results[1],results[2],]+list(results[0]) Result.append(varyList) return Result,tsum
• ## trunk/GSASIIphsGUI.py

 r953 resultsGrid.ForceRefresh() result = Results[r] for key,val in zip(result[4],result[3]): for key,val in zip(result[-1],result[4:-1]): vals = key.split(':') nObj,name = int(vals[0]),vals[1] wx.CallAfter(UpdateMCSA) G2plt.PlotStructure(G2frame,data) elif c == 1: if Results[r][1]: Results[r][1] = False else: Results[r][1] = True resultsTable.SetValue(r,c,Results[r][1]) resultsGrid.ForceRefresh() resultsSizer = wx.BoxSizer(wx.VERTICAL) resultVals = [] for result in Results: maxVary = max(maxVary,len(result[3])) resultVals.append(result[:3]+list(result[3])) maxVary = max(maxVary,len(result[-1])) resultVals.append(result[:-1]) rowLabels = [] for i in range(len(Results)): rowLabels.append(str(i)) colLabels = ['Select','Residual','Tmin',] for item in result[4]: colLabels.append(item) #            for i in range(maxVary): colLabels.append('variable:'+str(i)) Types = [wg.GRID_VALUE_BOOL,wg.GRID_VALUE_FLOAT+':10,4', colLabels = ['Select','Keep','Residual','Tmin',] for item in result[-1]: colLabels.append(item)   #from last result from for loop above Types = [wg.GRID_VALUE_BOOL,wg.GRID_VALUE_BOOL,wg.GRID_VALUE_FLOAT+':10,4', wg.GRID_VALUE_FLOAT+':10,4',]+maxVary*[wg.GRID_VALUE_FLOAT+':10,5',] resultsTable = G2gd.Table(resultVals,rowLabels=rowLabels,colLabels=colLabels,types=Types) for r in range(resultsGrid.GetNumberRows()): for c in range(resultsGrid.GetNumberCols()): if c == 0: if c in [0,1]: resultsGrid.SetReadOnly(r,c,isReadOnly=False) else: if G2frame.dataFrame.PhaseUserSize is None: mainSizer.FitInside(G2frame.dataFrame) Size = mainSizer.GetMinSize() Size = mainSizer.Fit() Size[0] += 40 Size[1] = max(Size[1],290) + 35 Size[1] = max(Size[1],350) + 35 MCSA.SetSize(Size) MCSA.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1) phaseName = generalData['Name'] MCSAdata = data['MCSA'] saveResult = [] for result in MCSAdata['Results']: if result[1]:       #keep? saveResult.append(result) MCSAdata['Results'] = saveResult covData = {} if 'PWDR' in reflName: pgbar.SetPosition(wx.Point(screenSize[2]-Size[0]-305,screenSize[1]+5)) pgbar.SetSize(Size) time1 = time.time() try: tsf = 0. for i in range(mcsaControls['Cycles']): MCSAdata['Results'].append(G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar)) print ' MC/SA runs completed: ',i Result,tsum = G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar) MCSAdata['Results'].append(Result) print ' MC/SA runs completed: %d residual: %.3f%%'%(i,100*Result[2]) tsf += tsum print ' MC/SA run time: %.2f'%(time.time()-time1) print ' Structure factor time: %.2f'%(tsf) finally: pgbar.Destroy() MCSAdata['Results'] = G2mth.sortArray(MCSAdata['Results'],2,reverse=False) UpdateMCSA() G2plt.PlotStructure(G2frame,data)