removed "PF" as a global variable and set it equal to something more shocked for Apollo samples (for a test)

src/globals/module_globals.f90

@@ -60,7 +60,6 @@ module module_globals
 real(DP),parameter  :: DISEJB = 100.0_DP       ! The extent of discontinuous ejecta in the unit of crater radii. It is used in ejecta_table_define.f90
 real(DP),parameter :: RCONT = 2.25267_DP       ! Coefficient of continuous ejecta size power law from Moore et al. (1974) - scaled from km to m
 real(DP),parameter :: EXPCONT = 1.006_DP       ! Exponent of continuous ejecta size power law from Moore et al. (1974) 
-real(DP),parameter :: PF     = 5.0e9          ! The shock pressure exceeding the hungoit elastic limit of common geological materials in Pa. (5 GPa)
 real(DP),parameter :: RAD_GP = 1.0_DP         ! The maximum radial position of producing impact glass spherules within a transient crater (unit of crater radii, crater%rad)
 
 
@@ -163,6 +162,7 @@ module module_globals
    integer(I4B) :: vnum      ! size of velocity distribution array
    integer(I4B) :: rcnum     ! size of real crater list array for quasi-MC
    integer(I4B) :: rccount   ! quasimc crater for iterating through the list of craters
+   integer(I4B) :: local     ! 
    integer(I4B) :: nqmc      ! current Quasi-MC crater for regolith distribution
    logical      :: currentqmc ! is the current crater a Quasi-MC crater?
    real(DP)     :: vescsq    ! Escape velocity at target

src/regolith/regolith_shock_damage_zone.f90

@@ -68,10 +68,11 @@ subroutine regolith_shock_damage_zone(crater,rm,eradi,depthb,xsfints)
    real(DP),parameter          :: c_dunite   = 6500.0
    real(DP),parameter          :: s_dunite   = 0.9
    real(DP),parameter          :: n          = -2.85
+   real(DP)                    :: pf         = 25.0e9    ! The shock pressure exceeding the hungoit elastic limit of common geological materials in Pa is 5 GPa.
    real(DP)                    :: pmax, up
    real(DP)                    :: rsh, xshints0
    real(DP)                    :: q1, q2, q3, thetaq 
-
+   
    ! Calculate the radius of shock pressure decay zone by a given shock 
    ! P(r) = Pmax * (r/r_p)**(n), where n is negative and we use -3 (fast decay
    ! that may be accounted for by acoustic fluidization).
@@ -80,7 +81,7 @@ subroutine regolith_shock_damage_zone(crater,rm,eradi,depthb,xsfints)
    ! The maximum peak shock pressure is estimated from planar impact 
    up       = crater%impvel * 0.5 
    pmax     = rho_dunite * (c_dunite + s_dunite * up) * up
-   rsh      = (crater%imp/2.0) * (PF / pmax)**(1.0/n)
+   rsh      = (crater%imp/2.0) * (pf / pmax)**(1.0/n)
    xshints0 = sqrt(rsh**2 - (crater%imp / 2.0)**2)
 
    if (eradi <= xshints0) then