GEANT321/gcons/ghmix.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1995/10/24 10:20:18  cernlib
   6 * Geant
   7 *
   8 *
   9 #include "geant321/pilot.h"
  10 *CMZ :  3.21/02 29/03/94  15.41.20  by  S.Giani
  11 *-- Author :
  12       SUBROUTINE GHMIX(A, W, N, AEFF)
  13 C
  14 C     ******************************************************************
  15 C     *                                                                *
  16 C     *   Works out an effective atomic weight AEFF for a material     *
  17 C     * with N elements of atomic weight A(I) in proportion W(I) by    *
  18 C     * weight.  The criterion is that the hadronic interaction        *
  19 C     * length of a 5 GeV/c pion is correct.  Errors on the calculated *
  20 C     * hadronic interaction length for other momenta and other        *
  21 C     * particles in GEANT version 3.04 are less than 1% in most cases.*
  22 C     * For details see Memorandum OPAL/0037N/JA/md, ref. Hadron       *
  23 C     * Milestone 84/003, Calculation of Hadronic Interaction Lengths  *
  24 C     * for Mixtures.                                                  *
  25 C.    *                                                                *
  26 C.    *    ==>Called by : GSMIXT                                       *
  27 C.    *       Author    J.Allison  *********                           *
  28 C     *                                                                *
  29 C     ******************************************************************
  30 C
  31 #include "geant321/gcunit.inc"
  32       DIMENSION A(N), W(N)
  33 C.
  34 C.    ------------------------------------------------------------------
  35 C.
  36 C         Work out PINT which is proportional to the interaction
  37 C         probability.  Also work out the mean atomic weight, i.e. that
  38 C         weighted by proportion by numbers, as starting value for
  39 C         iterative method of finding AEFF.
  40 C
  41       PINT   = 0.
  42       AINV   = 0.
  43       WTOT   = 0.
  44       DO 10 I = 1, N
  45          PINT   = PINT + W(I) * GHSIGM(5., 8, A(I)) / A(I)
  46          AINV   = AINV + W(I) / A(I)
  47          WTOT   = WTOT + W(I)
  48   10  CONTINUE
  49 C
  50       IF ( ABS ( WTOT - 1. ) .GT. 0.01 ) GO TO 98
  51 C
  52 C         Work out AEFF which gives PINT for 5 GeV/c pion.
  53 C         (This is a short Newton's method loop.)
  54 C
  55       AEFF   = 1. / AINV
  56       PNEW   = GHSIGM ( 5., 8, AEFF ) / AEFF
  57       DA     = 1.
  58       ITER   = 0
  59    20 CONTINUE
  60       ITER   = ITER + 1
  61       AEFF   = AEFF + DA
  62       POLD   = PNEW
  63       DAOLD  = DA
  64       PNEW   = GHSIGM ( 5., 8, AEFF ) / AEFF
  65       DP     = PNEW - POLD
  66       DA     = (PINT - PNEW ) * DAOLD / DP
  67       IF ( ( ITER .GT. 1 ) .AND. ( ABS ( DA ) .GT. ABS ( DAOLD ) ) )
  68      +                            GO TO 97
  69       IF ( ABS ( DA ) .GT. 0.01 ) GO TO 20
  70 C
  71       RETURN
  72 C
  73 C            Error conditions.
  74 C
  75   97  CONTINUE
  76       WRITE (CHMAIL,197)
  77       CALL GMAIL(0,0)
  78  197  FORMAT ( ' ***** GHMIX : ',
  79      +'HADRONIC INTERACTION MIXTURE ROUTINE NOT CONVERGING')
  80       RETURN
  81 C
  82   98  CONTINUE
  83       WRITE (CHMAIL,198) WTOT
  84       CALL GMAIL(0,0)
  85  198  FORMAT ( ' ***** GHMIX : ',
  86      +'FRACTIONS BY WEIGHT OF MIXTURES DO NOT ADD UP TO 1',F10.4)
  87       END