select_parents Subroutine

private subroutine select_parents(me, jfit, imom, idad)

Selects two parents from the population, using roulette wheel algorithm with the relative fitnesses of the phenotypes as the "hit" probabilities.

Reference

  • Davis 1991, chap. 1.

History

  • Jacob Williams : 3/10/2015 : rewrote this routine to return both parents, and also protect against the loop exiting without selecting a parent.

Type Bound

pikaia_class

Arguments

Type IntentOptional Attributes Name
class(pikaia_class), intent(inout) :: me
integer, intent(in), dimension(me%np) :: jfit
integer, intent(out) :: imom
integer, intent(out) :: idad

Calls

proc~~select_parents~~CallsGraph proc~select_parents pikaia_module::pikaia_class%select_parents proc~urand pikaia_module::pikaia_class%urand proc~select_parents->proc~urand genrand64_real1 genrand64_real1 proc~urand->genrand64_real1

Called by

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

Source Code

    subroutine select_parents(me,jfit,imom,idad)

    implicit none

    class(pikaia_class),intent(inout)   :: me
    integer,dimension(me%np),intent(in) :: jfit
    integer,intent(out)                 :: imom
    integer,intent(out)                 :: idad

    integer  :: np1,i,j
    real(wp) :: dice,rtfit
    integer,dimension(2) :: parents

    !initialize:
    np1 = me%np+1
    parents = -99

    !get two (unequal) parents:
    do j=1,2
        main: do
            dice = me%urand()*me%np*np1
            rtfit = 0.0_wp
            do i=1,me%np
                rtfit = rtfit+np1+me%fdif*(np1-2*jfit(i))
                if (rtfit>=dice) then
                    parents(j) = i
                    if (parents(1)/=parents(2)) exit main
                end if
            end do
        end do main
    end do

    imom = parents(1)
    idad = parents(2)

    end subroutine select_parents