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