mutate Subroutine

private subroutine mutate(me, gn)

Introduces random mutation in a genotype. Mutations occur at rate pmut at all gene loci.

Input

  • imut=1 Uniform mutation, constant rate
  • imut=2 Uniform mutation, variable rate based on fitness
  • imut=3 Uniform mutation, variable rate based on distance
  • imut=4 Uniform or creep mutation, constant rate
  • imut=5 Uniform or creep mutation, variable rate based on fitness
  • imut=6 Uniform or creep mutation, variable rate based on distance

Type Bound

pikaia_class

Arguments

Type IntentOptional Attributes Name
class(pikaia_class), intent(inout) :: me
integer, intent(inout), dimension(me%n*me%nd) :: gn

Calls

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

Called by

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

Source Code

    subroutine mutate(me,gn)

    implicit none

    class(pikaia_class),intent(inout)           :: me
    integer,dimension(me%n*me%nd),intent(inout) :: gn

    integer :: i,j,k,l,ist,inc,loc
    logical :: fix

    !Decide which type of mutation is to occur
    if (me%imut>=4 .and. me%urand()<=0.5_wp) then

        !CREEP MUTATION OPERATOR
        !Subject each locus to random +/- 1 increment at the rate pmut
        do i=1,me%n
            do j=1,me%nd

                if (me%urand()<me%pmut) then

                    !Construct integer
                    loc=(i-1)*me%nd+j
                    inc=nint( me%urand() )*2-1
                    ist=(i-1)*me%nd+1
                    gn(loc)=gn(loc)+inc

                    !This is where we carry over the one (up to two digits)
                    !first take care of decrement below 0 case
                    if (inc<0 .and. gn(loc)<0) then
                        if (j==1) then
                            gn(loc)=0
                        else
                            fix = .true.
                            do k=loc,ist+1,-1
                                gn(k)=9
                                gn(k-1)=gn(k-1)-1
                                if ( gn(k-1)>=0 ) then
                                    fix = .false.
                                    exit
                                end if
                            end do
                            if (fix) then
                                !we popped under 0.00000 lower bound; fix it up
                                if ( gn(ist)<0) then
                                    do l=ist,loc
                                        gn(l)=0
                                    end do
                                end if
                            end if
                        end if
                    end if

                    if (inc>0 .and. gn(loc)>9) then
                        if (j==1) then
                            gn(loc)=9
                        else
                            fix = .true.
                            do k=loc,ist+1,-1
                                gn(k)=0
                                gn(k-1)=gn(k-1)+1
                                if ( gn(k-1)<=9 ) then
                                    fix = .false.
                                    exit
                                end if
                            end do
                            if (fix) then
                                !we popped over 9.99999 upper bound; fix it up
                                if ( gn(ist)>9 ) then
                                    do l=ist,loc
                                        gn(l)=9
                                    end do
                                end if
                            end if
                        end if
                    end if

                end if

            end do
        end do

    else

        !UNIFORM MUTATION OPERATOR
        !Subject each locus to random mutation at the rate pmut
        do i=1,me%n*me%nd
            if (me%urand()<me%pmut) then
                gn(i)=int(me%urand()*10.0_wp)
            end if
        end do

    end if

    end subroutine mutate