ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_optimization.F90
Go to the documentation of this file.
1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
3!!!! !!!!
4!!!! ParaMonte: Parallel Monte Carlo and Machine Learning Library. !!!!
5!!!! !!!!
6!!!! Copyright (C) 2012-present, The Computational Data Science Lab !!!!
7!!!! !!!!
8!!!! This file is part of the ParaMonte library. !!!!
9!!!! !!!!
10!!!! LICENSE !!!!
11!!!! !!!!
12!!!! https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md !!!!
13!!!! !!!!
14!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
15!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
16
54! <li> [getMinSimplex](@ref pm_optimization::getMinSimplex) and [setMinSimplex](@ref pm_optimization::setMinSimplex) perform unconstrained derivative-free minimization of **multivariate** functions.<br>
55! The **Nelder–Mead method** (also **downhill simplex method**, **amoeba method**, or **polytope method**)
56! is a numerical method used to find the minimum or maximum of an objective function in a **multidimensional** space.<br>
57! It is a direct search method (based on function comparison) and is often applied to **nonlinear optimization problems** for which **derivatives are unknown**.<br>
58! The Nelder–Mead technique is a **heuristic search method** that can converge to **non-stationary points**.<br>
59! The Nelder–Mead technique was proposed by [John Nelder](https://en.wikipedia.org/wiki/John_Nelder) and [Roger Mead](https://en.wikipedia.org/wiki/Roger_Mead) in 1965.<br>
60! According to Numerical Recipes by Press et al. (1992):<br>
61! > The Nelder–Mead method crawls downhill in a straightforward fashion that makes almost no special assumptions about the target function.<br>
62! > This can be extremely slow, but it can also, in some cases, be extremely robust.<br>
63! > The storage requirement for this method is of order \f$\ms{ndim}^2\f$, where \f$\ms{ndim}\f$ is the number of dimensions of the function support.<br>
64! > As such, the Nedler-Mead Simplex method is most useful when the minimization tasks are minimal and isolated.<br>
78
79!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
80
82
83 use pm_kind, only: SK, IK, LK
84
85 implicit none
86
87 character(*, SK), parameter :: MODULE_NAME = "@pm_optimization"
88
89!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
90
151 interface isBracketMax
152
153 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
154
155#if RK5_ENABLED
156 pure elemental module function isBracketMax_RK5(xmax, xlow, xupp, fmax, flow, fupp) result(bracketed)
157#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
158 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMax_RK5
159#endif
160 use pm_kind, only: RKG => RK5
161 real(RKG) , intent(in) :: xmax, xlow, xupp, fmax, flow, fupp
162 logical(LK) :: bracketed
163 end function
164#endif
165
166#if RK4_ENABLED
167 pure elemental module function isBracketMax_RK4(xmax, xlow, xupp, fmax, flow, fupp) result(bracketed)
168#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
169 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMax_RK4
170#endif
171 use pm_kind, only: RKG => RK4
172 real(RKG) , intent(in) :: xmax, xlow, xupp, fmax, flow, fupp
173 logical(LK) :: bracketed
174 end function
175#endif
176
177#if RK3_ENABLED
178 pure elemental module function isBracketMax_RK3(xmax, xlow, xupp, fmax, flow, fupp) result(bracketed)
179#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
180 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMax_RK3
181#endif
182 use pm_kind, only: RKG => RK3
183 real(RKG) , intent(in) :: xmax, xlow, xupp, fmax, flow, fupp
184 logical(LK) :: bracketed
185 end function
186#endif
187
188#if RK2_ENABLED
189 pure elemental module function isBracketMax_RK2(xmax, xlow, xupp, fmax, flow, fupp) result(bracketed)
190#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
191 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMax_RK2
192#endif
193 use pm_kind, only: RKG => RK2
194 real(RKG) , intent(in) :: xmax, xlow, xupp, fmax, flow, fupp
195 logical(LK) :: bracketed
196 end function
197#endif
198
199#if RK1_ENABLED
200 pure elemental module function isBracketMax_RK1(xmax, xlow, xupp, fmax, flow, fupp) result(bracketed)
201#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
202 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMax_RK1
203#endif
204 use pm_kind, only: RKG => RK1
205 real(RKG) , intent(in) :: xmax, xlow, xupp, fmax, flow, fupp
206 logical(LK) :: bracketed
207 end function
208#endif
209
210 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
211
212 end interface
213
214!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
215
276 interface isBracketMin
277
278 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
279
280#if RK5_ENABLED
281 pure elemental module function isBracketMin_RK5(xmin, xlow, xupp, fmin, flow, fupp) result(bracketed)
282#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
283 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMin_RK5
284#endif
285 use pm_kind, only: RKG => RK5
286 real(RKG) , intent(in) :: xmin, xlow, xupp, fmin, flow, fupp
287 logical(LK) :: bracketed
288 end function
289#endif
290
291#if RK4_ENABLED
292 pure elemental module function isBracketMin_RK4(xmin, xlow, xupp, fmin, flow, fupp) result(bracketed)
293#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
294 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMin_RK4
295#endif
296 use pm_kind, only: RKG => RK4
297 real(RKG) , intent(in) :: xmin, xlow, xupp, fmin, flow, fupp
298 logical(LK) :: bracketed
299 end function
300#endif
301
302#if RK3_ENABLED
303 pure elemental module function isBracketMin_RK3(xmin, xlow, xupp, fmin, flow, fupp) result(bracketed)
304#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
305 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMin_RK3
306#endif
307 use pm_kind, only: RKG => RK3
308 real(RKG) , intent(in) :: xmin, xlow, xupp, fmin, flow, fupp
309 logical(LK) :: bracketed
310 end function
311#endif
312
313#if RK2_ENABLED
314 pure elemental module function isBracketMin_RK2(xmin, xlow, xupp, fmin, flow, fupp) result(bracketed)
315#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
316 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMin_RK2
317#endif
318 use pm_kind, only: RKG => RK2
319 real(RKG) , intent(in) :: xmin, xlow, xupp, fmin, flow, fupp
320 logical(LK) :: bracketed
321 end function
322#endif
323
324#if RK1_ENABLED
325 pure elemental module function isBracketMin_RK1(xmin, xlow, xupp, fmin, flow, fupp) result(bracketed)
326#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
327 !DEC$ ATTRIBUTES DLLEXPORT :: isBracketMin_RK1
328#endif
329 use pm_kind, only: RKG => RK1
330 real(RKG) , intent(in) :: xmin, xlow, xupp, fmin, flow, fupp
331 logical(LK) :: bracketed
332 end function
333#endif
334
335 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
336
337 end interface
338
339!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
340
434
435 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
436
437#if RK5_ENABLED
438 recursive module subroutine setBracketMax_RK5(getFunc, niter, xmax, xlow, xupp, fmax, flow, fupp)
439#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
440 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMax_RK5
441#endif
442 use pm_kind, only: RKG => RK5
443 procedure(real(RKG)) :: getFunc
444 integer(IK) , intent(inout) :: niter
445 real(RKG) , intent(inout) :: xlow, xupp
446 real(RKG) , intent(out) :: xmax, fmax
447 real(RKG) , intent(out) , optional :: flow, fupp
448 end subroutine
449#endif
450
451#if RK4_ENABLED
452 recursive module subroutine setBracketMax_RK4(getFunc, niter, xmax, xlow, xupp, fmax, flow, fupp)
453#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
454 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMax_RK4
455#endif
456 use pm_kind, only: RKG => RK4
457 procedure(real(RKG)) :: getFunc
458 integer(IK) , intent(inout) :: niter
459 real(RKG) , intent(inout) :: xlow, xupp
460 real(RKG) , intent(out) :: xmax, fmax
461 real(RKG) , intent(out) , optional :: flow, fupp
462 end subroutine
463#endif
464
465#if RK3_ENABLED
466 recursive module subroutine setBracketMax_RK3(getFunc, niter, xmax, xlow, xupp, fmax, flow, fupp)
467#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
468 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMax_RK3
469#endif
470 use pm_kind, only: RKG => RK3
471 procedure(real(RKG)) :: getFunc
472 integer(IK) , intent(inout) :: niter
473 real(RKG) , intent(inout) :: xlow, xupp
474 real(RKG) , intent(out) :: xmax, fmax
475 real(RKG) , intent(out) , optional :: flow, fupp
476 end subroutine
477#endif
478
479#if RK2_ENABLED
480 recursive module subroutine setBracketMax_RK2(getFunc, niter, xmax, xlow, xupp, fmax, flow, fupp)
481#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
482 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMax_RK2
483#endif
484 use pm_kind, only: RKG => RK2
485 procedure(real(RKG)) :: getFunc
486 integer(IK) , intent(inout) :: niter
487 real(RKG) , intent(inout) :: xlow, xupp
488 real(RKG) , intent(out) :: xmax, fmax
489 real(RKG) , intent(out) , optional :: flow, fupp
490 end subroutine
491#endif
492
493#if RK1_ENABLED
494 recursive module subroutine setBracketMax_RK1(getFunc, niter, xmax, xlow, xupp, fmax, flow, fupp)
495#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
496 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMax_RK1
497#endif
498 use pm_kind, only: RKG => RK1
499 procedure(real(RKG)) :: getFunc
500 integer(IK) , intent(inout) :: niter
501 real(RKG) , intent(inout) :: xlow, xupp
502 real(RKG) , intent(out) :: xmax, fmax
503 real(RKG) , intent(out) , optional :: flow, fupp
504 end subroutine
505#endif
506
507 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
508
509 end interface
510
511!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
512
606
607 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
608
609#if RK5_ENABLED
610 recursive module subroutine setBracketMin_RK5(getFunc, niter, xmin, xlow, xupp, fmin, flow, fupp)
611#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
612 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMin_RK5
613#endif
614 use pm_kind, only: RKG => RK5
615 procedure(real(RKG)) :: getFunc
616 integer(IK) , intent(inout) :: niter
617 real(RKG) , intent(inout) :: xlow, xupp
618 real(RKG) , intent(out) :: xmin, fmin
619 real(RKG) , intent(out) , optional :: flow, fupp
620 end subroutine
621#endif
622
623#if RK4_ENABLED
624 recursive module subroutine setBracketMin_RK4(getFunc, niter, xmin, xlow, xupp, fmin, flow, fupp)
625#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
626 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMin_RK4
627#endif
628 use pm_kind, only: RKG => RK4
629 procedure(real(RKG)) :: getFunc
630 integer(IK) , intent(inout) :: niter
631 real(RKG) , intent(inout) :: xlow, xupp
632 real(RKG) , intent(out) :: xmin, fmin
633 real(RKG) , intent(out) , optional :: flow, fupp
634 end subroutine
635#endif
636
637#if RK3_ENABLED
638 recursive module subroutine setBracketMin_RK3(getFunc, niter, xmin, xlow, xupp, fmin, flow, fupp)
639#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
640 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMin_RK3
641#endif
642 use pm_kind, only: RKG => RK3
643 procedure(real(RKG)) :: getFunc
644 integer(IK) , intent(inout) :: niter
645 real(RKG) , intent(inout) :: xlow, xupp
646 real(RKG) , intent(out) :: xmin, fmin
647 real(RKG) , intent(out) , optional :: flow, fupp
648 end subroutine
649#endif
650
651#if RK2_ENABLED
652 recursive module subroutine setBracketMin_RK2(getFunc, niter, xmin, xlow, xupp, fmin, flow, fupp)
653#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
654 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMin_RK2
655#endif
656 use pm_kind, only: RKG => RK2
657 procedure(real(RKG)) :: getFunc
658 integer(IK) , intent(inout) :: niter
659 real(RKG) , intent(inout) :: xlow, xupp
660 real(RKG) , intent(out) :: xmin, fmin
661 real(RKG) , intent(out) , optional :: flow, fupp
662 end subroutine
663#endif
664
665#if RK1_ENABLED
666 recursive module subroutine setBracketMin_RK1(getFunc, niter, xmin, xlow, xupp, fmin, flow, fupp)
667#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
668 !DEC$ ATTRIBUTES DLLEXPORT :: setBracketMin_RK1
669#endif
670 use pm_kind, only: RKG => RK1
671 procedure(real(RKG)) :: getFunc
672 integer(IK) , intent(inout) :: niter
673 real(RKG) , intent(inout) :: xlow, xupp
674 real(RKG) , intent(out) :: xmin, fmin
675 real(RKG) , intent(out) , optional :: flow, fupp
676 end subroutine
677#endif
678
679 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
680
681 end interface
682
683!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
684
768 interface getMinBrent
769
770 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
771
772#if RK5_ENABLED
773 recursive module function getMinBrent_RK5(getFunc, xlow, xupp, fmin, tol, niter) result(xmin)
774#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
775 !DEC$ ATTRIBUTES DLLEXPORT :: getMinBrent_RK5
776#endif
777 use pm_kind, only: RKG => RK5
778 procedure(real(RKG)) :: getFunc
779 integer(IK) , intent(inout) , optional :: niter
780 real(RKG) , intent(in) , optional :: xlow, xupp, tol
781 real(RKG) , intent(out) , optional :: fmin
782 real(RKG) :: xmin
783 end function
784#endif
785
786#if RK4_ENABLED
787 recursive module function getMinBrent_RK4(getFunc, xlow, xupp, fmin, tol, niter) result(xmin)
788#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
789 !DEC$ ATTRIBUTES DLLEXPORT :: getMinBrent_RK4
790#endif
791 use pm_kind, only: RKG => RK4
792 procedure(real(RKG)) :: getFunc
793 integer(IK) , intent(inout) , optional :: niter
794 real(RKG) , intent(in) , optional :: xlow, xupp, tol
795 real(RKG) , intent(out) , optional :: fmin
796 real(RKG) :: xmin
797 end function
798#endif
799
800#if RK3_ENABLED
801 recursive module function getMinBrent_RK3(getFunc, xlow, xupp, fmin, tol, niter) result(xmin)
802#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
803 !DEC$ ATTRIBUTES DLLEXPORT :: getMinBrent_RK3
804#endif
805 use pm_kind, only: RKG => RK3
806 procedure(real(RKG)) :: getFunc
807 integer(IK) , intent(inout) , optional :: niter
808 real(RKG) , intent(in) , optional :: xlow, xupp, tol
809 real(RKG) , intent(out) , optional :: fmin
810 real(RKG) :: xmin
811 end function
812#endif
813
814#if RK2_ENABLED
815 recursive module function getMinBrent_RK2(getFunc, xlow, xupp, fmin, tol, niter) result(xmin)
816#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
817 !DEC$ ATTRIBUTES DLLEXPORT :: getMinBrent_RK2
818#endif
819 use pm_kind, only: RKG => RK2
820 procedure(real(RKG)) :: getFunc
821 integer(IK) , intent(inout) , optional :: niter
822 real(RKG) , intent(in) , optional :: xlow, xupp, tol
823 real(RKG) , intent(out) , optional :: fmin
824 real(RKG) :: xmin
825 end function
826#endif
827
828#if RK1_ENABLED
829 recursive module function getMinBrent_RK1(getFunc, xlow, xupp, fmin, tol, niter) result(xmin)
830#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
831 !DEC$ ATTRIBUTES DLLEXPORT :: getMinBrent_RK1
832#endif
833 use pm_kind, only: RKG => RK1
834 procedure(real(RKG)) :: getFunc
835 integer(IK) , intent(inout) , optional :: niter
836 real(RKG) , intent(in) , optional :: xlow, xupp, tol
837 real(RKG) , intent(out) , optional :: fmin
838 real(RKG) :: xmin
839 end function
840#endif
841
842 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
843
844 end interface
845
846!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
847
935 interface setMinBrent
936
937 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
938
939#if RK5_ENABLED
940 recursive module subroutine setMinBrent_RK5(getFunc, xmin, xlow, xupp, fmin, tol, niter)
941#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
942 !DEC$ ATTRIBUTES DLLEXPORT :: setMinBrent_RK5
943#endif
944 use pm_kind, only: RKG => RK5
945 procedure(real(RKG)) :: getFunc
946 real(RKG) , value :: xlow, xupp, tol
947 real(RKG) , intent(inout) :: fmin, xmin
948 integer(IK) , intent(inout) :: niter
949 end subroutine
950#endif
951
952#if RK4_ENABLED
953 recursive module subroutine setMinBrent_RK4(getFunc, xmin, xlow, xupp, fmin, tol, niter)
954#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
955 !DEC$ ATTRIBUTES DLLEXPORT :: setMinBrent_RK4
956#endif
957 use pm_kind, only: RKG => RK4
958 procedure(real(RKG)) :: getFunc
959 real(RKG) , value :: xlow, xupp, tol
960 real(RKG) , intent(inout) :: fmin, xmin
961 integer(IK) , intent(inout) :: niter
962 end subroutine
963#endif
964
965#if RK3_ENABLED
966 recursive module subroutine setMinBrent_RK3(getFunc, xmin, xlow, xupp, fmin, tol, niter)
967#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
968 !DEC$ ATTRIBUTES DLLEXPORT :: setMinBrent_RK3
969#endif
970 use pm_kind, only: RKG => RK3
971 procedure(real(RKG)) :: getFunc
972 real(RKG) , value :: xlow, xupp, tol
973 real(RKG) , intent(inout) :: fmin, xmin
974 integer(IK) , intent(inout) :: niter
975 end subroutine
976#endif
977
978#if RK2_ENABLED
979 recursive module subroutine setMinBrent_RK2(getFunc, xmin, xlow, xupp, fmin, tol, niter)
980#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
981 !DEC$ ATTRIBUTES DLLEXPORT :: setMinBrent_RK2
982#endif
983 use pm_kind, only: RKG => RK2
984 procedure(real(RKG)) :: getFunc
985 real(RKG) , value :: xlow, xupp, tol
986 real(RKG) , intent(inout) :: fmin, xmin
987 integer(IK) , intent(inout) :: niter
988 end subroutine
989#endif
990
991#if RK1_ENABLED
992 recursive module subroutine setMinBrent_RK1(getFunc, xmin, xlow, xupp, fmin, tol, niter)
993#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
994 !DEC$ ATTRIBUTES DLLEXPORT :: setMinBrent_RK1
995#endif
996 use pm_kind, only: RKG => RK1
997 procedure(real(RKG)) :: getFunc
998 real(RKG) , value :: xlow, xupp, tol
999 real(RKG) , intent(inout) :: fmin, xmin
1000 integer(IK) , intent(inout) :: niter
1001 end subroutine
1002#endif
1003
1004 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1005
1006 end interface
1007
1008!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1009
1090
1091 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1092
1093#if RK5_ENABLED
1094 recursive module function isFailedMinPowell_RK5(getFunc, xmin, fmin, tol, niter) result(failed)
1095#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1096 !DEC$ ATTRIBUTES DLLEXPORT :: isFailedMinPowell_RK5
1097#endif
1098 use pm_kind, only: RKG => RK5
1099 procedure(real(RKG)) :: getFunc
1100 real(RKG) , intent(inout) , contiguous :: xmin(:)
1101 real(RKG) , intent(inout) , optional :: fmin
1102 real(RKG) , intent(in) , optional :: tol
1103 integer(IK) , intent(in) , optional :: niter
1104 logical(LK) :: failed
1105 end function
1106#endif
1107
1108#if RK4_ENABLED
1109 recursive module function isFailedMinPowell_RK4(getFunc, xmin, fmin, tol, niter) result(failed)
1110#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1111 !DEC$ ATTRIBUTES DLLEXPORT :: isFailedMinPowell_RK4
1112#endif
1113 use pm_kind, only: RKG => RK4
1114 procedure(real(RKG)) :: getFunc
1115 real(RKG) , intent(inout) , contiguous :: xmin(:)
1116 real(RKG) , intent(inout) , optional :: fmin
1117 real(RKG) , intent(in) , optional :: tol
1118 integer(IK) , intent(in) , optional :: niter
1119 logical(LK) :: failed
1120 end function
1121#endif
1122
1123#if RK3_ENABLED
1124 recursive module function isFailedMinPowell_RK3(getFunc, xmin, fmin, tol, niter) result(failed)
1125#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1126 !DEC$ ATTRIBUTES DLLEXPORT :: isFailedMinPowell_RK3
1127#endif
1128 use pm_kind, only: RKG => RK3
1129 procedure(real(RKG)) :: getFunc
1130 real(RKG) , intent(inout) , contiguous :: xmin(:)
1131 real(RKG) , intent(inout) , optional :: fmin
1132 real(RKG) , intent(in) , optional :: tol
1133 integer(IK) , intent(in) , optional :: niter
1134 logical(LK) :: failed
1135 end function
1136#endif
1137
1138#if RK2_ENABLED
1139 recursive module function isFailedMinPowell_RK2(getFunc, xmin, fmin, tol, niter) result(failed)
1140#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1141 !DEC$ ATTRIBUTES DLLEXPORT :: isFailedMinPowell_RK2
1142#endif
1143 use pm_kind, only: RKG => RK2
1144 procedure(real(RKG)) :: getFunc
1145 real(RKG) , intent(inout) , contiguous :: xmin(:)
1146 real(RKG) , intent(inout) , optional :: fmin
1147 real(RKG) , intent(in) , optional :: tol
1148 integer(IK) , intent(in) , optional :: niter
1149 logical(LK) :: failed
1150 end function
1151#endif
1152
1153#if RK1_ENABLED
1154 recursive module function isFailedMinPowell_RK1(getFunc, xmin, fmin, tol, niter) result(failed)
1155#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1156 !DEC$ ATTRIBUTES DLLEXPORT :: isFailedMinPowell_RK1
1157#endif
1158 use pm_kind, only: RKG => RK1
1159 procedure(real(RKG)) :: getFunc
1160 real(RKG) , intent(inout) , contiguous :: xmin(:)
1161 real(RKG) , intent(inout) , optional :: fmin
1162 real(RKG) , intent(in) , optional :: tol
1163 integer(IK) , intent(in) , optional :: niter
1164 logical(LK) :: failed
1165 end function
1166#endif
1167
1168 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1169
1170 end interface
1171
1172!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1173
1258
1259 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1260
1261#if RK5_ENABLED
1262 recursive module subroutine setMinPowell_RK5(getFunc, xmin, fmin, dirset, tol, niter)
1263#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1264 !DEC$ ATTRIBUTES DLLEXPORT :: setMinPowell_RK5
1265#endif
1266 use pm_kind, only: RKG => RK5
1267 procedure(real(RKG)) :: getFunc
1268 real(RKG) , intent(inout) , contiguous :: xmin(:), dirset(:,:)
1269 real(RKG) , intent(inout) :: fmin
1270 real(RKG) , intent(in) :: tol
1271 integer(IK) , intent(inout) :: niter
1272 end subroutine
1273#endif
1274
1275#if RK4_ENABLED
1276 recursive module subroutine setMinPowell_RK4(getFunc, xmin, fmin, dirset, tol, niter)
1277#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1278 !DEC$ ATTRIBUTES DLLEXPORT :: setMinPowell_RK4
1279#endif
1280 use pm_kind, only: RKG => RK4
1281 procedure(real(RKG)) :: getFunc
1282 real(RKG) , intent(inout) , contiguous :: xmin(:), dirset(:,:)
1283 real(RKG) , intent(inout) :: fmin
1284 real(RKG) , intent(in) :: tol
1285 integer(IK) , intent(inout) :: niter
1286 end subroutine
1287#endif
1288
1289#if RK3_ENABLED
1290 recursive module subroutine setMinPowell_RK3(getFunc, xmin, fmin, dirset, tol, niter)
1291#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1292 !DEC$ ATTRIBUTES DLLEXPORT :: setMinPowell_RK3
1293#endif
1294 use pm_kind, only: RKG => RK3
1295 procedure(real(RKG)) :: getFunc
1296 real(RKG) , intent(inout) , contiguous :: xmin(:), dirset(:,:)
1297 real(RKG) , intent(inout) :: fmin
1298 real(RKG) , intent(in) :: tol
1299 integer(IK) , intent(inout) :: niter
1300 end subroutine
1301#endif
1302
1303#if RK2_ENABLED
1304 recursive module subroutine setMinPowell_RK2(getFunc, xmin, fmin, dirset, tol, niter)
1305#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1306 !DEC$ ATTRIBUTES DLLEXPORT :: setMinPowell_RK2
1307#endif
1308 use pm_kind, only: RKG => RK2
1309 procedure(real(RKG)) :: getFunc
1310 real(RKG) , intent(inout) , contiguous :: xmin(:), dirset(:,:)
1311 real(RKG) , intent(inout) :: fmin
1312 real(RKG) , intent(in) :: tol
1313 integer(IK) , intent(inout) :: niter
1314 end subroutine
1315#endif
1316
1317#if RK1_ENABLED
1318 recursive module subroutine setMinPowell_RK1(getFunc, xmin, fmin, dirset, tol, niter)
1319#if __INTEL_COMPILER && DLL_ENABLED && (_WIN32 || _WIN64)
1320 !DEC$ ATTRIBUTES DLLEXPORT :: setMinPowell_RK1
1321#endif
1322 use pm_kind, only: RKG => RK1
1323 procedure(real(RKG)) :: getFunc
1324 real(RKG) , intent(inout) , contiguous :: xmin(:), dirset(:,:)
1325 real(RKG) , intent(inout) :: fmin
1326 real(RKG) , intent(in) :: tol
1327 integer(IK) , intent(inout) :: niter
1328 end subroutine
1329#endif
1330
1331 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1332
1333 end interface
1334
1335!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1336
1337end module pm_optimization ! LCOV_EXCL_LINE
Generate and return the minimum value and the corresponding abscissa xmin of the input 1-dimensional ...
Generate and return .true. if and only if a concave quadratic curve can fit the specified input tripl...
Generate and return .true. if and only if a convex quadratic curve can fit the specified input triple...
Generate and return .true. if and only if the algorithm fails to find the minimum value and the corre...
Refine an initial input interval such that the final returned interval is guaranteed to contain the m...
Refine an initial input interval such that the final returned interval is guaranteed to contain the m...
Compute and return the minimum value and the corresponding abscissa xmin of the input 1-dimensional f...
Compute and return the minimum value and the corresponding abscissa xmin(1:ndim) of the input arbitra...
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter RK5
Definition: pm_kind.F90:478
integer, parameter RK4
Definition: pm_kind.F90:489
integer, parameter RK2
Definition: pm_kind.F90:511
integer, parameter RK3
Definition: pm_kind.F90:500
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RK1
Definition: pm_kind.F90:522
This module contains procedures, generic interfaces, and types for numerical optimizations of mathema...
character(*, SK), parameter MODULE_NAME