我编写了一个Fortran代码，该代码可以计算给定列表的ith-
排列{1,2,3,...,n}，而无需计算所有其他列表，这是n!我需要的以便找到TSP的ith-路径（旅行商问题）。

当n!大，代码给了我一些错误，我测试的第i个置换发现的是不是确切值。对于n =
10，根本没有问题，但是对于n=20，代码崩溃或找到错误的值。我认为这是由于Fortran进行大数运算（大数加法）而导致的错误。

我使用Visual Fortran
Ultimate2013。在附件中，您找到了我用于目标的子例程。WeightAdjMatRete是网络的每对结之间的距离矩阵。

    ! Fattoriale
RECURSIVE FUNCTION factorial(n) RESULT(n_factorial)
IMPLICIT NONE
REAL, INTENT(IN) :: n
REAL :: n_factorial
IF(n>0) THEN
    n_factorial=n*factorial(n-1)
ELSE
    n_factorial=1.
ENDIF
ENDFUNCTION factorial

! ith-permutazione di una lista
SUBROUTINE ith_permutazione(lista_iniziale,n,i,ith_permutation)
IMPLICIT NONE
INTEGER :: k,n
REAL :: j,f
REAL, INTENT(IN) :: i
INTEGER, DIMENSION(1:n), INTENT(IN) :: lista_iniziale
INTEGER, DIMENSION(1:n) :: lista_lavoro
INTEGER, DIMENSION(1:n), INTENT(OUT) :: ith_permutation
lista_lavoro=lista_iniziale
j=i
DO k=1,n
    f=factorial(REAL(n-k))
    ith_permutation(k)=lista_lavoro(FLOOR(j/f)+1)
    lista_lavoro=PACK(lista_lavoro,MASK=lista_lavoro/=ith_permutation(k))
    j=MOD(j,f)
ENDDO
ENDSUBROUTINE ith_permutazione

! Funzione modulo, adattata
PURE FUNCTION mood(k,modulo) RESULT(ris)
IMPLICIT NONE
INTEGER, INTENT(IN) :: k,modulo
INTEGER :: ris
IF(MOD(k,modulo)/=0) THEN
    ris=MOD(k,modulo)
ELSE
    ris=modulo
ENDIF
ENDFUNCTION mood

! Funzione quoziente, adattata
PURE FUNCTION quoziente(a,p) RESULT(ris)
IMPLICIT NONE
INTEGER, INTENT(IN) :: a,p
INTEGER :: ris
IF(MOD(a,p)/=0) THEN
    ris=(a/p)+1
ELSE
    ris=a/p
ENDIF
ENDFUNCTION quoziente

! Vettori contenenti tutti i payoff percepiti dagli agenti allo state vector attuale e quelli ad ogni sua singola permutazione
SUBROUTINE tuttipayoff(n,m,nodi,nodi_rete,sigma,bvector,MatVecSomma,VecPos,lista_iniziale,ith_permutation,lunghezze_percorso,WeightAdjMatRete,array_perceived_payoff_old,array_perceived_payoff_neg)
IMPLICIT NONE
INTEGER, INTENT(IN) :: n,m,nodi,nodi_rete
INTEGER, DIMENSION(1:nodi), INTENT(IN) :: sigma
INTEGER, DIMENSION(1:nodi), INTENT(OUT) :: bvector
REAL, DIMENSION(1:m,1:n), INTENT(OUT) :: MatVecSomma
REAL, DIMENSION(1:m), INTENT(OUT) :: VecPos
INTEGER, DIMENSION(1:nodi_rete), INTENT(IN) :: lista_iniziale
INTEGER, DIMENSION(1:nodi_rete), INTENT(OUT) :: ith_permutation
REAL, DIMENSION(1:nodi_rete), INTENT(OUT) :: lunghezze_percorso
REAL, DIMENSION(1:nodi_rete,1:nodi_rete), INTENT(IN) :: WeightAdjMatRete
REAL, DIMENSION(1:nodi), INTENT(OUT) :: array_perceived_payoff_old,array_perceived_payoff_neg
INTEGER :: i,j,k
bvector=sigma
FORALL(i=1:nodi,bvector(i)==-1)
    bvector(i)=0
ENDFORALL
FORALL(i=1:m,j=1:n)
    MatVecSomma(i,j)=bvector(m*(j-1)+i)*(2.**REAL(n-j))
ENDFORALL
FORALL(i=1:m)
    VecPos(i)=1.+SUM(MatVecSomma(i,:))
ENDFORALL
DO k=1,nodi
    IF(VecPos(mood(k,m))<=factorial(REAL(nodi_rete))) THEN
        CALL ith_permutazione(lista_iniziale,nodi_rete,VecPos(mood(k,m))-1.,ith_permutation)
        FORALL(i=1:(nodi_rete-1))
            lunghezze_percorso(i)=WeightAdjMatRete(ith_permutation(i),ith_permutation(i+1))
        ENDFORALL
        lunghezze_percorso(nodi_rete)=WeightAdjMatRete(ith_permutation(nodi_rete),ith_permutation(1))
        array_perceived_payoff_old(k)=(1./SUM(lunghezze_percorso))
    ELSE
        array_perceived_payoff_old(k)=0.
    ENDIF
    IF(VecPos(mood(k,m))-SIGN(1,sigma(m*(quoziente(k,m)-1)+mood(k,m)))*2**(n-quoziente(k,m))<=factorial(REAL(nodi_rete))) THEN
        CALL ith_permutazione(lista_iniziale,nodi_rete,VecPos(mood(k,m))-SIGN(1,sigma(m*(quoziente(k,m)-1)+mood(k,m)))*2**(n-quoziente(k,m))-1.,ith_permutation)
        FORALL(i=1:(nodi_rete-1))
            lunghezze_percorso(i)=WeightAdjMatRete(ith_permutation(i),ith_permutation(i+1))
        ENDFORALL
        lunghezze_percorso(nodi_rete)=WeightAdjMatRete(ith_permutation(nodi_rete),ith_permutation(1))
        array_perceived_payoff_neg(k)=(1./SUM(lunghezze_percorso))
    ELSE
        array_perceived_payoff_neg(k)=0.
    ENDIF
ENDDO
ENDSUBROUTINE tuttipayoff