File: /nfs/home/0/users/jenkins/mfix.git/model/check_data/check_bc_inflow.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     ! Subroutine: CHECK_BC_INFLOW                                          !
4     ! Author: J.Musser                                    Date: 01-Mar-14  !
5     !                                                                      !
6     ! Purpose: Provided a detailed error message on common inflow BC       !
7     !                                                                      !
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
9           SUBROUTINE CHECK_BC_INFLOW(M_TOT, SKIP, BCV)
10     
11     ! Modules 
12     !---------------------------------------------------------------------//
13           use bc, only: bc_x_g, bc_x_s
14           use bc, only: bc_t_g, bc_t_s, bc_theta_m
15           use bc, only: bc_k_turb_g, bc_e_turb_g
16           use bc, only: bc_scalar
17           use param, only: dim_m
18           use param1, only: undefined, one, zero 
19           use physprop, only: inert_species
20           use physprop, only: mu_g0
21           use physprop, only: mw_avg, nmax
22           use physprop, only: ro_g0, x_s0
23           use run, only: energy_eq, granular_energy, k_epsilon
24           use run, only: solids_model, solve_ros, species_eq
25           use scalars, only: nscalar
26           use toleranc, only: compare
27     
28           use error_manager
29           IMPLICIT NONE
30     
31     ! Dummy arguments
32     !---------------------------------------------------------------------//
33           INTEGER, INTENT(in) :: BCV
34           INTEGER, INTENT(in) :: M_TOT
35           LOGICAL, INTENT(in) :: SKIP(DIM_M)
36     
37     ! Local variables
38     !---------------------------------------------------------------------//
39     ! loop/variable indices
40           INTEGER :: M, N
41           DOUBLE PRECISION SUM_X
42     ! Index of inert species
43           INTEGER :: INERT
44     !---------------------------------------------------------------------//
45     
46           CALL INIT_ERR_MSG("CHECK_BC_INFLOW")
47     
48     ! Check temperature dependency.
49           IF((ENERGY_EQ .OR. RO_G0==UNDEFINED .OR.MU_G0==UNDEFINED) .AND. &
50              BC_T_G(BCV)==UNDEFINED) THEN
51              WRITE(ERR_MSG, 1000) trim(iVar('BC_T_g',BCV))
52              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
53           ENDIF
54     
55     ! Sum together defiend gas phase species mass fractions.
56           SUM_X = ZERO
57           DO N = 1, NMAX(0)
58              IF(BC_X_G(BCV,N) /= UNDEFINED) THEN
59                 SUM_X = SUM_X + BC_X_G(BCV,N)
60              ELSE
61                 BC_X_G(BCV,N) = ZERO
62              ENDIF
63           ENDDO
64     
65     ! Enforce that the species mass fractions must sum to one.
66           IF(.NOT.COMPARE(ONE,SUM_X)) THEN
67     
68              IF(SPECIES_EQ(0)) THEN
69                 WRITE(ERR_MSG, 1110) BCV
70                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
71      1110 FORMAT('Error 1110: BC_X_g(',I3,',:) do NOT sum to ONE and the ',&
72              'gas phase',/'species equations are solved. Please correct ', &
73              'the mfix.dat file.')
74     
75              ELSEIF(RO_G0 == UNDEFINED .AND. MW_AVG == UNDEFINED) THEN
76                 WRITE(ERR_MSG, 1111) BCV
77                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
78      1111 FORMAT('Error 1111: BC_X_g(',I3,',:) do NOT sum to ONE and the ',&
79              'gas phase',/'is compressible and MW_AVG is UNDEFINED.',/     &
80              'Please correct the mfix.dat the mfix.dat file.')
81     
82              ELSEIF(.NOT.COMPARE(SUM_X,ZERO)) THEN
83                 WRITE(ERR_MSG, 1112) BCV
84                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
85      1112 FORMAT('Error 1112: BC_X_g(',I3,',:) do not sum to ONE or ZERO ',&
86              'and they',/'are not needed. Please correct the mfix.dat ',   &
87              'the mfix.dat file.')
88     
89              ELSE
90                 BC_X_G(BCV,:) = ZERO
91                 BC_X_G(BCV,1) = ONE
92              ENDIF
93           ENDIF
94     
95     
96     ! Verify that species mass fractions are defined for mass flow BCs whe
97     ! using variable solids density. Needed to calculation RO_s
98           DO M = 1, M_TOT
99     
100     ! If this phase is not present, clear out x_s for the BC and
101     ! cycle the solids loop. No need to continue checks.
102              IF(SKIP(M)) THEN
103                IF(SPECIES_EQ(M))THEN
104                    BC_X_S(BCV,M,:) = ZERO
105                    BC_X_S(BCV,M,1) = ONE
106                 ENDIF
107                 CYCLE
108              ENDIF
109     
110     ! Sum together defined species mass fractions.
111              SUM_X = ZERO
112              DO N = 1, NMAX(M)
113                 IF(BC_X_S(BCV,M,N) /= UNDEFINED) THEN
114                    SUM_X = SUM_X + BC_X_S(BCV,M,N)
115                 ELSE
116                    BC_X_S(BCV,M,N) = ZERO
117                 ENDIF
118              ENDDO
119     
120     ! Enforce that the species mass fractions must sum to one.
121              IF(.NOT.COMPARE(ONE,SUM_X)) THEN
122                 IF(SPECIES_EQ(M)) THEN
123                    WRITE(ERR_MSG, 1210) BCV, M
124                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
125      1210 FORMAT('Error 1210: BC_X_s(',I3,',',I2,':) do NOT sum to ONE ',  &
126              'and the solids phase',/'species equations are solved. ',     &
127              'Please correct the mfix.dat file.')
128     
129                 ELSEIF(SOLVE_ROS(M)) THEN
130                    WRITE(ERR_MSG, 1211) BCV, M
131                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
132      1211 FORMAT('Error 1211: BC_X_s(',I3,',',I2,':) do NOT sum to ONE ',  &
133              'and the solids phase',/'density is calculated. Please ',     &
134              'correct the mfix.dat file.')
135     
136                 ELSEIF(.NOT.COMPARE(SUM_X,ZERO)) THEN
137                    WRITE(ERR_MSG, 1212) BCV, M
138                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
139      1212 FORMAT('Error 1212: BC_X_s(',I3,',',I2,':) do NOT sum to ONE ',  &
140              'or ZERO and',/'they are not needed. Please correct the ',    &
141              'mfix.dat file.')
142     
143                 ELSE
144                    BC_X_S(BCV,M,:) = ZERO
145                    BC_X_S(BCV,M,1) = ONE
146                 ENDIF
147              ENDIF
148     
149     ! Set the solids density for the BC region.
150              IF(SOLVE_ROs(M)) THEN
151     ! Verify that the species mass fraction for the inert material is not
152     ! zero in the IC region when the solids is present.
153                 INERT = INERT_SPECIES(M)
154                 IF(BC_X_S(BCV,M,INERT) == ZERO) THEN
155                    WRITE(ERR_MSG,1213) M, BCV
156                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
157      1213 FORMAT('Error 1213: No inert species for phase ',I2,' in BC ',   &
158              'region',I3,'.',/'Unable to calculate solids phase density. ',&
159              'Please refer to the Readme',/' file for required variable ', &
160              'soilds density model input parameters and',/' make the ',   &
161              'necessary corrections to the mfix.dat file.')
162     
163                 ENDIF
164              ENDIF
165           ENDDO
166     
167     
168           DO M = 1, M_TOT
169     ! Check solids phase temperature dependency.
170              IF(ENERGY_EQ .AND. BC_T_S(BCV,M)==UNDEFINED) THEN
171                 IF(SKIP(M)) THEN
172                    BC_T_S(BCV,M) = BC_T_G(BCV)
173                 ELSE
174                    WRITE(ERR_MSG, 1000) trim(iVar('BC_T_s',BCV,M))
175                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
176                 ENDIF
177              ENDIF
178     
179     ! Check granular energy dependency
180              IF(GRANULAR_ENERGY) THEN
181                 IF(BC_THETA_M(BCV,M) == UNDEFINED) THEN
182                    IF(SKIP(M) .OR. SOLIDS_MODEL(M) /= 'TFM') THEN
183                       BC_THETA_M(BCV,M) = ZERO
184                    ELSE
185                       WRITE(ERR_MSG,1000) trim(iVar('BC_Theta_m',BCV,M))
186                       CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
187                    ENDIF
188                 ENDIF
189              ENDIF
190           ENDDO
191     
192     ! Check K-Epsilon BCs.
193           IF(K_Epsilon) THEN
194              IF(BC_K_Turb_G(BCV) == UNDEFINED) THEN
195                 WRITE(ERR_MSG, 1000) trim(iVar('BC_K_Turb_G',BCV))
196                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
197     
198              ELSEIF(BC_E_Turb_G(BCV) == UNDEFINED) THEN
199                 WRITE(ERR_MSG, 1000) trim(iVar('BC_E_Turb_G',BCV))
200                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
201              ENDIF
202           ENDIF
203     
204     ! Check scalar equation BCs.
205           DO N = 1, NScalar
206              IF(BC_Scalar(BCV,N) == UNDEFINED) THEN
207                 WRITE(ERR_MSG, 1001) trim(iVar('BC_Scalar',BCV,N))
208                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
209              ENDIF
210           ENDDO
211     
212           CALL FINL_ERR_MSG
213     
214           RETURN
215     
216      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
217              'correct the mfix.dat file.')
218     
219      1001 FORMAT('Error 1001: Illegal or unknown input: ',A,' = ',A,/   &
220              'Please correct the mfix.dat file.')
221     
222           END SUBROUTINE CHECK_BC_INFLOW
223     
224     
225     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
226     !                                                                      !
227     ! Subroutine: CHECK_BC_MASS_INFLOW                                     !
228     ! Author: J.Musser                                    Date: 01-Mar-14  !
229     !                                                                      !
230     ! Purpose: Provided a detailed error message on BC                     !
231     !                                                                      !
232     ! Comments:                                                            !
233     !     The velocities at the inflow face are fixed and the momentum     !
234     !     equations are not solved in the inflow cells. Since the flow is  !
235     !     into the domain all other scalars that are used need to be       !
236     !     specified (e.g., mass fractions, void fraction, etc.,)           !
237     !                                                                      !
238     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
239           SUBROUTINE CHECK_BC_MASS_INFLOW(M_TOT, SKIP, BCV)
240     
241     ! Modules 
242     !---------------------------------------------------------------------//
243           use bc, only: bc_ep_g, bc_p_g
244           use bc, only: bc_rop_s, bc_ep_s, bc_x_s
245           use bc, only: bc_massflow_g
246           use eos, only: eoss
247           use param, only: dim_m
248           use param1, only: undefined, one, zero 
249           use physprop, only: ro_g0, ro_s0
250           use physprop, only: inert_species, x_s0, base_ros
251           use run, only: solve_ros
252           use toleranc, only: compare
253           use error_manager
254           IMPLICIT NONE
255     
256     ! Dummy arguments
257     !---------------------------------------------------------------------//
258           INTEGER, INTENT(in) :: BCV
259           INTEGER, INTENT(in) :: M_TOT
260           LOGICAL, INTENT(in) :: SKIP(DIM_M)
261     
262     ! Local variables
263     !---------------------------------------------------------------------//
264     ! loop/variable indices
265           INTEGER :: M
266           DOUBLE PRECISION :: SUM_EP
267     ! Solids phase density in BC region.
268           DOUBLE PRECISION :: BC_ROs(DIM_M)
269     ! Index of inert species
270           INTEGER :: INERT
271     
272     !---------------------------------------------------------------------//
273     
274           CALL INIT_ERR_MSG("CHECK_BC_MASS_INFLOW")
275     
276     ! Check gas phase volume fraction.
277           IF(BC_EP_G(BCV) == UNDEFINED) THEN
278              WRITE(ERR_MSG, 1000) trim(iVar('BC_EP_g',BCV))
279              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
280           ENDIF
281     
282     ! Verify compressible boundary condition variables.
283           IF(RO_G0 == UNDEFINED) THEN
284              IF(BC_P_G(BCV) == UNDEFINED) THEN
285                 IF(BC_MASSFLOW_G(BCV) /= UNDEFINED .AND.                   &
286                    BC_MASSFLOW_G(BCV) /= ZERO) THEN
287                    WRITE(ERR_MSG, 1100) trim(iVar('BC_P_g',BCV))
288                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
289                 ENDIF
290      1100 FORMAT('Error 1100: ',A,' must be specified for compressible ',  &
291              'flows',/'when specifying BC_MASSFLOW_g to make the ',        &
292              'conversion to velocity.',/'Please correct the mfix.dat file.')
293     
294              ELSEIF(BC_P_G(BCV) <= ZERO) THEN
295                 WRITE(ERR_MSG, 1101) BCV, trim(iVal(BC_P_G(BCV)))
296                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
297              ENDIF
298      1101 FORMAT('Error 1101: Pressure must be greater than zero for ',    &
299              'compressible flow',/' >>>  BC_P_g(',I3,') = ',A,/'Please ',  &
300              'correct the mfix.dat file.')
301           ENDIF
302     
303     ! Calculate the solids volume fraction from the gas phase if there is
304     ! only one solids phase.
305           IF(M_TOT == 1 .AND. BC_EP_S(BCV,1) == UNDEFINED) THEN
306              BC_EP_S(BCV,1) = ONE - BC_EP_g(BCV)
307           ENDIF
308     
309     ! Bulk density or solids volume fraction must be explicitly defined
310     ! if there are more than one solids phase.
311           IF(M_TOT > 1 .AND. .NOT.COMPARE(BC_EP_g(BCV),ONE)) THEN
312              DO M = 1, M_TOT
313                 IF(BC_ROP_S(BCV,M) == UNDEFINED .AND. &
314                    BC_EP_S(BCV,M) == UNDEFINED) THEN
315                    WRITE(ERR_MSG, 1200) M, BCV, 'BC_ROP_s and BC_EP_s'
316                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
317                 ENDIF
318              ENDDO
319           ENDIF
320      1200 FORMAT('Error 1200: Insufficient solids phase ',I2,' data ',     &
321              'for BC',I3,'. ',/A,' not specified.',/'Please correct the ', &
322              'mfix.dat file.')
323     
324     ! Initialize the sum of the total volume fraction.
325           SUM_EP = BC_EP_G(BCV)
326           DO M = 1, M_TOT
327     
328     ! If this phase is not present, clear out EPs and ROPs for the BC and
329     ! cycle the solids loop. No need to continue checks.
330              IF(SKIP(M)) THEN
331                 BC_EP_S(BCV,M)  = ZERO
332                 BC_ROP_S(BCV,M) = ZERO
333                 CYCLE
334              ENDIF
335     
336     ! Set the solids density for the BC region.
337              IF(.NOT.SOLVE_ROs(M)) THEN
338                 BC_ROs(M) = RO_s0(M)
339              ELSE  
340     ! presence of non-zero inert species is checked by bc_inflow
341                 INERT = INERT_SPECIES(M)
342                 BC_ROs(M) = EOSS(BASE_ROs(M), X_s0(M,INERT),&
343                    BC_X_S(BCV,M,INERT))
344              ENDIF
345     
346     ! If both input parameters are defined. Make sure they are equivalent.
347              IF(BC_ROP_S(BCV,M) /= UNDEFINED .AND.                         &
348                 BC_EP_S(BCV,M) /= UNDEFINED) THEN
349     
350                 IF(.NOT.COMPARE(BC_EP_S(BCV,M)*BC_ROs(M),                  &
351                    BC_ROP_S(BCV,M))) THEN
352                    WRITE(ERR_MSG,1214) BCV
353                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
354                 ENDIF
355      1214 FORMAT('Error 1214: Illegal initial condition region : ',I3,/    &
356              'BC_EP_s and BC_ROP_s are inconsistent. Please correct the ',/&
357              'mfix.dat file.')
358     
359     ! Compute BC_EP_s from BC_ROP_s
360              ELSEIF(BC_EP_S(BCV,M) == UNDEFINED) THEN
361                 BC_EP_S(BCV,M) = BC_ROP_S(BCV,M) / BC_ROs(M)
362     
363     ! Compute BC_ROP_s from BC_EP_s and BC_ROs
364              ELSEIF(BC_ROP_S(BCV,M) == UNDEFINED) THEN
365                 BC_ROP_S(BCV,M) = BC_EP_S(BCV,M) * BC_ROs(M)
366     
367              ENDIF
368     ! Add this phase to the total volume fraction.
369              SUM_EP = SUM_EP + BC_EP_S(BCV,M)
370           ENDDO
371     
372     ! Verify that the volume fractions sum to one.
373           IF(.NOT.COMPARE(SUM_EP,ONE)) THEN
374              WRITE(ERR_MSG,1215) BCV, trim(iVal(SUM_EP))
375              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
376           ENDIF
377      1215 FORMAT('Error 1215: Illegal boundary condition region: ',I3,'. ',&
378              'Sum of volume',/'fractions does NOT equal ONE. (SUM = ',A,   &
379              ')',/'Please correct the mfix.dat file.')
380     
381           CALL FINL_ERR_MSG
382     
383           RETURN
384     
385      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
386              'correct the mfix.dat file.')
387     
388      1001 FORMAT('Error 1001: Illegal or unknown input: ',A,' = ',A,/   &
389              'Please correct the mfix.dat file.')
390     
391           END SUBROUTINE CHECK_BC_MASS_INFLOW
392     
393     
394     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
395     !                                                                      !
396     ! Subroutine: CHECK_BC_P_INFLOW                                        !
397     ! Author: J.Musser                                    Date: 01-Mar-14  !
398     !                                                                      !
399     ! Purpose: Provided detailed error message on bc                       !
400     !                                                                      !
401     ! Comments:                                                            !
402     !     Unlike the MI boundary, for the PI boundary the velocities at    !
403     !     the inflow face are calculated by solving the momentum eqns      !
404     !     and are not fixed. In this way, the PI is similar to the PO      !
405     !     except that the flow is into the domain and hence all other      !
406     !     scalars (e.g., mass fractions, void fraction, temperature,       !
407     !     etc.,) at the inflow cells need to be specified. To satisfy      !
408     !     the error routines at the start of the simulation, both the      !
409     !     tangential and normal components at the inflow also need to      !
410     !     be specified. The velocities values essentially serve as IC.     !
411     !                                                                      !
412     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
413           SUBROUTINE CHECK_BC_P_INFLOW(M_TOT, SKIP, BCV)
414     
415     ! Modules
416     !---------------------------------------------------------------------//
417           use bc, only: bc_p_g, bc_rop_s
418           use bc, only: bc_u_g, bc_v_g, bc_w_g
419           use bc, only: bc_u_s, bc_v_s, bc_w_s
420           use geometry, only: no_i, no_j, no_k
421           use param, only: dim_m
422           use param1, only: undefined, one, zero 
423           use physprop, only: ro_g0
424           use error_manager
425           IMPLICIT NONE
426     
427     ! Dummy arguments
428     !---------------------------------------------------------------------//
429           INTEGER, INTENT(in) :: BCV
430           INTEGER, INTENT(in) :: M_TOT
431           LOGICAL, INTENT(in) :: SKIP(DIM_M)
432     
433     ! Local variables
434     !---------------------------------------------------------------------//
435           INTEGER :: M
436     
437     !---------------------------------------------------------------------//
438     
439           CALL INIT_ERR_MSG("CHECK_BC_P_INFLOW")
440     
441     ! Remove checks on bc_ep_g/bc_rop_s; using routine check_bc_outflow
442     
443           IF (BC_P_G(BCV) == UNDEFINED) THEN
444              WRITE(ERR_MSG,1000) 'BC_P_g', BCV
445              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
446      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
447              'correct the mfix.dat file.')
448     
449           ELSEIF (BC_P_G(BCV)<=ZERO .AND. RO_G0==UNDEFINED) THEN
450              WRITE(ERR_MSG, 1101) BCV, trim(iVal(BC_P_G(BCV)))
451              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
452      1101 FORMAT('Error 1101: Pressure must be greater than zero for ',    &
453              'compressible flow',/' >>>  BC_P_g(',I3,') = ',A,/'Please ',  &
454              'correct the mfix.dat file.')
455           ENDIF
456     
457     ! Check that velocities are also specified. These are essentially used
458     ! as initial conditions for the boundary region. If they are not
459     ! specified then a default value is set here otherwise check_data_20
460     ! will complain and cause MFIX to exit.
461           IF(BC_U_G(BCV) == UNDEFINED) THEN
462              BC_U_G(BCV) = ZERO
463              IF(.NOT.NO_I) WRITE(ERR_MSG, 1300) trim(iVar('BC_U_g',BCV))
464           ENDIF
465     
466           IF(BC_V_G(BCV) == UNDEFINED) THEN
467              BC_V_G(BCV) = ZERO
468              IF(.NOT.NO_J) WRITE(ERR_MSG, 1300) trim(iVar('BC_V_g',BCV))
469           ENDIF
470     
471           IF(BC_W_G(BCV) == UNDEFINED) THEN
472              BC_W_G(BCV) = ZERO
473              IF(.NOT.NO_K) WRITE(ERR_MSG, 1300) trim(iVar('BC_W_g',BCV))
474           ENDIF
475     
476           DO M = 1, M_TOT
477              IF (SKIP(M)) THEN
478                 BC_U_S(BCV,M) = ZERO
479                 BC_V_S(BCV,M) = ZERO
480                 BC_W_S(BCV,M) = ZERO
481              ELSE
482                 IF(BC_U_S(BCV,M) == UNDEFINED) THEN
483                    BC_U_S(BCV,M) = ZERO
484                    IF(BC_ROP_S(BCV,M) /= ZERO .AND. .NOT.NO_I) &
485                       WRITE(ERR_MSG, 1300) trim(iVar('BC_U_s',BCV,M))
486                 ENDIF
487     
488                 IF(BC_V_S(BCV,M) == UNDEFINED) THEN
489                    BC_V_S(BCV,M) = ZERO
490                    IF(BC_ROP_S(BCV,M) /= ZERO .AND. .NOT.NO_J) &
491                       WRITE(ERR_MSG, 1300) trim(iVar('BC_V_s',BCV,M))
492                 ENDIF
493     
494                 IF(BC_W_S(BCV,M) == UNDEFINED) THEN
495                    BC_W_S(BCV,M) = ZERO
496                    IF(BC_ROP_S(BCV,M) /= ZERO .AND. .NOT.NO_K) &
497                       WRITE(ERR_MSG, 1300) trim(iVar('BC_W_s',BCV,M))
498                 ENDIF
499              ENDIF
500           ENDDO
501     
502      1300 FORMAT('Warning 1300: ',A,' was undefined. This variable was ', &
503              'set ',/ 'to zero to be used as the inital value in the BC ',&
504              'region.')
505     
506           CALL FINL_ERR_MSG
507     
508           RETURN
509           END SUBROUTINE CHECK_BC_P_INFLOW
510