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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Subroutine: DES_READ_RESTART                                        !
4     !  Purpose : Reads either single restart file or multiple restart      !
5     !  fles (based on bdist_io) flag.                                      !
6     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
7           SUBROUTINE READ_RES0_DES
8     
9           use cdist
10           use compar
11           use des_allocate
12           use des_bc
13           use des_rxns
14           use des_thermo
15           use desmpi
16           use discretelement
17           use error_manager
18           use machine
19           use mfix_pic, only: MPPIC, DES_STAT_WT
20           use mpi_utility
21           use param1
22           use read_res1_des
23           use run
24     
25           implicit none
26     
27           INTEGER :: LC1, LC2
28           INTEGER :: lDIMN, lNEXT_REC
29     
30           DOUBLE PRECISION :: VERSION
31     
32           lDIMN = merge(2,3,NO_K)
33     
34           CALL INIT_READ_RES_DES(trim(RUN_NAME), VERSION, lNEXT_REC)
35     
36           CALL READ_RES_DES(lNEXT_REC, VTP_FINDEX)
37           CALL READ_RES_DES(lNEXT_REC, TECPLOT_FINDEX)
38           CALL READ_RES_DES(lNEXT_REC, DTSOLID)
39     
40     ! Position data is read and used to setup pARRAY reads.
41           CALL READ_PAR_POS(lNEXT_REC)
42     
43           CALL READ_RES_pARRAY(lNEXT_REC, iGLOBAL_ID)
44     
45           DO LC1 = 2, 4
46              CALL READ_RES_pARRAY(lNEXT_REC, PEA(:,LC1))
47           ENDDO
48     
49           DO LC1 = 1, lDIMN
50              CALL READ_RES_pARRAY(lNEXT_REC, DES_VEL_NEW(LC1,:))
51           ENDDO
52     
53           DO LC1 = 1, merge(1,3,NO_K)
54              CALL READ_RES_pARRAY(lNEXT_REC, OMEGA_NEW(LC1,:))
55           ENDDO
56     
57           CALL READ_RES_pARRAY(lNEXT_REC, DES_RADIUS)
58           CALL READ_RES_pARRAY(lNEXT_REC, RO_SOL)
59     
60           IF(MPPIC) CALL READ_RES_pARRAY(lNEXT_REC, DES_STAT_WT)
61           IF(ENERGY_EQ) CALL READ_RES_pARRAY(lNEXT_REC, DES_T_s_NEW)
62     
63           IF(ANY_SPECIES_EQ) THEN
64              DO LC1=1, DIMENSION_N_S
65                 CALL READ_RES_pARRAY(lNEXT_REC, DES_X_s(:,LC1))
66              ENDDO
67           ENDIF
68     
69           IF(VERSION >= 1.1) THEN
70              CALL READ_RES_DES(lNEXT_REC, DES_USR_VAR_SIZE)
71              DO LC1=1, DES_USR_VAR_SIZE
72                 CALL READ_RES_pARRAY(lNEXT_REC, DES_USR_VAR(LC1,:))
73              ENDDO
74           ENDIF
75     
76     ! RES2 does not need the collision of BC information.
77           IF(RUN_TYPE == 'RESTART_2') RETURN
78     
79     ! Collision/neighbor data is read and used to setup cARRAY reads.
80           CALL READ_PAR_COL(lNEXT_REC)
81     
82           CALL READ_RES_cARRAY(lNEXT_REC, PV_PAIR(:))
83           DO LC1=1, lDIMN
84              CALL READ_RES_cARRAY(lNEXT_REC, PFN_PAIR(LC1,:))
85              CALL READ_RES_cARRAY(lNEXT_REC, PFT_PAIR(LC1,:))
86           ENDDO
87     
88     ! Save the number of BCMI's read from input file, then read the
89     ! value from the restart file.
90           CALL READ_RES_DES(lNEXT_REC, DEM_BCMI)
91     
92     ! Allocation of MIs is done here to ignore changes to the mfix.dat
93     ! file during RES1.
94           IF(DEM_BCMI > 0) CALL ALLOCATE_DEM_MI
95     
96     ! Only save the number of mass inflows for RESTART_1. This allows
97     ! for mass inflows to be added/removed with RESTART_2.
98     ! Todo: Prune entering/exiting flagged particles for RESTART_2.
99           DO LC1=1, DEM_BCMI
100              CALL READ_RES_DES(lNEXT_REC, DEM_MI_TIME(LC1))
101              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%VACANCY)
102              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%OCCUPANTS)
103              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%WINDOW)
104              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%OFFSET)
105              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%L)
106     
107              LC2 = DEM_MI(LC1)%OCCUPANTS
108     
109              allocate(DEM_MI(LC1)%W(LC2))
110              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%W(:))
111              allocate(DEM_MI(LC1)%H(LC2))
112              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%H(:))
113              allocate(DEM_MI(LC1)%P(LC2))
114              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%P(:))
115              allocate(DEM_MI(LC1)%Q(LC2))
116              CALL READ_RES_DES(lNEXT_REC, DEM_MI(LC1)%Q(:))
117           ENDDO
118     
119           CALL FINL_READ_RES_DES
120     
121     
122           WRITE(ERR_MSG,"('DES restart file read at Time = ',g12.5)") TIME
123           CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
124     
125           RETURN
126           END SUBROUTINE READ_RES0_DES
127