[SOLVED] CS计算机代考程序代写 python Computational Methods

Computational Methods
for Physicists and Materials Engineers
5
Systems of Nonlinear Equations

xi are unknowns to be determined Aij and bi are coefficients (knowns)
Matrix form:
Recap: system of linear equations
N linear equations: AxAxAxb
11 1 12 2 1N N 1 AxAxAxb
21 1 22 2 2N N 2 
AxAxAxb N1 1 N2 2 NN N N
AAAxb  11 12 1N   1   1 
21 22 2N 2 2 AAAxb
 A A Ax b
 N1 N2 NN   N   N  A xb

System of nonlinear equations
From linear to nonlinear
Ax  b
(I  I  A)x  b
If g is a linear operator
x  (I  A) x  b  f(x)
linear constant A(αx)  αAx
operator vector
However, in general, f cannot be written in this form, e.g. f(x)  (xT x)a
f(x)  sin(x)
If so, the system of equations is nonlinear
System of nonlinear equations can only be solved by iterative method
g(αx)  αg(x) For a matrix A

Root finding problem: Single‐variate single equation
x4 5×3 2×2 x10 Single‐variate nonlinear equation
-4 -2 0
sin(x)  1 x
f 0 -20
System of nonlinear equations: The product of two numbers is 10 and the sum of their squares is 36; what are the numbers?
x
xy  10
x2 y2 36
6 y0 -6
Problems
40 f 0 -40
2
x
10 20
-6 0 6
x

R
If there are two phases, we have g(x) for each phase:
Problems
Calculation of phase diagrams
Consider an A‐B binary system. Molar Gibbs free energy is a function of
the composition of B, x. E.g. from regular solution model, gα(x)(1x)gα xgα (1x)xωα RT (1x)ln(1x)xlnx
gαA, g αB ωα
: molar Gibbs free energies for pure A and B in α structure : the parameter which describes interaction energy in α
: gas constant (8.314 J∙K‒1∙mol‒1)
AB

 gβ(x)(1x)gβ xgβ (1x)xωβ RT (1x)ln(1x)xlnx
gα(x)(1x)gα xgα (1x)xωα RT (1x)ln(1x)xlnx
AB AB

Both gα and gβ vary with temperature T. Phase diagram can be produced by the common tangent construction at each T

T
Calculation of phase diagrams
g
At each temperature T, the slope of gα at x1 is same as the slope of gβ at x2:
x
Problems
gα(x)(1x)gα xgα (1x)xωα RT (1x)ln(1x)xlnx
AB
gβ(x)(1x)gβ xgβ (1x)xωβ RT (1x)ln(1x)xlnx
AB
Both gα and gβ vary with temperature T. Find the tangent line common for two curves. The tangent points x1 and x2 correspond to the phase boundaries at this T
gα  gβ
x xx1 x xx2



T
Calculation of phase diagrams
g
gβ
The two lines should be the same line: f1(x) = f2(x)
x
xx 1
f (x)gα(x ) f (x)gβ(x )
(xx ) (xx )
Problems
gα(x)(1x)gα xgα (1x)xωα RT (1x)ln(1x)xlnx
AB
gβ(x)(1x)gβ xgβ (1x)xωβ RT (1x)ln(1x)xlnx
AB
The tangent lines at x1 of gα and at x2 of gβ are
gα
1 1 x 1
xx1
2 2 x 2 x 2
gα gα x gα(x )
gα
xx2 xx1
x gβ(x ) x 1 1 x 2 2
gα x
2121
xx1
xx1
(x x )gβ(x )gα(x )


(xx )

T
Calculation of phase diagrams
g
x
In reality, the function g(x) is more complicated than the regular solution model. Numerical method is required
AB
gα(x)gα gα (12x)ωα RTln BA
x 1x
gβ(x)gβ gβ (12x)ωβ RTln BA
x 1x
Problems
gα(x)(1x)gα xgα (1x)xωα RT (1x)ln(1x)xlnx
AB
gβ(x)(1x)gβ xgβ (1x)xωβ RT (1x)ln(1x)xlnx
For each T, solve the system of nonlinear equations: gα(x )gβ(x )
12
gα(x )(x x )gβ(x )gα(x ) 12121



System of nonlinear equations
Fixed‐point problem:
Find x such that f(x) = x
Root‐finding problem:
Find x such that f(x) = 0

Basics of iterative method: Banach fixed point theorem
An operator f is called contraction if there exists a constant q  [0, 1) s.t.
f(x)f(y) q xy forallx,y Banach fixed point theorem:
Let f be a contraction. There exists a unique element x* (fixed point) s.t. f(x*)x*
The sequence constructed by
xn1 :f(xn)
converges to the unique fixed point x* with arbitrary initial guess x0
Priori error estimate: Posteriori error estimate:
qn
xn x* 1q x1 x0
xnx* q xnxn1 1q

Banach fixed point theorem for 1 equation w/ 1 variable
An operator f is called contraction if there exists a constant q  [0, 1) s.t. f(x)f(y)qxy forallx,y
Banach fixed point theorem:
Let f be a contraction. There exists a unique solution x* (fixed point) s.t.
f (x*)  x * The sequence constructed by
Find the intersection between 2 curves:
y = x and y = f(x)
xn1 : f(xn)
converges to the unique fixed point x* with arbitrary initial guess x0
Priori error estimate: Posteriori error estimate:
qn
xn x 1q x1 x0
xnx q xnxn1 1q

Requirement for f: There exists a constant q  [0, 1) s.t. f(x)f(y)qxy forallx,y
What does this mean?
Recap of calculus
The mean value theorem: f(x) is continuous on [a, b] and differentiable on (a, b). There exists ξ  (a, b) s.t.
f(ξ) f(b) f(a) ba
f(x) f(b)
f(a)
aξ ξ’b
x

Requirement for f: There exists a constant q  [0, 1) s.t. f(x)f(y)qxy forallx,y
What does this mean?
There exists ξ  (x, y) s.t.
f(ξ) f(x)f(y)
xy
f(x) f(y)  f(ξ)(xy)  f(ξ) xy sup f(ξ) xy
f (x)  f (y)  q x  y So, convergence condition can be restated as:
sup f(x) 1 xD
In whole range slope of curve always < 1ξDqyy=xExample f (x)  xxy = f(x)yy=xFixed pointy = f(x)x*xExample f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solution Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx* x0xInitial guessy = f(x)f(x0)y = f(x)yy=xExample f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx* x0x f(x0)y = f(x)yy=xExample f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx* x1 x0 = f(x0)xf(x1)Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx* x1 x0xy = f(x) f(x1)Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx*x2 x1 x0 = f(x1)xy = f(x)f(x2)Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx*x2 x1 x0xy = f(x) f(x2)Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx*x3 x2 x1 x0 = f(x2)xy = f(x)Example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx*x3 x2 x1 x0xlimxn x* ny = f(x)Homework problemTry this f(x) with |f’(x)| > 1
y
y=x
x0 Try this initial guess
x
Example f (x)  x
y = f(x)

y
y=x
Another example f (x)  x
y = f(x) x

y
y=x
Another example f (x)  x
Since slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionFixed pointx*y = f(x) x yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx* x0 Initialy = f(x) xguess f(x0)y = f(x) xyy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx* x0 f(x0)y = f(x) xyy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x* x0 = f(x0) f(x1)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x* x0y = f(x) xf(x1)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x* x2 x0 = f(x1)y = f(x) x f(x2)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x* x2 x0y = f(x) x f(x2)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x3x* x2 x0 = f(x2)y = f(x) xf(x3)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx1 x3x* x2 x0y = f(x) x f(x3)Another example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionyy=xx x x*x x x 13420y = f(x) x= f(x3)yy=xAnother example f (x)  xSince slope |f’(x)| < 1 for all x, there is a unique solution. We can use successive approximation: xn+1 := f(xn) to find this solutionx x x*x x x 13420y = f(x) xlimxn x* n Homework problem Try this f(x) with |f’(x)| > 1
y
y=x
Another example f (x)  x
x0 Try this initial guess
y = f(x) x

If the convergence condition is not satisfied, how to rearrange the equation s.t. the condition is satisfied?
y
y = f(x) 45o
y=x
f(x)x withsupf(x)1 xD
y  f(x) x g(y) x g(x) yx 
α2 > 45o
x

x=y x
45o
α2 < 45o yx = g(y) If the convergence condition is not satisfied, how to rearrange the equation s.t. the condition is satisfied?g(x)x withsupg(x)1 xDmaxiter = 1000tol = 1e-10x = 10.0for n in range(1, maxiter+1):# max allowed iterations# tolerance# initial guessxprev = xx = f(x)Python : Fixed‐point iterationif abs(x – xprev) < tol:print(f”Solution is x = {x}”)breakx f(x)0.4arctan(x) f(x) 0.4×2 1Always |f’(x)| ≤ 0.4. So, xn+1 := f(xn) converges to the fixed point x*import numpy as npdef f(x):return 0.4 * np.arctan(x)if n == maxiter:print(f”Convergence is not achieved in {n} steps!”) import scipy.optimizeimport numpy as npdef f(x):return 0.4 * np.arctan(x)Python : Fixed‐point iterationx f(x)0.4arctan(x)print(scipy.optimize.fixed_point(f, [10.0], xtol=1e-10, maxiter=1000))# scipy.optimize.fixed_point(function, [initial guess], tolerance, maximum allowed steps)Python : Fixed‐point iterationx  f (x)  10arctan(x) f(x) 10×2 1About x = 0, |f’(x)| > 1. So, xn+1 := f(xn) won’t find x* = 0
f f(x)0 ‐20
20
‐20 ‐10 0 10 20
x
g(f)
f 10arctan(x)xtan(0.1f)g(f) f g(f)
Aboutf=0,|g’(f)|<1.So,fn+1 :=g(fn)canfindf*=0f  x  Banach fixed point theorem for a system of equations w/ N variables x  f(x ,x ,,x )Jacobian matrix of f is1 112 N 2212Nxx , f(x)f(x ,x ,,x ) x  f (x ,x ,,x ) N N12 Nf ff 111x1 x2 xN  f f f  f 222 f(x)xx1 x2 xN  ff fMMMxx x 12N Banach fixed point theorem for a system of equations w/ N variables Let f be continuous and differentiable with the property:There exists a unique solution x* s.t. f(x*)x*The sequence constructed byxn1 :f(xn)sup f(x) 1 xDconverges to the unique fixed point x* with arbitrary initial guess x0Priori error estimate: Posteriori error estimate:qnxn x* 1q x1 x0xnx* q xnxn1 1qExample Solve the system of nonlinear equations:x1 0.5cosx1 0.5sinx2 x2 0.5sinx1 0.5cosx2f f1 1f(x)x1 x2 0.5sinx1 0.5cosx2   0.5cosx 0.5sinx  x1 x2  f f2 2 1 2 (0.5sin x )2  (0.5cos x )2  (0.5cos x )2  (0.5sin x )2 1/2  0.5 21122f(x)So, we can construct the successive approximation:s.t.xn1 :f(xn) x1,n1 0.5cosx1,n 0.5sinx2,nlimxn x* x2,n1 0.5sinx1,n 0.5cosx2,n n import numpy as npdef f(x1, x2):f1 = 0.5*np.cos(x1) – 0.5*np.sin(x2) f2 = 0.5*np.sin(x1) + 0.5*np.cos(x2) return (f1, f2)maxiter = 10000 # max allowed iterations tol = 1e-10 # tolerancex1_arr = np.array([]); x2_arr = np.array([])x1 = 0; x2 = 0x1_arr = np.append(x1_arr, x1)x2_arr = np.append(x2_arr, x2)for n in range(1, maxiter+1):# initial guessExamplex1_prev = x1; x2_prev = x2x1, x2 = f(x1, x2)x1_arr = np.append(x1_arr, x1)x2_arr = np.append(x2_arr, x2)if np.linalg.norm([x1 – x1_prev, x2 – x2_prev], 2) < tol:breakif n == maxiter:print(f”Solution is (x1, x2) = ({x1}, {x2}).”)print(f”Convergence is not achieved in {n} steps!”) import scipy.optimizeimport numpy as npExample Use SciPydef f_arr(x_arr):return np.array([0.5*np.cos(x_arr[0]) – 0.5*np.sin(x_arr[1]),0.5*np.sin(x_arr[0]) + 0.5*np.cos(x_arr[1])])# note that function used by scipy.optimize.fixed_point must return numpy arrayprint(scipy.optimize.fixed_point(f_arr, [0., 0.], xtol=1e-10, maxiter=1000))x2Example Solve the system of nonlinear equations:0.8 0.6 0.4 0.2Different initial guesses converge to the same fixed pointx1 0.5cosx1 0.5sinx2 x2 0.5sinx1 0.5cosx20.00.0 0.2x 0.4 0.6 1 For example,f(x)0 f(x)x2 10Root finding : Newton’s methodSolve the system of nonlinear equations: f(x)0This is a root‐finding problemAgain, we start from 1D. Solve a nonlinear equation:f (x)  sin(x)cos(2x) f(x)cos(x)x3If g is a linear operatorg(αx)  αg(x) Find x s.t. f(x) = 0Newton’s method for 1 variablef(x)xxto 1st order,f(x) f(x0) f(x0)(xx0) f(x)Newton’s method for 1 variableLet x0 be approximation to the solution. By Taylor expansion about x = x0xx0x Newton’s method for 1 variableLet x0 be approximation to the solution. By Taylor expansion about x = x0to 1st order,f(x) f(x0) f(x0)(xx0)f(x )f(x )(x x )0x x  f(x0)001010f (x0)f(x)xx1 x0xto 1st order,f(x) f(x1) f(x1)(xx1) f(x)Newton’s method for 1 variableLet x0 be approximation to the solution. By Taylor expansion about x = x1xx1 x0x Newton’s method for 1 variableLet x0 be approximation to the solution. By Taylor expansion about x = x1to 1st order,f(x) f(x1) f(x1)(xx1)f(x )f(x )(x x )0x x  f(x1)112121f (x1)f(x)xx2 x1 x0x Newton’s method for 1 variablef(x)xx n1 nf(xn) f(xn)limxn x nxx2 x1 x0ximport numpy as np fromsympyimport*f (x ) xn1 xn   nNewton’s method for 1 variableExamplef(x)cos(x)x3 0 f(x)sin(x)3x2x = symbols(‘x’, real=True)expr = cos(x) – x**3dexpr = diff(expr, x)f = lambdify([x], expr, “numpy”)fprime = lambdify([x], dexpr, “numpy”)# convert sympy expression to python function (see Lecture 2) maxiter = 10000tol = 1e-5x = 0.5 # if initial guess is x = 0, won’t converge for n in range(1, maxiter+1):xprev = xx = x – f(x) / fprime(x)if abs(x – xprev) < tol:print(f”Solution is x = {x}”)breakif n == maxiter:print(f”Convergence is not achieved in {n} steps!”)f (xn ) import numpy as npimport scipy.optimizeNewton’s method for 1 variableExamplef(x)cos(x)x3 0 f(x)sin(x)3x2x = symbols(‘x’, real=True)expr = cos(x) – x**3dexpr = diff(expr, x)f = lambdify([x], expr, “numpy”) fprime = lambdify([x], dexpr, “numpy”)Use SciPyprint(scipy.optimize.newton(f, 0.5, fprime=fprime, tol=1e-5, maxiter=10000))# scipy.optimize.newton(function, initial guess, derivative of function, tolerance, maximum allowed steps)# it’s possible that there are multiple solutions. try different initial guess# if fprime is not provided, secant method will be appliedNewton’s method for N variablesLet x0 be approximation to the solution x. By Taylor expansion of themulti‐variate function about x = x0 to 1st order, f(x)  f(x0 )  f(x0 )(x  x0 )Recall: f’(x) denotes Jacobian matrixf(x0)f(x0)(x1 x0)0x1 x0 [f(x0)]1f(x0)xn1 xn [f(xn)]1f(xn)In practice, NOT calculate the inverse of Jacobian. At Step (n), we solvethe system of linear equations:f(xn ) (xn  xn1 )  f(xn ) Gaussian, LU, QRAxbIterative:Direct:Jacobi, Gauss‐Seidel, SOR, GMRES, CG Problems of the original Newton’s methodTwo issues(i) Initial guess does not guarantee convergence→ Damped Newton method(ii) Computation of the Jacobian [f’(xn)]‒1 is costly → Quasi‐Newton methodConvergence of Newton’s method xn1 xn [f(xn)]1f(xn)g(xn)This successive approximation is equivalent to that used to solving thefixed point problem:The convergence condition isx  g(x) sup g(x) 1xDTheorem: Let f be twice continuously differentiable. Assume that x* is a zero of f s.t. f’(x*) is nonsingular. Then, the Newton’s method is locally convergent: there exists a neighborhood B of x* s.t. Newton iterations converge to x* for all initial guess x0  BNewtons’ method is local. What about fixed point method? n1 n n n1 In damped Newton method,Damping factor λn > 0
Damped Newton method
Correction to the guess at previous step:
try  1
In Newton’s method,
Δxnxn1xn[f(xn)] f(xn) x x Δx xtry
try  1
xn1 xn [f(xn)] f(xn)
x xλΔx n1 n n n
Initially, λn = 0 (start from n = 0). When ‖f(xn+1)‖ > ‖f(xn)‖, reject this xn+1 and reduce λn (usually halved) until ‖f(xn+1)‖ < ‖f(xn)‖f(x)Convergence of Newton’s methodx1x2 x0xf(x)Convergence of Newton’s methodx1 x1xBy Newton’s method, when the initial guess x0 is sufficiently close to the zero point, iteration converges; otherwise, not necessarily convergesx2 x0 x0 x1tryxf(x) f(x0)Convergence of Newton’s methodf(x1try)|f(x1try)| > |f(x0)|
In damped Newton method, reject this trial
x0

f(x) f(x0)
x1try half f(x1try)
x1try
x
Convergence of Newton’s method
λ0 := λ0/2 i.e. half x1try – x0 |f(x1try)| < |f(x0)| Accept this trial, let x1 = x1tryx2 x0 xn1 xn Anf(xn) An = [f’(x0)]‒1 which keeps same at all stepsChord methodFor 1 variable:xn1 xn [f(x0)]1f(xn) xn1 xn [f(x0)]1 f(xn)Quasi‐Newton method xn1 xn [f(xn)]1f(xn)Computation of Jacobian f’(xn) at each step is costly Replace [f’(xn)]‒1 by some approximate matrix An:For 1 variable: replace derivative by difference xn1 xn [f(xn)]1 f(xn)Quasi‐Newton methodxn1 xn [f(xn)]1f(xn) Computation of Jacobian f’(xn) at each step is costlyReplace [f’(xn)]‒1 by some approximate matrix An: xn1 xn Anf(xn)Secant method f(x)f(x )  f(x) n n1 Jnxxn n n11 xn1xnJn f(xn) Quasi‐Newton methodBut Jn which satisfies these equations is not unique Jn1 is known from the previous stepBroyden’s methodSecant method for N variables: replace derivative by difference f(xn)(xn xn1)f(xn)f(xn1)JΔx Δf nnnJΔx Δf nnnFind the Jn which is the closest to Jn1 Recall: Frobenius normA A2min JnJn1 Fik  FJnΔxnΔfn i,k1ΔfJ Δx JJn n1 nΔxTn n1 ΔxT Δx n nn1xn1 xn Jn f(xn)n 1/2Damping factor λn > 0 1
1
1 λ Δx  J Δf T 1
Quasi‐Newton method
Good Broyden’s method
1
1 1 Δx J Δf T1
Jn Jn1 n n1 nΔxnJn1 ΔxT J1 Δf
xn1 xn Jn f(xn)
Since xn+1 is updated by Jn , it is convenient to directly update Jn
nn1 n 1
1
Damped good Broyden’s method
1
try 1
xn1 xn Jn f(xn)
Δx xtry x J1f(x ) n n1 n n n
x xλΔx n1 n n n
JnJn1n nn1nΔxnJn1
ΔxT J1 Δf nn1 n

Fixed‐point problem:
Equivalent root‐finding problem:
import scipy.optimize
import numpy as np
def f_fp(x):
# function for fixed-point problem
return 0.4 * np.arctan(x)
def f_rf(x):
# function for root-finding problem
return 0.4 * np.arctan(x) – x
Python : Root finding
x  0.4arctan(x) 0.4arctan(x) x  0
print(scipy.optimize.fixed_point(f_fp, [10], xtol=1e-10, maxiter=1000))
sol = scipy.optimize.root(f_rf, [10], method=’hybr’)
# scipy.optimize.root(function, initial guess, method)
# ‘hybr’: hybrid method (default); ‘broyden1′: good broyden print(sol.x)
# return .x: solution; .success: if exit with success

Fixed‐point problem: Root‐finding problem:
x1 0.5cosx1 0.5sinx2 x2 0.5sinx1 0.5cosx2
def f_arr_fp(x_arr):
# function for fixed-point problem
Python : Root finding
0.5cosx1 0.5sinx2 x1 0 0.5sinx1 0.5cosx2 x2 0
return np.array(
[0.5*np.cos(x_arr[0]) – 0.5*np.sin(x_arr[1]),
0.5*np.sin(x_arr[0]) + 0.5*np.cos(x_arr[1])])
def f_arr_rf(x_arr):
# function for root finding problem
return np.array(
[0.5*np.cos(x_arr[0]) – 0.5*np.sin(x_arr[1]) – x_arr[0],
0.5*np.sin(x_arr[0]) + 0.5*np.cos(x_arr[1]) – x_arr[1]])
print(scipy.optimize.fixed_point(f_arr_fp, [0, 0], xtol=1e-10, maxiter=1000))
sol = scipy.optimize.root(f_arr_rf, [0, 0], method=’hybr’) print(sol.x)

System of nonlinear equations
Iterative method
Fixed‐point problem
f(x)x
Recursion formula according to Banach theorem
Root‐finding problem
f(x)0
(i) Newton’s method
(ii) Damped Newton method: convergence (iii) Quasi‐Newton method: speed
Chord method
Secant method
Broyden’s method
Good Broyden’s method Damped good Broyden’s method