Hello.
        Here is the message I received from Berkley concerning the correction
to the Monte Carlo of the stopping power.

Note that the subroutine has a new argument (material_type). The old version
always assumed that the first collision was in the stopper, never in the
target.

        If you have any question, don't hesitate to contact me. I will be at
Risoe until mid-august.

                                                  Jules Gascon
                                                  GASCON@DKNBI51

******************************************************************************

Hi
 
Your message to Jose Bacelar regarding his lifetime program was forwarded
to him while he was visiting here. Since then we have also been examining
the program.
 
We think that we agree that there is a problem with the subroutine
get_coll_epslon, which calculates the energy for the next nuclear scattering.
We also do not fully understand what Jose is doing, particularly in the
nested do loops. We have rewritten the subroutine, in a way which we think
is correct, and have included this below.
 
        Essentially this replaces the loops over div and n with a single
do loop over div = 1, 24 although it converges (very well) after only about
16-20. The k_e(2)*sqrt(epslon) comes from converting
                dn --> drho --> depslon
using deps/drho = k_e * sqrt(epslon).
We used a better five term expansion in sqrt(esplon) and thus the k_e has
disappeared (and been replaced by a_dedx, b_dedx etc). This expansion also
has to be included in the summation to calculate INT.
 
The first collision energies came out to be much difference (larger) with
the new version, but the resulting lifetimes did not change. We think that
this is because after the first few collisions the energy intervals are
so small that the precise value of the energy does not matter.
 
From your note you imply that you get a significant difference in the
lifetimes. Could you bitnet us a copy of your program and also a copy of
Bernard Hass' program so that we can compare.
 
Many thanks,
        Con Beausang (BEAUSANG@LBL)
        Dick Diamond (RMD@LBL)
 
 
        subroutine get_coll_epslon(material_type,epslon,coll_epslon)
        implicit none
        include 'libdir:dechist.cmn'
 
        integer div,material_type
        real delta_epslon,bl_int,int,previous_int
        real coll_energy, dedx, drhodeps, sinangle
        real dedxtest
 
        coll_int = -log ( 1.0 - ran(seed2)) / mu(material_type)
 
        sinangle = sin(angle_cm_min)
 
        previous_int = 0.0
        delta_epslon =  (epslon - fx_abs_min_limit) / 2.
        coll_epslon = epslon -delta_epslon/2.
 
        do div = 1, 24
        int = previous_int
        epslon_min = coll_epslon * sinangle
        if (epslon_min .lt. fx_abs_min_limit)
        1                       epslon_min = fx_abs_min_limit
        call get_blaug_int(coll_epslon,epslon_min,bl_int)
        coll_energy = coll_epslon / energy_epslon_conv(material_type)
        dedx   = a_dedx(material_type)
        1       + b_dedx(material_type) * coll_energy**0.5
        1       + c_dedx(material_type) * coll_energy
        1       + d_dedx(material_type) * coll_energy ** 1.5
        1       + e_dedx(material_type) * coll_energy ** 2.
c
c       In case the fit goes negative dedxtest !!!!
c
        if (dedx.le.0.025) dedx = 0.025
 
        drhodeps = x_rho_conv(material_type) /
        1       energy_epslon_conv(material_type) / dedx
        int = int + bl_int * drhodeps * delta_epslon
 
        if (int .gt. coll_int)goto 10
 
        previous_int = int
        coll_epslon = coll_epslon - 0.75 * delta_epslon
        goto 11
 
10      coll_epslon = coll_epslon + 0.25 * delta_epslon
11      delta_epslon =  (delta_epslon) / 2.
        end do
 
        end
