diff --git a/src/simulation.py b/src/simulation.py index c46f9d3..54e2ecc 100644 --- a/src/simulation.py +++ b/src/simulation.py @@ -174,12 +174,12 @@ def simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, l run_results = [res[0] for res in successful_results] loss_rates = [res[1] for res in successful_results] + n = len(run_results) # reject if not enough successful run - if len(run_results) >= min_runs: + if n >= min_runs: # statistics mean_rt = np.mean(run_results) std_dev = np.std(run_results, ddof=1) - n = len(run_results) # confidence interval t_value = t.ppf((1 + confidence_level)/2, n-1) @@ -196,8 +196,8 @@ def simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, l losses.append(mean_loss) print(f"C={c}, λ={lambda_val:.2f}, Mean RT={mean_rt:.2f} ± {ci:.2f}, Mean Loss Rate={mean_loss:.2%}") - elif len(run_results) > 0: - print(f"λ={lambda_val:.2f} skipped - only {len(run_results)} successful run(s)") + elif n > 0: + print(f"λ={lambda_val:.2f} skipped - only {n} successful run(s)") continue else: print(f"Stopped at λ={lambda_val:.2f} - no successful run") @@ -213,6 +213,85 @@ def simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, l loss_data[c] = (np.array(lambda_points), np.array(losses)) + # plot curves + if len(lambda_vals)>1: + ax_rt.set_xlim(left=0) + ax_rt.set_xlabel('Arrival Rate (λ)') + ax_rt.set_ylabel('Mean Response Time') + ax_rt.set_title(f'Mean Response Time vs Arrival Rate ({num_runs} runs, {int(confidence_level*100)}% CI)') + ax_rt.legend() + ax_rt.grid(True) + + return loss_data, mean_at_lambda1 + +def simulation_wrapper2(simulation_time, num_runs, min_runs, confidence_level, lambda_vals, max_loss_value=0.1): + fig, (ax_rt, ax_loss) = plt.subplots(2,1, figsize=(12,10), gridspec_kw={'height_ratios': [4, 3], 'hspace': 0.4}, sharex=False) + loss_data, mean_at_lambda1 = simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, lambda_vals, ax_rt) + + if len(lambda_vals) > 1: + for c, (lams, losses) in loss_data.items(): + lambda_vals2 = [] + if len(lams)==0: + print(f"C={c}: no data at all.") + continue + pos_idx = np.where(losses > 0)[0] + if pos_idx.size > 0: + first_lambda = lams[pos_idx[0]] + else: + first_lambda = None + last_lambda = lams[-1] + + print(f"C={c:>1} first λ with loss > 0: " + f"{first_lambda if first_lambda is not None else 'never'}; " + f" last λ with data: {last_lambda}") + + lambda_vals2.extend([last_lambda + i/1000 for i in range(-100,51)]) + + + + lambda_points = [] + losses = [] + ci_lower = [] + ci_upper = [] + print(f"\nProcessing C={c}") + + with Pool() as pool: + for lambda_val in lambda_vals2: + args_list = [(c, lambda_val, simulation_time, max_loss_value) for _ in range(num_runs)] + results = pool.map(run_single_simulation, args_list) + + successful_results = [res for res in results if res is not None] + loss_rates = [res[1] for res in successful_results] + + n = len(loss_rates) + + if n >= min_runs: + # statistics + mean_loss = np.mean(loss_rates) + std_dev = np.std(loss_rates, ddof=1) + + # confidence interval + t_value = t.ppf((1 + confidence_level)/2, n-1) + ci = t_value * std_dev / np.sqrt(n) + + # store results + lambda_points.append(lambda_val) + losses.append(mean_loss) + ci_lower.append(mean_loss - ci) + ci_upper.append(mean_loss + ci) + + print(f"C={c}, λ={lambda_val:.4f}, Mean Loss Rate={mean_loss:.2%} ± {ci:.2f}") + elif n > 0: + print(f"λ={lambda_val:.4f} skipped - only {n} successful run(s)") + continue + else: + print(f"Stopped at λ={lambda_val:.4f} - no successful run") + break + + ax_loss.plot(lambda_points, losses, label=f'C={c}') + ax_loss.fill_between(lambda_points, ci_lower, ci_upper, alpha=0.2) + + # determine optimal C for lamba = 1 if mean_at_lambda1: sorted_C = sorted(mean_at_lambda1.items(), key=lambda item: item[1][0]) @@ -224,91 +303,15 @@ def simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, l else: print("\nNo valid λ=1 data for any C.") + #plot curves + if len(lambda_vals) > 1: + ax_loss.set_xlim(left=2) + ax_loss.set_xlabel('Arrival Rate (λ)') + ax_loss.set_ylabel('Mean Loss Rate') + ax_loss.set_title(f'Mean Loss Rate vs Arrival Rate ({num_runs} runs, {int(confidence_level*100)}% CI)') + ax_loss.legend() + ax_loss.grid(True) - # plot curves - if len(lambda_vals)>1: - ax_rt.set_xlim(left=0) - ax_rt.set_xlabel('Arrival Rate (λ)') - ax_rt.set_ylabel('Mean Response Time') - ax_rt.set_title(f'Mean Response Time vs Arrival Rate ({num_runs} runs, {int(confidence_level*100)}% CI)') - ax_rt.legend() - ax_rt.grid(True) - - return loss_data - -def simulation_wrapper2(simulation_time, num_runs, min_runs, confidence_level, lambda_vals, max_loss_value=0.1): - fig, (ax_rt, ax_loss) = plt.subplots(2,1, figsize=(12,10), gridspec_kw={'height_ratios': [4, 3], 'hspace': 0.4}, sharex=False) - loss_data = simulation_wrapper1(simulation_time, num_runs, min_runs, confidence_level, lambda_vals, ax_rt) - - for c, (lams, losses) in loss_data.items(): - lambda_vals2 = [] - if len(lams)==0: - print(f"C={c}: no data at all.") - continue - pos_idx = np.where(losses > 0)[0] - if pos_idx.size > 0: - first_lambda = lams[pos_idx[0]] - else: - first_lambda = None - last_lambda = lams[-1] - - print(f"C={c:>1} first λ with loss > 0: " - f"{first_lambda if first_lambda is not None else 'never'}; " - f" last λ with data: {last_lambda}") - - lambda_vals2.extend([last_lambda + i/1000 for i in range(-100,51)]) - - - - lambda_points = [] - losses = [] - ci_lower = [] - ci_upper = [] - print(f"\nProcessing C={c}") - - with Pool() as pool: - for lambda_val in lambda_vals2: - args_list = [(c, lambda_val, simulation_time, max_loss_value) for _ in range(num_runs)] - results = pool.map(run_single_simulation, args_list) - - successful_results = [res for res in results if res is not None] - loss_rates = [res[1] for res in successful_results] - - n = len(loss_rates) - - if n >= min_runs: - # statistics - mean_loss = np.mean(loss_rates) - std_dev = np.std(loss_rates, ddof=1) - - # confidence interval - t_value = t.ppf((1 + confidence_level)/2, n-1) - ci = t_value * std_dev / np.sqrt(n) - - # store results - lambda_points.append(lambda_val) - losses.append(mean_loss) - ci_lower.append(mean_loss - ci) - ci_upper.append(mean_loss + ci) - - print(f"C={c}, λ={lambda_val:.4f}, Mean Loss Rate={mean_loss:.2%} ± {ci:.2f}") - elif n > 0: - print(f"λ={lambda_val:.4f} skipped - only {n} successful run(s)") - continue - else: - print(f"Stopped at λ={lambda_val:.4f} - no successful run") - break - - ax_loss.plot(lambda_points, losses, label=f'C={c}') - ax_loss.fill_between(lambda_points, ci_lower, ci_upper, alpha=0.2) - - ax_loss.set_xlim(left=2) - ax_loss.set_xlabel('Arrival Rate (λ)') - ax_loss.set_ylabel('Mean Loss Rate') - ax_loss.set_title(f'Mean Loss Rate vs Arrival Rate ({num_runs} runs, {int(confidence_level*100)}% CI)') - ax_loss.legend() - ax_loss.grid(True) - - plt.show() + plt.show()