File: /nfs/home/0/users/jenkins/mfix.git/model/des/read_particle_input.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     ! Subroutine: READ_PAR_INPUT                                           !
3     !                                                                      !
4     ! Purpose: Read the particle input and broadcasts the particle data to !
5     ! respective processors.                                               !
6     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
7           SUBROUTINE READ_PAR_INPUT
8     
9     !-----------------------------------------------
10     ! Modules
11     !-----------------------------------------------
12           USE discretelement
13           use funits
14           use compar
15           use desmpi
16           use cdist
17     
18           use mpi_init_des, only: des_scatter_particle
19           use mpi_utility
20           use geometry, only: NO_K
21     
22           use error_manager
23     
24           implicit none
25     !-----------------------------------------------
26     ! Local variables
27     !-----------------------------------------------
28     ! indices
29           integer :: i,j,k
30     ! index of particle
31           INTEGER :: lcurpar
32     ! local unit
33           INTEGER, PARAMETER :: lunit=10
34     ! local filename
35           character(30) lfilename
36     ! IO Status:
37           INTEGER :: IOS
38     ! Flag to indicate if file exists.
39           LOGICAL :: lEXISTS
40     ! Read dimension: 2D vs 3D data
41           integer :: RDMN
42     !-----------------------------------------------
43     
44     
45           CALL INIT_ERR_MSG("READ_PAR_INPUT")
46     
47     
48           IOS = 0
49           RDMN = merge(2,3,NO_K)
50     
51     ! Setup the file name based on distributed or serial IO.
52           IF(bDIST_IO) THEN
53              lFILENAME = ''
54              WRITE(lFILENAME,'("particle_input_",I4.4,".dat")') myPE
55           ELSE
56              lFILENAME= "particle_input.dat"
57           ENDIF
58     
59     ! Check the the file exists and open it.
60           IF(bDIST_IO .OR. myPE == PE_IO) THEN
61              INQUIRE(FILE=lFILENAME, EXIST=lEXISTS)
62              IF(.NOT.LEXISTS) THEN
63                 WRITE(ERR_MSG, 1100)
64                 CALL FLUSH_ERR_MSG
65                 IOS = 1
66              ELSE
67                 OPEN(UNIT=lUNIT, FILE=lFILENAME, FORM="FORMATTED")
68              ENDIF
69           ENDIF
70     
71     ! Collect the error message and quit.
72           CALL GLOBAL_ALL_SUM(IOS)
73           IF(IOS /= 0) CALL MFIX_EXIT(myPE)
74     
75      1100 FORMAT('Error 1100: FATAL - DEM particle input file not found!')
76     
77     
78     ! Read the file
79     !----------------------------------------------------------------->>>
80     ! In distributed IO the first line of the file will be number of
81     ! particles in that processor
82           IF (bdist_io) then
83              read(lunit,*) pip
84              DO lcurpar = 1,pip
85                 pea(lcurpar,1) = .true.
86                 read (lunit,*) (des_pos_new(k,lcurpar),k=1,RDMN),&
87                    des_radius(lcurpar), ro_sol(lcurpar),&
88                    (des_vel_new(k,lcurpar),k=1,RDMN)
89              ENDDO
90     
91     
92     ! Serial IO (not bDIST_IO)
93           ELSE
94     !----------------------------------------------------------------->>>
95     
96     ! Read into temporary variable and scatter
97              IF (myPE .eq. PE_IO) THEN
98     
99     ! Allocate and initialize temporary variables.
100                 ALLOCATE (dpar_pos(particles,3)); dpar_pos=0.0
101                 ALLOCATE (dpar_vel(particles,3)); dpar_vel=0.0
102                 ALLOCATE (dpar_rad(particles));   dpar_rad=0.0
103                 ALLOCATE (dpar_den(particles));   dpar_den = 0.0
104     ! Loop through the input file.
105                 DO lcurpar = 1, particles
106                    read (lunit,*,IOSTAT=IOS)                               &
107                    (dpar_pos(lcurpar,k),k=1,RDMN),dpar_rad(lcurpar),       &
108                    dpar_den(lcurpar),(dpar_vel(lcurpar,k),k=1,RDMN)
109     
110     ! Report read errors.
111                    IF(IOS > 0) THEN
112                       WRITE(ERR_MSG,1200)
113                       CALL FLUSH_ERR_MSG
114                       EXIT
115      1200 FORMAT('Error 1200: Error reported when reading particle input ',&
116              'file.',/'A common error is 2D input for 3D cases.')
117     
118     ! Report End-of-File errors.
119                    ELSEIF(IOS < 0) THEN
120                       WRITE(ERR_MSG,1201) &
121                          trim(iVal(lcurpar)), trim(iVal(Particles))
122                       CALL FLUSH_ERR_MSG
123                       EXIT
124      1201 FORMAT('Error 1201: Error reported when reading particle input ',&
125              'file.',/'End-of-File found for particle ',A,' and ',A,1X,    &
126              'entries are expected.')
127     
128                    ENDIF
129     
130                 ENDDO
131     
132              ENDIF
133     
134              CALL GLOBAL_ALL_SUM(IOS)
135              IF(IOS /= 0) CALL MFIX_EXIT(myPE)
136     
137              CALL DES_SCATTER_PARTICLE
138     
139              IF(myPE == PE_IO) &
140                 deallocate (dpar_pos,dpar_vel,dpar_rad,dpar_den)
141     
142           ENDIF   ! end if/else bdist_io
143     !-----------------------------------------------------------------<<<
144     
145           IF(bDIST_IO .OR. myPE == PE_IO) CLOSE(lUNIT)
146     
147     
148           CALL FINL_ERR_MSG()
149     
150           RETURN
151     
152       999 IF(dmp_log)write(unit_log,"(/1X,70('*')//,A,/10X,A,/1X,70('*')/)")&
153              ' From: read_par_input -',&
154              ' particle_input.dat file is missing.  Terminating run.'
155           CALL MFIX_EXIT(myPE)
156     
157           END SUBROUTINE READ_PAR_INPUT
158     
159