1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv 2 ! 3 ! Subroutine: CFRELVEL 4 ! Purpose: Calculate the normal and tangential components of the 5 ! relative velocity between contacting particles 6 ! 7 ! Author: Jay Boyalakuntla Date: 12-Jun-04 8 ! Reviewer: Rahul Garg Date: 01-Aug-07 9 ! 10 ! Comments: Relative (translational) velocity required for eqn 6 11 ! from the following paper: 12 ! Tsuji Y., Kawaguchi T., and Tanak T., "Lagrangian numerical 13 ! simulation of plug glow of cohesionless particles in a 14 ! horizontal pipe", Powder technology, 71, 239-250, 1992 15 ! 16 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 17 18 SUBROUTINE CFRELVEL(L, II, VRN, VSLIP, NORM, DIST_LI) 19 20 !----------------------------------------------- 21 ! Modules 22 !----------------------------------------------- 23 USE discretelement, only: DES_VEL_NEW, DES_RADIUS, OMEGA_NEW, DES_CROSSPRDCT 24 IMPLICIT NONE 25 !----------------------------------------------- 26 ! Dummy arguments 27 !----------------------------------------------- 28 ! indices of particle-particle contact pair 29 INTEGER, INTENT(IN) :: L, II 30 ! distance between particle centers 31 DOUBLE PRECISION, INTENT(IN) :: DIST_LI 32 ! unit normal vector along the line of contact pointing from 33 ! particle L to particle II 34 DOUBLE PRECISION, INTENT(IN) :: NORM(3) 35 ! slip velocity at point of contact 36 DOUBLE PRECISION, INTENT(OUT) :: VSLIP(3) 37 ! normal component of relative contact velocity (scalar) 38 DOUBLE PRECISION, INTENT(OUT) :: VRN 39 !----------------------------------------------- 40 ! Local variables 41 !----------------------------------------------- 42 ! translational relative velocity 43 DOUBLE PRECISION :: VRELTRANS(3) 44 ! rotational velocity at point of contact 45 DOUBLE PRECISION :: V_ROT(3), OMEGA_SUM(3) 46 ! distance from the contact point to the particle centers 47 DOUBLE PRECISION :: DIST_CL, DIST_CI 48 !----------------------------------------------- 49 50 ! translational relative velocity 51 VRELTRANS(:) = (DES_VEL_NEW(L,:) - DES_VEL_NEW(II,:)) 52 53 ! calculate the distance from the particle center to the contact point, 54 ! which is taken as the radical line 55 ! dist_ci+dist_cl=dist_li; dist_ci^2+a^2=ri^2; dist_cl^2+a^2=rl^2 56 DIST_CL = (DIST_LI**2 + DES_RADIUS(L)**2 - DES_RADIUS(II)**2)/& 57 (2.d0*DIST_LI) 58 DIST_CI = DIST_LI - DIST_CL 59 60 OMEGA_SUM(:) = OMEGA_NEW(L,:)*DIST_CL + & 61 OMEGA_NEW(II,:)*DIST_CI 62 63 ! calculate the rotational relative velocity 64 V_ROT = DES_CROSSPRDCT(OMEGA_SUM, NORM) 65 66 ! total relative velocity 67 VRELTRANS(:) = VRELTRANS(:) + V_ROT(:) 68 69 ! normal component of relative velocity (scalar) 70 VRN = DOT_PRODUCT(VRELTRANS,NORM) 71 72 ! slip velocity of the contact point 73 ! Equation (8) in Tsuji et al. 1992 74 VSLIP(:) = VRELTRANS(:) - VRN*NORM(:) 75 76 RETURN 77 END SUBROUTINE CFRELVEL 78