Changeset 975 for trunk/GSASIImath.py

Timestamp:

Jun 29, 2013 11:39:57 AM (10 years ago)

Author:

vondreele

Message:

a few mods to HessianLSQ

File:

• trunk/GSASIImath.py (modified) (6 diffs)

trunk/GSASIImath.py

@@ -42 +42 @@
 atan2d = lambda y,x: 180.*np.arctan2(y,x)/np.pi
 
-def HessianLSQ(func,x0,Hess,args=(),ftol=1.49012e-8,xtol=1.49012e-8, maxcyc=0):
+def HessianLSQ(func,x0,Hess,args=(),ftol=1.49012e-8,xtol=1.49012e-8, maxcyc=0,Print=False):
     
     """
@@ -66 +66 @@
     :param int maxcyc: The maximum number of cycles of refinement to execute, if -1 refine 
         until other limits are met (ftol, xtol)
+    :param bool Print: True for printing results (residuals & times) by cycle
 
     :returns: (x,cov_x,infodict) where
@@ -100 +101 @@
     lamMax = lam
     nfev = 0
+    if Print:
+        print ' Hessian refinement on %d variables:'%(n)
     while icycle < maxcyc:
-        lamMax = max(lamMax,lam)
+        time0 = time.time()
         M = func(x0,*args)
         nfev += 1
@@ -112 +115 @@
         Anorm = np.outer(Adiag,Adiag)
         Yvec /= Adiag
-        Amat /= Anorm        
+        Amat /= Anorm
         while True:
             Lam = np.eye(Amat.shape[0])*lam
@@ -135 +138 @@
                 print 'ouch #3 chisq1 ',chisq1,' stuck > chisq0 ',chisq0
                 break
+        lamMax = max(lamMax,lam)
+        if Print:
+            print ' Cycle: %d, Time: %.2fs, Chi**2: %.3g, Lambda: %.3g'%(icycle,time.time()-time0,chisq1,lam)
         if (chisq0-chisq1)/chisq0 < ftol:
             break
@@ -141 +147 @@
     nfev += 1
     Yvec,Amat = Hess(x0,*args)
+    Lam = np.eye(Amat.shape[0])*lam
+    Amatlam = Amat*(One+Lam)
     try:
-        Bmat = nl.inv(Amat)
+        Bmat = nl.inv(Amatlam)
         return [x0,Bmat,{'num cyc':icycle,'fvec':M,'nfev':nfev,'lamMax':lamMax,'psing':[]}]
     except nl.LinAlgError: