I made a damage simulator utility program. I don't really know how to upload it, so I'll just include it on this page.

I guess I did figure out how to upload it https://crawl.develz.org/mantis/view.php?id=2964

I think this is portable, but I've only used it in Visual Studio C++ 2010.

// Simulation.cpp : Defines the entry point for the console application.
//

#include <iostream>
#include <iomanip>
#include <string>
#include <cstdlib>
#include <time.h>

// Set up some globals, cause that's simpler for simulation purposes
enum SIM_TYPE
{
  UNARMED = 0,
  MELEE,
  THROWN_DART,
  THROWN_JAVELIN,
  SLING,
  BOW,
  CROSSBOW,

  SIM_TYPE_COUNT
} simulation_type = UNARMED;
const char simulation_names[7][15] = {
  "Unarmed",
  "Melee",
  "Thrown Dart",
  "Thrown Javelin",
  "Sling",
  "Bow",
  "Crossbow"
};
unsigned int weapon_skill = 0;
unsigned int fighting_skill = 0;
unsigned int strength = 5;
unsigned int dexterity = 5;
unsigned int weapon_damage = 0;
unsigned int weapon_bonus = 0;
unsigned int weapon_strength_weight = 5;
unsigned int ammo_damage = 0;
unsigned int ammo_bonus = 0;
unsigned int slaying_damage_bonus = 0;

enum BRAND_TYPES
{
  NONE = 0,
  VORPAL,
  ELEMENTAL,
  ELECTROCUTION,
  HOLY,
  SLAYING,
  STEEL,

  BRAND_TYPES_COUNT
} weapon_brand = NONE;
const char brand_names[7][14] = {
  "none",
  "vorpal",
  "elemental",
  "electrocution",
  "holy",
  "slaying",
  "steel"
};
BRAND_TYPES ammo_brand = NONE;
unsigned int monster_ac = 0;

// Output array for histogram generation
const unsigned int MAX_TRACKED_DAMAGE = 200;
unsigned int damage_distribution[MAX_TRACKED_DAMAGE];
unsigned int maximum_simulated_damage = 0;

// Random function, use standard c library for now
unsigned int random2(unsigned int max)
{
  if (max <= 1)
    return (0);

  unsigned int partn = RAND_MAX / max;

  while (true)
  {
    unsigned int bits = rand();
    unsigned int val  = bits / partn;

    if (val < max)
      return ((int)val);
  }
}

unsigned int random2avg(unsigned int max, unsigned int rolls)
{
    unsigned int sum = random2(max);

    for (unsigned int i = 0; i < (rolls - 1); i++)
        sum += random2(max + 1);

    return (sum / rolls);
}

// Function prototypes
unsigned int weighted_attribute();

// helper function
inline void set_to_input(unsigned int & value, unsigned int maximum = 1000)
{
  std::string input;
  std::getline(std::cin, input);

  if ( input.length() > 0 )
  {
    unsigned int input_val = atoi(input.c_str());
    value = std::min(maximum, input_val);
  }
}

void get_parameters()
{
  std::string input;

  // Get the type of simulation we are doing
  std::cout << "Choose weapon type (enter to keep current choice):" << std::endl;
  for (int count = 0; count < SIM_TYPE_COUNT; ++count)
  {
    std::cout << "  " << count + 1 << " - " << simulation_names[count] << std::endl;
  }
  std::cout << "Current Choice: " << simulation_type + 1 << " - " << simulation_names[simulation_type] << std::endl;

  std::getline(std::cin, input);
  if ( input.length() > 0 && input[0] >= '1' && input[0] <= '8' )
  {
    simulation_type = static_cast<SIM_TYPE>(input[0] - '1');
  }

  // Get skills
  if ( simulation_type < THROWN_DART )
  {
    std::cout << "Enter Fighting skill level (currently " << fighting_skill << ")" << std::endl;
    set_to_input(fighting_skill, 27);
  }

  std::cout << "Enter appropriate weapon skill (currently " << weapon_skill << ")" << std::endl;
  set_to_input(weapon_skill, 27);

  // Get stats
  std::cout << "Enter strength (currently " << strength << ")" << std::endl;
  set_to_input(strength);

  std::cout << "Enter dexterity (currently " << dexterity << ")" << std::endl;
  set_to_input(dexterity);

  // Get weapon parameters
  // Infer some default values based on the simulation type for ease of use
  switch (simulation_type)
  {
  case UNARMED:
    // Base punch damage, not counting mutations
    weapon_damage = 3;
    break;
  case MELEE:
    // No good default
    break;
  case THROWN_DART:
    weapon_damage = 0;
    ammo_damage = 2;
    break;
  case THROWN_JAVELIN:
    weapon_damage = 0;
    ammo_damage = 10;
    break;
  case SLING:
    // the default will be bullets
    weapon_damage = 6;
    ammo_damage = 6;
    break;
  case BOW:
    // Default to longbow
    weapon_damage = 6;
    ammo_damage = 7;
    break;
  case CROSSBOW:
    weapon_damage = 5;
    ammo_damage = 9;
    break;
  default:
    break;
  }

  if (simulation_type != THROWN_DART && simulation_type != THROWN_JAVELIN)
  {
    std::cout << "Enter weapon damage (currently " << weapon_damage << ")" << std::endl;
    set_to_input(weapon_damage);

    std::cout << "Enter weapon damage bonus (currently " << weapon_bonus << ")" << std::endl;
    set_to_input(weapon_bonus);
  }

  switch (simulation_type)
  {
  case UNARMED:
    weapon_strength_weight = 4;
    std::cout << "Your weighted attribute is " << weighted_attribute() << std::endl;
    break;
  case MELEE:
    std::cout << "Enter weapon strength weight (currently " << weapon_strength_weight << ")" << std::endl;
    set_to_input(weapon_strength_weight, 10);
    std::cout << "Your weighted attribute is " << weighted_attribute() << std::endl;
    break;
  case THROWN_DART:
    // doesn't use stat weighting
    break;
  case THROWN_JAVELIN:
    // doesn't use standard stat weighting
    break;
  case SLING:
    // Hrm, stat weighting only used for calculating delay, so ignoring it here
    break;
  case BOW:
    // once again, stat weighting only used for delay
    break;
  case CROSSBOW:
    // once again, stat weighting only used for delay
    break;
  }

  if (simulation_type == SLING)
  {
    std::cout << "Enter ammo damage, 4 is sling stone, 6 is sling bullet (currently " << ammo_damage << ")" << std::endl;
    set_to_input(ammo_damage);
  }

  if (simulation_type >= THROWN_DART)
  {
    std::cout << "Enter ammo bonus (currently " << ammo_bonus << ")" << std::endl;
    set_to_input(ammo_bonus);
  }

  if (simulation_type != THROWN_DART && simulation_type != THROWN_JAVELIN)
  {
    // Get the multiplier for the damage brand bonus
    std::cout << "Enter weapon brand" << std::endl;
    for (int count = 0; count < BRAND_TYPES_COUNT; ++count)
    {
      std::cout << "  " << count + 1 << " - " << brand_names[count] << std::endl;
    }
    std::cout << "Current brand: " << weapon_brand + 1 << " - " << brand_names[weapon_brand] << std::endl;
    unsigned int brand_index = weapon_brand + 1;
    set_to_input(brand_index,BRAND_TYPES_COUNT+1);
    --brand_index;
    if (brand_index >= BRAND_TYPES_COUNT)
      brand_index = 0;
    weapon_brand = static_cast<BRAND_TYPES>(brand_index);
  }

  if (simulation_type >= THROWN_DART)
  {
    // Get the brand for the ammunition too
    // Get the multiplier for the damage brand bonus
    std::cout << "Enter ammo brand" << std::endl;
    for (int count = 0; count < BRAND_TYPES_COUNT; ++count)
    {
      std::cout << "  " << count + 1 << " - " << brand_names[count] << std::endl;
    }
    std::cout << "Current ammo brand: " << ammo_brand + 1 << " - " << brand_names[ammo_brand] << std::endl;
    unsigned int brand_index = ammo_brand + 1;
    set_to_input(brand_index,BRAND_TYPES_COUNT+1);
    --brand_index;
    if (brand_index >= BRAND_TYPES_COUNT)
      brand_index = 0;
    ammo_brand = static_cast<BRAND_TYPES>(brand_index);

    // Note that the ammo brand overrides the bow brand if the ammo is branded
    // This only happens for thrown/fired items thanks to the parent if check
    // Oh right, also vorpal still works on the launcher
    if ( ammo_brand != NONE )
    {
      if ( weapon_brand != VORPAL )
        weapon_brand = NONE;
    }
  }

  // Get the total slaying bonus
  std::cout << "Enter the total slaying damage bonus (currently " << slaying_damage_bonus << ")" << std::endl;
  set_to_input(slaying_damage_bonus);

  // Get the monster's AC for damage reduction
  std::cout << "Enter target AC (currently " << monster_ac << ")" << std::endl;
  set_to_input(monster_ac);
}

unsigned int weighted_attribute()
{
  // Note that this form of the calculation seems to be a hell of a lot simpler than the actual code
  // return (strength * strength_weight + dexterity * (10-strength_weight))/20;
  // But, the actual code rounds str + dex / 2 off with integer math, so we should use the more complicated form
  signed int modified = ((signed int)dexterity - (signed int)strength)/2;
  return (signed int)strength + modified*(10-(signed int)weapon_strength_weight)/10;
}

unsigned int sim_melee_attack()
{
  unsigned int potential_damage = weapon_damage;

  if ( simulation_type == UNARMED )
  {
    // You need to account for extra damage from a damaging form or mutations when setting the sim parameters
    // So set weapon_damage = 3 + mutation modifiers or = form damage + str/dex modifiers
    // Yes, for Dragon/Stone/Bladehand forms you need to do the base damage math for me, sorry
    potential_damage += weapon_skill;
  }

  // Modify damage randomly based on stat bonus
  unsigned int attribute = weighted_attribute();
  potential_damage = potential_damage
    * (78
      + (attribute > 11 ? (random2(attribute - 11) * 2): 0)
      - (attribute < 9 ? (random2(9 - attribute) * 3) : 0) )
    / 78;

  // Add slaying bonus to potential damage
  potential_damage += slaying_damage_bonus;

  // randomize the damage, adding 0-2 randomly for great justice!
  if ( potential_damage <= 0 )
    potential_damage = 0;
  else
    potential_damage = random2(3) + random2(potential_damage);

  if ( simulation_type == MELEE )
  {
    // Apply a random multiplier based on weapon skill
    potential_damage = potential_damage * (25 + random2(weapon_skill + 1)) / 25;
  }

  // Apply a random multiplier based on fighting skill
  potential_damage = potential_damage * (30 + random2(fighting_skill + 1)) / 30;

  // NOT SIMULATED +1d10 from might

  if ( simulation_type == MELEE )
  {
    // Add weapon damage bonus
    potential_damage += random2(weapon_bonus + 1);

    // NOT SIMULATED dwarven, orcish, demonspawn racial damage bonuses
  }

  // NOT SIMULATED stabbing damage bonuses
  // If you're a solid stabbity, you probably your entire effective damage from that
  // This is an assumption, I haven't really studied the stabbity stabbity code-ity

  // Reduce damage based on monster AC
  unsigned int ac_reduction = random2(monster_ac + 1);
  if ( potential_damage > ac_reduction )
    potential_damage -= ac_reduction;
  else
    potential_damage = 0;

  // NOT SIMULATED 1/3rd damage based on target petrification...if you petrify things,
  // hopefully you'll figure this out for yourself

  // NOT SIMULATED auxiliary unarmed attacks
  // these probably add considerable damage at the cost of extra training in unarmed skill
  // and equippment considerations

  // Simulate brand damage
  if ( potential_damage > 0 )
  {
    switch (weapon_brand)
    {
    case VORPAL:
      potential_damage += 1 + random2(potential_damage) / 4;
      break;
    case ELEMENTAL:
      // Assume not weak or resistant
      potential_damage += 1 + random2(potential_damage) / 2;
      break;
    case ELECTROCUTION:
      // Assuming they're not flying or resistant
      if ( random2(3) == 0 )
      {
        potential_damage += 10 + random2(15);
      }
      break;
    case HOLY:
      // The only reason we would calculate this is if we are fighting that type of enemy
      // So assume it applies
      // Oh come on, this random distribution may or may not be similar to the slaying one
      // but probably is slightly different
      potential_damage += 1 + (random2(potential_damage * 15) / 10);
      break;
    case SLAYING:
      // The only reason we would calculate this is if we are fighting that type of enemy
      // So assume it applies
      potential_damage += 1 + random2(3*potential_damage/2);
      break;
    default:
      break;
    }
  }

  return potential_damage;
}

unsigned int sim_thrown_attack()
{
  // NOT SIMULATED: dwarven/orcish thrown weapons

  // crazy negative numbers
  signed int extra_damage = (10 * (weapon_skill / 2 + strength - 10)) / 12;
  // the comment in the source seems to say that the extra_damage is a multiplier,
  // and it's trying to reduce the multiplier to 1/3rd of the damage bonus,
  // BUT I don't think it is a multiplier, and this only makes sense for base
  // damage 10 javelins
  extra_damage = (extra_damage * (3 * weapon_damage + ammo_bonus)) / 30;

  if ( simulation_type == THROWN_DART )
  {
    // Lets add in our throwing skill again
    extra_damage += weapon_skill * 3 / 5;
  }

  if ( simulation_type == THROWN_JAVELIN )
  {
    extra_damage += weapon_skill * 3 / 5;

    // multiply the extra damage by a random factor determined by strength
    signed int multiplier = (15 + (strength - 15) / 2) * 100 / 15;
    multiplier = std::min(std::max(90,multiplier),160);

    if (multiplier > 100)
      extra_damage = extra_damage * (100 + random2avg(multiplier - 100, 2)) / 100;
    else if (multiplier < 100)
      extra_damage = extra_damage * (100 - random2avg(100 - multiplier, 2)) / 100;

    // multiply the extra damage by a random factor determined by dexterity
    multiplier = (20 + (dexterity - 20) / 2) * 100 / 20;
    multiplier = std::min(std::max(100,multiplier),150);

    if (multiplier > 100)
      extra_damage = extra_damage * (100 + random2avg(multiplier - 100, 2)) / 100;
    else if (multiplier < 100)
      extra_damage = extra_damage * (100 - random2avg(100 - multiplier, 2)) / 100;

    // confusing enough yet?
  }

  // now we decide our effective die size by using our base damage and
  // adding a random amount of our extra_damage, effectively dividing that by half
  unsigned int potential_damage = ammo_damage;
  
  // NOT SIMULATED a negative extra damage penalty, which might happen for low strength
  // low skill characters
  if (extra_damage > 0)
  {
    potential_damage += random2(extra_damage + 1);
  }

  // Apply our brand damage bonus here (basically only steel darts/javelins)
  // Elemental damage is applied in a later step
  if ( ammo_brand == STEEL )
    potential_damage = potential_damage * 150 / 100;

  // Randomize our slaying bonus, becuase effectively dividing it by half for ranged
  // attacks is cool, right?
  potential_damage += random2(slaying_damage_bonus + 1);

  // Add in our enchanted bonus, can you believe this wasn't randomly added in?
  potential_damage += ammo_bonus;

  // And here is where it gets weird, if our damage die size is > 24, we shrink it
  // and add more damage dice
  // On the plus side, if you can get your damage above the 24 threshold, it turns
  // into a much more normalized distribution
  // But because of all the randomization, it seems hard to do
  unsigned int die_count = 1;
  while ( potential_damage > 24 )
  {
    die_count *= 2;
    potential_damage /= 2;
  }

  // And roll the dice
  unsigned int damage = die_count;

  for (unsigned int i = 0; i < die_count; i++)
    damage += random2(potential_damage);
  
  // Reduce damage based on monster AC
  unsigned int ac_reduction = random2(monster_ac + 1);
  if ( damage > ac_reduction )
    damage -= ac_reduction;
  else
    damage = 0;

  // Simulate brand damage
  if ( damage > 0 )
  {
    switch (ammo_brand)
    {
    case ELEMENTAL:
      // Assume not weak or resistant
      damage += 1 + random2(damage) / 2;
      break;
    default:
      break;
    }
  }

  return damage;
}

unsigned int sim_launched_attack()
{
  // Lets start of randomizing our random numbers
  unsigned int base_damage = weapon_damage + random2(ammo_damage + 1);

  unsigned int multiplier = 100;

  if ( simulation_type == SLING )
  {
    // And slings remove randomization, cause their numbers are tiny to begin with
    // I think this brings them much closer to bows and crossbows
    base_damage = ammo_damage;
  }

  // NOT SIMULATED: extra damage for racial launcher
  // I think you can just set the weapon damage to 4/7/6 for bow/longbow/crossbow to simulate this

  // idea: demonic bow/sling/crossbow created from the bones/flesh/teeth of the damned
  // cause demonspawn deserve a racial damage bonus?

  // Hmm, randomizing the bonus from the ammo effectively divides it by half
  unsigned int extra_damage = weapon_bonus + random2(ammo_bonus + 1);

  // And we randomize the bonus again!  on average halving weapon and quartering ammo bonus
  extra_damage = random2(extra_damage + 1);

  switch (simulation_type)
  {
  case SLING:
    // add a strength bonus which maxes at sling damage bonus + 1
    {
      unsigned int strength_bonus = (((10 * (strength - 10)) / 9) * (2 * base_damage + ammo_bonus)) / 20;
      extra_damage += std::min(weapon_bonus + 1, strength_bonus);
    }
    // Modify the dice multiplier by a random factor based on skill
    multiplier = multiplier * (14 + random2(weapon_skill + 1)) / 14;
    break;
  case BOW:
    // add a strength bonus which maxes at sling damage bonus + 1, this requires much less strength than slings
    {
      unsigned int strength_bonus = (((10 * (strength - 10)) / 4) * (2 * base_damage + ammo_bonus)) / 20;
      extra_damage += std::min(weapon_bonus + 1, strength_bonus);
    }
    // Modify the dice multiplier by a random factor based on skill
    // Not quite as good a multiplier as slings get, maybe to balance the higher base damage
    multiplier = multiplier * (17 + random2(weapon_skill + 1)) / 17;
    break;
  case CROSSBOW:
    // No stat bonus to damage, instead they'll get the full launcher bonus again later on
    // Modify the dice multiplier by a random factor based on skill
    // even worse multiplier than bows
    multiplier = multiplier * (22 + random2(weapon_skill + 1)) / 22;
    break;
  default:
    break;
  }
  
  // now we decide our effective die size by using our base damage and
  // adding a random amount of our extra_damage, effectively dividing that by half
  unsigned int potential_damage = base_damage + random2(extra_damage + 1);

  // Apply our steel or vorpal brand damage bonus here
  // Elemental damage is applied later
  if ( weapon_brand == VORPAL )
    multiplier = multiplier * 130 / 100;

  if ( ammo_brand == STEEL )
    multiplier = multiplier * 150 / 100;

  // scale the damage by our multiplier
  potential_damage = potential_damage * multiplier / 100;

  // Randomize our slaying bonus, becuase effectively dividing it by half for ranged
  // attacks is cool, right?
  potential_damage += random2(slaying_damage_bonus + 1);

  // Add in our enchanted bonus, can you believe this wasn't randomly added in?
  potential_damage += ammo_bonus;

  if ( simulation_type == CROSSBOW )
  {
    // This wasn't part of the extra damage for crossbows
    // Interestingly that means the damage bonus on crossbows isn't re-randomized
    potential_damage += weapon_bonus;
  }

  // And here is where it gets weird, if our damage die size is > 24, we shrink it
  // and add more damage dice
  // On the plus side, if you can get your damage above the 24 threshold, it turns
  // into a much more normalized distribution
  // But because of all the randomization, it seems hard to do
  unsigned int die_count = 1;
  while ( potential_damage > 24 )
  {
    die_count *= 2;
    potential_damage /= 2;
  }

  // And roll the dice
  unsigned int damage = die_count;

  for (unsigned int i = 0; i < die_count; i++)
    damage += random2(potential_damage);
  
  // Reduce damage based on monster AC
  unsigned int ac_reduction = random2(monster_ac + 1);
  if ( damage > ac_reduction )
    damage -= ac_reduction;
  else
    damage = 0;

  // Simulate brand damage
  if ( damage > 0 )
  {
    // The weapon_brand is used unless the ammo is branded
    switch ((ammo_brand == NONE) ? weapon_brand : ammo_brand)
    {
    case ELEMENTAL:
      // Assume not weak or resistant
      damage += 1 + random2(damage) / 2;
      break;
    case HOLY:
      // The only reason we would calculate this is if we are fighting that type of enemy
      // So assume it applies
      // Oh come on, this random distribution may or may not be similar to the slaying one
      // but probably is slightly different
      damage += 1 + (random2(damage * 15) / 10);
      break;
      // TODO: implement double damage for silver
    //case SILVER:
    //  break;
    default:
      break;
    }
  }

  return damage;
}

void run_simulation()
{
  // Clear the output histogram
  for (int i = 0; i < MAX_TRACKED_DAMAGE; ++i)
    damage_distribution[i] = 0;
  maximum_simulated_damage = 0;

  double average_damage = 0;

  // Run 1000 damage calculations and record the result to the damage histogram
  for (int i = 0; i < 1000; ++i)
  {
    unsigned int damage = 0;
    switch (simulation_type)
    {
    case UNARMED:
    case MELEE:
      damage = sim_melee_attack();
      break;
    case THROWN_DART:
      damage = sim_thrown_attack();
      break;
    case THROWN_JAVELIN:
      damage = sim_thrown_attack();
      break;
    case SLING:
    case BOW:
    case CROSSBOW:
      damage = sim_launched_attack();
      break;
    }

    if ( damage >= MAX_TRACKED_DAMAGE )
    {
      average_damage += damage - (MAX_TRACKED_DAMAGE + 1);
      damage = MAX_TRACKED_DAMAGE - 1;
    }

    // This should set the upper bound for our histogram printout
    maximum_simulated_damage = std::max(damage,maximum_simulated_damage);

    ++damage_distribution[damage];
  }

  ++maximum_simulated_damage;
  maximum_simulated_damage = std::min(maximum_simulated_damage, MAX_TRACKED_DAMAGE);

  bool set_median_damage = false;
  unsigned int median_damage = 0;
  unsigned int hit_count = 0;
  // output the resulting histogram and statistical values
  std::cout << "Simulated Damage Histogram" << std::endl
            << "Damage Count" << std::endl;
  for (unsigned int i = 0; i < maximum_simulated_damage; ++i)
  {
    std::cout << std::setw(6) << i << " " << damage_distribution[i] << std::endl;
    average_damage += i * damage_distribution[i];
    hit_count += damage_distribution[i];
    if ( !set_median_damage && hit_count > 500 )
    {
      set_median_damage = true;
      median_damage = i;
    }
  }

  average_damage /= 1000.0;
  std::cout << std::endl
    << "Mean Damage: " << average_damage << " Median Damage: " << std::setprecision(3) << median_damage << std::endl;

  std::cout << "Chance of zero damage: " << (double)damage_distribution[0] / 10.0 << "%" << std::endl;
}

int main(int argc, char* argv[])
{
  srand( (unsigned int)time(NULL) );
  for (int i = 0; i < 10000; ++i)
    rand();

  std::cout << "Welcome to the Dungeon Crawl Stone Soup damage simulator for version 0.7.1" << std::endl;

  // Main input loop
  for (;;)
  {
    // Output main menu
    std::cout << std::endl
              << "Choose an option, invalid input exits." << std::endl
              << "1 - Set simulation parameters" << std::endl
              << "2 - Run simulation" << std::endl;

    std::string input;

    std::getline(std::cin,input);

    if (input.length() > 0 && input[0] == '1')
    {
      get_parameters();
    }
    else if (input.length() > 0 && input[0] == '2')
    {
      run_simulation();
    }
    else
    {
      break;
    }
  }

  return 0;
}