#include "cfx5ext.h"
dllexport(pt_injection)
      SUBROUTINE PT_INJECTION(NLOC,NRET,NARG,RET,ARG,CRESLT,
     &                        CZ,DZ,IZ,LZ,RZ)
C
C=======================================================================
C User routine: User defined particle injection
C=======================================================================
C
C ------------------------------
C        Argument list
C ------------------------------
C
      INTEGER NARG, NRET, NLOC
C
      REAL ARG(NLOC,NARG), RET(NLOC,NRET)
C
      CHARACTER*(4) CRESLT
C
      INTEGER IZ(*)
      CHARACTER CZ(*)*(1)
      DOUBLE PRECISION DZ(*)
      LOGICAL LZ(*)
      REAL RZ(*)
C
C=======================================================================
C
C     Return variables:
C     -----------------
C
C     Particle position             : CRD_PT = RET(1:NLOC,1:3)
C     Particle diameter             : DIAM_PT = RET(1:NLOC,4)
C     Particle velocity             : VEL_PT = RET(1:NLOC,5:7)
C     Particle temperature          : TEMP_PT = RET(1:NLOC,8)
C     Particle mass flow rate       : MRATE_PT = RET(1:NLOC,9)
C=======================================================================
C
C-----------------------------------------------------------------------
C     Calculate the particle quantities
C-----------------------------------------------------------------------
C
      CALL MYINJECTION_SUB (NLOC,RET(1,1),RET(1,4),RET(1,5),RET(1,8),
     &                      RET(1,9))
C
      END

      SUBROUTINE MYINJECTION_SUB (NPART,CRD_PT,DIAM_PT,VEL_PT,TEMP_PT,
     &                            MRATE_PT)
C
C=======================================================================
C     Calculate the particle positions, velocities and diameters
C=======================================================================
C
C ------------------------------
C     Preprocessor includes
C ------------------------------
C
#include "cfd_sysdep.h"
#include "cfd_constants.h"
C
C ------------------------------
C        Argument list
C ------------------------------
C
      INTEGER            NPART
      REAL               CRD_PT(3,NPART), VEL_PT(3,NPART), 
     &                   DIAM_PT(NPART), TEMP_PT(NPART), MRATE_PT(NPART)
C
C ------------------------------
C        Local variables
C ------------------------------
C
      INTEGER            ISEED, IPART, ILOC, NPART_INJECT, NPART_DIM,
     &                   IPART_S, IPART_F
      REAL               CRD_POS(3,4), TEMP(4), CONE_ANGLE(4), 
     &                   VEL_MAG(3), MRATE(4), DIAM(4),
     &                   DALP, BETA, VEL_Z, VEL_XY, ALP
C
C ------------------------------
C    Executable statements
C ------------------------------
C
C-----------------------------------------------------------------------
C     Inject particles at four different positions
C-----------------------------------------------------------------------
C
C---- Number of particles per injection location
C
      NPART_INJECT = NPART/4
C
C---- 1. position
C
      ILOC = 1
C
      CRD_POS(1,ILOC)  = 0.25
      CRD_POS(2,ILOC)  = 0.25
      CRD_POS(3,ILOC)  = 1.75
C
      TEMP(ILOC)       = 310.0
      CONE_ANGLE(ILOC) = 10.0
      VEL_MAG(ILOC)    = 2.1
      MRATE(ILOC)      = 1.1E-4
      DIAM(ILOC)       = 1.E-4
C
C---- 2. position
C
      ILOC = 2
C
      CRD_POS(1,ILOC)  = 0.25
      CRD_POS(2,ILOC)  = 0.75
      CRD_POS(3,ILOC)  = 1.75
C
      TEMP(ILOC)       = 320.0
      CONE_ANGLE(ILOC) = 15.0
      VEL_MAG(ILOC)    = 2.2
      MRATE(ILOC)      = 1.2E-4
      DIAM(ILOC)       = 2.E-4
C
C---- 3. position
C
      ILOC = 3
C
      CRD_POS(1,ILOC)  = 0.75
      CRD_POS(2,ILOC)  = 0.25
      CRD_POS(3,ILOC)  = 1.75
C
      TEMP(ILOC)       = 330.0
      CONE_ANGLE(ILOC) = 20.0
      VEL_MAG(ILOC)    = 2.3
      MRATE(ILOC)      = 1.3E-4
      DIAM(ILOC)       = 3.E-4
C
C---- 4. position
C
      ILOC = 4
C
      CRD_POS(1,ILOC)  = 0.75
      CRD_POS(2,ILOC)  = 0.75
      CRD_POS(3,ILOC)  = 1.75
C
      TEMP(ILOC)       = 340.0
      CONE_ANGLE(ILOC) = 25.0
      VEL_MAG(ILOC)    = 2.4
      MRATE(ILOC)      = 1.4E-4
      DIAM(ILOC)       = 4.E-4
C
C-----------------------------------------------------------------------
C     Loop over all four injection locations
C-----------------------------------------------------------------------
C
      DO ILOC=1,4
C
C---- Calculate first and last particle number
C
        IPART_S   = (ILOC-1)*NPART_INJECT + 1
        IPART_F   = MIN(NPART,ILOC*NPART_INJECT)
        NPART_DIM = IPART_F - IPART_S + 1
C
C---- Set the particle injection position, temperature, diameter and
C     mass flow rate
C
        DO IPART=IPART_S,IPART_F
C
          CRD_PT(1,IPART) = CRD_POS(1,ILOC)
          CRD_PT(2,IPART) = CRD_POS(2,ILOC)
          CRD_PT(3,IPART) = CRD_POS(3,ILOC)
C
          TEMP_PT(IPART)  = TEMP(ILOC)
C
          DIAM_PT(IPART)  = DIAM(ILOC)
C
          MRATE_PT(IPART) = MRATE(ILOC)
C
        ENDDO
C
C---- Set the particle injection velocities
C
        DALP   = 2.0*PI/FLOAT(NPART_DIM)
        BETA   = PI/180.*CONE_ANGLE(ILOC)
        VEL_Z  = VEL_MAG(ILOC)*COS(BETA)
        VEL_XY = VEL_MAG(ILOC)*SIN(BETA)
C
        DO IPART=IPART_S,IPART_F
C
          ALP = FLOAT(IPART)*DALP
C
          VEL_PT(1,IPART) = VEL_XY*COS(ALP)
          VEL_PT(2,IPART) = VEL_XY*SIN(ALP)
          VEL_PT(3,IPART) = - VEL_Z
C
        ENDDO
C
      ENDDO
C
      END
