From fb41fe86c010dd2cb7ae11ea876ee20c8051b5f6 Mon Sep 17 00:00:00 2001 From: vijay gopal chilkuri Date: Wed, 29 Apr 2020 16:05:02 +0200 Subject: [PATCH 01/20] Added a small example for tJ with 4 sites. --- examples/tJ_4.inp | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) create mode 100644 examples/tJ_4.inp diff --git a/examples/tJ_4.inp b/examples/tJ_4.inp new file mode 100644 index 0000000..49976a9 --- /dev/null +++ b/examples/tJ_4.inp @@ -0,0 +1,30 @@ +4 +2 +1 +1 +0 +false +1,2,3 +2,3,4 +1,1,1 +1.0,-1.0,0.0000 +1.0,-1.0,0.0000 +-1.00,0.0,0.00 +0.0,0.0,0.0,0.0 +2 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +0 +0 +0 From 6a379ba977d086ca5916189878c7a78045a0e311 Mon Sep 17 00:00:00 2001 From: vijay gopal chilkuri Date: Wed, 29 Apr 2020 16:06:13 +0200 Subject: [PATCH 02/20] Added a sample output for the small tJ ham with 4 sites and 1 hole. --- examples/tJ_4.out | 85 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 85 insertions(+) create mode 100644 examples/tJ_4.out diff --git a/examples/tJ_4.out b/examples/tJ_4.out new file mode 100644 index 0000000..b7e7de8 --- /dev/null +++ b/examples/tJ_4.out @@ -0,0 +1,85 @@ + +1-D t-J Eigenproblem, n=12 + + start: 0 end: 12 + HAMILTONIEN t-J + Le nombre de trou est : 1 + Famille 1 : F + + + LECTURE DES ATOMES, DES LIAISONS, DES INTEGRALES + + + ================ CLUSTER 1 ================== + + Liaisons entre les atomes 3 + Les atomes 1 et 2 forment la liaison 1 qui est de type 1 + Les atomes 2 et 3 forment la liaison 2 qui est de type 1 + Les atomes 3 et 4 forment la liaison 3 qui est de type 1 + ============================================= + Le nombre total d atomes est 4 + ============================================= + Nombre de J differents 3 + type de liaison 1 1 + type de liaison 2 1 + type de liaison 3 1 + type de J 1.0000000000000000 + type de J -1.0000000000000000 + type de J 0.0000000000000000 + Parametres : Jz 1 = 1.0000000000000000 + Parametres : Jxy 1 = 1.0000000000000000 + Parametre : t 1 = -1.0000000000000000 + Parametres : Jz 2 = 1.0000000000000000 + Parametres : Jxy 2 = 1.0000000000000000 + Parametre : t 2 = -1.0000000000000000 + Parametres : Jz 3 = 1.0000000000000000 + Parametres : Jxy 3 = 1.0000000000000000 + Parametre : t 3 = -1.0000000000000000 + coucoudslect3 + coucou + Parametres pour le t-J + xj1 = 1.0000000000000000E-002 + xj2 = 0.0000000000000000 + xt1 = -0.20000000000000001 + xt2 = 0.0000000000000000 + xv1 = 0.0000000000000000 + xv2 = 0.0000000000000000 + xv3 = 0.0000000000000000 + xbj = 0.0000000000000000 + xbt = 0.0000000000000000 + xeneparJ = 0.0000000000000000 + xeneperpJ = 0.0000000000000000 + xeneparT = 0.0000000000000000 + xeneperpT = 0.0000000000000000 + xenediagJ = 0.0000000000000000 + xenediagT = 0.0000000000000000 + xspar = -0.0000000000000000 + xsperp = -0.0000000000000000 + Le systeme comporte 0 plaquettes. + Spin total 0 + Nombre de vecteurs demande 8 + Nombre maximal d iterations de Davidson 280 + Vecteur calcule le plus bas 1 + Variable Nes4 (vecteurs d essai) 0 + Nombre de determinants en donnees 100 + Variable Ysuiv (suivre le vecteur initial) F + Seuil au dela duquel seront ecrits les vecteurs 9.9999999392252903E-009 + Option d ecriture des determinants sur FIL2 T + =======nombre de centres de spin alpha===== 2 + nt2= 3 + nt1 4 + nt2= 3 nbeta= 1 + Time used to build the matrix: 0.007272 +time = 0.007272 mpiid = 0 + Time used to assemble the matrix: 0.000024 + Time used: 0.000336 + Number of iterations of the method: 1 + Solution method: krylovschur + + Number of converged eigenvalues: 2 + Stopping condition: tol=1e-09, maxit=10000000 + k ||Ax-kx||/||kx|| + ----------------- ----------------- ------------------ + -3.855773 6.4802e-16 0.500000 0.000000 0.000000 1.000000 + -3.414214 1.59311e-15 0.500000 0.000000 0.000000 1.000000 + From 850fc2d06a674834bc931aac047c3c086cc25e29 Mon Sep 17 00:00:00 2001 From: vijay gopal chilkuri Date: Wed, 29 Apr 2020 16:10:18 +0200 Subject: [PATCH 03/20] Added example output and input for tJ ham 4 sites cycl case. --- examples/tJ_cycl_4.inp | 30 +++++++++++++ examples/tJ_cycl_4.out | 100 +++++++++++++++++++++++++++++++++++++++++ 2 files changed, 130 insertions(+) create mode 100644 examples/tJ_cycl_4.inp create mode 100644 examples/tJ_cycl_4.out diff --git a/examples/tJ_cycl_4.inp b/examples/tJ_cycl_4.inp new file mode 100644 index 0000000..4258df4 --- /dev/null +++ b/examples/tJ_cycl_4.inp @@ -0,0 +1,30 @@ +4 +2 +1 +1 +0 +false +1,2,3,1 +2,3,4,4 +1,1,1,1 +1.0,-1.0,0.0000 +1.0,-1.0,0.0000 +-1.00,0.0,0.00 +1.0,1.0,1.0,1.0 +12 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +0 +0 +0 diff --git a/examples/tJ_cycl_4.out b/examples/tJ_cycl_4.out new file mode 100644 index 0000000..b470a50 --- /dev/null +++ b/examples/tJ_cycl_4.out @@ -0,0 +1,100 @@ + +1-D t-J Eigenproblem, n=12 + + start: 0 end: 12 + HAMILTONIEN t-J + Le nombre de trou est : 1 + Famille 1 : F + + + LECTURE DES ATOMES, DES LIAISONS, DES INTEGRALES + + + ================ CLUSTER 1 ================== + + Liaisons entre les atomes 4 + Les atomes 1 et 2 forment la liaison 1 qui est de type 1 + Les atomes 2 et 3 forment la liaison 2 qui est de type 1 + Les atomes 3 et 4 forment la liaison 3 qui est de type 1 + Les atomes 1 et 4 forment la liaison 4 qui est de type 1 + ============================================= + Le nombre total d atomes est 4 + ============================================= + Nombre de J differents 3 + type de liaison 1 1 + type de liaison 2 1 + type de liaison 3 1 + type de liaison 4 1 + type de J 1.0000000000000000 + type de J -1.0000000000000000 + type de J 0.0000000000000000 + Parametres : Jz 1 = 1.0000000000000000 + Parametres : Jxy 1 = 1.0000000000000000 + Parametre : t 1 = -1.0000000000000000 + Parametres : Jz 2 = 1.0000000000000000 + Parametres : Jxy 2 = 1.0000000000000000 + Parametre : t 2 = -1.0000000000000000 + Parametres : Jz 3 = 1.0000000000000000 + Parametres : Jxy 3 = 1.0000000000000000 + Parametre : t 3 = -1.0000000000000000 + Parametres : Jz 4 = 1.0000000000000000 + Parametres : Jxy 4 = 1.0000000000000000 + Parametre : t 4 = -1.0000000000000000 + coucoudslect3 + coucou + Parametres pour le t-J + xj1 = 1.0000000000000000E-002 + xj2 = 0.0000000000000000 + xt1 = -0.20000000000000001 + xt2 = 0.0000000000000000 + xv1 = 0.0000000000000000 + xv2 = 0.0000000000000000 + xv3 = 0.0000000000000000 + xbj = 0.0000000000000000 + xbt = 0.0000000000000000 + xeneparJ = 0.0000000000000000 + xeneperpJ = 0.0000000000000000 + xeneparT = 0.0000000000000000 + xeneperpT = 0.0000000000000000 + xenediagJ = 0.0000000000000000 + xenediagT = 0.0000000000000000 + xspar = -0.0000000000000000 + xsperp = -0.0000000000000000 + Le systeme comporte 0 plaquettes. + Spin total 0 + Nombre de vecteurs demande 8 + Nombre maximal d iterations de Davidson 280 + Vecteur calcule le plus bas 1 + Variable Nes4 (vecteurs d essai) 0 + Nombre de determinants en donnees 100 + Variable Ysuiv (suivre le vecteur initial) F + Seuil au dela duquel seront ecrits les vecteurs 9.9999999392252903E-009 + Option d ecriture des determinants sur FIL2 T + =======nombre de centres de spin alpha===== 2 + nt2= 3 + nt1 4 + nt2= 3 nbeta= 1 + Time used to build the matrix: 0.007798 +time = 0.007798 mpiid = 0 + Time used to assemble the matrix: 0.000026 + Time used: 0.000336 + Number of iterations of the method: 1 + Solution method: krylovschur + + Number of converged eigenvalues: 12 + Stopping condition: tol=1e-09, maxit=10000000 + k ||Ax-kx||/||kx|| + ----------------- ----------------- ------------------ + -1.000000 8.49044e-16 0.500000 0.000000 0.000000 1.000000 + -1.000000 9.30854e-16 0.500000 0.000000 0.000000 1.000000 + -1.000000 2.25866e-15 0.500000 0.000000 0.000000 1.000000 + 1.000000 1.44702e-15 0.510159 0.000000 0.000000 1.000000 + 1.000000 1.19807e-15 0.605082 0.000000 0.000000 1.000000 + 1.000000 6.83261e-16 1.438653 0.000000 0.000000 1.000000 + 3.000000 7.89504e-16 0.822335 0.000000 0.000000 1.000000 + 3.000000 5.25436e-16 1.007741 0.000000 0.000000 1.000000 + 3.000000 9.298e-17 1.374264 0.000000 0.000000 1.000000 + 3.000000 1.62369e-16 0.876840 0.000000 0.000000 1.000000 + 3.000000 2.12039e-16 0.752834 0.000000 0.000000 1.000000 + 5.000000 4.04585e-16 1.500000 0.000000 0.000000 1.000000 + From 6a456a13ba59c8d285d2a41c7675f12857f4b923 Mon Sep 17 00:00:00 2001 From: vijay gopal chilkuri Date: Wed, 29 Apr 2020 17:28:58 +0200 Subject: [PATCH 04/20] Changes in .gitignore for MacOS. --- .gitignore | 1 + 1 file changed, 1 insertion(+) diff --git a/.gitignore b/.gitignore index 5b2aa47..94f6eea 100644 --- a/.gitignore +++ b/.gitignore @@ -2,3 +2,4 @@ */obj/* *mod.mod* irpf90.a +*DS_Store* From 0f145b96418842c6d8ab2b87e0a03b90140f5476 Mon Sep 17 00:00:00 2001 From: vijay gopal chilkuri Date: Fri, 1 May 2020 13:16:05 +0200 Subject: [PATCH 05/20] Added bigger examples and improved notebook functions. --- examples/tJ_4.inp | 2 +- examples/tJ_4x5.inp | 30 +++++++++++++++++ examples/tJ_4x6.inp | 30 +++++++++++++++++ examples/tJ_cycl_4.inp | 2 +- examples/tJ_cycl_4.out | 8 ++--- notebooks/t_J_Model.ipynb | 68 +++++++++++++++++++++------------------ src/ex1.c | 20 +++++++----- 7 files changed, 115 insertions(+), 45 deletions(-) create mode 100644 examples/tJ_4x5.inp create mode 100644 examples/tJ_4x6.inp diff --git a/examples/tJ_4.inp b/examples/tJ_4.inp index 49976a9..3d35ce0 100644 --- a/examples/tJ_4.inp +++ b/examples/tJ_4.inp @@ -11,7 +11,7 @@ false 1.0,-1.0,0.0000 -1.00,0.0,0.00 0.0,0.0,0.0,0.0 -2 +12 1 1 1 diff --git a/examples/tJ_4x5.inp b/examples/tJ_4x5.inp new file mode 100644 index 0000000..ff58fe9 --- /dev/null +++ b/examples/tJ_4x5.inp @@ -0,0 +1,30 @@ +20 +26 +4 +1 +0 +false +1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20 +2, 3, 4, 1, 6, 7, 8, 5, 10, 11, 12, 9, 14, 15, 16, 13, 18, 19, 20, 17, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 1, 2, 3, 4 +1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 +1.0,-1.0,0.00 +1.0,-1.0,0.00 +-1.000,0.0,0.00 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +2 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +0 +0 +0 diff --git a/examples/tJ_4x6.inp b/examples/tJ_4x6.inp new file mode 100644 index 0000000..214788e --- /dev/null +++ b/examples/tJ_4x6.inp @@ -0,0 +1,30 @@ +24 +7 +4 +1 +0 +false +1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24 +2, 3, 4, 1, 6, 7, 8, 5, 10, 11, 12, 9, 14, 15, 16, 13, 18, 19, 20, 17, 22, 23, 24, 21, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 1, 2, 3, 4 +1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 +1.0,-1.0,0.00 +1.0,-1.0,0.00 +-1.000,0.0,0.00 +0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 +2 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +1 +0 +0 +0 diff --git a/examples/tJ_cycl_4.inp b/examples/tJ_cycl_4.inp index 4258df4..69f61d8 100644 --- a/examples/tJ_cycl_4.inp +++ b/examples/tJ_cycl_4.inp @@ -11,7 +11,7 @@ false 1.0,-1.0,0.0000 -1.00,0.0,0.00 1.0,1.0,1.0,1.0 -12 +120 1 1 1 diff --git a/examples/tJ_cycl_4.out b/examples/tJ_cycl_4.out index b470a50..198417b 100644 --- a/examples/tJ_cycl_4.out +++ b/examples/tJ_cycl_4.out @@ -74,10 +74,10 @@ nt2= 3 nt1 4 nt2= 3 nbeta= 1 - Time used to build the matrix: 0.007798 -time = 0.007798 mpiid = 0 - Time used to assemble the matrix: 0.000026 - Time used: 0.000336 + Time used to build the matrix: 0.009882 +time = 0.009882 mpiid = 0 + Time used to assemble the matrix: 0.000258 + Time used: 0.004857 Number of iterations of the method: 1 Solution method: krylovschur diff --git a/notebooks/t_J_Model.ipynb b/notebooks/t_J_Model.ipynb index d10d2f1..60cf8bb 100644 --- a/notebooks/t_J_Model.ipynb +++ b/notebooks/t_J_Model.ipynb @@ -160,12 +160,12 @@ }, { "cell_type": "code", - "execution_count": 150, + "execution_count": 128, "metadata": {}, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb4AAAEuCAYAAADx63eqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3deVwU9f8H8NeyKLAsopSComKaeCCHJ9nhkeKRllmUmlZmh6mVlaX2zbOsTDPFNE1Ns9O0TMsDw7zKr+ZtWWpGoKIoKvd97Of3xycQZJdjf+zM7Hdez8djH7EzQ/v2zey8Zz6fz3zGIIQQICIi0gkXtQMgIiJSEgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgufk1i1apXaITillStXqh2CU2Le7LNixQq1Q6AqMAghhNpBUOUMBgP4p6o+5s0+zJt9mDfnwCs+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFeOMGTNmqB0EVa5+/fro0qWL2mE4nZSUFPTv31/tMJyOEAI9e/ZUOwynw7w5B97OQP/TcnNz4e7urnYYRKQhLHxERKQr7OMjIiJdYeEjIiJdYeEjIiJdYeEjIiJdYeHToMzMTEybNg1BQUHw9vZG/fr1cdttt+GTTz5ROzRNKywsxEcffYR+/fohJCQEISEh6N+/P5YuXYqCggK1w3NKzzzzjNohaFZRURE++ugjTJ06FXv37i2zbtasWSpFRVXBUZ0aNGjQIAwePBi9e/fG2rVrkZWVhaFDh2LWrFnw9/fH22+/rXaImjRs2DDUrVsXjz/+OBo3bgwASEhIwOrVq5GcnIyvv/5a5Qi1KTk52epyIQRCQ0ORkJCgcETO4amnnkJ2dja6dOmCzz77DN27d8f7778PAOjQoQOOHDmicoRkCwufBoWGhuL48eMl7zt37oyDBw/CYrGgbdu2OHXqlIrRaVdgYCD++uuvaq/TO6PRiICAgDKP0yl+vM6FCxeQn5+vYnTaFRISgt9++w2AbG0YO3Ysrl69iq+++gq33XYbjh49qnKEZAubOjXI09MTv/zyCwDg+++/h4+PDwDAxcWFz/qqgI+PD9atWweLxVKyzGKx4Ouvv0a9evVUjEzbmjdvjl27diEuLq7k9c8//yAuLg6+vr5qh6dZpU8IXF1dsWzZMoSGhuLuu+9GZmamipFRZVj4NGjJkiV4+eWXUa9ePcyZMwcLFy4EAFy5cgXjxo1TOTrtWrNmDb755hv4+fkhMDAQgYGB8PPzw/r167FmzRq1w9OsF198ESkpKVbXTZw4UeFonEenTp0QHR1dZtn06dPxxBNPID4+Xp2gqErY1KlRJ0+exIULF3DbbbfBbDaXLI+Ojka/fv1UjEzbfv31VxgMBrRo0QKnTp3Cvn370LZtW9xzzz1qh6ZpBw4cgMFgQOfOnfHnn38iOjoarVu3Zt4qwbw5JxY+DVq4cCE+/PBDtG7dGseOHUNUVBQGDRoEgJ3mFZk5cya2bt2KwsJCRERE4MCBA+jRowdiYmLQt29fvP7662qHqEk35u3XX39Fz549mbdKMG9OTJDmtGvXTmRkZAghhIiLixMdO3YUCxYsEEIIERYWpmZomtauXTtRWFgosrKyhJeXl0hLSxNCCJGdnS2Cg4NVjk67mDf7MG/Oy1XtwkvlWSyWkubNZs2aYdeuXYiMjMTZs2c5uKUCrq6uMBqNMJlMaNGiBerUqQMA8PDwgIsLu7NtYd7sw7w5L/51NMjX1xfHjh0reW82m7Fp0yZcvXoVv//+u4qRaVvt2rWRnZ0NADh8+HDJ8rS0NB6IKsC82Yd5c17s49OghIQEuLq6ws/Pr9y6vXv34o477lAhKu3Ly8uDm5tbueVXr15FYmIigoODVYhK+5g3+zBvzouFj4iIdIXX404iNTVV7RCcUm5urtohOKWMjAy1Q3BKzJtzYOFzEj179lQ7BKfUsGFDtUNwSiEhIWqH4JSYN+fApk4nUTx3IlUP82Yf5s0+zJtz4BUfERHpCgsfERHpCgsfERHpCgsfERHpCgsfERHpCgufRqWkAPPmAS1aAHLazjzUrQv06AFs3gwUFakcoFZdvgy8+SYQEAB4eiIXAHx8gH79gB07AI64syohAXjtNcDfHzCZACAXPj7AAw8A//0v02ZTbCwwfjzQsCHg4YEcALj5ZmD4cIBPUdEs3s6gMenpwNixwLffAi4uwL9TAZZhNgMeHsCsWcAzzygfoyZduSKTsXUrYDAA1m5cN5sBb2/gvfeAoUOVj1GDzp0DnnoK2LNHvs/LK7veYJCF0NcXWLQI6N9f+Rg16fRpmbhDh+RZaEFB2fVGI+DmBjRrBixdCtx1lyphknUsfBqSmCi/HwkJ5Q9A1phMwKhRwMKF8gClW//8A9x5J3D1avkDkDUmEzBpEjBtmuNj07Dff5ctCGlpVWtB8PCQ5wxjxzo8NG3btw/o2xfIzKzapbCHB7ByJU+2NISFTyMyM4GOHeUxvLCw6r/n6Qm8/DLwxhuOi03Trl4FQkOBS5cAi6Xqv2cyAXPmAOPGOS42DTt3DggLk03q1eHhAaxaBQwZ4pi4NO/kSSA8HKju1GQeHsDGjUBEhGPiomphH59GTJ0KnD0LjB4NHDwoW+pWrbq+vk0buTw5Wb5iYuSyrCx5Fq7bpxU9/7xs5hwzxnriSps6VZ6h9+ol25BffVVeXuvQyJGyWX3cOOtpCwiQqcrIuP6aMgXIyQGeeEJeJerSkCHyLNVW4gBZ5BYvlvtlaiqwe7dM3EMPVa0phxyOhU8DcnOBFSvkd+LiRdl3t3Jl2W0uXgQiI+U4jZtvBr7/HlizRq4rKADmz1cm1lGjRqFBgwZo165dybKpU6ciJCQEYWFh6NOnDy5evKhMMMnJwIYNMgG2EleseXN54Ckdm8UCfPihIqFay1uxefPmwWAw4OrVq4rEEh8vW+uKiipPW926gJeXfM2aJZcZDMDq1YqEajVvM2bMgL+/P8LCwhAWFoYtW7YoE8zRo3IwixAVJ27ZMvlFbdNG/vell+TyoiJg/XpFQrW1v33wwQdo3bo1goKCMHHiREVi0STFn/lO5Xz6qRBmsxDyGyVfb74pxKpVZZcVv4xGIcaOFSIr6/oyDw8h0tIcH+vu3bvF4cOHRVBQUMmytFIfHBUVJUaPHu34QIQQYu5cIUymqiVu61Yh+vcXIi5OiF69ri/39hYiP9/hoVrLmxBCnDt3TvTp00c0bdpUXLlyxeFxCCHEK68IUbt2xWkLCJDbGo3W98HGjYWwWBwfq7W8TZ8+XcydO9fxH36jESPKJ+TGxLVqJb+IXl7WExcaqkio1vK2Y8cO0atXL5GbmyuEEOLy5cuKxKJFvOLTgGXLZOtJVaSkyCvEDz4A3n77+nJXV0CJE99u3brBx8enzLI6deqU/JyVlQWDUiNtVqywPuz1RpGR8nJ669by64SQTVEOZi1vAPDSSy9hzpw5yuUM8motP79q2549C5w/Ly9sbrrp+vLUVOD4ccfEV5qtvClOCGDduspHAXXpIpM2c6Zs6vztN3lPSLHTpxVpXreWtyVLlmDy5MklD89t0KCBw+PQKhY+DUhMrPq29erJEfnPPSdbXorl58tb2NTy+uuvo0mTJvjiiy/whlIjbarSNGg2yzOE8eOtrxdCtcRt3LgR/v7+CA0NVfRzq/Jox6tXgU6dZF9fx46yqfOLL66vNxrV3d8WLVqEkJAQjBo1CinVHaFjj5ycqo06a9wYCA6WnaCNGskv6urVQOvWcn3t2kBSkmNjteGvv/7Czz//jPDwcHTv3h0HDx5UJQ4t0GThK+7z0MsrNja+WvnJzpa3Bn36KVC/vlyWl5ePF198tUbjmjt3bpVjeuutt3D+/HkMHz4cixYtshJzNvz8/Go0vtRr1yoPbMYM4LPP5Fm4FRkZGXhixIgajWtNcedrBbKzs/H2229XepJgsVhgNBprNL6CgsrvXcjKAg4flhc4SUny+N23b/FkCkBaWir69bu3RuM6XsVLyDFjxiA2NhbHjh1Dw4YNMWHCBKvb1WRs3p6eyKvKPR85OfIsdNYs2fe8Zw+wcyfQpw8AIDU9Hbd17FijsZ0/f75KeSssLERycjL279+PuXPn4uGHHy73CKWMjIwaja1FixaafDivJgvfhAkTIITQzSskpFm1c+TiIkfk+/vL9yZTbSxfPrdG43r11VerHdfw4cPx7bfflltuMplw6dKlGo2vbvE/viK9egEvvCAvqxMTgSZNgLVrgX879r3q1MGq776r0biGVuF+rdjYWMTFxSE0NBTNmjVDQkICOnTogEuXLpXZzsXFBUVFRTUan9lsrNofsxQhiuOR//X2rouff/6hRuOq6pWvr68vjEYjXFxc8PTTT+PAgQM2Yq652NKKiuDmUoXD5W+/2U4egLpmM/afOlWjsTVp0qRKeWvcuDEeeOABGAwGdOnSBS4uLuUGVHl5edVobLGxsfDy8qpSfErSZOHTm759ZQsIcH3CB6Ox7M+9e8v7rlxcZLPT++/L/r6TJ+XvCQHccYc68Z85c6bk540bN6J1cbOOo919t0wOYDtxvXoB7drJ5IWFydF4o0fL4eaAPDsPD1cm3lKCg4ORlJSE+Ph4xMfHo3Hjxjhy5Aj8/Pwc/tldu17/2VbaunQBAgPlCE4fHzlJws6d8hYIQKZN4RbaEoml+ga+++47qyNla5yLC1D66eq2Erdnj7xJ8rXX5Pvbbwd69gS2bZO/5+oqRxir4P7778fOnTsByGbP/Px83HzzzarEojpBqouPF8LdXQhAiOnTy6+fPl2IyEghTp4UIiNDiKQkITZtEiI4WJQMFuvUSZlYhw4dKvz8/ISrq6vw9/cXK1asEA888IAICgoSwcHBYuDAgSIhIUGZYI4evT6q01bibhxVV3pUp8EgxIABioRqLW+lBQQEKDaqMybm+ihiW2kbOlSIf/4RIjNTiIsXhVi9WghfX1EyqnjUKEVCtZq3ESNGiHbt2ong4GBx7733iosXLyoTzFdfXR+tWdH+1ratEP/9r0zeH38Icf/9crmbmxD/+Y8ioVrLW15enhg+fLgICgoS7du3Fz/99JMisWgRZ27RiF695BzK9jCb5T20kZE1G5NTCA4GTpyw73fNZmDTJqB795qNSeMsFtnia+/tlh4ewIED8kJaV/LzZad68WVvdbm7A2fOyAEwpCo2dWrEnDnFs+JXT61ach7cQYNqPCTnsGCBPBJXl5sb0L490K1bzcekcS4uQFSUfWnz8JATVeuu6AGyP2L2bPu+qCYT8PjjLHoaYZwxY8YMtYMgOfI5KEjOyFLVuTpr15az5u/dC5S6lU5fmjeXSdi+veqJc3MDbrkF2LVL/qxDbdvKArhvX9Xm9QZk0QsOlhfJrq6OjU+zOncGrl2Tg1iqmjiTSc4+/8UX10cHkar4V9CQ+++Xhc9slpNP22IwyG3atZP38un4PlTpqafkvR0eHhWfjbu4yMR27Srb6jQ42kxJr70mL5jd3Su++jMaZVr79pX3+uv0XOG6efPkkz3c3StOhqurTOzQofJswVj90bTkGOzj06DMTODzz2XzZ1KS/P6kpaXB29sbubmyP/DVV2XXlEHPjyO6UWqq7Ox87z3ZD2M0luQNeXnAwIHAK6/IIYtMXIkrV4Dly+V8r3l5MjXp6TJveXlyitOXX5aDYqmUCxeAJUvkCOF/nwySlp4u97eCAuDRR4EXX7x+8zppBgufhgkhp4VKSADuvfchbNu2DsHB8mHPVAEh5N3XiYmIvO8+fBMTI4/aeh26XUVFRfKBA0lJwKBBkdi+/Rt06iRnCqIKFBQAv/4KXL2KBwYPxvodO+TJVUXNNqQqFj4nYTAYwD9V9TFv9mHe7MO8OQf28RERka6w8BERka6w8BERka6w8BERka6w8BERka6w8DmJf/75R+0QnBLzZh/mzT7Hjh1TOwSqAt7OQEREusIrPiIi0hUWPiIi0hUWPiIi0hUWPiIi0hUWPvqfdu3aNbVDICKNYeHToOjo6JKf09LS8OSTTyIkJASPPPIILl++rGJk2jZ58mRcvXoVAHDo0CE0b94c4eHhCAgIwO7du1WOTrs6dOiAWbNmITY2Vu1QnMqhQ4fQs2dPjBgxAufPn0dERAS8vb3RuXNnHD16VO3wqAIsfBr0n//8p+TnCRMmoGHDhvjhhx/QuXNnjB49WsXItG3z5s24+d9HD7366qv4+uuv8ffffyMmJgYTJkxQOTrtSklJQWpqKnr27IkuXbpg/vz5uHjxotphad7YsWMxceJEDBgwALfffjtGjx6NtLQ0zJ49G2PHjlU7PKoAC5/GHTp0CLNmzUJAQABeeuklxMfHqx2SZhUWFqKwsBAAkJOTg86dOwMAAgMDkZeXp2ZomlavXj289957OHfuHObNm4czZ86gQ4cO6NmzJ5YtW6Z2eJpVUFCA/v37Y9iwYTAYDIiMjAQA9OrVC7m5uSpHRxVh4dOgpKQkvP/++5g3bx7S09PLPN/L8u+Tnqm8sWPH4p577sGOHTvQr18/jB8/Hrt378b06dMRxseHV8ldd92FDz/8EBcuXMCkSZOwb98+tUPSLHd3d/z4449Yt24dDAYDNmzYAADYvXs3jEajytFRRVzVDoDKe/rpp5GRkQEAePzxx3H16lXUr18fly5d4gG8As8//zyCg4OxZMkS/PXXXygsLMSZM2dw//33Y8qUKWqHp1mBgYHllhmNRvTr1w/9+vVTISLnsHTpUkycOBEuLi7Ytm0blixZgpEjR8Lf3x/Lly9XOzyqAKcs06hTp07hwoULCA8Ph9lsLlkeHR3Ng1EFmDf7MG/2OXnyJC5evMi8ORtBmrNw4UIRGBgoBg0aJAICAsSGDRtK1rVv317FyLQtKiqKebMD9zf7REVFiVatWjFvTohNnRq0bNkyHD58GGazGfHx8YiMjER8fDzGjx9fpr+Pylq+fDnzZgfub/ZZvnw5Dh06xLw5IRY+DbJYLCXNJs2aNcOuXbsQGRmJs2fP8gtVAebNPsybfZg358VRnRrk6+tb5rleZrMZmzZtwtWrV/H777+rGJm2MW/2Yd7sw7w5Lw5u0aCEhAS4urrCz8+v3Lq9e/fijjvuUCEq7WPe7MO82Yd5c14sfEREpCts6nQS2dnZaofglJg3+zBv9rly5YraIVAV8IrPSRgMBnaY24F5sw/zZh/mzTnwio+IiHSFhY+IiHSFhY+IiHSFhY+IiHSFhY+IiHSFhU/rrlwBjh9HOACcPAnwAZdVcukSICfV6ILTp4H8fLUjchIXLgBHj6IzAJw5A/z7YF+yTQjg3Dng8GEA6ITYWICPzdQ2Fj4tKioCtmwBevQAmjQBunVDNACEhwM33ww8/7w8KFEZhYXA+vUyTc2aAd27A8A2dOoE1K8PTJokD1B0g7w84MsvgdBQ4NZbgR498CMAtG8P+PoCM2bIMwkqIzsb+PhjoFUroHVr4O67ASAGISGAvz8wZw5w7ZraUZI1vI9Pa44fB+65B0hPBzIzrW9TqxZgNAL9+skDloeHsjFq0N69wP33y2P4v8/wLad2bcDFBRgyBFi+XKZR97ZuBYYOlZctthLn7i7XjxkDzJsnk6hza9YATz0FGAy2v6YeHjJt//kPMGWK3Ja0gYVPS/buBfr2BbKyqra9u7s83fzlF6DUQzD1ZvNm4OGH5Rl4VZhMQKdOQEyMLIa69fnnwDPPADk5VdveZAIiIoBvv5UnXjq1YIEsZlVNm6enPLdYvpzFTytY+LTizBmgY0fbZ922uLsDXbsC27fr8kz8yBHgrruqXvSKeXgAAwcCa9c6Ji7N27kTGDCg6kfvYiYT8MQTwKJFjolL49avB0aMqH7aPD1lU/vUqY6Ji6pHf0dKrZo0yfqVXuvWwE8/Aampsjjef3/Z9bm5wMGD8vJFBVFRUWjXrh2CgoKwYMECxT//uedk0Rs3TqYhNxdYtarsNnffLccFZWUBO3YATZvKA9fmzcUDEhxv1KhRaNCgAdq1a1eyLDk5GREREWjZsiUiIiKQkpKiTDBCXL/Ss5W4WrWAdeuAuDi5vewwvd6xFRenSKjW8rZu3ToEBQXBxcUFhw4dUiQOQHa9P/tsxWkLDwd+/FH27SUlyRMrPz+57739tnJ9ftby9uqrr6J169YICQnB4MGDkZqaqkwwGsTCpwVJSbKv5cahYEYjsHEjsGkT4OMjD1affw60bFl2u8xM4N13lYv3XydOnMDy5ctx4MABHD9+HJs2bcLff/+t2Of/9Rdw9Kj8+eJFYNYsYOXKstvcdJM8S586Vabw0CHg66/lutxc2WWlhJEjRyI6OrrMstmzZ6NXr144c+YMevXqhdmzZysTzK+/AomJ8mdbiQNkE/qIEde3LWaxAAsXOj5OWM9bu3btsH79enTr1k2RGIpt3Xp9ULWttNWrByxbJgdXBQTIBpziwmgwACtWKBOrtbxFRETgxIkT+O233xAYGIh33nlHmWC0SJD63nxTCHd3IeS59fVXUJAQGRlll23bJsQbb5Tf1t1diLg4RcNeu3atGDVqVMn7N954Q7z77ruKff64cULUqlU2DW++KcSqVdffP/20EHv3Xn9vMgmRnS1Eq1bX05acrEy8cXFxIigoqOR9YGCguHjxohBCiIsXL4rAwEBlAnnwQSFcXCpOXOnX+fNCdO9edpnZLEROjiLh3pi3Yt27dxcHDx5UJAYhhLjrrvKpqShtgBDt2wuRnn79fYMGQhQVKROvrbwJIcT69evFI488okwgGsQrPi1Yv77q9+cZDECp5osSRqNsx1NQu3bt8PPPP+PatWvIzs7Gli1bcP78ecU+/4cfgIKCircJCpIDZYtlZwOxsXI5IAe3/Pyz42KsyOXLl9GwYUMAgJ+fHy5fvqzMB8fE/P9vNHNxkR2sOiGEHHtWXd26AX/8cf19Vpbc/9S2cuVK9O/fX+0wVMPCpwW2+nZOn5bNoK++Cri6yhF13bvLAQY3ys+3/f9xkDZt2mDSpEno06cP+vXrh7CwMBgVHO2Xnl75NmYzkJZWdllaGuDlJX+2WBRPm1UGgwEGpYb8VXXUcEUMBm0kTiEZGdUfyBocDEybJr++xVxd1U/bW2+9BVdXVwwfPlzdQFSkycI3e/bskgOBHl6x8fHWE1FYKAezDBggbyCeMEH2licklNs0r6AAL73yisNifPbZZ60+nPTJJ5/E4cOHsWfPHtSrVw+BgYEl6woKChyat9TU5Er3pcxMoE6dssvq1Lk+eDYzMx0jRz5ao3GtWbOm0rgAwNfXF4n/9p8lJiaiQYMG5baxWCwwGo01Gl9BUVGV4qtIaloa+g8cWKNxHS99aV4DajI2b28zCiprXiilRQvZJzh+vOwqLZaWlobw8E41Glt1Wlk++eQTbNq0CV988QUMhvInWhkZGTUaW7NmzZBp60ZHFWmy8E2ePBlCCN28WnTpYjsZv/8uZ3C5+WZ5w3rz5sCBA+U2czObMf/TTx0W49KlS2GycqWZlJQEADh37hzWr1+PRx55pGRdrVq1HJq35s19Kt2X/vhDTkhSzGSSB6Xi5qc6depgy5bPajSuoUOHVhoXANx3331YvXo1AGD16tUYNGhQuW1cXFxQVFRUo/HVqlevSvFVpK6XF7YeOFCjcYWW/kPVgJqMzWLJhKtr1WY8aNpU3l305ptyLFppHh7eiIs7VKOxNWnSpEpxRUdHY86cOfj++++tfpcBwMvLq0Zji4+Ph1mD9xhrsvDpzmOPyRt9rAkOBtzc5I1nEyYADRsCn3xSfruCAlkYFfbggw+ibdu2uPfee7F48WLUrVtXsc9+5BF5GyMgm6Hc3OR/S//83XeyS/SBB+SyadOA336TrciAHKJePFLfkYYNG4auXbvi9OnTaNy4MT7++GNMnjwZMTExaNmyJbZv347Jkyc7PhAAiIyUbW6A7cQBsgPUza38z4D8uUMHh4dqLW/fffcdGjdujH379mHAgAHo27evw+MwGOSESsUXSbbS1qiR7GpftAj46KPy/5+GDeVoT0ezlrfnnnsOGRkZiIiIQFhYGJ599lnHB6JVgtSXliaEh4ewOixszhw57DAjQ4gtW4Ro0aL8Ni4uQjzwgNr/CsVduCCEm5tMwfTp5ddPny7X9eolxMmTcjTnzp1CBATI5bVqyZGhuvPHH9f3t4oSZ22UcECA/N033nBsjBq0Z48Qnp4Vp23aNPlzRkbZFyB/d+lSJSMmWzhzi1aMGgV89pl9s+GbTLJtpWvXmo9L4+65B9i2zb5Biu7ucsRnqW5J/ejUyf67993dgfh4OYG1jgghexpsdclXxmQCLl/W9eyCmsGmTq2YPVvebW2lw7lCJpOcqPK22xwTl8YtXlx+8EpVmEzACy/otOgB8q5qW83rFTGZ5BQkOit6gPxqfvGFfXPCm0zA0qUselrBwqcVDRoAe/bI4lfc/1IZkwno1UvXs9/ecou82PX2rvpUpSaTfEKDUhOlaFJwsJwRqDrFz2SSc3W99JLj4tK4228HvvqqesXPZAJmzgQefdRxcVH1sPBpSWCgfHpqeLj8Ztl6bo6np1z/wgvAhg1VL5T/ozp2lFORhYbKtNi638pslqmbNk1ON6nTc4XrevSQd2W3aiUTY+vMwWyWZxbvvy8fMqdzgwbJOQACAmRqbO1HXl7yPHbVKuCVV5SNkSqhdicj2XDqlBBjxsgecaNR5ABCGAxCNGsmxIcfynmQqJzjx4UYOVIOenF1FQLIEQaDEK1by6mlsrPVjlCjDhwQ4uGHhahdWwhXV5FbvL+FhQmxZo0QeXlqR6g5FosQu3cLMXCg3Nfk9Hm5wsVFiK5dhdiwQYiCArWjJGs4uMUZ5OTgJpMJ1woKdH91V1VCyFn0PT19UFiYrOfHxyHxVtwAAB50SURBVFWPEEBWFup5eSGlsFDXz92rDiGKJ0uoi6KiVD0+IcypsPA5CYPBAP6pqo95sw/zZh/mzTnwvISIiHSFhY+IiHSFhY+IiHSFhY+IiHSFhc9J+PhU/iQCKo95s48LhyXahXlzDhzV6SRycnLgYc9cSTrHvNnHYrHwIG4H5s05sPAREZGu8NSEiIh0hYWPiIh0hYWPiIh0hYWPiIh0hYWPiKiGJCUlqR0CVQELn5Pp37+/2iFoVnp6Ol577TU8+uij+PLLL8usGzt2rEpRad+lS5cwZswYjBs3DteuXcOMGTMQHByMhx9+GImJiWqHp1nJycllXteuXUOXLl2QkpKC5ORktcOjCvB2Bg06cuSI1eVCCAwcOJAHIxsefPBBtGzZErfddhtWrlyJWrVq4csvv4Sbmxs6dOhgM696169fPwwYMABZWVn48ssvMXz4cDzyyCPYsGEDtm/fjo0bN6odoia5uLggICCgzLKEhAQ0btwYBoMB//zzj0qRUWVY+DTIaDSie/fuVh9vsn//fuTk5KgQlfaFhYXh2LFjJe/feustbNmyBd9//z0iIiJY+Gxo3749jh49CgBo2rQpzp07V7LuxpzSdfPmzUNMTAzmzp2L4OBgAMAtt9yCuLg4lSOjyvCpphrUpk0bfPTRR2jZsmW5dU2aNFEhIueQl5dXZuaM119/Hf7+/ujWrRsyMzNVjk67LBZLyc+PPfZYmXVFRUVKh+M0JkyYgCFDhuCll15CkyZNMHPmTBgMBrXDoipgH58GzZgxo8zBqLQPPvhA4Wicx7333osdO3aUWTZy5EjMmzcPtWvXVikq7Rs0aFDJicGsWbNKlv/9999o1aqVWmE5hcaNG2PdunXo0aMHIiIikJ2drXZIVAVs6tSoU6dO4cKFCwgPD4fZbC5ZHh0djX79+qkYmbbZytvWrVs5MKgC3N/sUzpvRqMRsbGxaNeuHfOmdYI0JyoqSgQGBopBgwaJgIAAsWHDhpJ17du3VzEybVu4cCHzZgfmzT78njov9vFp0PLly3H48GGYzWbEx8cjMjIS8fHxGD9+vNUBLyQtW7aMebMD82Yffk+dFwufBlkslpLmpmbNmmHXrl2IjIzE2bNn+YWqAPNmH+bNPsyb8+LgFg3y9fUtM4TcbDZj06ZNuHr1Kn7//XcVI9M25s0+zJt9mDfnxcEtGpSQkABXV1f4+fmVW7d3717ccccdKkSlfcybfZg3+zBvzouFj4iIdIVNnU4iJiZG7RCcEvNmnzVr1qgdglNi3pwDr/ichMFgYIe5HZg3+zBv9mHenAOv+IiISFdY+IiISFdY+IiISFdY+IiISFdY+IiISFc4ZZmG5eQAmzcD588DwGtYvBjo2BEIDwf42K8KZGQAmzYBiYmYDABLlgBduwJhYWpHpmkpKTJtSUkAMAkffQR06wa0aaN2ZBqXlARs2QJcuYKJALB8OdCrF9C8udqRkQ28nUGDYmOBhQuBlStlgcvLA/LzC+Hu7gqjEfDzAyZOBIYPBzw91Y5WQ/78E3j/feDLLwFXVyAvDwX5+ajl4SET2by5TNxDDwHu7mpHqxlHjgDz5gHr15ekDQUFBTCZakEIICgImDQJGDQIqFVL7Wg1Qgjgv/8F3nsPiI4GjEYgPx8FBQWoZTIBFgvQubPc3/r3l+tJM1j4NGb5cmD8eKCwECgosL2dpydQpw6waxcQGKhYeNokBDB7NvDmm0B+PlDRU8PNZsDXVyaucWPFQtQii0Uel5csAXJz5XtbzGZ53vDTT8DNNysXoyYVFgJPPQWsWyebZSo6hJrNQPv2sunGy0u5GKlCLHwaMn8+MGUKUNWHOBsMsvgdOKDz4jd5MvDBB1VPnNEI+PgAR48C/v6OjU2jhACefBL4+uuqp612baBhQ3mF6OPj2Pg0y2IB7r9fngFUNXFubsCttwK//somGo3g4BaN+PFH4PXXq/5dAuTBKz0d6NEDyMpyWGja9tVX1St6gLwiTEkB7r5bnr3rUFRU9YoeIC+mExOBe+6p+CLnf9rrr1ev6AGy7Tg2FhgyxHFxUbWw8GnEa6/JVpMb7dwpl2dkyNepU2XXFxc/taYIbNasGYKDgxEWFoZOnTop++FCyKs9WwehIUNkv19mJvD338Cdd15fV1gIXLwoByWoIDU1FZGRkWjdujXatGmDffv2KfbZBQXAzJnl01a8jxW/CgtlX3Np+fnAiROAguGWOH36NMLCwkpederUwYIFC5QLIDNTnjFY298CAmRzZnKyPDv44IOy/Xq5ubJg3vgFVsj8+fMRFBSEdu3aYdiwYcjNzVUlDs1Q7mHvZMuJE0J4eAghj+RlXzt3CvHkk9bXlX7deqsQFovysQcEBIgrV64o/8FCCLFrlxBms/WE9O4tRHy8EOHhQhgMQjRqJF83bnf77aqE/thjj4nly5cLIYTIy8sTKSkpin32N9/YTlvxy9NTiIwMIe66q/w6g0GIwYMVC9eqwsJC4evrK+Lj45X70CVLZGKsJWzzZiFWrRLCzU0IX18hfvtNiOefL7uNq6sQzz6rXLz/SkhIEM2aNRPZ2dlCCCEeeughsWrVKsXj0BJe8WnAggXyTPr/IzEROHiwZuJxGu+9Z7uNd+ZM4I03ZL+KEPLq7uLF8tsdOSKboRSUlpaGPXv24MknnwQA1K5dG3Xr1lXs8+fMkRcvFXnwQTlK/+efy68TAti6FbhyxTHxVcVPP/2EFi1aICAgQLkPnTvX9v52yy3A2rWyWfPyZTnSMyio7DaFhcCnn1avmbSGFBYWIicnB4WFhcjOzkajRo0Uj0FLWPg0YN++igcivvOOPMj88gvQvbv1bSwWOVZDaQaDAX369EHHjh2xbNkyZT/88GHrnU0uLkCnTkD9+sCZM/JGyA8+sH4LQ+3aQKmnaCshLi4O9evXxxNPPIH27dvjqaeeQpaCnbQnTlS+zeOPy2O0LW5ushVZLWvWrMGwYcOU+8CiIiAuzvb6BQuAoUMBDw+gUSN5C0N0dPntjEYgPt5hYVrj7++PV155BU2bNkXDhg3h7e2NPn36KBqD1rDwaUB6uu11kybJYeT+/sCyZcAPP1i/L7agAEhLc1yMtvzyyy84cuQItm7disWLF2PPnj3KfbitYuHrKwtaZCRw113yxvX27eWQ2RsVFSmeuMLCQhw5cgRjxozB0aNH4enpidmzZyv2+db6kktr2lSeYK1eXfF2auxvAJCfn4/vv/8eDz30kHIfmp5e8U2Me/bIK7z0dODCBeDQIWDDhvLbubgonriUlBRs3LgRcXFxuHjxIrKysvD5558rGoPmqN3WeiOLxSJmzpwpAOjodarSPrzi19atQjz3nLV12QIY57AYR48eLbKysir8202fPl3MnTu35H1+fr5D83bZVpLq1pUBPPbY9WUPPCDEkSPltk0DxBAHxjhlyhRRVFRUJk+JiYkiICCg5P2ePXvEPffcU/L+2rVrDt7fcircx15/XXafVrwvpgigl8NiXLJkic39bMOGDSIiIqLc8tjYWIfF4waIQlvJMBhkf/J//iNE7dpC+PgIsWGDEO++W27bVECEOPBv+/HHH5fLy9q1a8WoUaNK3q9evVqMGTOm5H10dLTD4mnUqJHIzMy0fdBQieYKnx5FRFSt6AFCbNlSvs8cEMLLS4j165WNOzMzU6Snp5f83LVrV7F161blAggNtZ2oc+eEePTR6+8HD7Za+ISXlxA//6xczP+68847xalTp4QQ8oThlVdeUeyz/f0r3sdOnxbiiScq3sZkEuLkScVCLmPIkCFi5cqVyn+wrRFBN90k19epc33ZoEFC/P57+W3d3IS4fFnRsPfv3y/atm0rsrKyhMViEY899phYuHChojFoDQufBmzcKI+/N35HvL2F6NNHfleMRiEeeUSIzEwhWra0fvzOzVU27tjYWBESEiJCQkJE27ZtxaxZs5QNYMUK2wejmTOFOHBAiPr15RXgnj1CvPFG+e38/IS44YpMCUePHhUdO3YUwcHBYtCgQSI5OVmxz37rLdujiLt2lftYZaM+27ZVLNwyMjMzhY+Pj0hNTVX+w194QYhatawnJDZWiEmT5BfV21uehX7xRfntundXPm4hxLRp00SrVq1EUFCQGDFihMhV+mChMSx8GlBYKE8ab/yO3HyzPHanpwuRkiLEvn1ylP6N29WuLcSECWr/K1SQlSUvPawdiFxdhVi8WCYuMVGIqCh5BlF6Gw8PIUo1zepFUpIQ7u7W07Z0qRCfflpx0TObhVi9Wu1/hQr+/tt24kJD5b1HyclCXLkixNdfC9GgQfmz082b1f5XkBCCU5ZpxNtvA2+9Zd9IZw8POcKuWbMaD0v7nn8eWLFC3iBcXZ6ewLlzupx/6+GHgY0b7buNxtsbuHRJp/N89+gB7N1r34w/jRrJ/Y0TVquOozo1YuJEOQK/ugcTkwn48EOdFj1A3pTWsqUcxVkdJpOc7kaHRQ+QI4QbNZJPY6gOk0k+ukiXRQ+4vs+4VPPQ6ekpZwli0dMEFj6NcHWVMx6Fh8sruKrw8ADefRcYOdKhoWmbh4ec161Nm6ofjT085FXiwIGOjU3D6taVFy4BAfKevMoYDPLY/c03ZWd+0x0/P5k4P7+qPaPJxUXOJL9tGxAa6vj4qEpY+DTEbAa2bwemTpWPfrH2FBNXV3l879JFnnk/95zycWrOTTcB+/cDEybII7q1xNWuLRPXvTuwYweg5M3PGtWokZy4ZswYue+ZzeW3cXOTr3795OPn+vdXPk7NufVWOenByJHyJMraExc8PGTi7r9f3tN3xx2Kh0m2sY9Po4qK5LRQixcDZ88CJ0/+g5CQ5rj9dvm8vtat1Y5QowoK5I3DS5cCFy7gn9On0TwsTD6J4bnn5NRSVE5urny83IoVsv/ur7/+Qfv2zXHPPbIw6vTpTZXLzAS++EJOc3PlCv4+cwa3duoEDB4MPP20nD2INIeFz0kYDAbwT1V9zJt9mDf7MG/OgU2dRESkKyx8RESkKyx8RESkKyx8RESkKyx8RESkKyx8TiIyMlLtEJxS79691Q7BKXF/sw/3N+fA2xmcRFpaGry9vdUOw+mkpqaibt26aofhdLi/2Yf7m3Ng4SMiIl1hUycREekKCx8REekKCx8REekKCx8REekKCx8RkR3S0tIwefJktG7dGj4+PrjpppvQpk0bTJ48GampqWqHRxVg4dOgDh06YNasWYiNjVU7FKdy6NAh9OzZEyNGjMD58+cREREBb29vdO7cGUePHlU7PM3KzMzEtGnTEBQUBG9vb9SvXx+33XYbPvnkE7VD07SHH34Y9erVw65du5CcnIxr165h586dqFevHh5++GG1w6MK8HYGDbrlllvw4IMPYu3atfDz88OwYcMwZMgQNGrUSO3QNK1Lly6YOXMmUlNTMXHiRMyfPx+RkZH46aefMGXKFOzbt0/tEDVp0KBBGDx4MHr37o21a9ciKysLQ4cOxaxZs+Dv74+3335b7RA1qVWrVjh9+nS115H6WPg0qEOHDjhy5AgA4Oeff8ZXX32F9evXo02bNhg2bBieeeYZlSPUpvbt25dc2TVt2hTnzp2zuo7KCg0NxfHjx0ved+7cGQcPHoTFYkHbtm1x6tQpFaPTrj59+qB37954/PHH4evrCwC4fPkyPvnkE8TExGD79u0qR0i2sKlT4+666y58+OGHuHDhAiZNmsSrlgq4u7vjxx9/xLp162AwGLBhwwYAwO7du2E0GlWOTrs8PT3xyy+/AAC+//57+Pj4AABcXFz4UNUKfP3117h27Rq6d++OevXqwcfHBz169EBycjLWrl2rdnhUEUGaM2TIELVDcErHjh0Tffr0Ef369RMnT54UL7zwgvD29hZt27YVe/fuVTs8zTp+/Ljo3LmzqFu3rrjjjjvE6dOnhRBCJCUliaioKJWj07aTJ0+KmJgYkZGRUWb51q1bVYqIqoJNnU5m1apVeOKJJ9QOw+kwb/Zh3mxbuHAhFi9ejDZt2uDYsWOIiorCoEGDAJTtriDtYeFzMjf2XVHVMG/2Yd5sCw4Oxr59+2A2mxEfH4/IyEg8+uijGD9+PPuUNc5V7QCovJCQEKvLhRC4fPmywtE4D+bNPsybfSwWC8xmMwCgWbNm2LVrFyIjI3H27Fn2jWocC58GXb58Gdu2bUO9evXKLBdC4Pbbb1cpKu1j3uzDvNnH19cXx44dQ1hYGADAbDZj06ZNGDVqFH7//XeVo6OKsPBp0MCBA5GZmVnyhSqtR48eygfkJJg3+zBv9vn000/h6lr2EOrq6opPP/0Uo0ePVikqqgr28RERka7wPj4nkZOTo3YITik9PV3tEJwS9zf7cH9zDix8TuKVV15ROwSnNGbMGLVDcErc3+zD/c05sKnTSRgMBo4UswPzZh/mzT7Mm3PgFR8REekKCx8REekKCx8REekKCx8REekKCx8REekKC59WCQHs3w888QTQrRsOAEBEBPDGG8ClS2pHp1lCADt3AiNGAHfdBQC/om9fYO5c4No1taPTMIsF2LoVGDIEuPNO7AeA/v2BRYsA3ptmU2EhsH498OCDgJzdbR8GDgQ+/hjIzlY7OrKFtzNojRDAJ58Ab78NJCYCOTnyoFTM3V1uExEBzJwJdOigWqhaUlgILFkCzJkDpKYCWVkyTcVMJpnG++4DZswA2rRRLVRtyc0FoqKA99+XR+rMzLLrPT1l4oYMAaZNA265RZ04NSYjQ55MLVok972MjLLrzWaZtpEjgddfBxo1UiVMsoGFT0vy8oBhw4Bt2yo/XTQYAA8PYOVKeVDSscxM4N57gQMHKk+bi4tM27ffAn37KhOfZiUnyxOokyflCVZFjEZZBLduLb600a2LF4EePYDz5+V5Q0Vq1QLq1AF27ABsPASDVMDCpxUWC/DAA8CPP1Z+ECrNZAK+/BL49wGYelNQANx9N3DoUOUHodJMJnkM79bNcbFpWnY2EB4O/PUXkJ9f9d/z9AR++QWwMqG1HqSkAO3bAxcuyCu9qvL2lvvorbc6LjaqOvbxacWiRUBMTPmiN24ccPCgPKqvWlX+97KzgUcekaehOjR9OnDkSNmiV7s2sGIFEB8vu6eOHgX69Sv7e9nZ8ipRt91Xzz0H/P13+aL32WdyX0pLA06fBp58suz6rCygT5/qFcv/IY8+KnsgbBW9W2+VX+HPPiu7PCNDtjDwMkMbWPi0wGIBZs+23k538SIwa5Zs0qzo95cscVx8NowaNQoNGjRAu3btFP9sQLYML1pUPm2urrIZqnt3eaY9ZQqwdi0QEFB2u6Ki8gcoJZw/fx49e/ZE27ZtERQUhKioKGUDSE0FvvrK+iXyO+8AzZrJxN13n9z3buxHzs0FNmxQJNSyH5uLLl26IDQ0FEFBQZg+fbqin3/+PPDTTxXX/MWL5XnqjSwWIClJDrxSQ1FREdq3b4+BAweqE4DGsPBpwfbt5XvHi333HbBxY8VDEnNzZQUoKHBMfDaMHDkS0dHRin5mad98Y315drYc93P2rDzD3rwZiIsDOnYsu11WlhygoPRZuKurK+bNm4c///wT+/fvx+LFi/Hnn38qF8Ann8jOTmv+/PP6kV0I+WrRouw2GRnAu+86NERr3NzcsGPHDhw/fhzHjh1DdHQ09u/fr9jnL15c8b4yZIg8p/jpJ+vrMzPl/qaGqKgotOGIrhIsfFqwcGH50XTVVVQk+wcV1K1bN/j4+Cj6maUtWGD7fKG0Bg2AwEDgjz/Kr7t2TQ6KUVLDhg3R4d+rKC8vL7Rp0wYXLlxQLoCFCyseBbR4sTwrOH1atutt2VJ+m5MngdhYx8VohcFggNlsBgAUFBSgoKAABoNBsc//6CPZymCNl5e80+jllyv+f+zcqfxtNQkJCdi8eTOeeuopZT9Yw1j4tODvv////4+CAtmppSPnzlW+jasr8MUXwOrV8jh+I4NB3bTFx8fj6NGjCA8PV+5DK+sPHjdOHsnvvFPepGbtaF+7tiqJKyoqQlhYGBo0aICIiAjF8pafX3F/8Jtvynv3Kjt/cXOrfJua9uKLL2LOnDlwsXWVr0OazMT8+fNhMBh08zpn7YhcTQW5uZj43HMOi/HZZ59FdjXvyC0+I3fUKykprcLPNxhkH15+vhzLYU1GRiaGDh3lsBinTp0KS+n7MEvJzMzEgw8+iAULFqBOnToly5OTkx2aN4uty5bSLBZg716gcWPAyjPmUtPScF/v3g6LcebMmVbDMhqNOHbsGBISEnDgwAGcOHGiZN3Ro0cdFo+bWz1YLNbzFhoK9O4NzJ9feVrT01MRGtrVYXGuvGEswKZNm9CgQQN0vLGdv5Rt27Y5LJ6WLVsioyrNMkoTpL62bYt7U2y/3nxTiFWrbK/38BBiyRLFQ4+LixNBQUGKf64QQvj5VZyylSuF2LFDCHd329vUqSPEunXKx56fny/69Okj5s2bp/yHe3hUvr8Vv5YvF2LBgvLLvb2F2LlT+dhLmTlzppg7d64in1VQIITBYD1F48cLkZkpRGKifGVkCJGdLcThw+W39fIS4sQJRUIWQggxefJk4e/vLwICAoSvr6/w8PAQw4cPVy4AjdLkFZ/udOhge7CB0SjbR4zGsj9b265tW8fGqTEVDSZdskTOznLvvRXf31dYCAQF1XxsFRFC4Mknn0SbNm3wcmWdQo7QsqX15fXryxEanp5yf+zTR06oYG20Rl4e0KqVY+O8wZUrV5CamgoAyMnJQUxMDFq3bq3IZ7u6Av7+1tctWybH/4SFydfSpXJAlbUJEgoKyo8udqR33nkHCQkJiI+Px5o1a3D33Xfj888/Vy4AjWLh04IXX5RTkVkzZYo8cr/2mryJKDdXLrtRvXrFk1MqZtiwYejatStOnz6Nxo0b4+OPP1b08ydMkFND3ahpU+DZZ+VB6NIlOQAmI0Pe7nijVq2Un75s7969+Oyzz7Bjxw6EhYUhLCwMW6wNIHGUV1+1njghZLNmQoK8U/u99+S++cMP5bft3h1o2NDxsZaSmJiInj17IiQkBJ07d0ZERISiw/NfeklOfHCjnBzg8uXrr8xM+TW9erXsdkYj8NBD1lNPyuLMLVrRtq0cKWcPT095/9Xzz9dsTBpnscg5EC9ftu/3vbyA5ct1OONbbq4c6mpv34vZLO/j69WrZuPSuJQUub9VZ4ag0jw8gP/+V7eT3mgKr/i04q23rJ9OVoWbG/DYYzUbjxNwcZEzt9iTNhcXOYfi4ME1H5fmubvLy2V7EufqKi+p77675uPSuHr15MNS7EmbmxvQqROLnlaw8GnF4MGyWam63yqzWd4A7+3tmLg07tlnZRNmddJmMMirvV275Kh8XZo6VU5Q7eFR9d8xGuXRPyZGJlGHoqJkl7ytnglratUC/PzkPBSkDSx8WjJrFjB5sjyKV3bPjbu7PAjt3i1nzdUpg0HeWPz00zJtlR2PTSbZyvfrrzqfMNjFBVi3Tk6M7ulZ+fYmE9CkCXD4sK6fsVOrlpwnomfPqqXN01P2Ix88KL+upA0sfFpiMMgz8V275JMt3d3LXsoUX6rcfLN8yNfp03weH+QxfMECIDoaGDBANiuVPiM3GOSFccOGcnaNU6cUH5CoTbVqyRsd16+XTZfu7jJ5xVxc5JG7WTM50OXECVn8dM7DA9i0Cfj8c6BrV5m20i0HxU9wat1aziR48KAcMEvawcEtWnbtmpxd+dw5fDB7Np6fPl0WugEDrN/SQADkLFvr1skZMubMWYg33ngB4eHyJmNOXlGBs2dlEUxMxMK5c/HCm2/K0Zt33qnbps2qOHNGjvW5fBmYN28+3nrrJUREAJ07qx0Z2cLC5yQMBgP4p6o+5s0+zJt9mDfnwPNfIiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+IiLSFRY+J1FUVKR2CE6JebNPYWGh2iE4JebNORiEEELtIIiIiJTCKz4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItIVFj4iItKV/wOWKFdG5mwvRQAAAABJRU5ErkJggg==\n", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb4AAAEuCAYAAADx63eqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nOydd3gUVffHv7vZtE0DFBIISegEUulRlNASqjQRadIRxYpI8JWuiAKvNFGKSBGQJgi8NKWDkRK60gktgUBIQuqmbfb8/rim7yabPNmdyW/P53n2YXdm1v36zcyce889c6+CiAgMwzAMYyEopRbAMAzDMOaEAx/DMAxjUXDgYxiGYSwKDnwMwzCMRcGBj2EYhrEoOPAxDMMwFgUHPoZhGMai4MDHMAzDWBQc+BiGYRiLggMfwzAMY1Fw4GMYhmEsCg58DMMwjEXBgY9hGIaxKDjwMQzDMBYFBz6GYRjGouDAxzAMw1gUHPgYhmEYi4IDH8MwDGNRcOBjGIZhLAoOfAzDMIxFwYGPYRiGsSg48DEMwzAWBQc+hmEYxqLgwMcwDMNYFBz4GIZhGIuCAx/DMAxjUXDgYxiGYSwKDnwMwzCMRcGBr5Kg0WikllApYd/KB/tWPti3ygEHvkqCVquVWkKlhH0rH+xb+WDfKgcKIiKpRTAMwzCMueAeH8MwDGNRcOBjGIZhLAoOfAzDMIxFwYGPYRiGsSg48MmQK1eu5L3Pzs7G7Nmz0atXL3z++edcLl0CS5cuRVxcHADgzp07aNeuHapUqYI2bdrg77//llidfOnXrx82bNiA1NRUqaVUKu7evYtRo0Zh6tSpSE1NxdixY+Hr64s33ngD9+/fl1oeUwIc+GTIiBEj8t5/9tlnuHPnDiZOnIj09HS888470gmTOcuWLcOLL74IAPjoo48wYcIEJCYmYu7cuexbCZw5cwY7d+6Ep6cnBgwYgN9++w1ZWVlSy5I9I0aMQKtWreDo6IigoCB4e3tj//796Nq1K0aNGiW1PKYE+HEGGdKsWTNcvHgRABAYGIiIiAhYW1uDiBAQEFCoR8jk07hxY9y8eRMA0KpVK0REROTt8/f3Z98MkHu+JScnY9euXdi0aRMiIiLQs2dPDBo0CKGhoVJLlCUFr1NPT088fPhQ7z5GfnCPT4YkJSXht99+w/bt25GZmQlra2sAgEKhgEKhkFidfOnfvz9GjBiBu3fvom/fvli0aBEePHiANWvWwNPTU2p5siX3nHJ2dsZbb72Fffv24caNG2jTpg2++eYbidXJF6VSiVu3biEiIgIajQbnzp0DINLsOTk5EqtjSoJ7fDJkxIgRhQLcN998A1dXVzx58gRDhgzB4cOHJVQnb9asWYPly5cjMjISmZmZ8PDwQJ8+fTB58mS4uLhILU+WtGvXDidOnJBaRqXj8OHDGD9+PJRKJX788UcsXLgQV65cQVJSElauXIk+ffpILZExAAe+SsKwYcPw888/Sy2j0vHWW29h/fr1UsuodPD5Vj569uyJ3bt3Q6nkZJqcUUktgClOr169im07cuQIEhMTAQC7d+82t6RKgSHfcrezb/op6hsR4ejRo3y+lYK+8+3YsWN5PT32Tb5w4JMhUVFR8PHxwZgxY6BQKEBEiIiIwMSJE6WWJmuio6PRtGlT9q2M6Dvfzp07x76VAl+nlRhiZEdOTg4tWLCAOnfuTBcvXiQiorp160qsSv6wb+WDfSsf7Fvlhcf4ZEx0dDQmTJgAV1dX7N69u1C5NGMY9q18sG/lg32rfHCqU8bUrl0b27Ztw969e+Hs7Cy1nEoD+1Y+2Lfywb5VPrjHxzAMw1gUXHNbCcjKAh49SoFOJ7WSSkZmJlIePQK4bVcmMjKAx49T2Laykp6OxOhoPt8qARz4ZEpkJPDxx4CLC2BvD3h62sPaGmjUCFi1CkhLk1qhTLl2DXj7bcDREVCrYe/hAahUgL8/sHEjkJkptUJZcuECMHSoONccHIDate2hUgGtWwPbtwPZ2VIrlCFEwF9/Af36Aba2gJMTHD09AWtrIDgY2LsX4BlcZAmnOmXGkyfAoEHA6dPimtF3w3FwENfcxInAzJkAPysL4N494M03gX/+EaZptcWPcXQEFAph2oQJ4r2F888/4ny7e1e0CfTdpx0dRdthwQJg5Ejza5QlZ84AQ4aIC1aj0d/Lc3QULYnly0VwZGQDBz4ZcecO0LYtkJCg/75dFAcHoEsXYOtWwMrK9Ppky6VLQIcOQHIyjMoHq9XiprVihUUHvz//BLp1E9kDY+4CarXIQnz1lem1yZq9e4EBA0TAMwZ7e2DOHGEeIwusZs6cOVNqEQzw7BnQpg0QG2vcvRsQHZsHD8RLzyQSlsGDB8DLLwOJicaPrWRnAzdvAunpQMeOptUnU65dA9q3B8qyBF92NnD+PODkJM5Vi+TMGaBnT+ODHiBasSdPAvXqAX5+ptPGGA0nyWTCf/4DxMWJ4YFVq4D790UH5uJFoGtXcYy1NbBtm8jqEYlhBI0G+OUXcT2ag6ioKHTo0AFNmzaFj48PFi9eDADYtm0bfHx8oFQq82apNwvvvSeqMQyZ1qYN8McfQHy8aFVs3Qq4uYluzrffim62GTDk26RJk+Dt7Q1/f3/07ds3b5owUzNypAhkhmxr0gSIiBDZh4QE4OBBsU2jASZPFlaaA0O+TZs2Df7+/ggMDERoaCgeP35sejFEwODBIpAZMq4g06aJ73TqJIwbO9Ysg/OGPMvl22+/hUKhyFu02SKR7tl5JpfkZCJ7eyKASK0mmjGDyMuLSKEg6tFD7PfyIrK2JvroI6K2bYkePyYKDhbfUSqJ+vc3j9bHjx/T+fPn/9WdTA0bNqSrV6/StWvX6MaNGxQcHEwRERHmEkNka1uyaV27CnOcnITJP/1EtH+/MM7amuiDD8wkVb9vv//+O2VnZxMRUVhYGIWFhZlcy7VrwoqSbHNxEf/mnl8ffEB0+bL4bG9PNHu2yWUSkWHfkpKS8o5ZvHgxjRs3zvRi/vyTyMGhZONEqCOqV4/oyhWiR4+IOnUS2xwdiX780eQyDXlGRPTw4UMKDQ0lT09Pevbsmcm1yBUOfDLghx/E9ZR7zRR9Xb5M1K9f4W1RUfmBDxD3fynO4169etEff/yR99msgW/GDCI7O+NNA4iaNRM3qdzPjo5EGo159BagqG9ERDt27KDBgweb/LfHjSNSqYy3zcqKaPx4orS0/G0vvkik1ZpcajH0+TZnzhx65513TP/jffqIQGeMcfv3E3XrRnTvXn7gA4jq1yfS6UyvtQAFPXv99dfp0qVL5OXlZdGBj1OdMmDdOsMZkBo1xCMMV6+W/N+wtgb27694bSVx//59XLx4EW2kGvDZsEGkOYtSkmnt2hXerlQCx46ZTKI+DPm2evVqdOvWzeS/v22b/uIpfbY9fy4s/u47UZ+RS2ameATCnBT1bcqUKfDw8MDGjRvxxRdfmPbHiYA9e/SPIxc1rn9/YZC+C/LxYzEubSYKerZr1y64u7sjICDAbL8vVzjwyQBD4yUqlXj0bN06UYtREllZYozQXKSmpuL111/HokWLpJumKSGh+LaSTPPzA6ZPByZNyt+m05nVOEO+ffXVV1CpVBgyZIjJNSQnF99myLaqVcWzpO+/L4ayclEqpT/fvvrqK0RFRWHIkCFYunSpaQUYapkWNc7RUbQQPvpI//HW1mYzrqBnKpUKc+bMMX0DobIgdZdTH/PnzycAFvS6WyxzolAQbdpEtHev/rRU0VQnkEnAxArVtXbtWr1/n6ysLAoNDaVvv/222D5Dqc6kpCRq3LhxhepLKItp9esTRUcTDR1aaHsyQMMr+O9ZNBVXmm9r1qyhoKAgSktLK/YdjUZD1apVq+DzLbtM51ruMXFxRNWr5257TkCPCtV19OjRMp9vREQPHjwgHx+fQtvS0tIqVJsTQJnGnG///S/RtGn5xxRJdT4HqHUFn2///PNPqZ5duXKFqlevTl5eXuTl5UVWVlbk4eFBMTExhb4XExNTodrq1q1LycnJev9uUiLLwGdpBAYWv9GsXk105IjhIayigU+tJlq1yvRadTodvfXWW/TRRx/p3W/WMT4PD+NM8/QUN6Bx44ob6exMtGuXyaUa8m3//v3UpEkTio2NNbmGXJycynau4d9xPo0m/1x1diYKDze9VkO+3bp1K+/9kiVL6PXXXze1EFHlU5pxFy+KwfaYGPHSaoni44nCwsR+BweimzdNLLXka5SILH6MjwOfDJg1q/C1s2wZ0alT+gtebGxEIUtUFFFIiHgPiO/fvWt6rSdPniQA5OfnRwEBARQQEEB79+6lHTt2kLu7O9nY2FCNGjUoNDTU9GLeeUdUZpZkWq1aRHfuEE2cqP+ObmcnujImxpBv9evXp9q1a+dtM0d1Yt+++TUahmzr3FkEOaVSBMrFi0WBYu75Zq6aIEO+9evXj3x8fMjPz4969uxJ0dHRphfz8sulX6TVqhG5uua/Hj4UVcW5x1WvbvKqIEOeFYQDHyM5MTH5gc/TU2xLTydKScl/DR6cnzkpipeX6P1ZHDduiNr6kkybPl3sK7g9JUWYqVQSDRggnX6JCA8X9+GSbOvfn+j6dfE5NpZozx4iPz9hm7U10YcfSvv/IAm7dolWQGkXacFXwVSnvT3RV19J+r/ACHjKMpnQuzfwv/9Bb9FYaTg6Alu2AN27V7wu2RMUVP6n99VqUdHZqlWFSpI7RECDBmJ+zvJgZyfm+Kxfv2J1yZ6cHDH5QXmLU2xtgagooHr1itXFlBmu6pQJ8+aJuTfLip0d0KyZ/okjLILvvhMBrKzY24uJTi0s6AFietIVK4QFZUWtBoYNs8CgB4gJcb//vvzGhYVx0JMJHPhkQuPGYu7bstzD7ezEDWjvXgteoaFVK2DTprLdjOztRWth0ybT6ZI5nTsDP/xQNtvUajEX+A8/mE6X7BkwAPjii7JdqGo18MYbwKxZptPFlAlLvV3KknbtxFy2rq5iImBDqFTihtWpk8jylXSsRdCrF3DggHjozNHR8HE2NqK10K8fcPSoSD1ZMCNGiBS5k1PJttnaCttGjQJ27bLwlUAA4NNPgWXLREArKU1jby+M+/RTYM0ai14JRG7wGJ8MyckRvbi5c8XsGDY2QGpqMhwdnaHVAm+9JZ6PbdJEaqUyIysL2LlTGHftGmBri+SUFDjnLmA4dizwwQdA3bpSK5UV6eliNpe5c8UE6DY2QHJyMhwcnKFSiXnA330XcHeXWqnMSE0VswfNmwc8fQpYWyM5OVmcb3Z2YhmisWPFzC6MrODAJ3NiYsQsR8+epcLd3RH165dvSMviiI4GnjxBanw8HGvXFjlhOzupVcmeBw/ETEKJiRq4u6vRsKGYbIQpASJRKRQXB01qKtTu7kDDhtw1ljEc+BiGYRiLgsf4GIZhGIuCAx/DMAxjUXDgYxiGYSwKDnwMwzCMRcGBr5Lx9ttvSy1BtuTk5GDFihWYNm0awsPDC+2bPXu2RKrkj0ajwbx58zB//nxkZGRg7dq16NWrF8LCwpCamiq1vEpFo0aNpJbAGAFXdcqQBH0LrAIgIgQEBCA6OtrMiioHY8aMgUajQevWrbF+/XoEBwdjwYIFAIDmzZvjgrmXDK8kDBgwAB4eHkhPT8fNmzfRpEkTvPnmm9i9ezeePHmC9evXSy1Rljg5OUHx70PpubdRjUYDtVoNhUKBZH0r/jKygAOfDLGysoKXlxcK/mkUCgWICI8ePUJWVpaE6uSLv78/rly5AgDQarUYP3484uLisGnTJgQFBeFiwSXEmTwCAwNx6dIlEBFq1qyJmJiYvPMtICAgz1OmMB9++CESExMxf/58uLq6AgDq1q2Le/fuSayMKQ2V1AKY4tSrVw+HDx+Gp6dnsX0eHh4SKKocFGwQqFQqrFy5ErNmzULHjh05ZWcECoUC3bt3z+vFKBSKvPdMcZYsWYLz589j0KBB6NOnD95//332q5LAY3wy5OOPP8bz58/17gsLCzOzmspDy5YtceDAgULbZsyYgZEjR+L+/fvSiKoEtGzZMq9hsHr16rztkZGRcLL4iWBLpkWLFjh06BAAIDg4GBkZGRIrYoyBU52VhGHDhuHnn3+WWkalg30rH7m+ERH3YowkJiYGvr6+iI+Pl1oKUwqc6pQhvXr1KvSZiHD06FEkJiYCAHbv3i2FLNnDvpUP9q18FPUNADIzM/O2s2/yhQOfDImKioKPjw/GjBmTV2Rw7tw5TJw4UWppsoZ9Kx/sW/mIjo5G06ZNC/kWERHBvlUGiJEdOTk5tGDBAurcuTNdvHiRiIjq1q0rsSr5w76VD/atfLBvlRce45Mx0dHRmDBhAlxdXbF79248fPhQakmVAvatfLBv5YN9q3xwqlPG1K5dG9u2bcPevXvh7OwstZxKA/tWPti38sG+VT64x8cwDMNYFPwcXyUhMzNTagmVEvatfPAD/+WDfasccOCTK3FxwNy5gJcXoFZDZW8PODkBbdsCu3YBOTlSK5QlMTHA9OmAuztgbw/Y26vg4gJ07gz88Qeg00mtUKY8eABMmgS4uQF2drBzdgaqVAFeew04cQLgxJBebt0C3n8fqFEDsLUFXFzsUK0a8OabQESE1OoYQ3CqU24kJgLvvCOCm0IBpKcXP8bJCbCxAb78Enj3XfNrlCFPnwJjxgAHDwrb9E2g4egorJs/HxgyxPwaZcn9+8Do0cBff4lWQdF5YBUKQK0GqlcHvvsO6NlTEply4/p1YNQo4PJlQKsFsrML71cqATs7wMMDWL4caN9eEpmMATjwyYnHj0WP7vHj4jcgfajVwPDhwPffixuUhXLnDvDKK0BCQvEbkD7UauDTT4FZs0yvTdZcvgx06AAkJxuXQbC3F1mIDz4wvTYZEx4OdOsGpKYa1xG2twd+/JEbW3LCaubMmTOlFsFA3Hxeegl4+NC4uzcgjrt2DdBogE6dTKtPpsTGAq1bi3+Nzf5mZ4s0lKMjEBRkWn2y5f594OWXgefPjU9jarXA8eNA3bqAn59J5cmVq1dF760sQ3laLbB/P9CiBdCggcmkMWWAx/jkwvTpQFQUMG6cuCtnZABr1uTvHzwYSEnJf6WliRuWtzewcKFovVsgH3wghkPHj9dvGwC88YZoHyQnixtX796irfDZZ8Jyi2TkSCApCXjvPcPGjR4N3L4tzrf9+4GaNYVxY8aIlLwFMnCgaDitWiXaDsnJwMWLQNeu+cd07ChSoWlpwJEjgKenGLF4802Aa63kAQc+OZCeDvz0k7gqHj8GZs8GCsySDwD45RcxQJX7Gj8eiIwELlwQadGFC80iddSoUahRowZ8fX3ztr355psIDAxEYGAg6tSpg8DAQLNoiY8Hdu8WLWpDttWqBWzYAHzyCeDsLOo3fvlFDFkRiSyxOdDn26VLlxAUFITAwEC0bNkSZ8+eNY+Ye/eA06fFmJ4h44KDgTlzRCuhWjXxnU2bxD6FAli71ixS9fl2+fJlvPTSS/Dz88Nrr71mtgVfz58H7t4FrKxEgyk4GHBxAaZOBbZuFXVoL7wA7NgBTJsmbDt3DtiyRXxfpwN+/dX0OqOiotChQwc0bdoUPj4+WLx4MQCxwHVISAgaNmyIkJAQgyvAWATSTRrD5LF2LZGjI5G4F4vXl18SrVlTeFvB15EjRNOn53+2syNKTDS51OPHj9P58+fJx8dH7/5PPvmEZs2aZXIdRERz5xLZ25dsW+vWRE+fFj4mNpYoKEi8d3Ehysw0vVZ9voWEhNC+ffuIiGjv3r0UHBxseiFERJ98QmRjU7Jx8+cTLV2a/7lmTfHdevXE51q1iHQ6k0vV51vLli3p2LFjRET0008/0dSpU02ug4ho8GAiKyv9l+Ply0T9+hGNHUsUHp6/Xa0m0miIGjcWn/38TK/z8ePHdP78eSIiSk5OpoYNG9LVq1dp0qRJ9PXXXxMR0ddff01hYWGmFyNTuMcnB378sWyDBp6eQLt2QMHldqytgX37Kl5bEdq1a4dq1arp3UdE2Lp1KwYNGmRyHYDopOgrei3IuXMi7fTaa6LSrndv0bHOXVScSAxbmRp9vikUirzeSlJSEmrVqmV6IQCwfr1xxVMFC6Zy3+f2vJKTgUuXKl5bEfT5duvWLbRr1w4AEBISgu3bt5tcBxGwfbv+ceQaNYBGjUQa3cen8KiDRiMSMz4+4vPt20B0tGm11qxZE82bNwcAODk5oUmTJnj06BF27dqF4cOHAwCGDx+OnTt3mlaIjOEpy+TAkydlO37YMODkSTHIkEtmpqjwkJCTJ0/C1dUVDRs2NMvvxcWVfoxOJ9oHv/wiysuzssSYn0aTv18q2xYtWoQuXbrg008/hU6nw19//WWeHzZmfO7AAWDzZlGLf/u2GIPW6URJLCDyfRIZ5+Pjg127dqFPnz7Ytm0boswwUJuerj/oqVTAxo3AunXAzZuiYOrZs8LHJCWJ0QlAPIUUGwvUrm1yyQCA+/fv4+LFi2jTpg2ePn2KmjVrAgDc3Nzw9OlT84iQIbLs8X377bdQKBQW87ofGVk2g4YNE1daAbKysjDp448rVNe6Ir9RGps2bTLY20tOToa3t3eF6ouPL/0G3qkTMG+eqMSzsRHjMqtWAQEBYn9qagqGDh1Zobo2b95slF/Lli3DwoULERUVhYULF2L06NHFjtHpdLCysqpQfTnGVA0fPgzMmCG6Offvi1dKSl53JTEpCa917Vqhui4bWaC1evVq/PDDD2jRogVSUlJgY2Oj97iK1Obg4AKtNrPIfz+/8/z++2JbaqoYSy6Is7OwDgCSkxPRokVQhWq7evWq3v//1NRUvP7661i0aFGxOURzv1uUJ0+eVKi2+vXrIyX3f15OSJ1rZYjI37/4oIGhMb6XXyZKTS0+JqhWE/34o1nk3rt3r9gYX3Z2NtWoUYOioqLMooGIyN29dNsmTiTasaPwMb/9JrYDRM7O4rM5KOqbs7Mz6f4dJ9PpdOTk5GQeIUXPHWPGlBs2FOddlSr5g6MnT5pFrr7zLZebN29Sq1atTK4hJ4dIqSxsyerVYqjdzi5/29ixRH/+WfiyTEvLH+NzdCS6ccPkcikrK4tCQ0Pp22+/zdvWqFEjevz4MRGJccBGjRqZXohMkWWPz+Lo2lV0RwCRQrK1Ff8WfJ/L8OGiFV50TJBIPMUtEYcOHYK3tzdqmyuHA9Gby7XGkG0REcCrr+b38AIDxefcMb7MTOme5atVqxaO/zvAeOTIEbOliNG2bf57Q8bZ2uYPTHl4ACtXAosX56dJs7LyTTUzsf+mWHU6HWbPno133nnH5L+pVBb+3122DGjSRIwdF5wl6LffxDBov37CwunTxbl286bYr1IB9eqZVisRYfTo0WjSpAk++eSTvO29evXKy+KsW7cOvXv3Nq0QOSN15GWI6MGD/GbjjBnF98+YIfbZ2hI9f07UsWPxFrkZWr1ERAMHDiQ3NzdSqVTk7u5Oq1atIiKi4cOH07Jly8yiIZeLF0WLujTb3nuP6PZtouRkoshIUdQIECkURD17mkerPt9OnjxJzZs3J39/f2rdujWdO3fOPGIOHcrv9RkyzsVFlCqmphLFxBDNmZPf5VGpiEaPNotUfb4tWrSIGjZsSA0bNqTJkyfn9ZpNzebNRE5ORJ6e4nN6OlFKSv5r8GBhT6dORNevi2rOo0eJvLzyL98pU0yv8+TJkwSA/Pz8KCAggAICAmjv3r0UFxdHHTt2pAYNGlCnTp0oPj7e9GJkCk9ZJhc6dxbjKuXB0VE8V/X66xUqqTLg7w/8/Xf5vuvoCOzdKwpkLQoi0Yt79Kh837e3F13p3B6hhZCVJSo4k5LK931bW1Hh6e5esbqYssOpTrkwf35+xVxZsLERuZNevSpeUyVgyRJxHy4rdnZiCqlXX614TbJHoSi/cfb2YqJqCwt6gLjU5s0r32WqVosJbzjoyQOeq1Mu1Kwp5j/ctUtMRWIMNjZiGZnw8Px6aQujTh0xO8sffxhvm52daCscPSpa4RZJkybi/AkPN35uWHt7MUi6a5cYrLJAWrQQw5yXLhlvm1otqop//lmMFTLSw38GOdGrF7Bnjwhijo6Gj1MoAAcHMdp+8SLw4ovm0yhDRo4Uz+mp1cIWQyiVYv8rrwBnzpRssUUwaZJYasjOruTen5WV2N+jh4W3FgTz54uVPezsxMsQ1tZi/9ChYmq9gjVqjLTwGJ8cSUsT8yLOnStWVlWpkJycLJ7FycwEQkLETeuVVyx6OaKiJCWJVvX8+WLRASsrMSOKi4sLMjOBPn3EnJ2tWkmtVGbEx4u5YhcsEE/2K5X551tWlpiZecIEi12RwRAxMcCKFSJrrNWKSzE5WZxvWi0wYgTw4YdiVhdGXnDgkzNEYh6k6GhkPH8OO1dXMbbi6iq1MllDJFJRT54AaWlavPiiCv7+YtJgpgR0OjHpeWwsMpKTYVezJtCsWfEnsplCaLViary4OCArS4vq1VVo0aJ8Y4GMeeDAxzAMw1gUPMbHMAzDWBQc+BiGYRiLggMfwzAMY1Fw4GMYhmEsCg58MqRfv37YsGEDUsuyOC2Du3fvYtSoUZg6dSpSU1MxduxY+Pr64o033sD9gmsXMoXQ6XRYvXo1evTogYCAADRv3hwDBw7EsWPHpJYma7RaLVasWIGuXbvC398f/v7+6NatG5YvX45sY59uZySBqzpliLu7O1566SUcOXIEnTt3xqBBg9CjRw+D644xgnbt2mHQoEFISkrChg0bMHLkSAwYMAB//PEHNm7ciCNHjkgtUZaMHDkSXl5e6Ny5M3799Vc4Ozvj1Vdfxdy5c9G7d2988MEHUkuUJYMGDUKVKlUwfPjwvFVJoqOjsW7dOiQkJGDLli0SK2QMwYFPhjRr1gwXL15EcnIydu3ahU2bNiEiIgI9e/bEoEGDEBoaKrVEWZLrGwB4enri4cOHevcxhfH398eV3HWaAAQFBeH06dPIzMxEYGAgrl+/LqE6+dKoUSPcunWrzPsY6eFUpwzJXRnZ2dkZb731Fvbt24cbN26gTZs2+OabbyRWJ1+USiVu3bqFs2fPQqPR4Ny5cwCAO3fuICcnR2J18sXa2roQWjkAACAASURBVBqRkZEAgAsXLuRlFmxtbfWu0s0IqlWrhm3btkGn0+Vt0+l02LJlC6pWrSqhMqY0LHOmWZnjqGcSyRdeeAHvvPOOWRbdrKzMmzcPr732GpRKJXbu3Imvv/4aV65cQVJSEn788Uep5cmW+fPno0OHDrC1tYVWq8XmzZsBAM+ePUPPnj0lVidfNm/ejMmTJ2P8+PGoWrUqiAhJSUno0KFDnoeMPOFUp0w5e/YsFAoFWrVqhWvXruHAgQPw9vZG9+7dpZYma86cOQOlUolWrVrh6tWr2L9/P5o2bcq+lcKpU6egUqn4fCsn8fHxAICPPvoIGzZskFgNUxoc+GTIrFmzsH//fmi1WoSEhODMmTPo0KEDDh48iC5dumDKlClSS5QlRX07e/Ys2rdvz76VAvtWPnrpWQPzyJEj6NixIwBg9+7d5pbEGAkHPhni5+eHS5cuITMzE25uboiOjoazszPS09PRpk2bQoUITD7sW/lg38pH8+bN0bRpU4wZMwYKhQJEhEGDBuWlOYODgyVWyBiCi1tkiEqlgpWVFdRqNerXry+WhwFgb28PJa9kaRD2rXywb+Xj3LlzaNGiBb766iu4uLigffv2sLe3R3BwMAc9mcPFLTLExsYGGo0GarUa58+fz9uelJTEN6ISYN/KB/tWPpRKJSZMmIA33ngDEyZMgKurK7RardSyGCPgVKcMyczMhK2eVa7j4uIQExMDP14QVC/sW/lg3yqGvXv3Ijw8HHPmzJFaClMKHPgYhmEYi4LzGMz/a3jORMac8PlWOeDAJ1PS04G1a4GWLYGaNQG1Oh2enkD//sCpUwD30w2QkgIsWwYEBABubsh2dgbq1AHeegvgKcsM8vw5sGAB4OMDuLqK861uXWDcOIBnLCuB2Fhg9mygcWOgRg1kOTsD9eoBH38M3L0rtTrGAJzqlBnp6cDnnwM//ggoFEDRBRoUCkCtBtzcgLlzgddfl0an7EhOBiZOBDZuBJRKIC2t8H4rK8DWFqhbV9zheb5TAMCzZ+IevWOHsE2jKbxfpQKsrUVAXLQIaNtWGp2yIyoK+PBDYP9+cVFmZBTeb20tzrlWrYAlS4DAQGl0MnrhwCcjEhKAjh2BmzeLX0f6UKuBsDBgxgzTa5M1jx8D7dqJm1FWVunH29sD8+cD771nem0yJjISePVVIC4OMCZDZ28P/PQTMGiQ6bXJmitXgA4dgKQkwJg5YNVq4NdfgW7dTK+NMQoOfDIhIwN4+WXg6lXj7t25qNUi0zJhgum0yZqkJKB5c+DhQ6AspeRqNbBiBTB0qOm0yZinT0Un5OnTsqXN7e2B7dst+B5+/744354/L9v31Grg4EFxkTOSw2N8MuHLL8VYytixQESECIRr1uTvb9JEbE9IEK+DB8U2jQaYMkX0Ei2STz4BHj0Sg1H6jCvItGniLt+pkzDu7bfFnd8CGT0aiI8Hxo/Xb5uXl7AqJSX/NXWqSMW/8UbxFLzF8OaborH13nuGzzd7e+D770UeOTEROH5cnG+9ehnXtWZMDgc+GZCVJa6TjAyRtZs9G1i9uvAxjx+LwpZq1YAXXwR27wZyJ4DXaoHFi82jddSoUahRowZ8fX3ztk2bNg3+/v4IDAxEaGgoHj9+bB4xycnApk1AZqZh43KpV0/csYtqW7nS9Dqh37dcvv32WygUCsTFxZlFy+PHwKFD4h5cmm1VqgBOTuI1e3b+9o0bzSJVr28zZ86Eu7s7AgMDERgYiH379plHzLVrwN9/AzpdycatXCku1CZNxL+56ZisLOB//zOLVEPn23fffQdvb2/4+PggLCzMLFpkCTGSs2ULkZMTkWhji9eXXxKtWVN4W+7Lyopo/HiitLT8bWo1UWqq6bUeP36czp8/Tz4+PnnbkpKS8t4vXryYxo0bZ3ohRETffUfk4GCccfv3E3XrRnTvHlGnTvnbX3iBSKs1uVR9vhERPXz4kEJDQ8nT05OePXtmch1ERFOnEtnalmybl5c41spK/zlYrx6RTmd6rfp8mzFjBs2fP9/0P16UMWOIVKqSjWvcmCgpqfgFnftq3dosUvX5duTIEerUqRNlZGQQEdHTp0/NokWOcI9PBixbJlJJxvD8uegZfvcdUHCCCCsrwBwN33bt2qFatWqFtuXO7QgAaWlp5lu8dPny4tWb+ujfX/QK9+8vvi8rCzh5suK1FUGfbwAwYcIEzJs3z6wLvq5aJewwhgcPRM3Q6tXACy/kb3/6VIxHmxpDvknCxo2ljyO3bi1MmzVLpDqvXAH69cvff/kyEBNjWp3Q79uyZcvw2Wef5c3SU6NGDZPrkCsc+GRAdLTxx1atCri4AO+/X/ixtKwss1xPBpkyZQo8PDywceNGfPHFF+b50SdPSj/G0VG0ED76yPAx5krNFmHXrl1wd3dHQECAWX/336XjSiQuTjxD6uUFtGghUp0F05sqlWS2AQCWLl0Kf39/jBo1Cs/LWmhSHtLTjas6q10b8PMT44C1aokLdd06wNtb7Le1Ne68NQG3bt3CyZMn0aZNGwQHByMiIkISHXJAloFv4cKFUCgUFvO6c+dBmfzRaERn5+efgerVxbbMzGx89FFYheqaP3++0Zq++uorREVFYciQIVi6dKkezRq4ublVqL4kY+7gM2cC69eLVrgeUlNSMGrIkArVZczq2xqNBnPmzCm1kaDT6WBlZVWh+rKzSy/BT0sDzp8X1fqxseL+3aWLaEcAQFJSIrp06VWhui5fvlyqLgB49913ERkZiUuXLqFmzZqYOHGi3uMqUlsVtRqZxjy6kBsgZ88Wg6gnTgBHj+Y9N5qYnIyXmjevUG1RUVFG+abVapGQkIDTp09j/vz5GDBgAKhISW9KSkqFamvQoAFSjE1nmRFZBr4JEyaAiCzm5evrVWaPlEpRIe3uLj6r1dZYsWJeheqaNGlSmXUNGTIE27dvL7ZdrVbjyZMnFarPpVat0gV16iQeNI6JES8PD2DrVvEAJABHZ2es3r69QnUNHDiwVFmRkZG4d+8eAgICUKdOHURHR6N58+Z4UqQ3oFQqkZOTU6H61Gor4/6YBSDK1SP+dXGpguPHd1eoLmN7vq6urrCysoJSqcTYsWNx9uxZA5orTluiVgtbhRHpaH1rF1J+cKni6IhT165VqDYPDw+jfKtduzb69esHhUKB1q1bQ6lUFiuocnJyqlBtd+7cgZOTk1H6zIksA5+l0b69mOgByJ9gxMqq8PvOncVzV0qlSDstWCDG+wpOJ9W6tSTycfv27bz3u3btgnduWsfUvPxy/p3YkHGdOgG+vsK8wECRnxs3TpTRAqJ13rKlefQWwM/PD7Gxsbh//z7u37+P2rVr48KFC3BzczP5b7dokf/ekG2tWwONGolJSapVE5OPHD0qCmkBMUYo1aINMQVy+r/99pveStkKx8oqP12Z+1mfcSdOiGdK//Mf8fnll8XD7r//nv/devVMr1cPffr0wdGjRwGItGdWVhZefPFFSbRIDjGSc/s2kZ0dEUA0Y0bx/TNmEPXvT3T9OlFKClFsLNGePUR+fpRXLObnZx6tAwcOJDc3N1KpVOTu7k6rVq2ifv36kY+PD/n5+VHPnj0pOjraPGJOn86v6jRkXNGquqJVnZ06mUWqPt8K4uXlZbaqzj178osODdk2cCDR3buiUvjxY6J164hcXcV3lEqiwYPNIlWvb0OHDiVfX1/y8/Oj1157jR4/fmweMWvXEjk6ln6+NW1K9NdfwryrV4n69BHbbWyIJkwwi1R9vmVmZtKQIUPIx8eHmjVrRocPHzaLFjnCM7fIhLZtgb/+Kt93HR3FmN+QIRWrSfYQicmBC/Q4y4Sjo5iGxMLm7czJEROfP3tWvu+r1aIQtnnzitUle9LTgRo1yv/0vp2deBawbt2K1cWUGU51yoS5c8WED2XFykpci/37V7wm2aNQAN9+Wz7jbGyAhg1FDtnCsLIS55taXfbv2tmJ7J3FBT1AnGdTppTPOHt7oG9fDnoywWrmzJkzpRbBAJ6eolDljz+Mn3LSyko83nDqlBiHsUgaNxZBLDzc+OmgbGzE8hbh4YCDg2n1yZRmzcRsWpcuGW+bra0Ynjp8WFhokbRtKzIMt28bb5y9vVgma9cucdEyksOBT0YEBor7+N694vooKQA6OooCxVOnxL8WTdu2YhG5Q4eEcSWVnTs6CpP//FPM/WbBhIQIu8LDRedZpzN8rKOjqAE6elQUV1ksCoWYczMpKf9BWkPGKRSiYRUSAuzZI1oOjCzgMT4ZEhsrpvtbvFhUzxEBqalpcHJygFYrqukmTwZeey2/GpSBmAnghx/EiwjQ6ZCWmgoHJyfROg8KAiZNEg+kccs7j7t3xUxAq1aJItmcHDEDj5OTA7KzRVHipEmi+tiYin6L4do1cZFu2CCe6M/JQWpaGhydnUW1cPfuwKefivOOjZMVHPhkTE6OmNg9OhpITs7Giy9ao1kz0WFhSiA7W3RNYmKQnZoK6+rVRX1+nTpSK5M1GRnAkSOi4ZWRocMLLyjx8sv5z4oyBkhNFcbFxUGn1UL5wgtifcjc2SUY2cGBj2EYhrEouKqTYRiGsSg48DEMwzAWBQc+hmEYxqLgwMcwDMNYFBz4ZIhGo8G8efMwf/58ZGRkYO3atejVqxfCwsKQWt7pkiyURo0aSS2BYRiZwYFPhowYMQJPnz7FvXv30KNHD5w7dw6TJk0CEeHdd9+VWp5scXJygrOzM5ydneHk5AQnJydERkbmbWf0c6XAUjrZ2dmYPXs2evXqhc8//xwajUZCZfJm6dKlecv63LlzB+3atUOVKlXQpk0b/P333xKrY0qCH2eQIYGBgbh06RKICDVr1kRMTAwUCkXemmVX9K35xeDDDz9EYmIi5s+fD1dXVwBA3bp1ce/ePYmVyZvmzZvjwoULAICJEyciPj4eI0eOxM6dOxEfH4+ff/5ZYoXyxMfHB1evXgUA9OjRA2PGjEHfvn1x7NgxTJkyBeHh4RIrZAyhkloAYxiFQoHu3btD8e+sD7mrGjP6WbJkCc6fP49BgwahT58+eP/999kvIyjY9j18+DAiIiJgbW2Ndu3aGb04rCWiLTCnYGxsLPr27QsAaN++vSxXHWfy4VSnDGnZsmXeWN7q1avztuem7RjDtGjRAocOHQIABAcHIyMjQ2JF8icpKQk7duzA9u3bkZmZCet/58HjhlbJ9O/fHyNGjMDdu3fRt29fLFq0CA8ePMCaNWvg6ekptTymBDjVKVPOnj0LhUKBVq1a4dq1azhw4AAaN25cqAfIFKegbydPnsTRo0fRsmVLdO/eXWppsmXkyJGFPn/zzTdwdXXFkydPMGTIEBw+fFgiZfJn7dq1WLZsGSIjI5GZmQkPDw/06dMHkydPhouLi9TyGANw4JMhs2bNwv79+6HVahESEoIzZ86gQ4cOOHjwILp06YIpU6ZILVGWFPXt7NmzaN++PftmBGfOnIFSqSzU0PL29uYGQykUbGhdvXoVBw4cQJMmTdg3mcOBT4b4+fnh0qVLyMzMhJubG6Kjo+Hs7Iz09HS0adOGi1sMwL6VD24wlA/2rfLCgU+GNGvWDBf/Xeur4Hsgv+KTKQ77Vj64wVA+2LfKCxe3yBAbG5u856fOnz+ftz0pKQlKJf/JDMG+lQ+VSgUrKyuo1WrUr18/75lHe3t79q0E2LfKC/91ZMiJEyegVqsBoNAFlJ2djXXr1kklS/awb+WDGwzlg32rvHCqk2EsnMzMTNja2hbbHhcXh5iYGPj5+UmgSv6wb5UXDnwMwzCMRcH9cbmSkwPs2wf06AH4+iLpxReB5s2BDz4Abt+WWp1syc4Gtm8HQkMBHx+gevVEtGwJhIUBDx9KrU7GZGYCGzcCHToATZsixc0NaNUKmDEDiImRWp1s0WiAn34C2rUDvL2BmjWTERQEzJ0LxMdLrY4xBPf45IZWCyxcCPz3v+KqKroag7U1YGUFNGsGfP01EBwsjU6ZkZkJzJkDfPedsLDojFE2NoBSCbz8MvDNN+KezkAY9cUXwIoV+Z8LYmcHEAEhIcI4Hx/za5Qh8fHAtGnAunXivCp6mdrbC9t69RKXab160uhk9MOBT06kpYkr5fRpEfRKw94eWLAAeOcd02uTMc+fA507A9evA+nppR+vVgNr1gADBphem6yJiQHatwcePBAth5JQKIRxO3aI7rQFc/eu6OHFxooMQ0kolYCjI3DgAPDSS+bRx5QOBz65oNWKG8qpU0BZ5pdUq4Fly4Bhw0ynTcakp4sbyvXrQFaW8d+ztwd+/RWw2Ak2kpJE6vzhQ3HuGYtaDRw8KLrOFsjTp0BAAPDsGaDTGf89R0fgr78ArneRBzzGJxfmzwfOnCke9Ly9gcOHgcREMbbXp0/h/RqN6PHdv282qXLis8+AmzeBsWOBiAhh35o1hY/p2FEExrQ04MgRwNNTBMwBA4CEBGl0S864cUB0tPhXn3HW1sC2bcC9eyJnl5tS12jEuLOFTv49aJBIc777rn7b2rQB/vhDHBMbC2zdCri5iVRot25lC5aM6eDAJwdyckTKsmh608oK2LUL2LMHqFYNePttYMMGoGHD4t///nvz6S3A4sWL4evrCx8fHyxatMisv51bWJCRATx+DMyeDRRYzAIA8MILIjs3bZqw8Nw5YMsWsY+oeJA0FaNGjUKNGjXg6+ubty0hIQEhISFo2LAhQkJC8Pz5c/OIiYsDdu4UXWRDxgHAn38CQ4cWL27RakVQNAP6fNu2bRt8fHygVCpx7tw5s+gAgMhIkZDRag3bVrUqsHIlUKcO4OUlhkxzz7GkJOD3382jVZ9vkyZNgre3N/z9/dG3b18kJiaaR4wcIUZ6du8mcnIiEvfi/JePD1FKSuFtv/9O9MUXxY91ciLKyDCr7L///pt8fHwoLS2NsrOzqVOnTnT79m2z/f5PPxE5Oha24csvidasyf88dixReHj+Z7WaSKMhatxYfK5Zkygnx/Rajx8/TufPnycfH5+8bZMmTaKvv/6aiIi+/vprCgsLM70Q8WNE9vYlG1fwFRVFFBxc/Nw0A/p8u3btGt24cYOCg4MpIiLCLDqIiD76iMja2njbAKJmzYiSk/M/t29vHq36fPv9998pOzubiIjCwsLMd77JEO7xyYGlS4tX0xlCoQAKtOIKsX9/xWkyguvXr6NNmzZQq9VQqVQIDg7Gjh07zPb7331XvJquKD4+wOXL+Z81GtFyzy1OTE0VrXhT065dO1SrVq3Qtl27dmH48OEAgOHDh2Pnzp2mFwKI7IAxVUAlce+eyDGbGH2+NWnSBI0bNzb5bxdl9erSi1mK0q4d8O8i7QDEufbsWcXq0v+7xX0LDQ2FSiXWHg8KCkJ0dLTphcgUDnxywND43M2bYqBg0iRApRIl5cHBosCgKFlZYszGjPj6+uLkyZOIj4+HRqPBvn37EBUVZbbff/So9GMcHUWKqSBJSUDB9XzNKLkQT58+Rc2aNQEAbm5uePr0qXl+uCLuvNbW0hknAZmZYoy4LPj5AdOni8s3F1tbkSaVmtWrV6Nbt25Sy5AMWQa+hQsX5q3+bAmvh7du6TdCqxXFLD16AE+eABMnitFyPQEuOzMTYR98YDKN77zzTt68hLk0adIEkydPRmhoKLp27YrAwEBYWVnla8rONqlvz54lFbWhGKmpwL9zB+fh7JzfwU5JScWgQaMqVNfmzZtL1VWU3O8WRafTwcrKqkL16Up7dMEIEpOS0CskpEJ1XS7YNa8AKlKbnV1V6HTG+1a/vkjAfPSRGCrNJTk5EYGBL1WotrI2Nr/66iuoVCoMGTKk2L6UlJQK1dawYUOkGJvNMiMqqQXoY8KECZgwYYLUMsyHr2/hfEhB/v5bPGuVS3i4eGq2CNb29pi3cCHmjRtnGo0GGD16NEaPHg0A+Pzzz1G7du18TdbWIBM+LVOrVumTily9CvybTQQgOsv16+fb7ezsiNWrV+P11/UUd5gYV1dXxMTEoGbNmoiJiUGNGjWKHaNUKpGTk1OxP+zgYNxzoiVQxcUFu3ftkvUEChV57uXkiE6uMf9JT0/g0CHgyy9FLVpBnJyq4PTpU2jatMKklYm1a9diz549OHz4sN6GlpOTk0mvWbkgyx6fxREUJFKZ+vDzE/kRe3vR46tZE1i7tvhxSqV4wMjMxMbGAgAePnyIHTt2YPDgwWb77ebNxZAnIApgbW3FvwXf//abaFf06ye2TZ8OXLmSPzyVnS2JbQCAXr165a0asW7dOvTu3ds8P1xw9hVDxgFiupvcSZgLvgdE7k+qu7cEWFmJSs2Cn/XZVquWeGRm6dL8yXAKotUW/u+YkwMHDmDevHnYvXt33iomFovExTUMEdHffxevsst9zZtHlJAgqjv37SOqX1//cfXrE+l0Zpf+yiuvUJMmTcjf358OHTpk1t8+ciS/qnPGjOL7Z8wQ+zp1Irp+XVRzHj1K5OWVb1tQkHm0Dhw4kNzc3EilUpG7uzutWrWK4uLiqGPHjtSgQQPq1KkTxcfHm0fM1q35VcQlGXfvXvF9Xl5ECgVRr16m1fgv+nzbsWMHubu7k42NDdWoUYNCQ0PNouX774kcHEq2bfp08T4lpfALIFKpiN5+2yxS9fpWv359ql27NgUEBFBAQACNGzfOPGJkCM/cIheaNwcKrBheJhwcxPyeY8dWrCaZQySelSpvjYWTE7B+PWCujpZsyM4GatQQkyKUBwcH8UBa27YVq0vmpKQArq7lL4i1twcuXBBzUjDSwqlOufD11+LKKCsKhbiDmzHFKBcUCjExdXmyNlZW4ibWo0fF65I91tbiif7yGGdjI1KcFjhlmZOTWBylPLbZ2YnhUA568sBq5syZM6UWwQBo0EAMFISHG/+wUG7QCw8XY38WiL+/6Lhcvmy8bVZWYoaNU6eAKlVMq0+2BAUBt26JafCMNc7aWsy/dfJk+e7+/w/o1EmcN48eGW+bra1YneHgQdFuYKSHA5+cePllcUc+elQEtZIm9nNwEMeGh1t8MzIkRFh1+rT4XJJtjo6iAOHUKcDDwzz6ZIlCIXK8z56JyuGcnJJLFh0dReMsPByoXt18OmWGQiHmeL11S0yEUNr83o6OYgWxo0eLP1bDSAeP8cmRW7eAJUtE9aZSCWRkICc7G1a2tqL6s0YN8VTs0KGFn8S2cK5cEUOdW7aIXl1mJpCdnQM7OysolWI8MCwMePPN8mWV/98SESHWf9y9W5xfGRnI0WphlWuStzcweTLQty93Wf6FCDhxQth26FD++abVaqFWq6DTiWH7sDCgZ8/8QllGHnDgkzMaDfC//wHR0dBpNFBWqSKaj23b5tfxM8VIThb38CdPgPj4VNSu7Yg2bYCWLaVWJnPi48WE6LGxyEpNhY2rq5hzy9AUeQwAcZ7t3Svm/s7IyIKrqw06dSo+lzwjHzjwMQzDMBYFV3UyDMMwFgUHPoZhGMai4MDHMAzDWBQc+BiGYRiLggMfwzAMY1Fw4JMhqampmD59Onx8fODi4oLq1asjKCgIa/WtysDkodVqsWLFCnTt2hX+/v7w9/dHt27dsHz5cmSXdelsBgDw9ttvSy1BtuTk5GDFihWYNm0awsPDC+2bPXu2RKoYY+DHGWRI79690bdvX3Tu3Blbt25FWloaBg4ciNmzZ8Pd3R1z5syRWqIsGTRoEKpUqYLhw4fnrQsYHR2NdevWISEhAVu2bJFYoTxJSEjQu52IEBAQgGg9Cx8zwJgxY6DRaNC6dWusX78ewcHBWLBgAQCgefPmuHDhgsQKGUNw4JMhAQEBhVajbtWqFSIiIqDT6dC0aVPcuHFDQnXypVGjRrhlYDX7kvZZOlZWVvDy8iq0AKlCoQAR4dGjR8jKypJQnXzx9/fHlStXAIhsw/jx4xEXF4dNmzYhKCgIF8u72gpjcjjVKUMcHBzw559/AgB2796NatWqARCrcXM7xTDVqlXDtm3boCswWadOp8OWLVtQtWpVCZXJm3r16uHYsWO4d+9e3uvu3bu4d+8eXF1dpZYnWwo2CFQqFVauXImAgAB07NgRqampEipjSoMDnwxZtmwZPvnkE1StWhXz5s3DkiVLAADPnj3De++9J7E6+bJ582b8+uuvcHNzQ6NGjdCoUSO4ublhx44d2Lx5s9TyZMvHH3+M58+f690XFhZmZjWVh5YtW+LAgQOFts2YMQMjR47E/fv3pRHFGAWnOmXK9evX8ejRIwQFBcHR0TFv+4EDB9C1a1cJlcmbM2fOQKFQoH79+rhx4wZOnTqFpk2bonv37lJLkzVnz56FQqFAq1atcO3aNRw4cADe3t7sWymwb5UTDnwyZMmSJfjhhx/g7e2NS5cuYfHixej97zLhPGhumFmzZmH//v3QarUICQnB2bNn0b59exw8eBBdunTBlClTpJYoS4r6dubMGXTo0IF9KwX2rRJDjOzw9fWllJQUIiK6d+8etWjRghYtWkRERIGBgVJKkzW+vr6k1WopLS2NnJycKCkpiYiINBoN+fn5SaxOvrBv5YN9q7yopA68THF0Ol1eerNOnTo4duwY+vfvjwcPHnBxSwmoVCpYWVlBrVajfv36cP535U97e3solTycbQj2rXywb5UX/uvIEFdXV1y6dCnvs6OjI/bs2YO4uDj8/fffEiqTNzY2NtBoNACA8+fP521PSkriG1EJsG/lg32rvPAYnwyJjo6GSqWCm5tbsX3h4eFo27atBKrkT2ZmJmxtbYttj4uLQ0xMDPz8/CRQJX/Yt/LBvlVeOPAxDMMwFgWP8cmYiAhg5UogMhKIi0uDu7sDXn0VGDsWqF5danUyhQj480/gp5+ABw8Q//gxXmjYEOjUCRg5EqhSRWqFskSnAw4fBtauBaKjgfj4JNSt64IePYChQ4ECT9QwBcnJAfbuBTZuBGJikJSQAJcGDYC+fYEBAwB73GSYJQAAGQhJREFUe6kVMnrgHp/MIALWrwfmzAGiooCMDHFTysXeXnzu3h2YMQMICJBOq6zIyRGthLlzgbg4QKMRZuaiVgvj+vUDpk8HGjeWTquMyMwEli4Fvv0WSEkBik444uAgbBs6FJgyBfDykkan7EhLAxYsAJYsESampBTen9tSGDUK+M9/AD3DFox0cOCTEVlZ4gazb5+4rkpCoRBBcN06oH9/8+iTLWlpQJ8+wF9/iYBXElZWgJ0dsHMn0LmzefTJlMREIDQU+OcfID295GNVKhEEf/8daNPGPPpky5MnQIcOwP37omVaEtbWgIsLcOwY4ONjDnWMEXDgkwlEwOuvAwcOlH4TKoi9PbB1K9Czp+m0yRqtVqQxz54t/SZUELVa3MVfecV02mRMejoQFATcuCEaXMbi6AiEhwP+/qbTJmsSE4HmzUU6Rqs17jsKBeDsDFy4ANSrZ1p9jFFwza1M+P57cR8uGvSOHhXbUlLEq+jCDOnpwJtvAjEx5tMqK2bOBM6d0x/03nwTuHZN5O/u3Ckc5DQaoEeP4ikqC+GDD4BbtwoHvdxzLPel1YpMXkFSU4GQEMBilzccNgx4/Lh40PPyEmN9CQniYvzuO5FdAESrNiUF6NKlcPqdkQwOfDJApxNjeoaydO+/Dzg5iZe3t/7vL19uWo2GqFOnDvz8/BAYGIiWLVua98czM8UNRp9xnTuL8b6RI4Vx7doBd+8WPiYnB9iwwTxai5CYmIj+/fvD29sbTZo0walTp8z220lJohajaFsh9xxzchJDUunpwLZtxb+fng7s2mUerQW5efMmAgMD817Ozs5YtGiR+QRERwN//CHOu6L88AMQGwvUrAkEBgLBwcD48fn7dTqRIj1+3Hx6C7Bw4UL4+PjA19cXgwYNQkZZsiP/H5FoxhimAAcPEjk6EonmYOHX0aNEo0fr31fwVbUqUXa2+bV7eXnRs2fPzP/DRES//GLYuPBwolGjSjeuTh0inc7s0ocNG0Y//vgjERFlZmbS8+fPzfbbixcTqdUl2zJsGFFkpOH9LVuaTa5etFotubq60v379833o//5D5GtrX5Drl0j6tYt//O8eUTLlxc+RqEg6t7dfHr/JTo6murUqUMajYaIiN544w1as2aN2XXICe7xyYDFi4tX0xXk66+BZ89ElX5wsP5jtFrRGLUoFi3Sb5xSCbRsKZ75uH1bjMd8950oailKXJx4bsSMJCUl4cSJExg9ejQAMQNIFTM+ZrF4cek1QMOHAz//bHj/1avFO9Dm5PDhw6hfvz68zFlmuny5/t4eIM7FgQPFoHutWkC3bmLAviBE4pkRAyvemxKtVov09HRotVpoNBrUqlXL7BrkBAc+GXD7tuF9kyeL8XB3d1Gt/7//6R8fz84G7t0znUZDKBQKhIaGokWLFli5cqV5f9zQmmeuroCNjSh3ffVVkXpq1gyYOrX4sQqF2Y27d+8eqlevjpEjR6JZs2YYM2YM0kor461AHj0qeb+np2hgrVtn+BgbG2nOt1w2b96MQYMGme8Hs7JEjtgQJ06Iqs3kZGHwuXOicrgotrYiZWpG3N3d8emnn8LT0xM1a9aEi4sLQkNDzapBdkjd5dTHggULCIAFvR6UmpHLfe3fT/T++/r2ZREwyWQax40bR2lpacX+VtHR0URE9PTpU/L396fjx4/n7cvKyjKpb4mGTKpSRQgYNix/W79+RBcuFDs2BaCRJtQ4depUysnJKeRZREQEWVlZ0enTp4mI6MMPP6SpU6fm7Y+Pjzfx+ZZT4jk2ZQrRsWOlnYvPCXjNZBqXLVtm8P6QmZlJL7zwAj158qTQ9sjISJPpqQJQhiEzFAqi+/eJPv+cyMaGqFo1op07iebOLXbsc4CCTPi3/emnn4r5lZCQQB06dKDY2FjKysqi3r170/r16/P2HzhwwGR6GjRoQMnJyQb/llIhy8BnaTRtalzQA4j27SP64IPi2+3tiUq4V5iFGTNm0Pz58833g25uho16+JDorbfyP/ftqzfwkbMz0bZt5tNMRDExMeTl5ZX3+cSJE9TdjGM/9vYln2M3bxKNHFnyMS4uYvxZCnbu3EkhISHm/dHsbBHg9JnxwgviGGfn/G29exP9/XfxY52ciP75x6zSt27dSqNGjcr7vG7dOnr33XfNqkFucKpTBjRrJoaliuLiIh4wtrUVldGDB4vixKJDB4DY36SJ6bUWJC0tDSn/Pg6QlpaGP/74A76+vuYTUNIDwWvWiJr96tXFNGUTJgB79hQ/TqsFmjY1nUY9uLm5wcPDAzdv3gQgxquamlFDgwaG9730kkir66vmLEhmJtCoUcXqMpZNmzaZN80JiCf4DY2LxceLAc933xUXoouLGCS9cqX4sdnZIpdsRjw9PXH69GloNBoQEQ4fPowm5r5ZyA2pIy9DFBFB5OBQvHH44otEZ88SJScTPX9OdOoUUefO+hudtWubvzgxMjKS/P39yd/fn5o2bUqzZ882r4B9+wxXdapURN9/L4yLiRGljPoq8iRa2PfixYvUokUL8vPzo969e1NCQoLZfvvnnw3btny52F9a5sHcHa5cUlNTqVq1apSYmGj+H//vfw2XwwYEiC5wQgLRs2dEW7YQ1ahR+BgrK6KhQ82vm4imT59OjRs3Jh8fHxo6dChlZGRIokMu8MwtMqFJk+IPpxuLWi0qPz/8sGI1yR6dTjw3FRtbvu87OYmKoYEDK1aXzMnIEB3hkiqJS8LREfjtNwuc8e35c9HrK+8zcGq1mPYmMLBidTFlhlOdMmH2bHFdlAdbW5FZsTiUSjHhtIND2b+rUIjA169fxeuSOXZ2IvNbnvNNpQJq1wY6dqx4XbKnalVxoZXHOBsbEfA46MkC7vHJiLAwMXVZac9YFcTBQUwG0aKF6XTJGiIxO8u2bcYblxv0IiKkG6iSGJ0O6NULOHLE+LlhlUqgWjXg0iUxDmiRZGUB7dsDFy8a3/NTqURm4tIlYSAjOdzjkxFz5wKTJokGpUJR8rF2dqJm4+hRCw56gDBq9WoR/Izp+dnbizzf6dMWG/QAEcR++w3o3ds429Rq0dM7d86Cgx4gem6HDokHHY01rlEj0cjioCcbOPDJCIVCzLl8+LBYZcfOrvg6lk5OIuMSFibGBFu1kkSqvFAqxaJy//ufmAjYzk7kfwvi5CQebJ8xA7h+3fwlsDLE2hr45RexukdwsLDNxiZ/v0IhxvM8PIBvvhHLF/F6fBDBbO9esWpvq1biIlUVWNNbqRRBsUEDMU3OuXPi3GNkA6c6ZcyzZ8DmzcCDB0B8fCZq1bJFy5bAa68Vvs6YIjx6JO7mUVHIfP4cth4eok6/Sxf9z40wAMRMLL/+KuxLTs6Eh4ct2rcXmb3SMhAWzfXrYpaWmBikJiXBsV49ca61acPGyRQOfAzDMIxFwc1fhmEYxqLgwMcwDMNYFBz4GIZhGIuCAx/DMAxjUXDgY/5fEx8fL7UEhmFkBgc+GXKgwPILSUlJGD16NPz9/TF48GA8ffpUQmXy5rPPPkNcXBwA4Ny5c6hXrx7atGkDLy8vHD9+XGJ18qV58+aYPXs2IiMjpZZSqfi/9u49KKrrjgP4d9lVIOzGaBW0MmKbBseQKGrR+EqyY3EIajSVNrXRSmwsldTBRzXWYXw0TmNUUJwYHIlxGhJTZYwmYAa1JdqWahQVa4JBYgV5aCwJPlYg7OP0j2N47bLg1r33Mvf7mdkR9l64P3/cvb97zj333KKiIlitVsyePRuVlZWIjY1Fr169EBMTg7Nnz6odHnnBwqdBK1eubP566dKlGDBgAHJzcxETE4OkpCQVI9O2gwcPom/fvgCAZcuWYc+ePfjyyy9x5MgRLF26VOXotKuurg43btyA1WrF6NGjsXnzZtTU1KgdluYlJydj+fLlmDJlCsaNG4ekpCTcvHkT69evR3JystrhkRcsfBpXVFSEdevWISIiAosXL0Z5ebnaIWmWw+GAw+EAADQ0NCDm7rQ2kZGR+Pbbb9UMTdN69+6NTZs24cqVK0hLS0NZWRlGjhwJq9WKHTt2qB2eZtntdjzzzDOYNWsWDAYDEhISAACTJk1Co69PcCBFsPBp0PXr15Geno60tDTcunULrecYcLlcKkambcnJyYiPj0dBQQHi4uKQkpKCY8eOYfXq1YjmrPhdMnHiRLz55puorq7GK6+8guPHj6sdkmYFBQXh8OHDyMnJgcFgwIEDBwAAx44dg9FoVDk68oYTX2nQ/Pnzm59sPnfuXNTW1qJfv364du0aD+BeLFy4EI8//jgyMzNx8eJFOBwOlJWVYcaMGUhNTVU7PM2K9DBZt9FoRFxcHOLi4lSIqHvYvn07li9fjoCAABw6dAiZmZlITEzEwIEDkZWVpXZ45AWnLNOoL774AtXV1RgzZgzMZnPz+/n5+TwYecG8+YZ5882FCxdQU1PDvHU3Kj35nbzYunWriIyMFNOnTxcRERHiwIEDzctGjBihYmTalpGRwbz5gPubbzIyMsSQIUOYt26IXZ0atGPHDpw+fRpmsxnl5eVISEhAeXk5UlJS2lzvo7aysrKYNx9wf/NNVlYWioqKmLduiIVPg1wuV3O3yeDBg3H06FEkJCSgoqKCHygvmDffMG++Yd66L47q1KCwsDAUFxc3f282m5GXl4fa2lqcP39exci0jXnzDfPmG+at++LgFg2qqqqCyWRC//793ZYVFhZi/PjxKkSlfcybb5g33zBv3RcLHxER6Qqv8WnZf/4D/PnPwKVLaKitRXB4ODBuHPCLXwAPPKB2dJp14QLw7rtAeTnwzTcNGDQoGE89BcycCQQGqh2dhhUXA7t3A5WVaLh1C8GDBgGxscCzzwImHio8EQL49FMgJweoqgLq6m4hMvJBxMcDcXFAAC8maRJbfFqUnw+89hpw8iTgcgFNTS3LzGb5aUtMBH7/e2DwYLWi1BQhgP37Zdo+/xyw24G7s5cBACwW+W9SErB4MfD976sTp+Y4ncBf/gKsXy9PtBob5T73HYtFFr2FC+Xr7lyoetfUBLzzDvD668DVq0B9vdwHv2M2y3PTJUuABQuABx9UL1Zyx8KnJS4XkJIC7NoF3Lnjfd0ePWTz5aOPAKtVmfg0ym4H5s6VqegsbT17ygPSkSPAj3+sTHya1dAAPPcc8M9/dp64oCB59D56FBg6VJHwtOrmTdma+/e/ZcHzJjgYCA0Fjh0DIiKUiY86x8KnFUIAv/2t7KPr7NPU2gMPyBbixIn+i03DXC7ZhXn48L2lzWyWx/vhw/0Xm6bZ7cCkScCpU7KV1xUGgyx+p08DDz/s3/g0qqEBeOIJoLQU6Oq850Yj8L3vyZ7kAQP8Gx91DXugtSI723PRe/nlloPTrl3uP1dfD0yZAuj0gaubNrkXvZ49gbfektf4bt0Czp6VZ+it2WzAT34iD2S6tGIFUFTkXvSys4GaGtmsKS0Ffv3rlmVCALdvy4LpdCobr0bMnw9cvNhx0fvRj+Q+lZ3d8p7TCXzzDRAf37Y7lNTDwqcFQgBr1nhustTUAOvWAW+/3fHPOxzAzp1+C68j8+bNQ2hoKB577DHFtw3I//brr7unzWQCKiuBp54CevUCUlOBvXvdu5oaG+X7SqusrITVasWjjz6KqKgoZGRkKBvAnTvA9u2eq/5rr8nrxr16yUEt69YBI0e2LHe55FG81cOSldLY2IjRo0dj+PDhiIqKwurVqxXd/vXrwL593hvI27bJ89T2HA5ZME+e9F983jidTowYMQJTp05VJwCNYeHTgn/9S36qPNm/H/jwQ+8tuoYGID297aAEBSQmJrZ5WrzS8vJkj1179fXA2rVARYU8pzh4ELh8GRg1qu16NpssnEozmUxIS0tDSUkJTpw4gW3btqGkpES5AHbvlt2WnpSUtAymEkK+2ndr3r6tSuICAwNRUFCAc+fOobi4GPn5+Thx4oRi2+/sgQvPPw/cuAH87W+elzc2yh4KNWRkZGCozq/NtsbCpwVbttzbBSpP6uuBgoL7E08XPfnkk+jTp4+i22wtPV0egzsTGgpERsrRnu1duQKcO3f/Y/NmwIABGHm3FWWxWDB06FBUV1crF0BamvfBLNu2yeWlpXLI4scfu69z6pRMnoIMBkPzFGF2ux12ux2Gjgq4H2zd2nFrz2IB/vhHOYqzIy6XPFm7edM/8XWkqqoKBw8exEsvvaTshjWMhU8LPv/8/+/8t9uBsrL7E083cfFi5+uYTMB778nbIUtL3ZcbjV37Pf5SXl6Os2fPYsyYMcpttKLC+/KXX5ZH8gkTgA8+8HxBKzAQuHTJP/F54XQ6ER0djdDQUMTGxiqWN7sd+O9/O17+6qvyakNn5y89eshueCUtWrQIGzZsQABvKmzGTGiBzfb//w6HQ47k0JHOBqYYDHKQQVMT8LvfeV7H6VQvbTabDTNnzsSWLVvwoJI3enVlOKLLBRQWAuHh8ka09oRQJXFGoxHFxcWoqqrCyZMn8dlnnymyXZtNFi1Phg+XA6U2b+789wQEKJu2vLw8hIaGYlT7fn69U/o5SJ1xOp1i1apVAoBuXiUtV1M6fr36qhC7dnW4/A4gFvgxxqSkJHHnzh23v9fly5dFVFSUx79lU1OTn3N31WvK3n5biIICIYKCvKX2hgB+5rcYU1NThdPp9JibyZMni7S0NLdlX3/9tV/zVt+V/e27V1aWEFu2uL1fBwirH2Ncs2ZNp8eKtWvXio0bNzZ/f+bMGT/mzSQAp8cUpaQIYbMJcfWqfN2+LUR9vRCnT3e0vz3mtzh37tzZJkcrVqwQAwcOFBERESIsLEwEBweLF154oc06+fn5founf//+Ho8batNc4dOln/5UCIPB84HHaBQiMFCIP/1JiHfekV8bje7rWSxCHDqkeOjeCp+/jR/f8fE6M1OI48eFCAnxflwPCRHizBll43a5XGLOnDkiJSVF2Q1/JzLSczL69RPi+edlUgIChJg8WR7Rp01zXzcoSIjyckXDvn79uqirqxNCCFFfXy8mTJggcnNzFdt+v36e0xYcLERYWMtr40YhcnKE6NvXfd3AQCHu/hcU98knn4gpU6aos3GNYVenFixaBISEeF6WmiqvqP/hD8CcOfLr1FT39YKD5f1VCpo1axbGjh2L0tJShIeHY6fCt1QsWdIyFVlrgwbJuQCio4Fr1+QAmNu3gV/+0n3d8HBgxAj/x9paYWEhsrOzUVBQgOjoaERHR+NjTwNI/GXJEs/7mxCyW1NOOimHIC5aBOTmuq87apTiU5FcvXoVVqsVw4YNQ0xMDGJjYxUdnr9woZzApr2GBuCrr1peNpv8mNbWtl3PYJC33D70kDLxUsc4c4sWCAH84AedDzroSHAwsGqVvClZRxwOOWKzrs63nzeb5Ui9F1+8v3Fpns0mE+fr3fsWixwxNG3a/Y1L4776St7i2NWJbtoLCZFT5Y0de1/DIh+wxacFBgOwerXvT1wwmdrOsKETJhOwbJnvaQsMlPde6Y7ZDPzmN/KE6V4ZDLLJEh9//+PSuLAwYMYMz62+zphM8nbIJ564/3HRvWOLTyuEkMVrz557n6szL0+3E1W7XMD06fKm4XtpwISEAH//e9tJSXTFbgeefho4c+be5uq0WORUZ4884tfwtKq+Hhg9Wt451PqhKd4EBAB9+si5OgcO9G981DVs8WmFwSCnhpgzp+Prfa2ZTHK9Awd0W/QAeVDZt082QLqSth495Gxcf/2rjoseIBNx+DAwfnzXEhcYKGdaLizUbdED5HnmP/4BDBvWtZ6GoCBZ7D79lEVPS1j4tMRoBDIzgfffl30iQUHuNw+FhMguqhdflKeQsbHqxKohPXvKB4Hu3CkHtAQHy1S2ZjbL14IFwPnz7HICIPel/HzgjTeAIUPk9+1vcrZY5JnCkiVyogWV5mXVkt695ZM9NmyQ1/zMZvcZ4CwWeZ6wcqWcGeiHP1QlVOoAuzq1rKxMPpHh0iXk7d2Lqb/6lbwyPnu2/LSRR+fPy4eEVlQAOTm5SEychqefBn7+c98ua+nGqVPNT2DP3bcP0+bNAyZPls/s69lT7eg0SQhZBPfskfPJ79//IebPn46pU+UIzvYnYKQNLHzdhMFgAP9U94558w3z5hvmrXtgVycREekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKCx8REekKC183kZGRoXYI3RLz5pvNmzerHUK3lJ6ernYI1AUGIYRQOwgiIiKlsMVHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES6wsJHRES68j/dTBXFRDjOGQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] @@ -181,7 +181,7 @@ "import matplotlib.pyplot as plt\n", "\n", "Lx = 4\n", - "Ly = 4\n", + "Ly = 6\n", "#fig = plt.figure()\n", "#draw2DLattice(Lx,Ly,fig)\n", "draw2DLattice(Lx,Ly)\n", @@ -251,16 +251,17 @@ }, { "cell_type": "code", - "execution_count": 139, + "execution_count": 130, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ - "[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16]\n", - "[2, 3, 4, 1, 6, 7, 8, 5, 10, 11, 12, 9, 14, 15, 16, 13, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 1, 2, 3, 4]\n", - "[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]\n" + "[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24]\n", + "[2, 3, 4, 1, 6, 7, 8, 5, 10, 11, 12, 9, 14, 15, 16, 13, 18, 19, 20, 17, 22, 23, 24, 21, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 1, 2, 3, 4]\n", + "[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]\n", + "48\n" ] } ], @@ -268,7 +269,8 @@ "l1, l2, ltype = getInput(Lx,Ly)\n", "print(l1)\n", "print(l2)\n", - "print(ltype)" + "print(ltype)\n", + "print(len(l1))" ] }, { @@ -294,12 +296,8 @@ }, { "cell_type": "code", - "execution_count": 149, - "metadata": { - "jupyter": { - "source_hidden": true - } - }, + "execution_count": 2, + "metadata": {}, "outputs": [], "source": [ "\"\"\"\n", @@ -382,12 +380,8 @@ }, { "cell_type": "code", - "execution_count": 132, - "metadata": { - "jupyter": { - "source_hidden": true - } - }, + "execution_count": 116, + "metadata": {}, "outputs": [], "source": [ "\"\"\"\n", @@ -404,7 +398,7 @@ " for i in range(1,n-1):\n", " labelscenter[(i,j)]=countlabels\n", " countlabels+=1\n", - "\n", + "# print(labelscenter)\n", " edgescenter = dict()\n", " countedge = 1\n", " for j in range(1,m-1):\n", @@ -416,15 +410,27 @@ " edgescenter[((i,j),(i,j+1))] = countedge\n", " countedge += 1\n", "\n", - "\n", + "# print(edgescenter)\n", " l1 = []\n", " l2 = []\n", " itype = [1 for x in range(0,len(edgescenter))]\n", - "\n", + " \n", " for x in edgescenter:\n", - " if Lx + 1 not in x[1]:\n", + "# print(x[0],\"\\t\",x[1])\n", + "# print(l1)\n", + "# print(l2)\n", + " if Lx + 1 not in x[1] or Ly + 1 not in x[1]:\n", + "# print(\"-->\",x[1])\n", " l1.append(labelscenter[x[0]])\n", - " l2.append(labelscenter[x[1]])\n", + "# l2.append(labelscenter[x[1]])\n", + " if x[1][0] > Lx and x[1][1] <= Ly:\n", + " l2.append(labelscenter[(x[1][0]-Lx,x[1][1])])\n", + " elif x[1][0] <= Lx and x[1][1] > Ly:\n", + " l2.append(labelscenter[(x[1][0],x[1][1]-Ly)])\n", + " elif x[1][0] > Lx and x[1][1] > Ly:\n", + " l2.append(labelscenter[(x[1][0]-Lx,x[1][1]-Ly)])\n", + " else:\n", + " l2.append(labelscenter[x[1]])\n", " else:\n", " l1.append(labelscenter[x[0]])\n", " if x[1][0] > Lx:\n", @@ -440,24 +446,24 @@ }, { "cell_type": "code", - "execution_count": 133, + "execution_count": 117, "metadata": {}, "outputs": [ { "data": { "text/plain": [ - "([1, 2, 3, 4, 5, 6, 7, 8, 9, 1, 2, 3, 4, 5, 6, 7, 8, 9],\n", - " [2, 3, 1, 5, 6, 4, 8, 9, 7, 4, 5, 6, 7, 8, 9, 1, 2, 3],\n", - " [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1])" + "([1, 2, 3, 4, 5, 6, 1, 2, 3, 4, 5, 6],\n", + " [2, 1, 4, 3, 6, 5, 3, 4, 5, 6, 1, 2],\n", + " [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1])" ] }, - "execution_count": 133, + "execution_count": 117, "metadata": {}, "output_type": "execute_result" } ], "source": [ - "getInput(3,3)" + "getInput(2,3)" ] }, { diff --git a/src/ex1.c b/src/ex1.c index 98f89f0..1113ced 100644 --- a/src/ex1.c +++ b/src/ex1.c @@ -35,20 +35,24 @@ int main(int argc,char **argv) PetscLogDouble t1,t2,tt1,tt2; PetscErrorCode ierr; int mpiid; + SlepcInitialize(&argc,&argv,(char*)0,NULL); + MPI_Comm_rank(MPI_COMM_WORLD,&mpiid); char const* const fileName = argv[1]; FILE* file = fopen(fileName, "r"); Data getdata; - PetscInt nlocal; +//PetscInt nlocal; /* gather the input data */ + if(mpiid==0)printf("Reading data\n"); Data_new(file, &getdata); getdata.n = get_ntot(getdata.FAM1, getdata.natom, getdata.isz, getdata.ntrou, getdata.fix_trou1, getdata.fix_trou2); + if(mpiid==0)printf("Done reading data\n"); - nlocal = getdata.n/getdata.npar; +//nlocal = getdata.n/getdata.npar; PetscScalar *valxr; - PetscInt indxr[nlocal]; +//PetscInt indxr[nlocal]; char filename[PETSC_MAX_PATH_LEN]="FIL666"; PetscViewer viewer; PetscBool ishermitian; @@ -77,7 +81,7 @@ int main(int argc,char **argv) PetscReal trace2rdm=0.0; PetscReal trace2rdmfin=0.0; IS from, to; /* index sets that define the scatter */ - PetscInt idx_to[nlocal], idx_from[nlocal]; +//PetscInt idx_to[nlocal], idx_from[nlocal]; PetscScalar *values; int ndim=(getdata.natom/2)*((getdata.natom/2)-1)/2; double a, b, c; @@ -85,16 +89,16 @@ int main(int argc,char **argv) double gamma_pfin = 0.0, gamma_mfin = 0.0; double nel, s2dens; double nelfin, s2densfin; - double densmat2[getdata.natom][getdata.natom][getdata.natom][getdata.natom]; - memset(densmat2, 0, sizeof(densmat2)); +//double densmat2[getdata.natom][getdata.natom][getdata.natom][getdata.natom]; +//memset(densmat2, 0, sizeof(densmat2)); - SlepcInitialize(&argc,&argv,(char*)0,NULL); + if(mpiid==0)printf("Initializing Slepc vars\n"); ierr = PetscPrintf(PETSC_COMM_WORLD,"\n1-D t-J Eigenproblem, n=%D\n\n",getdata.n);CHKERRQ(ierr); ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); ierr = MatCreateAIJ(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,getdata.n,getdata.n,10*getdata.natom,NULL,10*getdata.natom,NULL,&A);CHKERRQ(ierr); ierr = MatMPIAIJSetPreallocation(A,10*getdata.natom,NULL,10*getdata.natom,NULL);CHKERRQ(ierr); + if(mpiid==0)printf("Done Initializing Slepc\n"); - MPI_Comm_rank(MPI_COMM_WORLD,&mpiid); ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); ierr = PetscTime(&tt1);CHKERRQ(ierr); ierr = PetscPrintf(PETSC_COMM_WORLD," start: %d end: %d \n",Istart, Iend);CHKERRQ(ierr); From ef3a4bd50fee7f0ba150267ed1b570b356eb35a8 Mon Sep 17 00:00:00 2001 From: vijay Date: Tue, 5 May 2020 19:47:20 +0200 Subject: [PATCH 06/20] Update README.md --- README.md | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/README.md b/README.md index 9b744f9..0a779e7 100644 --- a/README.md +++ b/README.md @@ -7,6 +7,24 @@ Double Exchange Hamiltonian: Complete Version (under GNU GENERAL PUBLIC LICENSE v2) +This program can perform Exact diagonalization calculations of various types of +model Hamiltonians. It is especially optimized for the t-J (or Double Exchange) +type model Hamiltonians. The core feature which the program is specialized for +is the adressing of determinant in an efficient manner to quickly construct the +Hamiltonian non-zero matrix-elements. Once the Hamiltonian is constructed in +its sparse format, it is stored in distributed memory for all linear algebra +operations. + +The main work of diagonalizing the Hamiltonian is performed using PETSc and +SLEPc helper functions. These functions return the eigenvectors which are +not stored to disk by default due to their large size. + +This project also contains subroutines which analyze the wavefunction in +its distributed memory form and calculates the various observables. The +output of the program are the energies and the various observables such as +the total Spin, various Spin-Spin correlation functions, and one-and two-body +density matrices. + _Dependencies_ --------------- From ca4a0c97ef626f5ec641c311fe3baf685ac9934f Mon Sep 17 00:00:00 2001 From: v1j4y Date: Fri, 22 Apr 2022 15:17:20 +0200 Subject: [PATCH 07/20] Fix for heis. --- src/providdet.irp.f | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/src/providdet.irp.f b/src/providdet.irp.f index afac84c..8fa422e 100644 --- a/src/providdet.irp.f +++ b/src/providdet.irp.f @@ -36,6 +36,10 @@ endif deth(count,2)=i-1 enddo + else + count+=1 + deth(count,1)=0 + deth(count,2)=1 endif !C if det=0 then exit From 81f32fbe058eef5e96200267d8e899974527445f Mon Sep 17 00:00:00 2001 From: v1j4y Date: Fri, 22 Apr 2022 15:35:46 +0200 Subject: [PATCH 08/20] Update readme for ms. --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 0a779e7..31c6958 100644 --- a/README.md +++ b/README.md @@ -66,7 +66,7 @@ _Using DEHam_ 140 # The largest number of non-zero elements per row (Multiple of Ndet) 1 # The total number of processors used in parallel (Multiple of Ndet) 1 # The number of holes -0 # The isz (ms-1/2) value +0 # The isz (nalpha-nbeta-1) (odd spins) value or (nalpha-nbeta) (even spins) true # Restrict the hole to the 1'st (i.e. half of natom) Family of states. *false* for no restrictions 1,2,3,1,2,3,4,5,6,7 # The topology of the system is specified here 2,3,4,8,7,6,5,6,7,8 # first and second line contain the two sites linked From 60e398661b15adae8e3cf423a2eb43d998db3f51 Mon Sep 17 00:00:00 2001 From: Eloh5 <104210679+Eloh5@users.noreply.github.com> Date: Fri, 22 Apr 2022 15:57:02 +0200 Subject: [PATCH 09/20] Update README.md --- README.md | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/README.md b/README.md index 31c6958..d18dd06 100644 --- a/README.md +++ b/README.md @@ -66,14 +66,14 @@ _Using DEHam_ 140 # The largest number of non-zero elements per row (Multiple of Ndet) 1 # The total number of processors used in parallel (Multiple of Ndet) 1 # The number of holes -0 # The isz (nalpha-nbeta-1) (odd spins) value or (nalpha-nbeta) (even spins) +0 # The Sz (nalpha-nbeta-1) (odd spins) value or (nalpha-nbeta) (even spins) true # Restrict the hole to the 1'st (i.e. half of natom) Family of states. *false* for no restrictions 1,2,3,1,2,3,4,5,6,7 # The topology of the system is specified here 2,3,4,8,7,6,5,6,7,8 # first and second line contain the two sites linked -1,1,1,2,2,2,2,3,3,3 # third line contains the type of link (1 for t or J, 2 for K and 3 for none) -.1430,-0.20,0.0000 # The three types of links this line gives J, K -.1430,-0.20,0.0000 # --1.00,0.0,0.00 # This line gives t +1,1,1,2,2,2,2,3,3,3 # third line contains the type of link (1 for t or J, 2 for K and 3 for J in the second family of states) +.1430,-0.20,0.0000 # The three types of links this line gives J, K, J2 along x axis +.1430,-0.20,0.0000 # The three types of links this line gives J, K, J2 along y axis +-1.00,0.0,0.00 # This line gives t 0.,0.,0.,0.,0.,0.,0.,0.,0. # Energy of each orbital + one extra term 2 # The total number of roots 1 # I The position of the first From 7e27a9a755364b18edfc1d80610470dad3a6f7ce Mon Sep 17 00:00:00 2001 From: vijay Date: Sat, 23 Apr 2022 15:49:11 +0200 Subject: [PATCH 10/20] Update for new IRPF90 version --- Makefile | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Makefile b/Makefile index f34b0b9..fc9e2be 100644 --- a/Makefile +++ b/Makefile @@ -25,7 +25,7 @@ ${BIN_DIR}: directories: ${OBJ_DIR} ${LIB_DIR} ${BIN_DIR} ${LIB_DIR}/irpf90.a: directories - cd ${SRC_DIR} && irpf90 init && $(MAKE) irpf90.a && cp irpf90.a ../${LIB_DIR} + cd ${SRC_DIR} && irpf90 init && $(MAKE) && cp IRPF90_temp/irpf90.a ../${LIB_DIR} ${OBJ_DIR}/get_ntot.o: ${SRC_DIR}/get_ntot.c directories chkopts ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} From f09f08f802c943adf9b6c4b79075f8fd4b42b0a3 Mon Sep 17 00:00:00 2001 From: v1j4y Date: Sat, 23 Apr 2022 16:35:09 +0200 Subject: [PATCH 11/20] Removed borken gif link. --- README.md | 3 --- 1 file changed, 3 deletions(-) diff --git a/README.md b/README.md index d18dd06..aefa0b6 100644 --- a/README.md +++ b/README.md @@ -107,9 +107,6 @@ A 2D t-J model Hamiltonian description and setup for using DEHam to solve for fe is provided in the notbooks folder. Please have a look about the details of using DEHam to study t-J Hamiltonians. -![](https://raw.githubusercontent.com/v1j4y/DEHam/master/notebooks/graph.png) - - _Publications using this code_ ------------------------------- From 87bacf09db3b6a97f12b731014e6acf965ee1bbd Mon Sep 17 00:00:00 2001 From: v1j4y Date: Sat, 23 Apr 2022 16:40:03 +0200 Subject: [PATCH 12/20] Restore figure and remove asciicast. --- README.md | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index aefa0b6..7b25a24 100644 --- a/README.md +++ b/README.md @@ -51,8 +51,6 @@ export C_INCLUDE_PATH+=:$PETSC_DIR/include/:$SLEPC_DIR/include:$PETSC_DIR/arch-l make ex1 ``` -[![asciicast](https://asciinema.org/a/Ng3tSNDoWBkV5C9ZYvbxCW43B.png)](https://asciinema.org/a/Ng3tSNDoWBkV5C9ZYvbxCW43B) - _Using DEHam_ --------------- @@ -107,6 +105,9 @@ A 2D t-J model Hamiltonian description and setup for using DEHam to solve for fe is provided in the notbooks folder. Please have a look about the details of using DEHam to study t-J Hamiltonians. +![](https://raw.githubusercontent.com/v1j4y/DEHam/master/notebooks/graph.png) + + _Publications using this code_ ------------------------------- From ce1dc54b3dc0258fa804ed59fa5c5dc646fe7f2d Mon Sep 17 00:00:00 2001 From: v1j4y Date: Mon, 16 May 2022 18:38:53 +0200 Subject: [PATCH 13/20] Fixed bugs for bool type and precision. --- src/elem_diag.irp.f | 3 ++- src/ex1.c | 2 +- src/extra_diag.irp.f | 1 + src/read2.c | 12 +++++++----- 4 files changed, 11 insertions(+), 7 deletions(-) diff --git a/src/elem_diag.irp.f b/src/elem_diag.irp.f index 96fc01e..8fd254f 100644 --- a/src/elem_diag.irp.f +++ b/src/elem_diag.irp.f @@ -3,7 +3,7 @@ subroutine elem_diag(xmatd) implicit none integer :: i - real*8 :: xmatd + real*8,intent(inout) :: xmatd logical :: yw ! write(6,*)'in elem_diag' @@ -32,6 +32,7 @@ subroutine elem_diag(xmatd) if(deter(i).ne.3) xmatd = xmatd + E(i) enddo xmatd = xmatd - E(natom+1) + if(yw)write(6,*)'xmatd=',xmatd !-----stockage de l element diag diff --git a/src/ex1.c b/src/ex1.c index 1113ced..c27b3a6 100644 --- a/src/ex1.c +++ b/src/ex1.c @@ -195,7 +195,7 @@ int main(int argc,char **argv) /* Save eigenvectors, if == ested */ - EPSGetConverged(eps,&nconv); + ierr = EPSGetConverged(eps,&nconv);CHKERRQ(ierr); if (getdata.print_wf) { PetscViewerASCIIOpen(PETSC_COMM_WORLD,filename,&viewer); PetscViewerSetFormat(viewer,PETSC_VIEWER_ASCII_MATLAB); diff --git a/src/extra_diag.irp.f b/src/extra_diag.irp.f index 8bbfd04..c659e68 100644 --- a/src/extra_diag.irp.f +++ b/src/extra_diag.irp.f @@ -26,6 +26,7 @@ subroutine extra_diag(tistart) count1=0 count2=1 tistart2=tistart + xmat=0.0d0 do j=1,nrows diff --git a/src/read2.c b/src/read2.c index 8f8b040..28ece1a 100644 --- a/src/read2.c +++ b/src/read2.c @@ -32,7 +32,9 @@ void Data_new(FILE* file, Data* dat) { dat->isz=atol(line); break; case 6: - dat->FAM1 = to_bool(line); + //dat->FAM1 = to_bool(line); + dat->FAM1 = 0; + dat->FAM1 = line && strcmp(line,"true")==0; break; case 7: arrayIdx=0; @@ -115,7 +117,7 @@ void Data_new(FILE* file, Data* dat) { /** * Convert the next token to a float value */ - val = strtof(token, &unconverted); + val = strtod(token, &unconverted); if (!isspace(*unconverted) && *unconverted != 0) { /** @@ -139,7 +141,7 @@ void Data_new(FILE* file, Data* dat) { /** * Convert the next token to a float value */ - val = strtof(token, &unconverted); + val = strtod(token, &unconverted); if (!isspace(*unconverted) && *unconverted != 0) { /** @@ -163,7 +165,7 @@ void Data_new(FILE* file, Data* dat) { /** * Convert the next token to a float value */ - val = strtof(token, &unconverted); + val = strtod(token, &unconverted); if (!isspace(*unconverted) && *unconverted != 0) { /** @@ -187,7 +189,7 @@ void Data_new(FILE* file, Data* dat) { /** * Convert the next token to a float value */ - val = strtof(token, &unconverted); + val = strtod(token, &unconverted); if (!isspace(*unconverted) && *unconverted != 0) { /** From 51d2fb62eca01314fbf5271040ef81f269fcfeee Mon Sep 17 00:00:00 2001 From: v1j4y Date: Mon, 13 Jun 2022 17:51:24 +0200 Subject: [PATCH 14/20] Fixed bug in reading true or false for FAM1. --- src/read2.c | 24 +++++++++++++++++++++++- 1 file changed, 23 insertions(+), 1 deletion(-) diff --git a/src/read2.c b/src/read2.c index 28ece1a..53f952a 100644 --- a/src/read2.c +++ b/src/read2.c @@ -1,5 +1,26 @@ #include "read2.h" +/* Compare two strings s1 and s2, assuming s1 is terminated + * by \n or a NULL, and s2 is terminated by a NULL. A match + * returns 1, a non-match returns 0. + */ +int +strcmpst1nl (const char * s1, const char * s2) +{ + char s1c; + do + { + s1c = *s1; + if (s1c == '\n') + s1c = 0; + if (s1c != *s2) + return 0; + s1++; + s2++; + } while (s1c); /* already checked *s2 is equal */ + return 1; +} + void Data_new(FILE* file, Data* dat) { //ata* dat = (Data*)malloc(sizeof(Data)); @@ -34,7 +55,8 @@ void Data_new(FILE* file, Data* dat) { case 6: //dat->FAM1 = to_bool(line); dat->FAM1 = 0; - dat->FAM1 = line && strcmp(line,"true")==0; + //dat->FAM1 = line && strcmp(line,"true")==0; + dat->FAM1 = strcmpst1nl(line,"true"); break; case 7: arrayIdx=0; From 3cd7118fb2bb2ac2ece7210b7a5be630a547cf1f Mon Sep 17 00:00:00 2001 From: vijay Date: Thu, 12 Jan 2023 16:04:36 +0100 Subject: [PATCH 15/20] Update README.md --- README.md | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/README.md b/README.md index 7b25a24..bf01a4f 100644 --- a/README.md +++ b/README.md @@ -72,19 +72,19 @@ true # Restrict the hole to the 1'st (i.e. half of natom) Family of states. .1430,-0.20,0.0000 # The three types of links this line gives J, K, J2 along x axis .1430,-0.20,0.0000 # The three types of links this line gives J, K, J2 along y axis -1.00,0.0,0.00 # This line gives t -0.,0.,0.,0.,0.,0.,0.,0.,0. # Energy of each orbital + one extra term +0.,0.,0.,0.,0.,0.,0.,0.,0. # Energy of each orbital + *one extra term* 2 # The total number of roots -1 # I The position of the first -1 # I SBox -1 # I -1 # I +1 # I a1 The position of the first +1 # I a2 SBox (starts from 0 ) +1 # I b1 +1 # I b2 1 # II The positions of the second -1 # II SBox -1 # II +1 # II SBox (starts from 0) +1 # II 1 # II 1 # III 1 # III The positions of the third -1 # III SBox +1 # III SBox (starts from 0) 1 # III 1 # positio of the hole 0 # fix the position of the first hole during the CI From 00f5fd9bc0ef81f710be6ddc45b8709a89af3d7f Mon Sep 17 00:00:00 2001 From: vijay Date: Thu, 12 Jan 2023 16:11:20 +0100 Subject: [PATCH 16/20] Update README.md --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index bf01a4f..727d6be 100644 --- a/README.md +++ b/README.md @@ -108,7 +108,7 @@ t-J Hamiltonians. ![](https://raw.githubusercontent.com/v1j4y/DEHam/master/notebooks/graph.png) -_Publications using this code_ +_Publications using this code please cite the following publication_ ------------------------------- 1. High-Spin Chains and Crowns from Double-Exchange Mechanism [doi:10.3390/cryst6040039](http://www.dx.doi.org/10.3390/cryst6040039) From 161ebebdbe2fa12b0fb6bcef6bc13b83d933afdd Mon Sep 17 00:00:00 2001 From: vijay Date: Thu, 12 Jan 2023 16:12:39 +0100 Subject: [PATCH 17/20] Update README.md --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 727d6be..2cd6316 100644 --- a/README.md +++ b/README.md @@ -108,7 +108,7 @@ t-J Hamiltonians. ![](https://raw.githubusercontent.com/v1j4y/DEHam/master/notebooks/graph.png) -_Publications using this code please cite the following publication_ +_Research using this code please cite the following publication_ ------------------------------- 1. High-Spin Chains and Crowns from Double-Exchange Mechanism [doi:10.3390/cryst6040039](http://www.dx.doi.org/10.3390/cryst6040039) From 0398124bfe5cb8c1aa6c1255a476635831a30d7a Mon Sep 17 00:00:00 2001 From: v1j4y Date: Tue, 7 Feb 2023 22:00:00 +0100 Subject: [PATCH 18/20] Projection working for 6_2h. --- Makefile | 7 +++++-- src/ex1.c | 35 +++++++++++++++++++++++------------ src/get_dmat.c | 43 +++++++++++++++++++++++++++++++++---------- src/get_dmat.h | 3 ++- src/get_s2.c | 10 +++++----- src/get_s2.h | 7 +++++-- src/read2.c | 14 ++++++++++---- src/read2.h | 4 +++- 8 files changed, 86 insertions(+), 37 deletions(-) diff --git a/Makefile b/Makefile index fc9e2be..49fb7c4 100644 --- a/Makefile +++ b/Makefile @@ -45,12 +45,15 @@ ${OBJ_DIR}/get_s2_mov.o: ${SRC_DIR}/get_s2_mov.c directories chkopts ${OBJ_DIR}/get_dmat.o: ${SRC_DIR}/get_dmat.c directories chkopts ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} +${OBJ_DIR}/get_proj.o: ${SRC_DIR}/get_proj.c directories chkopts + ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} + ${OBJ_DIR}/get_val_iaa2.o: ${SRC_DIR}/get_val_iaa2.c directories chkopts ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} ${OBJ_DIR}/ex1.o: ${SRC_DIR}/ex1.c -${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} -${BIN_DIR}/ex1: ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${OBJ_DIR}/ex1.o ${SRC_DIR}/read2.h ${SRC_DIR}/get_ntot.h ${SRC_DIR}/stimsyr.h chkopts - -${CLINKER} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -o ${BIN_DIR}/ex1 ${OBJ_DIR}/ex1.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${SLEPC_EPS_LIB} +${BIN_DIR}/ex1: ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${OBJ_DIR}/ex1.o ${SRC_DIR}/read2.h ${SRC_DIR}/get_ntot.h ${SRC_DIR}/stimsyr.h chkopts + -${CLINKER} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -o ${BIN_DIR}/ex1 ${OBJ_DIR}/ex1.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${SLEPC_EPS_LIB} # ${RM} ex1.o read2.o diff --git a/src/ex1.c b/src/ex1.c index c27b3a6..ac16b83 100644 --- a/src/ex1.c +++ b/src/ex1.c @@ -5,6 +5,7 @@ #include "stimsyr.h" #include "get_s2.h" #include "get_ntot.h" +#include "get_proj.h" #undef __FUNCT__ #define __FUNCT__ "main" @@ -84,13 +85,18 @@ int main(int argc,char **argv) //PetscInt idx_to[nlocal], idx_from[nlocal]; PetscScalar *values; int ndim=(getdata.natom/2)*((getdata.natom/2)-1)/2; + int ndimdmat2=(getdata.natom)*(getdata.natom)*(getdata.natom)*(getdata.natom); double a, b, c; double gamma_p = 0.0, gamma_m = 0.0; double gamma_pfin = 0.0, gamma_mfin = 0.0; double nel, s2dens; double nelfin, s2densfin; -//double densmat2[getdata.natom][getdata.natom][getdata.natom][getdata.natom]; -//memset(densmat2, 0, sizeof(densmat2)); + int nstates = getdata.ntrou*3; + double projvec[getdata.nroots*nstates], weightproj=0.0; + //double densmat2[getdata.natom][getdata.natom][getdata.natom][getdata.natom]; + ////double *densmat2; + ////densmat2 = (double *)malloc(sizeof(double)*ndimdmat2); + //memset(densmat2, 0, sizeof(double)*ndimdmat2); if(mpiid==0)printf("Initializing Slepc vars\n"); ierr = PetscPrintf(PETSC_COMM_WORLD,"\n1-D t-J Eigenproblem, n=%D\n\n",getdata.n);CHKERRQ(ierr); @@ -250,13 +256,17 @@ int main(int argc,char **argv) get_s2(xr, &Istart, &Iend, values, &getdata.natom, &norm, &norm2, &norm3, &norm4, &xymat, &xymat2, &xymat3, &xymat4, &weight3, &getdata.s21a1, &getdata.s21a2, &getdata.s21b1, &getdata.s21b2, &getdata.s22a1, &getdata.s22a2, &getdata.s22b1, &getdata.s22b2, &getdata.s23a1, &getdata.s23a2, - &getdata.s23b1, &getdata.s23b2, &getdata.postrou, natomax); + &getdata.s23b1, &getdata.s23b2, &getdata.postrou1, &getdata.postrou2, &getdata.postrou3, natomax); // get_s2_cyclic(xr, &Istart, &Iend, values, &getdata.natom, &norm, &norm2, &norm3, &norm4, &xymat, &xymat2, &xymat3, &xymat4, // &getdata.s21a1, &getdata.s21a2, &getdata.s21b1, &getdata.s21b2, &getdata.s22a1, &getdata.s22a2, // &getdata.s22b1, &getdata.s22b2, &getdata.s23a1, &getdata.s23a2, // &getdata.s23b1, &getdata.s23b2, &getdata.postrou, natomax); -// get_1rdm(values, &Istart, &Iend, &getdata.natom, &trace1rdm, natomax); -// get_2rdm(values, &Istart, &Iend, &getdata.natom, &trace2rdm, densmat2, natomax); + get_1rdm(values, &Istart, &Iend, &getdata.natom, &trace1rdm, natomax); + //get_2rdm(values, &Istart, &Iend, &getdata.natom, &trace2rdm, densmat2, natomax); + //get_2rdm(values, &Istart, &Iend, &getdata.natom, &trace2rdm, natomax); + get_proj(values, &Istart, &Iend, &getdata.natom, i, projvec, natomax); + weightproj = 0.0; + for(ii=0;ii<6;++ii) weightproj += projvec[(i)*6 + ii]*projvec[(i)*6+ii]; // analyse_(valxr, (Iend-Istart), &Istart, &Iend, &xymat, &norm); VecRestoreArray(vec2,&values); ierr = VecRestoreArray(xr, &valxr);CHKERRQ(ierr); @@ -269,7 +279,7 @@ int main(int argc,char **argv) MPI_Reduce(&norm2, &normfin2, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD); MPI_Reduce(&norm3, &normfin3, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD); MPI_Reduce(&norm4, &normfin4, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD); -// MPI_Reduce(&trace1rdm, &trace1rdmfin, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD); + MPI_Reduce(&trace1rdm, &trace1rdmfin, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD); // printf("done calc densmat\n"); // for(ll=0;ll void get_1rdm(PetscScalar *, PetscInt *, PetscInt *, int *, PetscReal *, const int natomax); -void get_2rdm(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, PetscReal *trace2rdm, double densmat2[*natom][*natom][*natom][*natom], const int natomax); +//void get_2rdm(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, PetscReal *trace2rdm, double densmat2[*natom][*natom][*natom][*natom], const int natomax); +void get_2rdm(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, PetscReal *trace2rdm, const int natomax); diff --git a/src/get_s2.c b/src/get_s2.c index 7a82498..f42d01f 100644 --- a/src/get_s2.c +++ b/src/get_s2.c @@ -26,7 +26,7 @@ void get_s2(Vec xr, PetscInt *Istart, PetscInt *Iend, PetscScalar *valxr, int *n PetscReal *norm, PetscReal *norm2, PetscReal *norm3, PetscReal *norm4, PetscReal *xymat, PetscReal *xymat2, PetscReal *xymat3, PetscReal *xymat4, PetscReal *weight3, int *s21a1, int *s21a2, int *s21b1, int *s21b2, int *s22a1, int *s22a2, - int *s22b1, int *s22b2, int *s23a1, int *s23a2, int *s23b1, int *s23b2, int *postrou, const int natomax){ + int *s22b1, int *s22b2, int *s23a1, int *s23a2, int *s23b1, int *s23b2, int *postrou1, int *postrou2, int *postrou3, const int natomax){ long int iaa2, iaa; long int iii; int ideter[natomax]; @@ -91,9 +91,9 @@ void get_s2(Vec xr, PetscInt *Istart, PetscInt *Iend, PetscScalar *valxr, int *n pos3=kko; } } - if(ntrouboit1==1 && pos1 == *postrou)okboit1=1; - if(ntrouboit2==1 && pos2 == *postrou)okboit2=1; - if(ntrouboit3==1 && pos3 == *postrou)okboit3=1; + if(ntrouboit1==1 && pos1 == *postrou1)okboit1=1; + if(ntrouboit2==1 && pos2 == *postrou2)okboit2=1; + if(ntrouboit3==1 && pos3 == *postrou3)okboit3=1; if(okboit1){ *norm2=*norm2+valxr[iiii]*valxr[iiii]; } @@ -681,6 +681,6 @@ void get_s2(Vec xr, PetscInt *Istart, PetscInt *Iend, PetscScalar *valxr, int *n ierr = PetscTime(&tt2); //printf(" norm = %18f weight = %18f weight/N = %18f tmpwe = %18f\n", *norm2, *weight3, *weight3/(*norm2),tmpwe); -//printf(" norm = %18f %18f xymat = %18f %18f | %d %d %d %d %d\n", *norm, *norm3, *xymat, *xymat3, *s22a1, *s22a2, *s22b1, *s22b2, *postrou); +//printf(" norm = %18f %18f %18f %18f xymat = %18f %18f %18f %18f | %d %d %d %d %d\n", *norm, *norm2, *norm3, *norm4, *xymat, *xymat2, *xymat3, *xymat4, *s21a1, *s21a2, *s21b1, *s21b2, *postrou); //ierr = PetscPrintf(PETSC_COMM_WORLD," Time used for the s2 loop: %f\n",tt2-tt1);CHKERRQ(ierr); } diff --git a/src/get_s2.h b/src/get_s2.h index facd7a2..b253958 100644 --- a/src/get_s2.h +++ b/src/get_s2.h @@ -4,8 +4,11 @@ #include #include -void get_s2(Vec xr, PetscInt *Istart, PetscInt *Iend, PetscScalar *valxr, int *natom, PetscReal *norm, PetscReal *norm2, PetscReal *norm3, PetscReal *norm4, PetscReal *xymat, PetscReal *xymat2, PetscReal *xymat3, PetscReal *xymat4, PetscReal *weight3, - int *s21a1, int *s21a2, int *s21b1, int *s21b2, int *s22a1, int *s22a2, int *s22b1, int *s22b2, int *s23a1, int *s23a2, int *s23b1, int *s23b2, int *postrou, const int natomax); +void get_s2(Vec xr, PetscInt *Istart, PetscInt *Iend, PetscScalar *valxr, int *natom, + PetscReal *norm, PetscReal *norm2, PetscReal *norm3, PetscReal *norm4, + PetscReal *xymat, PetscReal *xymat2, PetscReal *xymat3, PetscReal *xymat4, PetscReal *weight3, + int *s21a1, int *s21a2, int *s21b1, int *s21b2, int *s22a1, int *s22a2, + int *s22b1, int *s22b2, int *s23a1, int *s23a2, int *s23b1, int *s23b2, int *postrou1, int *postrou2, int *postrou3, const int natomax); void get_s2_mov(Vec, PetscInt *, PetscInt *, PetscScalar *, int *, PetscReal *, PetscReal *,PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, PetscReal *, int *, diff --git a/src/read2.c b/src/read2.c index 53f952a..72775a6 100644 --- a/src/read2.c +++ b/src/read2.c @@ -34,7 +34,7 @@ void Data_new(FILE* file, Data* dat) { /* note that fgets doesn't strip the terminating \n, checking its presence would allow to handle lines longer than sizeof(line) */ - if (count != 30){ + if (count != 32){ count++; switch(count){ case 1: @@ -266,15 +266,21 @@ void Data_new(FILE* file, Data* dat) { dat->s23b2=atol(line); break; case 27: - dat->postrou=atol(line); + dat->postrou1=atol(line); break; case 28: - dat->fix_trou1=atol(line); + dat->postrou2=atol(line); break; case 29: + dat->postrou3=atol(line); + break; + case 30: + dat->fix_trou1=atol(line); + break; + case 31: dat->fix_trou2=atol(line); break; - case 30: + case 32: dat->print_wf = atol(line); break; default: diff --git a/src/read2.h b/src/read2.h index 9a1a6ff..3cbe1ad 100644 --- a/src/read2.h +++ b/src/read2.h @@ -34,7 +34,9 @@ typedef struct { int s23a2; int s23b1; int s23b2; - int postrou; + int postrou1; + int postrou2; + int postrou3; long int fix_trou1; long int fix_trou2; long int print_wf; From 6df166e7191632cdc295b4e64eeacc4e4b0f76f4 Mon Sep 17 00:00:00 2001 From: v1j4y Date: Thu, 23 Feb 2023 18:17:30 +0100 Subject: [PATCH 19/20] Working on 9_3h. --- Makefile | 9 +- src/ex1.c | 5 +- src/get_proj_9_3h.c | 224 ++++++++++++++++++++++++++++++++++++++++++++ src/get_proj_9_3h.h | 8 ++ 4 files changed, 241 insertions(+), 5 deletions(-) create mode 100644 src/get_proj_9_3h.c create mode 100644 src/get_proj_9_3h.h diff --git a/Makefile b/Makefile index 49fb7c4..7ad10c8 100644 --- a/Makefile +++ b/Makefile @@ -45,7 +45,10 @@ ${OBJ_DIR}/get_s2_mov.o: ${SRC_DIR}/get_s2_mov.c directories chkopts ${OBJ_DIR}/get_dmat.o: ${SRC_DIR}/get_dmat.c directories chkopts ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} -${OBJ_DIR}/get_proj.o: ${SRC_DIR}/get_proj.c directories chkopts +${OBJ_DIR}/get_proj_9_3h.o: ${SRC_DIR}/get_proj_9_3h.c directories chkopts + ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} + +${OBJ_DIR}/get_proj_9_3h.o: ${SRC_DIR}/get_proj_9_3h.c directories chkopts ${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} ${OBJ_DIR}/get_val_iaa2.o: ${SRC_DIR}/get_val_iaa2.c directories chkopts @@ -54,6 +57,6 @@ ${OBJ_DIR}/get_val_iaa2.o: ${SRC_DIR}/get_val_iaa2.c directories chkopts ${OBJ_DIR}/ex1.o: ${SRC_DIR}/ex1.c -${CC} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -c -o $@ $< ${SLEPC_EPS_LIB} -${BIN_DIR}/ex1: ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${OBJ_DIR}/ex1.o ${SRC_DIR}/read2.h ${SRC_DIR}/get_ntot.h ${SRC_DIR}/stimsyr.h chkopts - -${CLINKER} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -o ${BIN_DIR}/ex1 ${OBJ_DIR}/ex1.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${SLEPC_EPS_LIB} +${BIN_DIR}/ex1: ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj_9_3h.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${OBJ_DIR}/ex1.o ${SRC_DIR}/read2.h ${SRC_DIR}/get_ntot.h ${SRC_DIR}/stimsyr.h chkopts + -${CLINKER} ${SLEPC_INCLUDE} ${PETSC_CC_INCLUDES} -o ${BIN_DIR}/ex1 ${OBJ_DIR}/ex1.o ${OBJ_DIR}/read2.o ${OBJ_DIR}/get_ntot.o ${OBJ_DIR}/get_s2.o ${OBJ_DIR}/get_s2_mov.o ${OBJ_DIR}/get_s2_cyclic.o ${OBJ_DIR}/get_dmat.o ${OBJ_DIR}/get_proj_9_3h.o ${OBJ_DIR}/get_val_iaa2.o ${LIB_DIR}/irpf90.a ${SLEPC_EPS_LIB} # ${RM} ex1.o read2.o diff --git a/src/ex1.c b/src/ex1.c index ac16b83..6eaa86f 100644 --- a/src/ex1.c +++ b/src/ex1.c @@ -5,7 +5,8 @@ #include "stimsyr.h" #include "get_s2.h" #include "get_ntot.h" -#include "get_proj.h" +//#include "get_proj.h" +#include "get_proj_9_3h.h" #undef __FUNCT__ #define __FUNCT__ "main" @@ -264,7 +265,7 @@ int main(int argc,char **argv) get_1rdm(values, &Istart, &Iend, &getdata.natom, &trace1rdm, natomax); //get_2rdm(values, &Istart, &Iend, &getdata.natom, &trace2rdm, densmat2, natomax); //get_2rdm(values, &Istart, &Iend, &getdata.natom, &trace2rdm, natomax); - get_proj(values, &Istart, &Iend, &getdata.natom, i, projvec, natomax); + get_proj_9_3h(values, &Istart, &Iend, &getdata.natom, i, projvec, natomax); weightproj = 0.0; for(ii=0;ii<6;++ii) weightproj += projvec[(i)*6 + ii]*projvec[(i)*6+ii]; // analyse_(valxr, (Iend-Istart), &Istart, &Iend, &xymat, &norm); diff --git a/src/get_proj_9_3h.c b/src/get_proj_9_3h.c new file mode 100644 index 0000000..89e571c --- /dev/null +++ b/src/get_proj_9_3h.c @@ -0,0 +1,224 @@ +#include +#include +#include +#include +#include +#include +#include "get_dmat.h" + + +/* + * + *---------------------------------------- + * Calculate the projection for 6_2h + *---------------------------------------- + * Input + * ===== + * valxr = The full vector + * Istart = Local starting id + * Iend = Local ending id + * Output + * ===== + * projvec + */ +void get_proj_9_3h(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, int iroot, double *projvec, const int natomax){ + + int ideter[natomax]; + int ideter1[natomax]; + int kko,kok,kkio; + int mmo,mom,mmio; + int p, q, r; + long int ii; + PetscInt iiii; + long int iii; + long int iaa2, iaa; + long int nrow=-1, ncol=-1; + double ms1=0.0, ms2=0.0 ,rest=0.0, norm=0.0; + double sq2; + double normvec[21]; + int sumMs=0; + int ntrouGauch = 0; + int ntrouMilieu = 0; + int ntrouDroit = 0; + double projvecmat[3*3*3*21]; + for(ii=0;ii<3*3*3*21;++ii) projvecmat[ii] = 0.0; + sq2 = sqrt(2.0); + + for(ii=*Istart;ii<*Iend;ii++) { + iii = ii + 1; + iiii = ii; + getdet_(&iii, ideter); + + ms1 = 0.0; + ms2 = 0.0; + p = -1; + q = -1; + // Find if Ms=5/2 + sumMs = 0; + ntrouGauch = 0; + ntrouMilieu = 0; + ntrouDroit = 0; + for(kko=0;kko<3;++kko){ + if(ideter[kko]==3){ + ntrouGauch += 1; + p = kko; + } + } + for(kko=3;kko<6;++kko){ + if(ideter[kko]==3){ + ntrouMilieu += 1; + q = kko - 3; + } + } + for(kko=6;kko<9;++kko){ + if(ideter[kko]==3){ + ntrouDroit += 1; + r = kko - 3; + } + } + if(ntrouGauch == 1 && ntrouMilieu == 1 && ntrouDroit == 1){ + sumMs = 1; + } + else{ + sumMs = 0; + } + // First box ms + for(kko=0;kko<=2;++kko){ + if(sumMs==1){ + if(ideter[kko]==1){ + ms1 = ms1 + 0.5; + } + else if(ideter[kko]==2){ + ms1 = ms1 - 0.5; + } + } + } + for(kok= 11;kok>=9;--kok){ + if(sumMs==1){ + if(ideter[kok]==1){ + ms1 = ms1 + 0.5; + } + else if(ideter[kok]==2){ + ms1 = ms1 - 0.5; + } + } + } + // Second box ms + for(kko=0;kko<=2;++kko){ + if(sumMs==1){ + if(ideter[kko]==1){ + ms2 = ms2 + 0.5; + } + else if(ideter[kko]==2){ + ms2 = ms2 - 0.5; + } + } + } + for(kok= 11;kok>=9;--kok){ + if(sumMs==1){ + if(ideter[kok]==1){ + ms2 = ms2 + 0.5; + } + else if(ideter[kok]==2){ + ms2 = ms2 - 0.5; + } + } + } + if(ideter[1] ==3 && ideter[4] == 3){ + norm += valxr[iiii]*valxr[iiii]; + } + + if(fabs(ms-2.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 0] = valxr[iiii]; + normvec[0] += 1.0; + } + else if(fabs(ms+2.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 1] = valxr[iiii]; + normvec[1] += 1.0; + } + else if(fabs(ms-1.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 2] += valxr[iiii]; + normvec[2] += 1.0; + } + else if(fabs(ms+1.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 3] += valxr[iiii]; + normvec[3] += 1.0; + } + else if(fabs(ms-0.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 4] += valxr[iiii]; + normvec[4] += 1.0; + } + else if(fabs(ms+0.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 5] += valxr[iiii]; + normvec[5] += 1.0; + } + else{ + rest+=valxr[iiii]; + } + + ms = 0.0; + + } + + for(p=0;p<3;++p){ + for(q=0;q<3;++q) { + projvecmat[p*3*6 + q*6 + 0] = projvecmat[p*3*6 + q*6 + 0]/sqrt(1.0); + projvecmat[p*3*6 + q*6 + 1] = projvecmat[p*3*6 + q*6 + 1]/sqrt(1.0); + projvecmat[p*3*6 + q*6 + 2] = projvecmat[p*3*6 + q*6 + 2]/sqrt(25.0); + projvecmat[p*3*6 + q*6 + 3] = projvecmat[p*3*6 + q*6 + 3]/sqrt(25.0); + projvecmat[p*3*6 + q*6 + 4] = projvecmat[p*3*6 + q*6 + 4]/sqrt(100.0); + projvecmat[p*3*6 + q*6 + 5] = projvecmat[p*3*6 + q*6 + 5]/sqrt(100.0); + } + } + //projvec[(iroot)*6 + 0] = projvec[(iroot)*6 + 0]/sqrt(1.0); + //projvec[(iroot)*6 + 1] = projvec[(iroot)*6 + 1]/sqrt(1.0); + //projvec[(iroot)*6 + 2] = projvec[(iroot)*6 + 2]/sqrt(25.0); + //projvec[(iroot)*6 + 3] = projvec[(iroot)*6 + 3]/sqrt(25.0); + //projvec[(iroot)*6 + 4] = projvec[(iroot)*6 + 4]/sqrt(100.0); + //projvec[(iroot)*6 + 5] = projvec[(iroot)*6 + 5]/sqrt(100.0); + //for(ii=0;ii<6;++ii){ + // projvec[(iroot)*6 + ii] = \ + // projvecmat[1*3*6 + 1*6 + ii]; + //} + for(ii=0;ii<6;++ii){ + projvecmat[0*3*6 + 0*6 + ii] *= 1.0/4.0; // (1,1) = 1.0/4.0 + projvecmat[0*3*6 + 1*6 + ii] *= sq2/4.0; // (1,2) = 1.4/4.0 + projvecmat[0*3*6 + 2*6 + ii] *= 1.0/4.0; // (1,3) = 1.0/4.0 + + projvecmat[1*3*6 + 0*6 + ii] *= sq2/4.0; // (2,1) = 1.4/4.0 + projvecmat[1*3*6 + 1*6 + ii] *= 2.0/4.0; // (2,2) = 2.0/4.0 + projvecmat[1*3*6 + 2*6 + ii] *= sq2/4.0; // (2,3) = 1.4/4.0 + + projvecmat[2*3*6 + 0*6 + ii] *= 1.0/4.0; // (3,1) = 1.0/4.0 + projvecmat[2*3*6 + 1*6 + ii] *= sq2/4.0; // (3,2) = 1.4/4.0 + projvecmat[2*3*6 + 2*6 + ii] *= 1.0/4.0; // (3,3) = 1.0/4.0 + } + + for(ii=0;ii<6;++ii){ + projvec[(iroot)*6 + ii] = \ + projvecmat[0*3*6 + 0*6 + ii] \ + + projvecmat[0*3*6 + 1*6 + ii] \ + + projvecmat[0*3*6 + 2*6 + ii] \ + + projvecmat[1*3*6 + 0*6 + ii] \ + + projvecmat[1*3*6 + 1*6 + ii] \ + + projvecmat[1*3*6 + 2*6 + ii] \ + + projvecmat[2*3*6 + 0*6 + ii] \ + + projvecmat[2*3*6 + 1*6 + ii] \ + + projvecmat[2*3*6 + 2*6 + ii]; + } + //projvec[(iroot)*6 + 0] = projvec[(iroot)*6 *6+ 0]*6/sqrt(9.0); + //projvec[(iroot)*6 + 1] = projvec[(iroot)*6 + 1]/sqrt(9.0); + //projvec[(iroot)*6 + 2] = projvec[(iroot)*6 + 2]/sqrt(25.0*9); + //projvec[(iroot)*6 + 3] = projvec[(iroot)*6 + 3]/sqrt(25.0*9); + //projvec[(iroot)*6 + 4] = projvec[(iroot)*6 + 4]/sqrt(100.0*9); + //projvec[(iroot)*6 + 5] = projvec[(iroot)*6 + 5]/sqrt(100.0*9); + //printf(" norm = %4.4f projnorm = %4.4f\n",norm,projvec[5]*projvec[5] +projvec[4]*projvec[4] +projvec[3]*projvec[3] +projvec[2]*projvec[2] +projvec[1]*projvec[1] +projvec[0]*projvec[0]); + //printf(" rest=%4.4f norm1=%2.2f norm2=%2.2f norm3=%2.2f norm4=%2.2f norm5=%2.2f norm6=%2.2f\n",rest,normvec[0],normvec[1],normvec[2],normvec[3],normvec[4],normvec[5]); + //printf(" proj ( 5/2,-5/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0], projvecseparate[1], projvecseparate[2], projvecseparate[3], projvecseparate[4], projvecseparate[5], projvecseparate[6], projvecseparate[7], projvecseparate[8]); + //printf(" proj (-5/2, 5/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+9], projvecseparate[1+9], projvecseparate[2+9], projvecseparate[3+9], projvecseparate[4+9], projvecseparate[5+9], projvecseparate[6+9], projvecseparate[7+9], projvecseparate[8+9]); + //printf(" proj ( 3/2,-3/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+18], projvecseparate[1+18], projvecseparate[2+18], projvecseparate[3+18], projvecseparate[4+18], projvecseparate[5+18], projvecseparate[6+18], projvecseparate[7+18], projvecseparate[8+18]); + //printf(" proj (-3/2, 3/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+27], projvecseparate[1+27], projvecseparate[2+27], projvecseparate[3+27], projvecseparate[4+27], projvecseparate[5+27], projvecseparate[6+27], projvecseparate[7+27], projvecseparate[8+27]); + //printf(" proj ( 1/2,-1/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+36], projvecseparate[1+36], projvecseparate[2+36], projvecseparate[3+36], projvecseparate[4+36], projvecseparate[5+36], projvecseparate[6+36], projvecseparate[7+36], projvecseparate[8+36]); + //printf(" proj (-1/2, 1/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+45], projvecseparate[1+45], projvecseparate[2+45], projvecseparate[3+45], projvecseparate[4+45], projvecseparate[5+45], projvecseparate[6+45], projvecseparate[7+45], projvecseparate[8+45]); +} /** END **/ + diff --git a/src/get_proj_9_3h.h b/src/get_proj_9_3h.h new file mode 100644 index 0000000..2612138 --- /dev/null +++ b/src/get_proj_9_3h.h @@ -0,0 +1,8 @@ +#include +#include +#include +#include +#include +#include + +void get_proj_9_3h(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, int iroot, double *projvec, const int natomax); From a18bfd2654840483811c649181a25de588b8c3ca Mon Sep 17 00:00:00 2001 From: v1j4y Date: Thu, 23 Feb 2023 18:30:58 +0100 Subject: [PATCH 20/20] Added files. --- src/get_proj.c | 193 +++++++++++++++++++++++++++++++++++++++++++++++++ src/get_proj.h | 8 ++ 2 files changed, 201 insertions(+) create mode 100644 src/get_proj.c create mode 100644 src/get_proj.h diff --git a/src/get_proj.c b/src/get_proj.c new file mode 100644 index 0000000..7985c95 --- /dev/null +++ b/src/get_proj.c @@ -0,0 +1,193 @@ +#include +#include +#include +#include +#include +#include +#include "get_dmat.h" + + +/* + * + *---------------------------------------- + * Calculate the projection for 6_2h + *---------------------------------------- + * Input + * ===== + * valxr = The full vector + * Istart = Local starting id + * Iend = Local ending id + * Output + * ===== + * projvec + */ +void get_proj(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, int iroot, double *projvec, const int natomax){ + + int ideter[natomax]; + int ideter2[natomax]; + int kko,kok,kkio; + int mmo,mom,mmio; + int p, q; + long int ii; + PetscInt iiii; + long int iii; + long int iaa2, iaa; + long int nrow=-1, ncol=-1; + double ms=0.0, rest=0.0, norm=0.0; + double sq2; + double normvec[6]; + int sumMs=0; + int ntrouGauch = 0; + int ntrouDroit = 0; + double projvecmat[3*3*6]; + for(ii=0;ii<3*3*6;++ii) projvecmat[ii] = 0.0; + sq2 = sqrt(2.0); + + for(ii=*Istart;ii<*Iend;ii++) { + iii = ii + 1; + iiii = ii; + getdet_(&iii, ideter); + + ms = 0.0; + p = -1; + q = -1; + // Find if Ms=5/2 + sumMs = 0; + ntrouGauch = 0; + ntrouDroit = 0; + for(kko=0;kko<3;++kko){ + if(ideter[kko]==3){ + ntrouGauch += 1; + p = kko; + } + } + for(kko=3;kko<6;++kko){ + if(ideter[kko]==3){ + ntrouDroit += 1; + q = kko - 3; + } + } + if(ntrouGauch == 1 && ntrouDroit == 1){ + sumMs = 1; + } + else{ + sumMs = 0; + } + for(kko=0;kko<=2;++kko){ + if(sumMs==1){ + if(ideter[kko]==1){ + ms = ms + 0.5; + } + else if(ideter[kko]==2){ + ms = ms - 0.5; + } + } + } + for(kok= 11;kok>=9;--kok){ + if(sumMs==1){ + if(ideter[kok]==1){ + ms = ms + 0.5; + } + else if(ideter[kok]==2){ + ms = ms - 0.5; + } + } + } + if(ideter[1] ==3 && ideter[4] == 3){ + norm += valxr[iiii]*valxr[iiii]; + } + + if(fabs(ms-2.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 0] = valxr[iiii]; + normvec[0] += 1.0; + } + else if(fabs(ms+2.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 1] = valxr[iiii]; + normvec[1] += 1.0; + } + else if(fabs(ms-1.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 2] += valxr[iiii]; + normvec[2] += 1.0; + } + else if(fabs(ms+1.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 3] += valxr[iiii]; + normvec[3] += 1.0; + } + else if(fabs(ms-0.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 4] += valxr[iiii]; + normvec[4] += 1.0; + } + else if(fabs(ms+0.5) < 1e-10){ + projvecmat[6*3*p + 6*q + 5] += valxr[iiii]; + normvec[5] += 1.0; + } + else{ + rest+=valxr[iiii]; + } + + ms = 0.0; + + } + + for(p=0;p<3;++p){ + for(q=0;q<3;++q) { + projvecmat[p*3*6 + q*6 + 0] = projvecmat[p*3*6 + q*6 + 0]/sqrt(1.0); + projvecmat[p*3*6 + q*6 + 1] = projvecmat[p*3*6 + q*6 + 1]/sqrt(1.0); + projvecmat[p*3*6 + q*6 + 2] = projvecmat[p*3*6 + q*6 + 2]/sqrt(25.0); + projvecmat[p*3*6 + q*6 + 3] = projvecmat[p*3*6 + q*6 + 3]/sqrt(25.0); + projvecmat[p*3*6 + q*6 + 4] = projvecmat[p*3*6 + q*6 + 4]/sqrt(100.0); + projvecmat[p*3*6 + q*6 + 5] = projvecmat[p*3*6 + q*6 + 5]/sqrt(100.0); + } + } + //projvec[(iroot)*6 + 0] = projvec[(iroot)*6 + 0]/sqrt(1.0); + //projvec[(iroot)*6 + 1] = projvec[(iroot)*6 + 1]/sqrt(1.0); + //projvec[(iroot)*6 + 2] = projvec[(iroot)*6 + 2]/sqrt(25.0); + //projvec[(iroot)*6 + 3] = projvec[(iroot)*6 + 3]/sqrt(25.0); + //projvec[(iroot)*6 + 4] = projvec[(iroot)*6 + 4]/sqrt(100.0); + //projvec[(iroot)*6 + 5] = projvec[(iroot)*6 + 5]/sqrt(100.0); + //for(ii=0;ii<6;++ii){ + // projvec[(iroot)*6 + ii] = \ + // projvecmat[1*3*6 + 1*6 + ii]; + //} + for(ii=0;ii<6;++ii){ + projvecmat[0*3*6 + 0*6 + ii] *= 1.0/4.0; // (1,1) = 1.0/4.0 + projvecmat[0*3*6 + 1*6 + ii] *= sq2/4.0; // (1,2) = 1.4/4.0 + projvecmat[0*3*6 + 2*6 + ii] *= 1.0/4.0; // (1,3) = 1.0/4.0 + + projvecmat[1*3*6 + 0*6 + ii] *= sq2/4.0; // (2,1) = 1.4/4.0 + projvecmat[1*3*6 + 1*6 + ii] *= 2.0/4.0; // (2,2) = 2.0/4.0 + projvecmat[1*3*6 + 2*6 + ii] *= sq2/4.0; // (2,3) = 1.4/4.0 + + projvecmat[2*3*6 + 0*6 + ii] *= 1.0/4.0; // (3,1) = 1.0/4.0 + projvecmat[2*3*6 + 1*6 + ii] *= sq2/4.0; // (3,2) = 1.4/4.0 + projvecmat[2*3*6 + 2*6 + ii] *= 1.0/4.0; // (3,3) = 1.0/4.0 + } + + for(ii=0;ii<6;++ii){ + projvec[(iroot)*6 + ii] = \ + projvecmat[0*3*6 + 0*6 + ii] \ + + projvecmat[0*3*6 + 1*6 + ii] \ + + projvecmat[0*3*6 + 2*6 + ii] \ + + projvecmat[1*3*6 + 0*6 + ii] \ + + projvecmat[1*3*6 + 1*6 + ii] \ + + projvecmat[1*3*6 + 2*6 + ii] \ + + projvecmat[2*3*6 + 0*6 + ii] \ + + projvecmat[2*3*6 + 1*6 + ii] \ + + projvecmat[2*3*6 + 2*6 + ii]; + } + //projvec[(iroot)*6 + 0] = projvec[(iroot)*6 *6+ 0]*6/sqrt(9.0); + //projvec[(iroot)*6 + 1] = projvec[(iroot)*6 + 1]/sqrt(9.0); + //projvec[(iroot)*6 + 2] = projvec[(iroot)*6 + 2]/sqrt(25.0*9); + //projvec[(iroot)*6 + 3] = projvec[(iroot)*6 + 3]/sqrt(25.0*9); + //projvec[(iroot)*6 + 4] = projvec[(iroot)*6 + 4]/sqrt(100.0*9); + //projvec[(iroot)*6 + 5] = projvec[(iroot)*6 + 5]/sqrt(100.0*9); + //printf(" norm = %4.4f projnorm = %4.4f\n",norm,projvec[5]*projvec[5] +projvec[4]*projvec[4] +projvec[3]*projvec[3] +projvec[2]*projvec[2] +projvec[1]*projvec[1] +projvec[0]*projvec[0]); + //printf(" rest=%4.4f norm1=%2.2f norm2=%2.2f norm3=%2.2f norm4=%2.2f norm5=%2.2f norm6=%2.2f\n",rest,normvec[0],normvec[1],normvec[2],normvec[3],normvec[4],normvec[5]); + //printf(" proj ( 5/2,-5/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0], projvecseparate[1], projvecseparate[2], projvecseparate[3], projvecseparate[4], projvecseparate[5], projvecseparate[6], projvecseparate[7], projvecseparate[8]); + //printf(" proj (-5/2, 5/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+9], projvecseparate[1+9], projvecseparate[2+9], projvecseparate[3+9], projvecseparate[4+9], projvecseparate[5+9], projvecseparate[6+9], projvecseparate[7+9], projvecseparate[8+9]); + //printf(" proj ( 3/2,-3/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+18], projvecseparate[1+18], projvecseparate[2+18], projvecseparate[3+18], projvecseparate[4+18], projvecseparate[5+18], projvecseparate[6+18], projvecseparate[7+18], projvecseparate[8+18]); + //printf(" proj (-3/2, 3/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+27], projvecseparate[1+27], projvecseparate[2+27], projvecseparate[3+27], projvecseparate[4+27], projvecseparate[5+27], projvecseparate[6+27], projvecseparate[7+27], projvecseparate[8+27]); + //printf(" proj ( 1/2,-1/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+36], projvecseparate[1+36], projvecseparate[2+36], projvecseparate[3+36], projvecseparate[4+36], projvecseparate[5+36], projvecseparate[6+36], projvecseparate[7+36], projvecseparate[8+36]); + //printf(" proj (-1/2, 1/2) : 1,4 %2.3f 1,5 %2.3f 1,6 %2.3f 2,4 %2.3f 2,5 %2.3f 2,6 %2.3f 3,4 %2.3f 3,5 %2.3f 3,6 %2.3f \n",projvecseparate[0+45], projvecseparate[1+45], projvecseparate[2+45], projvecseparate[3+45], projvecseparate[4+45], projvecseparate[5+45], projvecseparate[6+45], projvecseparate[7+45], projvecseparate[8+45]); +} /** END **/ + diff --git a/src/get_proj.h b/src/get_proj.h new file mode 100644 index 0000000..8b226b5 --- /dev/null +++ b/src/get_proj.h @@ -0,0 +1,8 @@ +#include +#include +#include +#include +#include +#include + +void get_proj(PetscScalar *valxr, PetscInt *Istart, PetscInt *Iend, int *natom, int iroot, double *projvec, const int natomax);