File: RELATIVE:/../../../mfix.git/model/write_spx1.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: WRITE_SPX1                                             C
4     !  Purpose: write out the time-dependent restart records (REAL)        C
5     !                                                                      C
6     !  Author: P. Nicoletti                               Date: 13-DEC-91  C
7     !  Reviewer: P. Nicoletti, W. Rogers, M. Syamlal      Date: 24-JAN-92  C
8     !                                                                      C
9     !  Revision Number:                                                    C
10     !  Purpose:                                                            C
11     !  Author:                                            Date: dd-mmm-yy  C
12     !  Reviewer:                                          Date: dd-mmm-yy  C
13     !                                                                      C
14     !  Literature/Document References:                                     C
15     !                                                                      C
16     !  Variables referenced: TIME, NSTEP, EP_g, P_g, P_star, U_g           C
17     !                        V_g, W_g, U_s, V_s, W_s, ROP_s, T_g, T_s      C
18     !                        IJKMAX2, MMAX                                 C
19     !  Variables modified: None                                            C
20     !                                                                      C
21     !  Local variables:  LC, N, NEXT_REC, NUM_REC                          C
22     !                                                                      C
23     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
24     !
25           SUBROUTINE WRITE_SPX1(L, unit_add)
26     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
27     !...Switches: -xf
28     !
29     !-----------------------------------------------
30     !   M o d u l e s
31     !-----------------------------------------------
32           USE cdist
33           USE compar
34           USE cutcell
35           USE fldvar
36           USE funits
37           USE geometry
38           USE machine
39           USE mpi_utility
40           USE output
41           USE param
42           USE param1
43           USE physprop
44           USE run
45           USE rxns
46           USE scalars
47           use discretelement, only: PRINT_DES_DATA
48           use discretelement, only: DISCRETE_ELEMENT
49           use discretelement, only: DES_CONTINUUM_COUPLED
50           use discretelement, only: PARTICLES, NFACTOR
51     
52     !//       USE tmp_array
53           IMPLICIT NONE
54     !-----------------------------------------------
55     !   G l o b a l   P a r a m e t e r s
56     !-----------------------------------------------
57     !-----------------------------------------------
58     !   D u m m y   A r g u m e n t s
59     !-----------------------------------------------
60     !
61     !             flag whether to write a particular SPx file
62           INTEGER L
63     
64     !              offset for use in post_mfix
65           INTEGER  unit_add
66     !-----------------------------------------------
67     !   L o c a l   P a r a m e t e r s
68     !-----------------------------------------------
69     !-----------------------------------------------
70     !   L o c a l   V a r i a b l e s
71     !-----------------------------------------------
72     !
73     ! local variables
74     !
75     !//
76           double precision, allocatable :: array1(:)     !//
77           double precision, allocatable :: array2(:)     !//
78     
79     !             loop counters
80           INTEGER LC, N
81     !
82     !             Pointer to the next record
83           INTEGER NEXT_REC
84     !
85     !              Number of records written each time step
86           INTEGER  NUM_REC
87     
88           INTEGER  uspx   ! UNIT_SPX + offset from post_mfix
89           CHARACTER(LEN=50), DIMENSION(1) :: LINE   !error message
90           double precision, dimension(:), allocatable :: TMP_VAR
91     
92           allocate(TMP_VAR(DIMENSION_3))
93     
94     !-----------------------------------------------
95           uspx = UNIT_SPX + unit_add
96     
97     !
98           if (myPE .eq.PE_IO) then
99              allocate (array1(ijkmax2))   !//
100              allocate (array2(ijkmax3))   !//
101           else
102              allocate (array1(1))   !//
103              allocate (array2(1))   !//
104           end if
105     !
106     ! ".SP1" FILE         EP_g    [ ROP_g, RO_g  must be calculated ...
107     !                                        not written out ]
108     !
109           SELECT CASE (L)
110           CASE (1)
111              if (myPE.eq.PE_IO.or.bDist_IO) then
112                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
113                 NUM_REC = NEXT_REC
114                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
115                 NEXT_REC = NEXT_REC + 1
116              end if
117     
118     ! Explicitly coupled simulations do not need to rebin particles to
119     ! the fluid grid every time step. However, this implies that the
120     ! fluid cell information and interpolation weights become stale.
121              IF(DISCRETE_ELEMENT .AND. .NOT.DES_CONTINUUM_COUPLED) THEN
122     ! Bin particles to fluid grid.
123                 CALL PARTICLES_IN_CELL
124     ! Calculate interpolation weights
125                 CALL CALC_INTERP_WEIGHTS
126     ! Calculate mean fields (EPg).
127                 CALL COMP_MEAN_FIELDS
128              ENDIF
129     
130     
131              if (bDist_IO) then
132                 IF(RE_INDEXING) THEN
133                    CALL UNSHIFT_DP_ARRAY(EP_g,TMP_VAR)
134                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
135                 ELSE
136                    call OUT_BIN_R(uspx+L,EP_g,size(EP_g),NEXT_REC)
137                 ENDIF
138     !           call OUT_BIN_R(uspx+L,EP_g,size(EP_g),NEXT_REC)
139              else
140                 call gatherWriteSpx (EP_g,array2, array1, uspx+L, NEXT_REC)   !//
141              end if
142              if (myPE .eq. PE_IO.or.bDist_IO) then
143                 NUM_REC = NEXT_REC - NUM_REC
144                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
145                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
146              end if
147     
148     
149     ! remove this redundant call here to write_des_data in case of new
150     ! coupled runs that contain at least a particle and involve some initial
151     ! settling of the system.
152     ! the call made in des_time_march is a better call for capturing the
153     ! initial state of such a des continuum coupled system
154              IF(DISCRETE_ELEMENT.AND.PRINT_DES_DATA .AND. &
155                 .NOT.(TRIM(RUN_TYPE)=='NEW' .AND. PARTICLES /=0 .AND. &
156                       NFACTOR >0 .AND. TIME == ZERO)) THEN
157                    CALL WRITE_DES_DATA
158              ENDIF
159     
160     
161     
162     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
163     !
164     ! ".SP2" FILE         P_g , P_star
165     !
166           CASE (2)
167              if (myPE.eq.PE_IO.or.bDist_IO) then
168                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
169                 NUM_REC = NEXT_REC
170                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
171                 NEXT_REC = NEXT_REC + 1
172              end if
173              if (bDist_IO) then
174                 IF(RE_INDEXING) THEN
175                    CALL UNSHIFT_DP_ARRAY(P_g,TMP_VAR)
176                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
177                    CALL UNSHIFT_DP_ARRAY(P_star,TMP_VAR)
178                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
179                 ELSE
180                   call OUT_BIN_R(uspx+L,P_g,size(P_g),NEXT_REC)
181                   call OUT_BIN_R(uspx+L,P_star,size(P_star),NEXT_REC)
182                 ENDIF
183     !           call OUT_BIN_R(uspx+L,P_g,size(P_g),NEXT_REC)
184     !           call OUT_BIN_R(uspx+L,P_star,size(P_star),NEXT_REC)
185              else
186                call gatherWriteSpx (P_g,array2, array1, uspx+L, NEXT_REC)   !//
187                call gatherWriteSpx (P_star,array2, array1, uspx+L, NEXT_REC)   !//
188              end if
189              if (myPE.eq.PE_IO.or.bDist_IO) then
190                 NUM_REC = NEXT_REC - NUM_REC
191                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
192                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
193              end if
194     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
195     !
196     ! ".SP3" FILE         U_g , V_g , W_g
197     !
198           CASE (3)
199              if (myPE.eq.PE_IO.or.bDist_IO) then
200                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
201                 NUM_REC = NEXT_REC
202                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
203                 NEXT_REC = NEXT_REC + 1
204              end if
205              if (bDist_IO) then
206                 IF(RE_INDEXING) THEN
207                    CALL UNSHIFT_DP_ARRAY(U_g,TMP_VAR)
208                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
209                    CALL UNSHIFT_DP_ARRAY(V_g,TMP_VAR)
210                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
211                    CALL UNSHIFT_DP_ARRAY(W_g,TMP_VAR)
212                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
213                 ELSE
214                    call OUT_BIN_R(uspx+L,U_g,size(U_g),NEXT_REC)
215                    call OUT_BIN_R(uspx+L,V_g,size(V_g),NEXT_REC)
216                    call OUT_BIN_R(uspx+L,W_g,size(W_g),NEXT_REC)
217                 ENDIF
218     !           call OUT_BIN_R(uspx+L,U_g,size(U_g),NEXT_REC)
219     !           call OUT_BIN_R(uspx+L,V_g,size(V_g),NEXT_REC)
220     !           call OUT_BIN_R(uspx+L,W_g,size(W_g),NEXT_REC)
221              else
222                call gatherWriteSpx (U_g,array2, array1, uspx+L, NEXT_REC)   !//
223                call gatherWriteSpx (V_g,array2, array1, uspx+L, NEXT_REC)   !//
224                call gatherWriteSpx (W_g,array2, array1, uspx+L, NEXT_REC)   !//
225              end if
226              if (myPE.eq.PE_IO.or.bDist_IO) then
227                 NUM_REC = NEXT_REC - NUM_REC
228                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
229                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
230              end if
231     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
232     !
233     ! ".SP4" FILE         U_s , V_s , W_s
234     !
235           CASE (4)
236              if (myPE.eq.PE_IO.or.bDist_IO) then
237                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
238                 NUM_REC = NEXT_REC
239                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
240                 NEXT_REC = NEXT_REC + 1
241              end if
242              if (bDist_IO) then
243                 DO LC = 1, MMAX
244                    IF(RE_INDEXING) THEN
245                       CALL UNSHIFT_DP_ARRAY(U_s(:,LC),TMP_VAR)
246                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
247                       CALL UNSHIFT_DP_ARRAY(V_s(:,LC),TMP_VAR)
248                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
249                       CALL UNSHIFT_DP_ARRAY(W_s(:,LC),TMP_VAR)
250                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
251                    ELSE
252                       call OUT_BIN_R(uspx+L,U_s(:,LC),Size(U_s(:,LC)),NEXT_REC)
253                       call OUT_BIN_R(uspx+L,V_s(:,LC),Size(V_s(:,LC)),NEXT_REC)
254                       call OUT_BIN_R(uspx+L,W_s(:,LC),Size(W_s(:,LC)),NEXT_REC)
255                    ENDIF
256                 ENDDO
257     !        DO LC = 1, MMAX
258     !          call OUT_BIN_R(uspx+L,U_s(:,LC),Size(U_s(:,LC)),NEXT_REC)
259     !          call OUT_BIN_R(uspx+L,V_s(:,LC),Size(V_s(:,LC)),NEXT_REC)
260     !          call OUT_BIN_R(uspx+L,W_s(:,LC),Size(W_s(:,LC)),NEXT_REC)
261     !        END DO
262              else
263              DO LC = 1, MMAX
264                 call gatherWriteSpx (U_s(:,LC),array2, array1, uspx+L, NEXT_REC)
265                 call gatherWriteSpx (V_s(:,LC),array2, array1, uspx+L, NEXT_REC)
266                 call gatherWriteSpx (W_s(:,LC),array2, array1, uspx+L, NEXT_REC)
267              END DO
268              end if
269              if (myPE.eq.PE_IO.or.bDist_IO) then
270                 NUM_REC = NEXT_REC - NUM_REC
271                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
272                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
273              end if
274     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
275     !
276     ! ".SP5" FILE         ROP_s
277     !
278           CASE (5)
279              if (myPE.eq.PE_IO.or.bDist_IO) then
280                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
281                 NUM_REC = NEXT_REC
282                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
283                 NEXT_REC = NEXT_REC + 1
284              end if
285              if (bDist_IO) then
286                 DO LC = 1, MMAX
287                    IF(RE_INDEXING) THEN
288                       CALL UNSHIFT_DP_ARRAY(ROP_s(:,LC),TMP_VAR)
289                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
290                    ELSE
291                       call OUT_BIN_R(uspx+L,ROP_s(:,LC),size(ROP_s(:,LC)), NEXT_REC)
292                    ENDIF
293                 ENDDO
294     !          DO LC = 1, MMAX
295     !            call OUT_BIN_R(uspx+L,ROP_s(:,LC),size(ROP_s(:,LC)), NEXT_REC)
296     !         END DO
297              else
298              DO LC = 1, MMAX
299                 call gatherWriteSpx (ROP_s(:,LC),array2, array1, uspx+L, NEXT_REC)
300              END DO
301              end if
302              if (myPE.eq.PE_IO.or.bDist_IO) then
303                 NUM_REC = NEXT_REC - NUM_REC
304                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
305                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
306              end if
307     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
308     !
309     ! ".SP6" FILE         T_g  , T_s
310     !
311           CASE (6)
312              if (myPE.eq.PE_IO.or.bDist_IO) then
313                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
314                 NUM_REC = NEXT_REC
315                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
316                 NEXT_REC = NEXT_REC + 1
317              end if
318              if (bDist_IO) then
319                 IF(RE_INDEXING) THEN
320                    CALL UNSHIFT_DP_ARRAY(T_g,TMP_VAR)
321                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
322                    DO LC = 1, MMAX
323                       CALL UNSHIFT_DP_ARRAY(T_s(:,LC),TMP_VAR)
324                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR), NEXT_REC)
325                    END DO
326                 ELSE
327                    call OUT_BIN_R(uspx+L,T_g,size(T_g), NEXT_REC)
328                    DO LC = 1, MMAX
329                       call OUT_BIN_R(uspx+L,T_s(:,LC),size(T_s(:,LC)), NEXT_REC)
330                    END DO
331                 ENDIF
332     !           call OUT_BIN_R(uspx+L,T_g,size(T_g), NEXT_REC)
333     !          DO LC = 1, MMAX
334     !            call OUT_BIN_R(uspx+L,T_s(:,LC),size(T_s(:,LC)), NEXT_REC)
335     !          END DO
336              else
337              call gatherWriteSpx (T_g,array2, array1, uspx+L, NEXT_REC)   !//
338               DO LC = 1, MMAX
339                 call gatherWriteSpx (T_s(:,LC),array2, array1, uspx+L, NEXT_REC)
340               END DO
341              end if
342              if (myPE.eq.PE_IO.or.bDist_IO) then
343                 NUM_REC = NEXT_REC - NUM_REC
344                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
345                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
346              end if
347     !
348     ! ".SP7" FILE         X_g, X_s
349     !
350           CASE (7)
351              if (myPE.eq.PE_IO.or.bDist_IO) then
352                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
353                 NUM_REC = NEXT_REC
354                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
355                 NEXT_REC = NEXT_REC + 1
356              end if
357              if (bDist_IO) then
358                 IF(RE_INDEXING) THEN
359                    DO N = 1, NMAX(0)
360                       CALL UNSHIFT_DP_ARRAY(X_G(:,N),TMP_VAR)
361                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
362                    END DO
363                    DO LC = 1, MMAX
364                       DO N = 1, NMAX(LC)
365                          CALL UNSHIFT_DP_ARRAY(X_s(:,LC,N),TMP_VAR)
366                          call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR), NEXT_REC)
367                       ENDDO
368                    END DO
369                 ELSE
370                    DO N = 1, NMAX(0)
371                       call OUT_BIN_R(uspx+L,X_G(:,N),size(X_G(:,N)), NEXT_REC)
372                    END DO
373                    DO LC = 1, MMAX
374                       DO N = 1, NMAX(LC)
375                          call OUT_BIN_R(uspx+L,X_s(:,LC,N),size(X_s(:,LC,N)), NEXT_REC)
376                       END DO
377                    END DO
378                 ENDIF
379     
380     !           DO N = 1, NMAX(0)
381     !             call OUT_BIN_R(uspx+L,X_G(:,N),size(X_G(:,N)), NEXT_REC)
382     !           END DO
383     !           DO LC = 1, MMAX
384     !            DO N = 1, NMAX(LC)
385     !               call OUT_BIN_R(uspx+L,X_s(:,LC,N),size(X_s(:,LC,N)), NEXT_REC)
386     !            END DO
387     !           END DO
388     
389              else
390               DO N = 1, NMAX(0)
391                 call gatherWriteSpx (X_G(:,N),array2, array1, uspx+L, NEXT_REC)
392               END DO
393               DO LC = 1, MMAX
394                 DO N = 1, NMAX(LC)
395                    call gatherWriteSpx (X_s(:,LC,N),array2, array1, uspx+L, NEXT_REC)
396                 END DO
397               END DO
398              end if
399              if (myPE.eq.PE_IO.or.bDist_IO) then
400                 NUM_REC = NEXT_REC - NUM_REC
401                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
402                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
403              end if
404     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
405     !
406     ! ".SP8" FILE         THETA_m
407     !
408           CASE (8)
409              if (myPE.eq.PE_IO.or.bDist_IO) then
410                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
411                 NUM_REC = NEXT_REC
412                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
413                 NEXT_REC = NEXT_REC + 1
414              end if
415              if (bDist_IO) then
416                 IF(RE_INDEXING) THEN
417                    DO LC = 1, MMAX
418                       CALL UNSHIFT_DP_ARRAY(THETA_m(:,LC),TMP_VAR)
419                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
420                    ENDDO
421                 ELSE
422                    DO LC = 1, MMAX
423                       call OUT_BIN_R(uspx+L,THETA_m(:,LC),size(THETA_m(:,LC)), NEXT_REC)
424                    END DO
425                 ENDIF
426     !          DO LC = 1, MMAX
427     !            call OUT_BIN_R(uspx+L,THETA_m(:,LC),size(THETA_m(:,LC)), NEXT_REC)
428     !         END DO
429              else
430               DO LC = 1, MMAX
431                 call gatherWriteSpx (THETA_m(:,LC),array2, array1, uspx+L, NEXT_REC)
432               END DO
433              end if
434              if (myPE.eq.PE_IO .or. bDist_IO) then
435                 NUM_REC = NEXT_REC - NUM_REC
436                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
437                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
438              end if
439     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
440     !
441     ! ".SP9" FILE         Scalar
442     !
443           CASE (9)
444              if (myPE.eq.PE_IO.or.bDist_IO) then
445                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
446                 NUM_REC = NEXT_REC
447                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
448                 NEXT_REC = NEXT_REC + 1
449              end if
450              if (bDist_IO) then
451                 IF(RE_INDEXING) THEN
452                    DO LC = 1, Nscalar
453                       CALL UNSHIFT_DP_ARRAY(Scalar(:,LC),TMP_VAR)
454                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
455                    ENDDO
456                 ELSE
457                    DO LC = 1, Nscalar
458                      call OUT_BIN_R(uspx+L,Scalar(:,LC),size(Scalar(:,LC)), NEXT_REC)
459                    END DO
460                 ENDIF
461     !          DO LC = 1, Nscalar
462     !            call OUT_BIN_R(uspx+L,Scalar(:,LC),size(Scalar(:,LC)), NEXT_REC)
463     !          END DO
464              else
465               DO LC = 1, Nscalar
466                 call gatherWriteSpx (Scalar(:,LC),array2, array1, uspx+L, NEXT_REC)
467               END DO
468              end if
469              if (myPE.eq.PE_IO.or.bDist_IO) then
470                 NUM_REC = NEXT_REC - NUM_REC
471                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
472                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
473              end if
474     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
475     !
476           CASE (10)  ! Reaction rates
477     
478              if (myPE.eq.PE_IO.or.bDist_IO) then
479                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
480                 NUM_REC = NEXT_REC
481                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
482                 NEXT_REC = NEXT_REC + 1
483              end if
484              if (bDist_IO) then
485                 IF(RE_INDEXING) THEN
486                    DO LC = 1, nRR
487                       CALL UNSHIFT_DP_ARRAY(ReactionRates(:,LC),TMP_VAR)
488                       call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
489                    ENDDO
490                 ELSE
491                    DO LC = 1, nRR
492                       call OUT_BIN_R(uspx+L,ReactionRates(:,LC),size(ReactionRates(:,LC)), NEXT_REC)
493                    END DO
494                 ENDIF
495     !          DO LC = 1, nRR
496     !            call OUT_BIN_R(uspx+L,ReactionRates(:,LC),size(ReactionRates(:,LC)), NEXT_REC)
497     !          END DO
498              else
499               DO LC = 1, nRR
500                 call gatherWriteSpx (ReactionRates(:,LC),array2, array1, uspx+L, NEXT_REC)
501               END DO
502              end if
503              if (myPE.eq.PE_IO.or.bDist_IO) then
504                 NUM_REC = NEXT_REC - NUM_REC
505                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
506                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
507              end if
508     !
509     ! ".SP11" FILE         turbulence
510     !
511           CASE (11)
512              if (myPE.eq.PE_IO.or.bDist_IO) then
513                 READ (uspx + L, REC=3) NEXT_REC, NUM_REC
514                 NUM_REC = NEXT_REC
515                 WRITE (uspx + L, REC=NEXT_REC) REAL(TIME), NSTEP
516                 NEXT_REC = NEXT_REC + 1
517              end if
518              if (K_Epsilon) then
519                 if (bDist_IO) then
520                 IF(RE_INDEXING) THEN
521                    CALL UNSHIFT_DP_ARRAY(K_Turb_G,TMP_VAR)
522                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
523                    CALL UNSHIFT_DP_ARRAY(E_Turb_G,TMP_VAR)
524                    call OUT_BIN_R(uspx+L,TMP_VAR,size(TMP_VAR),NEXT_REC)
525                 ELSE
526                    call OUT_BIN_R(uspx+L,K_Turb_G,size(K_Turb_G), NEXT_REC)
527                    call OUT_BIN_R(uspx+L,E_Turb_G,size(E_Turb_G), NEXT_REC)
528                 ENDIF
529     !             call OUT_BIN_R(uspx+L,K_Turb_G,size(K_Turb_G), NEXT_REC)
530     !             call OUT_BIN_R(uspx+L,E_Turb_G,size(E_Turb_G), NEXT_REC)
531     
532                 else
533                  call gatherWriteSpx (K_Turb_G,array2, array1, uspx+L, NEXT_REC)
534                  call gatherWriteSpx (E_Turb_G,array2, array1, uspx+L, NEXT_REC)
535                 end if
536               end if
537     
538                if (myPE.eq.PE_IO.or.bDist_IO) then
539                 NUM_REC = NEXT_REC - NUM_REC
540                 WRITE (uspx + L, REC=3) NEXT_REC, NUM_REC
541                 if(unit_add == 0) CALL FLUSH_bin (uspx + L)
542                end if
543     !        call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
544     !
545     !
546           CASE DEFAULT
547                 LINE(1) = 'Unknown SPx file index'
548                 CALL WRITE_ERROR ('WRITE_SPX1', LINE, 1)
549                 CALL MFIX_EXIT(myPE)
550           END SELECT
551     
552     !//      call unlock_tmp_array
553     !
554           deallocate (array1)
555           deallocate (array2)
556           deallocate (TMP_VAR)
557     !
558           RETURN
559           END SUBROUTINE WRITE_SPX1
560     
561           subroutine gatherWriteSpx(VAR, array2, array1, uspxL, NEXT_REC)
562             USE geometry
563             USE compar           !//
564             USE mpi_utility      !//d pnicol : for gatherWriteSpx
565             USE sendrecv         !//d pnicol : for gatherWriteSpx
566             USE cutcell
567             USE in_binary_512
568             IMPLICIT NONE
569             integer uspxL, NEXT_REC
570             double precision, dimension(ijkmax2) :: array1
571             double precision, dimension(ijkmax3) :: array2
572             double precision, dimension(DIMENSION_3) :: VAR
573             double precision, dimension(:), allocatable :: TMP_VAR
574     
575             allocate(TMP_VAR(DIMENSION_3))
576     
577     !       call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
578             IF(RE_INDEXING) THEN
579                TMP_VAR = UNDEFINED
580                CALL UNSHIFT_DP_ARRAY(VAR,TMP_VAR)
581                CALL gather (TMP_VAR,array2,root)
582             ELSE
583                CALL gather (VAR,array2,root)
584             ENDIF
585     !        call gather (VAR,array2,root)  !//d pnicol
586     
587     !       call MPI_Barrier(MPI_COMM_WORLD,mpierr)  !//PAR_I/O enforce barrier here
588             if (myPE.eq.PE_IO) then
589                call convert_to_io_dp(array2,array1,ijkmax2)
590                CALL OUT_BIN_R (uspxL, array1, IJKMAX2, NEXT_REC)
591             end if
592     
593             deallocate(TMP_VAR)
594     
595           End subroutine gatherWriteSpx
596     
597     
598     
599           subroutine gatherWriteSpx_netcdf(VAR, arr1, arr2 , arr4d, ncid, varid , &
600                     nx,ny,nz,ijkmax2_use , ijkmax3_use)
601     
602     
603          USE geometry
604          use param, only: dimension_3
605          USE compar           !//
606          USE mpi_utility      !//d pnicol : for gatherWriteSpx
607          USE sendrecv         !//d pnicol : for gatherWriteSpx
608          USE MFIX_netcdf
609          USE in_binary_512
610     
611          IMPLICIT NONE
612     
613          integer          :: ncid , varid , nx,ny,nz , ijkmax2_use , ijkmax3_use
614          integer          :: ii , jj , kk , ijk
615     
616          double precision ::  arr1(ijkmax2_use)
617          double precision ::  arr2(ijkmax3_use)
618          double precision ::  arr4d(nx,ny,nz,1)
619          double precision ::  var(dimension_3)
620     
621          call gather(var,arr2,root)
622          if (myPE .eq. PE_IO) then
623             call convert_to_io_dp(arr2,arr1,ijkmax2_use)
624     
625             ijk = 0
626             do kk = 1,nz
627                do jj = 1,ny
628                   do ii = 1,nx
629                      ijk = ijk + 1
630                      arr4d(ii,jj,kk,1) = arr1(ijk)
631                   end do
632                end do
633             end do
634     
635             call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid,arr4d) )
636     
637          end if
638     
639     
640         End subroutine gatherWriteSpx_netcdf
641     
642     
643     
644           subroutine gatherWriteSpx_netcdf_int(VAR, arr1, arr2 , arr4d, ncid, &
645                     varid , nx,ny,nz,ijkmax2_use , ijkmax3_use)
646     
647     
648          USE geometry
649          use param, only: dimension_3
650          USE compar           !//
651          USE mpi_utility      !//d pnicol : for gatherWriteSpx
652          USE sendrecv         !//d pnicol : for gatherWriteSpx
653          USE MFIX_netcdf
654          USE in_binary_512i
655     
656          IMPLICIT NONE
657     
658          integer          :: ncid , varid , nx,ny,nz , ijkmax2_use , ijkmax3_use
659          integer          :: ii , jj , kk , ijk
660     
661          integer ::  arr1(ijkmax2_use)
662          integer ::  arr2(ijkmax3_use)
663          integer :: arr4d(nx,ny,nz,1)
664          integer ::   var(dimension_3)
665     
666          call gather(var,arr2,root)
667          if (myPE .eq. PE_IO) then
668             call convert_to_io_i(arr2,arr1,ijkmax2_use)
669     
670             ijk = 0
671             do kk = 1,nz
672                do jj = 1,ny
673                   do ii = 1,nx
674                      ijk = ijk + 1
675                      arr4d(ii,jj,kk,1) = arr1(ijk)
676                   end do
677                end do
678             end do
679     
680             call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid,arr4d) )
681     
682          end if
683     
684     
685         End subroutine gatherWriteSpx_netcdf_int
686     
687     
688     
689     
690             subroutine copy_d_to_r(darr,rarr,nx,ny,nz)
691             implicit none
692     
693             integer          :: nx , ny , nz
694             double precision :: darr(*)
695             real             :: rarr(nx,ny,*)
696             integer          :: i , j , k , ijk
697     
698     
699             ijk = 0
700     
701             do i = 1,nx
702                do j = 1,ny
703                   do k = 1,nz
704                      ijk = ijk + 1
705                      rarr(i,j,k) = real(darr(ijk))
706                   end do
707                end do
708             end do
709     
710             return
711             end subroutine copy_d_to_r
712     
713     
714     
715     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
716     !
717     !                          write_mesh_netcdf
718     !
719     !
720             SUBROUTINE write_mesh_netcdf
721     
722             USE param
723             USE param1
724             USE fldvar
725             USE geometry
726             USE physprop
727             USE run
728     !       USE funits
729             USE scalars
730     !       USE output
731             USE rxns
732             USE cdist
733             USE compar
734             USE mpi_utility
735             USE MFIX_netcdf
736     !       USE tmp_array
737     
738             implicit none
739     
740             integer   :: ncid    , x_dimid , y_dimid , z_dimid  , t_dimid
741             integer   :: varid_x , varid_y , varid_z , L , dimids(4)
742             integer   :: varid_flag , coords_dimid , varid_coords , coords
743     
744             character(LEN=80) :: fname
745     
746             double precision, dimension(:) , allocatable :: xloc
747             double precision, dimension(:) , allocatable :: yloc
748             double precision, dimension(:) , allocatable :: zloc
749     
750             integer, dimension(:) , allocatable :: arr1
751             integer, dimension(:) , allocatable :: arr2
752             integer, dimension(:,:,:,:) , allocatable :: arr4d
753     
754     
755             if (.not. MFIX_usingNETCDF()) return
756     
757             if (.not. bGlobalNetcdf) return  ! no netCDF writes asked for
758     
759             if (.not. bFirst_netcdf_write) return
760     
761             if (myPE.eq.PE_IO) then
762                     allocate ( arr1(ijkmax2))
763                     allocate ( arr2(ijkmax3))
764                     allocate ( arr4d(imax2,jmax2,kmax2,1))
765                     allocate ( xloc(imax2) )
766                     allocate ( yloc(jmax2) )
767                     allocate ( zloc(kmax2) )
768             else
769                     allocate ( arr1(1))
770                     allocate ( arr2(1))
771                     allocate ( arr4d(1,1,1,1))
772                     allocate ( xloc(1) )
773                     allocate ( yloc(1) )
774                     allocate ( zloc(1) )
775             end if
776     
777             if (myPE.eq.PE_IO) then
778                     xloc(1) = -dx(1)
779                     do L = 2,imax2
780                             xloc(L) = xloc(L-1) + dx(L)
781                     end do
782     
783                     yloc(1) = -dy(1)
784                     do L = 2,jmax2
785                             yloc(L) = yloc(L-1) + dy(L)
786                     end do
787     
788                     zloc(1) = -dz(1)
789                     do L = 2,kmax2
790                             zloc(L) = zloc(L-1) + dz(L)
791                     end do
792     
793                     fname = trim(run_name) // "_MESH.nc"
794                     call MFIX_check_netcdf( MFIX_nf90_create(fname, NF90_64BIT_OFFSET, ncid) )
795     
796                     call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "x"   , imax2   , x_dimid  ) )
797                     call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "y"   , jmax2   , y_dimid  ) )
798                     call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "z"   , kmax2   , z_dimid  ) )
799                     call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "coordinates"   , 1   , coords_dimid  ) )
800                     call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "t"   , 1       , t_dimid  ) )  ! 4
801     
802                     dimids =  (/ x_dimid , y_dimid, z_dimid , t_dimid/)
803     
804     
805                     ! The dimids array is used to pass the IDs of the dimensions of
806                     ! the variables. Note that in fortran arrays are stored in
807                     ! column-major format.
808                     call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "x"  , NF90_DOUBLE, x_dimid, varid_x  ) )
809                     call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "y"  , NF90_DOUBLE, y_dimid, varid_y  ) )
810                     call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "z"  , NF90_DOUBLE, z_dimid, varid_z  ) )
811                     call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "coordinates"  , NF90_INT, coords_dimid, varid_coords  ) )
812                     call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "flag"  , NF90_INT, dimids, varid_flag  ) )  ! 9
813     
814     
815                     call MFIX_check_netcdf( MFIX_nf90_enddef(ncid) )
816     
817                     coords = 0
818                     if (COORDINATES .eq. 'CYLINDRICAL') coords = 1
819     
820                     call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_coords,coords) )
821                     call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_x,xloc) )
822                     call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_y,yloc) )
823                     call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_z,zloc) )
824             end if
825     
826     
827             ! needs to be called by all processes
828     
829             call gatherWriteSpx_netcdf_int(flag, arr1, arr2 , arr4d, ncid, &
830                    varid_flag , imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
831     
832     
833     
834             if (myPE.eq.PE_IO) then
835                     call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
836             end if
837     
838             deallocate ( arr1  )
839             deallocate ( arr2  )
840             deallocate ( arr4d )
841             deallocate ( xloc  )
842             deallocate ( yloc  )
843             deallocate ( zloc  )
844     
845             return
846             end subroutine write_mesh_netcdf
847     
848             SUBROUTINE write_netcdf(L, unit_add, the_time)
849     
850             USE param
851             USE param1
852             USE fldvar
853             USE geometry
854             USE physprop
855             USE run
856     !       USE funits
857             USE scalars
858     !       USE output
859             USE rxns
860             USE cdist
861             USE compar
862             USE mpi_utility
863             USE MFIX_netcdf
864     !       USE tmp_array
865     
866     
867             implicit none
868     
869             integer :: L , unit_add , I , n , ii
870     
871             integer   :: ncid , x_dimid , y_dimid , z_dimid
872             integer   :: t_dimid
873             integer   :: dimids(4) , varid_epg , varid_pg
874     
875             integer   :: varid_pstar  , varid_ug , varid_vg , varid_wg
876             integer   :: varid_tg , varid_x , varid_y , varid_z , varid_t
877             integer   :: varid_coords , coords_dimid , coords
878     
879             integer   :: varid_us(20) , varid_vs(20) , varid_ws(20)  !! MMAX
880             integer   :: varid_rops(20)  , varid_ts(20) !! mmax
881             integer   :: varid_thetam(20) !! mmax
882     
883             integer   :: varid_xg(20)  ! nmax(0)
884             integer   :: varid_xs(20,20)  ! mmax , MAX(nmax(1:mmax))
885     
886             integer   :: varid_scalar(20)  ! nscalar
887             integer   :: varid_rr(20)      ! nRR
888     
889             integer   :: varid_kturbg , varid_eturbg
890     
891     
892             character(LEN=80) :: fname, var_name
893             character(LEN=9) :: fname_index
894     
895             double precision, dimension(:) , allocatable :: arr1
896             double precision, dimension(:) , allocatable :: arr2
897     
898             double precision, dimension(:,:,:,:) , allocatable :: arr4d
899     
900     
901             double precision, dimension(:) , allocatable :: xloc
902             double precision, dimension(:) , allocatable :: yloc
903             double precision, dimension(:) , allocatable :: zloc
904     
905             double precision :: the_time
906             logical          :: file_exists
907     
908     
909     ! bWrite_netcdf(1)  : EP_g
910     ! bWrite_netcdf(2)  : P_g
911     ! bWrite_netcdf(3)  : P_star
912     ! bWrite_netcdf(4)  : U_g / V_g / W_g
913     ! bWrite_netcdf(5)  : U_s / V_s / W_s
914     ! bWrite_netcdf(6)  : ROP_s
915     ! bWrite_netcdf(7)  : T_g
916     ! bWrite_netcdf(8)  : T_s
917     ! bWrite_netcdf(9)  : X_g
918     ! bWrite_netcdf(10) : X_s
919     ! bWrite_netcdf(11) : Theta_m
920     ! bWrite_netCDF(12) : Scalar
921     ! bWrite_netCDF(13) : ReactionRates
922     ! bWrite_netCDF(14) : k_turb_g , e_turb_g
923     
924             if (.not. MFIX_usingNETCDF()) return
925             if (.not. bGlobalNetcdf) return
926     
927             call write_mesh_netcdf
928     
929             if (myPE.eq.PE_IO .and. .not.bDist_IO) then
930                allocate (arr1(ijkmax2))
931                allocate (arr2(ijkmax3))
932                allocate (arr4d(imax2,jmax2,kmax2,1))
933                allocate ( xloc(imax2) )
934                allocate ( yloc(jmax2) )
935                allocate ( zloc(kmax2) )
936     
937     
938     
939     
940             else
941                allocate (arr1(1))
942                allocate (arr2(1))
943                allocate (arr4d(1,1,1,1))
944                allocate ( xloc(1) )
945                allocate ( yloc(1) )
946                allocate ( zloc(1) )
947            end if
948     
949     !       CALL MPI_BARRIER(MPI_COMM_WORLD,mpierr)
950             if (myPE .ne. PE_IO) goto 1234
951     
952     
953             xloc(1) = -dx(1)
954             do II = 2,imax2
955                xloc(II) = xloc(II-1) + dx(II)
956             end do
957     
958     
959             yloc(1) = -dy(1)
960             do II = 2,jmax2
961                yloc(II) = yloc(II-1) + dy(II)
962             end do
963     
964             zloc(1) = -dz(1)
965             do II = 2,kmax2
966                zloc(II) = zloc(II-1) + dz(II)
967             end do
968     
969     
970     
971     
972             if (bFirst_netcdf_write .and. MFIX_usingNETCDF()) then
973                bFirst_netcdf_write = .false.
974                fname = trim(run_name) // '_netcdf_index.txt'
975                inquire (file=fname,exist=file_exists)
976     
977                if (file_exists) then
978                    open (unit=11,file=fname,status='old')
979                    read (11,*) netCDF_file_index
980                    close (unit=11)
981                else
982                    netCDF_file_index = 0
983                end if
984             end if
985     
986             fname_index = '_xxxxx.nc'
987             write (fname_index(2:6),'(i5.5)') netCDF_file_index
988             fname = trim(run_name)// fname_index
989             call MFIX_check_netcdf( MFIX_nf90_create(fname, NF90_64BIT_OFFSET, ncid) )
990             netCDF_file_index = netCDF_file_index + 1
991     
992             if (MFIX_usingNETCDF()) then
993                     fname = trim(run_name) // '_netcdf_index.txt'
994                     open (unit=11,file=fname,status='unknown')
995                     write (11,*) netCDF_file_index
996                     close (unit=11)
997             end if
998     
999     
1000             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "x"   , imax2   , x_dimid  ) )  ! 1
1001             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "y"   , jmax2   , y_dimid  ) )  ! 2
1002             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "z"   , kmax2   , z_dimid  ) )  ! 3
1003             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "coordinates"   , 1   , coords_dimid  ) )  ! 3
1004             call MFIX_check_netcdf( MFIX_nf90_def_dim(ncid, "t"   , 1       , t_dimid  ) )  ! 4
1005     
1006     
1007             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "x"  , NF90_DOUBLE, x_dimid, varid_x  ) )  ! 5
1008             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "y"  , NF90_DOUBLE, y_dimid, varid_y  ) )  ! 6
1009             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "z"  , NF90_DOUBLE, z_dimid, varid_z  ) )  ! 7
1010             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "coordinates"  , NF90_INT, coords_dimid, varid_coords  ) )  ! 7
1011             call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "t"  , NF90_DOUBLE, t_dimid, varid_t  ) )  ! 8
1012     
1013             dimids =  (/ x_dimid , y_dimid, z_dimid , t_dimid/)
1014     
1015             if (bWrite_netcdf(1)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "EP_g"  , NF90_DOUBLE, dimids, varid_epg  ) )  ! 9
1016             if (bWrite_netcdf(2)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_g"   , NF90_DOUBLE, dimids, varid_pg  ) )   ! 10
1017     
1018             if (bWrite_netcdf(3)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "P_star", NF90_DOUBLE, dimids, varid_pstar) )
1019             if (bWrite_netcdf(4)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "U_g"   , NF90_DOUBLE, dimids, varid_ug   ) )
1020             if (bWrite_netcdf(4)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "V_g"   , NF90_DOUBLE, dimids, varid_vg   ) )
1021             if (bWrite_netcdf(4)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "W_g"   , NF90_DOUBLE, dimids, varid_wg   ) )
1022             if (bWrite_netcdf(7)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, "T_g"   , NF90_DOUBLE, dimids, varid_tg   ) )
1023             do i = 1,1   ! mmax
1024                var_name = 'U_s_xxx'
1025                write (var_name(5:7),'(i3.3)') I
1026                if (bWrite_netcdf(5)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_us(I)) )
1027     
1028                var_name = 'V_s_xxx'
1029                write (var_name(5:7),'(i3.3)') I
1030                if (bWrite_netcdf(5)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_vs(I)) )
1031     
1032                var_name = 'W_s_xxx'
1033                write (var_name(5:7),'(i3.3)') I
1034                if (bWrite_netcdf(5)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ws(I)) )
1035     
1036                var_name = 'ROP_s_xxx'
1037                write (var_name(7:10),'(i3.3)') I
1038                if (bWrite_netcdf(6)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rops(I)) )
1039     
1040                var_name = 'T_s_xxx'
1041                write (var_name(5:7),'(i3.3)') I
1042                if (bWrite_netcdf(8)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_ts(I)) )
1043     
1044                var_name = 'Theta_m_xxx'
1045                write (var_name(9:11),'(i3.3)') I
1046                if (bWrite_netcdf(11)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_thetam(I)) )
1047     
1048                DO N = 1, NMAX(i)
1049                   var_name = 'X_s_xxx_xxx'
1050                   write (var_name(5:7) ,'(i3.3)') I
1051                   write (var_name(9:11),'(i3.3)') n
1052                   if (bWrite_netcdf(10)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xs(I,n)) )
1053                END DO
1054     
1055     
1056             end do
1057     
1058             do i = 1,nmax(0)
1059                var_name = 'X_g_xxx'
1060                write (var_name(5:7),'(i3.3)') I
1061                if (bWrite_netcdf(9)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_xg(I)) )
1062             end do
1063     
1064             do i = 1,nscalar
1065                var_name = 'Scalar_xxx'
1066                write (var_name(8:10),'(i3.3)') I
1067                if (bWrite_netCDF(12)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_scalar(I)) )
1068             end do
1069     
1070             do i = 1,nRR
1071                var_name = 'RRates_xxx'
1072                write (var_name(8:10),'(i3.3)') I
1073                if (bWrite_netCDF(13)) call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, var_name, NF90_DOUBLE, dimids, varid_rr(I)) )
1074             end do
1075     
1076     
1077             if (bWrite_netcdf(14) .and. k_Epsilon) then
1078                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'k_turb_g', NF90_DOUBLE, dimids, varid_kturbg) )
1079                call MFIX_check_netcdf( MFIX_nf90_def_var(ncid, 'e_turb_g', NF90_DOUBLE, dimids, varid_eturbg) )
1080             end if
1081     
1082     
1083     
1084             call MFIX_check_netcdf( MFIX_nf90_enddef(ncid) ) ! 11
1085     
1086      1234   continue
1087                bFirst_netcdf_write = .false.
1088     
1089     !       CALL MPI_BARRIER(MPI_COMM_WORLD,mpierr)
1090     
1091     
1092             if (myPE .eq. PE_IO) then
1093                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_t,the_time) )   ! 12
1094                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_x,xloc) )
1095                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_y,yloc) )
1096                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_z,zloc) )
1097                coords = 0
1098                if (COORDINATES .eq. 'CYLINDRICAL') coords = 1
1099                call MFIX_check_netcdf( MFIX_nf90_put_var(ncid,varid_coords,coords) )
1100             end if
1101     
1102             if (bWrite_netcdf(1)) then
1103     
1104                 call gatherWriteSpx_netcdf(EP_g, arr1, arr2 , arr4d, ncid, varid_epg , &
1105                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1106     
1107             end if
1108     
1109             if (bWrite_netcdf(2)) then
1110     
1111                 call gatherWriteSpx_netcdf(P_g, arr1, arr2 , arr4d, ncid, varid_pg , &
1112                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1113     
1114             end if
1115     
1116     
1117     
1118     
1119             if (bWrite_netcdf(3)) then
1120     
1121                 call gatherWriteSpx_netcdf(P_star, arr1, arr2 , arr4d, ncid, varid_pstar , &
1122                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1123     
1124             end if
1125     
1126             if (bWrite_netcdf(4)) then
1127     
1128                 call gatherWriteSpx_netcdf(U_g, arr1, arr2 , arr4d, ncid, varid_ug , &
1129                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1130     
1131                 call gatherWriteSpx_netcdf(V_g, arr1, arr2 , arr4d, ncid, varid_vg , &
1132                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1133     
1134                 call gatherWriteSpx_netcdf(W_g, arr1, arr2 , arr4d, ncid, varid_wg , &
1135                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1136     
1137             end if
1138     
1139             if (bWrite_netcdf(7)) then
1140     
1141                 call gatherWriteSpx_netcdf(t_g, arr1, arr2 , arr4d, ncid, varid_tg , &
1142                     imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1143     
1144             end if
1145     
1146     
1147             do i = 1,1   ! mmax
1148     
1149                if (bWrite_netcdf(5)) then
1150     
1151                   call gatherWriteSpx_netcdf(u_s(:,i) , arr1, arr2 , arr4d, ncid, varid_us(i) , &
1152                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1153     
1154                   call gatherWriteSpx_netcdf(v_s(:,i) , arr1, arr2 , arr4d, ncid, varid_vs(i) , &
1155                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1156     
1157                   call gatherWriteSpx_netcdf(w_s(:,i) , arr1, arr2 , arr4d, ncid, varid_ws(i) , &
1158                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1159     
1160                end if
1161     
1162                if (bWrite_netcdf(6)) then
1163     
1164                   call gatherWriteSpx_netcdf(ROP_s(:,i) , arr1, arr2 , arr4d, ncid, varid_rops(i) , &
1165                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1166     
1167                end if
1168     
1169                if (bWrite_netcdf(8)) then
1170     
1171                   call gatherWriteSpx_netcdf(T_s(:,i) , arr1, arr2 , arr4d, ncid, varid_ts(i) , &
1172                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1173     
1174                end if
1175     
1176                if (bWrite_netcdf(11)) then
1177     
1178                   call gatherWriteSpx_netcdf(Theta_m(:,i) , arr1, arr2 , arr4d, ncid, varid_thetam(i) , &
1179                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1180     
1181                end if
1182     
1183                if (bWrite_netcdf(10)) then
1184                   do N = 1,nmax(i)
1185     
1186                      call gatherWriteSpx_netcdf(X_s(:,i,N) , arr1, arr2 , arr4d, ncid, varid_xs(i,N) , &
1187                           imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1188     
1189                   end do
1190                end if
1191     
1192     
1193             end do
1194     
1195             if (bWrite_netcdf(9)) then
1196                do i = 1,nmax(0)
1197     
1198                    call gatherWriteSpx_netcdf(X_g(:,i) , arr1, arr2 , arr4d, ncid, varid_xg(i) , &
1199                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1200     
1201                end do
1202             end if
1203     
1204            if (bWrite_netcdf(12)) then
1205               do i = 1,nscalar
1206     
1207                  call gatherWriteSpx_netcdf(Scalar(:,i) , arr1, arr2 , arr4d, ncid, varid_scalar(i) , &
1208                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1209     
1210                end do
1211             end if
1212     
1213     
1214             if (bWrite_netcdf(13)) then
1215                do i = 1,nRR
1216     
1217                  call gatherWriteSpx_netcdf(ReactionRates(:,i) , arr1, arr2 , arr4d, ncid, varid_rr(i) , &
1218                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1219     
1220               end do
1221             end if
1222     
1223             if (bWrite_netcdf(14) .and. k_Epsilon) then
1224     
1225                call gatherWriteSpx_netcdf(k_turb_g , arr1, arr2 , arr4d, ncid, varid_kturbg , &
1226                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1227     
1228                call gatherWriteSpx_netcdf(e_turb_g , arr1, arr2 , arr4d, ncid, varid_eturbg , &
1229                        imax2,jmax2,kmax2,ijkmax2 , ijkmax3)
1230     
1231             end if
1232     
1233     
1234     
1235             ! Close the file. This frees up any internal netCDF resources
1236             ! associated with the file, and flushes any buffers.
1237     !       CALL MPI_BARRIER(MPI_COMM_WORLD,mpierr)
1238             if (myPE .eq. PE_IO) then
1239                call MFIX_check_netcdf( MFIX_nf90_close(ncid) )
1240             end if
1241     
1242     
1243             deallocate (arr1)
1244             deallocate (arr2)
1245             deallocate (arr4d)
1246             deallocate (xloc)
1247             deallocate (yloc)
1248             deallocate (zloc)
1249     
1250             return
1251             end subroutine write_netcdf
1252     
1253