Create the actual Adams-Moulton integrator: initialize the b coefficients.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| class(integrator_adams_moulton), | intent(inout) | :: | self | Integrator. |
||
| character(len=*), | intent(in) | :: | scheme | Selected scheme. |
||
| integer(kind=I_P), | intent(in), | optional | :: | iterations | Implicit iterations. |
|
| logical, | intent(in), | optional | :: | autoupdate | Enable cyclic autoupdate of previous time steps. |
|
| class(integrand_object), | intent(in), | optional | :: | U | Integrand molding prototype. |
|
| logical, | intent(in), | optional | :: | stop_on_fail | Stop execution if initialization fail. |
subroutine initialize(self, scheme, iterations, autoupdate, U, stop_on_fail)
!< Create the actual Adams-Moulton integrator: initialize the *b* coefficients.
class(integrator_adams_moulton), intent(inout) :: self !< Integrator.
character(*), intent(in) :: scheme !< Selected scheme.
integer(I_P), intent(in), optional :: iterations !< Implicit iterations.
logical, intent(in), optional :: autoupdate !< Enable cyclic autoupdate of previous time steps.
class(integrand_object), intent(in), optional :: U !< Integrand molding prototype.
logical, intent(in), optional :: stop_on_fail !< Stop execution if initialization fail.
if (self%is_supported(scheme=scheme)) then
call self%destroy
self%description_ = trim(adjustl(scheme))
select case(trim(adjustl(scheme)))
case('adams_moulton_0')
self%steps = 0 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 1.0_R_P
case('adams_moulton_1')
self%steps = 1 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 1.0_R_P/2.0_R_P
self%b(1) = 1.0_R_P/2.0_R_P
case('adams_moulton_2')
self%steps = 2 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -1.0_R_P/12.0_R_P
self%b(1) = 8.0_R_P/12.0_R_P
self%b(2) = 5.0_R_P/12.0_R_P
case('adams_moulton_3')
self%steps = 3 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 1.0_R_P/24.0_R_P
self%b(1) = -5.0_R_P/24.0_R_P
self%b(2) = 19.0_R_P/24.0_R_P
self%b(3) = 9.0_R_P/24.0_R_P
case('adams_moulton_4')
self%steps = 4 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -19.0_R_P/720.0_R_P
self%b(1) = 106.0_R_P/720.0_R_P
self%b(2) = -264.0_R_P/720.0_R_P
self%b(3) = 646.0_R_P/720.0_R_P
self%b(4) = 251.0_R_P/720.0_R_P
case('adams_moulton_5')
self%steps = 5 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 27.0_R_P/1440.0_R_P
self%b(1) = -173.0_R_P/1440.0_R_P
self%b(2) = 482.0_R_P/1440.0_R_P
self%b(3) = -798.0_R_P/1440.0_R_P
self%b(4) = 1427.0_R_P/1440.0_R_P
self%b(5) = 475.0_R_P/1440.0_R_P
case('adams_moulton_6')
self%steps = 6 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -863.0_R_P/60480.0_R_P
self%b(1) = 6312.0_R_P/60480.0_R_P
self%b(2) = -20211.0_R_P/60480.0_R_P
self%b(3) = 37504.0_R_P/60480.0_R_P
self%b(4) = -46461.0_R_P/60480.0_R_P
self%b(5) = 65112.0_R_P/60480.0_R_P
self%b(6) = 19087.0_R_P/60480.0_R_P
case('adams_moulton_7')
self%steps = 7 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 1375.0_R_P/120960.0_R_P
self%b(1) = -11351.0_R_P/120960.0_R_P
self%b(2) = 41499.0_R_P/120960.0_R_P
self%b(3) = -88547.0_R_P/120960.0_R_P
self%b(4) = 123133.0_R_P/120960.0_R_P
self%b(5) = -121797.0_R_P/120960.0_R_P
self%b(6) = 139849.0_R_P/120960.0_R_P
self%b(7) = 36799.0_R_P/120960.0_R_P
case('adams_moulton_8')
self%steps = 8 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -33953.0_R_P/3628800.0_R_P
self%b(1) = 312874.0_R_P/3628800.0_R_P
self%b(2) = -1291214.0_R_P/3628800.0_R_P
self%b(3) = 3146338.0_R_P/3628800.0_R_P
self%b(4) = -5033120.0_R_P/3628800.0_R_P
self%b(5) = 5595358.0_R_P/3628800.0_R_P
self%b(6) = -4604594.0_R_P/3628800.0_R_P
self%b(7) = 4467094.0_R_P/3628800.0_R_P
self%b(8) = 1070017.0_R_P/3628800.0_R_P
case('adams_moulton_9')
self%steps = 9 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 57281.0_R_P/7257600.0_R_P
self%b(1) = -583435.0_R_P/7257600.0_R_P
self%b(2) = 2687864.0_R_P/7257600.0_R_P
self%b(3) = -7394032.0_R_P/7257600.0_R_P
self%b(4) = 13510082.0_R_P/7257600.0_R_P
self%b(5) = -17283646.0_R_P/7257600.0_R_P
self%b(6) = 16002320.0_R_P/7257600.0_R_P
self%b(7) = -11271304.0_R_P/7257600.0_R_P
self%b(8) = 9449717.0_R_P/7257600.0_R_P
self%b(9) = 2082753.0_R_P/7257600.0_R_P
case('adams_moulton_10')
self%steps = 10 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -3250433.0_R_P/479001600.0_R_P
self%b(1) = 36284876.0_R_P/479001600.0_R_P
self%b(2) = -184776195.0_R_P/479001600.0_R_P
self%b(3) = 567450984.0_R_P/479001600.0_R_P
self%b(4) = -1170597042.0_R_P/479001600.0_R_P
self%b(5) = 1710774528.0_R_P/479001600.0_R_P
self%b(6) = -1823311566.0_R_P/479001600.0_R_P
self%b(7) = 1446205080.0_R_P/479001600.0_R_P
self%b(8) = -890175549.0_R_P/479001600.0_R_P
self%b(9) = 656185652.0_R_P/479001600.0_R_P
self%b(10) = 134211265.0_R_P/479001600.0_R_P
case('adams_moulton_11')
self%steps = 11 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 5675265.0_R_P/958003200.0_R_P
self%b(1) = -68928781.0_R_P/958003200.0_R_P
self%b(2) = 384709327.0_R_P/958003200.0_R_P
self%b(3) = -1305971115.0_R_P/958003200.0_R_P
self%b(4) = 3007739418.0_R_P/958003200.0_R_P
self%b(5) = -4963166514.0_R_P/958003200.0_R_P
self%b(6) = 6043521486.0_R_P/958003200.0_R_P
self%b(7) = -5519460582.0_R_P/958003200.0_R_P
self%b(8) = 3828828885.0_R_P/958003200.0_R_P
self%b(9) = -2092490673.0_R_P/958003200.0_R_P
self%b(10) = 1374799219.0_R_P/958003200.0_R_P
self%b(11) = 262747265.0_R_P/958003200.0_R_P
case('adams_moulton_12')
self%steps = 12 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -13695779093.0_R_P/2615348736000.0_R_P
self%b(1) = 179842822566.0_R_P/2615348736000.0_R_P
self%b(2) = -1092096992268.0_R_P/2615348736000.0_R_P
self%b(3) = 4063327863170.0_R_P/2615348736000.0_R_P
self%b(4) = -10344711794985.0_R_P/2615348736000.0_R_P
self%b(5) = 19058185652796.0_R_P/2615348736000.0_R_P
self%b(6) = -26204344465152.0_R_P/2615348736000.0_R_P
self%b(7) = 27345870698436.0_R_P/2615348736000.0_R_P
self%b(8) = -21847538039895.0_R_P/2615348736000.0_R_P
self%b(9) = 13465774256510.0_R_P/2615348736000.0_R_P
self%b(10) = -6616420957428.0_R_P/2615348736000.0_R_P
self%b(11) = 3917551216986.0_R_P/2615348736000.0_R_P
self%b(12) = 703604254357.0_R_P/2615348736000.0_R_P
case('adams_moulton_13')
self%steps = 13 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 24466579093.0_R_P/5230697472000.0_R_P
self%b(1) = -345457086395.0_R_P/5230697472000.0_R_P
self%b(2) = 2268078814386.0_R_P/5230697472000.0_R_P
self%b(3) = -9181635605134.0_R_P/5230697472000.0_R_P
self%b(4) = 25620259777835.0_R_P/5230697472000.0_R_P
self%b(5) = -52177910882661.0_R_P/5230697472000.0_R_P
self%b(6) = 80101021029180.0_R_P/5230697472000.0_R_P
self%b(7) = -94393338653892.0_R_P/5230697472000.0_R_P
self%b(8) = 86180228689563.0_R_P/5230697472000.0_R_P
self%b(9) = -61188680131285.0_R_P/5230697472000.0_R_P
self%b(10) = 33928990133618.0_R_P/5230697472000.0_R_P
self%b(11) = -15141235084110.0_R_P/5230697472000.0_R_P
self%b(12) = 8153167962181.0_R_P/5230697472000.0_R_P
self%b(13) = 1382741929621.0_R_P/5230697472000.0_R_P
case('adams_moulton_14')
self%steps = 14 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = -132282840127.0_R_P/31384184832000.0_R_P
self%b(1) = 1998759236336.0_R_P/31384184832000.0_R_P
self%b(2) = -14110480969927.0_R_P/31384184832000.0_R_P
self%b(3) = 61759426692544.0_R_P/31384184832000.0_R_P
self%b(4) = -187504936597931.0_R_P/31384184832000.0_R_P
self%b(5) = 418551804601264.0_R_P/31384184832000.0_R_P
self%b(6) = -710312834197347.0_R_P/31384184832000.0_R_P
self%b(7) = 934600833490944.0_R_P/31384184832000.0_R_P
self%b(8) = -963605400824733.0_R_P/31384184832000.0_R_P
self%b(9) = 781911618071632.0_R_P/31384184832000.0_R_P
self%b(10) = -499547203754837.0_R_P/31384184832000.0_R_P
self%b(11) = 251724894607936.0_R_P/31384184832000.0_R_P
self%b(12) = -102885148956217.0_R_P/31384184832000.0_R_P
self%b(13) = 50770967534864.0_R_P/31384184832000.0_R_P
self%b(14) = 8164168737599.0_R_P/31384184832000.0_R_P
case('adams_moulton_15')
self%steps = 15 ; allocate(self%b(0:self%steps)) ; self%b = 0.0_R_P
self%b(0) = 240208245823.0_R_P/62768369664000.0_R_P
self%b(1) = -3867689367599.0_R_P/62768369664000.0_R_P
self%b(2) = 29219384284087.0_R_P/62768369664000.0_R_P
self%b(3) = -137515713789319.0_R_P/62768369664000.0_R_P
self%b(4) = 451403108933483.0_R_P/62768369664000.0_R_P
self%b(5) = -1096355235402331.0_R_P/62768369664000.0_R_P
self%b(6) = 2039345879546643.0_R_P/62768369664000.0_R_P
self%b(7) = -2966365730265699.0_R_P/62768369664000.0_R_P
self%b(8) = 3414941728852893.0_R_P/62768369664000.0_R_P
self%b(9) = -3129453071993581.0_R_P/62768369664000.0_R_P
self%b(10) = 2285168598349733.0_R_P/62768369664000.0_R_P
self%b(11) = -1326978663058069.0_R_P/62768369664000.0_R_P
self%b(12) = 612744541065337.0_R_P/62768369664000.0_R_P
self%b(13) = -230992163723849.0_R_P/62768369664000.0_R_P
self%b(14) = 105145058757073.0_R_P/62768369664000.0_R_P
self%b(15) = 16088129229375.0_R_P/62768369664000.0_R_P
endselect
self%autoupdate = .true. ; if (present(autoupdate)) self%autoupdate = autoupdate
self%iterations = 1 ; if (present(iterations)) self%iterations = iterations
self%registers = self%steps
if (present(U)) call self%allocate_integrand_members(U=U)
else
call self%trigger_error(error=ERROR_UNSUPPORTED_SCHEME, &
error_message='"'//trim(adjustl(scheme))//'" unsupported scheme', &
is_severe=stop_on_fail)
endif
endsubroutine initialize