mardi 14 août 2018

Fortran - Random Number Array generator

I am trying to get a 3D matrix whose elements are 4D-vectors themselves choosen with some characteristics: 1) Their norm is 1 and 2) They are choosen randomly

I made this code in Fortran:

MODULE initial_conf
! Modulo en donde se define la configuración
! inicial con la que se va a trabajar
! Autor : Miguel Nava
IMPLICIT NONE

TYPE spin4d
  REAL, DIMENSION(4):: spin
END TYPE spin4d

CONTAINS

  REAL FUNCTION UniformRand(a,b)
    IMPLICIT NONE
    INTEGER, INTENT(IN)::a,b
    REAL:: number
    CALL RANDOM_SEED()
    CALL RANDOM_NUMBER(number)
    UniformRand=(b-a)*number + a
  END FUNCTION UniformRand

  !There is something strange here...
  !######
  SUBROUTINE RanUniVec(vector)
    IMPLICIT NONE
    REAL,DIMENSION(4)::vector
    REAL::number ! Estoy probando
    INTEGER::i
    LOGICAL::cond=.TRUE.
    DO WHILE (cond)
      DO i=1,4
        vector(i)=UniformRand(-1,1)
      END DO
      IF (NORM2(vector)<1) THEN
        cond=.FALSE.
      END IF
    END DO
    vector=vector/NORM2(vector)
  END SUBROUTINE RanUniVec

  SUBROUTINE InitialConf(Conf,m,n,l)
    IMPLICIT NONE
    TYPE (spin4d) , DIMENSION(m,n,l), INTENT(OUT) :: Conf
    INTEGER, INTENT(IN) :: m,n,l
    INTEGER::i,j,k
    DO k=1,l
      DO j=1,n
        DO i=1,m
          CALL RanUniVec(Conf(i,j,k)%spin)
        END DO
      END DO
    END DO
  END SUBROUTINE InitialConf
END MODULE initial_conf

PROGRAM test_random
USE initial_conf
IMPLICIT NONE

INTEGER,parameter::m=2,n=2,l=2
TYPE (spin4d) , DIMENSION(m,n,l) :: Conf
CALL InitialConf(Conf,m,n,l)
print 100,Conf
100 format(1x,"(",4F7.3,")")
end program test_random

For the moment the way I print the information is not relevant. What is relevant is the result of the last test, for example:

 ( -0.303  0.564 -0.466  0.610)
 (    NaN    NaN    NaN    NaN)
 (    NaN    NaN    NaN    NaN)
 (    NaN    NaN    NaN    NaN)
 (    Inf    NaN    Inf    NaN)
 (    NaN    NaN    NaN    NaN)
 (    NaN    NaN    NaN    NaN)
 (    NaN    NaN    NaN    NaN)

The first vector is what I was looking for, but the others of course are not. I still don't know where is the problem but I think it has to do with RanUniVec subroutine. What do you think?

And of course, I am new in Fortran.




Aucun commentaire:

Enregistrer un commentaire