adjmut Subroutine

private subroutine adjmut(me, oldph, fitns, ifit)

Dynamical adjustment of mutation rate:

  • imut=2 or imut=5 : adjustment based on fitness differential between best and median individuals
  • imut=3 or imut=6 : adjustment based on metric distance between best and median individuals

Type Bound

pikaia_class

Arguments

Type IntentOptional Attributes Name
class(pikaia_class), intent(inout) :: me
real(kind=wp), intent(in), dimension(me%n,me%np) :: oldph
real(kind=wp), intent(in), dimension(me%np) :: fitns
integer, intent(in), dimension(me%np) :: ifit

Called by

proc~~adjmut~~CalledByGraph proc~adjmut pikaia_module::pikaia_class%adjmut proc~pikaia pikaia_module::pikaia_class%pikaia proc~pikaia->proc~adjmut proc~solve_with_pikaia pikaia_module::pikaia_class%solve_with_pikaia proc~solve_with_pikaia->proc~pikaia

Source Code

    subroutine adjmut(me,oldph,fitns,ifit)

    implicit none

    class(pikaia_class),intent(inout)         :: me
    integer,dimension(me%np),intent(in)       :: ifit
    real(wp),dimension(me%n,me%np),intent(in) :: oldph
    real(wp),dimension(me%np),intent(in)      :: fitns

    integer  :: i
    real(wp) :: rdif

    real(wp),parameter :: rdiflo = 0.05_wp
    real(wp),parameter :: rdifhi = 0.25_wp
    real(wp),parameter :: delta  = 1.5_wp

    if (me%imut==2 .or. me%imut==5) then

        !Adjustment based on fitness differential
        rdif = abs(fitns(ifit(me%np)) - &
               fitns(ifit(me%np/2)))/(fitns(ifit(me%np)) + &
               fitns(ifit(me%np/2)))

    else if (me%imut==3 .or. me%imut==6) then

        !Adjustment based on normalized metric distance
        rdif=0.0_wp
        do i=1,me%n
            rdif=rdif+( oldph(i,ifit(me%np))-oldph(i,ifit(me%np/2)) )**2
        end do
        rdif=sqrt( rdif ) / real(me%n,wp)

    end if

    if (rdif<=rdiflo) then
        me%pmut=min(me%pmutmx,me%pmut*delta)
    else if (rdif>=rdifhi) then
        me%pmut=max(me%pmutmn,me%pmut/delta)
    end if

    end subroutine adjmut