Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1""" 

2Unified interfaces to root finding algorithms. 

3 

4Functions 

5--------- 

6- root : find a root of a vector function. 

7""" 

8__all__ = ['root'] 

9 

10import numpy as np 

11 

12 

13from warnings import warn 

14 

15from .optimize import MemoizeJac, OptimizeResult, _check_unknown_options 

16from .minpack import _root_hybr, leastsq 

17from ._spectral import _root_df_sane 

18from . import nonlin 

19 

20 

21def root(fun, x0, args=(), method='hybr', jac=None, tol=None, callback=None, 

22 options=None): 

23 """ 

24 Find a root of a vector function. 

25 

26 Parameters 

27 ---------- 

28 fun : callable 

29 A vector function to find a root of. 

30 x0 : ndarray 

31 Initial guess. 

32 args : tuple, optional 

33 Extra arguments passed to the objective function and its Jacobian. 

34 method : str, optional 

35 Type of solver. Should be one of 

36 

37 - 'hybr' :ref:`(see here) <optimize.root-hybr>` 

38 - 'lm' :ref:`(see here) <optimize.root-lm>` 

39 - 'broyden1' :ref:`(see here) <optimize.root-broyden1>` 

40 - 'broyden2' :ref:`(see here) <optimize.root-broyden2>` 

41 - 'anderson' :ref:`(see here) <optimize.root-anderson>` 

42 - 'linearmixing' :ref:`(see here) <optimize.root-linearmixing>` 

43 - 'diagbroyden' :ref:`(see here) <optimize.root-diagbroyden>` 

44 - 'excitingmixing' :ref:`(see here) <optimize.root-excitingmixing>` 

45 - 'krylov' :ref:`(see here) <optimize.root-krylov>` 

46 - 'df-sane' :ref:`(see here) <optimize.root-dfsane>` 

47 

48 jac : bool or callable, optional 

49 If `jac` is a Boolean and is True, `fun` is assumed to return the 

50 value of Jacobian along with the objective function. If False, the 

51 Jacobian will be estimated numerically. 

52 `jac` can also be a callable returning the Jacobian of `fun`. In 

53 this case, it must accept the same arguments as `fun`. 

54 tol : float, optional 

55 Tolerance for termination. For detailed control, use solver-specific 

56 options. 

57 callback : function, optional 

58 Optional callback function. It is called on every iteration as 

59 ``callback(x, f)`` where `x` is the current solution and `f` 

60 the corresponding residual. For all methods but 'hybr' and 'lm'. 

61 options : dict, optional 

62 A dictionary of solver options. E.g., `xtol` or `maxiter`, see 

63 :obj:`show_options()` for details. 

64 

65 Returns 

66 ------- 

67 sol : OptimizeResult 

68 The solution represented as a ``OptimizeResult`` object. 

69 Important attributes are: ``x`` the solution array, ``success`` a 

70 Boolean flag indicating if the algorithm exited successfully and 

71 ``message`` which describes the cause of the termination. See 

72 `OptimizeResult` for a description of other attributes. 

73 

74 See also 

75 -------- 

76 show_options : Additional options accepted by the solvers 

77 

78 Notes 

79 ----- 

80 This section describes the available solvers that can be selected by the 

81 'method' parameter. The default method is *hybr*. 

82 

83 Method *hybr* uses a modification of the Powell hybrid method as 

84 implemented in MINPACK [1]_. 

85 

86 Method *lm* solves the system of nonlinear equations in a least squares 

87 sense using a modification of the Levenberg-Marquardt algorithm as 

88 implemented in MINPACK [1]_. 

89 

90 Method *df-sane* is a derivative-free spectral method. [3]_ 

91 

92 Methods *broyden1*, *broyden2*, *anderson*, *linearmixing*, 

93 *diagbroyden*, *excitingmixing*, *krylov* are inexact Newton methods, 

94 with backtracking or full line searches [2]_. Each method corresponds 

95 to a particular Jacobian approximations. See `nonlin` for details. 

96 

97 - Method *broyden1* uses Broyden's first Jacobian approximation, it is 

98 known as Broyden's good method. 

99 - Method *broyden2* uses Broyden's second Jacobian approximation, it 

100 is known as Broyden's bad method. 

101 - Method *anderson* uses (extended) Anderson mixing. 

102 - Method *Krylov* uses Krylov approximation for inverse Jacobian. It 

103 is suitable for large-scale problem. 

104 - Method *diagbroyden* uses diagonal Broyden Jacobian approximation. 

105 - Method *linearmixing* uses a scalar Jacobian approximation. 

106 - Method *excitingmixing* uses a tuned diagonal Jacobian 

107 approximation. 

108 

109 .. warning:: 

110 

111 The algorithms implemented for methods *diagbroyden*, 

112 *linearmixing* and *excitingmixing* may be useful for specific 

113 problems, but whether they will work may depend strongly on the 

114 problem. 

115 

116 .. versionadded:: 0.11.0 

117 

118 References 

119 ---------- 

120 .. [1] More, Jorge J., Burton S. Garbow, and Kenneth E. Hillstrom. 

121 1980. User Guide for MINPACK-1. 

122 .. [2] C. T. Kelley. 1995. Iterative Methods for Linear and Nonlinear 

123 Equations. Society for Industrial and Applied Mathematics. 

124 <https://archive.siam.org/books/kelley/fr16/> 

125 .. [3] W. La Cruz, J.M. Martinez, M. Raydan. Math. Comp. 75, 1429 (2006). 

126 

127 Examples 

128 -------- 

129 The following functions define a system of nonlinear equations and its 

130 jacobian. 

131 

132 >>> def fun(x): 

133 ... return [x[0] + 0.5 * (x[0] - x[1])**3 - 1.0, 

134 ... 0.5 * (x[1] - x[0])**3 + x[1]] 

135 

136 >>> def jac(x): 

137 ... return np.array([[1 + 1.5 * (x[0] - x[1])**2, 

138 ... -1.5 * (x[0] - x[1])**2], 

139 ... [-1.5 * (x[1] - x[0])**2, 

140 ... 1 + 1.5 * (x[1] - x[0])**2]]) 

141 

142 A solution can be obtained as follows. 

143 

144 >>> from scipy import optimize 

145 >>> sol = optimize.root(fun, [0, 0], jac=jac, method='hybr') 

146 >>> sol.x 

147 array([ 0.8411639, 0.1588361]) 

148 

149 """ 

150 if not isinstance(args, tuple): 

151 args = (args,) 

152 

153 meth = method.lower() 

154 if options is None: 

155 options = {} 

156 

157 if callback is not None and meth in ('hybr', 'lm'): 

158 warn('Method %s does not accept callback.' % method, 

159 RuntimeWarning) 

160 

161 # fun also returns the Jacobian 

162 if not callable(jac) and meth in ('hybr', 'lm'): 

163 if bool(jac): 

164 fun = MemoizeJac(fun) 

165 jac = fun.derivative 

166 else: 

167 jac = None 

168 

169 # set default tolerances 

170 if tol is not None: 

171 options = dict(options) 

172 if meth in ('hybr', 'lm'): 

173 options.setdefault('xtol', tol) 

174 elif meth in ('df-sane',): 

175 options.setdefault('ftol', tol) 

176 elif meth in ('broyden1', 'broyden2', 'anderson', 'linearmixing', 

177 'diagbroyden', 'excitingmixing', 'krylov'): 

178 options.setdefault('xtol', tol) 

179 options.setdefault('xatol', np.inf) 

180 options.setdefault('ftol', np.inf) 

181 options.setdefault('fatol', np.inf) 

182 

183 if meth == 'hybr': 

184 sol = _root_hybr(fun, x0, args=args, jac=jac, **options) 

185 elif meth == 'lm': 

186 sol = _root_leastsq(fun, x0, args=args, jac=jac, **options) 

187 elif meth == 'df-sane': 

188 _warn_jac_unused(jac, method) 

189 sol = _root_df_sane(fun, x0, args=args, callback=callback, 

190 **options) 

191 elif meth in ('broyden1', 'broyden2', 'anderson', 'linearmixing', 

192 'diagbroyden', 'excitingmixing', 'krylov'): 

193 _warn_jac_unused(jac, method) 

194 sol = _root_nonlin_solve(fun, x0, args=args, jac=jac, 

195 _method=meth, _callback=callback, 

196 **options) 

197 else: 

198 raise ValueError('Unknown solver %s' % method) 

199 

200 return sol 

201 

202 

203def _warn_jac_unused(jac, method): 

204 if jac is not None: 

205 warn('Method %s does not use the jacobian (jac).' % (method,), 

206 RuntimeWarning) 

207 

208 

209def _root_leastsq(fun, x0, args=(), jac=None, 

210 col_deriv=0, xtol=1.49012e-08, ftol=1.49012e-08, 

211 gtol=0.0, maxiter=0, eps=0.0, factor=100, diag=None, 

212 **unknown_options): 

213 """ 

214 Solve for least squares with Levenberg-Marquardt 

215 

216 Options 

217 ------- 

218 col_deriv : bool 

219 non-zero to specify that the Jacobian function computes derivatives 

220 down the columns (faster, because there is no transpose operation). 

221 ftol : float 

222 Relative error desired in the sum of squares. 

223 xtol : float 

224 Relative error desired in the approximate solution. 

225 gtol : float 

226 Orthogonality desired between the function vector and the columns 

227 of the Jacobian. 

228 maxiter : int 

229 The maximum number of calls to the function. If zero, then 

230 100*(N+1) is the maximum where N is the number of elements in x0. 

231 epsfcn : float 

232 A suitable step length for the forward-difference approximation of 

233 the Jacobian (for Dfun=None). If epsfcn is less than the machine 

234 precision, it is assumed that the relative errors in the functions 

235 are of the order of the machine precision. 

236 factor : float 

237 A parameter determining the initial step bound 

238 (``factor * || diag * x||``). Should be in interval ``(0.1, 100)``. 

239 diag : sequence 

240 N positive entries that serve as a scale factors for the variables. 

241 """ 

242 

243 _check_unknown_options(unknown_options) 

244 x, cov_x, info, msg, ier = leastsq(fun, x0, args=args, Dfun=jac, 

245 full_output=True, 

246 col_deriv=col_deriv, xtol=xtol, 

247 ftol=ftol, gtol=gtol, 

248 maxfev=maxiter, epsfcn=eps, 

249 factor=factor, diag=diag) 

250 sol = OptimizeResult(x=x, message=msg, status=ier, 

251 success=ier in (1, 2, 3, 4), cov_x=cov_x, 

252 fun=info.pop('fvec')) 

253 sol.update(info) 

254 return sol 

255 

256 

257def _root_nonlin_solve(fun, x0, args=(), jac=None, 

258 _callback=None, _method=None, 

259 nit=None, disp=False, maxiter=None, 

260 ftol=None, fatol=None, xtol=None, xatol=None, 

261 tol_norm=None, line_search='armijo', jac_options=None, 

262 **unknown_options): 

263 _check_unknown_options(unknown_options) 

264 

265 f_tol = fatol 

266 f_rtol = ftol 

267 x_tol = xatol 

268 x_rtol = xtol 

269 verbose = disp 

270 if jac_options is None: 

271 jac_options = dict() 

272 

273 jacobian = {'broyden1': nonlin.BroydenFirst, 

274 'broyden2': nonlin.BroydenSecond, 

275 'anderson': nonlin.Anderson, 

276 'linearmixing': nonlin.LinearMixing, 

277 'diagbroyden': nonlin.DiagBroyden, 

278 'excitingmixing': nonlin.ExcitingMixing, 

279 'krylov': nonlin.KrylovJacobian 

280 }[_method] 

281 

282 if args: 

283 if jac: 

284 def f(x): 

285 return fun(x, *args)[0] 

286 else: 

287 def f(x): 

288 return fun(x, *args) 

289 else: 

290 f = fun 

291 

292 x, info = nonlin.nonlin_solve(f, x0, jacobian=jacobian(**jac_options), 

293 iter=nit, verbose=verbose, 

294 maxiter=maxiter, f_tol=f_tol, 

295 f_rtol=f_rtol, x_tol=x_tol, 

296 x_rtol=x_rtol, tol_norm=tol_norm, 

297 line_search=line_search, 

298 callback=_callback, full_output=True, 

299 raise_exception=False) 

300 sol = OptimizeResult(x=x) 

301 sol.update(info) 

302 return sol 

303 

304def _root_broyden1_doc(): 

305 """ 

306 Options 

307 ------- 

308 nit : int, optional 

309 Number of iterations to make. If omitted (default), make as many 

310 as required to meet tolerances. 

311 disp : bool, optional 

312 Print status to stdout on every iteration. 

313 maxiter : int, optional 

314 Maximum number of iterations to make. If more are needed to 

315 meet convergence, `NoConvergence` is raised. 

316 ftol : float, optional 

317 Relative tolerance for the residual. If omitted, not used. 

318 fatol : float, optional 

319 Absolute tolerance (in max-norm) for the residual. 

320 If omitted, default is 6e-6. 

321 xtol : float, optional 

322 Relative minimum step size. If omitted, not used. 

323 xatol : float, optional 

324 Absolute minimum step size, as determined from the Jacobian 

325 approximation. If the step size is smaller than this, optimization 

326 is terminated as successful. If omitted, not used. 

327 tol_norm : function(vector) -> scalar, optional 

328 Norm to use in convergence check. Default is the maximum norm. 

329 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

330 Which type of a line search to use to determine the step size in 

331 the direction given by the Jacobian approximation. Defaults to 

332 'armijo'. 

333 jac_options : dict, optional 

334 Options for the respective Jacobian approximation. 

335 alpha : float, optional 

336 Initial guess for the Jacobian is (-1/alpha). 

337 reduction_method : str or tuple, optional 

338 Method used in ensuring that the rank of the Broyden 

339 matrix stays low. Can either be a string giving the 

340 name of the method, or a tuple of the form ``(method, 

341 param1, param2, ...)`` that gives the name of the 

342 method and values for additional parameters. 

343 

344 Methods available: 

345 

346 - ``restart`` 

347 Drop all matrix columns. Has no 

348 extra parameters. 

349 - ``simple`` 

350 Drop oldest matrix column. Has no 

351 extra parameters. 

352 - ``svd`` 

353 Keep only the most significant SVD 

354 components. 

355 

356 Extra parameters: 

357 

358 - ``to_retain`` 

359 Number of SVD components to 

360 retain when rank reduction is done. 

361 Default is ``max_rank - 2``. 

362 max_rank : int, optional 

363 Maximum rank for the Broyden matrix. 

364 Default is infinity (i.e., no rank reduction). 

365 """ 

366 pass 

367 

368def _root_broyden2_doc(): 

369 """ 

370 Options 

371 ------- 

372 nit : int, optional 

373 Number of iterations to make. If omitted (default), make as many 

374 as required to meet tolerances. 

375 disp : bool, optional 

376 Print status to stdout on every iteration. 

377 maxiter : int, optional 

378 Maximum number of iterations to make. If more are needed to 

379 meet convergence, `NoConvergence` is raised. 

380 ftol : float, optional 

381 Relative tolerance for the residual. If omitted, not used. 

382 fatol : float, optional 

383 Absolute tolerance (in max-norm) for the residual. 

384 If omitted, default is 6e-6. 

385 xtol : float, optional 

386 Relative minimum step size. If omitted, not used. 

387 xatol : float, optional 

388 Absolute minimum step size, as determined from the Jacobian 

389 approximation. If the step size is smaller than this, optimization 

390 is terminated as successful. If omitted, not used. 

391 tol_norm : function(vector) -> scalar, optional 

392 Norm to use in convergence check. Default is the maximum norm. 

393 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

394 Which type of a line search to use to determine the step size in 

395 the direction given by the Jacobian approximation. Defaults to 

396 'armijo'. 

397 jac_options : dict, optional 

398 Options for the respective Jacobian approximation. 

399 

400 alpha : float, optional 

401 Initial guess for the Jacobian is (-1/alpha). 

402 reduction_method : str or tuple, optional 

403 Method used in ensuring that the rank of the Broyden 

404 matrix stays low. Can either be a string giving the 

405 name of the method, or a tuple of the form ``(method, 

406 param1, param2, ...)`` that gives the name of the 

407 method and values for additional parameters. 

408 

409 Methods available: 

410 

411 - ``restart`` 

412 Drop all matrix columns. Has no 

413 extra parameters. 

414 - ``simple`` 

415 Drop oldest matrix column. Has no 

416 extra parameters. 

417 - ``svd`` 

418 Keep only the most significant SVD 

419 components. 

420 

421 Extra parameters: 

422 

423 - ``to_retain`` 

424 Number of SVD components to 

425 retain when rank reduction is done. 

426 Default is ``max_rank - 2``. 

427 max_rank : int, optional 

428 Maximum rank for the Broyden matrix. 

429 Default is infinity (i.e., no rank reduction). 

430 """ 

431 pass 

432 

433def _root_anderson_doc(): 

434 """ 

435 Options 

436 ------- 

437 nit : int, optional 

438 Number of iterations to make. If omitted (default), make as many 

439 as required to meet tolerances. 

440 disp : bool, optional 

441 Print status to stdout on every iteration. 

442 maxiter : int, optional 

443 Maximum number of iterations to make. If more are needed to 

444 meet convergence, `NoConvergence` is raised. 

445 ftol : float, optional 

446 Relative tolerance for the residual. If omitted, not used. 

447 fatol : float, optional 

448 Absolute tolerance (in max-norm) for the residual. 

449 If omitted, default is 6e-6. 

450 xtol : float, optional 

451 Relative minimum step size. If omitted, not used. 

452 xatol : float, optional 

453 Absolute minimum step size, as determined from the Jacobian 

454 approximation. If the step size is smaller than this, optimization 

455 is terminated as successful. If omitted, not used. 

456 tol_norm : function(vector) -> scalar, optional 

457 Norm to use in convergence check. Default is the maximum norm. 

458 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

459 Which type of a line search to use to determine the step size in 

460 the direction given by the Jacobian approximation. Defaults to 

461 'armijo'. 

462 jac_options : dict, optional 

463 Options for the respective Jacobian approximation. 

464 

465 alpha : float, optional 

466 Initial guess for the Jacobian is (-1/alpha). 

467 M : float, optional 

468 Number of previous vectors to retain. Defaults to 5. 

469 w0 : float, optional 

470 Regularization parameter for numerical stability. 

471 Compared to unity, good values of the order of 0.01. 

472 """ 

473 pass 

474 

475def _root_linearmixing_doc(): 

476 """ 

477 Options 

478 ------- 

479 nit : int, optional 

480 Number of iterations to make. If omitted (default), make as many 

481 as required to meet tolerances. 

482 disp : bool, optional 

483 Print status to stdout on every iteration. 

484 maxiter : int, optional 

485 Maximum number of iterations to make. If more are needed to 

486 meet convergence, ``NoConvergence`` is raised. 

487 ftol : float, optional 

488 Relative tolerance for the residual. If omitted, not used. 

489 fatol : float, optional 

490 Absolute tolerance (in max-norm) for the residual. 

491 If omitted, default is 6e-6. 

492 xtol : float, optional 

493 Relative minimum step size. If omitted, not used. 

494 xatol : float, optional 

495 Absolute minimum step size, as determined from the Jacobian 

496 approximation. If the step size is smaller than this, optimization 

497 is terminated as successful. If omitted, not used. 

498 tol_norm : function(vector) -> scalar, optional 

499 Norm to use in convergence check. Default is the maximum norm. 

500 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

501 Which type of a line search to use to determine the step size in 

502 the direction given by the Jacobian approximation. Defaults to 

503 'armijo'. 

504 jac_options : dict, optional 

505 Options for the respective Jacobian approximation. 

506 

507 alpha : float, optional 

508 initial guess for the jacobian is (-1/alpha). 

509 """ 

510 pass 

511 

512def _root_diagbroyden_doc(): 

513 """ 

514 Options 

515 ------- 

516 nit : int, optional 

517 Number of iterations to make. If omitted (default), make as many 

518 as required to meet tolerances. 

519 disp : bool, optional 

520 Print status to stdout on every iteration. 

521 maxiter : int, optional 

522 Maximum number of iterations to make. If more are needed to 

523 meet convergence, `NoConvergence` is raised. 

524 ftol : float, optional 

525 Relative tolerance for the residual. If omitted, not used. 

526 fatol : float, optional 

527 Absolute tolerance (in max-norm) for the residual. 

528 If omitted, default is 6e-6. 

529 xtol : float, optional 

530 Relative minimum step size. If omitted, not used. 

531 xatol : float, optional 

532 Absolute minimum step size, as determined from the Jacobian 

533 approximation. If the step size is smaller than this, optimization 

534 is terminated as successful. If omitted, not used. 

535 tol_norm : function(vector) -> scalar, optional 

536 Norm to use in convergence check. Default is the maximum norm. 

537 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

538 Which type of a line search to use to determine the step size in 

539 the direction given by the Jacobian approximation. Defaults to 

540 'armijo'. 

541 jac_options : dict, optional 

542 Options for the respective Jacobian approximation. 

543 

544 alpha : float, optional 

545 initial guess for the jacobian is (-1/alpha). 

546 """ 

547 pass 

548 

549def _root_excitingmixing_doc(): 

550 """ 

551 Options 

552 ------- 

553 nit : int, optional 

554 Number of iterations to make. If omitted (default), make as many 

555 as required to meet tolerances. 

556 disp : bool, optional 

557 Print status to stdout on every iteration. 

558 maxiter : int, optional 

559 Maximum number of iterations to make. If more are needed to 

560 meet convergence, `NoConvergence` is raised. 

561 ftol : float, optional 

562 Relative tolerance for the residual. If omitted, not used. 

563 fatol : float, optional 

564 Absolute tolerance (in max-norm) for the residual. 

565 If omitted, default is 6e-6. 

566 xtol : float, optional 

567 Relative minimum step size. If omitted, not used. 

568 xatol : float, optional 

569 Absolute minimum step size, as determined from the Jacobian 

570 approximation. If the step size is smaller than this, optimization 

571 is terminated as successful. If omitted, not used. 

572 tol_norm : function(vector) -> scalar, optional 

573 Norm to use in convergence check. Default is the maximum norm. 

574 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

575 Which type of a line search to use to determine the step size in 

576 the direction given by the Jacobian approximation. Defaults to 

577 'armijo'. 

578 jac_options : dict, optional 

579 Options for the respective Jacobian approximation. 

580 

581 alpha : float, optional 

582 Initial Jacobian approximation is (-1/alpha). 

583 alphamax : float, optional 

584 The entries of the diagonal Jacobian are kept in the range 

585 ``[alpha, alphamax]``. 

586 """ 

587 pass 

588 

589def _root_krylov_doc(): 

590 """ 

591 Options 

592 ------- 

593 nit : int, optional 

594 Number of iterations to make. If omitted (default), make as many 

595 as required to meet tolerances. 

596 disp : bool, optional 

597 Print status to stdout on every iteration. 

598 maxiter : int, optional 

599 Maximum number of iterations to make. If more are needed to 

600 meet convergence, `NoConvergence` is raised. 

601 ftol : float, optional 

602 Relative tolerance for the residual. If omitted, not used. 

603 fatol : float, optional 

604 Absolute tolerance (in max-norm) for the residual. 

605 If omitted, default is 6e-6. 

606 xtol : float, optional 

607 Relative minimum step size. If omitted, not used. 

608 xatol : float, optional 

609 Absolute minimum step size, as determined from the Jacobian 

610 approximation. If the step size is smaller than this, optimization 

611 is terminated as successful. If omitted, not used. 

612 tol_norm : function(vector) -> scalar, optional 

613 Norm to use in convergence check. Default is the maximum norm. 

614 line_search : {None, 'armijo' (default), 'wolfe'}, optional 

615 Which type of a line search to use to determine the step size in 

616 the direction given by the Jacobian approximation. Defaults to 

617 'armijo'. 

618 jac_options : dict, optional 

619 Options for the respective Jacobian approximation. 

620 

621 rdiff : float, optional 

622 Relative step size to use in numerical differentiation. 

623 method : {'lgmres', 'gmres', 'bicgstab', 'cgs', 'minres'} or function 

624 Krylov method to use to approximate the Jacobian. 

625 Can be a string, or a function implementing the same 

626 interface as the iterative solvers in 

627 `scipy.sparse.linalg`. 

628 

629 The default is `scipy.sparse.linalg.lgmres`. 

630 inner_M : LinearOperator or InverseJacobian 

631 Preconditioner for the inner Krylov iteration. 

632 Note that you can use also inverse Jacobians as (adaptive) 

633 preconditioners. For example, 

634 

635 >>> jac = BroydenFirst() 

636 >>> kjac = KrylovJacobian(inner_M=jac.inverse). 

637 

638 If the preconditioner has a method named 'update', it will 

639 be called as ``update(x, f)`` after each nonlinear step, 

640 with ``x`` giving the current point, and ``f`` the current 

641 function value. 

642 inner_tol, inner_maxiter, ... 

643 Parameters to pass on to the "inner" Krylov solver. 

644 See `scipy.sparse.linalg.gmres` for details. 

645 outer_k : int, optional 

646 Size of the subspace kept across LGMRES nonlinear 

647 iterations. 

648 

649 See `scipy.sparse.linalg.lgmres` for details. 

650 """ 

651 pass