The ParaMonte Documentation Website
Current view: top level - kernel - SpecDRAM_mod.f90 (source / functions) Hit Total Coverage
Test: ParaMonte 1.5.1 :: MPI Parallel Kernel - Code Coverage Report Lines: 84 84 100.0 %
Date: 2021-01-08 13:07:16 Functions: 6 6 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       2             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       3             : !!!!
       4             : !!!!   MIT License
       5             : !!!!
       6             : !!!!   ParaMonte: plain powerful parallel Monte Carlo library.
       7             : !!!!
       8             : !!!!   Copyright (C) 2012-present, The Computational Data Science Lab
       9             : !!!!
      10             : !!!!   This file is part of the ParaMonte library.
      11             : !!!!
      12             : !!!!   Permission is hereby granted, free of charge, to any person obtaining a 
      13             : !!!!   copy of this software and associated documentation files (the "Software"), 
      14             : !!!!   to deal in the Software without restriction, including without limitation 
      15             : !!!!   the rights to use, copy, modify, merge, publish, distribute, sublicense, 
      16             : !!!!   and/or sell copies of the Software, and to permit persons to whom the 
      17             : !!!!   Software is furnished to do so, subject to the following conditions:
      18             : !!!!
      19             : !!!!   The above copyright notice and this permission notice shall be 
      20             : !!!!   included in all copies or substantial portions of the Software.
      21             : !!!!
      22             : !!!!   THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, 
      23             : !!!!   EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF 
      24             : !!!!   MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. 
      25             : !!!!   IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, 
      26             : !!!!   DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR 
      27             : !!!!   OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE 
      28             : !!!!   OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
      29             : !!!!
      30             : !!!!   ACKNOWLEDGMENT
      31             : !!!!
      32             : !!!!   ParaMonte is an honor-ware and its currency is acknowledgment and citations.
      33             : !!!!   As per the ParaMonte library license agreement terms, if you use any parts of 
      34             : !!!!   this library for any purposes, kindly acknowledge the use of ParaMonte in your 
      35             : !!!!   work (education/research/industry/development/...) by citing the ParaMonte 
      36             : !!!!   library as described on this page:
      37             : !!!!
      38             : !!!!       https://github.com/cdslaborg/paramonte/blob/main/ACKNOWLEDGMENT.md
      39             : !!!!
      40             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      41             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      42             : 
      43             : !> \brief
      44             : !> This module contains the classes and procedures for setting up the attributes of samplers of class [ParaDRAM_type](@ref paradram_mod::paradram_type).
      45             : !> For more information, see the description of this attributes in the body of the corresponding modules.
      46             : !> \author Amir Shahmoradi
      47             : 
      48             : module SpecDRAM_mod
      49             : 
      50             :     ! ParaDRAM Spec variable types
      51             :     use SpecDRAM_AdaptiveUpdateCount_mod                , only: AdaptiveUpdateCount_type
      52             :     use SpecDRAM_AdaptiveUpdatePeriod_mod               , only: AdaptiveUpdatePeriod_type
      53             :     use SpecDRAM_GreedyAdaptationCount_mod              , only: GreedyAdaptationCount_type
      54             :     use SpecDRAM_DelayedRejectionCount_mod              , only: DelayedRejectionCount_type
      55             :     use SpecDRAM_BurninAdaptationMeasure_mod            , only: BurninAdaptationMeasure_type
      56             :     use SpecDRAM_DelayedRejectionScaleFactorVec_mod     , only: DelayedRejectionScaleFactorVec_type
      57             : 
      58             :     ! ParaDRAM namelist variables
      59             :     use SpecDRAM_AdaptiveUpdateCount_mod                , only: adaptiveUpdateCount
      60             :     use SpecDRAM_AdaptiveUpdatePeriod_mod               , only: adaptiveUpdatePeriod
      61             :     use SpecDRAM_GreedyAdaptationCount_mod              , only: greedyAdaptationCount
      62             :     use SpecDRAM_DelayedRejectionCount_mod              , only: delayedRejectionCount
      63             :     use SpecDRAM_BurninAdaptationMeasure_mod            , only: burninAdaptationMeasure
      64             :     use SpecDRAM_DelayedRejectionScaleFactorVec_mod     , only: delayedRejectionScaleFactorVec
      65             : 
      66             :     implicit none
      67             : 
      68             :     type                                                :: SpecDRAM_type
      69             :         type(AdaptiveUpdateCount_type)                  :: AdaptiveUpdateCount
      70             :         type(AdaptiveUpdatePeriod_type)                 :: AdaptiveUpdatePeriod
      71             :         type(GreedyAdaptationCount_type)                :: GreedyAdaptationCount
      72             :         type(DelayedRejectionCount_type)                :: DelayedRejectionCount
      73             :         type(BurninAdaptationMeasure_type)              :: BurninAdaptationMeasure
      74             :         type(DelayedRejectionScaleFactorVec_type)       :: DelayedRejectionScaleFactorVec
      75             :     contains
      76             :         procedure, pass                                 :: nullifyNameListVar
      77             :         procedure, pass                                 :: setFromInputFile
      78             :         procedure, pass                                 :: setFromInputArgs
      79             :         procedure, pass                                 :: reportValues
      80             :         procedure, pass                                 :: checkForSanity
      81             :     end type SpecDRAM_type
      82             : 
      83             :     interface SpecDRAM_type
      84             :         module procedure                                :: constructSpecDRAM
      85             :     end interface SpecDRAM_type
      86             : 
      87             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      88             : 
      89             : contains
      90             : 
      91             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      92             : 
      93        1047 :     function constructSpecDRAM  ( nd &
      94             :                                 , methodName &
      95             :                                 !, chainSizeDef
      96             :                                 ) result(SpecDRAM)
      97             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
      98             :         !DEC$ ATTRIBUTES DLLEXPORT :: constructSpecDRAM
      99             : #endif
     100             :         use Constants_mod, only: IK
     101             :         implicit none
     102             :         integer(IK), intent(in)     :: nd
     103             :         character(*), intent(in)    :: methodName
     104             :        !integer(IK), intent(in)     :: chainSizeDef
     105             :         type(SpecDRAM_type)         :: SpecDRAM
     106        1047 :         SpecDRAM%AdaptiveUpdatePeriod                   = AdaptiveUpdatePeriod_type             (nd,methodName)
     107             :         ! ATTN: AdaptiveUpdateCount has to be constructed after AdaptiveUpdatePeriod. It depends on it.
     108        1047 :         SpecDRAM%AdaptiveUpdateCount                    = AdaptiveUpdateCount_type              (methodName) ! ,chainSizeDef,SpecDRAM%AdaptiveUpdatePeriod%def)
     109        1047 :         SpecDRAM%GreedyAdaptationCount                  = GreedyAdaptationCount_type            (methodName)
     110        1047 :         SpecDRAM%DelayedRejectionCount                  = DelayedRejectionCount_type            (methodName)
     111        1047 :         SpecDRAM%BurninAdaptationMeasure                = BurninAdaptationMeasure_type          (methodName)
     112        1047 :         SpecDRAM%DelayedRejectionScaleFactorVec         = DelayedRejectionScaleFactorVec_type   (nd,methodName)
     113        1047 :     end function constructSpecDRAM
     114             : 
     115             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     116             : 
     117        1047 :     subroutine nullifyNameListVar( SpecDRAM , nd )
     118             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
     119             :         !DEC$ ATTRIBUTES DLLEXPORT :: nullifyNameListVar
     120             : #endif
     121        1047 :         use Constants_mod, only: IK
     122             :         implicit none
     123             :         class(SpecDRAM_type), intent(in)    :: SpecDRAM
     124             :         integer(IK), intent(in)             :: nd
     125             :         ! nullify SpecDRAM variables to be read form the input namelist file
     126        1047 :         call SpecDRAM%AdaptiveUpdateCount           %nullifyNameListVar()
     127        1047 :         call SpecDRAM%AdaptiveUpdatePeriod          %nullifyNameListVar()
     128        1047 :         call SpecDRAM%GreedyAdaptationCount         %nullifyNameListVar()
     129        1047 :         call SpecDRAM%DelayedRejectionCount         %nullifyNameListVar()
     130        1047 :         call SpecDRAM%BurninAdaptationMeasure       %nullifyNameListVar()
     131        1047 :         call SpecDRAM%DelayedRejectionScaleFactorVec%nullifyNameListVar()
     132        1047 :     end subroutine nullifyNameListVar
     133             : 
     134             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     135             : 
     136        2094 :     subroutine setFromInputFile( SpecDRAM, Err )
     137             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
     138             :         !DEC$ ATTRIBUTES DLLEXPORT :: setFromInputFile
     139             : #endif
     140        1047 :         use Constants_mod, only: IK, RK
     141             :         use Err_mod, only: Err_type
     142             :         implicit none
     143             :         class(SpecDRAM_type), intent(inout) :: SpecDRAM
     144             :         type(Err_type), intent(out)         :: Err
     145        1047 :         Err%occurred = .false.
     146        1047 :         Err%msg = ""
     147        1047 :         call SpecDRAM%AdaptiveUpdateCount           %set(adaptiveUpdateCount)
     148        1047 :         call SpecDRAM%AdaptiveUpdatePeriod          %set(adaptiveUpdatePeriod)
     149        1047 :         call SpecDRAM%GreedyAdaptationCount         %set(greedyAdaptationCount)
     150        1047 :         call SpecDRAM%DelayedRejectionCount         %set(delayedRejectionCount)
     151        1047 :         call SpecDRAM%BurninAdaptationMeasure       %set(burninAdaptationMeasure)
     152        1047 :         call SpecDRAM%DelayedRejectionScaleFactorVec%set(delayedRejectionCount, delayedRejectionScaleFactorVec)
     153        1047 :     end subroutine setFromInputFile
     154             : 
     155             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     156             : 
     157        1041 :     subroutine setFromInputArgs ( SpecDRAM &
     158             :                                 ! input arguments to the specific ParaDRAM routine
     159             :                                 , adaptiveUpdateCount               &
     160             :                                 , adaptiveUpdatePeriod              &
     161             :                                 , greedyAdaptationCount             &
     162             :                                 , delayedRejectionCount             &
     163             :                                 , burninAdaptationMeasure           &
     164        1041 :                                 , delayedRejectionScaleFactorVec    &
     165             :                                 )
     166             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
     167             :         !DEC$ ATTRIBUTES DLLEXPORT :: setFromInputArgs
     168             : #endif
     169             : 
     170        1047 :         use Constants_mod, only: IK, RK
     171             :         implicit none
     172             :         class(SpecDRAM_type), intent(inout) :: SpecDRAM
     173             : 
     174             :         ! ParaDRAM variables
     175             :         integer(IK) , intent(in), optional  :: adaptiveUpdateCount
     176             :         integer(IK) , intent(in), optional  :: adaptiveUpdatePeriod
     177             :         integer(IK) , intent(in), optional  :: greedyAdaptationCount
     178             :         integer(IK) , intent(in), optional  :: delayedRejectionCount
     179             :         real(RK)    , intent(in), optional  :: burninAdaptationMeasure
     180             :         real(RK)    , intent(in), optional  :: delayedRejectionScaleFactorVec(:)
     181             : 
     182        1041 :         if (present(adaptiveUpdateCount))               call SpecDRAM%AdaptiveUpdateCount           %set(adaptiveUpdateCount)
     183        1041 :         if (present(adaptiveUpdatePeriod))              call SpecDRAM%AdaptiveUpdatePeriod          %set(adaptiveUpdatePeriod)
     184        1041 :         if (present(greedyAdaptationCount))             call SpecDRAM%GreedyAdaptationCount         %set(greedyAdaptationCount)
     185        1041 :         if (present(delayedRejectionCount))             call SpecDRAM%DelayedRejectionCount         %set(delayedRejectionCount)
     186        1041 :         if (present(burninAdaptationMeasure))           call SpecDRAM%BurninAdaptationMeasure       %set(burninAdaptationMeasure)
     187        1041 :         call SpecDRAM%delayedRejectionScaleFactorVec%set(SpecDRAM%DelayedRejectionCount%val, delayedRejectionScaleFactorVec)
     188             : 
     189        1041 :     end subroutine setFromInputArgs
     190             : 
     191             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     192             : 
     193         969 :     subroutine reportValues ( SpecDRAM              &
     194             :                             , prefix                &
     195             :                             , outputUnit            &
     196             :                             , isLeaderImage         &
     197             :                             , splashModeRequested   &
     198             :                             )
     199             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
     200             :         !DEC$ ATTRIBUTES DLLEXPORT :: reportValues
     201             : #endif
     202        1041 :         use Decoration_mod, only: GENERIC_OUTPUT_FORMAT
     203             :         use Decoration_mod, only: GENERIC_TABBED_FORMAT
     204             :         use Constants_mod, only: IK, RK, UNDEFINED
     205             :         use Err_mod, only: note
     206             :         implicit none
     207             :         class(SpecDRAM_type), intent(in)    :: SpecDRAM
     208             :         character(*), intent(in)            :: prefix
     209             :         integer(IK), intent(in)             :: outputUnit
     210             :         logical, intent(in)                 :: isLeaderImage, splashModeRequested
     211             :         integer(IK)                         :: i
     212             : 
     213         969 :         if (isLeaderImage) then
     214             : 
     215         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     216         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "adaptiveUpdatePeriod"
     217         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     218         359 :             write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%AdaptiveUpdatePeriod%val
     219         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%AdaptiveUpdatePeriod%desc )
     220             : 
     221             : 
     222         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     223         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "adaptiveUpdateCount"
     224         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     225         359 :             write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%AdaptiveUpdateCount%val
     226         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%AdaptiveUpdateCount%desc )
     227             : 
     228             : 
     229         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     230         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "greedyAdaptationCount"
     231         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     232         359 :             write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%GreedyAdaptationCount%val
     233         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%GreedyAdaptationCount%desc )
     234             : 
     235             : 
     236         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     237         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "burninAdaptationMeasure"
     238         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     239         359 :             write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%BurninAdaptationMeasure%val
     240         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%BurninAdaptationMeasure%desc )
     241             : 
     242             : 
     243         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     244         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "delayedRejectionCount"
     245         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     246         359 :             write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%DelayedRejectionCount%val
     247         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%DelayedRejectionCount%desc )
     248             : 
     249             : 
     250         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     251         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT) "delayedRejectionScaleFactorVec"
     252         359 :             write(outputUnit,GENERIC_OUTPUT_FORMAT)
     253         359 :             if ( size(SpecDRAM%DelayedRejectionScaleFactorVec%Val) == 0 ) then
     254         341 :                 write(outputUnit,GENERIC_TABBED_FORMAT) UNDEFINED
     255             :             else
     256        2068 :                 do i = 1, size(SpecDRAM%DelayedRejectionScaleFactorVec%Val)
     257        2068 :                     write(outputUnit,GENERIC_TABBED_FORMAT) SpecDRAM%DelayedRejectionScaleFactorVec%Val(i)
     258             :                 end do
     259             :             end if
     260         359 :             if (splashModeRequested) call note( prefix = prefix, outputUnit = outputUnit, newline = "\n", msg = SpecDRAM%DelayedRejectionScaleFactorVec%desc )
     261             : 
     262             :         end if
     263             : 
     264        1938 :     end subroutine reportValues
     265             : 
     266             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     267             : 
     268        1035 :     subroutine checkForSanity(SpecDRAM,Err,methodName,nd)
     269             : #if INTEL_COMPILER_ENABLED && defined DLL_ENABLED && (OS_IS_WINDOWS || defined OS_IS_DARWIN)
     270             :         !DEC$ ATTRIBUTES DLLEXPORT :: checkForSanity
     271             : #endif
     272         969 :         use Constants_mod, only: IK, RK
     273             :         use Err_mod, only: Err_type
     274             :         implicit none
     275             :         class(SpecDRAM_type), intent(inout) :: SpecDRAM
     276             :         integer(IK), intent(in)             :: nd
     277             :         character(*), intent(in)            :: methodName
     278             :         type(Err_type), intent(inout)       :: Err
     279        1035 :         call SpecDRAM%AdaptiveUpdateCount%checkForSanity            (Err,methodName)
     280        1035 :         call SpecDRAM%AdaptiveUpdatePeriod%checkForSanity           (Err,methodName)
     281        1035 :         call SpecDRAM%GreedyAdaptationCount%checkForSanity          (Err,methodName)
     282        1035 :         call SpecDRAM%DelayedRejectionCount%checkForSanity          (Err,methodName)
     283        1035 :         call SpecDRAM%BurninAdaptationMeasure%checkForSanity        (Err,methodName)
     284        1035 :         call SpecDRAM%DelayedRejectionScaleFactorVec%checkForSanity (Err,methodName,SpecDRAM%DelayedRejectionCount%val)
     285        2070 :     end subroutine checkForSanity
     286             : 
     287             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     288             : 
     289             : end module SpecDRAM_mod ! LCOV_EXCL_LINE

ParaMonte: Plain Powerful Parallel Monte Carlo Library 
The Computational Data Science Lab
© Copyright 2012 - 2021