From c56c5c8ec231ac4e6e63b539c075dd0c3c69f83d Mon Sep 17 00:00:00 2001 From: philippun Date: Thu, 10 Aug 2023 23:49:22 +0200 Subject: [PATCH] working BTCS 1D and 2D simulation with constant boundary --- examples/BTCS_1D_proto_example.cpp | 48 +++++++ examples/CMakeLists.txt | 6 +- proto/BTCS.ipynb | 222 ++++++++++++++--------------- src/BTCSv2.cpp | 55 +++++-- 4 files changed, 204 insertions(+), 127 deletions(-) create mode 100644 examples/BTCS_1D_proto_example.cpp diff --git a/examples/BTCS_1D_proto_example.cpp b/examples/BTCS_1D_proto_example.cpp new file mode 100644 index 0000000..fb04966 --- /dev/null +++ b/examples/BTCS_1D_proto_example.cpp @@ -0,0 +1,48 @@ +#include + +int main(int argc, char *argv[]) { + // ************** + // **** GRID **** + // ************** + + // create a linear grid with 20 cells + int cells = 20; + Grid grid = Grid(cells); + + MatrixXd concentrations = MatrixXd::Constant(1,20,0); + concentrations(0,0) = 2000; + // TODO add option to set concentrations with a vector in 1D case + grid.setConcentrations(concentrations); + + + // ****************** + // **** BOUNDARY **** + // ****************** + + // create a boundary with constant values + Boundary bc = Boundary(grid); + bc.setBoundarySideConstant(BC_SIDE_LEFT, 0); + bc.setBoundarySideConstant(BC_SIDE_RIGHT, 0); + + + // ************************ + // **** SIMULATION ENV **** + // ************************ + + // set up a simulation environment + Simulation simulation = Simulation(grid, bc, BTCS_APPROACH); // grid,boundary,simulation-approach + + // set the timestep of the simulation + simulation.setTimestep(0.1); // timestep + + // set the number of iterations + simulation.setIterations(100); + + // set kind of output [CSV_OUTPUT_OFF (default), CSV_OUTPUT_ON, CSV_OUTPUT_VERBOSE] + simulation.setOutputCSV(CSV_OUTPUT_VERBOSE); + + // **** RUN SIMULATION **** + + // run the simulation + simulation.run(); +} \ No newline at end of file diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index 97efe1c..8997c0e 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -1,18 +1,20 @@ add_executable(first_example first_example.cpp) add_executable(second_example second_example.cpp) add_executable(boundary_example1D boundary_example1D.cpp) +add_executable(FTCS_1D_proto_example FTCS_1D_proto_example.cpp) add_executable(FTCS_2D_proto_example FTCS_2D_proto_example.cpp) +add_executable(BTCS_1D_proto_example BTCS_1D_proto_example.cpp) add_executable(BTCS_2D_proto_example BTCS_2D_proto_example.cpp) add_executable(FTCS_2D_proto_example_mdl FTCS_2D_proto_example_mdl.cpp) -add_executable(FTCS_1D_proto_example FTCS_1D_proto_example.cpp) add_executable(reference-FTCS_2D_closed reference-FTCS_2D_closed.cpp) target_link_libraries(first_example tug) target_link_libraries(second_example tug) target_link_libraries(boundary_example1D tug) +target_link_libraries(FTCS_1D_proto_example tug) target_link_libraries(FTCS_2D_proto_example tug) +target_link_libraries(BTCS_1D_proto_example tug) target_link_libraries(BTCS_2D_proto_example tug) target_link_libraries(FTCS_2D_proto_example_mdl tug) -target_link_libraries(FTCS_1D_proto_example tug) target_link_libraries(reference-FTCS_2D_closed tug) # target_link_libraries(FTCS_2D_proto_example easy_profiler) diff --git a/proto/BTCS.ipynb b/proto/BTCS.ipynb index 8ee2b79..81b6e45 100644 --- a/proto/BTCS.ipynb +++ b/proto/BTCS.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 40, + "execution_count": 2, "metadata": {}, "outputs": [], "source": [ @@ -14,7 +14,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -70,7 +70,7 @@ }, { "cell_type": "code", - "execution_count": 41, + "execution_count": 4, "metadata": {}, "outputs": [], "source": [ @@ -84,7 +84,7 @@ }, { "cell_type": "code", - "execution_count": 217, + "execution_count": 5, "metadata": {}, "outputs": [], "source": [ @@ -94,7 +94,7 @@ }, { "cell_type": "code", - "execution_count": 260, + "execution_count": 6, "metadata": {}, "outputs": [], "source": [ @@ -102,7 +102,7 @@ "### CONSTANT ###\n", "################\n", "\n", - "# creates the coeffiecient matrix for a given row\n", + "# creates the coefficient matrix for a given row\n", "def create_coeff_matrix(alpha_x, row_index, s_x):\n", " # alias\n", " r = row_index\n", @@ -231,7 +231,7 @@ "### CLOSED ###\n", "##############\n", "\n", - "# creates the coeffiecient matrix for a given row\n", + "# creates the coefficient matrix for a given row\n", "def create_coeff_matrix(alpha_x, row_index, s_x):\n", " # alias\n", " r = row_index\n", @@ -339,7 +339,7 @@ }, { "cell_type": "code", - "execution_count": 280, + "execution_count": 7, "metadata": {}, "outputs": [ { @@ -411,106 +411,106 @@ " [1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 1.]]\n", "\n", "Test Concentrations t=t+100:\n", - "[[6.28591428e+01 5.97567801e+01 5.40111841e+01 4.64274274e+01\n", - " 3.79691932e+01 2.95580781e+01 2.19168076e+01 1.54897130e+01\n", - " 1.04428201e+01 6.72160139e+00 4.13431473e+00 2.43233389e+00\n", - " 1.37011777e+00 7.39680705e-01 3.83120860e-01 1.90615726e-01\n", - " 9.12968457e-02 4.24222400e-02 1.98821745e-02 1.12119045e-02]\n", - " [5.97567801e+01 5.68075320e+01 5.13455053e+01 4.41360389e+01\n", - " 3.60952541e+01 2.80992628e+01 2.08351211e+01 1.47252306e+01\n", - " 9.92742306e+00 6.38986213e+00 3.93026893e+00 2.31228800e+00\n", - " 1.30249670e+00 7.03174355e-01 3.64212236e-01 1.81208039e-01\n", - " 8.67909629e-02 4.03285243e-02 1.89009056e-02 1.06585500e-02]\n", - " [5.40111841e+01 5.13455053e+01 4.64086509e+01 3.98923724e+01\n", - " 3.26247065e+01 2.53975274e+01 1.88318306e+01 1.33094043e+01\n", - " 8.97290440e+00 5.77547886e+00 3.55237479e+00 2.08996222e+00\n", - " 1.17726205e+00 6.35564358e-01 3.29193342e-01 1.63784942e-01\n", - " 7.84460386e-02 3.64509491e-02 1.70835894e-02 9.63373367e-03]\n", - " [4.64274274e+01 4.41360389e+01 3.98923724e+01 3.42910501e+01\n", - " 2.80438435e+01 2.18314387e+01 1.61876371e+01 1.14406194e+01\n", - " 7.71301119e+00 4.96453891e+00 3.05358279e+00 1.79650883e+00\n", - " 1.01196167e+00 5.46324221e-01 2.82971023e-01 1.40787758e-01\n", - " 6.74313630e-02 3.13328401e-02 1.46848680e-02 8.28105286e-03]\n", - " [3.79691932e+01 3.60952541e+01 3.26247065e+01 2.80438435e+01\n", - " 2.29347645e+01 1.78541471e+01 1.32385436e+01 9.35634629e+00\n", - " 6.30784062e+00 4.06009008e+00 2.49727546e+00 1.46921755e+00\n", - " 8.27600633e-01 4.46793869e-01 2.31418841e-01 1.15138785e-01\n", - " 5.51465932e-02 2.56245656e-02 1.20095517e-02 6.77239541e-03]\n", - " [2.95580781e+01 2.80992628e+01 2.53975274e+01 2.18314387e+01\n", - " 1.78541471e+01 1.38990120e+01 1.03058789e+01 7.28368424e+00\n", - " 4.91049795e+00 3.16067973e+00 1.94406720e+00 1.14374953e+00\n", - " 6.44266629e-01 3.47817980e-01 1.80153846e-01 8.96326973e-02\n", - " 4.29302593e-02 1.99480907e-02 9.34913903e-03 5.27214237e-03]\n", - " [2.19168076e+01 2.08351211e+01 1.88318306e+01 1.61876371e+01\n", - " 1.32385436e+01 1.03058789e+01 7.64163230e+00 5.40072685e+00\n", - " 3.64104995e+00 2.34358978e+00 1.44149246e+00 8.48070648e-01\n", - " 4.77712648e-01 2.57901062e-01 1.33580985e-01 6.64611067e-02\n", - " 3.18320504e-02 1.47911669e-02 6.93222614e-03 3.90920308e-03]\n", - " [1.54897130e+01 1.47252306e+01 1.33094043e+01 1.14406194e+01\n", - " 9.35634629e+00 7.28368424e+00 5.40072685e+00 3.81696598e+00\n", - " 2.57331359e+00 1.65633306e+00 1.01877541e+00 5.99374289e-01\n", - " 3.37623616e-01 1.82271685e-01 9.44084172e-02 4.69714152e-02\n", - " 2.24973151e-02 1.04536635e-02 4.89935375e-03 2.76283093e-03]\n", - " [1.04428201e+01 9.92742306e+00 8.97290440e+00 7.71301119e+00\n", - " 6.30784062e+00 4.91049795e+00 3.64104995e+00 2.57331359e+00\n", - " 1.73487080e+00 1.11666292e+00 6.86835725e-01 4.04084817e-01\n", - " 2.27618334e-01 1.22883517e-01 6.36480556e-02 3.16670837e-02\n", - " 1.51671896e-02 7.04762745e-03 3.30303535e-03 1.86263918e-03]\n", - " [6.72160139e+00 6.38986213e+00 5.77547886e+00 4.96453891e+00\n", - " 4.06009008e+00 3.16067973e+00 2.34358978e+00 1.65633306e+00\n", - " 1.11666292e+00 7.18748669e-01 4.42087091e-01 2.60092297e-01\n", - " 1.46508290e-01 7.90949197e-02 4.09675600e-02 2.03827617e-02\n", - " 9.76247810e-03 4.53625956e-03 2.12602408e-03 1.19890202e-03]\n", - " [4.13431473e+00 3.93026893e+00 3.55237479e+00 3.05358279e+00\n", - " 2.49727546e+00 1.94406720e+00 1.44149246e+00 1.01877541e+00\n", - " 6.86835725e-01 4.42087091e-01 2.71918412e-01 1.59977266e-01\n", - " 9.01141476e-02 4.86496108e-02 2.51982789e-02 1.25370053e-02\n", - " 6.00469362e-03 2.79015723e-03 1.30767241e-03 7.37419252e-04]\n", - " [2.43233389e+00 2.31228800e+00 2.08996222e+00 1.79650883e+00\n", - " 1.46921755e+00 1.14374953e+00 8.48070648e-01 5.99374289e-01\n", - " 4.04084817e-01 2.60092297e-01 1.59977266e-01 9.41191353e-02\n", - " 5.30166932e-02 2.86219374e-02 1.48248578e-02 7.37587358e-03\n", - " 3.53273052e-03 1.64152815e-03 7.69340539e-04 4.33844533e-04]\n", - " [1.37011777e+00 1.30249670e+00 1.17726205e+00 1.01196167e+00\n", - " 8.27600633e-01 6.44266629e-01 4.77712648e-01 3.37623616e-01\n", - " 2.27618334e-01 1.46508290e-01 9.01141476e-02 5.30166932e-02\n", - " 2.98639564e-02 1.61225501e-02 8.35074541e-03 4.15478134e-03\n", - " 1.98996399e-03 9.24662067e-04 4.33364493e-04 2.44381787e-04]\n", - " [7.39680705e-01 7.03174355e-01 6.35564358e-01 5.46324221e-01\n", - " 4.46793869e-01 3.47817980e-01 2.57901062e-01 1.82271685e-01\n", - " 1.22883517e-01 7.90949197e-02 4.86496108e-02 2.86219374e-02\n", - " 1.61225501e-02 8.70402492e-03 4.50828782e-03 2.24302732e-03\n", - " 1.07431492e-03 4.99194087e-04 2.33958978e-04 1.31933544e-04]\n", - " [3.83120860e-01 3.64212236e-01 3.29193342e-01 2.82971023e-01\n", - " 2.31418841e-01 1.80153846e-01 1.33580985e-01 9.44084172e-02\n", - " 6.36480556e-02 4.09675600e-02 2.51982789e-02 1.48248578e-02\n", - " 8.35074541e-03 4.50828782e-03 2.33508742e-03 1.16178582e-03\n", - " 5.56446119e-04 2.58559763e-04 1.21180077e-04 6.83355566e-05]\n", - " [1.90615726e-01 1.81208039e-01 1.63784942e-01 1.40787758e-01\n", - " 1.15138785e-01 8.96326973e-02 6.64611067e-02 4.69714152e-02\n", - " 3.16670837e-02 2.03827617e-02 1.25370053e-02 7.37587358e-03\n", - " 4.15478134e-03 2.24302732e-03 1.16178582e-03 5.78028166e-04\n", - " 2.76850968e-04 1.28642322e-04 6.02912313e-05 3.39992757e-05]\n", - " [9.12968457e-02 8.67909629e-02 7.84460386e-02 6.74313630e-02\n", - " 5.51465932e-02 4.29302593e-02 3.18320504e-02 2.24973151e-02\n", - " 1.51671896e-02 9.76247810e-03 6.00469362e-03 3.53273052e-03\n", - " 1.98996399e-03 1.07431492e-03 5.56446119e-04 2.76850968e-04\n", - " 1.32599868e-04 6.16142144e-05 2.88769419e-05 1.62842106e-05]\n", - " [4.24222400e-02 4.03285243e-02 3.64509491e-02 3.13328401e-02\n", - " 2.56245656e-02 1.99480907e-02 1.47911669e-02 1.04536635e-02\n", - " 7.04762745e-03 4.53625956e-03 2.79015723e-03 1.64152815e-03\n", - " 9.24662067e-04 4.99194087e-04 2.58559763e-04 1.28642322e-04\n", - " 6.16142144e-05 2.86298280e-05 1.34180382e-05 7.56666548e-06]\n", - " [1.98821745e-02 1.89009056e-02 1.70835894e-02 1.46848680e-02\n", - " 1.20095517e-02 9.34913903e-03 6.93222614e-03 4.89935375e-03\n", - " 3.30303535e-03 2.12602408e-03 1.30767241e-03 7.69340539e-04\n", - " 4.33364493e-04 2.33958978e-04 1.21180077e-04 6.02912313e-05\n", - " 2.88769419e-05 1.34180382e-05 6.28867725e-06 3.54629466e-06]\n", - " [1.12119045e-02 1.06585500e-02 9.63373367e-03 8.28105286e-03\n", - " 6.77239541e-03 5.27214237e-03 3.90920308e-03 2.76283093e-03\n", - " 1.86263918e-03 1.19890202e-03 7.37419252e-04 4.33844533e-04\n", - " 2.44381787e-04 1.31933544e-04 6.83355566e-05 3.39992757e-05\n", - " 1.62842106e-05 7.56666548e-06 3.54629466e-06 1.99981734e-06]]\n" + "[[1.03406102e-02 2.94500761e-02 4.42425339e-02 5.30261971e-02\n", + " 5.54585202e-02 5.24218517e-02 4.55845468e-02 3.68489169e-02\n", + " 2.78821774e-02 1.98459123e-02 1.33384447e-02 8.49104197e-03\n", + " 5.13286742e-03 2.95307410e-03 1.62013743e-03 8.48919901e-04\n", + " 4.24937081e-04 2.01925178e-04 8.71721193e-05 2.40321537e-05]\n", + " [2.94500761e-02 8.38738682e-02 1.26002814e-01 1.51018702e-01\n", + " 1.57945963e-01 1.49297526e-01 1.29824869e-01 1.04945780e-01\n", + " 7.94084902e-02 5.65211929e-02 3.79879138e-02 2.41825024e-02\n", + " 1.46184154e-02 8.41036020e-03 4.61415423e-03 2.41772537e-03\n", + " 1.21022155e-03 5.75083262e-04 2.48266349e-04 6.84436157e-05]\n", + " [4.42425339e-02 1.26002814e-01 1.89292678e-01 2.26873779e-01\n", + " 2.37280529e-01 2.24288074e-01 1.95034511e-01 1.57658921e-01\n", + " 1.19294525e-01 8.49111829e-02 5.70688361e-02 3.63291144e-02\n", + " 2.19610889e-02 1.26347941e-02 6.93179449e-03 3.63212293e-03\n", + " 1.81810288e-03 8.63941426e-04 3.72967877e-04 1.02822111e-04]\n", + " [5.30261971e-02 1.51018702e-01 2.26873779e-01 2.71916020e-01\n", + " 2.84388867e-01 2.68816963e-01 2.33755563e-01 1.88959635e-01\n", + " 1.42978587e-01 1.01768970e-01 6.83989610e-02 4.35416919e-02\n", + " 2.63211197e-02 1.51432348e-02 8.30799388e-03 4.35322413e-03\n", + " 2.17905877e-03 1.03546348e-03 4.47014816e-04 1.23235833e-04]\n", + " [5.54585202e-02 1.57945963e-01 2.37280529e-01 2.84388867e-01\n", + " 2.97433845e-01 2.81147656e-01 2.44477981e-01 1.97627254e-01\n", + " 1.49537045e-01 1.06437135e-01 7.15364362e-02 4.55389587e-02\n", + " 2.75284751e-02 1.58378583e-02 8.68908336e-03 4.55290746e-03\n", + " 2.27901267e-03 1.08296042e-03 4.67519482e-04 1.28888687e-04]\n", + " [5.24218517e-02 1.49297526e-01 2.24288074e-01 2.68816963e-01\n", + " 2.81147656e-01 2.65753227e-01 2.31091425e-01 1.86806041e-01\n", + " 1.41349044e-01 1.00609098e-01 6.76194106e-02 4.30454425e-02\n", + " 2.60211350e-02 1.49706457e-02 8.21330676e-03 4.30360996e-03\n", + " 2.15422380e-03 1.02366219e-03 4.41920139e-04 1.21831301e-04]\n", + " [4.55845468e-02 1.29824869e-01 1.95034511e-01 2.33755563e-01\n", + " 2.44477981e-01 2.31091425e-01 2.00950511e-01 1.62441204e-01\n", + " 1.22913096e-01 8.74868017e-02 5.87999104e-02 3.74310888e-02\n", + " 2.26272367e-02 1.30180465e-02 7.14205726e-03 3.74229645e-03\n", + " 1.87325156e-03 8.90147442e-04 3.84281146e-04 1.05941024e-04]\n", + " [3.68489169e-02 1.04945780e-01 1.57658921e-01 1.88959635e-01\n", + " 1.97627254e-01 1.86806041e-01 1.62441204e-01 1.31311658e-01\n", + " 9.93585501e-02 7.07212008e-02 4.75317441e-02 3.02579531e-02\n", + " 1.82910487e-02 1.05233230e-02 5.77338359e-03 3.02513857e-03\n", + " 1.51426955e-03 7.19563348e-04 3.10639131e-04 8.56389336e-05]\n", + " [2.78821774e-02 7.94084902e-02 1.19294525e-01 1.42978587e-01\n", + " 1.49537045e-01 1.41349044e-01 1.22913096e-01 9.93585501e-02\n", + " 7.51808454e-02 5.35120496e-02 3.59654675e-02 2.28950452e-02\n", + " 1.38401426e-02 7.96259931e-03 4.36850033e-03 2.28900759e-03\n", + " 1.14579032e-03 5.44466287e-04 2.35048846e-04 6.47997321e-05]\n", + " [1.98459123e-02 5.65211929e-02 8.49111829e-02 1.01768970e-01\n", + " 1.06437135e-01 1.00609098e-01 8.74868017e-02 7.07212008e-02\n", + " 5.35120496e-02 3.80886839e-02 2.55994179e-02 1.62961827e-02\n", + " 9.85110493e-03 5.66760067e-03 3.10940114e-03 1.62926457e-03\n", + " 8.15548004e-04 3.87538964e-04 1.67302528e-04 4.61230046e-05]\n", + " [1.33384447e-02 3.79879138e-02 5.70688361e-02 6.83989610e-02\n", + " 7.15364362e-02 6.76194106e-02 5.87999104e-02 4.75317441e-02\n", + " 3.59654675e-02 2.55994179e-02 1.72053779e-02 1.09526702e-02\n", + " 6.62093111e-03 3.80919642e-03 2.08982961e-03 1.09502929e-03\n", + " 5.48130103e-04 2.60465075e-04 1.12444088e-04 3.09992877e-05]\n", + " [8.49104197e-03 2.41825024e-02 3.63291144e-02 4.35416919e-02\n", + " 4.55389587e-02 4.30454425e-02 3.74310888e-02 3.02579531e-02\n", + " 2.28950452e-02 1.62961827e-02 1.09526702e-02 6.97229584e-03\n", + " 4.21477956e-03 2.42487392e-03 1.33035233e-03 6.97078252e-04\n", + " 3.48930915e-04 1.65807929e-04 7.15801202e-05 1.97336541e-05]\n", + " [5.13286742e-03 1.46184154e-02 2.19610889e-02 2.63211197e-02\n", + " 2.75284751e-02 2.60211350e-02 2.26272367e-02 1.82910487e-02\n", + " 1.38401426e-02 9.85110493e-03 6.62093111e-03 4.21477956e-03\n", + " 2.54785040e-03 1.46584558e-03 8.04203081e-04 4.21386475e-04\n", + " 2.10930076e-04 1.00231528e-04 4.32704570e-05 1.19290695e-05]\n", + " [2.95307410e-03 8.41036020e-03 1.26347941e-02 1.51432348e-02\n", + " 1.58378583e-02 1.49706457e-02 1.30180465e-02 1.05233230e-02\n", + " 7.96259931e-03 5.66760067e-03 3.80919642e-03 2.42487392e-03\n", + " 1.46584558e-03 8.43339650e-04 4.62679258e-04 2.42434761e-04\n", + " 1.21353640e-04 5.76658437e-05 2.48946360e-05 6.86310854e-06]\n", + " [1.62013743e-03 4.61415423e-03 6.93179449e-03 8.30799388e-03\n", + " 8.68908336e-03 8.21330676e-03 7.14205726e-03 5.77338359e-03\n", + " 4.36850033e-03 3.10940114e-03 2.08982961e-03 1.33035233e-03\n", + " 8.04203081e-04 4.62679258e-04 2.53838528e-04 1.33006358e-04\n", + " 6.65779344e-05 3.16370632e-05 1.36578799e-05 3.76528954e-06]\n", + " [8.48919901e-04 2.41772537e-03 3.63212293e-03 4.35322413e-03\n", + " 4.55290746e-03 4.30360996e-03 3.74229645e-03 3.02513857e-03\n", + " 2.28900759e-03 1.62926457e-03 1.09502929e-03 6.97078252e-04\n", + " 4.21386475e-04 2.42434761e-04 1.33006358e-04 6.96926953e-05\n", + " 3.48855181e-05 1.65771940e-05 7.15645839e-06 1.97293709e-06]\n", + " [4.24937081e-04 1.21022155e-03 1.81810288e-03 2.17905877e-03\n", + " 2.27901267e-03 2.15422380e-03 1.87325156e-03 1.51426955e-03\n", + " 1.14579032e-03 8.15548004e-04 5.48130103e-04 3.48930915e-04\n", + " 2.10930076e-04 1.21353640e-04 6.65779344e-05 3.48855181e-05\n", + " 1.74623662e-05 8.29791415e-06 3.58225145e-06 9.87577426e-07]\n", + " [2.01925178e-04 5.75083262e-04 8.63941426e-04 1.03546348e-03\n", + " 1.08296042e-03 1.02366219e-03 8.90147442e-04 7.19563348e-04\n", + " 5.44466287e-04 3.87538964e-04 2.60465075e-04 1.65807929e-04\n", + " 1.00231528e-04 5.76658437e-05 3.16370632e-05 1.65771940e-05\n", + " 8.29791415e-06 3.94307268e-06 1.70224439e-06 4.69285352e-07]\n", + " [8.71721193e-05 2.48266349e-04 3.72967877e-04 4.47014816e-04\n", + " 4.67519482e-04 4.41920139e-04 3.84281146e-04 3.10639131e-04\n", + " 2.35048846e-04 1.67302528e-04 1.12444088e-04 7.15801202e-05\n", + " 4.32704570e-05 2.48946360e-05 1.36578799e-05 7.15645839e-06\n", + " 3.58225145e-06 1.70224439e-06 7.34867500e-07 2.02592857e-07]\n", + " [2.40321537e-05 6.84436157e-05 1.02822111e-04 1.23235833e-04\n", + " 1.28888687e-04 1.21831301e-04 1.05941024e-04 8.56389336e-05\n", + " 6.47997321e-05 4.61230046e-05 3.09992877e-05 1.97336541e-05\n", + " 1.19290695e-05 6.86310854e-06 3.76528954e-06 1.97293709e-06\n", + " 9.87577426e-07 4.69285352e-07 2.02592857e-07 5.58520626e-08]]\n" ] } ], @@ -540,19 +540,19 @@ "print(\"\\nTest Alpha:\")\n", "print(testAlpha)\n", "\n", - "testConcentrations = sim_closed(testConcentrations, testAlpha, testAlpha, iterations, delta_t, delta_x, delta_y, bc_left, bc_right, bc_top, bc_bottom)\n", + "testConcentrations = sim_constant(testConcentrations, testAlpha, testAlpha, iterations, delta_t, delta_x, delta_y, bc_left, bc_right, bc_top, bc_bottom)\n", "print(\"\\nTest Concentrations t=t+\" + str(iterations) + \":\")\n", "print(testConcentrations)" ] }, { "cell_type": "code", - "execution_count": 281, + "execution_count": 8, "metadata": {}, "outputs": [ { "data": { - "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfUAAAGdCAYAAAAL7+omAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA3wUlEQVR4nO3de3xU5Z3H8e8JkAlgEouFTCIBIxtBLqUUkBCVi5TQWKldtKBYLvVSXdCK1Bc1Ul+N3W6iVGlWUVxd5LIWsV0u8lqsEl6SIAUsSGIpUsQ1hVQz5gULCbdc59k/aEbHXGAyZ8J4zufd1/Oqc87zPPObk9HfPM95zjmWMcYIAAB85cVc7AAAAIA9SOoAADgESR0AAIcgqQMA4BAkdQAAHIKkDgCAQ5DUAQBwCJI6AAAO0fliB2AXv9+vTz/9VPHx8bIs62KHAwAIgTFGJ0+eVEpKimJiIjferKmpUV1dnS19xcbGKi4uzpa+7OKYpP7pp58qNTX1YocBAAhDeXm5evfuHZG+a2pqlNb3EvkqG23pz+v1qqysLKoSu2OSenx8vCTp8N4rlHBJeL/yyhtO2RGSJGn72T629LPl/wbZ0o8klZTb9C9MeTd7+pHU7e/2zK5cUtFgSz+S1PWzGlv66fR/p23pR5J00p6+zNmztvQjSaau3p6OGu35D60kGbv64i7aHaZB9dquNwL/LY+Euro6+SobVfZeXyXEh5cnqk/6lTb8sOrq6kjqkdA05Z5wSUzYf6z4Bvumfrp2sucQd6mNtaUfSYrpZtMX0MYvcqdYe5J65y72JfXONv3b0amTfTEpxp4EaiwbE6hdZ7ss+/69M7b1RVLvMP841B1x+jQhPvw8Ea0ck9QBALgQjcavxjB/rzUavz3B2IykDgBwFb+M/GHOwoTbPlIiNv/w/PPPKy0tTXFxcRo+fLjeeeedNusXFxdr+PDhiouL05VXXqkXXnghUqEBAFzMb9P/olFEkvprr72mefPmaeHChSopKdH111+v7OxsHTlypMX6ZWVluvHGG3X99derpKREjz76qH7yk59o7dq1kQgPAABHikhSX7x4se666y7dfffduvrqq1VQUKDU1FQtXbq0xfovvPCC+vTpo4KCAl199dW6++67deedd+qpp56KRHgAABdrNMaWEo1sT+p1dXV67733lJWVFbQ9KytLO3bsaLHNzp07m9WfNGmS9uzZo/r6llf71tbWqrq6OqgAAHA+TefUwy3RyPakfvToUTU2NiopKSloe1JSknw+X4ttfD5fi/UbGhp09OjRFtvk5+crMTExULjxDADA7SK2UO7L1xoaY9q8/rCl+i1tb5KTk6OqqqpAKS8vDzNiAIAb+GXUGGaJ1pG67Ze0ff3rX1enTp2ajcorKyubjcabeL3eFut37txZl112WYttPB6PPB6PPUEDAFyDS9pCEBsbq+HDh6uwsDBoe2FhoTIzM1tsM3r06Gb1N2/erBEjRqhLly52hwgAQIf65JNP9MMf/lCXXXaZunXrpm9+85t67733AvuNMcrNzVVKSoq6du2qcePGaf/+/SG/T0Sm3+fPn6///M//1Msvv6wDBw7ooYce0pEjR3TfffdJOjd1PnPmzED9++67T4cPH9b8+fN14MABvfzyy1q2bJkefvjhSIQHAHCxjl79fvz4cV177bXq0qWL/vCHP+iDDz7Q008/rUsvvTRQZ9GiRVq8eLGWLFmi3bt3y+v1auLEiTp58mRIny0id5SbNm2ajh07pl/+8peqqKjQ4MGD9cYbb6hv376SpIqKiqBr1tPS0vTGG2/ooYce0nPPPaeUlBQ988wzuuWWWyIRHgDAxfz/KOH2caGefPJJpaamavny5YFtV1xxReCfjTEqKCjQwoULNWXKFEnSypUrlZSUpNWrV+vee++94PeK2G1i58yZozlz5rS4b8WKFc22jR07Vnv37o1UOAAA2O7Ll1O3tN5r48aNmjRpkn7wgx+ouLhYl19+uebMmaN77rlH0rkbsPl8vqBLuz0ej8aOHasdO3aElNSd+ZgaAABaEe7K96YiSampqUGXV+fn5zd7v48//lhLly5Venq63nrrLd133336yU9+olWrVklSYKF4KJeCt8ZxD3QpbzgV9qNT/9ZwiU3RSBX1X7Oln2M13W3pR5Ia6jrZ0k9snS3dSJI61duzkjTGpn4kyQr3MU5NovTOU9HG+DlO6BiNRjY8pe3c/5eXlyshISGwvaWrsvx+v0aMGKG8vDxJ0rBhw7R//34tXbo0aH1ZqJeCt4SROgDAVfw2FUlKSEgIKi0l9eTkZA0cODBo29VXXx1YW+b1eiUppEvBW0NSBwAggq699lodPHgwaNuHH34YWDyelpYmr9cbdGl3XV2diouLW70UvDWOm34HAKAtfllqVGjT2i31caEeeughZWZmKi8vT1OnTtWf/vQnvfjii3rxxRclnZt2nzdvnvLy8pSenq709HTl5eWpW7dumj59ekhxkdQBAK7iN+dKuH1cqJEjR2r9+vXKycnRL3/5S6WlpamgoEB33HFHoM6CBQt09uxZzZkzR8ePH9eoUaO0efNmxcfHhxQXSR0AgAi76aabdNNNN7W637Is5ebmKjc3N6z3IakDAFyl0Ybp93DbRwpJHQDgKk5O6qx+BwDAIRipAwBcxW8s+U2Yq9/DbB8pJHUAgKsw/Q4AAKIeI3UAgKs0KkaNYY5pG22KxW4kdQCAqxgbzqkbzqkDAHDxcU4dAABEPUbqAABXaTQxajRhnlMP897xkUJSBwC4il+W/GFOVPsVnVmd6XcAABzCcSP17Wf7qGun8D5WRf3XbIpGer+6ty39fFqdYEs/kqTqLrZ00/mMfQtFOtXY86s3pt6+X89WvT0XrVgN9l38Yhpt6svvt6cfKfxnWAIdzMkL5RyX1AEAaIs959Sj88cs0+8AADgEI3UAgKucWygX5gNdmH4HAODi89twm1hWvwMAgIhipA4AcBUnL5QjqQMAXMWvGMfefIakDgBwlUZjqTHMp6yF2z5SOKcOAIBDMFIHALhKow2r3xuZfgcA4OLzmxj5w1wo54/ShXJMvwMA4BCM1AEArsL0OwAADuFX+KvXbXzOoa1sn37Pz8/XyJEjFR8fr169eun73/++Dh482GaboqIiWZbVrPz1r3+1OzwAABzL9pF6cXGx5s6dq5EjR6qhoUELFy5UVlaWPvjgA3Xv3r3NtgcPHlRCwufPDe/Zs6fd4QEAXM6em89E55I025P6m2++GfR6+fLl6tWrl9577z2NGTOmzba9evXSpZdeandIAAAE2HObWJck9S+rqqqSJPXo0eO8dYcNG6aamhoNHDhQP//5zzV+/PhW69bW1qq2tjbwurq6WpK05f8GqUttbFgxH6tpe0YhFJ9WJ5y/0gU4VWlfTJ7j9nwZY6vsWygSe9qeM1SdzjTY0o8kWXU29dXQaE8/ktRo05m8KL0cB0B4IvpTwxij+fPn67rrrtPgwYNbrZecnKwXX3xRa9eu1bp169S/f39NmDBB27Zta7VNfn6+EhMTAyU1NTUSHwEA4DBNz1MPt0SjiI7U77//fv35z3/W9u3b26zXv39/9e/fP/B69OjRKi8v11NPPdXqlH1OTo7mz58feF1dXU1iBwCcF9Pv7fDAAw9o48aN2rZtm3r37h1y+4yMDL3yyiut7vd4PPJ4POGECABwIXuuU3dJUjfG6IEHHtD69etVVFSktLS0dvVTUlKi5ORkm6MDAMC5bE/qc+fO1erVq/X6668rPj5ePp9PkpSYmKiuXbtKOjd1/sknn2jVqlWSpIKCAl1xxRUaNGiQ6urq9Morr2jt2rVau3at3eEBAFzObyz5w735TJQ+etX2pL506VJJ0rhx44K2L1++XLNnz5YkVVRU6MiRI4F9dXV1evjhh/XJJ5+oa9euGjRokDZt2qQbb7zR7vAAAC7nt2H63TXXqZsLuFRmxYoVQa8XLFigBQsW2B0KAACuwr3fAQCuYs+jV10yUgcAIJo1ylJjmNeZh9s+UqLzpwYAAAgZI3UAgKsw/Q4AgEM0Kvzpcxuf6GCr6PypAQAAQsZIHQDgKky/AwDgEDzQBQAAhzA2PDrVcEkbAADuk5ubK8uygorX6w3sN8YoNzdXKSkp6tq1q8aNG6f9+/e3671I6gAAV2mafg+3hGLQoEGqqKgIlH379gX2LVq0SIsXL9aSJUu0e/dueb1eTZw4USdPngz5szlu+r2kvLdiusWF1UdDXSebopFU3cWWbjzH7fv9FXfUnmkjzwm/Lf1IUpfqBlv66XS23pZ+JMmqqbOlH1NvX0xqtOdCGmNTP+c6s+l7YFc/knQBz6CAe12Mp7R17tw5aHTexBijgoICLVy4UFOmTJEkrVy5UklJSVq9erXuvffekN6HkToAAO1UXV0dVGpra1usd+jQIaWkpCgtLU233XabPv74Y0lSWVmZfD6fsrKyAnU9Ho/Gjh2rHTt2hBwPSR0A4CqN/3j0arhFklJTU5WYmBgo+fn5zd5v1KhRWrVqld566y299NJL8vl8yszM1LFjx+Tz+SRJSUlJQW2SkpIC+0LhuOl3AADaYuf0e3l5uRISEgLbPR5Ps7rZ2dmBfx4yZIhGjx6tfv36aeXKlcrIyJAkWVZwPMaYZtsuBCN1AADaKSEhIai0lNS/rHv37hoyZIgOHToUOM/+5VF5ZWVls9H7hSCpAwBcxa8YW0p71dbW6sCBA0pOTlZaWpq8Xq8KCwsD++vq6lRcXKzMzMyQ+2b6HQDgKo3GUmOY0++htH/44Yc1efJk9enTR5WVlfrVr36l6upqzZo1S5Zlad68ecrLy1N6errS09OVl5enbt26afr06SHHRVIHACCC/v73v+v222/X0aNH1bNnT2VkZGjXrl3q27evJGnBggU6e/as5syZo+PHj2vUqFHavHmz4uPjQ34vkjoAwFU6+jr1NWvWtLnfsizl5uYqNzc3rJgkkjoAwGWMDU9pMzzQBQCAi69RlhrDfCBLuO0jJTp/agAAgJAxUgcAuIrfhH7v9pb6iEYkdQCAq/htOKcebvtIic6oAABAyBipAwBcxS9L/jAXuoXbPlJI6gAAV+noO8p1JKbfAQBwCOeN1Mu7SXFxYXURW2dTLJI6n7Hn11xslX1LLT0n/Lb0E3e80ZZ+JKlLtT0HPeZUrS39SJJqbfoi1NXb048k09BgT0c2Lt010boMGGiFkxfKOS+pAwDQBr9suE1slJ5Tj86fGgAAIGSM1AEArmJsWP1uonSkTlIHALhKRz+lrSOR1AEAruLkhXK2R5WbmyvLsoKK1+tts01xcbGGDx+uuLg4XXnllXrhhRfsDgsAAMeLyEh90KBB2rJlS+B1p06dWq1bVlamG2+8Uffcc49eeeUV/fGPf9ScOXPUs2dP3XLLLZEIDwDgYky/h9pp587nHZ03eeGFF9SnTx8VFBRIkq6++mrt2bNHTz31FEkdAGA7J98mNiInBQ4dOqSUlBSlpaXptttu08cff9xq3Z07dyorKyto26RJk7Rnzx7V17d+047a2lpVV1cHFQAA3Mz2pD5q1CitWrVKb731ll566SX5fD5lZmbq2LFjLdb3+XxKSkoK2paUlKSGhgYdPXq01ffJz89XYmJioKSmptr6OQAAztQ0/R5uiUa2J/Xs7GzdcsstGjJkiL797W9r06ZNkqSVK1e22saygg+OMabF7V+Uk5OjqqqqQCkvL7chegCA0zk5qUf8krbu3btryJAhOnToUIv7vV6vfD5f0LbKykp17txZl112Wav9ejweeTweW2MFAOCrLOIX2tXW1urAgQNKTk5ucf/o0aNVWFgYtG3z5s0aMWKEunTpEunwAAAu4+SRuu1J/eGHH1ZxcbHKysr07rvv6tZbb1V1dbVmzZol6dy0+cyZMwP177vvPh0+fFjz58/XgQMH9PLLL2vZsmV6+OGH7Q4NAABHJ3Xbp9///ve/6/bbb9fRo0fVs2dPZWRkaNeuXerbt68kqaKiQkeOHAnUT0tL0xtvvKGHHnpIzz33nFJSUvTMM89wORsAACGyPamvWbOmzf0rVqxotm3s2LHau3ev3aEAANCMUfjXmRt7QrEd934HALgKd5QDAMAhSOpfId3+bqlTbHgHu1O9fRMrnWrs6Sv2tN+WfiSpS3WDTf3U2dKPJHU6VWtLP9aZGlv6kSRTa09MpsGe4y1JprHRpo7s+z7Z2heAsDguqQMA0BZG6gAAOISTk3p0PuUdAACEjJE6AMBVjLFkwhxph9s+UkjqAABX4XnqAAAg6jFSBwC4ipMXypHUAQCu4uRz6ky/AwDgEIzUAQCuwvQ7AAAO4eTpd5I6AMBVjA0j9WhN6pxTBwDAIRipAwBcxUgyYT5A075nedqLpA4AcBW/LFncUQ4AAEQzRuoAAFdh9ftXyCUVDercpSGsPmLq7TtbYldfnc6E95mC+jpbb0s/MadqbelHkqwzNbb0Y2rsi0l19hwn1dvUjyT57fk+GZv6sVW4JzmBC+Q3liyHXqfO9DsAAA5BUgcAuIox9pT2ys/Pl2VZmjdv3hdiMsrNzVVKSoq6du2qcePGaf/+/SH3TVIHALhK0zn1cEt77N69Wy+++KK+8Y1vBG1ftGiRFi9erCVLlmj37t3yer2aOHGiTp48GVL/JHUAADrAqVOndMcdd+ill17S1772tcB2Y4wKCgq0cOFCTZkyRYMHD9bKlSt15swZrV69OqT3IKkDAFzFzpF6dXV1UKmtbX2x7ty5c/Xd735X3/72t4O2l5WVyefzKSsrK7DN4/Fo7Nix2rFjR0ifjaQOAHCVpqe0hVskKTU1VYmJiYGSn5/f4nuuWbNGe/fubXG/z+eTJCUlJQVtT0pKCuy7UI67pA0AgLaEu9CtqQ9JKi8vV0JCQmC7x+NpVre8vFwPPvigNm/erLi4uFb7tKzg8/TGmGbbzoekDgBAOyUkJAQl9Za89957qqys1PDhwwPbGhsbtW3bNi1ZskQHDx6UdG7EnpycHKhTWVnZbPR+Pky/AwBc5dxIPdxz6hf+fhMmTNC+fftUWloaKCNGjNAdd9yh0tJSXXnllfJ6vSosLAy0qaurU3FxsTIzM0P6bIzUAQCu0tG3iY2Pj9fgwYODtnXv3l2XXXZZYPu8efOUl5en9PR0paenKy8vT926ddP06dNDioukDgDARbZgwQKdPXtWc+bM0fHjxzVq1Cht3rxZ8fHxIfVDUgcAuIpR+M9DD7d9UVFR0GvLspSbm6vc3Nyw+iWpAwBcxclPaWOhHAAADmF7Ur/iiitkWVazMnfu3BbrFxUVtVj/r3/9q92hAQDw+fx7uCUK2T79vnv3bjU2NgZe/+Uvf9HEiRP1gx/8oM12Bw8eDLrWr2fPnnaHBgCAZMP0u6J0+t32pP7lZPzEE0+oX79+Gjt2bJvtevXqpUsvvdTucAAACGLnHeWiTUTPqdfV1emVV17RnXfeed5b3Q0bNkzJycmaMGGCtm7dGsmwAABwpIiuft+wYYNOnDih2bNnt1onOTlZL774ooYPH67a2lr913/9lyZMmKCioiKNGTOm1Xa1tbVBT8Oprq6WJHX9rEadw/xUVqN9P8Gs+sbzV7qQfuoabOlHkqyaOns6qrWpH0mmjScbhaSu3p5+JJkGe465afTb0s+5vuz5PsnYF1PUDlmAVjh59XtEk/qyZcuUnZ2tlJSUVuv0799f/fv3D7wePXq0ysvL9dRTT7WZ1PPz8/X444/bGi8AwAWMFf458ShN6hGbfj98+LC2bNmiu+++O+S2GRkZOnToUJt1cnJyVFVVFSjl5eXtDRUAAEeI2Eh9+fLl6tWrl7773e+G3LakpCToSTUt8Xg8LT7iDgCAtjh5oVxEkrrf79fy5cs1a9Ysdf7SCe6cnBx98sknWrVqlSSpoKBAV1xxhQYNGhRYWLd27VqtXbs2EqEBANwuGu4TGyERSepbtmzRkSNHdOeddzbbV1FRoSNHjgRe19XV6eGHH9Ynn3yirl27atCgQdq0aZNuvPHGSIQGAIBjRSSpZ2VlybQyN7FixYqg1wsWLNCCBQsiEQYAAM2w+h0AACeJ0unzcPFAFwAAHIKROgDAVZh+BwDAKVj9DgCAU1j/KOH2EX04pw4AgEMwUgcAuAvT7wAAOISDkzrT7wAAOAQjdQCAuzj40askdQCAq/CUtq+QTv93Wp06NYTXiY1/Lauh0Z6O7OpHkqmvt6ejOpv6kWQawvybNbHrs0kyjX6b+rHvbydjT0wAnMlxSR0AgDY5eKEcSR0A4C4OPqfO6ncAAByCkToAwFUsc66E20c0IqkDANyFc+oAADgE59QBAEC0Y6QOAHAXpt8BAHAIByd1pt8BAHAIRuoAAHdx8EidpA4AcBdWvwMAgGjHSB0A4CrcUQ4AAKdw8Dl1pt8BAHAIkjoAAA7B9DsAwFUs2XBO3ZZI7Oe8pH7ytBRTf7GjCDCNjfZ01Oi3px9Jsikm09BgSz+SjcfJb9+JLttiMjb+7exiovSEINARuKQNAABEO+eN1AEAaIuDV7+T1AEA7uLgpM70OwAADkFSBwC4StMd5cItF2rp0qX6xje+oYSEBCUkJGj06NH6wx/+ENhvjFFubq5SUlLUtWtXjRs3Tvv372/XZws5qW/btk2TJ09WSkqKLMvShg0bgva3N7i1a9dq4MCB8ng8GjhwoNavXx9qaAAAnJ+xqVyg3r1764knntCePXu0Z88e3XDDDbr55psDuXHRokVavHixlixZot27d8vr9WrixIk6efJkyB8t5KR++vRpDR06VEuWLGlxf3uC27lzp6ZNm6YZM2bo/fff14wZMzR16lS9++67oYYHAEBUmTx5sm688UZdddVVuuqqq/Rv//ZvuuSSS7Rr1y4ZY1RQUKCFCxdqypQpGjx4sFauXKkzZ85o9erVIb9XyEk9Oztbv/rVrzRlypRm+9obXEFBgSZOnKicnBwNGDBAOTk5mjBhggoKCkINDwCAttk4Uq+urg4qtbW1bb51Y2Oj1qxZo9OnT2v06NEqKyuTz+dTVlZWoI7H49HYsWO1Y8eOkD+arefU2xvczp07g9pI0qRJk9psU1tb2+xgAgBwPnaeU09NTVViYmKg5Ofnt/ie+/bt0yWXXCKPx6P77rtP69ev18CBA+Xz+SRJSUlJQfWTkpIC+0Jh6yVtbQV3+PDhNtuF+oHy8/P1+OOPhxEtAADhKS8vV0JCQuC1x+NpsV7//v1VWlqqEydOaO3atZo1a5aKi4sD+y0r+A51xphm2y5ERFa/tye4UNvk5OSoqqoqUMrLy9sfMADAPZpuExtukQIr2ptKa0k9NjZW//RP/6QRI0YoPz9fQ4cO1b//+7/L6/VKUrNBbGVlZbPB7oWwNam3Nziv1xtyG4/H0+xgAgBwXh28+r3FEIxRbW2t0tLS5PV6VVhYGNhXV1en4uJiZWZmhtyvrUm9vcGNHj06qI0kbd68uV0fCACAtnT0deqPPvqo3nnnHf3tb3/Tvn37tHDhQhUVFemOO+6QZVmaN2+e8vLytH79ev3lL3/R7Nmz1a1bN02fPj3kzxbyOfVTp07po48+CrwuKytTaWmpevTooT59+gSCS09PV3p6uvLy8poFN3PmTF1++eWBBQUPPvigxowZoyeffFI333yzXn/9dW3ZskXbt28P+QMBABBNPvvsM82YMUMVFRVKTEzUN77xDb355puaOHGiJGnBggU6e/as5syZo+PHj2vUqFHavHmz4uPjQ36vkJP6nj17NH78+MDr+fPnS5JmzZqlFStWXFBwR44cUUzM55MEmZmZWrNmjX7+85/rscceU79+/fTaa69p1KhRIX8gAADa1MH3fl+2bFmb+y3LUm5urnJzc8OLSZJljDMerFxdXa3ExERN6HW3OsfEXuxwPsfz1C+sL56n3rGc8a89HKTB1KtIr6uqqipia6Sa8sSVj+WpU1xcWH011tTo4399NKLxtgf3fgcAwCEc9+hVc/asjGXTCMsOfptGaTaOrKJxVGzXaNZEYUy2YoQNhM/Bj151XFIHAKBNDk7qTL8DAOAQjNQBAK4S6nXmrfURjRipAwDgECR1AAAcgul3AIC7OHihHEkdAOAqTj6nTlIHALhPlCblcHFOHQAAh2CkDgBwF86pAwDgDE4+p870OwAADsFIHQDgLky/AwDgDEy/AwCAqMdIHQDgLky/AwDgEA5O6ky/AwDgEI4bqZu6ehnrYkfxBf4o/Dln/DZ1Y+NnsykmW5ko/NsBCJuTF8o5LqkDANAmB0+/k9QBAO7i4KTOOXUAAByCkToAwFU4pw4AgFMw/Q4AAKIdI3UAgKsw/Q4AgFMw/Q4AAKIdI3UAgLs4eKROUgcAuIr1jxJuH9GI6XcAAByCkToAwF2YfgcAwBmcfElbyNPv27Zt0+TJk5WSkiLLsrRhw4bAvvr6ev3sZz/TkCFD1L17d6WkpGjmzJn69NNP2+xzxYoVsiyrWampqQn5AwEA0CZjU4lCISf106dPa+jQoVqyZEmzfWfOnNHevXv12GOPae/evVq3bp0+/PBDfe973ztvvwkJCaqoqAgqcXFxoYYHAIBrhTz9np2drezs7Bb3JSYmqrCwMGjbs88+q2uuuUZHjhxRnz59Wu3Xsix5vd5QwwEAIHRROtIOV8TPqVdVVcmyLF166aVt1jt16pT69u2rxsZGffOb39S//uu/atiwYa3Wr62tVW1tbeB1dXX1uX9obJSs8Bb1G79D/9pNjP9iRxA5xuF/OwBh45x6O9XU1OiRRx7R9OnTlZCQ0Gq9AQMGaMWKFdq4caNeffVVxcXF6dprr9WhQ4dabZOfn6/ExMRASU1NjcRHAADgKyNiSb2+vl633Xab/H6/nn/++TbrZmRk6Ic//KGGDh2q66+/Xr/73e901VVX6dlnn221TU5OjqqqqgKlvLzc7o8AAHAiBy+Ui8j0e319vaZOnaqysjK9/fbbbY7SWxITE6ORI0e2OVL3eDzyeDzhhgoAcBmm30PQlNAPHTqkLVu26LLLLgu5D2OMSktLlZycbHd4AAA4Vsgj9VOnTumjjz4KvC4rK1Npaal69OihlJQU3Xrrrdq7d6/+53/+R42NjfL5fJKkHj16KDY2VpI0c+ZMXX755crPz5ckPf7448rIyFB6erqqq6v1zDPPqLS0VM8995wdnxEAgM9xR7nP7dmzR+PHjw+8nj9/viRp1qxZys3N1caNGyVJ3/zmN4Pabd26VePGjZMkHTlyRDExn08SnDhxQj/+8Y/l8/mUmJioYcOGadu2bbrmmmtCDQ8AgDY5efo95KQ+btw4mTYuG2prX5OioqKg17/5zW/0m9/8JtRQAADAF3DvdwCAuzh4+p1HrwIA3KWDL2nLz8/XyJEjFR8fr169eun73/++Dh48GBySMcrNzVVKSoq6du2qcePGaf/+/SF/NJI6AMBVms6ph1suVHFxsebOnatdu3apsLBQDQ0NysrK0unTpwN1Fi1apMWLF2vJkiXavXu3vF6vJk6cqJMnT4b02Zh+BwAggt58882g18uXL1evXr303nvvacyYMTLGqKCgQAsXLtSUKVMkSStXrlRSUpJWr16te++994Lfi5E6AMBdbJx+r66uDipffCZJa6qqqiSdu9RbOndpuM/nU1ZWVqCOx+PR2LFjtWPHjpA+GkkdAOAqljG2FElKTU0Neg5J0/1XWmOM0fz583Xddddp8ODBkhS4n0tSUlJQ3aSkpMC+C8X0OwAA7VReXh50K/Tz3b78/vvv15///Gdt37692T7LsoJeG2OabTsfkjoAwF1svKQtISHhgp9v8sADD2jjxo3atm2bevfuHdju9XolnRuxf/H26JWVlc1G7+fD9DsAwFU6evW7MUb333+/1q1bp7fffltpaWlB+9PS0uT1elVYWBjYVldXp+LiYmVmZob02RipAwAQQXPnztXq1av1+uuvKz4+PnCePDExUV27dpVlWZo3b57y8vKUnp6u9PR05eXlqVu3bpo+fXpI70VSBwC4SwffUW7p0qWSFHj+SZPly5dr9uzZkqQFCxbo7NmzmjNnjo4fP65Ro0Zp8+bNio+PDyksxyV109goY3FWoU0XcH9+AHCqjn6gy4U8E8WyLOXm5io3N7f9QYlz6gAAOIbjRuoAALTJwQ90IakDAFyF56kDAOAUDh6pc04dAACHYKQOAHCdaJ0+DxdJHQDgLsaEf2lvlF4azPQ7AAAOwUgdAOAqrH4HAMApWP0OAACiHSN1AICrWP5zJdw+ohFJHQDgLky/AwCAaMdIHQDgKqx+BwDAKRx88xmSOgDAVRipf5UYO1ZAAADw1eO8pA4AQFscvPqdpA4AcBUnT79zSRsAAA7BSB0A4C6sfgcAwBmYfv+Cbdu2afLkyUpJSZFlWdqwYUPQ/tmzZ8uyrKCSkZFx3n7Xrl2rgQMHyuPxaODAgVq/fn2ooQEA4GohJ/XTp09r6NChWrJkSat1vvOd76iioiJQ3njjjTb73Llzp6ZNm6YZM2bo/fff14wZMzR16lS9++67oYYHAEDbjE0lCoU8/Z6dna3s7Ow263g8Hnm93gvus6CgQBMnTlROTo4kKScnR8XFxSooKNCrr74aaogAALSK6fcQFRUVqVevXrrqqqt0zz33qLKyss36O3fuVFZWVtC2SZMmaceOHa22qa2tVXV1dVABAMDNbE/q2dnZ+u1vf6u3335bTz/9tHbv3q0bbrhBtbW1rbbx+XxKSkoK2paUlCSfz9dqm/z8fCUmJgZKamqqbZ8BAOBgfmNPiUK2r36fNm1a4J8HDx6sESNGqG/fvtq0aZOmTJnSajvLsoJeG2OabfuinJwczZ8/P/C6urqaxA4AOD/uKNd+ycnJ6tu3rw4dOtRqHa/X22xUXllZ2Wz0/kUej0cej8e2OAEA7mDJhnPqtkRiv4jfUe7YsWMqLy9XcnJyq3VGjx6twsLCoG2bN29WZmZmpMMDAMAxQh6pnzp1Sh999FHgdVlZmUpLS9WjRw/16NFDubm5uuWWW5ScnKy//e1vevTRR/X1r39d//zP/xxoM3PmTF1++eXKz8+XJD344IMaM2aMnnzySd188816/fXXtWXLFm3fvt2GjwgAwBdwR7nP7dmzR+PHjw+8bjqvPWvWLC1dulT79u3TqlWrdOLECSUnJ2v8+PF67bXXFB8fH2hz5MgRxcR8PkmQmZmpNWvW6Oc//7kee+wx9evXT6+99ppGjRoVzmcDAKAZJ1/SFnJSHzdunEwbv1Deeuut8/ZRVFTUbNutt96qW2+9NdRwAADAP3DvdwCAu7D6HQAAZ7CMkRXmOfFw20cKz1MHAMAhGKkDANzF/48Sbh9RiKQOAHAVpt8BAEDUY6QOAHAXVr8DAOAQ3FEOAABncPId5TinDgCAQzBSBwC4C9PvAAA4g+U/V8LtIxox/Q4AgEOQ1AEA7tI0/R5uCcG2bds0efJkpaSkyLIsbdiw4UshGeXm5iolJUVdu3bVuHHjtH///pA/GkkdAOAuxqYSgtOnT2vo0KFasmRJi/sXLVqkxYsXa8mSJdq9e7e8Xq8mTpyokydPhvQ+nFMHACDCsrOzlZ2d3eI+Y4wKCgq0cOFCTZkyRZK0cuVKJSUlafXq1br33nsv+H0YqQMAXKXp3u/hFkmqrq4OKrW1tSHHU1ZWJp/Pp6ysrMA2j8ejsWPHaseOHSH1RVIHALiLjefUU1NTlZiYGCj5+fkhh+Pz+SRJSUlJQduTkpIC+y4U0+8AALRTeXm5EhISAq89Hk+7+7IsK+i1MabZtvMhqQMA3MUo/Oeh/2OhXEJCQlBSbw+v1yvp3Ig9OTk5sL2ysrLZ6P18mH4HALiKnefU7ZCWliav16vCwsLAtrq6OhUXFyszMzOkvhipAwDcxciG28SGVv3UqVP66KOPAq/LyspUWlqqHj16qE+fPpo3b57y8vKUnp6u9PR05eXlqVu3bpo+fXpI70NSBwAgwvbs2aPx48cHXs+fP1+SNGvWLK1YsUILFizQ2bNnNWfOHB0/flyjRo3S5s2bFR8fH9L7kNQBAO5yER7oMm7cOJk22liWpdzcXOXm5oYVFkkdAOAufkmhLSpvuY8oxEI5AAAcgpE6AMBV7Fi9bufqdzuR1AEA7nIRzql3FKbfAQBwCEbqAAB3cfBInaQOAHAXByd1pt8BAHAIRuoAAHdx8HXqJHUAgKtwSRsAAE7BOfXPbdu2TZMnT1ZKSoosy9KGDRuC9luW1WL59a9/3WqfK1asaLFNTU1NyB8IAAC3Cnmkfvr0aQ0dOlQ/+tGPdMsttzTbX1FREfT6D3/4g+66664W635RQkKCDh48GLQtLi4u1PAAAGib30hWmCNtf3SO1ENO6tnZ2crOzm51v9frDXr9+uuva/z48bryyivb7NeyrGZtAQCwHdPv7fPZZ59p06ZNuuuuu85b99SpU+rbt6969+6tm266SSUlJW3Wr62tVXV1dVABAMDNIprUV65cqfj4eE2ZMqXNegMGDNCKFSu0ceNGvfrqq4qLi9O1116rQ4cOtdomPz9fiYmJgZKammp3+AAARzKfj9bbW+TCkfrLL7+sO+6447znxjMyMvTDH/5QQ4cO1fXXX6/f/e53uuqqq/Tss8+22iYnJ0dVVVWBUl5ebnf4AAAnCjeh2zF9HyERu6TtnXfe0cGDB/Xaa6+F3DYmJkYjR45sc6Tu8Xjk8XjCCREAAEeJ2Eh92bJlGj58uIYOHRpyW2OMSktLlZycHIHIAACu5jf2lCgU8kj91KlT+uijjwKvy8rKVFpaqh49eqhPnz6SpOrqav3+97/X008/3WIfM2fO1OWXX678/HxJ0uOPP66MjAylp6erurpazzzzjEpLS/Xcc8+15zMBANA64z9Xwu0jCoWc1Pfs2aPx48cHXs+fP1+SNGvWLK1YsUKStGbNGhljdPvtt7fYx5EjRxQT8/kkwYkTJ/TjH/9YPp9PiYmJGjZsmLZt26Zrrrkm1PAAAHAty5goPdsfourqaiUmJmqcblZnq8vFDgcAEIIGU68iva6qqiolJCRE5D2a8sS3U/9FnWPCW5PV4K/VlvKlEY23Pbj3OwDAXfw2XJLmlHPqAAB8pXFHOQAAEO0YqQMA3MXIhpG6LZHYjqQOAHAXpt8BAEC0Y6QOAHAXv19SmDeP8Tvk5jMAAHylMf0OAACiHSN1AIC7OHikTlIHALiLg+8ox/Q7AAAOwUgdAOAqxvhlwnx0arjtI4WkDgBwF2PCnz7nnDoAAFHA2HBOPUqTOufUAQBwCEbqAAB38fslK8xz4pxTBwAgCjD9DgAAoh0jdQCAqxi/XybM6XcuaQMAIBow/Q4AAKIdI3UAgLv4jWQ5c6ROUgcAuIsxksK9pC06kzrT7wAAOAQjdQCAqxi/kQlz+t1E6UidpA4AcBfjV/jT71zSBgDARefkkTrn1AEAcAjHjNSbfjU1qD7sewoAADpWg+oldcwIuMHUhj193hRvtHFMUj958qQkabveuMiRAADa6+TJk0pMTIxI37GxsfJ6vdrusydPeL1excbG2tKXXSwTrScGQuT3+/Xpp58qPj5elmW1WKe6ulqpqakqLy9XQkJCB0fYfsTd8b6qsRN3xyJu+xhjdPLkSaWkpCgmJnJnhmtqalRXV2dLX7GxsYqLi7OlL7s4ZqQeExOj3r17X1DdhISEqPkih4K4O95XNXbi7ljEbY9IjdC/KC4uLuoSsZ1YKAcAgEOQ1AEAcAhXJXWPx6Nf/OIX8ng8FzuUkBB3x/uqxk7cHYu4EW0cs1AOAAC3c9VIHQAAJyOpAwDgECR1AAAcgqQOAIBDOC6pP//880pLS1NcXJyGDx+ud955p836xcXFGj58uOLi4nTllVfqhRde6KBIz8nPz9fIkSMVHx+vXr166fvf/74OHjzYZpuioiJZltWs/PWvf+2gqKXc3Nxm7+/1ettsc7GPdZMrrriixeM3d+7cFutfrOO9bds2TZ48WSkpKbIsSxs2bAjab4xRbm6uUlJS1LVrV40bN0779+8/b79r167VwIED5fF4NHDgQK1fv77D4q6vr9fPfvYzDRkyRN27d1dKSopmzpypTz/9tM0+V6xY0eLfoKampkPilqTZs2c3e/+MjIzz9nsxj7ekFo+bZVn69a9/3WqfHXG8ERmOSuqvvfaa5s2bp4ULF6qkpETXX3+9srOzdeTIkRbrl5WV6cYbb9T111+vkpISPfroo/rJT36itWvXdljMxcXFmjt3rnbt2qXCwkI1NDQoKytLp0+fPm/bgwcPqqKiIlDS09M7IOLPDRo0KOj99+3b12rdaDjWTXbv3h0Ud2FhoSTpBz/4QZvtOvp4nz59WkOHDtWSJUta3L9o0SItXrxYS5Ys0e7du+X1ejVx4sTAcxBasnPnTk2bNk0zZszQ+++/rxkzZmjq1Kl69913OyTuM2fOaO/evXrssce0d+9erVu3Th9++KG+973vnbffhISEoONfUVFh653Bzne8Jek73/lO0Pu/8Ubb9xC/2MdbUrNj9vLLL8uyLN1yyy1t9hvp440IMQ5yzTXXmPvuuy9o24ABA8wjjzzSYv0FCxaYAQMGBG279957TUZGRsRiPJ/KykojyRQXF7daZ+vWrUaSOX78eMcF9iW/+MUvzNChQy+4fjQe6yYPPvig6devn/H7/S3uj4bjLcmsX78+8Nrv9xuv12ueeOKJwLaamhqTmJhoXnjhhVb7mTp1qvnOd74TtG3SpEnmtttusz1mY5rH3ZI//elPRpI5fPhwq3WWL19uEhMT7Q2uDS3FPWvWLHPzzTeH1E80Hu+bb77Z3HDDDW3W6ejjDfs4ZqReV1en9957T1lZWUHbs7KytGPHjhbb7Ny5s1n9SZMmac+ePaqvvziP1auqqpIk9ejR47x1hw0bpuTkZE2YMEFbt26NdGjNHDp0SCkpKUpLS9Ntt92mjz/+uNW60XispXPfm1deeUV33nlnqw8CanKxj/cXlZWVyefzBR1Tj8ejsWPHtvp9l1r/O7TVJtKqqqpkWZYuvfTSNuudOnVKffv2Ve/evXXTTTeppKSkYwL8gqKiIvXq1UtXXXWV7rnnHlVWVrZZP9qO92effaZNmzbprrvuOm/daDjeCJ1jkvrRo0fV2NiopKSkoO1JSUny+XwttvH5fC3Wb2ho0NGjRyMWa2uMMZo/f76uu+46DR48uNV6ycnJevHFF7V27VqtW7dO/fv314QJE7Rt27YOi3XUqFFatWqV3nrrLb300kvy+XzKzMzUsWPHWqwfbce6yYYNG3TixAnNnj271TrRcLy/rOk7Hcr3valdqG0iqaamRo888oimT5/e5oNFBgwYoBUrVmjjxo169dVXFRcXp2uvvVaHDh3qsFizs7P129/+Vm+//baefvpp7d69WzfccINqa2tbbRNtx3vlypWKj4/XlClT2qwXDccb7eOYp7Q1+fJoyxjT5gispfotbe8I999/v/785z9r+/btbdbr37+/+vfvH3g9evRolZeX66mnntKYMWMiHaakc/+BazJkyBCNHj1a/fr108qVKzV//vwW20TTsW6ybNkyZWdnKyUlpdU60XC8WxPq9729bSKhvr5et912m/x+v55//vk262ZkZAQtSrv22mv1rW99S88++6yeeeaZSIcqSZo2bVrgnwcPHqwRI0aob9++2rRpU5tJMlqOtyS9/PLLuuOOO857bjwajjfaxzEj9a9//evq1KlTs1/AlZWVzX4pN/F6vS3W79y5sy677LKIxdqSBx54QBs3btTWrVsv+BGyX5SRkXFRf0V3795dQ4YMaTWGaDrWTQ4fPqwtW7bo7rvvDrntxT7eTVcahPJ9b2oXaptIqK+v19SpU1VWVqbCwsKQH/8ZExOjkSNHXtS/QXJysvr27dtmDNFyvCXpnXfe0cGDB9v1fY+G440L45ikHhsbq+HDhwdWMjcpLCxUZmZmi21Gjx7drP7mzZs1YsQIdenSJWKxfpExRvfff7/WrVunt99+W2lpae3qp6SkRMnJyTZHd+Fqa2t14MCBVmOIhmP9ZcuXL1evXr303e9+N+S2F/t4p6Wlyev1Bh3Turo6FRcXt/p9l1r/O7TVxm5NCf3QoUPasmVLu37UGWNUWlp6Uf8Gx44dU3l5eZsxRMPxbrJs2TINHz5cQ4cODbltNBxvXKCLtUIvEtasWWO6dOlili1bZj744AMzb9480717d/O3v/3NGGPMI488YmbMmBGo//HHH5tu3bqZhx56yHzwwQdm2bJlpkuXLua///u/Oyzmf/mXfzGJiYmmqKjIVFRUBMqZM2cCdb4c929+8xuzfv168+GHH5q//OUv5pFHHjGSzNq1azss7p/+9KemqKjIfPzxx2bXrl3mpptuMvHx8VF9rL+osbHR9OnTx/zsZz9rti9ajvfJkydNSUmJKSkpMZLM4sWLTUlJSWCV+BNPPGESExPNunXrzL59+8ztt99ukpOTTXV1daCPGTNmBF398cc//tF06tTJPPHEE+bAgQPmiSeeMJ07dza7du3qkLjr6+vN9773PdO7d29TWloa9J2vra1tNe7c3Fzz5ptvmv/93/81JSUl5kc/+pHp3Lmzeffddzsk7pMnT5qf/vSnZseOHaasrMxs3brVjB492lx++eVRfbybVFVVmW7dupmlS5e22MfFON6IDEcldWOMee6550zfvn1NbGys+da3vhV0adisWbPM2LFjg+oXFRWZYcOGmdjYWHPFFVe0+qWPFEktluXLl7ca95NPPmn69etn4uLizNe+9jVz3XXXmU2bNnVo3NOmTTPJycmmS5cuJiUlxUyZMsXs37+/1ZiNufjH+oveeustI8kcPHiw2b5oOd5Nl9J9ucyaNcsYc+6ytl/84hfG6/Uaj8djxowZY/bt2xfUx9ixYwP1m/z+9783/fv3N126dDEDBgyw/cdJW3GXlZW1+p3funVrq3HPmzfP9OnTx8TGxpqePXuarKwss2PHjg6L+8yZMyYrK8v07NnTdOnSxfTp08fMmjXLHDlyJKiPaDveTf7jP/7DdO3a1Zw4caLFPi7G8UZk8OhVAAAcwjHn1AEAcDuSOgAADkFSBwDAIUjqAAA4BEkdAACHIKkDAOAQJHUAAByCpA4AgEOQ1AEAcAiSOgAADkFSBwDAIUjqAAA4xP8D7z4vfwrZU84AAAAASUVORK5CYII=", + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgIAAAGdCAYAAABgnRvHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA9oklEQVR4nO3dfVRU17k/8O8ZXmbUwqRGYSAiEleiol6K+AJYjMaIIYm1jVaStGAabOqNaTU0Kw0xuTG56xeSNjFE40ttVeJNgzQX3+6SVnFFRCvxRgJp2hpLr0S4hinFWxjRMAMz+/eHYeJxXmCcPTCZ8/1knbUyZ/Z55pnDCM/svc8+ihBCgIiIiDRJN9QJEBER0dBhIUBERKRhLASIiIg0jIUAERGRhrEQICIi0jAWAkRERBrGQoCIiEjDWAgQERFpWPhQJyCLw+HAZ599hqioKCiKMtTpEBGRD4QQuHTpEuLj46HTBe47and3N2w2m5RYkZGRMBgMUmINpZApBD777DMkJCQMdRpEROSHlpYWjBkzJiCxu7u7kZT4NZjb7FLimUwmNDU1feWLgZApBKKiogAAc4YvRbgS4VcsXcwoGSkBAOw3jZASx6GX96MSOjk9JopD3urUOluvnDhXeqTEAQDl8245gbrlfPsAAGGVk5OwyTtPokfOL1VhlxPnajCHpDhcgX2w9KIHJ1Dp/F0eCDabDeY2O5rqEhEd5V+vg+WSA0lp52Gz2VgIBIu+4YBwJQLhSqRfsXQ6vYyUAABKmJxYjnD/iptrBWUhYJdUCITJ61JUdJLen8ReTqHI+QMnJA6fCUXOz04oMruDJRUCYCEwaL441YMxtBsdpfO7EAglIVMIEBERDYRdOGD3s8azy+p1CgIsBIiISFMcEHD42dvj7/HBJGB9I5s3b0ZSUhIMBgPS0tJw/Phxr+2PHTuGtLQ0GAwG3Hrrrdi6dWugUiMiIg1zSPovVASkECgvL8eaNWuwdu1a1NfXIysrCzk5OWhubnbbvqmpCffccw+ysrJQX1+PZ555Bj/5yU9QUVERiPSIiIjoCwEpBNavX4+CggKsWLECkyZNQklJCRISErBlyxa37bdu3YqxY8eipKQEkyZNwooVK/DII4/g1VdfDUR6RESkYXYhpGyhQnohYLPZUFdXh+zsbNX+7OxsnDx50u0xtbW1Lu0XLlyI06dPo6fH/WVOVqsVFotFtREREfWnb46Av1uokF4ItLe3w263IzY2VrU/NjYWZrPZ7TFms9lt+97eXrS3t7s9pri4GEaj0blxMSEiIiLfBWyy4PXXggohvF4f6q69u/19ioqK0NnZ6dxaWlr8zJiIiLTAAQG7n1so9QhIv3xw1KhRCAsLc/n239bW5vKtv4/JZHLbPjw8HDfffLPbY/R6PfR6eQv/EBGRNvDyQTXpPQKRkZFIS0tDVVWVan9VVRUyMzPdHpORkeHS/vDhw5g+fToiIuStqEdERERqARkaKCwsxK9//Wvs2LEDZ86cwRNPPIHm5masXLkSwNVu/fz8fGf7lStX4vz58ygsLMSZM2ewY8cObN++HU8++WQg0iMiIg3jVQNqAVlZMDc3FxcvXsSLL76I1tZWTJkyBZWVlUhMTAQAtLa2qtYUSEpKQmVlJZ544gls2rQJ8fHx2LBhA5YsWRKI9IiISMMc8P9uFKGznBCgCBEaZY3FYoHRaMSdIx70/6ZDsaMlZQXYvy7p7oOGEL/pkFXW3Qfl3elPuSLr7oNWOXEAiO5gvPugpJsO8e6DmtYrelCN/ejs7ER0dHRAXqPv78QnZ2IR5edNhy5dcmDipL8HNN/BwnsNEBGRpvTN/Pc3RqgIuUJAFzPK79sI/98sk6RsgM7xcqZhWEfJ64hyRMqJpbPKm2IS2SEnlsH9shM3ZHibnG+ow9vk9VJEtF+REkfX2SUlDgCIy3JyElaJPSc2OedcSOz1Yi9F8LALSLj7oJxcgkHIFQJERETecI6AWsAWFCIiIqLgxx4BIiLSFAcU2OHfpGmHn8cHExYCRESkKQ5xdfM3Rqjg0AAREZGGsUeAiIg0xS5haMDf44MJCwEiItIUFgJqHBogIiLSMPYIEBGRpjiEAofw86oBP48PJiwEiIhIUzg0oMahASIiIg1jjwAREWmKHTrY/fweLPFemUOOhQAREWmKkDBHQHCOABER0VcT5wiocY4AERGRhrFHgIiINMUudLALP+cIhNC9BlgIEBGRpjigwOFnh7gDoVMJcGiAiIhIw0KuR8B+0wgoYXq/YnSOl1cfhaV1SIlz1y1NUuIAwOjILilxWruNUuIAQMM/4qXE+b+Wm6TEAYCer8n55+GI8O/zeK0RipwJSpFC3rcZxS7pQipZcQAISbEUiReJiVC63uwrjpMF1UKuECAiIvJGzhwBDg0QERFRCGCPABERacrVyYJ+3nSIQwNERERfTQ4JSwzzqgEiIiIKCewRICIiTeFkQTUWAkREpCkO6Lig0DVYCBARkabYhQK7n3cP9Pf4YMI5AkRERBrGHgEiItIUu4SrBuwcGiAiIvpqcggdHH5OFnSE0GRBDg0QERFpGHsEiIhIUzg0oMZCgIiINMUB/2f9O+SkEhSkDw0UFxdjxowZiIqKQkxMDL797W/j7NmzXo+prq6Goigu2yeffCI7PSIiIrqG9B6BY8eOYdWqVZgxYwZ6e3uxdu1aZGdn4y9/+QtGjBjh9dizZ88iOjra+Xj06NGy0yMiIo2Ts6BQ6Eyxk14I/P73v1c93rlzJ2JiYlBXV4c5c+Z4PTYmJgY33XST7JSIiIic5CwxzEJgwDo7OwEAI0eO7Ldtamoquru7kZycjGeffRbz5s3z2NZqtcJqtTofWywWAIBDHw5HeIRfOVtHyRv9ueuWJilxCkbXSIkDAPFh1v4bDcC53q9JiQMAIyNTpcT5fe8kKXEAoLvLKCVOZKe8Fcj0HWFS4oTr/fs3cq2wSEmxbPJ+HSm9vVLiCLtdShwiANi8eTN+8YtfoLW1FZMnT0ZJSQmysrLctt2zZw+2bNmChoYGWK1WTJ48GevWrcPChQudbUpLS/GDH/zA5djPP/8cBoNhwHkFtKQRQqCwsBDf/OY3MWXKFI/t4uLisG3bNlRUVGDPnj2YMGEC5s+fj5oaz3/8iouLYTQanVtCQkIg3gIREYUYBxQpmy/Ky8uxZs0arF27FvX19cjKykJOTg6am5vdtq+pqcGCBQtQWVmJuro6zJs3D4sWLUJ9fb2qXXR0NFpbW1WbL0UAEOAegccffxx//OMfceLECa/tJkyYgAkTJjgfZ2RkoKWlBa+++qrH4YSioiIUFhY6H1ssFhYDRETUr6EYGli/fj0KCgqwYsUKAEBJSQkOHTqELVu2oLi42KV9SUmJ6vFLL72E/fv347/+67+QmvplD6qiKDCZTL6/gWsErEfgxz/+MQ4cOICjR49izJgxPh+fnp6OxsZGj8/r9XpER0erNiIiov70rSPg7zZQNpsNdXV1yM7OVu3Pzs7GyZMnBxTD4XDg0qVLLsPsXV1dSExMxJgxY3Dfffe59BgMhPRCQAiBxx9/HHv27MF7772HpKSkG4pTX1+PuLg4ydkRERHJY7FYVNu1c9f6tLe3w263IzY2VrU/NjYWZrN5QK/z2muv4fLly1i2bJlz38SJE1FaWooDBw6grKwMBoMBs2fP9vol2h3pQwOrVq3CO++8g/379yMqKsr5Jo1GI4YNGwbgarf+hQsXsGvXLgBXu0DGjRuHyZMnw2az4e2330ZFRQUqKipkp0dERBrnEAoc/i4o9MXx1w9JP//881i3bp3bYxRF/ZpCCJd97pSVlWHdunXYv38/YmJinPvT09ORnp7ufDx79mxMmzYNGzduxIYNGwb6VuQXAlu2bAEAzJ07V7V/586dePjhhwEAra2tqgkSNpsNTz75JC5cuIBhw4Zh8uTJOHjwIO655x7Z6RERkcY5JCwx3LeOQEtLi2poWq/Xu7QdNWoUwsLCXL79t7W1ufQSXK+8vBwFBQV49913cdddd3ltq9PpMGPGjKHvERADuCNTaWmp6vFTTz2Fp556SnYqREREATWQOWqRkZFIS0tDVVUVvvOd7zj3V1VVYfHixR6PKysrwyOPPIKysjLce++9/eYihEBDQwOmTp068DcA3muAiIg0Rs5tiH07vrCwEHl5eZg+fToyMjKwbds2NDc3Y+XKlQBch8zLysqQn5+PN954A+np6c7ehGHDhsFovLrGyQsvvID09HTcdtttsFgs2LBhAxoaGrBp0yafcmMhQEREmmKHAruP6wC4i+GL3NxcXLx4ES+++CJaW1sxZcoUVFZWIjExEYDrkPkvf/lL9Pb2YtWqVVi1apVz//Lly5296h0dHXj00UdhNpthNBqRmpqKmpoazJw506fcWAgQERENgsceewyPPfaY2+euHzKvrq7uN97rr7+O119/3e+8WAgQEZGmDMXQQDBjIUBERJpih+9d++5ihIrQKWmIiIjIZ+wRICIiTeHQgBoLASIi0pShuOlQMGMhQEREmiJu4DbC7mKEitApaYiIiMhn7BEgIiJN4dCAWsgVAkKnQOj8vKtUpENSNsDoyC4pceLDXG9teaPGhH9NShyrkPPeACAm0iIlzgi9TUocALgS0f99MwZChMv7hSHCJHVHyooDAAO4expRMJF598FQEDolDREREfks5HoEiIiIvLFLuA2xv8cHExYCRESkKRwaUAudkoaIiIh8xh4BIiLSFAd0cPj5Pdjf44MJCwEiItIUu1Bg97Nr39/jg0nolDRERETkM/YIEBGRpnCyoBoLASIi0hQh4e6DgisLEhERfTXZocDu502D/D0+mIROSUNEREQ+Y48AERFpikP4P8bvkHMrkqDAQoCIiDTFIWGOgL/HB5PQeSdERETkM/YIEBGRpjigwOHnZD9/jw8mLASIiEhTuLKgGocGiIiINCzkegQUh4Di53ROnVVefdTabZQS51zv16TEAQCr6JIS5396vi4lDgC0dI+UEufS53opcQBAZ5VT8et65E0vVuySYsmKAwBCUqxQmoZNQY2TBdVCrhAgIiLyxgEJSwyH0ByB0ClpiIiIyGfsESAiIk0REq4aECHUI8BCgIiINIV3H1RjIUBERJrCyYJq0t/JunXroCiKajOZTF6POXbsGNLS0mAwGHDrrbdi69atstMiIiIiNwLSIzB58mQcOXLE+TgsLMxj26amJtxzzz344Q9/iLfffht/+MMf8Nhjj2H06NFYsmRJINIjIiIN49CAWkAKgfDw8H57Afps3boVY8eORUlJCQBg0qRJOH36NF599VUWAkREJB2XGFYLyCBHY2Mj4uPjkZSUhAceeADnzp3z2La2thbZ2dmqfQsXLsTp06fR09Pj8Tir1QqLxaLaiIiIyDfSC4FZs2Zh165dOHToEH71q1/BbDYjMzMTFy9edNvebDYjNjZWtS82Nha9vb1ob2/3+DrFxcUwGo3OLSEhQer7ICKi0NQ3NODvFiqkFwI5OTlYsmQJpk6dirvuugsHDx4EALz11lsej1EU9QkVXyxZev3+axUVFaGzs9O5tbS0SMieiIhCHQsBtYBfPjhixAhMnToVjY2Nbp83mUwwm82qfW1tbQgPD8fNN9/sMa5er4deL29deSIiIi0K+IWQVqsVZ86cQVxcnNvnMzIyUFVVpdp3+PBhTJ8+HREREYFOj4iINIY9AmrSC4Enn3wSx44dQ1NTE06dOoWlS5fCYrFg+fLlAK526efn5zvbr1y5EufPn0dhYSHOnDmDHTt2YPv27XjyySdlp0ZERMRC4DrShwb+93//Fw8++CDa29sxevRopKen4/3330diYiIAoLW1Fc3Nzc72SUlJqKysxBNPPIFNmzYhPj4eGzZs4KWDREREg0B6IbB7926vz5eWlrrsu+OOO/Dhhx/KToWIiMiFgP/rAAg5qQQF3muAiIg0hSsLqrEQICIiTWEhoBZyhYDO1gudvdevGJEd8uZQNvwjXkqckZGpUuIAQEyknFUYW7pHSokDAHXtchaEutI+XEocABjeKecfesRlh5Q4ABD+uV1KHJ3Nv38jKr1ycoKQd57gkBiLKMSFXCFARETkDXsE1FgIEBGRprAQUAv4gkJEREQUvNgjQEREmiKEAuHnN3p/jw8mLASIiEhTHFD8XkfA3+ODCYcGiIiINIyFABERacpQ3Wtg8+bNSEpKgsFgQFpaGo4fP+6x7Z49e7BgwQKMHj0a0dHRyMjIwKFDh1zaVVRUIDk5GXq9HsnJydi7d6/PebEQICIiTembI+Dv5ovy8nKsWbMGa9euRX19PbKyspCTk6O69861ampqsGDBAlRWVqKurg7z5s3DokWLUF9f72xTW1uL3Nxc5OXl4aOPPkJeXh6WLVuGU6dO+ZSbIoQIiSWTLRYLjEYj5k17GuFhBr9inb83WlJWwPDp7VLi3HlLo5Q4QGgvKPTZ+ZulxAGA4Z/KuQ12VIu8xW2+dsEmJU5k22UpcQBA6bgkJY64LC8n0W2VE6dX3sJLwi5r4aWQ+JXtolf0oBr70dnZiehoeb+Dr9X3d2Lm3tUIH6H3K1bvZSv++ztvDDjfWbNmYdq0adiyZYtz36RJk/Dtb38bxcXFA3rNyZMnIzc3F//2b/8GAMjNzYXFYsHvfvc7Z5u7774bX//611FWVjbg98IeASIi0hSZQwMWi0W1Wa2uRajNZkNdXR2ys7NV+7Ozs3Hy5MmB5exw4NKlSxg58ssvYLW1tS4xFy5cOOCYfVgIEBGRpsgcGkhISIDRaHRu7r7dt7e3w263IzY2VrU/NjYWZrN5QDm/9tpruHz5MpYtW+bcZzab/YrZh5cPEhGRpggJKwv2FQItLS2qoQG93vOQg6KoX1MI4bLPnbKyMqxbtw779+9HTEyMlJjXYiFARER0g6Kjo/udIzBq1CiEhYW5fFNva2tz+UZ/vfLychQUFODdd9/FXXfdpXrOZDLdUMzrcWiAiIg0ReDqnEu/Nh9eLzIyEmlpaaiqqlLtr6qqQmZmpsfjysrK8PDDD+Odd97Bvffe6/J8RkaGS8zDhw97jekOewSIiEhTHFCgDPLKgoWFhcjLy8P06dORkZGBbdu2obm5GStXrgQAFBUV4cKFC9i1axeAq0VAfn4+3njjDaSnpzu/+Q8bNgxGoxEAsHr1asyZMwevvPIKFi9ejP379+PIkSM4ceKET7mxR4CIiCjAcnNzUVJSghdffBHf+MY3UFNTg8rKSiQmJgIAWltbVWsK/PKXv0Rvby9WrVqFuLg457Z69Wpnm8zMTOzevRs7d+7Ev/zLv6C0tBTl5eWYNWuWT7mxR4CIiDRlqG469Nhjj+Gxxx5z+1xpaanqcXV19YBiLl26FEuXLvU5l2uFXCGgu9IDXZh/HR0GOWsAAQD+r+UmKXF+3ztJShwAGKGXsyjNpc/9W5DjWlfah0uJozfLWQQIAAwX5Szcou+Qt6BQeFePnEBWOZ8BABA2OTkJu7zzJCuWcITm4j1a5xAKFD8LAX+vOggmHBogIiLSsJDrESAiIvKmb+a/vzFCBQsBIiLSlKGaIxCsODRARESkYewRICIiTWGPgBoLASIi0hReNaDGQoCIiDSFkwXVOEeAiIhIw9gjQEREmnK1R8DfOQKSkgkCLASIiEhTOFlQjUMDREREGsYeASIi0hTxxeZvjFDBQoCIiDSFQwNqHBogIiLSMOmFwLhx46Aoisu2atUqt+2rq6vdtv/kk09kp0ZERPTl2IC/W4iQPjTwwQcfwG63Ox//6U9/woIFC/Dd737X63Fnz55FdHS08/Ho0aNlp0ZERARIGBpACA0NSC8Erv8D/vLLL2P8+PG44447vB4XExODm266SXY6REREKlxZUC2gcwRsNhvefvttPPLII1AU79VTamoq4uLiMH/+fBw9ejSQaREREdEXAnrVwL59+9DR0YGHH37YY5u4uDhs27YNaWlpsFqt+I//+A/Mnz8f1dXVmDNnjsfjrFYrrFar87HFYgEAKJ93Q9H5V6oNb7P332iAer4m5xR3dxmlxAGAKxFySlmdVV7X2PBOObEMF+WV6cP/4ZASJ7LDJiUOAOi6uqXEUT639t9ogESPpPfX0yMnDgAIOT87aXGA0PoK+RXHqwbUAloIbN++HTk5OYiPj/fYZsKECZgwYYLzcUZGBlpaWvDqq696LQSKi4vxwgsvSM2XiIg0QCj+j/GHUCEQsKGB8+fP48iRI1ixYoXPx6anp6OxsdFrm6KiInR2djq3lpaWG02ViIhIswLWI7Bz507ExMTg3nvv9fnY+vp6xMXFeW2j1+uh1+tvND0iItIoThZUC0gh4HA4sHPnTixfvhzh4eqXKCoqwoULF7Br1y4AQElJCcaNG4fJkyc7JxdWVFSgoqIiEKkREZHWcY1hlYAUAkeOHEFzczMeeeQRl+daW1vR3NzsfGyz2fDkk0/iwoULGDZsGCZPnoyDBw/innvuCURqREREdI2AFALZ2dkQHvpNSktLVY+feuopPPXUU4FIg4iIyAWvGlDjTYeIiEh7Qqhr31+86RAREZGGsUeAiIg0hUMDaiwEiIhIW3jVgAoLASIi0hjli83fGKGBcwSIiIg0jD0CRESkLRwaUGEhQERE2sJCQIVDA0RERBrGHgEiItIW3oZYhYUAERFpCu8+qBZ6hUC3ze8Bj+FtNjm5AHBEyLlVcmSnvOpThMsZEdL1yPuXEHHZISWOvkNOHACI7JDzOYj45+dS4gCAcqVbShxhtUqJAwCip1dOHLu8n51whNBvaaIAC71CgIiIyBtOFlRhIUBERNrCOQIqvGqAiIhIw9gjQEREmqKIq5u/MUIFCwEiItIWzhFQYSFARETawjkCKpwjQEREpGHsESAiIm3h0IAKCwEiItIWFgIqHBogIiLSMPYIEBGRtrBHQIWFABERaQuvGlDh0AAREZGGsRAgIiJN6VtZ0N/NV5s3b0ZSUhIMBgPS0tJw/Phxj21bW1vx0EMPYcKECdDpdFizZo1Lm9LSUiiK4rJ1d/t2l1IWAkREpC1C0uaD8vJyrFmzBmvXrkV9fT2ysrKQk5OD5uZmt+2tVitGjx6NtWvXIiUlxWPc6OhotLa2qjaDweBTbiwEiIiIAmz9+vUoKCjAihUrMGnSJJSUlCAhIQFbtmxx237cuHF44403kJ+fD6PR6DGuoigwmUyqzVcsBIiIiG6QxWJRbVar1aWNzWZDXV0dsrOzVfuzs7Nx8uRJv16/q6sLiYmJGDNmDO677z7U19f7HIOFABERaYoCCXMEvoiVkJAAo9Ho3IqLi11er729HXa7HbGxsar9sbGxMJvNN/w+Jk6ciNLSUhw4cABlZWUwGAyYPXs2GhsbfYoTcpcPCms3hOLwK0ZE+xVJ2QAjFDmXmOg7wqTEAQARJicnxS7vQtrwz+1y4nT1SIkDALou3ybceKJckRMHAMRlOZ9NYbVJiQMAordXUiD//t0GJJYIoYvF6UsSLx9saWlBdHS0c7der/d4iHLd3wMhhMs+X6SnpyM9Pd35ePbs2Zg2bRo2btyIDRs2DDhOyBUCREREgyU6OlpVCLgzatQohIWFuXz7b2trc+kl8IdOp8OMGTN87hHg0AAREWnLIF81EBkZibS0NFRVVan2V1VVITMz07/3cg0hBBoaGhAXF+fTcewRICIibRmCJYYLCwuRl5eH6dOnIyMjA9u2bUNzczNWrlwJACgqKsKFCxewa9cu5zENDQ0Ark4I/Mc//oGGhgZERkYiOTkZAPDCCy8gPT0dt912GywWCzZs2ICGhgZs2rTJp9xYCBAREQVYbm4uLl68iBdffBGtra2YMmUKKisrkZiYCODqAkLXrymQmprq/P+6ujq88847SExMxKeffgoA6OjowKOPPgqz2Qyj0YjU1FTU1NRg5syZPuWmCBEas2EsFguMRiPm35SHcCXSv2C3+H4dpie20SOkxOkdzsmCA4rDyYIDiyNzsqBNUiy7nM8AAAhZsULj1+NXQq/oQTX2o7Ozs98x9xvV93di3P/7f9D5uOjO9Rzd3fh07dqA5jtYfJ4jUFNTg0WLFiE+Ph6KomDfvn2q54UQWLduHeLj4zFs2DDMnTsXf/7zn/uNW1FRgeTkZOj1eiQnJ2Pv3r2+pkZERNS/IVhZMJj5XAhcvnwZKSkpePPNN90+//Of/xzr16/Hm2++iQ8++AAmkwkLFizApUuXPMasra1Fbm4u8vLy8NFHHyEvLw/Lli3DqVOnfE2PiIiIfODzHIGcnBzk5OS4fU4IgZKSEqxduxb3338/AOCtt95CbGws3nnnHfzoRz9ye1xJSQkWLFiAoqIiAFcnTRw7dgwlJSUoKyvzNUUiIiLPhmCyYDCTevlgU1MTzGazahlFvV6PO+64w+syirW1tS5LLy5cuNDrMVar1WVpRyIiov4M1d0Hg5XUQqBvsQRfl1E0m80+H1NcXKxa1jEhIcGPzImIiLQpIAsK3cgyir4eU1RUhM7OTufW0tJy4wkTEZF29C0x7O8WIqSuI9B3+0Oz2axa2ai/ZRRNJpPPSy/q9XqvazoTERG5xTkCKlJ7BJKSkmAymVTLKNpsNhw7dszrMooZGRkuSy8ePnxY6tKLREREAOcIXM/nHoGuri787W9/cz5uampCQ0MDRo4cibFjx2LNmjV46aWXcNttt+G2227DSy+9hOHDh+Ohhx5yHpOfn49bbrnFebvG1atXY86cOXjllVewePFi7N+/H0eOHMGJEyckvEUiIiLyxOdC4PTp05g3b57zcWFhIQBg+fLlKC0txVNPPYXPP/8cjz32GP75z39i1qxZOHz4MKKiopzHNDc3Q6f7sjMiMzMTu3fvxrPPPovnnnsO48ePR3l5OWbNmuXPeyMiInLFoQEVLjHsDpcYHhAuMTwwXGJ4gLjEsKYN5hLDtz73EsL8XGLY3t2Nc//+jDaXGCYiIqLQEXJ3HxS2Hoh+LlXsj66zS1I2QKSkbxTh+ggpcQAAknoEILFHQGfrlRNI4jdd5XOrlDjCKifO1ViS3l+PvJ4TWd/khYPfvmmQcGhAJeQKASIiIq9YCKhwaICIiEjD2CNARESaImMdgFBaR4A9AkRERBrGQoCIiEjDODRARETawsmCKiwEiIhIUzhHQI2FABERaU8I/SH3F+cIEBERaRh7BIiISFs4R0CFhQAREWkK5wiocWiAiIhIw9gjQERE2sKhARUWAkREpCkcGlDj0AAREZGGsUeAiIi0hUMDKiwEiIhIW1gIqHBogIiISMNCrkdA9NghlF7/Yly+IikbQLHbpcQJi4yQEgcAoChy4giJJXGvnPMkbD1S4gCA6LFJiuPf51EVq1dSLEmfSwAQDkmfA+GQEweQ+9mkkMPJgmohVwgQERF5xaEBFRYCRESkLSwEVDhHgIiISMPYI0BERJrCOQJqLASIiEhbODSgwqEBIiIiDWOPABERaQqHBtRYCBARkbZwaECFQwNEREQaxh4BIiLSFvYIqLAQICIiTVG+2PyNESo4NEBERKRh7BEgIiJt4dCACgsBIiLSFF4+qObz0EBNTQ0WLVqE+Ph4KIqCffv2OZ/r6enBz372M0ydOhUjRoxAfHw88vPz8dlnn3mNWVpaCkVRXLbu7m6f3xAREZFXQtIWInwuBC5fvoyUlBS8+eabLs9duXIFH374IZ577jl8+OGH2LNnD/7617/iW9/6Vr9xo6Oj0draqtoMBoOv6REREZEPfB4ayMnJQU5OjtvnjEYjqqqqVPs2btyImTNnorm5GWPHjvUYV1EUmEwmX9MhIiLyXQh9o/dXwOcIdHZ2QlEU3HTTTV7bdXV1ITExEXa7Hd/4xjfw7//+70hNTfXY3mq1wmq1Oh9bLBYAgLDbIRT/LoYQ18T1m90uJ44tCKdzOCT+SxIOOWHscuIAAHp6pISRmpOs8xSEPzsI/mamwcE5AmoBvXywu7sbTz/9NB566CFER0d7bDdx4kSUlpbiwIEDKCsrg8FgwOzZs9HY2OjxmOLiYhiNRueWkJAQiLdAREQU0gJWCPT09OCBBx6Aw+HA5s2bvbZNT0/H97//faSkpCArKwu//e1vcfvtt2Pjxo0ejykqKkJnZ6dza2lpkf0WiIgoFA3RZMHNmzcjKSkJBoMBaWlpOH78uMe2ra2teOihhzBhwgTodDqsWbPGbbuKigokJydDr9cjOTkZe/fu9TmvgBQCPT09WLZsGZqamlBVVeW1N8BtUjodZsyY4bVHQK/XIzo6WrURERH1p29owN/NF+Xl5VizZg3Wrl2L+vp6ZGVlIScnB83NzW7bW61WjB49GmvXrkVKSorbNrW1tcjNzUVeXh4++ugj5OXlYdmyZTh16pRPuUkvBPqKgMbGRhw5cgQ333yzzzGEEGhoaEBcXJzs9IiIiAbd+vXrUVBQgBUrVmDSpEkoKSlBQkICtmzZ4rb9uHHj8MYbbyA/Px9Go9Ftm5KSEixYsABFRUWYOHEiioqKMH/+fJSUlPiUm8+FQFdXFxoaGtDQ0AAAaGpqQkNDA5qbm9Hb24ulS5fi9OnT+M1vfgO73Q6z2Qyz2QybzeaMkZ+fj6KiIufjF154AYcOHcK5c+fQ0NCAgoICNDQ0YOXKlb6mR0RE5J3EoQGLxaLarG4mm9tsNtTV1SE7O1u1Pzs7GydPnrzht1FbW+sSc+HChT7H9Hkq+unTpzFv3jzn48LCQgDA8uXLsW7dOhw4cAAA8I1vfEN13NGjRzF37lwAQHNzM3S6L2uQjo4OPProozCbzTAajUhNTUVNTQ1mzpzpa3pEREReybxq4PqJ6s8//zzWrVun2tfe3g673Y7Y2FjV/tjYWJjN5hvOwWw2S4npcyEwd+5cCC+X+Xh7rk91dbXq8euvv47XX3/d11SIiIiGVEtLi2qOml6v99hWUdT3LBRCuOzzlYyYQXhxOhERUQBJvOnQQCarjxo1CmFhYS7f1Nva2ly+0fvCZDJJicnbEBMRkbYM8uWDkZGRSEtLc1l5t6qqCpmZmTf8NjIyMlxiHj582OeY7BEgIiJNGYqVBQsLC5GXl4fp06cjIyMD27ZtQ3Nzs3NSfFFRES5cuIBdu3Y5j+mblN/V1YV//OMfaGhoQGRkJJKTkwEAq1evxpw5c/DKK69g8eLF2L9/P44cOYITJ074lBsLASIiogDLzc3FxYsX8eKLL6K1tRVTpkxBZWUlEhMTAVxdQOj6NQWuXWa/rq4O77zzDhITE/Hpp58CADIzM7F79248++yzeO655zB+/HiUl5dj1qxZPuWmiIHM7vsKsFgsMBqNmKu7H+FKhF+xdAbPkz18pYRLqrUigrBmC8L16nmvgYGGCb6fHe81oG29ogfV2I/Ozs6ALRDX93ciJf8lhEX6d3dbu60bH+16JqD5DpYg/OtCREQUOIoQUPwsPP09PphwsiAREZGGsUeAiIi0ReLlg6GAhQAREWnKUFw1EMw4NEBERKRh7BEgIiJt4dCASugVAsIBwL/LmcQ1d0r0l7DbpcRRenulxAlajiC8fDCUL9WTKYRmT5M2cGhAjUMDREREGhZ6PQJERETecGhAhYUAERFpCocG1FgIEBGRtrBHQIVzBIiIiDSMPQJERKQ5odS17y8WAkREpC1C+H/ZawhdNsuhASIiIg1jjwAREWkKrxpQYyFARETawqsGVDg0QEREpGHsESAiIk1RHFc3f2OEChYCRESkLRwaUOHQABERkYaxR4CIiDSFVw2osRAgIiJt4YJCKiwEiIhIU9gjoBZ6hYDwfxaIcMj7CSuwS4kj7HLihDqZPzuIIJwWHELfQogoOIReIUBEROQNrxpQYSFARESawqEBNV4+SEREpGHsESAiIm3hVQMqLASIiEhTODSg5vPQQE1NDRYtWoT4+HgoioJ9+/apnn/44YehKIpqS09P7zduRUUFkpOTodfrkZycjL179/qaGhEREfnI50Lg8uXLSElJwZtvvumxzd13343W1lbnVllZ6TVmbW0tcnNzkZeXh48++gh5eXlYtmwZTp065Wt6RERE3glJW4jweWggJycHOTk5Xtvo9XqYTKYBxywpKcGCBQtQVFQEACgqKsKxY8dQUlKCsrIyX1MkIiLyiEMDagG5aqC6uhoxMTG4/fbb8cMf/hBtbW1e29fW1iI7O1u1b+HChTh58qTHY6xWKywWi2ojIiIi30gvBHJycvCb3/wG7733Hl577TV88MEHuPPOO2G1Wj0eYzabERsbq9oXGxsLs9ns8Zji4mIYjUbnlpCQIO09EBFRCHMIOVuIkH7VQG5urvP/p0yZgunTpyMxMREHDx7E/fff7/E4RVFUj4UQLvuuVVRUhMLCQudji8XCYoCIiPrHlQVVAn75YFxcHBITE9HY2Oixjclkcvn239bW5tJLcC29Xg+9Xi8tTyIi0gYFEuYISMkkOAR8ZcGLFy+ipaUFcXFxHttkZGSgqqpKte/w4cPIzMwMdHpERESa5nOPQFdXF/72t785Hzc1NaGhoQEjR47EyJEjsW7dOixZsgRxcXH49NNP8cwzz2DUqFH4zne+4zwmPz8ft9xyC4qLiwEAq1evxpw5c/DKK69g8eLF2L9/P44cOYITJ05IeItERETX4MqCKj4XAqdPn8a8efOcj/vG6ZcvX44tW7bg448/xq5du9DR0YG4uDjMmzcP5eXliIqKch7T3NwMne7LzojMzEzs3r0bzz77LJ577jmMHz8e5eXlmDVrlj/vjYiIyAUvH1TzuRCYO3cuhJdK6NChQ/3GqK6udtm3dOlSLF261Nd0iIiIyA+81wAREWkLrxpQYSFARESaoggBxc8xfn+PDyYsBNwRDnmh7NJC0WALoX/oRESesBAgIiJtcXyx+RsjRLAQICIiTeHQgFrAFxQiIiKi4MUeASIi0hZeNaDCQoCIiLSFKwuqsBAgIiJN4cqCapwjQEREpGEsBIiISFv6hgb83Xy0efNmJCUlwWAwIC0tDcePH/fa/tixY0hLS4PBYMCtt96KrVu3qp4vLS2FoiguW3d3t095sRAgIiJNURxyNl+Ul5djzZo1WLt2Lerr65GVlYWcnBw0Nze7bd/U1IR77rkHWVlZqK+vxzPPPIOf/OQnqKioULWLjo5Ga2urajMYDD7lxjkCREREAbZ+/XoUFBRgxYoVAICSkhIcOnQIW7ZsQXFxsUv7rVu3YuzYsSgpKQEATJo0CadPn8arr76KJUuWONspigKTyeRXbuwRICIibZE4NGCxWFSb1Wp1eTmbzYa6ujpkZ2er9mdnZ+PkyZNuU6ytrXVpv3DhQpw+fRo9PT3OfV1dXUhMTMSYMWNw3333ob6+3ufTwUKAiIi0RUjaACQkJMBoNDo3d9/u29vbYbfbERsbq9ofGxsLs9nsNkWz2ey2fW9vL9rb2wEAEydORGlpKQ4cOICysjIYDAbMnj0bjY2NPp0ODg0QERHdoJaWFkRHRzsf6/V6j20VRVE9FkK47Ouv/bX709PTkZ6e7nx+9uzZmDZtGjZu3IgNGzYM+D2wECAiIk2Rea+B6OhoVSHgzqhRoxAWFuby7b+trc3lW38fk8nktn14eDhuvvlmt8fodDrMmDHD5x4BDg0QEZG2DPLlg5GRkUhLS0NVVZVqf1VVFTIzM90ek5GR4dL+8OHDmD59OiIiIjy8LYGGhgbExcUNODeAhQAREVHAFRYW4te//jV27NiBM2fO4IknnkBzczNWrlwJACgqKkJ+fr6z/cqVK3H+/HkUFhbizJkz2LFjB7Zv344nn3zS2eaFF17AoUOHcO7cOTQ0NKCgoAANDQ3OmAPFoQEiItIWAcDHdQDcxvBBbm4uLl68iBdffBGtra2YMmUKKisrkZiYCABobW1VrSmQlJSEyspKPPHEE9i0aRPi4+OxYcMG1aWDHR0dePTRR2E2m2E0GpGamoqamhrMnDnTp9wUIfwcKAkSFosFRqMRc7EY4Yr7bhMiIgpOvaIH1diPzs7Ofsfcb1Tf34k7U59GeJhvi+5cr9fejffqXw5ovoOFPQJERKQtAhLuPiglk6DAOQJEREQaxh4BIiLSlhu8aZBLjBDBQoCIiLTFAcDzOj4DjxEiODRARESkYewRICIiTZG5smAoYCFARETawjkCKhwaICIi0jD2CBARkbawR0CFhQAREWkLCwEVDg0QERFpGHsEiIhIW7iOgAoLASIi0hRePqjGQoCIiLSFcwRUfJ4jUFNTg0WLFiE+Ph6KomDfvn2q5xVFcbv94he/8BiztLTU7THd3d0+vyEiIiIaOJ97BC5fvoyUlBT84Ac/wJIlS1yeb21tVT3+3e9+h4KCArdtrxUdHY2zZ8+q9hkM/t0vmoiIyIVDAIqf3+gdodMj4HMhkJOTg5ycHI/Pm0wm1eP9+/dj3rx5uPXWW73GVRTF5VgiIiLpODSgEtDLB//+97/j4MGDKCgo6LdtV1cXEhMTMWbMGNx3332or6/32t5qtcJisag2IiIi8k1AC4G33noLUVFRuP/++722mzhxIkpLS3HgwAGUlZXBYDBg9uzZaGxs9HhMcXExjEajc0tISJCdPhERhSTxZa/AjW5gj8CA7NixA9/73vf6HetPT0/H97//faSkpCArKwu//e1vcfvtt2Pjxo0ejykqKkJnZ6dza2lpkZ0+ERGFIn+LABlDC0EkYJcPHj9+HGfPnkV5ebnPx+p0OsyYMcNrj4Ber4der/cnRSIiIs0LWI/A9u3bkZaWhpSUFJ+PFUKgoaEBcXFxAciMiIg0zSHkbCHC5x6Brq4u/O1vf3M+bmpqQkNDA0aOHImxY8cCACwWC95991289tprbmPk5+fjlltuQXFxMQDghRdeQHp6Om677TZYLBZs2LABDQ0N2LRp0428JyIiIs+E4+rmb4wQ4XMhcPr0acybN8/5uLCwEACwfPlylJaWAgB2794NIQQefPBBtzGam5uh033ZGdHR0YFHH30UZrMZRqMRqampqKmpwcyZM31Nj4iIiHygCBEaMx4sFguMRiPmYjHClYihToeIiHzQK3pQjf3o7OxEdHR0QF6j7+/EXQn/inCdf3PMeh1WHGnZEtB8BwvvNUBERNrikHD5n5bnCBAREX2lcWVBlYCuI0BERETBjT0CRESkLQISegSkZBIUWAgQEZG2cGhAhUMDREREGsYeASIi0haHA4CfCwI5NLygEBER0VcahwZUODRARESkYewRICIibWGPgAoLASIi0hauLKjCoQEiIiINY48AERFpihAOCD9vI+zv8cGEhQAREWmLEP537XOOABER0VeUkDBHIIQKAc4RICIi0jD2CBARkbY4HIDi5xg/5wgQERF9RXFoQIVDA0RERBrGHgEiItIU4XBA+Dk0wMsHiYiIvqo4NKDCoQEiIiINY48AERFpi0MACnsE+rAQICIibRECgL+XD4ZOIcChASIiIg1jjwAREWmKcAgIP4cGRAj1CLAQICIibREO+D80wMsHiYiIvpLYI6DGOQJEREQaFjI9An3VWS96/F4ngoiIBlcvegAMzjftXmH1u2u/L99QEDKFwKVLlwAAJ1A5xJkQEdGNunTpEoxGY0BiR0ZGwmQy4YRZzt8Jk8mEyMhIKbGGkiJCZKDD4XDgs88+Q1RUFBRFcdvGYrEgISEBLS0tiI6OHuQMbxzzHnxf1dyZ9+Bi3vIIIXDp0iXEx8dDpwvcqHV3dzdsNpuUWJGRkTAYDFJiDaWQ6RHQ6XQYM2bMgNpGR0cHzYffF8x78H1Vc2feg4t5yxGonoBrGQyGkPjjLRMnCxIREWkYCwEiIiIN01QhoNfr8fzzz0Ov1w91Kj5h3oPvq5o78x5czJtCQchMFiQiIiLfaapHgIiIiNRYCBAREWkYCwEiIiINYyFARESkYSFXCGzevBlJSUkwGAxIS0vD8ePHvbY/duwY0tLSYDAYcOutt2Lr1q2DlOlVxcXFmDFjBqKiohATE4Nvf/vbOHv2rNdjqquroSiKy/bJJ58MUtbAunXrXF7fZDJ5PWaoz3WfcePGuT1/q1atctt+qM53TU0NFi1ahPj4eCiKgn379qmeF0Jg3bp1iI+Px7BhwzB37lz8+c9/7jduRUUFkpOTodfrkZycjL179w5a3j09PfjZz36GqVOnYsSIEYiPj0d+fj4+++wzrzFLS0vd/gy6u7sHJW8AePjhh11ePz09vd+4Q3m+Abg9b4qi4Be/+IXHmINxvil4hFQhUF5ejjVr1mDt2rWor69HVlYWcnJy0Nzc7LZ9U1MT7rnnHmRlZaG+vh7PPPMMfvKTn6CiomLQcj527BhWrVqF999/H1VVVejt7UV2djYuX77c77Fnz55Fa2urc7vtttsGIeMvTZ48WfX6H3/8sce2wXCu+3zwwQeqvKuqqgAA3/3ud70eN9jn+/Lly0hJScGbb77p9vmf//znWL9+Pd5880188MEHMJlMWLBggfO+G+7U1tYiNzcXeXl5+Oijj5CXl4dly5bh1KlTg5L3lStX8OGHH+K5557Dhx9+iD179uCvf/0rvvWtb/UbNzo6WnX+W1tbpa4Q19/5BoC7775b9fqVld7XrB/q8w3A5Zzt2LEDiqJgyZIlXuMG+nxTEBEhZObMmWLlypWqfRMnThRPP/202/ZPPfWUmDhxomrfj370I5Genh6wHPvT1tYmAIhjx455bHP06FEBQPzzn/8cvMSu8/zzz4uUlJQBtw/Gc91n9erVYvz48cLhcLh9PhjONwCxd+9e52OHwyFMJpN4+eWXnfu6u7uF0WgUW7du9Rhn2bJl4u6771btW7hwoXjggQek5yyEa97u/Pd//7cAIM6fP++xzc6dO4XRaJSbnBfu8l6+fLlYvHixT3GC8XwvXrxY3HnnnV7bDPb5pqEVMj0CNpsNdXV1yM7OVu3Pzs7GyZMn3R5TW1vr0n7hwoU4ffo0enqG5haTnZ2dAICRI0f22zY1NRVxcXGYP38+jh49GujUXDQ2NiI+Ph5JSUl44IEHcO7cOY9tg/FcA1c/N2+//TYeeeQRjzer6jPU5/taTU1NMJvNqnOq1+txxx13ePy8A55/Dt6OCbTOzk4oioKbbrrJa7uuri4kJiZizJgxuO+++1BfXz84CV6juroaMTExuP322/HDH/4QbW1tXtsH2/n++9//joMHD6KgoKDftsFwvmlwhEwh0N7eDrvdjtjYWNX+2NhYmM1mt8eYzWa37Xt7e9He3h6wXD0RQqCwsBDf/OY3MWXKFI/t4uLisG3bNlRUVGDPnj2YMGEC5s+fj5qamkHLddasWdi1axcOHTqEX/3qVzCbzcjMzMTFixfdtg+2c91n37596OjowMMPP+yxTTCc7+v1faZ9+bz3HefrMYHU3d2Np59+Gg899JDXm99MnDgRpaWlOHDgAMrKymAwGDB79mw0NjYOWq45OTn4zW9+g/feew+vvfYaPvjgA9x5552wWq0ejwm28/3WW28hKioK999/v9d2wXC+afCEzN0H+1z/rU4I4fWbnrv27vYPhscffxx//OMfceLECa/tJkyYgAkTJjgfZ2RkoKWlBa+++irmzJkT6DQBXP2l2Gfq1KnIyMjA+PHj8dZbb6GwsNDtMcF0rvts374dOTk5iI+P99gmGM63J75+3m/0mEDo6enBAw88AIfDgc2bN3ttm56erpqYN3v2bEybNg0bN27Ehg0bAp0qACA3N9f5/1OmTMH06dORmJiIgwcPev3DGiznGwB27NiB733ve/2O9QfD+abBEzI9AqNGjUJYWJhLpd3W1uZSkfcxmUxu24eHh+Pmm28OWK7u/PjHP8aBAwdw9OjRAd9O+Vrp6elDWq2PGDECU6dO9ZhDMJ3rPufPn8eRI0ewYsUKn48d6vPdd4WGL5/3vuN8PSYQenp6sGzZMjQ1NaGqqsrnW+HqdDrMmDFjSH8GcXFxSExM9JpDsJxvADh+/DjOnj17Q5/3YDjfFDghUwhERkYiLS3NOQO8T1VVFTIzM90ek5GR4dL+8OHDmD59OiIiIgKW67WEEHj88cexZ88evPfee0hKSrqhOPX19YiLi5Oc3cBZrVacOXPGYw7BcK6vt3PnTsTExODee+/1+dihPt9JSUkwmUyqc2qz2XDs2DGPn3fA88/B2zGy9RUBjY2NOHLkyA0VgkIINDQ0DOnP4OLFi2hpafGaQzCc7z7bt29HWloaUlJSfD42GM43BdBQzVIMhN27d4uIiAixfft28Ze//EWsWbNGjBgxQnz66adCCCGefvppkZeX52x/7tw5MXz4cPHEE0+Iv/zlL2L79u0iIiJC/Od//ueg5fyv//qvwmg0iurqatHa2urcrly54mxzfd6vv/662Lt3r/jrX/8q/vSnP4mnn35aABAVFRWDlvdPf/pTUV1dLc6dOyfef/99cd9994moqKigPtfXstvtYuzYseJnP/uZy3PBcr4vXbok6uvrRX19vQAg1q9fL+rr652z619++WVhNBrFnj17xMcffywefPBBERcXJywWizNGXl6e6qqZP/zhDyIsLEy8/PLL4syZM+Lll18W4eHh4v333x+UvHt6esS3vvUtMWbMGNHQ0KD6zFutVo95r1u3Tvz+978X//M//yPq6+vFD37wAxEeHi5OnTo1KHlfunRJ/PSnPxUnT54UTU1N4ujRoyIjI0PccsstQX2++3R2dorhw4eLLVu2uI0xFOebgkdIFQJCCLFp0yaRmJgoIiMjxbRp01SX4S1fvlzccccdqvbV1dUiNTVVREZGinHjxnn8hxIoANxuO3fu9Jj3K6+8IsaPHy8MBoP4+te/Lr75zW+KgwcPDmreubm5Ii4uTkRERIj4+Hhx//33iz//+c8ecxZi6M/1tQ4dOiQAiLNnz7o8Fyznu++yxeu35cuXCyGuXkL4/PPPC5PJJPR6vZgzZ474+OOPVTHuuOMOZ/s+7777rpgwYYKIiIgQEydOlF7QeMu7qanJ42f+6NGjHvNes2aNGDt2rIiMjBSjR48W2dnZ4uTJk4OW95UrV0R2drYYPXq0iIiIEGPHjhXLly8Xzc3NqhjBdr77/PKXvxTDhg0THR0dbmMMxfmm4MHbEBMREWlYyMwRICIiIt+xECAiItIwFgJEREQaxkKAiIhIw1gIEBERaRgLASIiIg1jIUBERKRhLASIiIg0jIUAERGRhrEQICIi0jAWAkRERBrGQoCIiEjD/j8DFlwiG/OW9QAAAABJRU5ErkJggg==", "text/plain": [ "
" ] diff --git a/src/BTCSv2.cpp b/src/BTCSv2.cpp index b344a38..5209379 100644 --- a/src/BTCSv2.cpp +++ b/src/BTCSv2.cpp @@ -6,11 +6,11 @@ using namespace Eigen; -static SparseMatrix createCoeffMatrix(MatrixXd &alpha, int dim, int rowIndex, double sx) { +static SparseMatrix createCoeffMatrix(MatrixXd &alpha, int numCols, int rowIndex, double sx) { // square matrix of column^2 dimension for the coefficients - SparseMatrix cm(dim, dim); - cm.reserve(VectorXi::Constant(dim, 3)); + SparseMatrix cm(numCols, numCols); + cm.reserve(VectorXi::Constant(numCols, 3)); // left column cm.insert(0,0) = 1 + sx * (calcAlphaIntercell(alpha(rowIndex,0), alpha(rowIndex,1)) @@ -18,13 +18,14 @@ static SparseMatrix createCoeffMatrix(MatrixXd &alpha, int dim, int rowI cm.insert(0,1) = -sx * calcAlphaIntercell(alpha(rowIndex,0), alpha(rowIndex,1)); // inner columns - int n = dim-1; + int n = numCols-1; for (int i = 1; i < n; i++) { cm.insert(i,i-1) = -sx * calcAlphaIntercell(alpha(rowIndex,i-1), alpha(rowIndex,i)); cm.insert(i,i) = 1 + sx * ( calcAlphaIntercell(alpha(rowIndex,i), alpha(rowIndex,i+1)) - + calcAlphaIntercell(alpha(rowIndex,i-1), alpha(i,1)) - ); + + calcAlphaIntercell(alpha(rowIndex,i-1), alpha(rowIndex,i)) + ) + ; cm.insert(i,i+1) = -sx * calcAlphaIntercell(alpha(rowIndex,i), alpha(rowIndex,i+1)); } @@ -96,21 +97,15 @@ static VectorXd createSolutionVector(MatrixXd &concentrations, MatrixXd &alphaX, } // first column -> additional fixed concentration change from perpendicular dimension - sv(0) += 2 * sx * alphaX(rowIndex,0) * bcLeft(0); + sv(0) += 2 * sx * alphaX(rowIndex,0) * bcLeft(rowIndex); // last column -> additional fixed concentration change from perpendicular dimension - sv(length-1) += 2 * sx * alphaX(rowIndex,length-1) * bcRight(length-1); + sv(length-1) += 2 * sx * alphaX(rowIndex,length-1) * bcRight(rowIndex); return sv; } -// BTCS solution for 1D grid -static void BTCS_1D(Grid &grid, Boundary &bc, double ×tep) { - // TODO -} - - static VectorXd solve(SparseMatrix &A, VectorXd &b) { SparseLU> solver; solver.analyzePattern(A); @@ -120,6 +115,38 @@ static VectorXd solve(SparseMatrix &A, VectorXd &b) { } +// BTCS solution for 1D grid +static void BTCS_1D(Grid &grid, Boundary &bc, double ×tep) { + int length = grid.getLength(); + double sx = timestep / (grid.getDeltaCol() * grid.getDeltaCol()); // TODO create method getDelta() for 1D case + + VectorXd concentrations_t1(length); + + SparseMatrix A; + VectorXd b(length); + + MatrixXd alpha = grid.getAlpha(); + VectorXd bcLeft = bc.getBoundarySideValues(BC_SIDE_LEFT); + VectorXd bcRight = bc.getBoundarySideValues(BC_SIDE_RIGHT); + + MatrixXd concentrations = grid.getConcentrations(); + A = createCoeffMatrix(alpha, length, 0, sx); // this is exactly same as in 2D + for (int i = 0; i < length; i++) { + b(i) = concentrations(0,i); // TODO check if this is really the only thing on right hand side of equation in 1D? + } + b(0) += 2 * sx * alpha(0,0) * bcLeft(0); + b(length-1) += 2 * sx * alpha(0,length-1) * bcRight(0); + + concentrations_t1 = solve(A, b); + + for (int j = 0; j < length; j++) { + concentrations(0,j) = concentrations_t1(j); + } + + grid.setConcentrations(concentrations); +} + + // BTCS solution for 2D grid static void BTCS_2D(Grid &grid, Boundary &bc, double ×tep) { int rowMax = grid.getRow();