import math,random,sys
def arg(i,d):
    return float(sys.argv[i]) if len(sys.argv)>i else d
N_EOL=int(arg(1,0))
VENT_MOYEN=arg(2,6.2)
PV_TRACKERS=arg(3,500.0)
H_TRACKERS=arg(4,2100.0)
PV_FLOTTANT=arg(5,24.0)
RACCORDEMENT=arg(6,320.0)
P_POMPES=arg(7,64.0)
P_PELTON=arg(8,48.0)
STOCK_MAX=arg(9,588.0)
PRIX_PV=arg(10,0.57)
LAT=25.0
H=8760
random.seed(2026)
P_EOL=6.0
FACTEUR_ROTOR=1.09
PERTES_EOL=0.08
SILLAGE=1.018
H_FLOTTANT=H_TRACKERS*0.75
RENDEMENT_POMPE=0.87
RENDEMENT_TURBINE=0.86
PLUIE_GWH=1.9
PRIX=[60]*6+[70]*3+[30]*7+[80]*2+[130]*4+[80]*2
def weibull_quantile(u,mean,k=2.0):
    c=mean/math.gamma(1+1/k)
    return c*(-math.log(1-u))**(1/k)
def normal_cdf(x):
    return 0.5*(1+math.erf(x/math.sqrt(2)))
vent=[]
z=0.0
rho=0.96
for h in range(H):
    z=rho*z+math.sqrt(1-rho*rho)*random.gauss(0,1)
    saison=1+0.15*math.cos(2*math.pi*(h/24-15)/365)
    u=min(max(normal_cdf(z),1e-6),1-1e-6)
    vent.append(weibull_quantile(u,VENT_MOYEN*saison))
COURBE=[(3,0.0),(4,0.05),(5,0.13),(6,0.24),(7,0.39),(8,0.56),(9,0.74),(10,0.89),(11,0.97),(11.5,1.0)]
def p_eolienne(v):
    if v<3 or v>=25: return 0.0
    if v>=11.5: return P_EOL
    for (v1,f1),(v2,f2) in zip(COURBE,COURBE[1:]):
        if v1<=v<v2: return P_EOL*(f1+(f2-f1)*(v-v1)/(v2-v1))
    return P_EOL
eol=[N_EOL*p_eolienne(v*FACTEUR_ROTOR)*(1-PERTES_EOL)*SILLAGE for v in vent]
def pv_forme(h):
    jour=h//24
    heure=h%24
    decl=23.45*math.sin(2*math.pi*(284+jour)/365)
    ha=15*(heure+0.5-12)
    sinalt=math.sin(math.radians(LAT))*math.sin(math.radians(decl))+math.cos(math.radians(LAT))*math.cos(math.radians(decl))*math.cos(math.radians(ha))
    return max(sinalt,0.0)
nuage=[]
c=0.8
for j in range(366):
    c=0.6*c+0.4*random.uniform(0.55,1.0)
    nuage.append(c)
brut=[pv_forme(h)*nuage[h//24] for h in range(H)]
somme=sum(brut)
pv=[min(b/somme*PV_TRACKERS*H_TRACKERS,0.85*PV_TRACKERS)+b/somme*PV_FLOTTANT*H_FLOTTANT for b in brut]
pluie_h=PLUIE_GWH*1000/H
stock=STOCK_MAX*0.5
vendu=turbine=pompe_e=ecrete=prod_pv=prod_eol=recette=0.0
heures_pelton=0
for h in range(H):
    heure=h%24
    prod=eol[h]+pv[h]
    prod_pv+=pv[h];prod_eol+=eol[h]
    stock=min(STOCK_MAX,stock+pluie_h)
    export=min(prod,RACCORDEMENT)
    surplus=prod-export
    pompage=0.0
    if surplus>0 and stock<STOCK_MAX:
        pompage=min(surplus,P_POMPES,(STOCK_MAX-stock)/RENDEMENT_POMPE)
        surplus-=pompage
    nuit=(heure>=23 or heure<6)
    creux=(heure>=5 and heure<17)
    if (nuit or creux) and stock<STOCK_MAX*0.95:
        dispo=min(export,P_POMPES-pompage,(STOCK_MAX-stock)/RENDEMENT_POMPE)
        if dispo>0:
            pompage+=dispo
            export-=dispo
    stock+=pompage*RENDEMENT_POMPE
    pompe_e+=pompage
    ecrete+=surplus
    turb=0.0
    if (heure>=17 or heure<5) and stock>0:
        turb=min(P_PELTON,stock*RENDEMENT_TURBINE,RACCORDEMENT-export)
        stock-=turb/RENDEMENT_TURBINE
        if turb>0: heures_pelton+=1
    turbine+=turb
    vendu+=export+turb
    recette+=(export+turb)*PRIX[heure]
eolien_m=N_EOL*P_EOL*229/240
pv_m=PV_TRACKERS*PRIX_PV+PV_FLOTTANT*(PRIX_PV+0.18)
gouttieres=PV_TRACKERS*8.3/500
couverture=2.0
machines_hydro=67.8*(P_PELTON+P_POMPES)/112
conduites=13.3*max(P_PELTON,P_POMPES)/64
bassins=19.1*STOCK_MAX/588
raccord=22*RACCORDEMENT/240
machines=eolien_m+pv_m+machines_hydro+raccord+0.5
civil=bassins+conduites+2.0+gouttieres+couverture
invest=machines+civil
opex=1.5*(STOCK_MAX/588)**0.5+2.8*N_EOL*P_EOL/240+0.006*(PV_TRACKERS+PV_FLOTTANT)+0.15
annuite_prive=machines*0.0578+civil*0.0442+opex
annuite_public=machines*0.0465+civil*0.0361+opex
print(f"Eolien : {N_EOL} x {P_EOL:.0f} MW a {VENT_MOYEN:.1f} m/s, PV trackers {PV_TRACKERS:.0f} MWc a {H_TRACKERS:.0f} h, PV flottant {PV_FLOTTANT:.0f} MWc, raccordement {RACCORDEMENT:.0f} MW")
print(f"Stock {STOCK_MAX:.0f} MWh, turbines {P_PELTON:.0f} MW, pompes {P_POMPES:.0f} MW")
print(f"Eolien net : {prod_eol/1000:.1f} GWh, PV : {prod_pv/1000:.1f} GWh, pompe : {pompe_e/1000:.1f} GWh, turbine : {turbine/1000:.1f} GWh, ecrete : {ecrete/1000:.1f} GWh")
print(f"Electricite vendue : {vendu/1000:.1f} GWh, dont pilotable : {100*turbine/vendu:.1f} %, turbines en marche : {heures_pelton} h")
print(f"Investissement : {invest:.1f} M euros (eolien {eolien_m:.1f}, PV {pv_m:.1f}, eau {gouttieres+couverture:.1f}, hydraulique {machines_hydro+conduites+bassins+2.0:.1f}, raccordement {raccord:.1f}), soit {invest/vendu*1000:.3f} M euros par GWh annuel")
print(f"LCOE prive / public : {annuite_prive*1e6/vendu:.1f} / {annuite_public*1e6/vendu:.1f} euros/MWh")
print(f"Valeur au profil de prix ensoleille (moyenne {sum(PRIX)/24:.1f}) : {recette/vendu:.1f} euros/MWh, recette {recette/1e6:.1f} M euros, temps de retour {invest/(recette/1e6-opex):.1f} ans")