import QuantLib as ql
import numpy as np
import matplotlib.pyplot as plt

# --- UNIT 1: Instrument Setup ---
def create_30y_swap(today, rate_quote):
    ql.Settings.instance().evaluationDate = today
    
    # We use a Relinkable handle to swap curves per path/time-step
    yield_handle = ql.RelinkableYieldTermStructureHandle()
    yield_handle.linkTo(ql.FlatForward(today, rate_quote, ql.Actual365Fixed()))
    
    calendar = ql.TARGET()
    settle_date = calendar.advance(today, 2, ql.Days)
    maturity_date = calendar.advance(settle_date, 30, ql.Years)
    
    index = ql.Euribor6M(yield_handle)
    
    fixed_schedule = ql.Schedule(settle_date, maturity_date, ql.Period(ql.Annual), 
                                 calendar, ql.ModifiedFollowing, ql.ModifiedFollowing, 
                                 ql.DateGeneration.Forward, False)
    float_schedule = ql.Schedule(settle_date, maturity_date, ql.Period(ql.Semiannual), 
                                 calendar, ql.ModifiedFollowing, ql.ModifiedFollowing, 
                                 ql.DateGeneration.Forward, False)
    
    swap = ql.VanillaSwap(ql.VanillaSwap.Payer, 1e6, fixed_schedule, rate_quote, 
                          ql.Thirty360(ql.Thirty360.BondBasis), float_schedule, 
                          index, 0.0, ql.Actual360())
    
    # Pre-calculate all fixing dates needed for the life of the swap
    fixing_dates = [index.fixingDate(d) for d in float_schedule]
    
    return swap, yield_handle, index, fixing_dates

# --- UNIT 2: The Simulation Loop ---
def run_simulation(n_paths=50):
    today = ql.Date(27, 1, 2026)
    swap, yield_handle, index, fixing_dates = create_30y_swap(today, 0.03)
    
    # HW Parameters: a=mean reversion, sigma=vol
    model = ql.HullWhite(yield_handle, 0.01, 0.01)
    process = ql.HullWhiteProcess(yield_handle, 0.01, 0.01)
    
    times = np.arange(0, 31, 1.0) # Annual buckets
    grid = ql.TimeGrid(times)
    
    rng = ql.GaussianRandomSequenceGenerator(ql.UniformRandomSequenceGenerator(
        len(grid)-1, ql.UniformRandomGenerator()))
    seq = ql.GaussianMultiPathGenerator(process, grid, rng, False)
    
    npv_matrix = np.zeros((n_paths, len(times)))

    for i in range(n_paths):
        path = seq.next().value()[0]
        # 1. Clear previous path's fixings to avoid data pollution
        ql.IndexManager.instance().clearHistories()
        
        for j, t in enumerate(times):
            if t >= 30: continue
            
            eval_date = ql.TARGET().advance(today, int(t), ql.Years)
            ql.Settings.instance().evaluationDate = eval_date
            
            # 2. Update the curve with simulated short rate rt
            rt = path[j]
            # Use pillars up to 35Y to avoid extrapolation crashes
            tenors = [0, 1, 2, 5, 10, 20, 30, 35]
            dates = [ql.TARGET().advance(eval_date, y, ql.Years) for y in tenors]
            discounts = [model.discountBond(t, t + y, rt) for y in tenors]
            
            sim_curve = ql.DiscountCurve(dates, discounts, ql.Actual365Fixed())
            sim_curve.enableExtrapolation()
            yield_handle.linkTo(sim_curve)
            
            # 3. MANUAL FIXING INJECTION
            # For every fixing date that has passed or is 'today'
            for fd in fixing_dates:
                if fd <= eval_date:
                    # Direct manual calculation to avoid index.fixing() error
                    # We calculate the forward rate for the 6M period starting at fd
                    start_date = fd
                    end_date = ql.TARGET().advance(start_date, 6, ql.Months)
                    
                    # Manual Euribor rate formula: (P1/P2 - 1) / dt
                    p1 = sim_curve.discount(start_date)
                    p2 = sim_curve.discount(end_date)
                    dt = ql.Actual360().yearFraction(start_date, end_date)
                    fwd_rate = (p1 / p2 - 1.0) / dt
                    
                    index.addFixing(fd, fwd_rate, True) # Force overwrite
            
            # 4. Valuation
            swap.setPricingEngine(ql.DiscountingSwapEngine(yield_handle))
            npv_matrix[i, j] = max(swap.NPV(), 0)
            
        ql.Settings.instance().evaluationDate = today

    return times, npv_matrix

# --- UNIT 3: Visualization ---
times, npv_matrix = run_simulation(n_paths=100)
ee = np.mean(npv_matrix, axis=0)

plt.figure(figsize=(10, 5))
plt.plot(times, ee, lw=2, label="Expected Exposure (EE)")
plt.fill_between(times, ee, alpha=0.3)
plt.title("30Y Swap EE Profile - Hull White")
plt.xlabel("Years"); plt.ylabel("Exposure"); plt.grid(True); plt.legend()
plt.show()