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