stdrep Subroutine

private subroutine stdrep(me, ph, fits, oldph, fitns, ifit, jfit, nnew)

Steady-state reproduction: insert offspring pair into population only if they are fit enough (replace-random if irep=2 or replace-worst if irep=3).

Type Bound

pikaia_class

Arguments

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

Calls

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

Called by

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

Source Code

    subroutine stdrep(me,ph,fits,oldph,fitns,ifit,jfit,nnew)

    implicit none

    class(pikaia_class),intent(inout)             :: me
    real(wp),dimension(me%n,2),intent(in)         :: ph
    real(wp),dimension(2),intent(in)              :: fits
    real(wp),dimension(me%n,me%np),intent(inout)  :: oldph
    real(wp),dimension(me%np),intent(inout)       :: fitns
    integer,dimension(me%np),intent(inout)        :: ifit
    integer,dimension(me%np),intent(inout)        :: jfit
    integer,intent(out)                           :: nnew

    integer  :: i,j,k,i1,if1
    real(wp) :: fit

    nnew = 0
    main_loop : do j=1,2

        !1. get offspring fitness
        fit = fits(j)

        !2. if fit enough, insert in population
        do i=me%np,1,-1

            if (fit>fitns(ifit(i))) then

                !make sure the phenotype is not already in the population
                if (i<me%np) then
                    if (all(oldph(1:me%n,ifit(i+1))==ph(1:me%n,j))) cycle main_loop
                end if

                !offspring is fit enough for insertion, and is unique

                !(i) insert phenotype at appropriate place in population
                if (me%irep==3) then
                    i1=1
                else if (me%ielite==0 .or. i==me%np) then
                    i1=int(me%urand()*me%np)+1
                else
                    i1=int(me%urand()*(me%np-1))+1
                end if
                if1 = ifit(i1)
                fitns(if1)=fit
                do k=1,me%n
                    oldph(k,if1)=ph(k,j)
                end do

                !(ii) shift and update ranking arrays
                if (i<i1) then

                    !shift up
                    jfit(if1)=me%np-i
                    do k=i1-1,i+1,-1
                        jfit(ifit(k))=jfit(ifit(k))-1
                        ifit(k+1)=ifit(k)
                    end do
                    ifit(i+1)=if1

                else

                    !shift down
                    jfit(if1)=me%np-i+1
                    do k=i1+1,i
                        jfit(ifit(k))=jfit(ifit(k))+1
                        ifit(k-1)=ifit(k)
                    end do
                    ifit(i)=if1

                end if

                nnew = nnew+1
                cycle main_loop

            end if

        end do

    end do main_loop

    end subroutine stdrep