File: RELATIVE:/../../../mfix.git/model/check_data/check_solids_dem.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  SUBROUTINE: CHECK_DES_SOLIDS                                        !
4     !  Author: J.Musser                                   Date: 02-FEB-14  !
5     !                                                                      !
6     !  Purpose:                                                            !
7     !                                                                      !
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
9           SUBROUTINE CHECK_SOLIDS_DEM
10     
11     ! Use the error manager for posting error messages.
12     !---------------------------------------------------------------------//
13           use error_manager
14     
15           implicit none
16     
17     ! Initialize the error manager.
18           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM")
19     
20     ! Particle-particle collision parameters.
21           CALL CHECK_SOLIDS_DEM_COLLISION
22     ! DES cohesion model parameters.
23           CALL CHECK_SOLIDS_DEM_COHESION
24     ! Particle-particle conduction model parameters.
25           CALL CHECK_SOLIDS_DEM_ENERGY
26     
27           CALL FINL_ERR_MSG
28     
29           RETURN
30     
31           END SUBROUTINE CHECK_SOLIDS_DEM
32     
33     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
34     !                                                                      !
35     !  SUBROUTINE CHECK_SOLIDS_DEM_ENERGY                                  !
36     !  Author: J.Musser                                   Date: 02-FEB-14  !
37     !                                                                      !
38     !  Purpose:                                                            !
39     !                                                                      !
40     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
41           SUBROUTINE CHECK_SOLIDS_DEM_ENERGY
42     
43           use des_thermo, only: DES_MIN_COND_DIST
44           use run, only: UNITS
45     
46     ! Global Parameters:
47     !---------------------------------------------------------------------//
48           use param1, only: UNDEFINED
49     
50           use error_manager
51     
52           IMPLICIT NONE
53     
54     
55     ! Initialize the error manager.
56           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM_ENERGY")
57     
58     ! Conduction Equations:
59     !'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
60     
61     ! Set the default value for the minimum distance separating particles'
62     ! surfaces.
63           IF(DES_MIN_COND_DIST == UNDEFINED)THEN
64              DES_MIN_COND_DIST = 1.0D-04 ! cm
65              IF (UNITS == 'SI') DES_MIN_COND_DIST = &
66                 DES_MIN_COND_DIST/100.0  ! m
67           ENDIF
68     
69           CALL FINL_ERR_MSG
70     
71           RETURN
72     
73           END SUBROUTINE CHECK_SOLIDS_DEM_ENERGY
74     
75     
76     
77     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
78     !                                                                      !
79     !  Module name: CHECK_SOLIDS_DEM_COHESION                              !
80     !  Author: J.Musser                                   Date: 11-DEC-13  !
81     !                                                                      !
82     !  Purpose: Check/set parameters for DES cohesion models.              !
83     !                                                                      !
84     !  Comments: Original code moved from CHECK_DES_DATA                   !
85     !                                                                      !
86     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
87           SUBROUTINE CHECK_SOLIDS_DEM_COHESION
88     
89     ! Global Variables:
90     !---------------------------------------------------------------------//
91     
92     ! Runtime Flag: Invoke a cohesion model for DES simulation
93           use discretelement, only: USE_COHESION
94     ! Largest discrete particle diameter.
95           use discretelement, only: MAX_RADIUS
96     
97     ! Runtime Flag: Invoke Square Well
98           use discretelement, only: SQUARE_WELL
99     ! Runtime Flag: Invoke Van der Waals model.
100           use discretelement, only: VAN_DER_WAALS
101     ! User specified parameter for increase neighbor search area.
102           use discretelement, only: FACTOR_RLM
103     
104     ! User specified parameters for Van der Waals model.
105           use discretelement, only: VDW_INNER_CUTOFF
106           use discretelement, only: VDW_OUTER_CUTOFF
107           use discretelement, only: HAMAKER_CONSTANT
108           use discretelement, only: WALL_VDW_INNER_CUTOFF
109           use discretelement, only: WALL_VDW_OUTER_CUTOFF
110           use discretelement, only: WALL_HAMAKER_CONSTANT
111           use discretelement, only: ASPERITIES
112           use discretelement, only: SURFACE_ENERGY
113           use discretelement, only: WALL_SURFACE_ENERGY
114     
115     
116     ! Global Parameters:
117     !---------------------------------------------------------------------//
118           use param1, only: UNDEFINED
119           use param1, only: ZERO
120           use constant, only: Pi
121     
122           use error_manager
123     
124           IMPLICIT NONE
125     
126     ! Local Variables:
127     !---------------------------------------------------------------------//
128     ! Neighborhood size for Van der Waals force.
129           DOUBLE PRECISION :: VDW_NEIGHBORHOOD
130     
131     
132     ! Override the following settings if cohesion not used.
133           IF(.NOT.USE_COHESION) THEN
134     !No more square well in the code
135              SQUARE_WELL = .FALSE.
136              VAN_DER_WAALS = .FALSE.
137              WALL_VDW_OUTER_CUTOFF = ZERO
138              RETURN
139           ENDIF
140     
141     
142     ! Initialize the error manager.
143           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM_COHESION")
144     
145     
146     ! Verify that only one cohesion model is specified.
147           IF (SQUARE_WELL .AND. VAN_DER_WAALS) THEN
148              WRITE(ERR_MSG,1002)
149              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
150     
151     ! Verify that at a cohesion model is specified.
152           ELSEIF(.NOT.SQUARE_WELL .AND. .NOT.VAN_DER_WAALS) THEN
153              WRITE(ERR_MSG,1003)
154              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
155           ENDIF
156     
157     
158     ! Van der Waals model checks.
159           IF (VAN_DER_WAALS) THEN
160     
161              IF (VDW_INNER_CUTOFF .EQ. UNDEFINED) THEN
162                 WRITE(ERR_MSG,1201) 'VDW_INNER_CUTOFF'
163                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
164              ENDIF
165     
166              IF(VDW_OUTER_CUTOFF .EQ. UNDEFINED) THEN
167                 WRITE(ERR_MSG,1201) 'VDW_OUTER_CUTOFF'
168                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
169              ENDIF
170     
171              IF(HAMAKER_CONSTANT .EQ. UNDEFINED) THEN
172                 WRITE(ERR_MSG,1201) 'HAMAKER_CONSTANT'
173                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
174              ENDIF
175     
176              IF (WALL_VDW_INNER_CUTOFF .EQ. UNDEFINED)THEN
177                 WRITE(ERR_MSG,1201) 'WALL_VDW_INNER_CUTOFF'
178                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
179              ENDIF
180     
181              IF (WALL_VDW_OUTER_CUTOFF .EQ. UNDEFINED)THEN
182                 WRITE(ERR_MSG,1201) 'WALL_VDW_OUTER_CUTOFF'
183                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
184              ENDIF
185     
186              IF(WALL_HAMAKER_CONSTANT .EQ. UNDEFINED) THEN
187                 WRITE(ERR_MSG,1201) 'WALL_HAMAKER_CONSTANT'
188                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
189              ENDIF
190     
191              VDW_NEIGHBORHOOD = 1.0d0 + (VDW_OUTER_CUTOFF/(2.d0*MAX_RADIUS))
192              IF (FACTOR_RLM < VDW_NEIGHBORHOOD) THEN
193                 WRITE(ERR_MSG,1202)
194                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
195              ENDIF
196     
197              IF (ASPERITIES < ZERO) THEN
198                 WRITE(ERR_MSG,1001) 'ASPERITIES', trim(iVal(ASPERITIES))
199                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
200              ENDIF
201     
202              SURFACE_ENERGY=HAMAKER_CONSTANT/&
203                 (24.d0*Pi*VDW_INNER_CUTOFF**2)
204     
205              WALL_SURFACE_ENERGY=WALL_HAMAKER_CONSTANT/&
206                 (24.d0*Pi*WALL_VDW_INNER_CUTOFF**2)
207     
208           ENDIF
209     
210     
211           CALL FINL_ERR_MSG
212     
213           RETURN
214     
215      1001 FORMAT('Error 1001: Illegal or unknown input: ',A, ' = ',A,/     &
216                  'Please correct the mfix.dat file.')
217     
218      1002 FORMAT('Error 1000: Cannot use SQUARE_WELL and VAN_DER_WAALS ',  &
219              'cohesion',/'models simultaneously.')
220     
221      1003 FORMAT('Error 1001: A cohesion model was not selected. Specify ',&
222              'one of the available models in the mfix.dat file.')
223     
224     
225     !<------------------- Van der Waals model messages. ----------------->!
226     
227      1201 FORMAT('Error 1201: Missing input data for Van der Waals ',      &
228              'cohesion model.',/'Input parameter ',A,' is UNDEFINED.')
229     
230      1202 FORMAT('Error 1202: VDW_OUTER_CUTOFF outside of the neighbor ',  &
231              'search distance.',/'Increase FACTOR_RLM to increase the ',   &
232              'search distance.')
233     
234           END SUBROUTINE CHECK_SOLIDS_DEM_COHESION
235     
236     
237     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
238     !                                                                      !
239     !  Subroutine: CHECK_SOLIDS_DEM_COLLISION                              !
240     !  Author: J.Musser                                   Date: 11-Dec-13  !
241     !                                                                      !
242     !  Purpose: Check user input data for DES collision calculations.      !
243     !                                                                      !
244     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
245           SUBROUTINE CHECK_SOLIDS_DEM_COLLISION
246     
247     
248     ! Global Variables:
249     !---------------------------------------------------------------------//
250     
251     ! User specified collision model
252           USE discretelement, only: DES_COLL_MODEL
253           USE discretelement, only: DES_COLL_MODEL_ENUM
254           USE discretelement, only: LSD
255           USE discretelement, only: HERTZIAN
256     ! Particle and wall friction coeff.
257           USE discretelement, only: MEW, MEW_W
258     ! Parameter constatns.
259           USE param1, only: ONE, ZERO, UNDEFINED
260     
261     !      USE mpi_utility
262           use error_manager
263     
264           IMPLICIT NONE
265     
266     ! Initialize the error manager.
267           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM_COLLISION")
268     
269     ! Check coefficient friction
270           IF(MEW == UNDEFINED) THEN
271              WRITE(ERR_MSG,1000) 'MEW'
272              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
273           ELSEIF (MEW < ZERO .OR. MEW_W > ONE) THEN
274              WRITE(ERR_MSG,1001) 'MEW', trim(iVal(MEW))
275              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
276           ENDIF
277     
278           IF(MEW_W == UNDEFINED) THEN
279              WRITE(ERR_MSG,1000) 'MEW_W'
280              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
281           ELSEIF(MEW_w < ZERO .OR. MEW_W > ONE) THEN
282              WRITE(ERR_MSG,1001) 'MEW_W', trim(iVal(MEW_W))
283              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
284           ENDIF
285     
286     ! Check collision model specific parameters.
287           SELECT CASE (trim(DES_COLL_MODEL))
288     ! Linear spring-dashpot model.
289           CASE('LSD')
290             DES_COLL_MODEL_ENUM = LSD
291             CALL CHECK_SOLIDS_DEM_COLL_LSD
292     ! Hertzian collision model.
293           CASE('HERTZIAN')
294              DES_COLL_MODEL_ENUM = HERTZIAN
295              CALL CHECK_SOLIDS_DEM_COLL_HERTZ
296     ! Unknown collision model.
297           CASE DEFAULT
298              WRITE(ERR_MSG,2000) TRIM(DES_COLL_MODEL)
299              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
300           END SELECT
301     
302      2000 FORMAT('Error 2000: Invalid particle-particle collision model:',&
303              A,/'Please correct the mfix.dat file.')
304     
305     
306           CALL FINL_ERR_MSG
307     
308           RETURN
309     
310      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
311              'correct the mfix.dat file.')
312     
313      1001 FORMAT('Error 1001: Illegal or unknown input: ',A,' = ',A,/      &
314              'Please correct the mfix.dat file.')
315     
316           END SUBROUTINE CHECK_SOLIDS_DEM_COLLISION
317     
318     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
319     !                                                                      !
320     !  Subroutine: CHECK_SOLIDS_DEM_COLL_LSD                               !
321     !  Author: J.Musser                                   Date: 11-Dec-13  !
322     !                                                                      !
323     !  Purpose: Check user input data for DES collision calculations.      !
324     !                                                                      !
325     !  References:                                                         !
326     !   - Schafer et al., J. Phys. I France, 1996, 6, 5-20 (see page 7&13) !
327     !   -  Van der Hoef et al., Advances in Chemical Engineering, 2006, 31,!
328     !      65-149 (pages 94-95)                                            !
329     !   - Silbert et al., Physical Review E, 2001, 64, 051302 1-14 (page 5)!
330     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
331           SUBROUTINE CHECK_SOLIDS_DEM_COLL_LSD
332     
333     
334     ! Global Variables:
335     !---------------------------------------------------------------------//
336     ! Number of discrete solids phases
337           USE discretelement, only: DES_MMAX, DES_D_P0, DES_RO_s
338     ! Particle and wall normal and tangential spring constants
339           USE discretelement, only: KN, KN_W
340           USE discretelement, only: KT, KT_W
341     ! Particle and wall tangential spring factor := KN/KT
342           USE discretelement, only: KT_FAC, KT_W_FAC
343     
344           use discretelement, only: DES_ETAN, DES_ETAN_WALL
345           use discretelement, only: DES_ETAT, DES_ETAT_WALL
346     
347     ! Coefficients of restitution: Normal and Tangential
348           USE discretelement, only: DES_EN_INPUT, DES_EN_WALL_INPUT
349           USE discretelement, only: DES_ET_INPUT, DES_ET_WALL_INPUT
350     
351     ! Tangential damping factors := ET/EN
352           USE discretelement, only: DES_ETAT_FAC, DES_ETAT_W_FAC
353     
354           use constant, only: PI
355     
356           use discretelement, only: DTSOLID
357     
358     ! Parameter constatns.
359           USE param1, only: ZERO, HALF, ONE, UNDEFINED
360     
361     !      USE mpi_utility
362           use error_manager
363     
364           IMPLICIT NONE
365     
366     
367     ! Local Variables:
368     !---------------------------------------------------------------------//
369     ! Loop index.
370           INTEGER :: M, L, LC
371     ! Flag to warn user.
372           LOGICAL :: FLAG_WARN
373     ! Collision length scale.
374           DOUBLE PRECISION :: TCOLL, TCOLL_TMP
375     ! Collision length scale.
376           DOUBLE PRECISION :: MASS_M, MASS_L, MASS_EFF
377     ! Alias for coefficient restitution
378           DOUBLE PRECISION :: EN
379     
380     ! Initialize the error manager.
381           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM_COLL_LSD")
382     
383     ! Initialize.
384           TCOLL = UNDEFINED
385     
386     ! Check for particle-particle normal spring constants.
387           IF(KN == UNDEFINED) THEN
388              WRITE(ERR_MSG, 1000) 'KN'
389              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
390           ENDIF
391     
392     ! Check for particle-particle tangential spring constant factors.
393           IF(KT_FAC == UNDEFINED) THEN
394              WRITE (ERR_MSG, 2100) 'KT_FAC'
395              CALL FLUSH_ERR_MSG()
396              KT_FAC = 2.0d0/7.0d0
397           ELSEIF(KT_FAC > ONE .OR. KT_FAC < ZERO) THEN
398              WRITE(ERR_MSG,1001) 'KT_FAC', trim(iVal(KT_FAC))
399              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
400           ENDIF
401     ! Calculate the particle-particle tangential spring factor.
402           KT = KT_FAC*KN
403     
404     ! Check for particle-wall normal spring constants.
405           IF(KN_W == UNDEFINED) THEN
406              WRITE(ERR_MSG, 1000) 'KN_W'
407              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
408           ENDIF
409     
410     ! Check for particle-wall tangential spring constant factors.
411           IF(KT_W_FAC == UNDEFINED) THEN
412              WRITE (ERR_MSG, 2100) 'KT_W_FAC'
413              CALL FLUSH_ERR_MSG()
414              KT_W_FAC = 2.0d0/7.0d0
415           ELSEIF(KT_W_FAC > ONE .OR. KT_W_FAC < ZERO) THEN
416              WRITE(ERR_MSG,1001) 'KT_W_FAC', trim(iVal(KT_W_FAC))
417              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
418           ENDIF
419     ! Calculate the particle-wall tangential spring factor.
420           KT_W = KT_W_FAC*KN_W
421     
422      2100 FORMAT('Warning 2100: Tangential spring factor ',A,' not ',      &
423              'specified in mfix.dat.',/'Setting to default: (2/7).')
424     
425     ! Check for particle-particle tangential damping coefficients
426           IF(DES_ETAT_FAC == UNDEFINED) THEN
427              WRITE (ERR_MSG, 2101) 'DES_ETAT_FAC'
428              CALL FLUSH_ERR_MSG
429              DES_ETAT_FAC = HALF
430           ELSEIF(DES_ETAT_FAC > ONE .OR. DES_ETAT_FAC < ZERO) THEN
431              WRITE(ERR_MSG,1001) 'DES_ETAT_FAC', iVal(DES_ETAT_FAC)
432              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
433           ENDIF
434     
435     ! Check for particle-wall tangential damping coefficients
436           IF(DES_ETAT_W_FAC == UNDEFINED) THEN
437              WRITE (ERR_MSG, 2101) 'DES_ETAT_W_FAC'
438              CALL FLUSH_ERR_MSG
439              DES_ETAT_W_FAC = HALF
440           ELSEIF(DES_ETAT_W_FAC > ONE .OR. DES_ETAT_W_FAC < ZERO) THEN
441              WRITE(ERR_MSG,1001) 'DES_ETAT_W_FAC', iVal(DES_ETAT_W_FAC)
442              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
443           ENDIF
444     
445      2101 FORMAT('Warning 2101: Tangential damping factor ',A,' not ',     &
446              'specified',/'in mfix.dat. Setting to default: (1/2).')
447     
448     
449           LC = 0
450           DO M = 1, DES_MMAX
451     
452     ! Calculate the mass of a phase M particle.
453              MASS_M = (PI/6.d0)*(DES_D_P0(M)**3)*DES_RO_S(M)
454     
455     ! Particle-Particle Collision Parameters ------------------------------>
456              DO L = M, DES_MMAX
457                 LC = LC+1
458     
459     ! Check particle-particle normal restitution coefficient
460                 IF(DES_EN_INPUT(LC) == UNDEFINED) THEN
461                    WRITE(ERR_MSG,1000) trim(iVar('DES_EN_INPUT',LC))
462                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
463                 ELSEIF(DES_EN_INPUT(LC) > ONE .OR.                         &
464                    DES_EN_INPUT(LC) < ZERO) THEN
465                    WRITE(ERR_MSG,1001) trim(iVar('DES_EN_INPUT',LC)),      &
466                       trim(iVal(DES_EN_INPUT(LC)))
467                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
468                 ENDIF
469                 EN = DES_EN_INPUT(LC)
470     
471     ! Calculate masses used for collision calculations.
472                 MASS_L = (PI/6.d0)*(DES_D_P0(L)**3)*DES_RO_S(L)
473                 MASS_EFF = MASS_M*MASS_L/(MASS_M+MASS_L)
474     
475     ! Calculate the M-L normal and tangential damping coefficients.
476                 IF(EN .NE. ZERO) THEN
477                    DES_ETAN(M,L) = 2.0D0*SQRT(KN*MASS_EFF) * ABS(LOG(EN))
478                    DES_ETAN(M,L) = DES_ETAN(M,L)/SQRT(PI*PI + (LOG(EN)**2))
479                 ELSE
480                    DES_ETAN(M,L) = 2.0D0*SQRT(KN*MASS_EFF)
481                 ENDIF
482                 DES_ETAT(M,L) = DES_ETAT_FAC*DES_ETAN(M,L)
483     
484     ! Store the entries in the symmetric matrix.
485                 DES_ETAN(L,M) = DES_ETAN(M,L)
486                 DES_ETAT(L,M) = DES_ETAT(M,L)
487     
488     ! Calculate the collision time scale.
489                 TCOLL_TMP = PI/SQRT(KN/MASS_EFF -                          &
490                    ((DES_ETAN(M,L)/MASS_EFF)**2)/4.d0)
491                 TCOLL = MIN(TCOLL_TMP, TCOLL)
492              ENDDO
493     
494     
495     ! Particle-Wall Collision Parameters ---------------------------------->
496     ! Check particle-wall normal restitution coefficient.
497              IF(DES_EN_WALL_INPUT(M) == UNDEFINED) THEN
498                 WRITE(ERR_MSG,1000) trim(iVar('DES_EN_WALL_INPUT',M))
499                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
500              ELSEIF(DES_EN_WALL_INPUT(M) > ONE .OR.                        &
501                 DES_EN_WALL_INPUT(M) < ZERO) THEN
502                 WRITE(ERR_MSG,1001) trim(iVar('DES_EN_WALL_INPUT',M)),     &
503                    trim(iVal(DES_EN_WALL_INPUT(M)))
504                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
505              ENDIF
506              EN = DES_EN_WALL_INPUT(M)
507     
508     ! Calculate masses used for collision calculations.
509              MASS_EFF = MASS_M
510     
511     ! Calculate the M-Wall normal and tangential damping coefficients.
512              IF(EN .NE. ZERO) THEN
513                 DES_ETAN_WALL(M) = 2.d0*SQRT(KN_W*MASS_EFF)*ABS(LOG(EN))
514                 DES_ETAN_WALL(M) = DES_ETAN_WALL(M)/SQRT(PI*PI+(LOG(EN))**2)
515              ELSE
516                 DES_ETAN_WALL(M) = 2.D0*SQRT(KN_W*MASS_EFF)
517              ENDIF
518              DES_ETAT_WALL(M) = DES_ETAT_W_FAC*DES_ETAN_WALL(M)
519     
520     ! Calculate the collision time scale.
521              TCOLL_TMP = PI/SQRT(KN_W/MASS_EFF -                           &
522                 ((DES_ETAN_WALL(M)/MASS_EFF)**2.d0)/4.d0)
523     !         TCOLL = MIN(TCOLL_TMP, TCOLL)
524           ENDDO
525     
526     ! if following are assigned warn user they are discarded
527           FLAG_WARN = .FALSE.
528           DO M = 1, DES_MMAX+DES_MMAX*(DES_MMAX-1)/2
529              IF(DES_ET_INPUT(M) .NE. UNDEFINED) FLAG_WARN = .TRUE.
530           ENDDO
531           IF (FLAG_WARN) THEN
532              WRITE(ERR_MSG,2102) 'DES_ET_INPUT'
533              CALL FLUSH_ERR_MSG
534           ENDIF
535     
536           FLAG_WARN = .FALSE.
537           DO M = 1, DES_MMAX
538              IF(DES_ET_WALL_INPUT(M) .NE. UNDEFINED) FLAG_WARN = .TRUE.
539           ENDDO
540           IF (FLAG_WARN)THEN
541              WRITE(ERR_MSG,2102) 'DES_ET_WALL_INPUT'
542              CALL FLUSH_ERR_MSG
543           ENDIF
544     
545      2102 FORMAT('Warning 2102: ',A,' values are not used ',/' with the',  &
546              ' linear spring-dashpot collision model.')
547     
548     ! Store the smalled calculated collision time scale. This value is used
549     ! in time-marching the DEM solids.
550           DTSOLID = TCOLL/50.d0
551     
552     
553           CALL FINL_ERR_MSG
554     
555           RETURN
556     
557      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
558              'correct the mfix.dat file.')
559     
560      1001 FORMAT('Error 1001: Illegal or unknown input: ',A,' = ',A,/      &
561              'Please correct the mfix.dat file.')
562     
563           END SUBROUTINE CHECK_SOLIDS_DEM_COLL_LSD
564     
565     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
566     !                                                                      !
567     !  Subroutine: CHECK_SOLIDS_DEM_COLL_HERTZ                             !
568     !  Author: J.Musser                                   Date: 11-Dec-13  !
569     !                                                                      !
570     !  Purpose: Check user input data for Hertzian collisions.             !
571     !                                                                      !
572     !  References:                                                         !
573     !   - Schafer et al., J. Phys. I France, 1996, 6, 5-20 (see page 7&13) !
574     !   -  Van der Hoef et al., Advances in Chemical Engineering, 2006, 31,!
575     !      65-149 (pages 94-95)                                            !
576     !   - Silbert et al., Physical Review E, 2001, 64, 051302 1-14 (page 5)!
577     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
578           SUBROUTINE CHECK_SOLIDS_DEM_COLL_HERTZ
579     
580     ! Global Variables:
581     !---------------------------------------------------------------------//
582     ! Number of discrete solids phases, diameters and densities
583           USE discretelement, only: DES_MMAX, DES_D_P0, DES_RO_s
584     ! User defined coefficients of restitution: Normal and Tangential
585           USE discretelement, only: DES_EN_INPUT, DES_EN_WALL_INPUT
586           USE discretelement, only: DES_ET_INPUT, DES_ET_WALL_INPUT
587     ! Damping coefficients: Normal and Tangential
588           use discretelement, only: DES_ETAN, DES_ETAN_WALL
589           use discretelement, only: DES_ETAT, DES_ETAT_WALL
590     ! Hertzian spring constants: Normal and Tangential
591           use discretelement, only: HERT_KN, HERT_Kwn
592           use discretelement, only: HERT_KT, HERT_Kwt
593     ! Tangential damping factors := ET/EN
594           USE discretelement, only: DES_ETAT_FAC, DES_ETAT_W_FAC
595     ! Particle and wall Young's modulus and Shear modulus
596           USE discretelement, only: E_YOUNG, Ew_YOUNG, G_MOD
597     ! Particle and wall Poisson ratio
598           USE discretelement, only: V_POISSON, Vw_POISSON
599     ! Solids time step-size.
600           use discretelement, only: DTSOLID
601     ! Particle and wall normal spring constants
602           USE discretelement, only: KN, KN_W
603     ! Particle and wall tangential spring factor := KN/KT
604           USE discretelement, only: KT_FAC, KT_W_FAC
605     ! The constant PI
606           use constant, only: PI
607     
608     ! Parameter constatns.
609           USE param1, only: ZERO, ONE, UNDEFINED
610     
611     !      USE mpi_utility
612           use error_manager
613     
614           IMPLICIT NONE
615     
616     ! Local Variables:
617     !---------------------------------------------------------------------//
618     ! Loop index.
619           INTEGER :: M, L, LC
620     ! Message for formatted output.
621           CHARACTER(len=64) :: MSG
622     ! Collision length scale.
623           DOUBLE PRECISION :: TCOLL, TCOLL_TMP
624     ! Particle and effective mass.
625           DOUBLE PRECISION :: MASS_M, MASS_L, MASS_EFF
626     ! Effective physical quantities. Radius, Youngs, Shear
627           DOUBLE PRECISION :: R_EFF, E_EFF, G_MOD_EFF, RED_MASS_EFF
628     ! Alias for coefficient restitution
629           DOUBLE PRECISION :: EN, ET
630     ! Shear modules for wall
631           DOUBLE PRECISION :: G_MOD_WALL
632     
633     ! Initialize the error manager.
634           CALL INIT_ERR_MSG("CHECK_SOLIDS_DEM_COLL_HERTZ")
635     
636     ! Initialize.
637           TCOLL = UNDEFINED
638     
639     ! check young's modulus and poisson ratio
640           IF(Ew_YOUNG == UNDEFINED ) THEN
641              MSG='Wall value for Youngs modulus'
642              WRITE(ERR_MSG,1002) 'Ew_YOUNG', MSG
643              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
644           ENDIF
645     
646           IF(Vw_POISSON == UNDEFINED) THEN
647              MSG='Wall value for Poissons ratio'
648              WRITE(ERR_MSG,1002) 'Vw_POISSON', MSG
649              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
650           ELSEIF (Vw_POISSON > 0.5d0 .OR. Vw_POISSON <= -ONE) THEN
651              WRITE(ERR_MSG,1001) 'Vw_POISSON',iVal(Vw_POISSON)
652              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
653           ENDIF
654     
655           G_MOD_WALL = 0.5d0*Ew_YOUNG/(1.d0+Vw_POISSON)
656     
657           DO M=1,DES_MMAX
658     
659              IF(E_YOUNG(M) == UNDEFINED) THEN
660                 MSG=''; WRITE(MSG,"('Phase ',I2,' Youngs modulus')") M
661                 WRITE(ERR_MSG,1002) 'E_YOUNG', MSG
662                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
663              ENDIF
664              IF(V_POISSON(M) == UNDEFINED) THEN
665                 MSG=''; WRITE(MSG,"('Phase ',I2,' Poissons ratio')") M
666                 WRITE(ERR_MSG,1002) 'V_POISSON', MSG
667                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
668              ELSEIF(V_POISSON(M) > 0.5d0 .OR.                              &
669                 V_POISSON(M) <= -ONE) THEN
670                 WRITE(ERR_MSG,1001) trim(iVar('V_POISSON',M)),             &
671                    iVal(V_POISSON(M))
672                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
673              ENDIF
674     ! Calculate the shear modulus for phase M.
675              G_MOD(M) = 0.5d0*E_YOUNG(M)/(1.d0+V_POISSON(M))
676           ENDDO
677     
678     
679           LC = 0
680           DO M=1,DES_MMAX
681     ! Calculate the mass of a phase M particle.
682              MASS_M = (PI/6.d0)*(DES_D_P0(M)**3)*DES_RO_S(M)
683     
684     ! Particle-Particle Collision Parameters ------------------------------>
685              DO L=M,DES_MMAX
686                 LC = LC+1
687     
688     ! Check particle-particle normal restitution coefficient
689                 IF(DES_EN_INPUT(LC) == UNDEFINED) THEN
690                    WRITE(ERR_MSG,1000) trim(iVar('DES_EN_INPUT',LC))
691                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
692     
693                 ELSEIF(DES_EN_INPUT(LC) > ONE .OR.                         &
694                    DES_EN_INPUT(LC) < ZERO) THEN
695                    WRITE(ERR_MSG,1001) trim(iVar('DES_EN_INPUT',LC)),      &
696                       trim(iVal(DES_EN_INPUT(LC)))
697                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
698                 ENDIF
699                 EN = DES_EN_INPUT(LC)
700     
701     ! Check particle-particle tangential restitution coefficient
702                 IF(DES_ET_INPUT(M) == UNDEFINED) THEN
703                    WRITE(ERR_MSG,1000) trim(iVar('DES_ET_INPUT',M))
704                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
705                 ELSEIF(DES_ET_INPUT(M) > ONE .OR.                          &
706                    DES_ET_INPUT(M) < ZERO) THEN
707                    WRITE(ERR_MSG,1001) trim(iVar('DES_ET_INPUT',M)),       &
708                       iVal(DES_ET_INPUT(M))
709                    CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
710                 ENDIF
711                 ET = DES_ET_INPUT(LC)
712     
713     
714     ! Calculate masses used for collision calculations.
715                 MASS_L = (PI/6.d0)*(DES_D_P0(L)**3)*DES_RO_S(L)
716                 MASS_EFF = (MASS_M*MASS_L)/(MASS_M+MASS_L)
717                 RED_MASS_EFF = (2.d0/7.d0)*MASS_EFF
718     ! Calculate the effective radius, Youngs modulus, and shear modulus.
719                 R_EFF = 0.5d0*(DES_D_P0(M)*DES_D_P0(L)/                    &
720                    (DES_D_P0(M) + DES_D_P0(L)))
721                 E_EFF = E_YOUNG(M)*E_YOUNG(L) /                            &
722                    (E_YOUNG(M)*(1.d0 - V_POISSON(L)**2) +                  &
723                     E_YOUNG(L)*(1.d0 - V_POISSON(M)**2))
724                 G_MOD_EFF = G_MOD(M)*G_MOD(L)/                             &
725                    (G_MOD(M)*(2.d0 - V_POISSON(L)) +                       &
726                     G_MOD(L)*(2.d0 - V_POISSON(M)))
727     
728     ! Calculate the spring properties and store in symmetric matrix format.
729                 HERT_KN(M,L)=(4.d0/3.d0)*SQRT(R_EFF)*E_EFF
730                 HERT_KT(M,L)= 8.d0*SQRT(R_EFF)*G_MOD_EFF
731     
732                 HERT_KN(L,M) = HERT_KN(M,L)
733                 HERT_KT(L,M) = HERT_KT(M,L)
734     
735     ! Calculate the normal coefficient.
736                 IF(EN .NE. ZERO) THEN
737                    DES_ETAN(M,L) = 2.d0*SQRT(HERT_KN(M,L)*MASS_EFF)*       &
738                       ABS(LOG(EN))
739                    DES_ETAN(M,L) = DES_ETAN(M,L)/                          &
740                       SQRT(PI*PI + (LOG(EN))**2)
741                 ELSE
742                    DES_ETAN(M,L) = 2.d0*SQRT(HERT_KN(M,L)*MASS_EFF)
743                 ENDIF
744                 DES_ETAN(L,M) = DES_ETAN(M,L)
745     
746     ! Calculate the tangential coefficients.
747                 IF(ET .NE. ZERO) THEN
748                    DES_ETAT(M,L) = 2.d0*SQRT(HERT_KT(M,L)*RED_MASS_EFF)*   &
749                       ABS(LOG(ET))
750                    DES_ETAT(M,L) = DES_ETAT(M,L)/ SQRT(PI*PI + (LOG(ET))**2)
751                 ELSE
752                    DES_ETAT(M,L) = 2.d0*SQRT(HERT_KT(M,L)*RED_MASS_EFF)
753                 ENDIF
754                 DES_ETAT(L,M) = DES_ETAT(M,L)
755     
756                 TCOLL_TMP = PI/SQRT(HERT_KN(M,L)/MASS_EFF -                &
757                    ((DES_ETAN(M,L)/MASS_EFF)**2)/4.d0)
758                 TCOLL = MIN(TCOLL_TMP, TCOLL)
759              ENDDO
760     
761     ! Particle-Wall Collision Parameters ---------------------------------->
762     ! Check particle-wall normal restitution coefficient.
763              IF(DES_EN_WALL_INPUT(M) == UNDEFINED) THEN
764                 WRITE(ERR_MSG,1000) trim(iVar('DES_EN_WALL_INPUT',M))
765                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
766              ELSEIF(DES_EN_WALL_INPUT(M) > ONE .OR.                        &
767                 DES_EN_WALL_INPUT(M) < ZERO) THEN
768                 WRITE(ERR_MSG,1001) trim(iVar('DES_EN_WALL_INPUT',M)),     &
769                    trim(iVal(DES_EN_WALL_INPUT(M)))
770                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
771              ENDIF
772              EN = DES_EN_WALL_INPUT(M)
773     
774     ! Check particle-wall tangential restitution coefficient
775              IF(DES_ET_WALL_INPUT(M) == UNDEFINED) THEN
776                 WRITE(ERR_MSG,1000) trim(iVar('DES_ET_WALL_INPUT',M))
777                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
778              ELSEIF(DES_ET_WALL_INPUT(M) > ONE .OR.                        &
779                 DES_ET_WALL_INPUT(M) < ZERO) THEN
780                 WRITE(ERR_MSG,1001) trim(iVar('DES_ET_WALL_INPUT',M)),     &
781                    trim(iVal(DES_ET_WALL_INPUT(M)))
782                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
783              ENDIF
784              ET = DES_ET_WALL_INPUT(M)
785     
786     ! Calculate masses used for collision calculations.
787              MASS_EFF = MASS_M
788              RED_MASS_EFF = (2.d0/7.d0)*MASS_EFF
789     ! Calculate the effective radius, Youngs modulus, and shear modulus.
790              R_EFF = 0.5d0*DES_D_P0(M)
791              E_EFF = E_YOUNG(M)*Ew_YOUNG /                                 &
792                 (E_YOUNG(M)*(1.d0-Vw_POISSON**2) +                         &
793                  Ew_YOUNG  *(1.d0-V_POISSON(M)**2))
794              G_MOD_EFF = G_MOD(M)*G_MOD_WALL /                             &
795                 (G_MOD(M)*(2.d0 - Vw_POISSON) +                            &
796                  G_MOD_WALL*(2.d0 - V_POISSON(M)))
797     
798     ! Calculate the spring properties.
799              HERT_Kwn(M) = (4.d0/3.d0)*SQRT(R_EFF)*E_EFF
800              HERT_Kwt(M) = 8.0*SQRT(R_EFF)*G_MOD_EFF
801     
802     ! Calculate the tangential coefficients.
803              IF(EN /= ZERO) THEN
804                 DES_ETAN_WALL(M) = 2.d0*SQRT(HERT_Kwn(M)*MASS_EFF)*&
805                    ABS(LOG(EN))
806                 DES_ETAN_WALL(M) = DES_ETAN_WALL(M)/&
807                    SQRT(PI*PI + (LOG(EN))**2)
808              ELSE
809                 DES_ETAN_WALL(M) = 2.d0*SQRT(HERT_Kwn(M)*MASS_EFF)
810              ENDIF
811     
812              IF(ET /= ZERO) THEN
813                 DES_ETAT_WALL(M) = 2.d0*SQRT(HERT_Kwt(M)*RED_MASS_EFF)*    &
814                     ABS(LOG(ET))
815                 DES_ETAT_WALL(M) = DES_ETAT_WALL(M)/SQRT(PI*PI+(LOG(ET))**2)
816              ELSE
817                 DES_ETAT_WALL(M) = 2.d0*SQRT(HERT_Kwt(M)*RED_MASS_EFF)
818              ENDIF
819     
820     ! Calculate the collision time scale.
821              TCOLL_TMP = PI/SQRT(HERT_Kwn(M)/MASS_EFF -                    &
822                 ((DES_ETAN_WALL(M)/MASS_EFF)**2)/4.d0)
823           ENDDO
824     
825     
826     ! If following are assigned warn user they are discarded.
827            IF(KN .NE. UNDEFINED) THEN
828               WRITE(ERR_MSG, 2200) 'KN'
829               CALL FLUSH_ERR_MSG
830            ENDIF
831            IF(KN_W .NE. UNDEFINED) THEN
832               WRITE(ERR_MSG, 2200) 'KN_W'
833               CALL FLUSH_ERR_MSG
834            ENDIF
835            IF(KT_FAC .NE. UNDEFINED) THEN
836               WRITE(ERR_MSG, 2200) 'KT_FAC'
837               CALL FLUSH_ERR_MSG
838            ENDIF
839            IF(KT_W_FAC .NE. UNDEFINED) THEN
840               WRITE(ERR_MSG, 2200) 'KT_W_FAC'
841               CALL FLUSH_ERR_MSG
842            ENDIF
843            IF(DES_ETAT_FAC .NE. UNDEFINED) THEN
844               WRITE(ERR_MSG, 2200) 'DES_ETAT_FAC'
845               CALL FLUSH_ERR_MSG
846            ENDIF
847            IF(DES_ETAT_W_FAC .NE. UNDEFINED) THEN
848               WRITE(ERR_MSG, 2200) 'DES_ETAT_W_FAC'
849               CALL FLUSH_ERR_MSG
850            ENDIF
851     
852      2200 FORMAT('Warning 2200: ',A,' values are not used ',/' with the',  &
853              ' linear spring-dashpot collision model.')
854     
855     
856     ! Store the smalled calculated collision time scale. This value is used
857     ! in time-marching the DEM solids.
858           DTSOLID = TCOLL/50.d0
859     
860     
861           CALL FINL_ERR_MSG
862     
863           RETURN
864     
865      1000 FORMAT('Error 1000: Required input not specified: ',A,/'Please ',&
866              'correct the mfix.dat file.')
867     
868      1001 FORMAT('Error 1001: Illegal or unknown input: ',A,' = ',A,/      &
869              'Please correct the mfix.dat file.')
870     
871      1002 FORMAT('Error 1002: Required input not specified: ',A,/          &
872              'Description:',A,/'Please correct the mfix.dat file.')
873     
874           END SUBROUTINE CHECK_SOLIDS_DEM_COLL_HERTZ
875     
876     
877     
878     
879     
880     
881     
882     
883     
884     
885