From d651ad346ba38a25556d274c4fde09e9cbdcc873 Mon Sep 17 00:00:00 2001 From: Tom Body Date: Thu, 16 Nov 2023 10:59:34 -0500 Subject: [PATCH] Implement model to calculate edge impurity conc. --- cfspopcon/algorithms/__init__.py | 18 +--- .../algorithms/edge_impurity_concentration.py | 71 ++++++++++++++ cfspopcon/algorithms/single_functions.py | 5 + .../algorithms/two_point_model_fixed_tet.py | 7 +- cfspopcon/formulas/__init__.py | 51 +++++----- .../scrape_off_layer_model/__init__.py | 3 + .../edge_impurity_concentration.py | 87 ++++++++++++++++++ cfspopcon/helpers.py | 2 + cfspopcon/named_options.py | 2 + cfspopcon/unit_handling/default_units.py | 6 +- example_cases/SPARC_PRD/input.yaml | 7 ++ tests/regression_results/PRD.json | 26 ++++++ tests/regression_results/SPARC_PRD_result.nc | Bin 2946296 -> 2956732 bytes 13 files changed, 239 insertions(+), 46 deletions(-) create mode 100644 cfspopcon/algorithms/edge_impurity_concentration.py create mode 100644 cfspopcon/formulas/scrape_off_layer_model/edge_impurity_concentration.py diff --git a/cfspopcon/algorithms/__init__.py b/cfspopcon/algorithms/__init__.py index 8a915b7d..d0db085e 100644 --- a/cfspopcon/algorithms/__init__.py +++ b/cfspopcon/algorithms/__init__.py @@ -6,6 +6,7 @@ from .beta import calc_beta from .composite_algorithm import predictive_popcon from .core_radiated_power import calc_core_radiated_power +from .edge_impurity_concentration import calc_edge_impurity_concentration from .extrinsic_core_radiator import calc_extrinsic_core_radiator from .fusion_gain import calc_fusion_gain from .geometry import calc_geometry @@ -40,6 +41,7 @@ Algorithms["two_point_model_fixed_tet"]: two_point_model_fixed_tet, Algorithms["calc_zeff_and_dilution_from_impurities"]: calc_zeff_and_dilution_from_impurities, Algorithms["use_LOC_tau_e_below_threshold"]: use_LOC_tau_e_below_threshold, + Algorithms["calc_edge_impurity_concentration"]: calc_edge_impurity_concentration, **SINGLE_FUNCTIONS, } @@ -53,22 +55,6 @@ def get_algorithm(algorithm: Union[Algorithms, str]) -> Union[Algorithm, Composi __all__ = [ - "calc_beta", - "calc_core_radiated_power", - "calc_extrinsic_core_radiator", - "calc_fusion_gain", - "calc_geometry", - "calc_heat_exhaust", - "calc_ohmic_power", - "calc_peaked_profiles", - "calc_plasma_current_from_q_star", - "calc_power_balance_from_tau_e", - "predictive_popcon", - "calc_q_star_from_plasma_current", - "two_point_model_fixed_fpow", - "two_point_model_fixed_qpart", - "two_point_model_fixed_tet", - "calc_zeff_and_dilution_from_impurities", "ALGORITHMS", "get_algorithm", ] diff --git a/cfspopcon/algorithms/edge_impurity_concentration.py b/cfspopcon/algorithms/edge_impurity_concentration.py new file mode 100644 index 00000000..8aea61e7 --- /dev/null +++ b/cfspopcon/algorithms/edge_impurity_concentration.py @@ -0,0 +1,71 @@ +"""Run the two point model with a fixed sheath entrance temperature.""" + + +from ..atomic_data import read_atomic_data +from ..formulas.scrape_off_layer_model import build_L_int_integrator, calc_required_edge_impurity_concentration +from ..named_options import Impurity +from ..unit_handling import Unitfull, convert_to_default_units, ureg +from .algorithm_class import Algorithm + +RETURN_KEYS = [ + "edge_impurity_concentration", +] + + +def run_calc_edge_impurity_concentration( + edge_impurity_species: Impurity, + q_parallel: Unitfull, + SOL_power_loss_fraction: Unitfull, + target_electron_temp: Unitfull, + upstream_electron_temp: Unitfull, + upstream_electron_density: Unitfull, + kappa_e0: Unitfull, + lengyel_overestimation_factor: Unitfull, + reference_electron_density: Unitfull = 1.0 * ureg.n20, + reference_ne_tau: Unitfull = 1.0 * ureg.n20 * ureg.ms, +) -> dict[str, Unitfull]: + """Calculate the impurity concentration required to cool the scrape-off-layer using the Lengyel model. + + Args: + edge_impurity_species: :term:`glossary link` + reference_electron_density: :term:`glossary link` + reference_ne_tau: :term:`glossary link` + q_parallel: :term:`glossary link` + SOL_power_loss_fraction: :term:`glossary link` + target_electron_temp: :term:`glossary link` + upstream_electron_temp: :term:`glossary link` + upstream_electron_density: :term:`glossary link` + kappa_e0: :term:`glossary link` + lengyel_overestimation_factor: :term:`glossary link` + + Returns: + :term:`edge_impurity_concentration` + """ + atomic_data = read_atomic_data() + + L_int_integrator = build_L_int_integrator( + atomic_data=atomic_data, + impurity_species=edge_impurity_species, + reference_electron_density=reference_electron_density, + reference_ne_tau=reference_ne_tau, + ) + + edge_impurity_concentration = calc_required_edge_impurity_concentration( + L_int_integrator=L_int_integrator, + q_parallel=q_parallel, + SOL_power_loss_fraction=SOL_power_loss_fraction, + target_electron_temp=target_electron_temp, + upstream_electron_temp=upstream_electron_temp, + upstream_electron_density=upstream_electron_density, + kappa_e0=kappa_e0, + lengyel_overestimation_factor=lengyel_overestimation_factor, + ) + + local_vars = locals() + return {key: convert_to_default_units(local_vars[key], key) for key in RETURN_KEYS} + + +calc_edge_impurity_concentration = Algorithm( + function=run_calc_edge_impurity_concentration, + return_keys=RETURN_KEYS, +) diff --git a/cfspopcon/algorithms/single_functions.py b/cfspopcon/algorithms/single_functions.py index 364b6f19..6de5cde7 100644 --- a/cfspopcon/algorithms/single_functions.py +++ b/cfspopcon/algorithms/single_functions.py @@ -71,5 +71,10 @@ return_keys=["plasma_stored_energy"], name="calc_plasma_stored_energy", ) +calc_upstream_electron_density = Algorithm.from_single_function( + lambda nesep_over_nebar, average_electron_density: nesep_over_nebar * average_electron_density, + return_keys=["upstream_electron_density"], + name="calc_upstream_electron_density", +) SINGLE_FUNCTIONS = {Algorithms[key]: val for key, val in locals().items() if isinstance(val, Algorithm)} diff --git a/cfspopcon/algorithms/two_point_model_fixed_tet.py b/cfspopcon/algorithms/two_point_model_fixed_tet.py index 587c2a80..3346d9d1 100644 --- a/cfspopcon/algorithms/two_point_model_fixed_tet.py +++ b/cfspopcon/algorithms/two_point_model_fixed_tet.py @@ -21,8 +21,7 @@ def run_two_point_model_fixed_tet( target_electron_temp: Unitfull, q_parallel: Unitfull, parallel_connection_length: Unitfull, - average_electron_density: Unitfull, - nesep_over_nebar: Unitfull, + upstream_electron_density: Unitfull, toroidal_flux_expansion: Unitfull, fuel_average_mass_number: Unitfull, kappa_e0: Unitfull, @@ -35,7 +34,7 @@ def run_two_point_model_fixed_tet( q_parallel: :term:`glossary link` parallel_connection_length: :term:`glossary link` average_electron_density: :term:`glossary link` - nesep_over_nebar: :term:`glossary link` + upstream_electron_density: :term:`glossary link` toroidal_flux_expansion: :term:`glossary link` fuel_average_mass_number: :term:`glossary link` kappa_e0: :term:`glossary link` @@ -48,7 +47,7 @@ def run_two_point_model_fixed_tet( target_electron_temp=target_electron_temp, parallel_heat_flux_density=q_parallel, parallel_connection_length=parallel_connection_length, - upstream_electron_density=nesep_over_nebar * average_electron_density, + upstream_electron_density=upstream_electron_density, toroidal_flux_expansion=toroidal_flux_expansion, fuel_average_mass_number=fuel_average_mass_number, kappa_e0=kappa_e0, diff --git a/cfspopcon/formulas/__init__.py b/cfspopcon/formulas/__init__.py index 68239238..d32b4c55 100644 --- a/cfspopcon/formulas/__init__.py +++ b/cfspopcon/formulas/__init__.py @@ -1,6 +1,6 @@ """Formulas for POPCONs analysis.""" -from . import energy_confinement_time_scalings, fusion_reaction_data, plasma_profile_data, radiated_power +from . import energy_confinement_time_scalings, fusion_reaction_data, plasma_profile_data, radiated_power, scrape_off_layer_model from .average_fuel_ion_mass import calc_fuel_average_mass_number from .beta import calc_beta_normalised, calc_beta_poloidal, calc_beta_toroidal, calc_beta_total from .confinement_regime_threshold_powers import ( @@ -41,54 +41,55 @@ ) __all__ = [ - "energy_confinement_time_scalings", - "fusion_reaction_data", - "calc_density_peaking", - "calc_effective_collisionality", - "plasma_profile_data", - "calc_fuel_average_mass_number", + "calc_1D_plasma_profiles", "calc_B_pol_omp", "calc_B_tor_omp", "calc_beta_normalised", "calc_beta_poloidal", "calc_beta_toroidal", - "calc_troyon_limit", - "calc_greenwald_density_limit", "calc_beta_total", "calc_bootstrap_fraction", + "calc_bremsstrahlung_radiation", + "calc_change_in_dilution", + "calc_change_in_zeff", + "calc_confinement_transition_threshold_power", "calc_current_relaxation_time", + "calc_density_peaking", + "calc_effective_collisionality", "calc_f_shaping", + "calc_fuel_average_mass_number", "calc_fusion_power", + "calc_greenwald_density_limit", "calc_greenwald_fraction", + "calc_impurity_charge_state", + "calc_impurity_radiated_power_mavrin_coronal", + "calc_impurity_radiated_power_mavrin_noncoronal", + "calc_impurity_radiated_power_post_and_jensen", + "calc_impurity_radiated_power_radas", + "calc_impurity_radiated_power", "calc_LH_transition_threshold_power", "calc_LI_transition_threshold_power", "calc_loop_voltage", - "calc_resistivity_trapped_enhancement", "calc_neoclassical_loop_resistivity", "calc_neutron_flux_to_walls", "calc_normalised_collisionality", - "calc_1D_plasma_profiles", "calc_ohmic_power", "calc_peak_pressure", "calc_plasma_current", "calc_plasma_surface_area", "calc_plasma_volume", + "calc_q_star", + "calc_resistivity_trapped_enhancement", "calc_rho_star", "calc_Spitzer_loop_resistivity", - "calc_q_star", - "calc_triple_product", - "calc_tau_e_and_P_in_from_scaling", - "thermal_calc_gain_factor", - "calc_confinement_transition_threshold_power", - "calc_change_in_zeff", - "calc_change_in_dilution", - "calc_bremsstrahlung_radiation", "calc_synchrotron_radiation", - "calc_impurity_radiated_power_post_and_jensen", - "calc_impurity_radiated_power_radas", - "calc_impurity_radiated_power", - "calc_impurity_radiated_power_mavrin_coronal", - "calc_impurity_radiated_power_mavrin_noncoronal", + "calc_tau_e_and_P_in_from_scaling", + "calc_triple_product", + "calc_troyon_limit", + "energy_confinement_time_scalings", + "fusion_reaction_data", + "plasma_profile_data", "radiated_power", - "calc_impurity_charge_state", + "scrape_off_layer_model", + "thermal_calc_gain_factor", ] diff --git a/cfspopcon/formulas/scrape_off_layer_model/__init__.py b/cfspopcon/formulas/scrape_off_layer_model/__init__.py index e58586a5..5293a841 100644 --- a/cfspopcon/formulas/scrape_off_layer_model/__init__.py +++ b/cfspopcon/formulas/scrape_off_layer_model/__init__.py @@ -2,6 +2,7 @@ These are mostly based on the two-point-model, from :cite:`stangeby_2018`. """ +from .edge_impurity_concentration import build_L_int_integrator, calc_required_edge_impurity_concentration from .lambda_q import calc_lambda_q from .parallel_heat_flux_density import calc_parallel_heat_flux_density from .solve_target_first_two_point_model import solve_target_first_two_point_model @@ -12,4 +13,6 @@ "calc_parallel_heat_flux_density", "calc_lambda_q", "solve_target_first_two_point_model", + "calc_required_edge_impurity_concentration", + "build_L_int_integrator", ] diff --git a/cfspopcon/formulas/scrape_off_layer_model/edge_impurity_concentration.py b/cfspopcon/formulas/scrape_off_layer_model/edge_impurity_concentration.py new file mode 100644 index 00000000..23c3b710 --- /dev/null +++ b/cfspopcon/formulas/scrape_off_layer_model/edge_impurity_concentration.py @@ -0,0 +1,87 @@ +"""Lengyel model to compute the edge impurity concentration.""" +from typing import Callable + +import numpy as np +import xarray as xr +from scipy.interpolate import InterpolatedUnivariateSpline + +from ...named_options import Impurity +from ...unit_handling import Unitfull, convert_units, magnitude, ureg, wraps_ufunc + + +def build_L_int_integrator( + atomic_data: dict[Impurity, Unitfull], + impurity_species: Impurity, + reference_electron_density: Unitfull, + reference_ne_tau: Unitfull, +) -> Callable[[Unitfull, Unitfull], Unitfull]: + r"""Build an interpolator to calculate L_{int}$ between arbitrary temperature points. + + $L_int = \int_a^b L_z(T_e) sqrt(T_e) dT_e$ where $L_z$ is a cooling curve for an impurity species. + This is used in the calculation of the radiated power associated with a given impurity. + + Args: + atomic_data: :term:`glossary link` + impurity_species: [] :term:`glossary link` + reference_electron_density: [n20] :term:`glossary link` + reference_ne_tau: [n20 * ms] :term:`glossary link` + + Returns: + :term:`L_int_interpolator` + """ + if isinstance(impurity_species, xr.DataArray): + impurity_species = impurity_species.item() + + Lz_curve = atomic_data[impurity_species].noncoronal_electron_emission_prefactor.sel( + dim_log_electron_density=np.log10(magnitude(convert_units(reference_electron_density, ureg.m**-3))), + dim_log_ne_tau=np.log10(magnitude(convert_units(reference_ne_tau, ureg.m**-3 * ureg.s))), + ) + + log_electron_temp = Lz_curve.dim_log_electron_temperature + electron_temp = np.power(10, log_electron_temp) + Lz_sqrt_Te = Lz_curve * np.sqrt(electron_temp) + + interpolator = InterpolatedUnivariateSpline(electron_temp, magnitude(Lz_sqrt_Te)) + def L_int(start_temp, stop_temp): + return interpolator.integral(start_temp, stop_temp) + + return wraps_ufunc( + input_units=dict(start_temp=ureg.eV, stop_temp=ureg.eV), return_units=dict(L_int=ureg.W * ureg.m**3 * ureg.eV**1.5) + )(L_int) + + +def calc_required_edge_impurity_concentration( + L_int_integrator: Callable[[Unitfull, Unitfull], Unitfull], + q_parallel: Unitfull, + SOL_power_loss_fraction: Unitfull, + target_electron_temp: Unitfull, + upstream_electron_temp: Unitfull, + upstream_electron_density: Unitfull, + kappa_e0: Unitfull, + lengyel_overestimation_factor: Unitfull, +) -> Unitfull: + """Calculate the relative concentration of an edge impurity required to achieve a given SOL power loss fraction. + + N.b. this function does not ensure consistency of the calculated impurity concentration + with the parallel temperature profile. You may wish to implement an iterative solver + to find a consistent set of L_parallel, T_t and T_u. + + Args: + L_int_integrator: :term:`glossary link` + q_parallel: :term:`glossary link` + SOL_power_loss_fraction: :term:`glossary link` + target_electron_temp: :term:`glossary link` + upstream_electron_temp: :term:`glossary link` + upstream_electron_density: :term:`glossary link` + kappa_e0: :term:`glossary link` + lengyel_overestimation_factor: :term:`glossary link` + + Returns: + :term:`impurity_concentration` + """ + L_int = L_int_integrator(target_electron_temp, upstream_electron_temp) + + numerator = q_parallel**2 - ((1.0 - SOL_power_loss_fraction) * q_parallel) ** 2 + denominator = 2.0 * kappa_e0 * (upstream_electron_density * upstream_electron_temp) ** 2 * L_int + + return numerator / denominator / lengyel_overestimation_factor diff --git a/cfspopcon/helpers.py b/cfspopcon/helpers.py index cd03cb95..f349affc 100644 --- a/cfspopcon/helpers.py +++ b/cfspopcon/helpers.py @@ -33,6 +33,8 @@ def convert_named_options(key: str, val: Any) -> Any: # noqa: PLR0911, PLR0912 return make_impurities_array(list(val.keys()), list(val.values())) elif key == "core_radiator": return Impurity[val] + elif key == "edge_impurity_species": + return Impurity[val] elif key == "lambda_q_scaling": return LambdaQScaling[val] elif key == "SOL_momentum_loss_function": diff --git a/cfspopcon/named_options.py b/cfspopcon/named_options.py index f905efc3..56609254 100644 --- a/cfspopcon/named_options.py +++ b/cfspopcon/named_options.py @@ -40,6 +40,8 @@ class Algorithms(Enum): calc_P_SOL = auto() use_LOC_tau_e_below_threshold = auto() calc_plasma_stored_energy = auto() + calc_edge_impurity_concentration = auto() + calc_upstream_electron_density = auto() class ProfileForm(Enum): diff --git a/cfspopcon/unit_handling/default_units.py b/cfspopcon/unit_handling/default_units.py index 6219eb58..6f348b3e 100644 --- a/cfspopcon/unit_handling/default_units.py +++ b/cfspopcon/unit_handling/default_units.py @@ -131,6 +131,10 @@ vertical_minor_radius="m", z_effective="", zeff_change_from_core_rad="", + edge_impurity_species=None, + lengyel_overestimation_factor="", + upstream_electron_density="n19", + edge_impurity_concentration="", ) @@ -252,7 +256,7 @@ def convert_to_default_units(value: Quantity, key: str) -> Quantity: def convert_to_default_units(value: Union[float, Quantity, xr.DataArray], key: str) -> Union[float, Quantity, xr.DataArray]: """Convert an array or scalar to default units.""" unit = DEFAULT_UNITS[key] - if unit is None or unit == "": + if unit is None: return value elif isinstance(value, (xr.DataArray, Quantity)): return convert_units(value, unit) diff --git a/example_cases/SPARC_PRD/input.yaml b/example_cases/SPARC_PRD/input.yaml index a34044ab..09ed7c94 100644 --- a/example_cases/SPARC_PRD/input.yaml +++ b/example_cases/SPARC_PRD/input.yaml @@ -42,7 +42,9 @@ algorithms: - calc_P_SOL - calc_average_total_pressure - calc_heat_exhaust + - calc_upstream_electron_density - two_point_model_fixed_tet + - calc_edge_impurity_concentration # Finally, we calculate several parameters which aren't used in other # calculations, but which are useful for characterizing operational # points. These can be used later when masking inaccessible operational @@ -166,3 +168,8 @@ fraction_of_P_SOL_to_divertor: 0.6 kappa_e0: 2600.0 # Calculate P_rad_SOL such that the target electron temperature in eV reaches this value (for FixedTargetElectronTemp) target_electron_temp: 25.0 + +# Which species should be injected to cool the scrape-off-layer? +edge_impurity_species: Argon +# Factor applied to the result of the Lengyel model to give an absolute low-Z impurity concentration. +lengyel_overestimation_factor: 4.3 diff --git a/tests/regression_results/PRD.json b/tests/regression_results/PRD.json index e32ad68e..b06cb83d 100644 --- a/tests/regression_results/PRD.json +++ b/tests/regression_results/PRD.json @@ -259,6 +259,18 @@ }, "data": 25.0 }, + "edge_impurity_species": { + "dims": [], + "attrs": {}, + "data": "Argon" + }, + "lengyel_overestimation_factor": { + "dims": [], + "attrs": { + "units": "dimensionless" + }, + "data": 4.3 + }, "average_electron_density": { "dims": [], "attrs": { @@ -999,6 +1011,13 @@ }, "data": 5941.561499571544 }, + "upstream_electron_density": { + "dims": [], + "attrs": { + "units": "_1e19_per_cubic_metre" + }, + "data": 7.5 + }, "upstream_electron_temp": { "dims": [], "attrs": { @@ -1034,6 +1053,13 @@ }, "data": 0.41293376247700014 }, + "edge_impurity_concentration": { + "dims": [], + "attrs": { + "units": "dimensionless" + }, + "data": 0.08444946529174713 + }, "greenwald_fraction": { "dims": [], "attrs": { diff --git a/tests/regression_results/SPARC_PRD_result.nc b/tests/regression_results/SPARC_PRD_result.nc index c00afd126bb0b98cb54caeee15726601060bbacb..3809c7d692bdb5248f1e75baf570ba43554a0f5d 100644 GIT binary patch delta 11328 zcmch72{cySA1<#U6`2x=A|atcn$?5SBpJ#SDheqgB=c)1r6e3OCZS}gkSX*<)T^k> znKIAw>?PGbzJLE(_uh5ax@+CF?q2ITzqR*y_VDcAIeR*X-uIF*KQDtZ;5(yp;Jhtm zH=VFy^rzgV2T}Yft#l%Qa+|W1P6#o1QGDn`43iI~o=&W1cAy+}Bx0CtDb{piIg1Y^ z#F2<%ai-+b2?18;C9RGO3=GT~XDlug0h~rm5-fCLp4FKm?MOH-a$cfAU$|>xYHv(% zEi$7xI1=aC{3(%iVv@~^Qb=+ZyHmOsnz4ddw#1Uc?nt;UxlP$ZC)k$Wrsz5n!Ambw zT=8r_ksGJBJTth`dJ{9J<1q!^J6oR46PLz{or6Zxp?N8Z9#<*=M=N*XyJnj@4cwv7v%~`T(ZhYOv&@~VxUy{)OS)frPu({uK=)|bzmH@$2}R>44g78d`P|7a(XBka7; zlk%%CQsQ!mnAJX%c+$sfY$+c}udZ>XjF4UxFb`UGa{`3!L@F| zt|M&L>62E5C|jo~Br7DpQ2zha6=&>hOwU=E+7aK@DUw75LPUI_paX=xxEg5_5UFmZCdCuyUIy7dd`@q~udLOg}2 zma-si2BEm&25Iw&k_`*##l$h`8;fLx%E*STA?6ng5TnwjY_dYtWL(gKXqFZw++}o0 zViQp(vru3sp&?85gwP;S{;!=RG&YhwA~a8wEm&b@;_xQ2U}1KmY!lfa(r1dB7t)1@ z{LN&3;dO-kmW6Z~;^~$Jo;)EVN0uk7Kt3|di2^y}b+W<-$jXH+O)r}}nOdH)b}+Ry zwX?UdGPbuM-D_ri&feOVXxny{%yf)cw{@ZLlSIr`r$w^DhGa%|B2b=}NRU5Hrk^LS z%P&-BNzmm>csSZ;En=)!rQ7PtRJF!ldlG0r|5InP%y56{!}=Vy;Q zAyD%W$Cc@|h_ZH#Et#moW7+9s4whg%%V5{pcV8HUak1CJ<*37Z#r`-dDO`+ZaS$mTQhupSdamab}D4>38c zKD1&c$74Yi89e^I{=MK*{CX38MJIUYzrO5GnM4rPZ^eryCS5#ur&ub*x(AQDdzB4j zm*BZ=N{nD^7&0^!XqFZZ$USt|as0|d6cA=1D|Wv_vCFXQkAw0kY1{C2anmY%dhep; z;-Z0y&>PF1^7f+If4*bymt=e`J#?opB@~UvHv~G`e@4qnV*l+j@caK0rNly_EW3D4Ja6ingv@f>t4U{ubj6 zXkYC!yxut(odl0oeap^reBQ)}fwZ1AOcjL~_Oe|PPG!Sr;gXaB zo7WiEch5L-)&i4Dv$dDzJ;GF7o260*ZAaOCk#B_E==`0D}-T} zeI6T24q-GuY|m%GLm1aq3|Ve@1(R$on#cBECl$YURjn-luTqBJ0XDR4XC^i~;wJwF zA=)!Rgy?E6rpCQl+Uz!-9t+fv%hy7xp9uG+fAo;~Fvxdy#20zz_ZOS?IHO?X zG_}eo&m6_KjR*8zYNEs^KSX$d@_!*{Y1NjK(xC38FVT=eDFbc3%Y&O7l%phLGSd{?J1}0 zFc5q#Wnb_?44c%>3zV(GXgYPp(q^HT7*`ASdnvUKlPqf7anCnksv-54%7&}|n_~OC z)nV6)cHp+`A;VO8HH2D9uXr(Oh1iuTzqkvRA~7dz9V@jKX{vhZ7d0}GWf;BZ^6&@b zFE^R1?z@jdg}oMzA4^eOvPtfv!7xgw>E4~MMU3!?!P@DBMvaD_O9!lItstY+FGw-UWEOMiu$>~8*!AKZP$m7Bk7sgE%#Ir zq!s1`I7$p5>%F+M5N!?eEzF(W)gGWw-=p-=VqX;3a+=g{hP=<$6t9hU82 zPdSH{=A|!-yqeI~*88&7OaL9f^t zUg;=e^i?Ihf>Hp+&F-DJROgM!tunczPPUj9HCFwkZS_a6%1+-{>fKvLa1AXm+jvv~ z0gGg2KE_-?BzFsA@nj(4`Gqf4#tI|(ncA6<11pd&rMbNz%>vm=I76pyc_H63`n`{S zCJNP**yPIQQJkvtw!X*%B{s>Mb}rLGS%_I>Dw8zI3E^ZJasaC0OMxw)Hfn0Ku8GzC zKwZNZ>Y^oSX=svoXUTt;oSojpUD$Gr6Kyp*t!F|-&@uh|Wnl7GbRSI5_GV#6Z;8Gn zjWY}b?x|nX4+LP?va|VAjtNFfg$j=DFUL5Qb#kUY0Fx(%tyu~uFs&NdH=z6Ok7CuP z?dH9o_QG99#>S*a9RV@a_#u<`G(=_>z2h^yig@E#re%8^k=$`W-FT}v(zp3ZJ}TLR zZ00Z~T}2z@U+6sJG`a?bJ9tD~Uu&T_;CbH)gQqAtY*cXN_7uuYU)R=b)I+(WlesAC z8&u9}SXs~xpyt&STh)pCsCy$naV_RCnwAM|YP{Hw7Hab74@NI>w0#kbQJUY6jw$uw zI-3e~s|Ft6lTAWz^~AI1{rAa{D_&gftRRN(4pAG|N@KKLkm{{PPEHS8vtO`0%Z&SQImoY1Ws_qyU@?HU^u&+cNcLVqd&3- zzRO!-yda04?fiL6K6%HXRuqHjG}@)Qt|NaGYc;1^?yVPwhf7Y~eUV!Tyvh)LkgF9@ z9^BeWXYU~XsmDnvoApSZ*VwhLk~)cW@n770;dRKWwc@;BEr9%0|19ROe1t;Aw+9=v zxKVt3w}e^CCD6W0Zhl+bhcd2_^m9tbP<{-7cilgu^3x;Rjg_6K2@EmXd-@jY9%xdQ z-AzX0^bx^4v&E$1jjQYv6VO)WpjmLQ8l5a!=IgY{fy<&^VG-5iBl?(i++jy9K z0cvi%PQA(Ng1Y;ACe7KdqjBoh{ljP0p(Tlr_0gLwv^BjP;ErfQCm$o{U`80aJ$EYV zYxtsXUH&3o(P9iTu~E6=7&0-kCOh0BITB+~)HET2F~J>vt1k8orr4U7X$>)AMw0u> zkq7dB6zeF4gp0F=!}H9&fVK7>2pnP0^Cw0SRh1DvVwsHi7pcj8{xwLh+i-ukR|3-J z?Y=swuR&Hw;L#ta`jGcdv*v`cFA6?U^%=H37eFyj6EW1k479xb&^qHYD6Lq;m_?UF z`9`VJiNADE`SRP^sfr-fxDPegE7zkgxbItsrU@FqN2laP-A2o+Lz4AVb~b#JZ*NPT_SvBF8JD5DjTCBJ6V=@q!%-KOxvJF98BL2MFPFTL zM9b$boO|QA&^9OH(`T%YPW;%+JZOpT{MB#61tZXByF7>z@*IPvHy5dSw_qefA@)X< z3&w7<+o?t>Vd9eiq2^|iqzcr< z{JMHD#Rg6MA$4K`q>p>=nwZ5#p`Am)-s7bYI!|#Q_I{#H*%{ov3@Zi@|Lh7n)=g{hq%i7nw|H^jo?zXkYt`*L(GCbUOHUzx1#| z52a=KgP+{!D>Ytgt09NMs@CY@H5M4TwMe5>>h zf-b?4T%la}NARhQiqwQQdw+PH2u(k;aW#VWuSgyZ+JmTfefiprvxv8-xq3x>F_P8K z*XYEO%j&R2)7vTwk;S+6?oC}TqMSc=(olL& z_4~l?CY1BP8C7wjp)w`&@^t|d)c6^PQMv+AmrrdI_dGO=rULZD!7&`KZreP+7 zQJrX7ux3BTW%KHd%_1-H2F4a{&rtcQ;*h1&XVg6Q zIIYLO2X&oq6vx)-pvhgbS5=AJM<^VX3;y;9?N}>0^Ida9bD3MOcS& zEz?k03U{+RTi+!9S$gX#Ad zemW@5VpjTa6NCMazlz%nwOSc}!0VyC_#*#y1Z}wR^I@1WqE1@P1#xXeJQM$8{v)kO zdc$AYtU?Z4I4i^FM?+>h1J%JTOB=Z#6sx*KN>HFW9sT%yCyK%bTL&FJgQh(tSv!3Y zrT&Xote|tETtZj6e`N_OKhHT0h0UUdwrj_-RyWl1#qKzj5sjuoHr7SQBG7ta?@g__ zO0++Is_W(EgDxSzo^Q0P=qazNX#OIG{`WVvjADLZn34LqV(%7maarV4^7gg`#-DA} z-O)|QWWixy2dM{`{v{N&e*1IGTCzXR+0OV!u_3=Df0%w2ykbpl&I*1<(5g?eVuGI# zC6`()8kmTX;QNd zRVWH>6c@dh44O`yQ{k~Bl!k2SG$_7JEqxT>+FM6#y-YpPW$N|;qWjS+wJn^f z@CpVt!%x0j0K@LB#?f9EFveaWboAUJOi)gk?RX)FskLqeKB?sVKVu!|RLtmU*fsF>K!Y-ZMpK&hL=3aOy#|tUu=l>Uzf^U zjibT2_E$7Feq1VmcfAO&1;tX$OFvRQPKbGAdla-7~GqYOv8jsR{MoBSXOhg z8MO1`z2{8KP@1z{>HeNDa_98=-q98+s+JT}hekSd$upnEsZ~W=Q13gYw}+pGZ~A-n z6p0gP{U9vWB(Vt{t4}oYhw`8+Rar5n`V4x#m#8UCe#bz>$kY9;FEKoSlI^OB62`Qf z6E<#T#e`1p-Zgg3OPfg4Z-1pE0-#mNnu6cDIt?hRNR&jkr$KF7e^WNs@s^$Kiaq%5`Q}?vAZgjyw z`B}f-yIU}_d&KJdsc#r_-pAT{+5rTx%$qXNQ;L(DA1Mo}dQizYSS`={6#LvQh;nzah?U$964^cqH){I2m4Bht&4W z5Jjp!GWA!#6MSoqT#t=DtpUTNVoeU^Wo9UfSiNJQde38jeJsnE5SFwAUKZT@yNxUmn2@0JTv8B`rtupexln<)(`Q0T z=DCpg-E*4%fIL#aW{N4w`XEzdAItiq?~&_buK4ELF%;}r8L+ryD>)r&*VcHD3lL9w zg6-J>lnw}!hb(VUZkN~e@uLo^l*RNcKI);?`~{%{FVqhZh1w5W@U7rTkOoo zMAgvxw>Q3EnkQ$}ylhq$9RO|M2&BC-20C-s1B7E(Xe9yeM~&Q6EF?4GojM6RPzr^r%a@;&b0 zP+t>_qQ}V%kJTT6=D%l0c0m+Me^%5MOpwcJFR@3=+g74VhavPv!aQm{f1N6NbPx?3 z=F1ru-^90HCD&Yv#?iLfouimT4joYebjfx1(S7)sL(cJo=-bIm4O70k34=D9XQe(0 zVPwjT-XZLZvF^=RB-v&$F+cyH)hhtgN+JG~A2Xn5-R3FnUxM$auX#8qS^QCK(8IB8 zbrJ`>4$#j1%zTYNHK&CP^ zFT11RF>)Qov#$m$;Ta1D0;T&*w^S_&>};U60jC!yopgg)M}ImuC}*1{svVS zFUv(v$D%gs>$SR?D`?o1$+F*w13T;QfetGoy13K~@rE_2HNB3nr>3#R_ zq0e-E%+CWp7<|0yT}-zUMyV?uesvzJ!uaC!i3(l@Op0fUT-auZX&c5(K?(_=56EoU zALotl;q?U#LMQ(!&OJUT8rTl6&8MVXo>ZpzAnmN=BbF*|}X zDW52s9{JHJCi3mmn*FGvo|0H(Q-|7|byo(1n9-nRaq*5Q4b9qJd+zb>LYt}D&K$QG zbTnO*KhG(S?%;r3D<1KoFKB~0qoX(m>uH>A{I@W=b8(GnSuVz9#tv@jSWU&`ao*K; zwuWQ+0sV(rM-67z>|=DM%;0>rx{1HZO`w*y=ASd=h1`+*g zrR4#UP+t#n+KsxQSJK6NH3J9`n?HsIfkm>43)1xwNYDDvq|w(02;2^Uf+IK z1I?a>AzOLa(MGKvnSASmj(MYUF`+SZ=SH|S#wMfh4MQorWIDM=+z=AaQBUreF6`vK zbsOXAy}Z*Z+L(0Ut$lG%8PlcCgCkVM9L#F{qm^kTkNIVH!*pC(|0>oP+iBBa058_a zbd!;*2;BQob=7eNM0}dQe5Otnam@4FuA76A7~^)m)bk)xQ}2Fod(wl<9TnCosxIWf zb*N|YV@nhq4q3M2u?LDu&aPB#B6r6hcJ`U=bU>MI*C6#R1W;bkDc1kb9aJS<<|&LU zLhVRH6xSpj4PgZj4Bo9q^D`!U*+CApMVK`;e~>|^(9!v~`Z4Hky)CF!coBUA@iEpJ z?=Ykcd&Unl7`4C0we~GD#%+VM?)rvd@Io^s=bF;bJS zGdoT3A#=N7&zsp&>7Pl{Sf69 zAyXbk{HQ9RI;Sg5LX3vk|uYgfY1)n6mQTca@=!DKY~wDYYo^FVkY3(u+IGC(jo-%QB-x(*CHa!Lia+-Ma1^h zxcrk>hr}D@28^M%k(zMyoF2m`WNvGl8moVS-0QlzO>v^+UNJGr!b=-P&4=y@_mzP* zaCN(k+bGI>x7B=6|Ag}PPy>9HLRIbO+|ftR@OAx*+mhWppfzFQ=`9PJvfXxSVcZrOW?Y(`*oD%w_G+|nSu%yla}XEJ5bFm ze#P9s9$!rkXIi^yqERx*;L-OjXfeE2NN@H-J9BVU=Cd{Ej8rmI6;k#?kHb4Fn?ZH- zM}A0bw* z{E#r}ip@S08Wn&3HnOFQIW9xEX&Hu#ONE7QDkdxXuU$-kO-c30vix zG%Yc1&d@yp!@NY_@|vo7`b{095jB>#$RaUkNs(W`T%DqqCDy8=^nrk|CMhS8w$9V{ z0@gE0`h%Zw)bVj^QZy*BRTH7>e8Ac$ElSAR2yF@+bM;VS%oe4m0tXyT5+ijPdX}rY z7~PTx)z8z4K(RhTUrP8JBJ`ucdP9^-5>|VJ+63bE1U)BV*ym|dVEXJm%1GE6qm&m| zXq=~?Y)jB=5#jcTdK3BfB+Z2|+!3SCi`ec+&<`Q# zoi{Z7aS`dxMNR)R4_nu)ra#$0y6Zzt|7#7_Zmy3dR=VHQ^uMQ9U2kjp({W6>u4?+p zH56UE{bU@j9+gv&J+E=banz%pFCg2ias|KkzQCEoW%s+B8`yHI{vJ-<6PyK%dsKfP zhaUc#lLL(QsrnGbK0emT5yJhd{|JSCRUaedo#ZTG+qTu>Q6;R4^+7{Rb#Wf;4D zHFy-m7x_TOtJuD%`eV=s)c)Jp7*PEwI4-Gmajaiba~Algv%`P+CVyUm;~BLli?#pr zE^LE*X(NaDpxV2I$AkP*jcr`N!k033v42I?1q1@B-p58j%|3uz@Ju^!gT+|0xJSpQ}z=-s@tc;8*v37B=$XFId^_3vhPLQbFA Kbmv;D@AiLQqqG75