#! /usr/bin/env python3
#
def duel_many_test ( p ):

#*****************************************************************************80
#
## duel_many_test() tests duel_simulation() many times, and averages results.
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    15 January 2023
#
#  Author:
#
#    John Burkardt
#
#  Input:
#
#    real P(2), the probabilities that player I will hit an opponent
#    in a single shot.  Each P(I) should be between 0 and 1.
#
  import numpy as np

  print ( '' )
  print ( 'duel_many_test():' )
  print ( '  duel_many() simulates many duels between combatants' )
  print ( '  who fire alternately at each other, with known probabilities' )
  print ( '  of hitting the other.' )
  
  duel_num = 1000
 
  print ( '' )
  print ( '  In this example, the one shot probabilities are:' )
  print ( '  Player[0]: %g' % ( p[0] ) )
  print ( '  Player[1]: %g' % ( p[1] ) )
  print ( '' )
  print ( '  The number of duels to fight is %d' % ( duel_num ) )

  s, shot_average = duel_simulation ( p, duel_num )

  print ( '' )
  print ( '  Result of the duels:' )
  for player in range ( 0, 2 ):
    print ( '  Player %d winning percentage is %g' % ( player, s[player] ) )

  print ( '' )
  print ( '  The average number of shots fired was %g' % ( shot_average ) )

  return

def duel_once ( p ):

#*****************************************************************************80
#
## duel_once() returns the outcome of a single duel.
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    15 January 2023
#
#  Author:
#
#    John Burkardt
#
#  Reference:
#
#    Martin Shubik,
#    "Does the Fittest Necessarily Survive?",
#    in Readings in Game Theory and Political Behavior,
#    edited by Martin Shubik,
#    Doubleday, 1954,
#    LC: H61.S53.
#
#  Input:
#
#    real P(2), the probabilities that player I will hit an opponent
#    in a single shot.  Each P(I) should be between 0 and 1.
#
#  Output:
#
#    integer SURVIVOR, the survivor of the duel.
#
  from numpy.random import default_rng
  import numpy as np

  rng = default_rng ( )

  shot_num = 0
 
  while ( True ):

    for player in range ( 0, 2 ):
      shot_num = shot_num + 1
      hit = ( rng.random ( ) <= p[player] )
      if ( hit ):
        survivor = player
        return survivor, shot_num

def duel_once_test ( p ):

#*****************************************************************************80
#
## duel_once_test() tests duel_once().
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    15 January 2023
#
#  Author:
#
#    John Burkardt
#
#  Input:
#
#    real P(2), the probabilities that player I will hit an opponent
#    in a single shot.  Each P(I) should be between 0 and 1.
#
  import numpy as np

  print ( '' )
  print ( 'duel_once_test():' )
  print ( '  duel_once() simulates a single duel between combatants' )
  print ( '  who fire alternately at each other, with known probabilities' )
  print ( '  of hitting the other.' )
  
  print ( '' )
  print ( '  In this example, the one shot probabilities are:' )
  print ( '  Player[0]: %g' % ( p[0] ) )
  print ( '  Player[1]: %g' % ( p[1] ) )

  survivor, shot_num = duel_once ( p )

  print ( '' )
  print ( '  Result of this duel:' )
  print ( '  Player %d survived after %d turns' % ( survivor, shot_num ) )

  return

def duel_simulation ( p, duel_num ):

#*****************************************************************************80
#
## duel_simulation() simulates a duel.
#
#  Discussion:
#
#    Player 1 fires at player 2, and hits with a probability of P(1).
#    If Player 2 misses, then Player 2 fires at Player 1, hitting with
#    a probability of P(2).
#
#    The duel continues with alternating shots until only one player 
#    survives.
#
#    The simulation is intended to estimate the probabilities that a
#    player will survive, and the number of turns required.
#
#    The exact probability that player 1 will survive is
#
#      P(1) / ( P(1) + P(2) - P(1) * P(2) )
#
#    Player 2's chance is
#
#     P(2) * ( 1 - P(1) ) / ( P(1) + P(2) - P(1) * P(2) )
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    16 August 2017
#
#  Author:
#
#    John Burkardt
#
#  Reference:
#
#    Paul Nahin,
#    Duelling Idiots and Other Probability Puzzlers,
#    Princeton University Press, 2000,
#    ISBN13: 978-0691009797,
#    LC: QA273.N29.
#
#    Martin Shubik,
#    "Does the Fittest Necessarily Survive?",
#    in Readings in Game Theory and Political Behavior,
#    edited by Martin Shubik,
#    Doubleday, 1954,
#    LC: H61.S53.
#
#  Input:
#
#    real P(2), the probabilities that player I will hit an opponent
#    in a single shot.  Each P(I) should be between 0 and 1.
#
#    integer DUEL_NUM, the number of duels to run.
#
#  Output:
#
#    real S(2), the estimated probablity that player I survives.
#
#    real shot_AVERAGE, the average number of turns.
#
  import numpy as np

  s = np.zeros ( 2 )
  shot_average = 0

  for duel in range ( 0, duel_num ):
    survivor, shot_num = duel_once ( p )
    s[survivor] = s[survivor] + 1
    shot_average = shot_average + shot_num

  s = s / float ( duel_num )
  shot_average = shot_average / float ( duel_num )

  return s, shot_average

def duel_simulation_test ( ):

#*****************************************************************************80
#
## duel_simulation_test() tests duel_simulation().
#
#  Licensing:
#
#    This code is distributed under the MIT license.
#
#  Modified:
#
#    23 March 2021
#
#  Author:
#
#    John Burkardt
#
  import numpy as np
  import platform

  print ( '' )
  print ( 'duel_simulation_test():' )
  print ( '  Python version: %s' % ( platform.python_version ( ) ) )
  print ( '  Test duel_simulation()' )

  p = np.array ( [ 0.10, 0.15 ] )
  duel_once_test ( p )

  p = np.array ( [ 0.15, 0.10 ] )
  duel_once_test ( p )

  p = np.array ( [ 0.10, 0.15 ] )
  duel_many_test ( p )

  p = np.array ( [ 0.15, 0.10 ] )
  duel_many_test ( p )

  duel_num = 1000
  p = np.array ( [ 0.10, 0.15 ] )
  s, shot_average = duel_simulation ( p, duel_num )
  print ( '' )
  print ( '  Player 0 survival probability = ', s[0] )
  print ( '  Player 1 survival probability = ', s[1] )
  print ( '  Average number of shots fired = ', shot_average )

  shot_num = 1000
  p = np.array ( [ 0.15, 0.10 ] )
  s, shot_average = duel_simulation ( p, shot_num )
  print ( '' )
  print ( '  Player 0 survival probability = ', s[0] )
  print ( '  Player 1 survival probability = ', s[1] )
  print ( '  Average number of shots fired = ', shot_average )
#
#  Terminate.
#
  print ( '' )
  print ( 'duel_simulation_test():' )
  print ( '  Normal end of execution.' )

  return

def timestamp ( ):

#*****************************************************************************80
#
## timestamp() prints the date as a timestamp.
#
#  Licensing:
#
#    This code is distributed under the MIT license. 
#
#  Modified:
#
#    06 April 2013
#
#  Author:
#
#    John Burkardt
#
  import time

  t = time.time ( )
  print ( time.ctime ( t ) )

  return None

if ( __name__ == '__main__' ):
  timestamp ( )
  duel_simulation_test ( )
  timestamp ( )