* *----------------------------------------------------------------------* * * * Copyright (C) 2020: CERN * * All Rights Reserved. * * * * Source routine or FLUKA 4: * * * * Created on 24 September 2020 by David Horvath & Roberto Versaci * * ELI ERIC * * * * Modified on 23 Jan 2024 by David Horvath * * ELI ERIC * * * * This is a simplified user written source routine utilizing a * * separate source routine library. * * * * It is intended as an alternative new-user-friendly version of the * * source.f routine. Existing FLUKA 4 source routines remain * * compatible. * * * * Note that the beam card still has some meaning - in the scoring the * * maximum momentum used in deciding the binning is taken from the * * beam momentum. Other beam card parameters are obsolete. * * * * Output variables: * * * * nomore = if > 0 the run will be terminated * * * *----------------------------------------------------------------------* * * Quick start guide: * ------------------ * * This user source routine template aims to modernize the legacy routine * by implementing modern Fortran conventions and to provide built in * sampling functions. * * The users only need to change / add code between the BEGINNING and END * marks, one section for declaration of user variables, and one for * assigning values to the beam parameters. * * By default there is no user variable defined, and all code lines for * parameter assignment are commented out. These comments start with the * symbol: '*'. To enable one, the '*' needs to be deleted. * * (Note: In Fortran each code line should start in column 7 or further in.) * * Every beam parameter has a default value based on the FLUKA input file. * A parameter assignment should only be used if the default value has to be * modified. * * There are three ways to assign a value to a parameter: * * 1. Direct assignment: A parameter is equal to a value. For example: * * momentum_energy = 0.1D0 * * If the parameter defined as a double precision, then the assigned * value should be represented as double precision as well so as to * not loose numerical precision. To do this a 'D' exponential mark * must be used. * * 2. Using a sampling function: A parameter is assigned to a value, * which is calculated by a separate function. For example: * * coordinate_x = sample_flat_distribution( [min], [max] ) * * The parameters with the '[' and ']' brackets need to be replaced * with numbers, or user variables containing the desired values, like: * * coordinate_x = sample_flat_distribution( -1.0D0, 1.0D0 ) * * 3. Using a sampling subroutine: They are similar to function, but they * are not returning values directly; instead they modify the variables * in their argument list. For example: * * call sample_annular_distribution( [rmin], [rmax], coordinate_x, coordinate_y ) * * The example above has two input parameters with brackets, and two * output parameters (without bracket). The input parameters have to be * provided, as for functions. The output parameter names usually don't * need to be changed, but there are cases, where a subset of possible * output parameters has to be selected. * * For further details see the FLUKA manual. * *----------------------------------------------------------------------* module source_variables implicit none integer, save :: particle_code integer, save :: heavyion_atomic_number, heavyion_mass_number, heavyion_isomer double precision, save :: momentum_energy, particle_weight logical, save :: energy_logical_flag double precision, save :: divergence_x, divergence_y logical, save :: gaussian_divergence_logical_flag double precision, save :: coordinate_x, coordinate_y, coordinate_z integer, save :: direction_flag double precision, save :: direction_cosx, direction_cosy, direction_cosz double precision, save :: polarization_cosx, polarization_cosy, polarization_cosz double precision, save :: particle_age double precision, save :: kshort_component double precision, save :: delayed_radioactive_decay end module source_variables include 'source_library.inc' *=== Source ===========================================================* subroutine SOURCE ( nomore ) use source_library use source_variables implicit none double precision xdummy integer :: nomore, debug_lines = 100 logical, save :: first_run = .true., debug_logical_flag = .false. type(phase_space) phase_space_entry ! Function declarations double precision sample_flat_momentum_energy double precision sample_gaussian_momentum_energy double precision sample_maxwell_boltzmann_energy double precision sample_histogram_momentum_energy double precision sample_spectrum_momentum_energy double precision sample_discrete_momentum_energy double precision sample_flat_distribution double precision sample_gaussian_distribution double precision sample_histogram_file double precision sample_spectrum_file double precision sample_discrete_file double precision FLRNDM ! ==================================== ! BEGINNING of user declared variables ! ==================================== double precision proton_ratio double precision random_number ! ============================== ! END of user declared variables ! ============================== nomore = 0 call initialization() if ( first_run ) then ! ================================== ! BEGINNING of custom initialization ! ================================== ! ============================ ! END of custom initialization ! ============================ first_run = .false. end if ! ============================== ! BEGINNING of customizable code ! ============================== ! Exercise 1A divergence_x = 0.4D0 ! Exercise 1B divergence_y = WHASOU(1) ! Exercise 2 coordinate_x = sample_flat_distribution( -5.0D0, 5.0D0 ) coordinate_y = sample_gaussian_distribution( 0.0D0, 4.0D0 ) ! Exercise 3 momentum_energy = sample_histogram_momentum_energy( "histogram.txt", "MeV" ) energy_logical_flag = .true. ! Exercise 4 proton_ratio = 0.25D0 ! Uniform random number between [0,1) random_number = FLRNDM(xdummy) if ( random_number < proton_ratio ) then ! Proton particle particle_code = 1 else ! Electron particle particle_code = 3 end if ! Enable debug output debug_logical_flag = .true. ! ============================================== ! END of customizable code - Do not change below ! ============================================== if ( nomore .eq. 0 ) then call set_primary() if ( debug_logical_flag ) call print_primary( debug_lines ) end if return *=== End of subroutine Source =========================================* end