From 9a7e0e15f69ce7e2aa2b7ff9df998b2b5b937edd Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Tue, 17 Mar 2026 04:49:42 +0900 Subject: [PATCH 01/17] Add tetrapod model implementation in C and Python --- sasmodels/models/tetrapod.c | 60 +++++++++++++++++++++++++++++ sasmodels/models/tetrapod.py | 74 ++++++++++++++++++++++++++++++++++++ 2 files changed, 134 insertions(+) create mode 100644 sasmodels/models/tetrapod.c create mode 100644 sasmodels/models/tetrapod.py diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c new file mode 100644 index 00000000..622e4b83 --- /dev/null +++ b/sasmodels/models/tetrapod.c @@ -0,0 +1,60 @@ +const double A = + 109.5 / 2.0 * M_PI / 180.0; // half of the angle between arms in radians + +static double u_n(int n, double theta, double alpha) { + const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2}; + return cos(A) * cos(theta + phi[n] * 2) + + sin(A) * sin(theta + phi[n] * 2) * cos(alpha - phi[n]); +} + +static double A_n(double q, double u, double L, double R) { + double quL2 = q * u * L * 0.5; + + double term1 = (fabs(quL2) > 1e-12) ? sin(quL2) / quL2 : 1.0; + + double mu = sqrt(fmax(0.0, 1.0 - u * u)); + double qmuR = q * mu * R; + + double term2 = (fabs(qmuR) > 1e-12) ? 2.0 * sas_J1(qmuR) / qmuR : 1.0; + + return term1 * term2; +} + +static double form_volume(double L, double R) { + // V = 4 * \pi * R^2 * L + return 4.0 * M_PI * R * R * L; +} + +static double Iq(double q, double L, double R, double sld_particle, + double sld_solvent) { + double contrast = sld_particle - sld_solvent; + double total = 0.0; + + for (int dtheta = 0; dtheta < GAUSS_N; dtheta++) { + double theta = + M_PI_2 * (GAUSS_Z[dtheta] + 1.0); // map from [-1, 1] to [0, pi] + double w_theta = + GAUSS_W[dtheta] * M_PI_2; // adjust weight for the new range + + double integral_alpha = 0.0; + for (int dalpha = 0; dalpha < GAUSS_N; dalpha++) { + double alpha = + M_PI * (GAUSS_Z[dalpha] + 1.0); // map from [-1, 1] to [0, pi] + double w_alpha = + GAUSS_W[dalpha] * M_PI; // adjust weight for the new range + + double sum_arms = 0.0; + for (int n = 0; n < 4; n++) { + for (int m = 0; m < 4; m++) { + double u = u_n(n, theta, alpha); + double v = u_n(m, theta, alpha); + sum_arms += A_n(q, u, L, R) * A_n(q, v, L, R) * + cos(q * (u - v) * L / 2.0) * sin(theta); + } + } + integral_alpha += sum_arms * w_alpha; + } + total += integral_alpha * w_theta; + } + return contrast * contrast * total * 2 * M_PI; +} \ No newline at end of file diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py new file mode 100644 index 00000000..d0c51b94 --- /dev/null +++ b/sasmodels/models/tetrapod.py @@ -0,0 +1,74 @@ +r""" +Definition +---------- + +Calculates the scattering from a tetrapod-shaped structure. A tetrapod consists +of four cylindrical arms radiating from a central point, oriented along the +(1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions. + +The scattering intensity is calculated as: + +.. math:: + + I(q) = (\Delta \rho)^2 \left[ 4 I_{\text{arm}}(q) + 6 I_{\text{corr}}(q) \right] + +where $\Delta\rho$ is the scattering length density contrast between the +tetrapod and the solvent. + +The arm scattering $I_{\text{arm}}(q)$ is obtained by integrating over all +orientations of each cylindrical arm: + +.. math:: + + I_{\text{arm}}(q) = \int_0^{\pi/4} F_{\text{cylinder}}(q, \theta, L, R)^2 \sin(\theta) d\theta + +where $F_{\text{cylinder}}$ is the form factor of a cylinder. + +The correlation term $I_{\text{corr}}(q)$ accounts for coherent scattering +between different arms. + +Geometry +-------- + +The four arms are positioned along the following directions: +- Arm 1: (1, 1, 1) +- Arm 2: (-1, -1, 1) +- Arm 3: (-1, 1, -1) +- Arm 4: (1, -1, -1) + +Each arm has length $L$ and radius $R$. + +References +---------- + +#. J. S. Pedersen, *J. Appl. Cryst.*, 33 (2000) 488-494 + +Authorship and Verification +---------------------------- + +* **Author:** NIST IGOR/DANSE **Date:** pre 2010 +* **Last Modified by:** Paul Butler **Date:** March 20, 2016 +""" + +import numpy as np +from numpy import cos, inf, pi, sin + +name = "tetrapod" +title = "Tetrapod with four cylindrical arms" +description = """ + Calculates the scattering from a tetrapod structure with four cylindrical + arms radiating from a central point. +""" +category = "shape:cylinder" + +# [ "name", "units", default, [lower, upper], "type", "description"], +parameters = [ + ["length", "Ang", 400, [0, inf], "volume", "Cylindrical arm length"], + ["radius", "Ang", 30, [0, inf], "volume", "Cylindrical arm radius"], + ["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Tetrapod scattering length density"], + ["sld_solvent", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Solvent scattering length density"], +] + +source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] +valid = "radius >= 0.0 && length >= 0.0" +have_Fq = False From 40769a5b2de74053d110c1571b9de735827262c4 Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Tue, 17 Mar 2026 05:37:05 +0900 Subject: [PATCH 02/17] Fix incorrect references to cited papers. --- sasmodels/models/tetrapod.py | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index d0c51b94..f2658268 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -41,13 +41,13 @@ References ---------- -#. J. S. Pedersen, *J. Appl. Cryst.*, 33 (2000) 488-494 +#. Seoki Kyoo Seo *Korean J. Chem. Eng.* 34(2017) 1192-1198 Authorship and Verification ---------------------------- -* **Author:** NIST IGOR/DANSE **Date:** pre 2010 -* **Last Modified by:** Paul Butler **Date:** March 20, 2016 +* **Author:** Yuhei Yamada +* **Last Modified by:** """ import numpy as np @@ -70,5 +70,4 @@ ] source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] -valid = "radius >= 0.0 && length >= 0.0" have_Fq = False From c995274dbea32d1c2f195b7eb7e5ce886190b73f Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Mon, 16 Mar 2026 20:43:08 +0000 Subject: [PATCH 03/17] Initial plan From 4399102559eccf0b2e925dd66ca8423682106ece Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Mon, 16 Mar 2026 21:03:40 +0000 Subject: [PATCH 04/17] Fix u_n formula in tetrapod.c and update docstring Co-authored-by: IndigoCarmine <79552713+IndigoCarmine@users.noreply.github.com> --- sasmodels/models/tetrapod.c | 7 +++-- sasmodels/models/tetrapod.py | 54 +++++++++++++++++++++++++----------- 2 files changed, 42 insertions(+), 19 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index 622e4b83..f2bf7c22 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -3,8 +3,9 @@ const double A = static double u_n(int n, double theta, double alpha) { const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2}; - return cos(A) * cos(theta + phi[n] * 2) + - sin(A) * sin(theta + phi[n] * 2) * cos(alpha - phi[n]); + const double sign[4] = {1.0, -1.0, 1.0, -1.0}; + return sign[n] * cos(A) * cos(theta) + + sin(A) * sin(theta) * cos(alpha - phi[n]); } static double A_n(double q, double u, double L, double R) { @@ -39,7 +40,7 @@ static double Iq(double q, double L, double R, double sld_particle, double integral_alpha = 0.0; for (int dalpha = 0; dalpha < GAUSS_N; dalpha++) { double alpha = - M_PI * (GAUSS_Z[dalpha] + 1.0); // map from [-1, 1] to [0, pi] + M_PI * (GAUSS_Z[dalpha] + 1.0); // map from [-1, 1] to [0, 2*pi] double w_alpha = GAUSS_W[dalpha] * M_PI; // adjust weight for the new range diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index f2658268..54d0cc10 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -6,35 +6,52 @@ of four cylindrical arms radiating from a central point, oriented along the (1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions. -The scattering intensity is calculated as: +The scattering intensity is calculated as a powder average over all +orientations: .. math:: - I(q) = (\Delta \rho)^2 \left[ 4 I_{\text{arm}}(q) + 6 I_{\text{corr}}(q) \right] + I(q) = \frac{(\Delta \rho)^2}{4\pi} + \int_0^{\pi} \int_0^{2\pi} + \left|\sum_{n=1}^{4} F_n(q, \theta, \varphi)\right|^2 + \sin\theta \, d\theta \, d\varphi -where $\Delta\rho$ is the scattering length density contrast between the -tetrapod and the solvent. - -The arm scattering $I_{\text{arm}}(q)$ is obtained by integrating over all -orientations of each cylindrical arm: +where $\Delta\rho$ is the SLD contrast and $F_n$ is the form factor amplitude +of the $n$-th cylindrical arm: .. math:: - I_{\text{arm}}(q) = \int_0^{\pi/4} F_{\text{cylinder}}(q, \theta, L, R)^2 \sin(\theta) d\theta + F_n(q, \theta, \varphi) = \text{sinc}\!\left(\frac{q u_n L}{2}\right) + \cdot \frac{2 J_1(q \mu_n R)}{q \mu_n R} + +with $u_n = \hat{q} \cdot \hat{a}_n$ the projection of $\hat{q}$ onto the +arm axis, $\mu_n = \sqrt{1 - u_n^2}$, $L$ the arm length, and $R$ the arm +radius. Expanding the squared modulus into a double sum and exploiting the +reality of $F_n$ gives: -where $F_{\text{cylinder}}$ is the form factor of a cylinder. +.. math:: -The correlation term $I_{\text{corr}}(q)$ accounts for coherent scattering -between different arms. + I(q) = \frac{(\Delta \rho)^2}{4\pi} + \int \sum_{n=1}^{4}\sum_{m=1}^{4} + F_n F_m \cos\!\left(\frac{q(u_n-u_m)L}{2}\right) + \sin\theta \, d\theta \, d\varphi + +The cosine factor is the interference term between the centres of arms $n$ +and $m$, which are displaced by $\tfrac{L}{2}\hat{a}_n$ from the junction. Geometry -------- -The four arms are positioned along the following directions: -- Arm 1: (1, 1, 1) -- Arm 2: (-1, -1, 1) -- Arm 3: (-1, 1, -1) -- Arm 4: (1, -1, -1) +The four arms are oriented along tetrahedral directions. With +$A = 109.5°/2$ (the half-angle between arms), the arm unit vectors and the +corresponding projections $u_n$ are + +.. math:: + + u_n = s_n \cos A \cos\theta + \sin A \sin\theta \cos(\varphi - \varphi_n) + +where $(s_n, \varphi_n) = (+1,\ 0),\ (-1,\ \pi/2),\ (+1,\ \pi),\ (-1,\ 3\pi/2)$ +for $n = 1, 2, 3, 4$ respectively. Each arm has length $L$ and radius $R$. @@ -71,3 +88,8 @@ source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] have_Fq = False + +tests = [ + # Default parameters: length=400, radius=30, sld=4, sld_solvent=1 + [{}, 0.1, 1.0026479478514564e-3], +] From 03c63fbbc53501098c52d5d2353f6c3c754bffe2 Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 18 Mar 2026 12:24:34 +0900 Subject: [PATCH 05/17] Refactor Fq_n function in tetrapod.c and remove text generated through hallucination --- sasmodels/models/tetrapod.c | 12 ++++-------- sasmodels/models/tetrapod.py | 7 +------ 2 files changed, 5 insertions(+), 14 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index f2bf7c22..4ea80584 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -8,17 +8,13 @@ static double u_n(int n, double theta, double alpha) { sin(A) * sin(theta) * cos(alpha - phi[n]); } -static double A_n(double q, double u, double L, double R) { +static double Fq_n(double q, double u, double L, double R) { double quL2 = q * u * L * 0.5; - double term1 = (fabs(quL2) > 1e-12) ? sin(quL2) / quL2 : 1.0; - double mu = sqrt(fmax(0.0, 1.0 - u * u)); double qmuR = q * mu * R; - double term2 = (fabs(qmuR) > 1e-12) ? 2.0 * sas_J1(qmuR) / qmuR : 1.0; - - return term1 * term2; + return sas_sinx_x(quL2) * sas_2J1x_x(qmuR); } static double form_volume(double L, double R) { @@ -49,7 +45,7 @@ static double Iq(double q, double L, double R, double sld_particle, for (int m = 0; m < 4; m++) { double u = u_n(n, theta, alpha); double v = u_n(m, theta, alpha); - sum_arms += A_n(q, u, L, R) * A_n(q, v, L, R) * + sum_arms += Fq_n(q, u, L, R) * Fq_n(q, v, L, R) * cos(q * (u - v) * L / 2.0) * sin(theta); } } @@ -57,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle, } total += integral_alpha * w_theta; } - return contrast * contrast * total * 2 * M_PI; + return contrast * contrast * total * 2 / M_PI; } \ No newline at end of file diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index 54d0cc10..0dff0b0b 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -6,7 +6,7 @@ of four cylindrical arms radiating from a central point, oriented along the (1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions. -The scattering intensity is calculated as a powder average over all +The scattering intensity is calculated as an average over all orientations: .. math:: @@ -88,8 +88,3 @@ source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] have_Fq = False - -tests = [ - # Default parameters: length=400, radius=30, sld=4, sld_solvent=1 - [{}, 0.1, 1.0026479478514564e-3], -] From 7037bea3789b90ca531953fa71cb9cac18acbfd2 Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 18 Mar 2026 12:26:16 +0900 Subject: [PATCH 06/17] Remove unused numpy imports in tetrapod.py --- sasmodels/models/tetrapod.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index 0dff0b0b..fa9e3a80 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -67,8 +67,7 @@ * **Last Modified by:** """ -import numpy as np -from numpy import cos, inf, pi, sin +from numpy import inf name = "tetrapod" title = "Tetrapod with four cylindrical arms" From 6a2524746f6d4d6ce82bbd9eacc8f2d38d4fd91e Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sat, 28 Mar 2026 02:16:06 +0900 Subject: [PATCH 07/17] Update Iq function to include form_volume in return value and enhance author attribution in tetrapod.py --- sasmodels/models/tetrapod.c | 2 +- sasmodels/models/tetrapod.py | 3 ++- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index 4ea80584..4bcc14fb 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -53,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle, } total += integral_alpha * w_theta; } - return contrast * contrast * total * 2 / M_PI; + return contrast * contrast * total * 2 / M_PI * form_volume(L, R); } \ No newline at end of file diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index fa9e3a80..eae3521f 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -63,7 +63,7 @@ Authorship and Verification ---------------------------- -* **Author:** Yuhei Yamada +* **Author:** Yuhei Yamada (Github user name: Indigo Carmine, https://orcid.org/0009-0003-9780-4135) * **Last Modified by:** """ @@ -87,3 +87,4 @@ source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] have_Fq = False +opencl = True From f7de9f6e56d82ef82227a87dc522fe66ef64132e Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sat, 28 Mar 2026 02:35:24 +0900 Subject: [PATCH 08/17] Fix Iq function to remove form_volume calculation and adjust return value --- sasmodels/models/tetrapod.c | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index 4bcc14fb..aeeb8c07 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -53,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle, } total += integral_alpha * w_theta; } - return contrast * contrast * total * 2 / M_PI * form_volume(L, R); + return 1e-4 * contrast * contrast * total * 2 / M_PI; } \ No newline at end of file From 26b9812b552c099aaaab25091d7ac79d3c171d2c Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sat, 28 Mar 2026 02:48:58 +0900 Subject: [PATCH 09/17] Refactor angle calculation for clarity in tetrapod.c --- sasmodels/models/tetrapod.c | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index aeeb8c07..b0189fd9 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -1,5 +1,5 @@ -const double A = - 109.5 / 2.0 * M_PI / 180.0; // half of the angle between arms in radians +// half of the angle between arms in radians +const double A = acos(-1 / 3) / 2.0 * M_PI / 180.0; static double u_n(int n, double theta, double alpha) { const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2}; From 0e0450c6f2f148ce73cbfb539556ba37327ac3ca Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 17 Jun 2026 11:36:23 +0900 Subject: [PATCH 10/17] Add test values and fix OpenCL compatibility issues. --- sasmodels/models/tetrapod.c | 16 +++++++++------- sasmodels/models/tetrapod.py | 15 +++++++++++---- 2 files changed, 20 insertions(+), 11 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index b0189fd9..89b1e7bb 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -1,11 +1,11 @@ -// half of the angle between arms in radians -const double A = acos(-1 / 3) / 2.0 * M_PI / 180.0; +// acos(-1.0/3.0)/2.0 = half of tetrahedral angle ~54.7356 deg +#define TETRAHEDRAL_HALF_ANGLE acos(-1.0 / 3.0) / 2.0 static double u_n(int n, double theta, double alpha) { const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2}; const double sign[4] = {1.0, -1.0, 1.0, -1.0}; - return sign[n] * cos(A) * cos(theta) + - sin(A) * sin(theta) * cos(alpha - phi[n]); + return sign[n] * cos(TETRAHEDRAL_HALF_ANGLE) * cos(theta) + + sin(TETRAHEDRAL_HALF_ANGLE) * sin(theta) * cos(alpha - phi[n]); } static double Fq_n(double q, double u, double L, double R) { @@ -25,6 +25,7 @@ static double form_volume(double L, double R) { static double Iq(double q, double L, double R, double sld_particle, double sld_solvent) { double contrast = sld_particle - sld_solvent; + const double V_arm = M_PI * R * R * L; double total = 0.0; for (int dtheta = 0; dtheta < GAUSS_N; dtheta++) { @@ -45,13 +46,14 @@ static double Iq(double q, double L, double R, double sld_particle, for (int m = 0; m < 4; m++) { double u = u_n(n, theta, alpha); double v = u_n(m, theta, alpha); - sum_arms += Fq_n(q, u, L, R) * Fq_n(q, v, L, R) * - cos(q * (u - v) * L / 2.0) * sin(theta); + double arm = Fq_n(q, u, L, R) * Fq_n(q, v, L, R) * + cos(q * (u - v) * L / 2.0) * sin(theta); + sum_arms += arm; } } integral_alpha += sum_arms * w_alpha; } total += integral_alpha * w_theta; } - return 1e-4 * contrast * contrast * total * 2 / M_PI; + return 1e-4 * contrast * contrast * V_arm * V_arm * total / (4.0 * M_PI); } \ No newline at end of file diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index eae3521f..3a4f311a 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -79,12 +79,19 @@ # [ "name", "units", default, [lower, upper], "type", "description"], parameters = [ - ["length", "Ang", 400, [0, inf], "volume", "Cylindrical arm length"], - ["radius", "Ang", 30, [0, inf], "volume", "Cylindrical arm radius"], - ["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Tetrapod scattering length density"], - ["sld_solvent", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Solvent scattering length density"], + ["length", "Ang", 1000, [0, inf], "volume", "Cylindrical arm length"], + ["radius", "Ang", 50, [0, inf], "volume", "Cylindrical arm radius"], + ["sld", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Tetrapod scattering length density"], + ["sld_solvent", "1e-6/Ang^2", 0, [-inf, inf], "sld", "Solvent scattering length density"], ] source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"] have_Fq = False opencl = True + + +test = [ + [{}, 1.099643429303014e-03, 2.746377685546875e03], + [{}, 1.033727919136565e-01, 4.374285638332367e-01], + [{}, 3.145051218117226e-01, 3.527995198965073e-03], +] From 5faa83f95e1f58a4d1566870fdc2d41efaa64acc Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 17 Jun 2026 12:56:06 +0900 Subject: [PATCH 11/17] Fix tetrapod model descriptions --- sasmodels/models/tetrapod.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index 3a4f311a..ed4ebbe1 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -43,7 +43,7 @@ -------- The four arms are oriented along tetrahedral directions. With -$A = 109.5°/2$ (the half-angle between arms), the arm unit vectors and the +$A = 109.5 /2$ (the half-angle between arms), the arm unit vectors and the corresponding projections $u_n$ are .. math:: From 3ba6a039392b9ea87965671d0773b160aa060380 Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 17 Jun 2026 13:25:20 +0900 Subject: [PATCH 12/17] Update Tetrapod model to support CoreShellModel --- sasmodels/models/tetrapod.c | 36 ++++++++++++++------------ sasmodels/models/tetrapod.py | 50 ++++++++++++++++++++++-------------- 2 files changed, 51 insertions(+), 35 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index 89b1e7bb..f23d943d 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -8,24 +8,28 @@ static double u_n(int n, double theta, double alpha) { sin(TETRAHEDRAL_HALF_ANGLE) * sin(theta) * cos(alpha - phi[n]); } -static double Fq_n(double q, double u, double L, double R) { +// L: arm length, R: outer radius, t: shell thickness, R_c = R - t: core radius +static double Fq_n(double q, double u, double L, double R, double t, + double contrast_core, double contrast_shell) { + double R_c = R - t; double quL2 = q * u * L * 0.5; - double mu = sqrt(fmax(0.0, 1.0 - u * u)); - double qmuR = q * mu * R; - - return sas_sinx_x(quL2) * sas_2J1x_x(qmuR); + double V_c = M_PI * R_c * R_c * L; + double V = M_PI * R * R * L; + return sas_sinx_x(quL2) * + (contrast_core * V_c * sas_2J1x_x(q * mu * R_c) + + contrast_shell * V * sas_2J1x_x(q * mu * R)); } -static double form_volume(double L, double R) { - // V = 4 * \pi * R^2 * L +static double form_volume(double L, double R, double t) { + // V = 4 * pi * R^2 * L (outer volume of 4 arms) return 4.0 * M_PI * R * R * L; } -static double Iq(double q, double L, double R, double sld_particle, - double sld_solvent) { - double contrast = sld_particle - sld_solvent; - const double V_arm = M_PI * R * R * L; +static double Iq(double q, double L, double R, double t, + double sld_core, double sld_shell, double sld_solvent) { + double contrast_core = sld_core - sld_shell; + double contrast_shell = sld_shell - sld_solvent; double total = 0.0; for (int dtheta = 0; dtheta < GAUSS_N; dtheta++) { @@ -46,14 +50,14 @@ static double Iq(double q, double L, double R, double sld_particle, for (int m = 0; m < 4; m++) { double u = u_n(n, theta, alpha); double v = u_n(m, theta, alpha); - double arm = Fq_n(q, u, L, R) * Fq_n(q, v, L, R) * - cos(q * (u - v) * L / 2.0) * sin(theta); - sum_arms += arm; + double Fn = Fq_n(q, u, L, R, t, contrast_core, contrast_shell); + double Fm = Fq_n(q, v, L, R, t, contrast_core, contrast_shell); + sum_arms += Fn * Fm * cos(q * (u - v) * L / 2.0) * sin(theta); } } integral_alpha += sum_arms * w_alpha; } total += integral_alpha * w_theta; } - return 1e-4 * contrast * contrast * V_arm * V_arm * total / (4.0 * M_PI); -} \ No newline at end of file + return 1e-4 * total / (4.0 * M_PI); +} diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index ed4ebbe1..f0e1a4e8 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -16,22 +16,26 @@ \left|\sum_{n=1}^{4} F_n(q, \theta, \varphi)\right|^2 \sin\theta \, d\theta \, d\varphi -where $\Delta\rho$ is the SLD contrast and $F_n$ is the form factor amplitude -of the $n$-th cylindrical arm: +where $F_n$ is the core-shell form factor amplitude of the $n$-th arm: .. math:: F_n(q, \theta, \varphi) = \text{sinc}\!\left(\frac{q u_n L}{2}\right) - \cdot \frac{2 J_1(q \mu_n R)}{q \mu_n R} - -with $u_n = \hat{q} \cdot \hat{a}_n$ the projection of $\hat{q}$ onto the -arm axis, $\mu_n = \sqrt{1 - u_n^2}$, $L$ the arm length, and $R$ the arm -radius. Expanding the squared modulus into a double sum and exploiting the -reality of $F_n$ gives: + \left[ + (\rho_\text{core} - \rho_\text{shell})\, V_\text{core}\, + \frac{2 J_1(q \mu_n R_\text{core})}{q \mu_n R_\text{core}} + + (\rho_\text{shell} - \rho_\text{solvent})\, V_\text{shell}\, + \frac{2 J_1(q \mu_n R)}{q \mu_n R} + \right] + +with $u_n = \hat{q} \cdot \hat{a}_n$, $\mu_n = \sqrt{1 - u_n^2}$, $L$ the arm +length, $R$ the arm outer radius, $R_\text{core} = R - t$ the core radius +($t$ = shell thickness), and $V = \pi R^2 L$, $V_\text{core} = \pi R_\text{core}^2 L$. +Expanding the squared modulus into a double sum gives: .. math:: - I(q) = \frac{(\Delta \rho)^2}{4\pi} + I(q) = \frac{1}{4\pi} \int \sum_{n=1}^{4}\sum_{m=1}^{4} F_n F_m \cos\!\left(\frac{q(u_n-u_m)L}{2}\right) \sin\theta \, d\theta \, d\varphi @@ -53,7 +57,7 @@ where $(s_n, \varphi_n) = (+1,\ 0),\ (-1,\ \pi/2),\ (+1,\ \pi),\ (-1,\ 3\pi/2)$ for $n = 1, 2, 3, 4$ respectively. -Each arm has length $L$ and radius $R$. +Each arm has length $L$, outer radius $R$, and shell thickness $t$. References ---------- @@ -70,18 +74,21 @@ from numpy import inf name = "tetrapod" -title = "Tetrapod with four cylindrical arms" +title = "Core-shell tetrapod with four cylindrical arms" description = """ - Calculates the scattering from a tetrapod structure with four cylindrical - arms radiating from a central point. + Calculates the scattering from a core-shell tetrapod structure with four + cylindrical arms radiating from a central point. Each arm has a cylindrical + core and a coaxial shell. """ category = "shape:cylinder" # [ "name", "units", default, [lower, upper], "type", "description"], parameters = [ - ["length", "Ang", 1000, [0, inf], "volume", "Cylindrical arm length"], - ["radius", "Ang", 50, [0, inf], "volume", "Cylindrical arm radius"], - ["sld", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Tetrapod scattering length density"], + ["length", "Ang", 1000, [0, inf], "volume", "Arm length"], + ["radius", "Ang", 50, [0, inf], "volume", "Arm outer radius"], + ["thickness", "Ang", 10, [0, inf], "volume", "Arm shell thickness"], + ["sld_core", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Arm core scattering length density"], + ["sld_shell", "1e-6/Ang^2", 0.5, [-inf, inf], "sld", "Arm shell scattering length density"], ["sld_solvent", "1e-6/Ang^2", 0, [-inf, inf], "sld", "Solvent scattering length density"], ] @@ -91,7 +98,12 @@ test = [ - [{}, 1.099643429303014e-03, 2.746377685546875e03], - [{}, 1.033727919136565e-01, 4.374285638332367e-01], - [{}, 3.145051218117226e-01, 3.527995198965073e-03], + # thickness=0 reduces to uniform cylinder: contrast = sld_core - sld_solvent = 1 + [{"thickness": 0}, 1.099643429303014e-03, 2.746377685546875e+03], + [{"thickness": 0}, 1.033727919136565e-01, 4.374285638332367e-01], + [{"thickness": 0}, 3.145051218117226e-01, 3.527995198965073e-03], + # default core-shell parameters + [{}, 1.099643429303014e-03, 1.846798618097820e+03], + [{}, 1.033727919136565e-01, 1.811197691988146e-01], + [{}, 3.145051218117226e-01, 1.080396557066592e-03], ] From 08584bf66bc755ee06f9dded46573c36362c8feb Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 17 Jun 2026 13:45:59 +0900 Subject: [PATCH 13/17] add radius_effective_modes and the function --- sasmodels/models/tetrapod.c | 38 ++++++++++++++++++++++++------------ sasmodels/models/tetrapod.py | 10 ++++++++-- 2 files changed, 33 insertions(+), 15 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index f23d943d..58301ed6 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -16,18 +16,27 @@ static double Fq_n(double q, double u, double L, double R, double t, double mu = sqrt(fmax(0.0, 1.0 - u * u)); double V_c = M_PI * R_c * R_c * L; double V = M_PI * R * R * L; - return sas_sinx_x(quL2) * - (contrast_core * V_c * sas_2J1x_x(q * mu * R_c) + - contrast_shell * V * sas_2J1x_x(q * mu * R)); + return sas_sinx_x(quL2) * (contrast_core * V_c * sas_2J1x_x(q * mu * R_c) + + contrast_shell * V * sas_2J1x_x(q * mu * R)); } static double form_volume(double L, double R, double t) { - // V = 4 * pi * R^2 * L (outer volume of 4 arms) - return 4.0 * M_PI * R * R * L; + // V = 4 * pi * (R + t)^2 * L (outer volume of 4 arms) + return 4.0 * M_PI * (R + t) * (R + t) * L; } -static double Iq(double q, double L, double R, double t, - double sld_core, double sld_shell, double sld_solvent) { +static double radius_effective(int mode, double L, double R, double t) { + switch (mode) { + default: + case 1: // equivalent volume sphere + return cbrt(form_volume(L, R, t) / M_4PI_3); + case 2: // length of tetrapod arms (L) + return L; + } +} + +static double Iq(double q, double L, double R, double t, double sld_core, + double sld_shell, double sld_solvent) { double contrast_core = sld_core - sld_shell; double contrast_shell = sld_shell - sld_solvent; double total = 0.0; @@ -45,16 +54,19 @@ static double Iq(double q, double L, double R, double t, double w_alpha = GAUSS_W[dalpha] * M_PI; // adjust weight for the new range + double u[4], F[4]; + for (int n = 0; n < 4; n++) { + u[n] = u_n(n, theta, alpha); + F[n] = Fq_n(q, u[n], L, R, t, contrast_core, contrast_shell); + } double sum_arms = 0.0; for (int n = 0; n < 4; n++) { - for (int m = 0; m < 4; m++) { - double u = u_n(n, theta, alpha); - double v = u_n(m, theta, alpha); - double Fn = Fq_n(q, u, L, R, t, contrast_core, contrast_shell); - double Fm = Fq_n(q, v, L, R, t, contrast_core, contrast_shell); - sum_arms += Fn * Fm * cos(q * (u - v) * L / 2.0) * sin(theta); + sum_arms += F[n] * F[n]; + for (int m = n + 1; m < 4; m++) { + sum_arms += 2.0 * F[n] * F[m] * cos(q * (u[n] - u[m]) * L / 2.0); } } + sum_arms *= sin(theta); integral_alpha += sum_arms * w_alpha; } total += integral_alpha * w_theta; diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index f0e1a4e8..2bb158f9 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -97,13 +97,19 @@ opencl = True +radius_effective_modes = [ + "equivalent volume sphere", + "length of tetrapod arms (L)", +] + + test = [ # thickness=0 reduces to uniform cylinder: contrast = sld_core - sld_solvent = 1 - [{"thickness": 0}, 1.099643429303014e-03, 2.746377685546875e+03], + [{"thickness": 0}, 1.099643429303014e-03, 2.746377685546875e03], [{"thickness": 0}, 1.033727919136565e-01, 4.374285638332367e-01], [{"thickness": 0}, 3.145051218117226e-01, 3.527995198965073e-03], # default core-shell parameters - [{}, 1.099643429303014e-03, 1.846798618097820e+03], + [{}, 1.099643429303014e-03, 1.846798618097820e03], [{}, 1.033727919136565e-01, 1.811197691988146e-01], [{}, 3.145051218117226e-01, 1.080396557066592e-03], ] From a1461ec05ece31f59c8eb01a8c43e96a6c6f191b Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Wed, 24 Jun 2026 15:47:51 +0900 Subject: [PATCH 14/17] add DOI, adn rename core_radius value from radius. --- sasmodels/models/tetrapod.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index 2bb158f9..2aca63bc 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -62,7 +62,7 @@ References ---------- -#. Seoki Kyoo Seo *Korean J. Chem. Eng.* 34(2017) 1192-1198 +#. Seoki Kyoo Seo *Korean J. Chem. Eng.* 34(2017) 1192-1198 DOI:10.1007/s11814-016-0341-x Authorship and Verification ---------------------------- @@ -85,7 +85,7 @@ # [ "name", "units", default, [lower, upper], "type", "description"], parameters = [ ["length", "Ang", 1000, [0, inf], "volume", "Arm length"], - ["radius", "Ang", 50, [0, inf], "volume", "Arm outer radius"], + ["core_radius", "Ang", 50, [0, inf], "volume", "Arm core radius"], ["thickness", "Ang", 10, [0, inf], "volume", "Arm shell thickness"], ["sld_core", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Arm core scattering length density"], ["sld_shell", "1e-6/Ang^2", 0.5, [-inf, inf], "sld", "Arm shell scattering length density"], From 8ccde0b64a5c294fcbd6e88daf21a24b49ff7012 Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sun, 5 Jul 2026 21:45:56 +0900 Subject: [PATCH 15/17] add note about overwlap. --- sasmodels/models/tetrapod.py | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index 2aca63bc..a3796653 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -59,6 +59,21 @@ Each arm has length $L$, outer radius $R$, and shell thickness $t$. +.. note:: + + Each of the four arms is modelled as a complete cylinder of length $L$ + extending from the central junction, so the arms overlap near the origin. + This overlap is neglected in two ways, following the treatment used in the + reference. First, the particle volume is overestimated (the effect being + largest when the arms are short and wide), which affects the calculated + intensity through the volume normalisation. Second, the scattering amplitude + is the plain sum of the four cylinder amplitudes, so the overlapping region + near the origin is counted more than once; because this region is small + compared with the arms, the resulting artefacts are expected to appear + mainly in the high-$q$ range. Consequently the model is only valid for + long-arm tetrapods, i.e. for arm lengths much larger than the arm width, + $L \gg R = R_\text{core} + t$. + References ---------- From fa04b9a3efd220e8ca9dd8e492b270f150f208dd Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sun, 5 Jul 2026 22:26:40 +0900 Subject: [PATCH 16/17] fix value name "radius" and define precaluculated values. --- sasmodels/models/tetrapod.c | 22 ++++++++++++---------- sasmodels/models/tetrapod.py | 15 ++++++++------- 2 files changed, 20 insertions(+), 17 deletions(-) diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c index 58301ed6..36494528 100644 --- a/sasmodels/models/tetrapod.c +++ b/sasmodels/models/tetrapod.c @@ -1,23 +1,25 @@ -// acos(-1.0/3.0)/2.0 = half of tetrahedral angle ~54.7356 deg -#define TETRAHEDRAL_HALF_ANGLE acos(-1.0 / 3.0) / 2.0 +// Half of the tetrahedral angle acos(-1/3)/2 ~54.7356 deg. +// Precompute its cosine and sine: cos = sqrt(1/3), sin = sqrt(2/3). +#define COS_HALF_ANGLE 0.57735026918962573 // sqrt(1/3) +#define SIN_HALF_ANGLE 0.81649658092772603 // sqrt(2/3) static double u_n(int n, double theta, double alpha) { const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2}; const double sign[4] = {1.0, -1.0, 1.0, -1.0}; - return sign[n] * cos(TETRAHEDRAL_HALF_ANGLE) * cos(theta) + - sin(TETRAHEDRAL_HALF_ANGLE) * sin(theta) * cos(alpha - phi[n]); + return sign[n] * COS_HALF_ANGLE * cos(theta) + + SIN_HALF_ANGLE * sin(theta) * cos(alpha - phi[n]); } -// L: arm length, R: outer radius, t: shell thickness, R_c = R - t: core radius +// L: arm length, R: core radius, t: shell thickness, R_o = R + t: outer radius static double Fq_n(double q, double u, double L, double R, double t, double contrast_core, double contrast_shell) { - double R_c = R - t; + double R_o = R + t; double quL2 = q * u * L * 0.5; double mu = sqrt(fmax(0.0, 1.0 - u * u)); - double V_c = M_PI * R_c * R_c * L; - double V = M_PI * R * R * L; - return sas_sinx_x(quL2) * (contrast_core * V_c * sas_2J1x_x(q * mu * R_c) + - contrast_shell * V * sas_2J1x_x(q * mu * R)); + double V_c = M_PI * R * R * L; + double V_o = M_PI * R_o * R_o * L; + return sas_sinx_x(quL2) * (contrast_core * V_c * sas_2J1x_x(q * mu * R) + + contrast_shell * V_o * sas_2J1x_x(q * mu * R_o)); } static double form_volume(double L, double R, double t) { diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py index a3796653..857fdaa4 100644 --- a/sasmodels/models/tetrapod.py +++ b/sasmodels/models/tetrapod.py @@ -23,14 +23,14 @@ F_n(q, \theta, \varphi) = \text{sinc}\!\left(\frac{q u_n L}{2}\right) \left[ (\rho_\text{core} - \rho_\text{shell})\, V_\text{core}\, - \frac{2 J_1(q \mu_n R_\text{core})}{q \mu_n R_\text{core}} - + (\rho_\text{shell} - \rho_\text{solvent})\, V_\text{shell}\, \frac{2 J_1(q \mu_n R)}{q \mu_n R} + + (\rho_\text{shell} - \rho_\text{solvent})\, V_\text{outer}\, + \frac{2 J_1(q \mu_n (R+t))}{q \mu_n (R+t)} \right] with $u_n = \hat{q} \cdot \hat{a}_n$, $\mu_n = \sqrt{1 - u_n^2}$, $L$ the arm -length, $R$ the arm outer radius, $R_\text{core} = R - t$ the core radius -($t$ = shell thickness), and $V = \pi R^2 L$, $V_\text{core} = \pi R_\text{core}^2 L$. +length, $R$ the arm core radius, $R + t$ the outer radius ($t$ = shell +thickness), and $V_\text{core} = \pi R^2 L$, $V_\text{outer} = \pi (R+t)^2 L$. Expanding the squared modulus into a double sum gives: .. math:: @@ -57,7 +57,8 @@ where $(s_n, \varphi_n) = (+1,\ 0),\ (-1,\ \pi/2),\ (+1,\ \pi),\ (-1,\ 3\pi/2)$ for $n = 1, 2, 3, 4$ respectively. -Each arm has length $L$, outer radius $R$, and shell thickness $t$. +Each arm has length $L$, core radius $R$, shell thickness $t$, and hence outer +radius $R + t$. .. note:: @@ -72,7 +73,7 @@ compared with the arms, the resulting artefacts are expected to appear mainly in the high-$q$ range. Consequently the model is only valid for long-arm tetrapods, i.e. for arm lengths much larger than the arm width, - $L \gg R = R_\text{core} + t$. + $L \gg R + t$. References ---------- @@ -100,7 +101,7 @@ # [ "name", "units", default, [lower, upper], "type", "description"], parameters = [ ["length", "Ang", 1000, [0, inf], "volume", "Arm length"], - ["core_radius", "Ang", 50, [0, inf], "volume", "Arm core radius"], + ["radius", "Ang", 50, [0, inf], "volume", "Arm core radius"], ["thickness", "Ang", 10, [0, inf], "volume", "Arm shell thickness"], ["sld_core", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Arm core scattering length density"], ["sld_shell", "1e-6/Ang^2", 0.5, [-inf, inf], "sld", "Arm shell scattering length density"], From e1f36e598c1a1260cae1a774fe1044464b19234d Mon Sep 17 00:00:00 2001 From: IndigoCarmine Date: Sun, 5 Jul 2026 23:26:23 +0900 Subject: [PATCH 17/17] add a tetrapod figure --- sasmodels/models/img/tetrapod.png | Bin 0 -> 29342 bytes sasmodels/models/tetrapod.py | 7 +++++++ 2 files changed, 7 insertions(+) create mode 100644 sasmodels/models/img/tetrapod.png diff --git a/sasmodels/models/img/tetrapod.png b/sasmodels/models/img/tetrapod.png new file mode 100644 index 0000000000000000000000000000000000000000..06f0d90b4c7c3b3932feb794da4c2c3e2ccf49a8 GIT binary patch literal 29342 zcmYIv1yqz>v^LVA!jJ;e!hi@!$Iuee(hbre-3^k`NJxitNOzZXcMT2FF!TTe%s=|w z|K7EjHS5fJIqyE_?0oif!apfWzrZHLMnXb*A@fl}1qtby2jWA9iH`V%JVzdi_hTa3%lMTKem7EyBIMhC?Y! zaBr6%c$uK&w-pWyp4TD3L`4oz|7KJ1_q*$m8qe`cre4|rAt^93J${32JiQ>*KhpNG<&R{%b8hYi#$<)+>8Qqc56@huXO|Q&Wu$pNCD^S|n8W8>#Nd@;YpU+x&a1Q8_0iX)w@ja`XU%5Mpy(njU zpox@*NHZY?v{lmHn+H%@en;lAq?xotMSfxV0&`O588QV5m~I#-L5(IZx_s!kHZc=)T&Ax z7nvSa=WCB1QiPVnZ8SNav*ROgUFb=y7rRpH;p477Dl$+8?^bjJwF$lqbzbLPxAk$Q zL5h|jFmw)H!eyBx6@4crn>S{T<4q&_lTb}Z)#?=rPHgx=rxmAw!sN2euOEfn-p`c$ zR`>9xG2sCy7?zFT(b3>y;j;~_e-e6E!am$LKJi%LYsU?QwN3WQeV9H%f72K(bMtnK zKI9cPd0$X`S}7gVJaz<*aeu76^6{}s{Lgw(>yY_R8W9}rVbX@zXu)Af9($nmAv6-b zu8D`y2`+agyj$M?-THjGZbs^uNNpB~H!XFwn4(YT8qQcBy(T92f+H6n zppzb|xjn|5sG_Fc%}2%g%Bp`TulblL`vc}sg?_8YY#wD_=nL|u?HG_{y-svHHfjea z4KCu$X!6hXQ89DkR1$izw9}@Bm+N$dbWZqd#s%Wn%(uAr0F|raHRd=6$DaACFX}+n z`*nN!brzNu6g6*XVhB|jgWo$Q)y+Tes(f6oWT%X6hIGJgBRg+#>wt_& zQw-OpT9dy)*U~{Tz_#ZPWaMxrC3QBycBo8kYqm|Dqc%nj1VP#oepPmm(V(ixu3ct9 zKN0p_9U%w!2T`8|bYX2~zir4-U%c*-g!a5M*{gQ|+lF}JI>)__wRV?Ii404`!fd+T zisbjXeGhK1-Bz&oBH&&3;gpH@fZ<_2fKdZ*Q58{(dRSe ztC?LpTc1(=E}9szw>)U?1G<7}{;ri4liVs&VE*Qq$YNo3cook_sTYmmOJi2YSW@=K zmK5U4ns-6m1RjI`4A8=V(`bico09VlTpS)cXi>pqLgX@F3drD82 zH`U2z+bXPMy9xcXrKOWP82|fRu#mbBz&&X+k$YS6RCS_xKx5@dp>T%n4+p>p;amjO z{>>H)I1q4xmXSf5vP34;?Cp;*zyrHppV!AOH-aPy__Y7qqJB47{FIk6J&=u^2iHbd zx>6{R-P3RDw~_&KUfolAcP>`OJFyLv00L220#EYP8!~y{shb6U)M!S7S>kBV4$^L= z-nMrmc-|X&z@4O8rzw3cujf@VZqWo&BS|Hrv%uML4%AadhNHP!|T}+y54^ z7$ajRtj<&uvzGa@^?ABksZ1JgtwxNm5=?*wRa{X(0C5`svtVPV(v){(B0A9H6rIV; zuijZck9v<%8aKeBkc@m8^#6S>=Vk;rLSW<;cq)CKu%J+xz2>&{yLBu!cH;w`^+pc< zKhOF_tP;rU%x$)!nzg@h)DFmA@q$M9#|)og zeR_hr$M>diat~7*MA|?^>)ffANyPS4c8(Le!`=6 zV-$WD%STDl(Sj;UQ}mv@n6cipa11D^2M*`F`}rz?Fxvi378{ZQV+kt>Su#T^@)z5(WptoZXU}{wv*pt0{$INCHB?Bf^s9A;}Be# z({kkgNv<*Xx>uLGdFXT}$~Mtlaibt$!(GG$@2k7-22ueBn7nd+>}04M?TZ6jpMK|u z^D*Y)(_8hc>9Ly_Y6lHfFk?p(c{gvYkTS`Ta04Bw_$$9Gcw;>=i;NauPxBW;?4GcrgA0m1au3hgb*PNa_987z~CM3 z_0}wrhowxUhdFz_;mhnAbed7gP2id@oBiz%{cC+u=XmklAB6(<+!4BdX z#G=zd6nwu^E>PXNVg(}H0WBiU`zsTLfGmY=+M9!BfG0b*so`H|Gtx|dHiy@KB(x9m zLsNhE)0%hgAwi=J7_cmH{AjGIxDRmHykl{*;fXFJv9+ysRSm_959!#=%NCJU?wLZf zwVxT2>{X-VYCI7qtCp6VJrf>o_5Zfj--s}pbZ;uYAk2ed>_yf7v2=jf>g1QcRiS># zTR%D*)nn}jJqLsYxIdhwH@>w_@^+%sLn^0q6tLV^)A0fl890*N(rh@}ZE~F&kL_wQ zkuW>pB=#0LhUa6t_kvs83(8hTwwd9tQU(
G-rVtj)xS{CwR{F={ckS$&kbT9NA zODsw5`w>E7!gY35dF((8J$6CQS6X@Bkr`OHMi$TB7#{z-qhx$A&bRQ56E$h?3U;>(-r0TX zKMIha%a1^s?5DD>h-!X^wm~p~qv2Tyfe-Vcs&YGczCC$W5hIep;fwYvvf=fVdndW? zg1z`$Hu0{`Ct$7sz!uK&9&m03gbB}*Jg1wV@OXGF^F|v);Cwf;Ktt!muGNz0jFW`-hDWkLz%;KW~0t`~xBW#e>38Q(yDfP44XRP_tFc6p4zoyCO8)8ehj zHr_c3%_z%-II>0Cs8670Y;voM7|d5UrZ=_}bcYF%5!sHaFWgku&Uj{@V!1_ZH(aRm z+qHPtKrkOe%pbjJ^Dt0i0g|^X!q)8y0c4ObDyJm`Ub;0^?v^Pg5QaY z>CK>_F=3^fUopEO5}GEZAb%X+E5r0`finkHlB=qzRsyfh`MZI*){oL($D-W_G}}+! zcvY^G7#W&40NQGiPy#SU4he33{n4a5O_g4`(>5neNUgHgI2{!>R9BYdRY-Ot$Cjc}FLFl^(mjx8^}3a@}&m!@*fu+3p<$gh=0NMGioBSJ6;t zrb_eMw&7y2H(Vl~8=Dd-s(l~Jsk&QRL6F5_+gD0wkD|H4g)3-1?bZBtum8oyUq^%5 zx9BTS5it+|8r&k>Ehy$r@>ORrT^yn)4zIuAzPTpJTHlm5sDAO)AKHSW)yVLx(SL=$ zw(Y2h^^XYrOtEHkB7lcDHhh=P^rTf6zNW*&Q0l(LUhRqqVTTQd!Cw|1(2;Qo#+^6^ zg3=#M`n>fb{V%ng61>+svRiF1J@x*Js#RVw84pM+kF2WYl^xgDbR|-=etHsSY}ijS zVFdpC%_nXX*3Wc4XkTReTQaB5Gdx@!HHa6CI$ww6I4Y;$hc-2QkNU|g=*i4bA&6YD zsV-66cl4=sjo@%GL#j?Ls5172chuZKxO>?hPc+Exp(~fz69(1|2JZg14Ieh<;TojjMX8J-Xq$$0+DZ-`OrXsn4AbdiEnRd8dPr|GkUM z#T>cT!+iY}r1{5Yw?1Lyx5+!a#2ZBih`1`G=$E^RVWwh?LJx6!x77?nm)n9;g4 z11WIlFK8(0@C_MV8y(){PI;E^QRCj zsnWQEWaXFCkh}qu`$xLh22o4*@o*8x2K29u9oOX#+agbN&^2#+XLm$#Wc)|^Y@h=f zL-)CwX#5%a?vpeV3eR!p{cGH2>#U09iO-e_e11!9Zd`Hiq12UBU zaa@OoCw_MBD?%A1_P(XB%Wn`GlXZ>#1gnZkJc66oXDFqkZGH_koZi+JWkqilVXtXN zlNQG2P1?DKPHf-$f?KP`!qvHx zq4A8h)?A={dUalj71^r}!xh^4px4p$JCMN1y>Tyi84cs9dPm{Ok*D77!mo$&*RHrR z0wLhj)@bO$g@p7S9b+}>p_GtlHz#~)UU4TGvJNodlf6qeq`N=yj)O7N&Bp

5k1$ z3p>fWGaJ-Kb3@n3J(#X;T2wHc=}isUKKh*q?{>P--Ss@1VH~NeQWB1Ati;ThD=HIN zkl?zC^B2Q))}e00&rHF+S({irfO9&uvY@AB_tGA zR_g~=USK}U4V*7RFNUskPgj5@J)gR=`~1)o1P)Y~G6$6|w&!vt#2w4-6wnmW#21}9 z49+(E@=CVtQVk}@wkwIRhx&R_tp8^DIvN#S6aHuszYS*Xl-{^+aiO&fH?Lq>*^<6_ zEf9DM;q5UplUz}_X#l^1m{amQCgZj49fuz_RI(vY>j>hLlIA2EJ84oEP&ZeLKAk^4 ze|Z{EnnoOYU)wxk0g_&t;Is0mXG-d$y`Nm?=`87D>T2F8%$aM5P}C25{J45;eod~u zx`WkmwuW+i&IfozoL+)i2zFcOZ-0yARjjS>oXD!QO9AnQI1qm5hg-viu2UbKpVsV+ z;t((=01y0UQ3`6OX{>a+D%oZ5Uc&S^{XWg4lI}nO?eKRi=>*l8)2=Vg8S1*QrgXek zPS>Ku7x)|l1)IYcfKv+=7?*kl!bS9nOQPKX(r_I@k^XeZx)%CR1~cr3Gs^vSOZQHu zL&Yv3F=)tYb_Bk#H;a#%AMKVODwXzmO!;@vxgssZ(R5c!w}{wp2=r33&BkzOjKyBd%izJx+~wrE@f4n*4}N|30v6C%(_KK?$AGX&x;L- z`>Dvy>X(gK>)+PAJwoBV#=grPtb{>kaJIRo-QUrHzzj0Dm{m7oA5e$=iX~`WQ~b<^u5o74-5EugO~L4){P*C$X2|VhHf6A=i}U__VRE(p`@Bob`>=ztY%L%M-b>-(gfEGNXp6&9c;BJ*?u8 zQBGXB3uynXjkMaE7SaofiLL^c5Zw??6HbrHS4Y=uy^2SjMx}K&E=2Bdbc{S=wR0XB zC+@7LF)G_kpwpLz*u{TPxKG9)r2v@<888!`*75~Zv%SyEjHx@2KKlQf` zNHC{uyEZYrlB82m=^xPi;iBrh>yTvaTXILMz%wuweWe!@6JkvI(zSW9y+Qk}eYmqQHIx0T&f}%Zp*Oklfpn^+}4N zHS!3TFZu}RDB~UTr;XUf!K4F%ej#FG0p}cIiG`-O8~3n{&vscAoMw6-N-r0Mm<#7Z zMzOtkaLTngC(LQmgX%1ve5QFO>b-Wk#qqJa?G0)I*o2J;CNdQ6j_%^!saSza{H_9p z2e0TuM3bmHJidRan{5*1X0_DNlJ_(++zB{;_Q^=}yXrU4cwKbcqsavB8y~s;*saye z5xuHs<=3;HsdSDe6Sy6w|L&RoUSZ)r3&Grd9FYN`|DG~=vo}R@J`r>g3ICtydH(n$3gmJI~9D^G2Bo(BO{>GYj%t{+0 zc?8wYQ*Tdq906qn?vR^Bq4}<<9)fwPdHJ-|kAerUjFfCqtj#&ib=+!h#xnCVj)>8b zR)k+Ue(+Zlfj&VzX|w%KRX^&?_si%$?4TXpeH*afO>UUra1b*0jZ0=AkRt5btfoNf z9dg_5QTO-byfb@_^eXf)zUbu?o_piC|2~DkLJrjbcmGk{qXqi&)H6}ahpjEzR<{#l zN5|?rHR|_$!}uQM@~9$2&zj93y**hrW7p!M1)GaqVWnpr7P-BcK8d0$9IEk4Y9G_z zY8QzrKHORFCg0{vJTfD}|5z`a4+&u$+GaJM2pZIxFmu$BMAb#hs~FTdh9=wXzfWLD zXyNIRI8}TTj92AuIgGW?L&HGQpmX-Pcji@(oF{Yx!pz|6`6ao1)aIXY*BMT5U9vu4 zt6Q`+n=;iB1Kmg_xW)C||8Tt0oPg`IfihNv9KfWJ(naz!X|%@Q7>jwnV{xJ9G+S>v z5vWy1A6js^u3+|&pw9wtX+P6kkON8lpgQGM@2U`bZwF0FX&cxVUdxrVx_6yTzS=53OX^Xa0w6dkR>Td4DWr8q zd&fDz4-grJLbTL;X=xJCUrbmV5f4VTmC{O@*<996szIBCpE^UB-4mHTtD_qozd&B9TxKZ?0z4NOsH?(C2lF)SHQ{A5Mc zTmP_;Uu?^k@sW-0Wzd#n8VLGa&9-{J=s|H;$;s^bixVUwzdq-38 z0#>&EH=;o=PH5YE7$k~#&h9FcYs6eOe;Psxq;lH?v&~Glf$%mWyT)sdLI|+!05*_U6qT%;ss&HM_id zMZmu1H6Iffa=B2t6d3-c+fZ>ue?rqXIhF(}gIXC&;8CE4fUISzFHnWSQ+Uov+Bc|E zlKZFpfbJXD-G{V@sP=_B`ylXmyM;-qA%yj1)}O7;N(ZoqEwa;L?4s=vGjUm?JLBfz{%t~V0tUa26gM^{Brmeb}E=M{g2?Y<7X84tG2lVoqj z#n>jl`TMJ<8;wD#jl6kG>ZMsda*txnfa#w3yz`LdY>={M13WOZ)eI`oSPm7I4f6XA zc<^3_8Dw2%?9I-oYwFLfym5%BMXY!UY5L~1FLf}?o}Y-&VM=3>ofj}(qI<+N*U)pK zsOHYWn_IgZ{;QtUYoc9)g$`TssvZh$7pa_3irXu{6zH7QM6zGI^-;w|6Tl%JEBOuU zK3pHA#8Dj(y<9wtp)kB{hfeFHely@6%6_bSJbsaM@qJ-VZ1BO)n54IV4ex0MPu;?0TXq5fnm?QncI~V(-agEg) z;#m%yu!6SvcUSo@fDSx&J;T;A>Rmcep87_Kj0ko?fRdpJ+zjk83`it&7P2RKblSLZ z!&}6W5V`!=e7&yzJJ%7pQ;&Fz+3-%bVV zZ)&~{E0%Iy)_mpQ{qSrayaauV0}=leSX1fj50(rZmV(6U2l3e_M6U^Lij23q<6fO- ztGiJ*m_?x!AsP8oig-a;r|XGS(}pONP90fpvm{OY-3#U%0Y2JkFTWtVTEO4Kgtu$k zMKjN=A1v3hM6X-=Yml;FTKf_jdg>qQGAChyCrQqX$BD7Pk&xBGK^iuOSUGRNS& zhVAwV7>BjydPj1blRxb!iAfB+_GI6vj;i=7K-4?m$_BNO2)D4ZfS`%UJt%~{T3g4Q^Cq9*zyeVH2} z_r6>6APux{O&C(iY-0+ceQ+|BT$4KcLDI+kLTL#}scF_FUO;h>{8nipdoJ|m$gsuO zf322ayPMq60zZ6V>lCZ)EzpcRk-$r7qbPXI>r6Jo0Z}14A(@vU zxlN|Gdb8`n6Xb9b4|yIwN{F`PJ^+$Br)YJ z4>6)STJrf1cALD>Re>e1P+F z2ijee#&Pd)EFWR6tt{NKNl$ zf5d)e&O*21X#Snj0&{>tAkLSbFbWsefN{F^>j2dNym*y9p1e88$VOSRIoEndnueko z+bf1Z&4?1QMtxMOene65oLP4xkFw0SL&(!-vzCv_Pxh655j}`MK)QR^^(@-}jm@*J z`bz|R;s?6PyZB7q^YL#Q$G%35#=@t+7wtqSJ)*1UC)ngn#$4z;1~@k`vT@7jO+_3U z;n9z7*%PbvT*9o#Gfgcvp~@ROc)KK@`Zd@);|T6TH(wzBw)eY4+v zyMqv(Poo3B%oosWi>M{VNEE2nA_^hDTL*{l^2<>TF6BK2aH5d*nxi*txkL;u3MJ_t z`tWz$_zlw}9=CJ$Ho;M8x$iTxS<^xcJ<$Y@!oN#2oR6)V-NN>MHR)GqdS;W_`ama} z?yGXtP3b%)g4_pUz}> zZ{0hj?(w&q9{fUkdjQD%iwephL_x7>!j7-M31HjL8k-`7GgUWE)v`S$R#jj?nvkE4 ztT(xp0$(hhlXAKKF4b*tZqzrrJZyw+uqXvVJFwQSj8jh1aAY!h_P`u=R1eorA%nh@ z-$FfJ*!X24EV^_`XE7fLZEa#)C4lMd=D(PyL?!7BWca}ZhSfkq4gG~4j^#GI-MSiHUZG+#|4n?Hp6|S8HrjhF zxJP?iz`V)3hQ#(rIg`Syr?v+l3;jyPdzU@J+B(1yNCDyE_CoX%8WvNK|442OVBIB# z!4iBGd;#xoqBQ1ew7C$jy$Cr+VSM3mW1qwb04Wv(RPP<2PVnmkX!C0C$=N_nUwXuVOJ4+{k8~LM z6?of(Bl^ta$d#K(Fr-mWdP}iv6|Ilp!Yh~N*lj(bzG_<}P0=vm8>uSmqKj$ysJ*f- z*SOnxuurLTg77o$iL zI)&*s>L=DRRv?y*$_dY3C!kc8Br@?Z{B@{=wfCr$)M)|b@Nthk#rM8bF0|Srd%bH| zHNXy63}W}A6Z^f)mErW2d!1WyBD5mZmGVngKG&Voa36{xU7mJq^G9@#jsTCUqOvP%a~$cY_eU#|GL|XcZ*Ba1F9|Rx z%}0NN&A+}p32wOvS$VOu<>Lvw*VxD`Pz7^JZPar>RtawzBk;xo^pFEEx(|swzOE5` zhE0ezSkdk)nr40m(YOX4Ni2ugM)oEV*tr|4Y9nx!wa3U%SnJ5+3kbhISm!|~QtIQP z2a(?$|ekokrryo$IL3i7d!Ym$&Pil0JHiKWGfxhiS00hVeBQMs6R z_R&M+m_8U9@Tr;R=d2GzT9;g{J&`-UH2wEi3RC5D4oWMJ35^|?9$=mjU;bzI!kX|G zjcfYf%LWI87SCrKo??Tg<~Rf1*IPzO<`pNB{}_2Pua~qSU#_n;y*TI=v2n5Q_1)^W z7DuJxyDsw!81TB-gY((k0TL>GG8|yEC&7!lhSCrOXlzu8*1^uw*m2i_Xgy2l^KzIA zfHt?*1Gy(&9(=cWE>;k7XO*a6?=lcxYpR!6aS`6co9Avu>t`{+Xy)TvPT#|fubE`_ zm-;Zq?{f7~>r$_2BRR)e`FBk<@C${fx3bS;8JsJ zte*&p^dry#rBN|Dv{ZYZF(ed18q%Tgdi~l%ZJbEUQg3MtcnnNxT z$&Kfqoi@s;nM0gt4RpDG5IFDbA5DkMGuVMVi+7?ZiahSSt1t&IxB5-dn~5= zAvHXDTiV`~0ho9F-aveyOe1}lTSXEbo-}8<&T#u#lIZ?EWL?cdchchrqD6;4{4jz{ zUrerwxVi8)Tg&pkRS^0YrKj5g*<`gJ=j+i$v80pWT1x~8P7!Ud!Rsmu3FJ@E0OuMN z>?3A&4s471LAa$k}2plDs@$YgF-7MiD0WBy;(;4qZAlw*z zGHf}a<%HEMJEFr^3~Kv|q`;f!<&SzJ{1g~#(-{duyF#xuLjwfHy-u(c81nu+Fp&zvu>bj|Ybmk{o*Jq# zI;bd{=VCQO@0;Gg+ozeuR#_^uUD`lvm@JG`P{1N%`!^#XDy#51E5a!keUYQ{b&HV@ zBS-N09BuCr&Fb?ZbJY`Nx8Kx@zmUq`cANp(a)EXxt7(A-x-%#FrZQ3)qj-|lq+Gf! z{VMJ{+MAVxr7DDfzP7#MBqO={77-G#vEOFe(e{=3X_X}0cqaU$^jF1JmiOrU>S z)(;gKH}t0Z)|P1HI5UenFqVq$ zxkb1BVP(c2;0x4YkK1~oys5#G#KG?fnkLbIRMX^hj`MMaSd9YN^!t5vMr&^ z?}JrH_{LH9#*6tKwyw z%6*fIwx-QHJU-HYeA6&y^^N0Xfe1|m7-a*fc4Rr*J>Q|-5QMe_sHbxDiVv4Ojr`vm z!sZKsziKt~z5&gcB{WLKH7f`blJQIQix?JA=Jj8#=zm3zFTp1In<)kUTLmp%!_+Ln z-51(LhK1E1@2g_Vzb`zsu6x(pRsCs_nnygND)x_0KBZ*@#dw-D9`;LJd>u#kGx*<}G)Uq7=ADZ&RwS9Qks0yI#N5(RB(STOsaf*L`>z~n zF@O1`if9NS0th6q*_lZ7BI|2y<;zJ=Bu^5gv@Di9;yP=+5(}s0Uo1%fwI15Pe9xYg ztTH=!m)c`8+45q|VJT>TQfmIEmgz2ymDQskQkrBmKtS8P0D@?$i1pvpbk?|jDzW{Z z=M}#b5F+HOnORS0i7o6RZ(gjAd1EG`A7Bao-8$NT8yz(k)JShpHRWBPhlK*XsCXSv zYc+cF<78(?mNaD{6-SCQ5cjWT{ZlWbXRiE)<@H0w)vl)U^&iPom9Ew5OTI5F%)F@C zbu}m7j2_TYC(=1W>_@^s!fl{s2!#ze_Nv3-LZP(+))>Cb9_UfSEyEm|cR^ zs$Z!g66lWLU|JsbtNDY_md$1sKdx(W4M?9^sBDbsN94Y93Js#)jeY$Y?%TT9xW|eZ zaZY+{lmVZ&P+(@@rA@%k$@d{}a`P6n@+4Y0dMny6lW;<0nh^4(O224;EZu`y$&lP2 z&XB`v1V2A1loTg2a6kkiawb%Hi-z>_HD8zKUgiL74Jz!n1@!A#56Ae{K2uy8NK2Iv zE4`CPledZg`w986&olH9t+uVMYr2QMKZqXv5e;4L@_7>Rt5|eGPetV&HFYeVgnE$8 ziqo@+<_8SVi&|Eob$z-&ze5t|@QG$Yh z^+)fU7S_y0$%z?YN_qBvexvo;^c^c!)@{4U$V|n964Hu4{LyB$<%h`QjL08p|Euys z6*Kjp!|Z1Kxvz0w;gY<%iH*Z08Iq2C4wEXAenbm>6b>pdd7d=tBux1|c3@|l-Aq|! zX_RW=T-N%~>^P-%d46rm>w>>dKZo5GfZ6Q{X5C^Dfx>J$_iGL+nH@C#zt9geywy$6 zFJ7zh8ihkDGWa;Oz%%1!%q#)!ns-Kem>2j`+hN^WHGG3wI>DK--T)GVn>WKgSwh~Y zdjmPA_2*QzAx0LqUGmzHpJ04MSn|o~|wItmz^J{Gyhy3MUDjB`N9iWq@9A5NPeGBjgo9|D}rA`xaSyQrudSi=!JBKzrW@msOY;0}Dkg8X+U zm-JF+axlsQDfOP+i%5aWP8FV>N>!0Zbo1?niE?o7L??8nY_+Q3t{!ciIwy`Bk54QZ z5A!W5c3Icm5?95H>gwIl5hgb`_r+d)s3Ok%L|<99+<~k7_FEHb!vHVHQScZ6?Q@9Xn zpvyi(e5Dg;htpy!S%i}LvJNzPDB4=CK>k5aK|**YUq+HFIp0x|jEUBYAzU1LEPDMG z83tYX<8KzgDgV#+7Z1PIXCaLvS>_U4#`M7-OAqU4cz+n9)Ig_3aiE+2!Ai|+-mD~l z)vnsDVGLGzigFgvscURCn>&+Mu$q=LcQ`%XQm1RXTyDc*#&n_~UPLrwB_lvu$0=4s zzh5`HD%Oa$npBt|R#b4%z_m(Trz!@TspFJvEwHpmD?_tFD8FPHQ@k2rAmUBvn^3;M zk8e2}r`iZV5u0>8DQWU{wsevHy%KPSWLdyh1J&rbUV`@n*Dq~u*6*Ni;&_13w#}s@ z1Wn^gJWzejnUMY{ix?KV{3u2XPES z@%eITqu&{Xrt=ljOy9`Jr!nNqh^O_NTFEmLnPSzc7$u`=W0$KX(0?yj+~d}!{YW4~ z_95*HNGu+E3pHMX0Ejn$H_$>1@dcrD!< zuC0hq;93+;ErVK;ZY&con)I>8i+sk?CW*bWJVYPRmsT9JOw*7{DPuGvVYIpbOii7q zbGnQ4Zv;l#$&++jzh&w;ucq2fr}9F_;eE2;DKgxQ1NZB>2p_M|F4j=uTPc!tzlm_4 zAIV3Lw)V-i`^s&V!0N@TiR#QgImq{+UnmURwoj_-xp&02^=D!oog-v0nd&kzs=%w?tf&O@a3@%yx3 zQjaNKG&+u&`rNG>ibcdZ4`s-rbMU21Y*>DQW`d`*qO++Y=~VcBZiJ(Iho^+(;{=hfRwqPIlI zdxSd7F{iDmx=(ldKG$%;WN%^UWA8UOuI}#38Z#6dJ!Y*=1|z#dg-Twz#`dhT_PD0@ zfVI#Qn^hj^cU(OwtT#fL19r8>hs5eI)8A>`+*NNw-Z~HGikF~%6wRkLR`rpW_$o5g z|HdY!ZRiLfzC0u`BaWGm_t9#vC*Ld05-hpS@NL}|xYEGdxZcLzv-@%mZe4474D^-` z{QS5}A+|C>+rtLn!p~!9(%W7)#ueHlW(c=6hLU z;io{yXWPh9XwUnrnh0y_5~~;OIXTlWLNyKhkL7Vmr8RbJPmBcN6lT5pAjitd(w{2) zvOio;(n~qkTF*DR&FHO_-L`B|4HA<3<+3r`sM(@KAd%j2O2e;Hpn)o>#ZPngF!ij< z8eSg%_;8TIGU;{4(b0X7Zxu}jp`1A}aBw0jsWnJ6h zEmYMel2n>kBVHTGyY7Y~sDSasmes&p6;psR@RKx9eC<##bMbW zQf(eJ2>#f?73XE}k?z>UdVMFbCXuJMUoFoi4n(_kz~;qL={LgmJ+$S_f&x)4F-IuF zF=3DgWDI*4-(4JBu93TNL8FdSCaJ$aP)22oku&#{_L`9gfj+%gV7DK2Tw+pOIiKf8 z9dPpv@wSIIC_LwFW)T09k4G&ZR63GXV&l11!_gza%9W@0{x$G5bPJ?kPPndT6}=)i->s2t{1i9fttyHvkfzgW0} z`o$=6C^bfpX~%w4o%!^|=lVIIFYBkSa&nf3T*5hfs$TwS*(=~Nk4pBiykk)$U76WmX^yAciHq+ne|na>b>t!kMXd_>+?1Xmu5@dNX7=< z*Nm;?WO@*-a@+00{HoORwQ^nq8SJA*k`S(`mU0;p)?kuu>8Sod!Q`E01}fdo_OB|^ z`Q_RIwd7izWB zL5jF``8(R$+^k+wTr^dx4ZU_!>S(T1`$@fz)SMooo4 z0r@d__FMLBP%>db_*yg?$~LWWqTGcF$GA}jPQOWa(|6Htb+vIz*C4G{D~ux^{jPf! zk}AJTH}l>3ga(}vA7DUjCZhDMMK9-;lZKcfmau!4sYxf@IDMq%q8q5eU{wH+iH~FG zC_LNthAw$;@b6i2VTpyp#6%a?qkHQmfCaV-FbU6^AqzLH=q{FuPB7roa2@Uc!-<$ zf*5v>ullMhIc6wlY~_|&#FBbVwZ+LA!?vCu`24j6ghgzluEx#^IU0^f@UwEb0$PaN z>J*e+4k8-5on69K4PltSRepv@U>#=V_Ln+RM2<4?FtOW6d6br!#V5~~lH%OwrTB>$ z`>{(hcXM|BGA=wi&XOVjS|V{E#qCJ z2hW_d+u7x}_jXTF?s+S3mzzPZT7j`3z1!A)@O5GI8;aJ0!dNA=nHa^as-0(`9jH~f z^8Vm@a|(9=EB9T?w`reLDO7&UBxIpJ_q|R^DN)kT)zZ{F;~_P16XVzB3#nm!V%6zV z)3j(Yw6T2s(;n3+-_2@7ep-O@MH}MI=K05`DY?Uvnw&63_yB9EOi|1t=)2^;*M3Sk zv~*0Sb~e6i)ogG5`80&rI%d{t0kHSQIqvEuh!7uP_RZnK)_S8vLXTijp1@ZFZLOjG z`9r&0mH|ZTUejg`kuc<$on$%6Lh+WQ#%`Q)ICIKR3%37QH$D{}osZ%P?qtdqkNEge zp7ex{`aC|oY1v^!Pk zpP$$*ECJ)ty$n{>+s}Jhm=Ec$?k*?Pqie@|DB!Db6<=0N*lxKYV?6+!;j&JR`G{k& zT2k;FD6xa*=R;dyn^AMxZ5~}kG-5ocz<={*_4uvDGt6~X9mIPo|kbZ!8rsn zG`Vg^Do&Wo{0B(x_sqRIv^xGe;N@e<0=*9Qxc~INJ4*Zj}^Yx~?OQ8}5-y z$9ga*jgmKo29t~X>8G9`R9}2C?iv5U-5wo7;tIb@x{GvA10RX$(fs`mbZZMTYiK%o zjW*2s56UD1k$-K^`)2$11;jDt`Bl6MunG>-iI(_T$W~otL7d!dO!!yEKjZv`!p?n% zp?`e7^Vb;r8awKOH{z$b@jX}%>Y}&enyO8Vh%K5()z8kW&)Yr^f#XA@XxB)|A~uA$7jfxe85L4g2_(~H$o7}U!AXEPTBHf9qU<%|oxQnH|QpHukn}bF1n8gS6Po+6fF1S@(z# zdS|v6wU z4hj&}cRb4`8#*+<-}mq9oD);5`^(;378^PEpB(Kt8(k=@fm)h<;?p1Ku7j^C<6f5z zw~|QjVePs(f|thx{+VahFr^yS%>6t^o%ZgL;Z^);EJantybs;B21YYmQP^v5v0QOd zXFK?drUPio!X?XiZ8tjz`dX)1HdSs_kT#5NwLcRKszPclcEhscah^yEWNH7ciCdPv;2kw?`;4r zX==^>cf6jQBrw7RLXrcYn7f=rMO~ z)u1rWORG#9`!Q9%M{60@pGN7;Lla`}QZ~V0@SIHBUX5p<&fbZbjqMGk=qprw(hb3q zD?;Q?I~CU3mw=PirW74hoqc=-O7QmraYnPOanfDzx>?r4Q>M~zX5Wh0j~ zr3zll{z~h&OP50j%(e!Up?dgjfv!v8ewo^;5Lp~;-N9?Y?%&Lr(-I3WHK7laqwf&fro3+2yKP%#uQJFN@ z{2AYHAdtvGQUt|sJ}gjq)?&}e?z&czWLoSy)JV!jiBmj6*k7+N=QtNGE1Z^UG*dl) zx^eD$KjRl_==rYXja*!aAbw9w=q#PJulMU5<^a#4IAr1-`(S8)lqS&JA^XJo=#uTc3 z85pLt?gpD0h>&bjE>hO8lbGBs17*+M4KtoOisy~@{D{WUnY$Z=)W&l+B*wm$RIDt4 z2x&+Jy&9yRdp+0Ud9ny0jP>01Pg9G0^cvhG-u9|HB=0e$p1hE7rpx7=SpP`h+mb4h zfDWLfvPE- zr*JwPiJnDJ`(Sb3-*+tX-BR1czKtepp&#BtO4Kx`Uj(r=~D9qAEm z=fm~gUVB{E^_fJv{_Aql8UHX5vGcVdHGx<)%Lw`|h^9qlFh-Fy}sgtF~?cdu3M=%4ckMvX29e*U1!3G3$}l1LsH1nvPm}GDmVG z%|WW%MM=diTa6JHi#}AHErqRL2cG#jB#Jt=vOtI~LsT<(EO^GQ4)a>~v4wX*w(>7B zX*D$e*5TdJ*xZGm2%949t);}iF%d{POPJ+&edXm_3#!BzjnFjpX+&0OgokDGGRNh> z%{U3z&fehMewCD=%eWag+A0G3D!2-uJ;p+7Z&ocf@mm@OuCPYgc&omRqDHbyi*&D- z7Gy_q#?dp#5ybYNa3@qz9E`lkUP(-!xHUJX6PEh|7IOFDFHy=P$}FxtL;N zWcMa%LiCVK%DUpxU6jqL?DC9p>F1|&nUwt(+qD>FL;8w2f()_v zszk5N?c^!20xC(w?7?30+8UbDc@diGM97F0;d0UlV{4k^Lm5ipSI&2IMx)`27js#i zd`^__m$jee+~FUz?;r0qt*o)_UAO0sC41ysn_IlIVfqE5D!MePZ*o z@sE0rq70#?tAccq+}TW3boAc%P*DH?G;;e*jrGs?`P66s1uC^a7mNAZ)~pbc%_e4J z$ITK7;~7M5Rjqh2)3|~xe{1q${h;;mLFKj68Rr3LlNIDeX75eHvk>b*p_@XZq`umZ z?Hm8vER2XtLe|kseCBbnL=l&{=7&ULZyhe`qC^7fciNjxK*`LKJvLg72_WV2RcuZJeA9Rt;xE!#mw z!m8%w;X6TY4XqJpNSSmu-L>vvgQbf{(6*z&(jAVmT?rBO+Uj^gxzc`1+27n~U?nc@ zFU(T;_J7PYVMWDU(mJX{#5r&i-fwxv2ni!6uh+)bzG6$~pu5C_`QBbb^byTKpcNT6 z=|#7#pRmH8Mb7$PG^sim1?6o*N_*WD7lH?wV> zq;4jun!#^oA7z7wxrE-Z7x;r1ZkwgH*>tEv){E1hMhy+OShHzd}8zS z{<~XEx1hydWQw$p!FDZ+t_~#CtVcC0LRd=G?s#P%V^YW90IErxKzPYrz=YGg7Q51z z?O2vd7yc|b(%ow%>V{9xLoBO_BEU78=C9@BViuvsYoCL*O~SplOGnj9vAF5vM6+kp zKC1jUz3b|XRWy+pa(1+*g4{)(475Yh+%~hEi!@|)l_GNp>3xp8G3gx`WI}lG%eT?R zkm&+nk0otz#4LapuA2&{wG`waeZ>wS<@d|LrkgjG1O+9gxM~1==Wa$Et*@=a&T&Xl zX2L<*tcbb74qfbhV+;NHcbwAxl4jiMVj}Yd^+dX(Y%^XC(n!|1Wfnm%HpqxrT2lM1 ze_A+nIh_IPhD^eC_@>if-hJ3|*8P@|jJgG;RZ=P~s)V7FEk~$CU`tU(VJxspO{N+} zM@yDsSfDkl^rZNN)}8T>A^OtEl=JFA4+8x9*U3!lMt)%n+@sx8&78-i7~9~ttSE?Q z0Ox9r)JKI@Z*VQ9)4E^=;;?M$z624c`BUb4XAVhqb%-C-q zp_(}(+35xD+jA{CUDsPn|A3VYyy|=bi~JJnhN0y(>v^8vysCC{5Pc%B{j((?Y31W;Ok)? zvXPpAwgxLP6m@#9 z&*@HlUVcy+DnRdk)-vT{v~(L{pxsMe7=)0HCM@~E;u<#2d+{8?Y4T*JcZ&2@^B#DP`VOU>fY!_QsM4EYB|;G z!Tj=RPZN%ur^ggcOWhEZ#v9F;IFir3m&tdWmVmA{GyV3YU#SA$J>mE~K3L#M>!b`s z1-wpkTZ>&)2AQR|r-HFPKyb&I&2icj7L)&Y3U`uC$yz6c=4&WYiOP9N&&a_=)o8*` zr}x{K!-3R76oHaOL^1Z-xmM6Vi&2WR#PJIh^3dcgko9B3Z3=x^Q|~oqVh)RI+v$w> zyOs+ETMKNTg|AwvG5I|Gf?sA~mb(W-@bDUc{IhJYf!JzFo{=O{Sl}H1M?qpIWNa7PCz^v zTT(r(asDi7o~sVl*zoJ01=;S5#)Ufwx@|hvoq&1OFt)z zbV+>rrZ4f+(n>Nc%9~N;^UqPWItRuDNgBA>aoaVL%fRCh= zRZqw^B%gla3cn5O#)-YfvGqQ~`<(8+Wq_Pd=wJb!8tdi1eC_?K2KRBYsmJcZY5)dK zkMS!ntdN%5+-%obWcb}MVSs1gNKaU2_;j6thlE3r0V@yC80@H~S`;P^V^&bo!S2Pt z?IZr;ER|r`vv-zMxw4wXU2~W1%29?|u~=EQpFb+nt}*t=Xv)Z);8KMvA|cBLUu+F< zteTh>L*b`u3$j|$nJZQu__`3Ov|FaHTpS8dpzHkD1({N+fpls+XAD4`D2pO@9fDt`L%kYtF_OaEf{w4 zl(GF^5?fX@Wpvl~IcdI2QQx5Jy`=3P3(52QIVL>|rC+NiY0ty~9j|UiqH8ff=>6pi z*!6^UDe5JIPZGBT)>^lM^;bRlfDW@hENkh-jG?;H7D5WTE040eNomoHf@jkSWOEa+ z!=kmwmb%r0-{Rb87kP zON$Fxr1O3-y#~{->2s6XLTsM5vX7dqA+#$|fd%*l3!8J8%7mW&YbrZSM-J1Y#6u@) z$ltiC)887+%gFs3rAU4XF9Cy7{P7lEPo@U)PK)^9W@HI1mc(y zz-24I7+w`vl5)OKsPhtW7TCAom+bwvF|zA{4-ZFnN^Oq#jp@$B0WA-vn|VOzLpA}C zbPv-^gv-y$z46711v~h zk6*4L9ZiH@o>o8L$g4$~^$je8h@Zz|UyL*J?veowd4jX+g~am>YtF5JV{lk*umdNAtYQa9(`_fcFU!b znr{T|hSPQujEwi$+(blUOEJqkfkCD=Ix(JKYA*W1VXN-eJqEBcO`dNiiV1jRk*_|I z;|ev?%mVKQKKKq~*Lvc8=YXz721q?!w}@ZbioyLtSTohzoUE1+A1+t{kW3lS1!Id8 zKA#+#WC*5+sfd5aY#g*eWIr?*G9T*0F}A}$edvon$I%k_d`p}m zVznl4OxJa8%m16*ihAjH*tWYxh^~PwGTMwr@1^IvT&-o(UBq?W*S974!)~>Pk2>5Lnt_6L;O@e3` zrmuk@bDzJe%}`F$5YmvOJJ}!&dyUlPyX}nbJ`TIZq3{~86<2Z`p581+*LnuM%fQOQ!a(48$ztoSJ z_E<&(6cpBb+hIckgeyPDsmw;jTfcQywG^_bH7W+{`z7S8kCade}1S^E<)L zf8Eu}Q+vxRFTna}MfKofsATq6QzkGVR|oo{Iuf=jqP1?I1ta!)T?#r(ASlNa0Z9t| zXY=Xi!dF(zFkAfN`mDsrA{{UQ;xAaxaf>Ilh2Y}~V9g|nzO>F?1~s_YE!7d3_sl~6 zcu!@wtjO85>|t{>?m>6Tgeoq_JA^tGDa8A?6OBx}bqNuJrMB4@z!1F|EO_PjA5!y7 z#y=Rz_}sS&vhXiltZpt)|BYE?Om zUk9iPF;VM%9oj&Y;KzSv`t||0sQ1aYZtu+`$euWx&UK@a#g@8XZWC(&@jtXD+eUvm z*l92AK2m8{kqE47U^+oBN`Yvw`%(g1%~MEodxNWwr#-VQT<_-#I(keb8TMdd(PXx8 z&cOETsER8J;hMBGfL2*wrFPyevzhNlZgtJ`kVO)ibJc%I)AL{O-k(Le2HmR2rs95CM$szM529kHXF78T%;9C?Abo}su zs<-6gNu!p8UvA(+04NnN-;$=aPsyJF$dO?%w$ZvYZ|G`J_&SdW<#MiWo}!+aPLd&` zu@czwd}cAI_~=3IX&yQYYREf#`)=tKF3FF4BRE|2A_^;0xD)p0a~U*oa_em7nLgD0 zY8me~{QaeWmd|-Sf69huE}^H}(RB$1vfq>W^yXjF36S8j%X06e*uf*#w6<-vB|)_K z{57OMgrdfEc-u7B80NfLMriC~bW1FTBQ}4*{QDdI3#%9d5+DFQK4Da%ss*?ww|J=W z#`4N8%UyO6dL_N5WorE`EK4_Lp3bj`4^;jll(0Zcp*EK*~DJqomsoY~s5l zBUFFb+?#=}Mey$%UC8g+(FIYaZ>h-ww1cRMN6&Qab@?QGv|#MYS_R6&DsAQkk=AZX zidsqy->dRPHt!taIC)_YWW@USfxrlYO%5odf}~Ocl(r(B@o=2M0KI|&R>qOjd%pTz z(^LU}Zt^GFF84%eD5jC+eOad6p4_Lr2YCF)_d|BxQKFIllM18_|{nczZ$3N^&kjv25 zyzi4HQ1kNw`X?p~;kDFoS-Vj!adELPjkcw!wMwxb0M^=P>~`LB`1-*NJbq{b(IVZi zOQPj(ja{+A3M%0Oa`H7cZhSvs6}MhsmkJ7fLPyzmmQ%Gat<%_cx#q%up5G9MV>;UG z$t~MG0O%{xbdC!f#=x6Q$6WDEMVAL*VN6GVs zi+3Vs&&xcsG=Br8U*|zlKVigr@WduWmk&=8J z`*FVd=8(l_aBC;x2k!UTX}g^Co_~ng0@PYi;bB~Q>4)rG$?Mq6Pl7{?s_@+vT_`HU zV6BwAl*EXH<~V8H1%!)ZYg{msSPo{^8vXqSv@zTsKp%J3&UfF+Va#E zLczg@?w3+0VTw#CzH_tzREL{y@wXqYdiLYzYoaB+(x|q68Rf||svk&qKi_nK=SZn= zqt#i_OveNHWwwXi4wd8hfd_fwBcJGhUH1E+%Adp&dWaBcP$5iPy;Rl<-g=9dF!1mz z+1y-|PMPMDIG^--s&6R3zCTM>`X@4vfA0ltD^+AsJ}>LQW26(Gs=j{uZ=4Oh zRCHeL9`H1Mq6&KxSmg$8qP1Juerh6BVGE8NPKw8NeT&qwk=aU{zrzE?9`?CZhDAGv zuXknAgJ29~e17PMU)$>K;r)8D>h9I$I&0pjs^Q`o!}rznboPLyg0Rh>%xNGxmH>wT z)6#@7vBTdHHH<&{B$+TL-r(AR774xmL)x4Dd~+fK{!TEp&R^t87l=ygEg@p&?mSIZ z&@BJ;=dnnJgzT+)auzF);~_&rw!jro5!y>Yl_5MPb33V{=3^HY9wtm z;_uwk8T#Gv-5$328<9er$@%XcWdB29>%)up-!|CR>pe@1x_(5imTjv5ZfvBxnX>j_WvfQf&kNubM)-wtd@lQ0VcA)rNZ?)2@7ME?3; zN`PT$i}2G`Q+0v1g+H{BH_vaMELj~WFZ2F;DKTYd`~|$Pi2&v&WYI9KU3^ ztfJ70&pP8)cZ)C|e{n&rDm;k-7`$kKH)B!mYgcWvTB;x*)_cQ#FUs4P@%VwRTQ5~@ z32i29$JWMBLH##BVH?};;E&*3*zs<3@7YNZmQrOL3-x9@n@f-q6CX=7 zJk|%%bHyd|#MfI2Ebnu)kLN6;(En}@4C}d_loJu4ock3w zi39}4JQx)r+3FpPb<~^a7V7f)&N^zM9#&?3yyLBF)L;Ua6LWlq04SwaE~ehP0U2>{ zB{_hK**rj}X83nc8j0S1Do0@ca_faKxGjZaM;e8UTDAVjHyxJ=0Cip0cw)7p9Ngj& zb9aSR0s+C{|9()u`Ye=ap`deC8dd)_0KBJE<(NoZUSp$gYHNqFl&73eGw1Q+V2HPo z`hqkkM}Xgx@h|7UZaPc+9udEImYMqF32uIv$D5$Tfya)qQi>gt{<-aIMxZ3C{T{%>y+rA)r9YMPOw0;wRltiE6{~EdB_l`NZS6Kn(?pe-7 ze!om8unb{aMQ-{Hm2X1l;^>1iVm~M&Ixek*-My-Yp_eNBc$HOU78qXl`>D2egd%5^ zjjwYii{j=0saXhb&Vo(c(&crbNl}wbK-_~mzAOy21n2w_rEp>Wo4~I91A8=)d05y^ z!|2@UlojpT55B+Slm39`fLN!!IiG&UifAZua>dC=6jROog3}{y>@zqUu(Rx;qvTzmnaX^_g62 z89Gig zJ2!ZvkRZ~Y=PzhsL$XWd&Zlepoc?3E!M5;=H<0W02RknL?~Z~mQiifepgFKr8HN%_ zf`!w|_Nu8om1ZAQU4ILWA*plU@oEq|WK$lyuRho|cHE(kbm`JiL&)ktZTaKJr=A_FZxgDzZJla=I`Qyo4hSf%wY%L2yKa{>8_pvxS}rrR9g z0eLI?U;qLzA+kO7p!XXIRJn5{WQTg(13G*e{UGg#AqYx&oo}<8#^}G5*8@S_lp;Ry zt-&q_;<;qIhNghcN_;r-vIxN`M^@_P;?9VQJ2KbS>h#ug>hB7B*A+@PWhmX`apGBb z0F(@b>ok+dx97c~+w1du)PLcydM#$u>h7{b^4~H%!Pjm=?-ou7rK~}^uF*pu0~*%K zrxA$KUpU5lX*TE$pf6ym0O*TO(T2{?=iG1Kt~!jZmfHc_nUq&n_w`ZJ;-3aFqHOPW z{*Qzi4?M|aTs(*VO5suUHF-M0Pn)Tz$IB9o5cBmAnCJc~iMf2&qI!+zVnx&VCUt`& z#d&@FEBh$)l{WlRAjO~|?~?He zpDKimM8H;CUe3Kx7;$Xhp!QmtF>d=BB!=h;xh{!F>%sLbOwQTGNXtBzD!&*rGM__q z#(FpYGc45fsYho^G$ZXV_T8h;B9wCE)juC~JihS1$2i-FzYUbs?Ajk=YuJsyoN0a3N^<`yn^csXi~2;UUU0oX zi|b3fP}gS}^|N~VZY5X!qX<^-s|`ns?!@oI4bR?2{LgxS*d3KzjcfM k)bsytU;Kamyc4Hn&RUNna4`dRN