penalisation_s.f90 Source File


Contents

Source Code


Source Code

submodule (penalisation_m) penalisation_s
    implicit none

contains

    module subroutine set_parameters_penalisation(self, dphi, filename, &
                pen_method_in, &
                penfuns_type_in, hermite_order_in, &
                charfun_parwidth_in, charfun_radlimwidth_in, &
                dirindfun_parwidth_in, dphi_max_in, charfun_at_target_in, &
                max_step_size_in, rho_min_in)
        class(penalisation_t), intent(inout) :: self
        real(GP), intent(in) :: dphi
        character(*), intent(in), optional :: filename
        integer, intent(in), optional :: pen_method_in
        integer, intent(in), optional :: penfuns_type_in
        integer, intent(in), optional :: hermite_order_in
        real(GP), intent(in), optional :: charfun_parwidth_in
        real(GP), intent(in), optional :: charfun_radlimwidth_in
        real(GP), intent(in), optional :: dirindfun_parwidth_in
        real(GP), intent(in), optional :: dphi_max_in
        real(GP), intent(in), optional :: charfun_at_target_in
        real(GP), intent(in), optional  :: max_step_size_in
        real(GP), intent(in), optional :: rho_min_in

        integer :: io_error
        character(len=256) :: io_errmsg 
        
        character(len=150) :: penfuns_typec, pen_methodc
        integer :: hermite_order
        real(GP) :: charfun_parwidth, charfun_radlimwidth, dirindfun_parwidth, &
                    dphi_max
        real(GP) :: dphi_max_sigmoid
        real(GP) :: charfun_at_target, max_step_size, dphi_max_default, &
                    rho_min
    
        namelist / params_penalisation / &
            pen_methodc, &
            penfuns_typec, &
            hermite_order, &
            charfun_parwidth, &
            charfun_radlimwidth, &
            dirindfun_parwidth, &
            dphi_max, &
            charfun_at_target, &
            max_step_size, &
            rho_min
       
        self%dphi = dphi
       
        if ( (.not.present(filename))       .and. &
           (present(pen_method_in))         .and. &
           (present(penfuns_type_in))       .and. & 
           (present(hermite_order_in))      .and. &
           (present(charfun_parwidth_in))   .and. &
           (present(charfun_radlimwidth_in)).and. &
           (present(dirindfun_parwidth_in)) .and. &
           (present(dphi_max_in))           .and. &
           (present(charfun_at_target_in))  .and. &
           (present(max_step_size_in))      .and. &
           (present(rho_min_in)))then  

            ! input via explicit setting       
            self%pen_method             = pen_method_in
            self%penfuns_type           = penfuns_type_in
            self%hermite_order          = hermite_order_in
            self%charfun_parwidth       = charfun_parwidth_in
            self%charfun_radlimwidth    = charfun_radlimwidth_in
            self%dirindfun_parwidth     = dirindfun_parwidth_in
            self%dphi_max               = dphi_max_in
            self%charfun_at_target      = charfun_at_target_in
            self%max_step_size          = max_step_size_in
            self%rho_min                = rho_min_in
            
        elseif ( (present(filename))               .and. &
           (.not.present(pen_method_in))           .and. & 
           (.not.present(penfuns_type_in))         .and. & 
           (.not.present(hermite_order_in))        .and. &
           (.not.present(charfun_parwidth_in))     .and. &
           (.not.present(charfun_radlimwidth_in))  .and. &
           (.not.present(dirindfun_parwidth_in))   .and. &
           (.not.present(dphi_max_in))             .and. &  
           (.not.present(charfun_at_target_in))    .and. &
           (.not.present(max_step_size_in))        .and. &
           (.not.present(rho_min_in))) then  
            ! Input via file 

            ! default parameters
            dphi_max_default = 4.0_GP*self%dphi

            pen_methodc         = 'VIA_TRACE'
            penfuns_typec       = 'HERMITE'
            hermite_order       = 3
            charfun_parwidth    = 0.0_GP
            charfun_radlimwidth = 0.0_GP
            dirindfun_parwidth  = 0
            dphi_max            = dphi_max_default
            charfun_at_target   = 0.5_GP
            max_step_size       = dphi_max_default
            rho_min             = 0.0

            open(unit = 20, file = filename, status = 'old', action = 'read',&
                iostat = io_error, iomsg = io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, GRILLIX_ERR_NAMELIST, __LINE__, __FILE__)
            endif

            read(20, nml = params_penalisation, iostat = io_error, iomsg = io_errmsg)
            if (io_error /= 0) then
                call handle_error(io_errmsg, GRILLIX_ERR_NAMELIST, __LINE__, __FILE__)
            endif          
               
            close(20)
            
            select case(pen_methodc)
                case('VIA_TRACE')
                    self%pen_method = PEN_METHOD_VIA_TRACE
                case('VIA_IN_VESSEL')
                    self%pen_method = PEN_METHOD_VIA_IN_VESSEL
                case('VIA_SUBTRACE')
                    self%pen_method = PEN_METHOD_VIA_SUBTRACE
                case('VIA_STABLE_TRACE')
                    self%pen_method = PEN_METHOD_VIA_STABLE_TRACE
                case('VIA_RHO')
                    self%pen_method = PEN_METHOD_VIA_RHO
                case('NONE')
                    self%pen_method = PEN_METHOD_NONE
                case default
                    call handle_error('Pen_method not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__, &
                                      error_info_t(pen_methodc))
            end select  
            
            select case(penfuns_typec)
                case('TANH')
                    self%penfuns_type = PENFUNS_TYPE_TANH
                case('HERMITE')
                    self%penfuns_type = PENFUNS_TYPE_HERMITE
                case default
                    call handle_error('Penfuns_type not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__, &
                                      error_info_t(penfuns_typec))
            end select
            
            self%hermite_order       = hermite_order
            self%charfun_parwidth    = charfun_parwidth
            self%charfun_radlimwidth = charfun_radlimwidth
            self%dirindfun_parwidth  = dirindfun_parwidth
            self%dphi_max            = dphi_max
            self%charfun_at_target   = charfun_at_target
            self%rho_min             = rho_min

            if(abs(max_step_size - dphi_max_default) > GP_EPS) then
                self%max_step_size = max_step_size
            else
                self%max_step_size = dphi_max
            endif

            
        else
            call handle_error('Input to subroutine not consistent', &
                              GRILLIX_ERR_OTHER, __LINE__, __FILE__)
        endif
        
        ! Determine maximum angle to be traced based on sigmoid step width
        select case(self%penfuns_type)
            case(PENFUNS_TYPE_TANH)
                dphi_max_sigmoid = 10.0_GP * max(self%charfun_parwidth, self%dirindfun_parwidth) ! after this distance from the target the tanh step would have dropped to ~1E-9 
            case(PENFUNS_TYPE_HERMITE)
                dphi_max_sigmoid = 1.2_GP * max(self%charfun_parwidth, self%dirindfun_parwidth)  ! after this distance from the target the hermite step dropped to 0
            case default
                call handle_error('Penfuns_type not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__, &
                                  error_info_t('penfuns_type:', [self%penfuns_type]))
        end select  
        self%dphi_max = max(self%dphi_max, dphi_max_sigmoid)
        
        ! Check that charfun_at_target input is within valid range
        if ( self%charfun_at_target <= 100*GP_EPS .or. self%charfun_at_target >= 1.0_GP - 100*GP_EPS) then
            call handle_error('Invalid input value, must be within open interval ]0, 1[', &
                                GRILLIX_ERR_OTHER, __LINE__, __FILE__, &
                                error_info_t('charfun_at_target', &
                                r_info=[self%charfun_at_target]))
        endif

    end subroutine 

    module subroutine display_penalisation(self)
        !! Displays information on penalisation
        class(penalisation_t), intent(inout) :: self
        !! Instance of the type

        character(len=:), allocatable :: penfuns_typec, pen_methodc

        select case(self%pen_method)
            case(PEN_METHOD_VIA_TRACE)
                pen_methodc = 'VIA_TRACE'
            case(PEN_METHOD_VIA_IN_VESSEL)
                pen_methodc = 'VIA_IN_VESSEL'
            case(PEN_METHOD_VIA_SUBTRACE)
                pen_methodc = 'VIA_SUBTRACE'
            case(PEN_METHOD_VIA_STABLE_TRACE)
                pen_methodc = 'VIA_STABLE_TRACE'
            case(PEN_METHOD_NONE)
                pen_methodc = 'NONE'
            case default
                pen_methodc = 'not set' 
        end select 

        select case(self%penfuns_type)
            case(PENFUNS_TYPE_TANH)
                penfuns_typec = 'TANH'
            case(PENFUNS_TYPE_HERMITE)
                penfuns_typec = 'HERMITE'
            case default
                penfuns_typec = 'not set' 
        end select 

        if (.not. is_rank_info_writer) then
            return
        endif

        write(get_stdout(),*)'Information on penalisation ---------------------------'
        write(get_stdout(),101)pen_methodc, self%pen_method, &
                penfuns_typec, self%penfuns_type,           &
                self%hermite_order,                         &
                self%charfun_parwidth,                      &
                self%charfun_radlimwidth,                   &
                self%dirindfun_parwidth,                    &
                self%dphi,                                  &
                self%dphi_max,                              &
                self%charfun_at_target,                     &
                self%max_step_size,                         &
                self%n_p_inds,                              &
                self%n_pb_inds
                                
101     FORMAT( 3X,'pen_method          = ',A20,' (',I5,')' /, &
                3X,'penfuns_type        = ',A20,' (',I5,')' /, &
                3X,'hermite_order       = ',I5              /, &
                3X,'charfun_parwidth    = ',ES14.6E3        /, &
                3X,'charfun_radlimwidth = ',ES14.6E3        /, &
                3X,'dirindfun_parwidth  = ',ES14.6E3        /, &
                3X,'dphi                = ',ES14.6E3        /, &
                3X,'dphi_max            = ',ES14.6E3        /, &
                3X,'charfun_at_target   = ',ES14.6E3        /, &
                3X,'max_step_size       = ',ES14.6E3        /, &
                3X,'n_p_inds            = ',I9              /, &
                3X,'n_pb_inds           = ',I9)
                write(get_stdout(),*)'-------------------------------------------------------'

    end subroutine

end submodule