File: N:\mfix\model\utilities_mod.f

1     MODULE utilities
2     
3       IMPLICIT NONE
4     
5     CONTAINS
6     
7     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
8     !                                                                      !
9     !  function: mfix_isnan                                                !
10     !  Purpose: check whether argument is NAN                              !
11     !                                                                      !
12     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
13           LOGICAL FUNCTION mfix_isnan(x)
14     
15     ! Dummy arguments
16     !---------------------------------------------------------------------//
17           double precision, intent(in) :: x
18     
19     ! Local variables
20     !---------------------------------------------------------------------//
21           CHARACTER(LEN=80) :: notnumber
22     !---------------------------------------------------------------------//
23     
24           mfix_isnan = .False.
25           WRITE(notnumber,*) x
26     ! To check for NaN's in x, see if x (a real number) contain a letter "N"
27     ! "n" or symbol "?", in which case it is a NaN (Not a Number)
28     
29           IF(INDEX(notnumber,'?') > 0 .OR.     &
30              INDEX(notnumber,'n') > 0 .OR.     &
31              INDEX(notnumber,'N') > 0 ) THEN
32             mfix_isnan = .TRUE.
33              RETURN
34           ENDIF
35     
36           RETURN
37           END FUNCTION mfix_isnan
38     
39     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
40     !                                                                      C
41     !  function: MAX_VEL_INLET                                             C
42     !  Purpose: Find maximum velocity at inlets.                           C
43     !                                                                      C
44     !  Author: S. Benyahia                                Date: 26-AUG-05  C
45     !                                                                      C
46     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
47           DOUBLE PRECISION FUNCTION MAX_VEL_INLET()
48     
49     ! Modules
50     !---------------------------------------------------------------------//
51           use bc, only: bc_defined, bc_type_enum
52           use bc, only: mass_inflow, p_inflow
53           use bc, only: bc_plane
54           use bc, only: bc_k_b, bc_k_t, bc_j_s, bc_j_n, bc_i_w, bc_i_e
55           use compar, only: dead_cell_at
56           use fldvar, only: u_g, v_g, w_g, u_s, v_s, w_s
57           use functions, only: funijk, is_on_mype_owns
58           use functions, only: im_of, jm_of, km_of
59           use mpi_utility, only: global_all_max
60           use param, only: dimension_bc
61           use param1, only: zero, small_number
62           use physprop, only: mmax
63           use run, only: units
64           use toleranc, only: max_allowed_vel, max_inlet_vel, max_inlet_vel_fac
65           IMPLICIT NONE
66     
67     ! Local variables
68     !---------------------------------------------------------------------//
69           INTEGER :: L, I, J, K, IJK, M
70           DOUBLE PRECISION :: maxVEL
71     !---------------------------------------------------------------------//
72     
73     ! initializing
74           maxVEL = ZERO
75     
76           DO L = 1, DIMENSION_BC
77     
78              IF(.NOT.(BC_DEFINED(L).AND.(BC_TYPE_ENUM(L)==MASS_INFLOW .OR. &
79                 BC_TYPE_ENUM(L) == P_INFLOW))) CYCLE
80     
81              DO K = BC_K_B(L), BC_K_T(L)
82              DO J = BC_J_S(L), BC_J_N(L)
83              DO I = BC_I_W(L), BC_I_E(L)
84     
85                 IF (.NOT.IS_ON_myPE_OWNS(I,J,K)) CYCLE
86                 IF (DEAD_CELL_AT(I,J,K)) CYCLE
87                 IJK = FUNIJK(I,J,K)
88     
89                 SELECT CASE (BC_PLANE(L))
90                 CASE ('E'); maxVEL = max(maxVEL,abs(U_G(IJK)))
91                 CASE ('N'); maxVEL = max(maxVEL,abs(V_G(IJK)))
92                 CASE ('T'); maxVEL = max(maxVEL,abs(W_G(IJK)))
93                 CASE ('W'); maxVEL = max(maxVEL,abs(U_G(IM_OF(IJK))))
94                 CASE ('S'); maxVEL = max(maxVEL,ABS(V_G(JM_OF(IJK))))
95                 CASE ('B'); maxVEL = max(maxVEL,abs(W_G(KM_OF(IJK))))
96                 END SELECT
97     
98                 DO M=1,MMAX
99                    SELECT CASE (BC_PLANE(L))
100                    CASE ('E'); maxVEL = max(maxVEL,abs(U_s(IJK,M)))
101                    CASE ('N'); maxVEL = max(maxVEL,abs(V_s(IJK,M)))
102                    CASE ('T'); maxVEL = max(maxVEL,abs(W_s(IJK,M)))
103                    CASE ('W'); maxVEL = max(maxVEL,abs(U_s(IM_OF(IJK),M)))
104                    CASE ('S'); maxVEL = max(maxVEL,abs(V_s(JM_OF(IJK),M)))
105                    CASE ('B'); maxVEL = max(maxVEL,abs(W_s(KM_OF(IJK),M)))
106                    END SELECT
107                 ENDDO
108     
109              ENDDO
110              ENDDO
111              ENDDO
112           ENDDO
113     
114           CALL GLOBAL_ALL_MAX(maxVEL, MAX_VEL_INLET)
115     
116     ! If no inlet velocity is specified, use an upper limit defined in
117     ! toleranc_mod.f
118           IF(MAX_VEL_INLET <= SMALL_NUMBER) THEN
119              MAX_VEL_INLET = MAX_ALLOWED_VEL
120              IF(UNITS == 'SI') MAX_INLET_VEL = 1D-2*MAX_ALLOWED_VEL
121           ELSE
122     ! Scale the value using a user defined scale factor
123              MAX_VEL_INLET = 100.0d0*MAX_INLET_VEL_FAC*MAX_VEL_INLET
124           ENDIF
125     
126           RETURN
127           END FUNCTION MAX_VEL_INLET
128     
129     
130     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
131     !                                                                      C
132     !  Function: CHECK_VEL_BOUND()                                         C
133     !  Purpose: Check velocities upper bound to be less than speed of      C
134     !           sound                                                      C
135     !                                                                      C
136     !  Author: S. Benyahia                                Date: 25-AUG-05  C
137     !                                                                      C
138     !                                                                      C
139     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
140           LOGICAL FUNCTION CHECK_VEL_BOUND ()
141     
142     ! Modules
143     !---------------------------------------------------------------------//
144           use compar, only: ijkstart3, ijkend3
145           use discretelement, only: des_continuum_coupled, des_continuum_hybrid
146           use fldvar, only: ep_g, u_g, v_g, w_g
147           use fldvar, only: ep_s, u_s, v_s, w_s
148           use functions, only: fluid_at 
149           use indices, only: i_of, j_of, k_of
150           use mpi_utility, only: global_all_or
151           use physprop, only: mmax
152           use toleranc, only: max_inlet_vel
153     
154           IMPLICIT NONE
155     
156     ! Local variables
157     !---------------------------------------------------------------------//
158           INTEGER :: M
159     ! Indices
160           INTEGER :: IJK
161           LOGICAL :: ALL_IS_ERROR
162     !---------------------------------------------------------------------//
163     
164     !!$omp   parallel do private(IJK)
165     ! initializing
166           CHECK_VEL_BOUND = .FALSE.
167           ALL_IS_ERROR    = .FALSE.
168     
169     LOOP_FLUID : DO IJK = IJKSTART3, IJKEND3
170     
171              IF (FLUID_AT(IJK)) THEN
172                 IF(ABS(U_G(IJK)) > MAX_INLET_VEL .OR. &
173                    ABS(V_G(IJK)) > MAX_INLET_VEL .OR. &
174                    ABS(W_G(IJK)) > MAX_INLET_VEL) THEN
175                    CHECK_VEL_BOUND = .TRUE.
176                    WRITE(*,1000) MAX_INLET_VEL, I_OF(IJK), J_OF(IJK), K_OF(IJK), &
177                                  EP_g(IJK), U_G(IJK), V_G(IJK), W_G(IJK)
178                    EXIT LOOP_FLUID
179                 ENDIF
180     
181                 IF (.NOT.DES_CONTINUUM_COUPLED .OR. DES_CONTINUUM_HYBRID) THEN
182                    DO M = 1, MMAX
183                      IF(ABS(U_S(IJK,M)) > MAX_INLET_VEL .OR. &
184                         ABS(V_S(IJK,M)) > MAX_INLET_VEL .OR. &
185                         ABS(W_S(IJK,M)) > MAX_INLET_VEL) THEN
186                        CHECK_VEL_BOUND = .TRUE.
187                        WRITE(*,1010) MAX_INLET_VEL, I_OF(IJK), J_OF(IJK), K_OF(IJK), M, &
188                                      EP_s(IJK, M), U_S(IJK,M), V_S(IJK,M), W_S(IJK,M)
189                        EXIT LOOP_FLUID
190                      ENDIF
191                    ENDDO
192                 ENDIF   ! end if(.not.des_continuum_coupled or des_continuum_hybrid)
193              ENDIF
194     
195           ENDDO LOOP_FLUID
196     
197           CALL GLOBAL_ALL_OR(CHECK_VEL_BOUND, ALL_IS_ERROR)
198           IF(ALL_IS_ERROR) CHECK_VEL_BOUND = .TRUE.
199     
200           RETURN
201      1000 FORMAT(1X,'Message from: CHECK_VEL_BOUND',/&
202                 'WARNING: velocity higher than maximum allowed velocity: ', &
203                 G12.5, '(to change this adjust the scale factor MAX_INLET_VEL_FAC)'/&
204                 'in this cell: ','I = ',I4,2X,' J = ',I4,2X,' K = ',I4, /&
205                 '  ','Epg = ', G12.5, 'Ug = ', G12.5, 'Vg = ', G12.5, 'Wg = ', G12.5)
206      1010 FORMAT(1X,'Message from: CHECK_VEL_BOUND',/&
207                 'WARNING: velocity higher than maximum allowed velocity: ', &
208                 G12.5,/&
209                 'in this cell: ','I = ',I4,2X,' J = ',I4,2X,' K = ',I4,' M = ',I4, /&
210                 '  ','Eps = ', G12.5,'Us = ', G12.5, 'Vs = ', G12.5, 'Ws = ', G12.5)
211     
212           END FUNCTION CHECK_VEL_BOUND
213     
214     
215     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
216     !                                                                      !
217     !  Function: SEEK_COMMENT                                              !
218     !  Purpose: determine if (and where) a comment character appears       !
219     !           in a data input line                                       !
220     !     The function SEEK_COMMENT returns the index to where a comment   !
221     !     character was found in the input data line.  Equals MAXCOL + 1   !
222     !     if no-comment characters in the line.                            !
223     !                                                                      !
224     !  Author: P.Nicoletti                                Date: 25-NOV-91  !
225     !                                                                      !
226     !                                                                      !
227     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
228           INTEGER FUNCTION SEEK_COMMENT (LINE, MAXCOL)
229           IMPLICIT NONE
230     
231     ! Dummy arguments
232     !---------------------------------------------------------------------//
233     ! input data line
234           CHARACTER(len=*) :: LINE
235     ! maximum column of input data line to search
236           INTEGER :: MAXCOL
237     
238     ! Local parameters
239     !---------------------------------------------------------------------//
240     ! the number of designated comment characters
241           INTEGER, PARAMETER :: DIM_COMMENT = 2
242     
243     ! Local variables
244     !---------------------------------------------------------------------//
245     ! loop indicies
246           INTEGER :: L, L2
247     ! the comment characters
248           CHARACTER, DIMENSION(DIM_COMMENT) :: COMMENT_CHAR
249     !---------------------------------------------------------------------//
250           DATA COMMENT_CHAR/'#', '!'/
251     
252           DO L = 1, MAXCOL
253              DO L2 = 1, DIM_COMMENT
254                 IF (LINE(L:L) == COMMENT_CHAR(L2)) THEN
255                    SEEK_COMMENT = L
256                    RETURN
257                 ENDIF
258              END DO
259           END DO
260           SEEK_COMMENT = MAXCOL + 1
261           RETURN
262           END FUNCTION SEEK_COMMENT
263     
264     
265     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
266     !                                                                      !
267     !  Function: SEEK_END                                                  !
268     !  Purpose: determine where trailing blanks begin in a line            !
269     !     The function SEEK_END returns the index to where the last        !
270     !     character was found in the input data line.  Equals MAXCOL       !
271     !     if no trailing blank characters in the line.                     !
272     !                                                                      !
273     !  Author: P.Nicoletti, M. Syamlal                    Date: 7-AUG-92   !
274     !                                                                      !
275     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
276           INTEGER FUNCTION SEEK_END (LINE, MAXCOL)
277           IMPLICIT NONE
278     
279     ! Dummy arguments
280     !---------------------------------------------------------------------//
281     ! maximum column of input data line to search
282           INTEGER :: MAXCOL
283     ! input data line
284           CHARACTER :: LINE*(*)
285     
286     ! Local variables
287     !---------------------------------------------------------------------//
288           INTEGER :: L
289     !---------------------------------------------------------------------//
290           SEEK_END = 0
291           DO L = 1, MAXCOL
292              IF (LINE(L:L) /= ' ') SEEK_END = L
293           END DO
294           RETURN
295           END FUNCTION SEEK_END
296     
297     
298     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
299     !                                                                      !
300     !  Function: LINE_TOO_BIG                                              !
301     !  Purpose: return an error condition if input data is located past    !
302     !           column MAXCOL in the data input file                       !
303     !     The function LINE_TOO_BIG returns a value greater than 0 to      !
304     !     indicate an error condition (data passed column MAXCOL in LINE)  !
305     !                                                                      !
306     !  Author: P.Nicoletti                                Date: 25-NOV-91  !
307     !                                                                      !
308     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
309           INTEGER FUNCTION LINE_TOO_BIG (LINE, LINE_LEN, MAXCOL)
310           IMPLICIT NONE
311     
312     ! Dummy arguments
313     !---------------------------------------------------------------------//
314     ! input data line
315           CHARACTER(LEN=*) :: LINE
316     ! length of input data line
317           INTEGER :: LINE_LEN
318     ! maximum column that non-blank charcaters are
319     ! are in the input data line
320           INTEGER :: MAXCOL
321     
322     ! Local variables
323     !---------------------------------------------------------------------//
324     ! Loop index
325           INTEGER :: L
326     !---------------------------------------------------------------------//
327           DO L = MAXCOL + 1, LINE_LEN
328              IF (LINE(L:L) /= ' ') THEN
329                 LINE_TOO_BIG = L
330                 RETURN
331              ENDIF
332           ENDDO
333           LINE_TOO_BIG = 0
334           RETURN
335           END FUNCTION LINE_TOO_BIG
336     
337     
338     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
339     !                                                                      !
340     ! Function: BLANK_LINE                                                 !
341     ! Author: P. Nicoletti                                Date: 25-NOV-91  !
342     !                                                                      !
343     ! Purpose: Return .TRUE. if a line contains no input or only spaces.   !
344     !                                                                      !
345     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
346           LOGICAL FUNCTION BLANK_LINE (line)
347           IMPLICIT NONE
348     
349     ! Dummy arguments
350     !---------------------------------------------------------------------//
351           CHARACTER :: LINE*(*)
352     
353     ! Local variables
354     !---------------------------------------------------------------------//
355           INTEGER :: L
356     !---------------------------------------------------------------------//
357           BLANK_LINE = .FALSE.
358           DO L=1, len(line)
359              IF(line(L:L)/=' ' .and. line(L:L)/='    ')RETURN
360           ENDDO
361           BLANK_LINE = .TRUE.
362           RETURN
363           END FUNCTION BLANK_LINE
364     
365     
366     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
367     !                                                                      !
368     !  Subroutine: bound_x                                                 !
369     !  Purpose: bound the values of x array                                !
370     !                                                                      !
371     !  Author: M. Syamlal                                 Date: 15-SEP-98  !
372     !                                                                      !
373     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
374           SUBROUTINE BOUND_X(ARRAY, IJKMAX2)
375     
376     ! Modules
377     !---------------------------------------------------------------------//
378           USE param, only: dimension_3
379           USE param1, only: one, zero
380           IMPLICIT NONE
381     
382     ! Dummy arguments
383     !---------------------------------------------------------------------//
384     ! Maximum dimension
385           INTEGER :: IJKMAX2
386     ! Array
387           DOUBLE PRECISION :: Array(DIMENSION_3)
388     !--------------------------------------------------------------------//
389           IF (IJKMAX2 > 0) THEN
390              ARRAY(:) = MIN(ONE,MAX(ZERO,ARRAY(:)))
391           ENDIF
392           RETURN
393           END SUBROUTINE BOUND_X
394     
395     END MODULE utilities
396