# Changeset 2838

Ignore:
Timestamp:
May 24, 2017 2:32:22 PM (5 years ago)
Message:

add a bit of printing to LS & show no. of penalty fxns

Location:
trunk
Files:
3 edited

Unmodified
Added
Removed
• ## trunk/GSASIImath.py

 r2775 Yvec /= Adiag Amat /= Anorm if Print: print 'initial chi^2 %.5g'%(chisq0) while True: Lam = np.eye(Amat.shape[0])*lam lam *= 10. if Print: print 'matrix modification needed; lambda now %.1e'%(lam) print 'new chi^2 %.5g; matrix modification needed; lambda now %.1e'%(chisq1,lam) else: x0 += Xvec Anorm = np.outer(Adiag,Adiag) Lam = np.eye(Amat.shape[0])*lam Amatlam = Amat/Anorm  #*(One+Lam)              #don't scale Amat to Marquardt array Amatlam = Amat/Anorm #    Amatlam = Amat*(One+Lam)                #scale Amat to Marquardt array? try: Bmat = nl.inv(Amatlam)/Anorm  #*(One+Lam)      #don't rescale Bmat to Marquardt array Bmat = nl.inv(Amatlam)/Anorm #        Bmat = Bmat*(One+Lam)               #rescale Bmat to Marquardt array? return [x0,Bmat,{'num cyc':icycle,'fvec':M,'nfev':nfev,'lamMax':lamMax,'psing':[], 'Converged': ifConverged, 'DelChi2':deltaChi2}] except nl.LinAlgError:
• ## trunk/GSASIIpwd.py

 r2808 def getFCJVoigt3(pos,sig,gam,shl,xdata): 'needs a doc string' Df = pyd.pypsvfcj(len(xdata),xdata-pos,pos,sig,gam,shl) #    Df = pyd.pypsvfcjo(len(xdata),xdata-pos,pos,sig,gam,shl)
• ## trunk/GSASIIstrMath.py

 r2837 negWt = {} pWsum = {} pWnum = {} for phase in Phases: pId = Phases[phase]['pId'] for name,rest in names: pWsum[name] = 0. pWnum[name] = 0 itemRest = phaseRest[name] if itemRest[rest] and itemRest['Use']: pWt.append(wt/esd**2) pWsum[name] += wt*((obs-calc)/esd)**2 pWnum[name] += 1 elif name in ['Torsion','Rama']: coeffDict = itemRest['Coeff'] pWt.append(wt/esd**2) pWsum[name] += wt*(restr/esd)**2 pWnum[name] += 1 elif name == 'ChemComp': for i,[indx,factors,obs,esd] in enumerate(itemRest[rest]): pWt.append(wt/esd**2) pWsum[name] += wt*((obs-calc)/esd)**2 pWnum[name] += 1 elif name == 'Texture': SHkeys = textureData['SH Coeff'][1].keys() pWt.append(wt/esd1**2) pWsum[name] += wt*(-Z1[ind[0]][ind[1]]/esd1)**2 pWnum[name] += 1 if ifesd2: Z2 = 1.-Z pWt.append(wt/esd2**2) pWsum[name] += wt*(Z2/esd2)**2 pWnum[name] += 1 for phase in Phases: cell = General['Cell'][1:7] pWsum[name] = 0.0 pWnum[name] = 0 for hist in Phases[phase]['Histograms']: if not Phases[phase]['Histograms'][hist]['Use']: pWt.append(wt) pWsum[name] += wt*(Y[ind])**2 pWnum[name] += 1 pWsum['PWLref'] = 0. pWnum[name] = 0 for item in varyList: if 'PWLref' in item and parmDict[item] < 0.: pWt.append(negWt[pId]) pWsum['PWLref'] += negWt[pId]*(parmDict[item])**2 pWnum[name] += 1 pVals = np.array(pVals) pWt = np.array(pWt)         #should this be np.sqrt? return pNames,pVals,pWt,pWsum return pNames,pVals,pWt,pWsum,pWnum def penaltyDeriv(pNames,pVal,HistoPhases,calcControls,parmDict,varyList): #reflection processing begins here - big arrays! iBeg = 0 time0 = time.time() while iBeg < nRef: iFin = min(iBeg+blkSize,nRef) refl.T[7] = np.copy(refl.T[9]) refl.T[10] = atan2d(fbs[0],fas[0])  #ignore f' & f" #        GSASIIpath.IPyBreak() #                refl.T[10] = atan2d(np.sum(fbs,axis=0),np.sum(fas,axis=0)) #include f' & f" iBeg += blkSize #    print ' %d sf time %.4f\r'%(nRef,time.time()-time0) #    print 'sf time %.4f, nref %d, blkSize %d'%(time.time()-time0,nRef,blkSize) def StructureFactorDerv2(refDict,G,hfx,pfx,SGData,calcControls,parmDict): if not SGData['SGInv'] and 'S' in calcControls[hfx+'histType'] and phfx+'Flack' in parmDict: Flack = 1.-2.*parmDict[phfx+'Flack'] #    time0 = time.time() time0 = time.time() #reflection processing begins here - big arrays! iBeg = 0 #        GSASIIpath.IPyBreak() iBeg += blkSize #    print ' %d derivative time %.4f\r'%(nRef,time.time()-time0) #    print 'derv time %.4f, nref %d, blkSize %d'%(time.time()-time0,nRef,blkSize) #loop over atoms - each dict entry is list of derivatives for all the reflections for i in range(len(Mdata)): else: StructureFactor2(refDict,G,hfx,pfx,SGData,calcControls,parmDict) #            print 'sf calc time: %.3fs'%(time.time()-time0) badPeak = False for iref,refl in enumerate(refDict['RefList']): for item in ['BabA','BabU']: if phfx+item in varylist: dMdvh[varylist.index(phfx+item)][iref] = w*dFdvDict[pfx+item][iref]*parmDict[phfx+'Scale']*ref[11+im] dMdvh[varylist.index(phfx+item)][iref] = w*dFdvDict[phfx+item][iref]*parmDict[phfx+'Scale']*ref[11+im] elif phfx+item in dependentVars: depDerivDict[phfx+item][iref] = w*dFdvDict[pfx+item][iref]*parmDict[phfx+'Scale']*ref[11+im] depDerivDict[phfx+item][iref] = w*dFdvDict[phfx+item][iref]*parmDict[phfx+'Scale']*ref[11+im] else:   #F refinement for iref,ref in enumerate(refDict['RefList']): GetFobsSq(Histograms,Phases,parmDict,calcControls) pNames,pVals,pWt,pWsum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) pNames,pVals,pWt,pWsum,pWnum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) if np.any(pVals): dpdv = penaltyDeriv(pNames,pVals,HistoPhases,calcControls,parmDict,varylist) continue        #skip non-histogram entries GetFobsSq(Histograms,Phases,parmDict,calcControls) pNames,pVals,pWt,pWsum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) pNames,pVals,pWt,pWsum,pWnum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) if np.any(pVals): dpdv = penaltyDeriv(pNames,pVals,HistoPhases,calcControls,parmDict,varylist) dlg.Destroy() raise G2obj.G2Exception('User abort')         #Abort!! pDict,pVals,pWt,pWsum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) pDict,pVals,pWt,pWsum,pWnum = penaltyFxn(HistoPhases,calcControls,parmDict,varylist) if len(pVals): pSum = np.sum(pWt*pVals**2) for name in pWsum: if pWsum[name]: print '  Penalty function for %8s = %12.5g'%(name,pWsum[name]) print '  Penalty function for %5d %8ss = %12.5g'%(pWnum[name],name,pWsum[name]) print 'Total penalty function: %12.5g on %d terms'%(pSum,len(pVals)) Nobs += len(pVals)
Note: See TracChangeset for help on using the changeset viewer.