Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(RootFinderDataType), | intent(inout) | :: | RootFinderData |
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
REAL(r64) FUNCTION BrentMethod( RootFinderData )
! FUNCTION INFORMATION:
! AUTHOR Dimitri Curtil (LBNL)
! DATE WRITTEN March 2006
! MODIFIED
! RE-ENGINEERED na
! PURPOSE OF THIS FUNCTION:
! This function computes the next iterate using the Brent's method.
! If new iterate does not lie within the lower and upper points then
! the secant method is used instead.
!
! Convergence rate is at best quadratic.
!
! PRECONDITION:
! Lower and upper points must be defined and distinct.
!
! POSTCONDITION:
! - LowerPoint%X < XCandidate < UpperPoint%X
! - RootFinderData%CurrentMethodType update with current solution method.
!
! METHODOLOGY EMPLOYED:
! Inverse quadratic interpolation using the last 3 best iterates.
! The next root estimate is x = B + P/Q whereby B is the current best estimate
! of the root.
!
! REFERENCES:
! na
! USE STATEMENTS:
! na
IMPLICIT NONE ! Enforce explicit typing of all variables in this routine
! FUNCTION ARGUMENT DEFINITIONS:
TYPE(RootFinderDataType), INTENT(INOUT) :: RootFinderData ! Data used by root finding algorithm
! FUNCTION PARAMETER DEFINITIONS:
! na
! INTERFACE BLOCK SPECIFICATIONS
! na
! DERIVED TYPE DEFINITIONS
! na
! FUNCTION LOCAL VARIABLE DECLARATIONS:
REAL(r64) :: XCandidate
REAL(r64) :: A, FA, B, FB, C, FC
REAL(r64) :: R, S, T, P, Q
! FLOW:
! Only attempt Brent's method if enough history points are available
! and if the root finder is converging (not diverging) since the previous
! iterate.
!
! We assume that;
! - the root is bracketed between the lower and upper points (see AdvanceRootFinder() ).
! - there are at least 3 history points
IF ( RootFinderData%NumHistory == 3 ) THEN
A = RootFinderData%History(2)%X
FA = RootFinderData%History(2)%Y
B = RootFinderData%History(1)%X
FB = RootFinderData%History(1)%Y
C = RootFinderData%History(3)%X
FC = RootFinderData%History(3)%Y
! Should never happen if CheckRootFinderConvergence() is invoked prior to this subroutine
IF (FC == 0.0d0) THEN
BrentMethod = C
RETURN
! Should never happen if CheckRootFinderConvergence() is invoked prior to this subroutine
ELSE IF (FA == 0.0d0) THEN
BrentMethod = A
RETURN
ELSE
R = FB/FC
S = FB/FA
T = FA/FC
P = S*(T*(R-T)*(C-B)-(1.0d0-R)*(B-A))
Q = (T-1.0d0)*(R-1.0d0)*(S-1.0d0)
! Only accept correction if it is small enough (75% of previous increment)
IF ( ABS(P) <= 0.75d0*ABS(Q*RootFinderData%Increment%X) ) THEN
RootFinderData%CurrentMethodType = iMethodBrent
XCandidate = B + P/Q
! Check that new candidate is within range and brackets
IF ( .NOT.CheckRootFinderCandidate(RootFinderData, XCandidate) ) THEN
! Recovery method
XCandidate = FalsePositionMethod( RootFinderData )
END IF
ELSE
! Recover from bad correction with bisection
! Biscetion produced the best numerical performance in testing compared to
! - Secant
! - False position (very slow recovery)
XCandidate = BisectionMethod( RootFinderData )
END IF
END IF
ELSE
! Not enough history to try Brent's method yet: use Secant's method
XCandidate = SecantMethod( RootFinderData )
END IF
BrentMethod = XCandidate
RETURN
END FUNCTION BrentMethod