Skip to content

Instantly share code, notes, and snippets.

@jthorton
Last active October 28, 2021 13:40
Show Gist options
  • Save jthorton/20d23fe0c682413158f40a714ee183bf to your computer and use it in GitHub Desktop.
Save jthorton/20d23fe0c682413158f40a714ee183bf to your computer and use it in GitHub Desktop.
How to train your force field follow along notebook.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"id": "28ae309a",
"metadata": {},
"source": [
"General transferable force fields, such as [Parsley](https://pubs.acs.org/doi/full/10.1021/acs.jctc.1c00571) or [Sage](https://openforcefield.org/force-fields/force-fields/),\n",
"allow us to rapidly parameterize small molecules covering vast amounts of chemical space for use in molecular dynamics simulations. However, with every force field comes the inherent assumption that\n",
"parameters fit to a set of representative molecules in a training set are transferable to any new molecules with similar chemistry.\n",
"Despite potential accuracy problems caused by poor transferability or discrete atom types the use \n",
"of force fields in fields like drug discovery has been a massive success. Within the Open Force Field, we\n",
"use a unique method based on [SMARTS](https://www.daylight.com/dayhtml/doc/theory/theory.smarts.html) patterns to link\n",
"the force field parameters to chemical substructures, avoiding the use of atom types altogether, resulting in a more\n",
"[compact and manageable force field](https://doi.org/10.1021/acs.jctc.8b00640). While we have seen the number of unique\n",
"torsion parameters grow between versions of the force field to improve parameter performance the number of parameters\n",
"has remained very low when compared to other state-of-the-art force fields like OPLS3.\n",
"\n",
"\n",
"| Force Field | Number of unique torsion parameters |\n",
"|---|:----------------:|\n",
"|OpenFF-1.0.0| 157|\n",
"|OpenFF-1.2.0| 163|\n",
"|OpenFF-1.3.0| 167|\n",
"|OpenFF-2.0.0| 167|\n",
"|OPLS3|48,142|\n",
"|OPLSe|146,669|\n",
"\n",
"\n",
"This explosion in torsion parameters allows for high accuracy over a large area of chemical space with the types often becoming\n",
"very specific to certain chemistry. However, due to improved computer power and the reduced computational cost of generating\n",
"high accuracy reference data via the use of machine learned potentials such as [ANI](https://pubs.rsc.org/en/content/articlelanding/2017/SC/C6SC05720A),\n",
"a better solution may be to generate bespoke parameters on the fly. Enter :tada:[BespokeFit](https://github.com/openforcefield/bespoke-fit) :tada:\n",
", BespokeFit is a python package which allows users to easily generate bespoke force field parameters for\n",
"their molecules under study which compliment the base OpenFF force field. Here we can think of BespokeFit as our general fitting\n",
"scheme which can be applied to new molecules to generate bespoke parameters on the fly for new unseen chemistry.\n",
"\n",
"You may have seen [previously](https://openforcefield.org/community/news/science-updates/ff-training-example-2021-07-01/)\n",
"how BespokeFit, QCSubmit and QCFractal can be combined to fit general force fields and this time we are going to show\n",
"its other side and walk through an application to a simple molecule.\n",
"\n",
"This will have the following structure:\n",
"- Building the general BespokeFit optimisation workflow\n",
"- Building and inspecting the molecule specific optimisation schema made from the general workflow\n",
"- Setup and run the BespokeExecutor\n",
"- Inspecting the fitted parameters\n",
"\n",
"## Building the general workflow\n",
"\n",
"BespokeFit aims to provide a reproducible parameter optimization workflow for SMIRNOFF based force fields.\n",
"As such normal BespokeFit execution starts with a general fitting workflow. This captures every process in the\n",
"workflow along with any adjustable settings such as how the reference data should be generated. The default workflow\n",
"is designed to optimize bespoke torsion parameters at the same level of theory as that used in the mainline openff force fields.\n",
"Here however, we will be starting with a blank workflow and building it up as we go.\n",
"\n",
"First, we start with our optimization engine, which is an easy choice as we currently only support the fantastic\n",
"[ForceBalance](https://github.com/leeping/forcebalance) package (watch this space - optimizers are coming!).\n",
"Let's start by creating the ForceBalance optimization schema."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "63eeb695",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'type': 'ForceBalance',\n",
" 'max_iterations': 10,\n",
" 'job_type': 'optimize',\n",
" 'penalty_type': 'L2',\n",
" 'step_convergence_threshold': 0.01,\n",
" 'objective_convergence_threshold': 0.01,\n",
" 'gradient_convergence_threshold': 0.01,\n",
" 'n_criteria': 2,\n",
" 'eigenvalue_lower_bound': 0.01,\n",
" 'finite_difference_h': 0.01,\n",
" 'penalty_additive': 1.0,\n",
" 'initial_trust_radius': -0.25,\n",
" 'minimum_trust_radius': 0.05,\n",
" 'error_tolerance': 1.0,\n",
" 'adaptive_factor': 0.2,\n",
" 'adaptive_damping': 1.0,\n",
" 'normalize_weights': False,\n",
" 'extras': {}}"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from openff.bespokefit.schema.optimizers import ForceBalanceSchema\n",
"\n",
"fb = ForceBalanceSchema()\n",
"fb.dict()"
]
},
{
"cell_type": "markdown",
"id": "e2e7bdd3",
"metadata": {},
"source": [
"As you can tell there are many options here and in some cases, it might not be clear what a valid input is. For example,\n",
"what other penalty types could we use? :sparkles:[Pydantic](https://github.com/samuelcolvin/pydantic):sparkles: to the rescue, pydantic\n",
"allows us to \"_define how data should be in pure, canonical python_\" and has run time validation to ensure our data is always correct.\n",
"What's more we even get a `schema` method for each model for free which fully describes the model with a short description for each field\n",
"and information on the acceptable inputs, it basically comes with its own documentation:book:. As such pydantic is used\n",
"extensively throughout BespokeFit to try and catch possible input errors ahead of run time, there is nothing worse than queuing up\n",
"a calculation for it to be instantly returned with an error on line one:man_facepalming:. Now the schema is\n",
"validated during assignment as we see below. First lets start by inspecting the optimizer schema using:"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "4e19a4aa",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'title': 'ForceBalanceSchema',\n",
" 'description': 'A class containing the main ForceBalance optimizer settings to use during an\\noptimization.\\n\\nPriors and target definitions are stored separately as part of an\\n``OptimizationSchema``.',\n",
" 'type': 'object',\n",
" 'properties': {'type': {'title': 'Type',\n",
" 'default': 'ForceBalance',\n",
" 'enum': ['ForceBalance'],\n",
" 'type': 'string'},\n",
" 'max_iterations': {'title': 'Max Iterations',\n",
" 'description': 'The maximum number of optimization iterations to perform.',\n",
" 'default': 10,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'integer'},\n",
" 'job_type': {'title': 'Job Type',\n",
" 'description': 'The mode to run force balance in.',\n",
" 'default': 'optimize',\n",
" 'enum': ['optimize'],\n",
" 'type': 'string'},\n",
" 'penalty_type': {'title': 'Penalty Type',\n",
" 'description': 'The penalty type.',\n",
" 'default': 'L2',\n",
" 'enum': ['L1', 'L2'],\n",
" 'type': 'string'},\n",
" 'step_convergence_threshold': {'title': 'Step Convergence Threshold',\n",
" 'description': 'The step size convergence criterion.',\n",
" 'default': 0.01,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'objective_convergence_threshold': {'title': 'Objective Convergence Threshold',\n",
" 'description': 'The objective function convergence criterion.',\n",
" 'default': 0.01,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'gradient_convergence_threshold': {'title': 'Gradient Convergence Threshold',\n",
" 'description': 'The gradient norm convergence criterion.',\n",
" 'default': 0.01,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'n_criteria': {'title': 'N Criteria',\n",
" 'description': 'The number of convergence thresholds that must be met for convergence.',\n",
" 'default': 2,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'integer'},\n",
" 'eigenvalue_lower_bound': {'title': 'Eigenvalue Lower Bound',\n",
" 'description': 'The minimum eigenvalue for applying steepest descent correction.',\n",
" 'default': 0.01,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'finite_difference_h': {'title': 'Finite Difference H',\n",
" 'description': 'The step size for finite difference derivatives in many functions.',\n",
" 'default': 0.01,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'penalty_additive': {'title': 'Penalty Additive',\n",
" 'description': 'The factor for the multiplicative penalty function in the objective function.',\n",
" 'default': 1.0,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'initial_trust_radius': {'title': 'Initial Trust Radius',\n",
" 'description': \"The initial value of the optimizers adaptive trust radius which 'adapts' (i.e. increases or decreases) based on whether the last step was a good or bad step.\",\n",
" 'default': -0.25,\n",
" 'type': 'number'},\n",
" 'minimum_trust_radius': {'title': 'Minimum Trust Radius',\n",
" 'description': 'The minimum value of the optimizers adaptive trust radius.',\n",
" 'default': 0.05,\n",
" 'type': 'number'},\n",
" 'error_tolerance': {'title': 'Error Tolerance',\n",
" 'description': 'Steps that increase the objective function by more than this will be rejected.',\n",
" 'default': 1.0,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'adaptive_factor': {'title': 'Adaptive Factor',\n",
" 'description': 'The amount to change the step size by in the event of a good / bad step.',\n",
" 'default': 0.2,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'adaptive_damping': {'title': 'Adaptive Damping',\n",
" 'description': 'A damping factor that restraints the trust radius to trust0.',\n",
" 'default': 1.0,\n",
" 'exclusiveMinimum': 0,\n",
" 'type': 'number'},\n",
" 'normalize_weights': {'title': 'Normalize Weights',\n",
" 'description': 'Whether to normalize the weights for the fitting targets',\n",
" 'default': False,\n",
" 'type': 'boolean'},\n",
" 'extras': {'title': 'Extras',\n",
" 'description': 'Extra settings (mostly logging settings) to include in the ForceBalance input file.',\n",
" 'default': {},\n",
" 'type': 'object',\n",
" 'additionalProperties': {'type': 'string'}}}}"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"fb.schema()"
]
},
{
"cell_type": "markdown",
"id": "dc8c9bdc",
"metadata": {},
"source": [
"We can now see that `penalty_type` only accepts the values L2 and L1 so lets try something else and see what happens.\n"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "dfd324e1",
"metadata": {},
"outputs": [
{
"ename": "ValidationError",
"evalue": "1 validation error for ForceBalanceSchema\npenalty_type\n unexpected value; permitted: 'L1', 'L2' (type=value_error.const; given=L3; permitted=('L1', 'L2'))",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mValidationError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m/var/folders/9q/nm__l0v13fggc72v94p0hybw0000gq/T/ipykernel_60581/1577149853.py\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mfb\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mpenalty_type\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m\"L3\"\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;32m~/miniconda3/envs/bespokefit_new/lib/python3.9/site-packages/pydantic/main.cpython-39-darwin.so\u001b[0m in \u001b[0;36mpydantic.main.BaseModel.__setattr__\u001b[0;34m()\u001b[0m\n",
"\u001b[0;31mValidationError\u001b[0m: 1 validation error for ForceBalanceSchema\npenalty_type\n unexpected value; permitted: 'L1', 'L2' (type=value_error.const; given=L3; permitted=('L1', 'L2'))"
]
}
],
"source": [
"fb.penalty_type = \"L3\""
]
},
{
"cell_type": "markdown",
"id": "e6190353",
"metadata": {},
"source": [
"Nice, we get an informative error message before even running the workflow!\n"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "1c087d13",
"metadata": {},
"outputs": [],
"source": [
"fb.penalty_type = \"L1\""
]
},
{
"cell_type": "markdown",
"id": "f8fc6983",
"metadata": {},
"source": [
"Now we can set it to L1 and start to build up our workflow and see what other pieces we might need.\n",
"\n",
"{{< note >}}\n",
"By default the workflow model comes ready to fit bespoke torsion parameters using Open Force Field best practices and\n",
"default settings. Meaning only users who want absolute control over every setting (or developers showing of their work :grin:)\n",
"would need to do this manual set up .\n",
"{{< /note >}}\n",
"\n",
"\n",
"Lets start by creating our `BespokeWorkflowFactory` which acts as a fitting template for future optimizations\n",
"and add our `ForceBalanceSchema` optimizer settings. We can also remove all other defaults as we will go through them in\n",
"the next sections.\n"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "69875e28",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n",
"Warning: Unable to load toolkit 'AmberTools'. \n"
]
},
{
"data": {
"text/plain": [
"{'initial_force_field': 'openff_unconstrained-1.3.0.offxml',\n",
" 'optimizer': {'type': 'ForceBalance',\n",
" 'max_iterations': 10,\n",
" 'job_type': 'optimize',\n",
" 'penalty_type': 'L1',\n",
" 'step_convergence_threshold': 0.01,\n",
" 'objective_convergence_threshold': 0.01,\n",
" 'gradient_convergence_threshold': 0.01,\n",
" 'n_criteria': 2,\n",
" 'eigenvalue_lower_bound': 0.01,\n",
" 'finite_difference_h': 0.01,\n",
" 'penalty_additive': 1.0,\n",
" 'initial_trust_radius': -0.25,\n",
" 'minimum_trust_radius': 0.05,\n",
" 'error_tolerance': 1.0,\n",
" 'adaptive_factor': 0.2,\n",
" 'adaptive_damping': 1.0,\n",
" 'normalize_weights': False,\n",
" 'extras': {}},\n",
" 'target_templates': [],\n",
" 'parameter_hyperparameters': [],\n",
" 'target_torsion_smirks': ['[!#1]~[!$(*#*)&!D1:1]-,=;!@[!$(*#*)&!D1:2]~[!#1]'],\n",
" 'expand_torsion_terms': True,\n",
" 'generate_bespoke_terms': True,\n",
" 'fragmentation_engine': None,\n",
" 'default_qc_specs': []}"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from openff.bespokefit.workflows import BespokeWorkflowFactory\n",
"\n",
"workflow = BespokeWorkflowFactory(\n",
" fragmentation_engine=None, \n",
" optimizer=fb, \n",
" parameter_hyperparameters=[], \n",
" target_templates=[], \n",
" default_qc_specs=[]\n",
")\n",
"workflow.dict()"
]
},
{
"cell_type": "markdown",
"id": "ef169733",
"metadata": {},
"source": [
"### Stage 1 Fragmentation\n",
"\n",
"BespokeFit makes extensive use of the fantastic [openff-fragmenter](https://github.com/openforcefield/openff-fragmenter)\n",
"package where possible to reduce the cost of QM torsion drives as they can quickly become very expensive for large\n",
"drug-like molecules. Here we will add the WBO based fragmentation scheme to the pipeline with default settings. You can find\n",
"out more about how fragmenter works in a pre-print [here](https://www.biorxiv.org/content/10.1101/2020.08.27.270934v1)."
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "47583377",
"metadata": {},
"outputs": [],
"source": [
"from openff.fragmenter.fragment import WBOFragmenter\n",
"\n",
"fragmenter = WBOFragmenter()\n",
"workflow.fragmentation_engine = fragmenter"
]
},
{
"cell_type": "markdown",
"id": "45516d4d",
"metadata": {},
"source": [
"By default, fragmenter will not fragment terminal rotatable bonds (such as methyl groups) as these are assumed to be\n",
"trivial and well described by our chosen force field. However, to keep things simple in our example we will be working with\n",
"`BrCO` and so we need to change this behaviour. Luckily BespokeFit/fragmenter allow us to overwrite this behaviour\n",
"by defining a SMARTS pattern (chemical substructure) which will be used to select the rotatable bonds."
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "4a6f879d",
"metadata": {},
"outputs": [],
"source": [
"workflow.target_torsion_smirks = [\"[*]~[!$(*#*)&!D1:1]-,=;!@[!$(*#*)&!D1:2]~[*]\"]\n",
"\n"
]
},
{
"cell_type": "markdown",
"id": "f5aaf609",
"metadata": {},
"source": [
"Here we have taken the default SMARTS pattern `[!#1]~[!$(*#*)&!D1:1]-,=;!@[!$(*#*)&!D1:2]~[!#1]` used to find torsions\n",
"and allowed the connecting atoms to the central bond to be hydrogens. This will ensure our simple molecule is processed\n",
"by the workflow, as without this change the workflow would reject the molecule as with no `rotatable` bonds there would be\n",
"no work to do :man_shrugging:.\n",
"\n",
"\n",
"### Stage 2 Reference Generation\n",
"\n",
"As you may have guessed we are going to need some target reference data to optimize our force field parameters against\n",
"and that's where the target templates come in. These `target schemas` have two main functions:\n",
"- Describe the contribution to the total error function in the parameter optimization\n",
"- Generate specific reference data tasks that need to be computed in order to train with the target\n",
"\n",
"In this case we will be using the `TorsionProfileTargetSchema` which requires a torsiondrive as input data. That is\n",
"a series of constrained geometry optimizations around the targeted rotatable bond, usually ranging from -180 -> 180 degrees\n",
"in 15 degree increments. There is slightly more to a torsiondrive than that and you can read more about how it works in the\n",
"[paper](https://doi.org/10.1063/5.0009232), but for now the main point is that we should get back a series of geometries\n",
"and energies to fit to. The actual torsion profile target then contributes the average RMSE between the QM and MM\n",
"energies at each geometry to the objective function and has a small number of settings exposed which we can see here."
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "493a380f",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{'weight': 1.0,\n",
" 'reference_data': None,\n",
" 'extras': {},\n",
" 'type': 'TorsionProfile',\n",
" 'attenuate_weights': True,\n",
" 'energy_denominator': 1.0,\n",
" 'energy_cutoff': 10.0}"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from openff.bespokefit.schema.targets import TorsionProfileTargetSchema\n",
"\n",
"target = TorsionProfileTargetSchema()\n",
"workflow.target_templates = [target, ]\n",
"target.dict()"
]
},
{
"cell_type": "markdown",
"id": "2e8df449",
"metadata": {},
"source": [
"You will also note that we passed the target in a list as the optimizer can use multiple targets to construct the objective\n",
"function, you can check out [Simons blog](https://openforcefield.org/community/news/science-updates/ff-training-example-2021-07-01/)\n",
"to see this in action.\n",
"\n",
"{{< warning >}}\n",
"All target templates passed to the optimizer will be exercised at the same time as multi-stage fits are not yet supported.\n",
"{{< /warning >}}\n",
"\n",
"We now need to tell BespokeFit what program, method and basis should be used to actually compute the reference tasks. As\n",
"we use the game changing [QCEngine](https://github.com/MolSSI/QCEngine) to power these calculations we have a wide range\n",
"of possibilities with one standard interface. BespokeFit defines the target specification using the `QCSpec` class from\n",
"QCSubmit which has built in validation to again catch any errors early in the process. Choosing a new QC task\n",
"specification for all tasks in the BespokeFit workflow is then as simple as:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "fb61da9a",
"metadata": {},
"outputs": [],
"source": [
"from openff.qcsubmit.common_structures import QCSpec\n",
"\n",
"xtb_spec = QCSpec(\n",
" method=\"gfn2xtb\", \n",
" basis=None, \n",
" program=\"xtb\", \n",
" spec_name=\"xtb\", \n",
" spec_description=\"gfn2xtb\"\n",
")\n",
"\n",
"workflow.default_qc_specs = [xtb_spec]"
]
},
{
"cell_type": "markdown",
"id": "31fd8d4e",
"metadata": {},
"source": [
"To keep things moving quickly we have chosen to use a semi-empirical QC method, but by default BespokeFit will use the\n",
"same level of theory as that used to fit the main line force fields.\n",
"\n",
"\n",
"### Stage 3 Parameter Optimization\n",
"\n",
"Finally, we have the parameter optimization stage and as we have already set up our optimizer using `ForceBalance` the last step\n",
"is to define some parameter hyperparameters. These tell BespokeFit which parameters we would like to fit and allows us to\n",
"define a prior to the parameter."
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "36422a09",
"metadata": {},
"outputs": [],
"source": [
"from openff.bespokefit.schema.smirnoff import ProperTorsionHyperparameters\n",
"\n",
"prior = ProperTorsionHyperparameters()\n",
"workflow.parameter_hyperparameters = [prior]"
]
},
{
"cell_type": "markdown",
"id": "9bc8da1d",
"metadata": {},
"source": [
"BespokeFit now knows we intend to optimize the torsion parameters to the target and with these final two settings we\n",
"can generate SMARTS patterns specific to the molecules that pass through the workflow and to fully expand\n",
"the torsion parameters to use all available k values."
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "ae6d5d03",
"metadata": {},
"outputs": [],
"source": [
"workflow.generate_bespoke_terms = True\n",
"workflow.expand_torsion_terms = True"
]
},
{
"cell_type": "markdown",
"id": "6d4b7bb3",
"metadata": {},
"source": [
"All that is left to do now is to save the workflow to file for later use, this is now our general workflow template which\n",
"we can apply to all molecules that pass into the BespokeFit pipeline. The serialized workflow is also a fantastic provenance\n",
"store, making it easy to tell exactly what settings were used to compute the parameters in a project."
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "7bee7f91",
"metadata": {},
"outputs": [],
"source": [
"# write to json\n",
"workflow.export_factory(\"workflow.json\")\n",
"# or yaml\n",
"workflow.export_factory(\"workflow.yaml\")"
]
},
{
"cell_type": "markdown",
"id": "43f67feb",
"metadata": {},
"source": [
"## Building the Molecule Specific Schema\n",
"\n",
"Now we can use our general fitting schema to build a molecule specific optimization schema.\n",
"This schema fully defines the optimization protocol that should be applied to the molecule including information\n",
"about the reference data generation tasks which BespokeFit will automatically perform locally for us. In this\n",
"example we will be using `BrCO`, so lets create the molecule using the `openff-toolkit`, save it to file for later and\n",
"create a `BespokeOptimizationSchema` for it using our workflow:"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "9927fd31",
"metadata": {},
"outputs": [],
"source": [
"from openff.toolkit.topology import Molecule\n",
"\n",
"target_molecule = Molecule.from_smiles(\"BrCO\")\n",
"target_molecule.generate_conformers(n_conformers=1)\n",
"# write to file for later\n",
"target_molecule.to_file(file_path=\"BrCO.sdf\", file_format=\"sdf\")\n",
"# make the specific schema\n",
"schema = workflow.optimization_schema_from_molecule(molecule=target_molecule)"
]
},
{
"cell_type": "markdown",
"id": "d111bcbd",
"metadata": {},
"source": [
"If we now inspect the schema we find two main changes, i) the input molecule has been inserted, ii) some bespoke SMARTS patterns\n",
"have been created for the rotatable torsion identified in the molecule. These patterns and the initial force field values\n",
"can be easily viewed via:"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "769ba025",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#35H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#1H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.6985935181583, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=0, unit=kilocalorie/mole)},\n",
" ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#1H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#35H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.6985935181583, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=0, unit=kilocalorie/mole)}}"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"schema.initial_parameter_values"
]
},
{
"cell_type": "markdown",
"id": "d3d9f94b",
"metadata": {},
"source": [
"You may also notice that as requested BespokeFit has expanded the k terms for each torsion parameter with each new term\n",
"set to zero by default. The non-zero values are then taken from the base force field to form our initial parameter values.\n",
"We can also use the `openff-toolkit` to query the molecule and identify the atoms matching these new patterns and\n",
"understand how the parameters will be applied to the molecule by visualizing the substructures."
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "c8cfc73f",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAADICAIAAAC7/QjhAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dd1xUZ/Y/8M90OtKrRjGiogIqiAUVG2U1RrNqojFm0zXRmLbJZrPFzf5MW5NoNIlrsu5XY0nQxGxs2EFFEVAEMiCgiErvMAxl2v39MQSGYYY6zDBzz/uVP3bvvTNzhjyZc8vznMNhGAaEEEIIW3FNHQAhhBBiSpQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgirUSIkhBDCapQICSGEsBolQkIIIaxGiZAQQgir8U0dACEWTaHAvXt48ADl5WhpAY8HV1f4+sLPD1ZWpg6OEAIAHIZhTB0DIZZIKsXFi0hLA5cLuRwqVfsukQhKJUaOxLx5cHMzXYiEEIASISED4sYNxMVBpYJSqfcYDgc8HkJCsGABuPSQghCToURIiEExDI4ehVgMmaxHxwsEcHfHU09BJBrgyAghutF5KFauXBkeHn78+HGdewsKCsLDw8PDwyUSiZEDI2bpxAn8+itkspWHD4fv3n08N1fnUQW1teG7d4fv3i1paEBZGfbu7erakRAykGiyDFJTU2/fvl1aWqpzr1QqTUxMBCCXy40bFzFDmZlIT4dcDiC1uPh2dXVpQ4POA6UyWeL9+wDkKhUUClRUIC4OCxcaNVpCCAC6IiTEYFpacPw4+nbCJJcjPR0lJYaOiRDSPUqEhBhIUlKHqaG9JZfj9GnDRUMI6SlKhIQYAsMgObmPl4NtCgtBj6IJMTpKhIQYQnk5FIr+vgmXi7w8Q0RDCOkFSoSEGEJxMfq/EkkmQ0GBAYIhhPQGzRpttWXLlv3793feLpVKjR8MMT81NTrvi265cmV/Zmbn7VJ9qwyrqgwbFyGkW5QIW926devWrVumjoKYrZYWnZtvVVbeqqzsxfvQKh1CjI4SYauPP/54xYoVnbfn5eVFRkYaPx5iZvTUhfl4wYIV48Z13p5XVRX53Xc6XiAUGjYuQki3KBG2cnFxGT58eOftdGuU9IiTE4TCzmXVXKythw8Z0vlwvbdGnZ0NHhohpGs0WYYQQ/DxMcCbCATQdTZGCBlQlAgJMQQ3N/D7fX+FYTBqlCGiIYT0AiVCQgyBw0FYGASCfr3J0KGwtzdQQISQnqJESIiBhIX1t61gQICBQiGE9AIlQkIMRCTCokX9uihMTu5XtVJCSJ/QrFGEhYX5+Ph4eXnp3Gtrazt79mwAgn7e9SJsMH487t9HejpksjBfXx8HBy89tzpthcLZw4cDEGheRFZU4Pp1hIYaJVZCSCvqUN8jEomEw+HY2dmZOpABI5OhqgoSCRQKcDiwtYWjIxwdTR2WGWIYnDiBjIyedKiXtLRwOBw7zbWD1tbYsAHW1gMYISGkI0qE3Tt58uRzzz23evXqTz75xNSxGFpJCdLTkZ2NhgYIBOBwoFKBywXDQKkEl4thwxAYiLFjDTAlklVu3sSJE2CYLipxn8zLe+6XX1YHBn6yYEGHHWFhiI4e8AgJIb+hRNi9tLS0kJAQPp+fmZnp7+9v6nAMJD8fp0+juhpKZTfPpYRCcDiYOhUzZvR3ViSrNDbi4kWkpYHDgVze+Y+cVlISsmsXn8vNfPllfxeX9h1cLl56Ce7uRo2WEBbjbdq0ydQxDHZeXl73799PTU29f//+ypUrTR1OvzU34/BhXLqE+nqoVN33TFAqoVSiuBipqfD2hq46KUQHgQAPP4zp0zF0KJycYGODpibN+6Ve9vb36+pSi4vv19WtnDCh/YUMg6oqBAWZIGZCWImuCHukvLzc39+/rq7u5MmT0WZ926q6uubrrflltzUvUPhcnrPI0dvGlcfldfNygQBz52Lq1IEN0lJJJNixQzMXlkul/tu31zU3n1y9OvrhhzscvHIlzOr2Q03pvfyblzS38IUiZ68R3g8H8gRUQJUMapQIe+qTTz555513xo4dm56ebq4zSGtrsWuXuER8+PaZzjut+aJZXpOmegR28yYCAebMwbRpAxKhxbt4ERcuaG74JDHxnTNnxrq5pa9dK+BpnIg4OeGVV8Dr7tRk0BBf/uXwx2s7b7e2HzLr8demPvqi8UMipIfo1mhPhYaGxsbG5uTkuLq6hoWFmTqc3lMq8e23aGioaKrOqsm34gkjvEP8HHz9HHw9bFxkKnmdrOFOfaGAKxhm59nV+6hUuH8fvr5wcjJW6BbExweZmWhubtsQ6uMTKxbnVFa62tqG+fq2H9ncDJEIQ4eaIMg+qbifk5V4zMrWIWLVW37Bs/yCZ3mMCJA1SevKC+/ciBeIrIcF0LIQMkjRgvqeEgqF//rXvwBs2rSpslcd5gaJhAQ0NLQ9ERTxRTO8gtX/RA6dtjZgWZCLP4BLpWlKprs13XI5Dh/uYj4k0YvPR8c5okIe71+RkQA2xcdXNjZ2OFj9r8ysiGwdZix7Rf1P5HN/X7v9XNDc5QAuHfpCqaBWi2SQokTYC4sXL46KiqqpqTG/y+imJiQlddH0lcPhhLqPA9CiaKluqet8gFKllKk0Xi6XIyVlAAJlgYAArRYTi0ePjnr44Zqmpk3x8R2OlMm07qOaHQ6HE7rwDwBapPXVJXc7H6CUy2TN1OmMmBglwt75/PPPBQLBzp07MzIyTB1Lb/QgWtVvF4ICDg9Avazhi8yDvxTEN8gbY++c3py2+8Mbu/93N771aLkcSUkDFa3Fi44Gh6O54fOoKAGPtzM1NaOsrMORaWkoLjZqbIamUrbeORAIrQDUVxZ/8cK0X754o6GmPPbD5zcvG/nh8lH/2/q6SWMkbEeJsHfGjh27du1apVL5+utm9Z+uWNzF5aDa9YpsAC5WQxxF9gCUjKqmpT6/vmhv7vFbNQV2AmshV+Ao/K22DsOgqQk1NQMct4Xy8MCkSZobxrq5rQ0JUapUr8fFdTiSYRAX1/0Sl0Hsetw+AC4+Ix3dfQEoFXL1/NK9f3n81tU4Oyc3oZWto5shujkS0ldULqTX/vGPfxw8ePD8+fNHjhxZunSpqcPpgepqFBVpbZMp5erMB6BR0ZRTd7+oocyGb/XYiLkctF+s1MkaPGxc1o9/3NnKEYBSpWzdweGAy0VZGU2Z6aO5cyEWa86a+UdExMHMzPN37x7Jzl46dmz7kQ8eQCzG+PEmCLL3ZE1SdeYD0FhflZN8pijnho2D82Nv7uBw2k+76yqKPEYErN95ydl7BAClvPtydIQMHEqEvebk5LRp06b169e/9dZbMTExVlZWpo5IP7kciYm4fLlzWZMmRfOxexc1twi5gkeHR3jbumkducxvvjoLAuiw0FClgpSe7vSVjQ1mzcLp020bnKytN0VErD9x4q3Tp2NGjbLSrGl35gxGjzaLsj5NkppjX76tuUVoZfvoxs+9R2nXB1j29k51FgRACw2JadGt0b5Yu3bthAkT8vPzt27daupY9GAYpKfjiy+QkAClsvN+K55wvk/YfJ+web5hYR4TRjr4ylWKg7fjjty9oFJ1ON6er78AdD/b77FcWBjcOpx2rA0JmeDhkV9Ts1Xr+Wt9Pa5cMWpsfWVl6zD/6ffmP/3evKf/HLb4+ZETZ8tbmg7+8+kjn72q6jhr1N6JasiRwYKuCPuCx+Nt3bp13rx5mzdvXrNmjbe3t6kj6qioCHFxKCzs4hD18okOL2oo25t7PKMq18vGdarHBH0vbMfjUXuKfuFyER2N775r28DjcrdGR8/bs2fzxYtrgoK8NVs4Xb6M4ODB/wdXL5/Q3FKUc2PvXx7PuHDYa+SEqY++YKrACOkCndH30dy5c5csWdLQ0PDee++ZOhYNEgmOHMG333adBXXysfOY4DIKQFZNfo9eIJdjsJ0BmB0/P4wapblh7ogRS8aMaZDJ3jt3rsORCgXOnjVqbAbiM3rShIjHAGQlHjN1LIToRomw7z799FMrK6s9e/YkJyebOpbfHgfu2KFvpUSBpHhf7vEWZVezEuz41gAaFU09+kQvLwzm56PmIjpaq47ap1FRVnz+nps3k7WmOP36K+7dM2psPVCQeWXf31a2NEq6OMZuiBuAxvoqYwVFSO9QIuw7Pz+/1157jWGYjRs3mrhka24uvvoKZ8/qbAZbK5McunNmT87RO/WF18oyu3ibQmkZgCEih+4/USDA7Nl9DZdocHbWaknv5+T02tSpDLDx5EntcTWYllLUlj849NGLe/687E5awrVfvu3iyMKc6wCGeAwzVmiE9A4lwn557733vLy8kpKSDhw4YJoISkqwezcOHkRtbeedcpUivvj6l7/+kFWTL+DyZ3uHTPfU3dxHyaguldy4U18IYILTyG4+lMuFpye0uiWQPouIgJ2d5ob3Zs3ysrdPKiw8kNnxxKW0FDdvGjU2XeQtTfEHPv1y3aysxGMCkfXslW9Of2ydziOVCvml2G130hIATJi1xLhhEtJTNFmmX+zs7DZv3vzss8++8847S5YssbW1Nd5nNzUhIQHJyTovERiGya7JP110ra5FAiDAeWSkT5h6pbxag0y6K/sn9f9WqZS1sgb1XdOxTn6Brt11/xEKsXy5ob4HgUiEiAgca3+EZicUbp4799n//e+ds2eXjBljK9RYXXDuHAICIBKZIE71uEo8dvq//6wrLwSHExD+SOQzf1WvlFdrqC7b9XprnzKVQl5bXqi+azp2+sLAuctMEjMh3aJE2F9PP/30zp07k5OTP/744/fff98YH6lSISUFFy6gpUXn/hJpxckHVx40lALwsnWLHjq9c0MJJaMqkVa0/V8ul+dt6z7ZdcxEtzGaC+p1EAqxZg00JzSS/ps0Cdevo6SkbcPTwcE7U1OTi4o+Tkx8f86c9iOlUly6hPnzjR9jye2Mk7v++iA7BYDXw4HRL7w/LGCK1jFKhbzkdvtTai5f4D0qeHLUkxMjV2ouqCdkUKF+hAZw9erVGTNmiESi7Ozs4R3rKRtefj7i4lBRoXNng0waX3L9RsUtBoyd0DbCa/KkTolNyajkSu1ya0KegNvxd4oB06KQARDxhBx1YUw+H1ZWeOopuNMKsAHw4AF279bccPXBgxn/+Y+Iz89ev374kCHtO3g8rFsHFxejhdZQXRZ/8NMbpw4wjMrO2SNi5ZuTolZpJTalQi5vbtR6odDalsvrcLbNMKoWqQSAyNaBw+nylIsQY6FEaBirV6/ev3//8uXLY2NjB+ozqqpw6hTy8nTuVDKq1HLxheLUFqWMx+GGuAXM8QkV8fQU7HByglQKhum2AGkrDgc8HsaNQ0yMqW7KscKhQ8jK0tyw+qef9mdkLB83LlbrXvTo0XjiCSNEpFTIU0/subD/Xy2NEh5fEBKzZs7qt0U2dD+AWBRKhIZRVFQ0evRoqVQaHx8/2+DTKZubcfkykpJ01ogBkFt7L+7BlZqWegD+jg9FD5vupG/mp7095s1DYCAUCiQn4+pVKBSQyzvXYGslFEKlwqhRiIigC8EBV1+PHTs0z06K6utH79ghlcni//CH2Vo3G1avxsjupjX1T27ymbhv/lZTeg+Af+iC6Bffd/J8aEA/kRCToERoMP/85z//9re/BQcHp6am8jquDOs7hkFGBk6fRqP2TSe1yubaUw+u3K57AMDVakjU0OkPO+rpac7nIywMs2ZBc+YFw6CwELdvIz8fVVVoaWnNiEIh7O3h7Y1RozBqFK0XNJ4LF3CxQw3YfyYk/O3ChWBPz9QXX+Rp1rRzdcW6dQNU5a6y8Papb/9++/oFAK6+D0c9/4+HJ8/p9lWEmClKhAbT1NQUEBBQUFDw73//+8UXXzTAOxYUIC4OWg3q2j5O0ZJQnJpSkaViVNZ80WyvyaHu47j65iP4+yMmBprPmcjgJJfjyy9R194buUkuD/jyy4La2n8/8siLkyd3ODgmBlO0p6v0U5OkNuHgpykn9qiUCmv7IbOfeCN04R+0nvMRYmEoERpSbGzs448/7ubmlpubO6Q/WaeuDufP66sRo2JUaZU554uSGxXNXA53ouvouT5TbPh6Lto8PREdjYfojpb5yMjAkSOaG2LF4scPHXKztc3dsGGI5tW5lRU2bICNjUE+VqVUpJ35/vx3HzXWV3N5/IkLnpj71J9sHJwN8uaEDGa8TZs2mToGyzFu3Lj4+PisrCy5XB4VFdWXt5DLcekSfvxRcya9prv1RT/cOZ1WeUuuUoxw8HliZOREt7ECrq4TdmtrzJ+PxYvpQtDMuLvj7l3Ni8Jx7u7xBQVZFRVypTJKs46B+vlux2qlfXM3/fIPm59NO3NQ3tI0Iij8iT//Z2LkKoFIf+MRQiwIXREa2M2bN0NCQrhcbmZm5ujRo3vxSvXjwLNn0dCgc39uVVVKVert+gIAzlaOUb7T/Ifouc7jchEaijlzaIanuSopwTffaJZKuFlaGrJrF5fDyVy3brSra/uRHA7Wru3PPKbc3NyU2E9uXzsBwNl7RNRzm/ynLOhH6ISYn75eETIMJBLU1KClBXw+DDU3xPx5enoWFhampKTcvXt31apV7TskEjQ2orlZ95+ruBiHDiE5WWex0AaZ7INLl5786ScbnpX7EGam16RlI+a5WetpDe/nhyeeQFAQ+Gb4XIfGlZq9PWpqNB8Pe9rZFdbXpxQV3a2tXTWhY5OssjIMH96HP1dDQ8MHH3zw5JNP2tg7uQsaZi5/ddkfv3Ib1l1dIUIsTi+vCGUyZGTg5k2UloLLBZcLhoFCARsbjBqFkBDqywOgvLx89OjRtbW1J375JcbXF2lprfc51RP8FArY2rb+uby8IJHg7Fn9jwOZ3Wlp7507Vy6VcjmcdaEhHy+YYyvQ80zIxQWRkfA3wx8yGledSaXYvl2zeFC5VDp6+/ba5uYTTz4Zo3U7lM+HStVhXHVJpVLt3r37vffeKy8v53K569au/XjzJtshbl2/ihBL1eNEyDBISkJ8PADtqxaGAYcDDgcCAdzd8eij0Lx1w0qfbtny1h//OMbNLePVVwU6F/9xueDxYGsLqVTfqvbkoqKNJ08mFRYCmOLjsy0mZqqvr84jYWWF8HBMnWp+l1A0rrpw+TI6diX89MqVt06fHuPqmrFunUDnv2suF3w+PDyweLG+P1dycvLGjRuTkpIATJkyZdu2bVOnTh2A6AkxGz27NSqV4v/+D9nZkMl0rOluq5OkVEIiQVoabG27PSe1ZE1NoTk5sUlJOZWVTkLhtKG6FvYxDFQqNDfrXMleWF+//sSJV0+cKKyv93Fw2PG7322PiRmqszs5h4PAQKxciZEjB2hJ2QCicdU1Hx/8+iua2ttDhnp7x2Zl5VRWOllb6x1X+v9chYWF69evf/XVVwsLC318fHbs2LF9+/ahOt+HEDbpQSJsaMA336C2NqOwUFxRoWIYZ2sdc8kYhrlQUHC3psZFJBLdvw8OB8NY2X5MKsWuXbyamhGOjgcyM68WFj47caKdUE+ps04a5fLPrl59/NCh1OJiG4HgjenTY5cvD/Xx0V2VcfhwPPEEQkLQ4/cfRGhcdYvLhYMDxOK2DTwud8SQIT0aVyoV7t4FoF4509jY+Nlnnz3++OOpqak2NjZvvPFGbGxsaGgoVfskBN0nQqUS//kP6uqgUj37v//99fx5a4Fgvp9f5wPlKtWoL77Ym56+ZMwYXzs7FBbC2Zl1RbmUSvz3v6ipgUrl7+KSXFQkLi9vkMkW9ey53dGcnEcOHDiclSVTKhf5+x978sllAQFCnXfAHBwQE4OoKK0+dmaDxlUPubmhsBA1NW0bejGuVCoUFsLJ6ei1a4888sjhw4dlMtmiRYuOHTu2bNkyoTmePBEyMLq7mXbunPrXqtdvLJfjl1/0rQSwWAkJqK5uu8v3WVSUgMfbdf369eLirl93o6Rk5u7diw8eLKitneTldfGZZ46uWjVc5/o/dWv4DRsQpLvFrnmgcdVzUVGtd4l/e5zfi3F1//7MRx5ZvHhxQUHBpEmTLl68ePTo0QFvkEKIuekyEdbVITW1pw0KOlMqcfZsH19rjiQSJCVBLm97uDXG1fWV0FAVw7wWF6dvUlJlY+PGkyenfPPN5fv3XWxstkZHJ7/wwkx9hWACA7FhAyIizHJpRBsaV71iZdU6ovo2rgoKXOztt27dmpycPHPmTKNFTYgZ6TIRXr2qr91BjyiVEIv1VYu2QFevdr7E+XtEhJut7eX793/MztbaJVcqtyUljdy27Ytr17gczqthYXdefXXj1Kk8nXNevL3x7LNYutQSOuLSuOqVq1fR6Ule78bVxo0bn3/eYIXgCbE4+hOhutBJH25eaeJw0Ok/VIuVnt75932IlZW6vfgbp041alwDnc3PD9q587W4uPqWlvl+fulr126LiXHU2eTB3h6LFuH552EZs/toXPVW/8eVtbVWm0NCiCb9ibCmBgpFf99eLkdubn/fxCzU1Oi71/fCpElBnp4P6uo+v3oVQE5l5cL9+xfs3ZtdUTHa1fX4k0+eWbNmrJuutcw8HsLC8MormDy58zWBuaJx1SsGGVdyOXJyjBk1IeZF/6Om8nKd69IK6+sv3bvXebtC3zm+ni5ClqaiAjyezt8sHpf7eVTU3D17Nl+6dKem5rv0dIVK5WRt/c6MGa9Pm6Z7UigAf39ER8NJTx0180XjqlcMNa4qKowUMCFmSH8ibGzUef9qf0bGfj31wHSTSltLhFi2ykqdZULVZj70ULCnZ3pZ2X/T0vhc7itTpvwjIsJFX/ccV1dERUGzyYAloXHVK4YaV83NAxgkIWau15MPx7u7T9JV3UPFMPt0/pApFPjqK0RHY+TIPsRnBuRyJCbi0iV9z73O5ue/FhcnLi/ncjgMsGfpUu2iyZp4PLz8suX/vndC40qbYccVIUQ//YnQxkbnLayF/v4fzZ/febtMqdT9gwWgshL79llgk3SGQXY2Tp/WbB2n6XZ19Z/PnTskFgMY6ewc5OHxU3b2Z1evPjF+PFdfqrO1tfAsSOOqWwMxrnTOwyKEAOgqEbq793dqn5bcXNy5g5AQzJ1rliXBtBQXIy4ODx7o3CmVyf515crHly83KxS2QuFb06f/KTxcxTApO3ZcLy7el5GxRt9yeJ2zZiwJjauu0bgixOj0J0InJ/D5fVv1XNrQIJXJRjo7a+9QKnHtGrKyMHs2Jk0y10sfiQQJCbhxA7rWMjMM811GxjtnzpQ2NHCAp4KCPlmwwPO3Qmj/b+7cp48cefvMmSVjxjh07porEKBXvXzNEY0rfWhcEWIi+pdPqNsa9GkR7ttnzoz98suNJ0/Wa3RTayeR4NgxfPutvtPewUv9g/vll7h+XeevVUpR0Yzdu58+cqS0oSHUxyfxuef2Ll3qqVEO9KnAwBnDhpU1NHx8+bKO92cYBAQMXPiDAo2rzmhcEWJSXVaWmTatD519VACfy1UolV9cuzZmx47/u3lTpbMKVHEx/vtf/Pyz2dSNzM7Gjh2Ii4OuX+Gi+vo1R46EffPN1QcPfBwc9ixdeu355zs3yuFwONuio7kczpYrV/Kqqjrs4/EwZgxsbQfuGwwWNK40DfS4AuDnx4pxRUhfdfl75OiIkBAIBL19x92PPpry4ovhw4aVSCTP/PyzupCmjkMZBunp+OILxMcbYJH1wFFPyoiNRW1t551NcvnHly+P2bHju/R0K4HgnfDw7FdeWRMUpK/BzWRv79WBgTKl8k9aBTO5XERGDkT4gw6NKzXjjCtgUP8RCBkEujsxnzcPjo59OH+f7O198ZlnYpcvf2jIkOvFxbN2715x6NA9Xf/BQy5HQgK2b0d6em8/ZcA1NSEuDl99hTt3dO4/mpMT8OWXfzp7Vt0TR/zyyx/Nn2/f+SFNR58sWOAgEv2UnX2m7W0FAixebAl1RHuIxpVxxpVafj5u3zZU7IRYHo6+6vXt1A1UpdKMoqLKpqaHHB11zFb4rYEqgBBvb83H9Y1y+SeJiZ8kJjbJ5TYCwYawsL/MmqW3oejw4YiOhodHX7+O4ahUSEvD+fP6ijunlZS8Fhd38d49ABO9vLZGR8/S1zJCl80XL/7l/Plx7u43167li0SYORNs6wxA40oXQ44rzfMMV1esW9eHMw9C2KAHiRCAVIr9+1FV1UWRi64V1tf/+dy5fenpDODr4LB53rynAgN13+ThcDBhAiIjTflUIz8fp06hvFznzqrGxvcTEr5MSVGqVC42Nn+dNWv9lCm6W0bo16JQjP/qq9vV1V8tXrzu73/HpEmGiNvc0LjSYOBxtXDhutDQDvuiojB1ap9jJ8SC9SwRAmAYJCUhPh6A9s+WutIVhwOBANbWaGzUNzn+WmHhxri4a4WFAMJ8fbdFR4f5+ur+OCsrhIdj6tS+TS/su+pqnDunr1S/XKn8KiXl7/Hxdc3NAh5vXUjI+3Pm6G4ZAcDREY2N4HB0/8pzuT/eurXswAFnJ6fcvDwXFxfDfQezQuPKsOMK+DEra1lsrLO1de6GDR3KrVlZYcMG6CvsRwiL8TZt2tSjAzkcDB2KsDDY2aGxEVIp+HwIBK0/KLa2GDcOv/sdFixAcDCkUp01kX0dHJ6dOHGks/PVwsKcysrdaWl3amqmDx2q446WQoH8fIjFcHKCcZKETIZLl/Djj/pO2JGTU3X6dGRycrNCscjf/38rV64ODLTS2SDXzg6RkVi6FNOmwcEBjY1obASXCz4fPB4YBnZ2GDcu4A9/uJKRIRaLW1paoqOjB/TLDV40rvo/rjqeywa4uV158EBcXt6iVEZrVqxVKCCTwd/fgF+OEMvQ4ytCLQwDiQSNjRAIYBAaCjYAABEFSURBVG+vo6JHURHi4lBYqPPVDTLZlitXPrp8uUWhsBMK35w+/d3wcJG+rut+foiOHsDSGOoOeWfOQCrVfUBFBeLi1PMa9gUEuE2aFKWvIjaXi9BQzJmDzvMaJJLWyXt2dm0TJsVicXBwMIC0tLTx48cb4suYORpXOnU9rkpK8P33mulQXF4evHMngLS1a8e7u7cfzOHgpZcGxbNSQgaTHl8RauFwIBLBzg42NrrvMjk4YOJEODmhsLDzPRwhjxcxfPjygID8mhpxRUVCQUFFWdkikUj3SXpNDa5fR1MThg6Fvh+1PisqQmwsUlJ033ZrakJCAn7+GVVVsLJCRERgZOTDrq6638rfH6tWYcIE3UGKRLC2hrW15p/L3d29pKQkOTk5Ly9vzZo1BvlC5o3GVWfdjisXF9TVobS0bZu7rW2JRJJcVJRXVaVddK2yEsHBff9ehFiivl4R9py6iH5ior7FTGfz898+ffpgS8vomprWk3TNc1hN1taYPRtTphimhlZ9Pc6dg76CzurT+dOnIZW2FkPpYp5FPxonVVdX+/v7V1VVHT16dNGiRX14B5aicaVJKsX27ZpL8qubmvy3b69qbDy6atUirduhK1Zg7Niefh1CWGDgE6Fal5MFGKWSk5yMhAQ0N4PHQ2goIiL01sv39ER0NHozp1ybXI7kZFy8qHey4t27iItrfRw1YgSiouDpqftIKytERCA0tD8T07du3fr666+PHDlSLBaLulsrRjqgcdUmMREdl9JvTUp6PS5upLOz+OWXO9wfHjIEr7xi+ItgQsyWsRKhmvqHQN+sgaYmxMcjJQUqVftJur4fAn3Nd5qaUF8PlQpWVjqWbHfX4AZ1dTh/vnUJtqMj5s6FvnL+XC4mTsTcuf2fhqdQKIKDg8Vi8ZYtW958881+vhsb0bgCoFTi66+hUV9NoVIF79wpLi/fEhn55vTpHQ6eO5d161YJ0c+4iRDdLyhGRQVOnWothOHm1tWtIT4fYWGYNQtCIYqLkZKC3Fy0tLSe6jIMFAp4eCA4GMHBEApRUoK4OOgsygVAJsOVK7h8GQoFBALMmIHwcL1nzSNGdHWrrffOnj27YMECe3v73NxcT31XCaQLNK4A5OTg++81N5zNz1+wd6+9SJS7YYNmkW4IhVi/nkWVjAjpktEToVpzc/tJuk45OYiLQ00N8NtJupOT7iPt7CASQSKBXK6zcj8EAnA4cHdHUZHuAxgGWVntp/PjxiEyEo6Ouj/O2Rnz5g1ELf9FixYdP378hRde2LVrl8HfnC1oXO3bp1W2bdGBA8dzc1+YPHnXI490ODIoCEuW9PfjCLEIJkqEapWV7SfpnSmVSEnBhQtoaQGP19p5VesRmnrNdX8UF+PkydbGPd7eiIlBp9L+rQQCTJ/e1el8/9y+fXv8+PFyufzatWshISED8RFsweZxVVmJr7/WPA+4XV09/quv5ErltRdeCPH27nDw88/Dx8cAH0qImevr8gmDsLFBYCC8vVFUhOZm7b1cLnx9ERwMmQzFxSgsRHo6RCJ4ebX/SHX+tZLLUV8PiQQqFYTCrn7OJBKcPo3jx1FXB3t7REVh0SK9J+yBgVi1Cv7+A1et0dnZuba29sqVK2Kx+Nlnn9XXZIB0j83jysYGjY0oKmrb4GxtXdvcrF5i/+zEiR3GVUmJGfcxJsRwTHpF2Ealaj9J10nr/Do6GsOGdTiAYZCZiZQUFBW1nw7b2cHfHzNnat/+6sk1QRsfH0RHQ1/JLoOSSCT+/v6lpaU//PDDihUrjPCJFo6d46q5Gdu3az4rlbS0+G/fXtrQ8MPy5SvGjetw8NKlCAw0fAyEmJXBkQjVGhoQH48bN3r9xKWlBbGxrY9GbGzg4QGhEPX1KC0Fw4DPx2OPtT996flTInt7zJuHwEBjnjLv2rXrpZdeGjp06K1bt2yoLKRBsHBcJSfj5EnNDbuuX3/p6NGhjo631q+30ewEaWeH9ev1ZmtC2MGkt0a1CIXw98eoUSgvR3299l71xITJkwGguBhlZbCzw7BhYBgcOoS8PAiFeOQRLF2KiRMxYQJCQhAUhMpKVFUhOxsjRsDRETIZ9u1DXR08PbFsGWbOhLW1jkgEAsyciWXL4O1t5BtHEydOPHHixK1bt0Qi0ezZs4350RaLhePK2xvZ2Zql3SZ6eZ3Iy7tVWSni82cPH95+pEwGDgcjRgxgMIQMeoPpilBTbi5OnNC7KqumBomJiI4Gn4+sLMTGAsCaNfDz0z5SqcTevbh3D25uePllcDgQi9HUhEmTer2SzFgSExNnzpxpZWWVnZ39UH/Wd5PO2DOu7t7F3r2aGxLv35+5e7eVQJD9yisPaYbB4+Hll6GrGSQhLDGYrgg1ubhg8mSIRCgshFKpvdfaun1+wfHjqK1FQADCw3W8D5cLb2+kpqKxEQ89BCcnuLvrPR/39MTy5QgP11t8xCiGDRuWlZWVnp5eXl7++9//3oSRWCD2jCsnJ5SWaq6vH+bomFVRkV5aWi6V/l5znYa60LnWs0NC2GQQd6xWrz5+5ZWuHubL5a0zHbpo3eDh0VpuPz+/dUvni2Bra0RH48UXtedKmMiWLVtsbGwOHjx46dIlU8dicdgzrqKitAqXb4mMtBEIDmZmXrp3r8OR2dnt34IQ9hnEiVDNwQFLl+pd8FRR0Xpe7+XV1Zuol0+1lefXPG3nchEWho0bERY2eOaRDx069M0332QY5rXXXlPpWxtO+oMN48rJSasl/VBHxzenT2eA1+LiVFppOy5ObxUCQizdoE+Eaj4+eO45LFmiXae/bY64vvr9mns7F9/y88PatYiOHoSz5t59991hw4bduHFjz549po7Fcln8uJo1C5qV1YB3w8OHOTreKCnZc/NmhyMrKnDjhlFjI2TQMJNECIDDQVAQXn21QwGOtnPYrtcjq28QaT4T4nAQGoqnnhrAvqz9Y21t/cEHHwB499136zvPdSSGYtnjSijEvHmaG6wFgg/mzQPw7rlz9VrLK8+fR1OTMaMjZJAwn0SoJhR2mJjedsatb8W05l7NqQoCgXHWyPfHqlWrwsPDy8rKPvzwQ1PHYukseFwFBWnd/l01YUL4sGFlDQ0faj2BbmrCxYtGjY2QwcHcEiHQYZ1y2yxwjdlxOlRWar+Qy9W73nnQ4HA427Zt43K5n332WV5enqnDsXSWOq44HERHd9zA2RYTw+VwPrt6NU/rCyYno6LCqOERMgiYYSL082ufC+fo2NpKpqBA7/FtMwA1T9Xlcr09UQeTSZMmrVmzRiaTvf3226aOxdJZ8Ljy9cWECZobJnl5rQkKkimVb5850+FIlQpxcUaNjZBBwAwT4ZgxHZ7cqCfBX78OuVz38WlpaGmBQNChx42PDzQLTQ1iH330kYODw88//3zq1ClTx2LRLHtcLVigFdhH8+c7iEQ/37p1SqtNR34+6PYDYRkzTIQeHh3Kc4SFwcoKdXX45Rcd878LC3HuHABMmdLa8pthIBRCq2H3IObh4fHuu+8CeOONNxQKhanDsVyWPa7s7TFjhuYGDzu7d2fOBPDGqVMKrS8YF6ej2gAhlmuwVpbpmpMTbt1q/XkSieDsjOxslJUhJwd8Png8yGQoL0dSEk6cgFwOX18sXdp6vs/hYMgQxMQMnlWD3ZoyZcoPP/yQk5Pj4eExZcoUU4djuSx7XPn6IjNTsy/VFB+fH8TinMpKD1vbKZoTapqaYGWlt4EiIRZnsNYa7daBA7h7F21XSLdv49gx1NZqH6aeHP+730EobN3C5+OZZ6DVoXTQO3LkyGOPPebk5JSXl+fi4mLqcCyXZY8rsRiHD2tuOJKd/dgPPzhZW+dt2OCi2e1EJML69VprEAmxVOZ5RQhg1ChkZEAub61r5eyMkBB4esLaGnZ2cHSEjw8CAxETg8mT2ydBCASYNw9jx5ow8L4ZO3bs1atXxWJxU1NTTEyMqcOxXJY9rtzdUVCgWXN8rJvb1QcPxOXlTQpFzKhR7UcqlZDJ4O9vgiAJMTqzvSIEIJHgP/+BVIoePjkTCBAejlmzBjisgZKVlRUUFMQwTFpa2oSOkwCJIVn2uCotxa5dmmVRsyoqgr7+mgHSXnppgrp6apuhQ+HhAW9vDB8+6JaFEGI4ZntFCEAkQnAwSkshkXRTJpHPB5+PJUsQEmKs4AzPzc2trKwsOTk5Nzd3zZo1pg7Hcln2uLKzQ309SkraNrjZ2pY1NCQXFeVWV68JCupwcH09iotx9y6SkyEWw8EBrq7GDpiQgWfOV4Rtbt/GmTOoqYFKpV3vSiiESoWgIMyZA/Nv+F5VVeXv7z979ux9+/ZR//oBZ6njSirF9u2aRXOqGhv9t2+fPXz4vsces+li+YdQCHd3LF8OBwdjxEmIsVhEIlSrqEBeHu7fR20tFArY2MDDA35+ePjhQbq0q09KSkq8um6JQAzLIsfVlSvouJS+RCLxUtcQ6BqXC6EQa9Z005eDELNiQYmQHXbt2nXgwIFp06bpK0D61FNPPXjw4I9//OPChQuNHBsxG0olvv5as4DcruvXD2RmTvP1/XD+fJ2veOqnnx7U1/9x+vSF/v4QifD883SblFgMM1xQz2537txJSEjIyMjQd8C1a9cSEhJKNB4CEaKNx0NkpOaGO9XVCQUFGWVl+l5xragooaCgpKEBAFpasH9/TycTETLoUSIkhJX8/bvpttg1qZRaVRCLQYmQEFaqrOymyVTX5HIkJWnWqSHEfFEiJISVUlK6WRzSLQ4HmZkGioYQU6JESAgrZWX1NxHKZPj1VwNFQ4gpUSIkhH2am9HUhP7PGC8tNUQ0hJgY39QBkL5ITk5esGCBzl1FRUVGDoaYn9pa8Pmdey0lFxUt2LtX5yuK6ut1bFUoIJO1Vx4nxDxRIjRLlZWVZ8+eNXUUxGzp6TZc2dh4Nj+/F+/D40Eup0RIzB0lQrMUERHx7bff6tw1f/78goIC44ZDzI2emjgRw4d/u3ixzl3z9+4t6NyOSqk04/I6hPyGEqFZsrGxGTlypM5dAvphIt1ydIRCAYbRaiNsIxCMdHbW+QpBW88pTTweXQ4SC0CTZQhhH2trWFtrZcG+8PQ0RDSEmBglQkJYacwYcPv3n79AgPHjDRQNIaZEiZAQVpoypb+JEEBgoCFCIcTEKBESwkpubhg5su+5UCBAWBisrAwaEyGmQYmQELZauLCPcz45HNjYYPZsQwdEiGnQrFEzM2LEiBkzZowbN07fASEhIe7u7p40i4F0y94eK1bg++8hl49wcpoxbNg4d3d9x4Z4e7vb2nra2QGAQIAnnwSffj2IhaDGvISwW04OfvxR3xJ7bRxOa4d6b+8BDosQ46FESAjrlZXh++8hlXaTDoVCuLpixQo4OhorMkKMgRIhIQRQqZCWhosX0dwMlapD93n1c0RHR8ydizFjDLD6kJBBhhIhIURDaSkKClBSgoYGcDhwcIC3N/z8oKfiDCEWgBIhIYQQVqPlE4QQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYTVKhIQQQliNEiEhhBBWo0RICCGE1SgREkIIYbX/DyFOkY2prPTwAAAAlHpUWHRyZGtpdFBLTCByZGtpdCAyMDIxLjA5LjExAAB4nHu/b+09BiDgZUAANiBmBeIGRmUGBSDNyAamWDjAFBMjgwZIkBGd5mZgZGBkYmBkZmBkYWBiZRBhEA8CScCNvb75rj0Dg4M9iANk2wEpFQg7dj9QfD+I/S1nuv1Dt2VgtqbPon1ANlh9+gSuA0BKFcQWAwCe8xnC5LL9iAAAAOV6VFh0TU9MIHJka2l0IDIwMjEuMDkuMTEAAHicfVFJDsIwDLznFfMBIjtu0uZIFwFCtBIU/sAV8X/hgEqKgMY5jJ3xeIlBOsd2f73jfVxrDEALN8aIixCROSAB1N1m16MZ1/UUaYZzP54Q4DVD7ZO5HofDFGHUN7AtxZFnkC29CmuGJXqBiejQaLQQIZbZ8xdPMGBFtgqBnE/or2KBrVbm6qm4YutiJKl+EL0SVUg4xKrEAjEkorMcPIViqcmubz+28NpLPfRt3ksyl6dXB5KHTG6RJ0muz/2yMsK82Fw6+dNnKTYPtiBiN3mcahEAAAB4elRYdFNNSUxFUyByZGtpdCAyMDIxLjA5LjExAAB4nD2LsQ3AIAwEV0kJkrFsDAZElzR0GQCxRGqGD6FI8zr9/fc27sv0NuyO8zmmcR5ZIykQVEeYVclHWJQiVMIgQiyfZOS8GByjL4XCngtryQn+ijGJpwj7bucLnfgYD0jjddUAAACVelRYdHJka2l0UEtMMSByZGtpdCAyMDIxLjA5LjExAAB4nHu/b+09BiDgZUAANiBmBeIGRmUGBSDNyAamWDjAFBMjgwZIkBGd5mZgZGBkYmBkZmBkYWBiZRBhEA8CScCNvb75rj0Dg4M9iANk2wEpFQg7dj9QfD+I/S1nuv1Dt2VgtqbPon1ANlh9+gSuA0BKFcQWAwCe8xnCh5RHQwAAAOZ6VFh0TU9MMSByZGtpdCAyMDIxLjA5LjExAAB4nH1RSQ7CMAy85xXzASI7btLmSBcBQrQSFP7AFfF/4YBKioDGOYyd8XiJQTrHdn+9431cawxACzfGiIsQkTkgAdTdZtejGdf1FGmGcz+eEOA1Q+2TuR6HwxRh1DewLcWRZ5AtvQprhiV6gYno0Gi0ECGW2fMXTzBgRbYKgZxP6K9iga1W5uqpuGLrYiSpfhC9ElVIOMSqxAIxJKKzHDyFYqnJrm8/tvDaSz30bd5LMpenVweSh0xukSdJrs/9sjLCvNhcOvnTZyk2D7YgYjf33Q8TAAAAeXpUWHRTTUlMRVMxIHJka2l0IDIwMjEuMDkuMTEAAHicPYuxDcAgDARXSQmSsWwMBkSXNHQZALFEaoYPoUjzOv399zbuy/Q27I7zOaZxHlkjKRBUR5hVyUdYlCJUwiBCLJ9k5LwYHKMvhcKeC2vJCf6KMYmnCPtu5wud+BgPlKiVqAAAAABJRU5ErkJggg==\n",
"text/plain": [
"<IPython.core.display.Image object>"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from rdkit.Chem import Draw, rdDepictor\n",
"\n",
"parameters = schema.parameters\n",
"rd_mol = target_molecule.to_rdkit()\n",
"rdDepictor.Compute2DCoords(rd_mol)\n",
"rd_mols = []\n",
"atoms = []\n",
"for parameter in parameters:\n",
" matches = target_molecule.chemical_environment_matches(query=parameter.smirks)\n",
" flat_matches = [atom for match in matches for atom in match]\n",
" rd_mols.append(rd_mol)\n",
" atoms.append(flat_matches)\n",
"\n",
"Draw.MolsToGridImage(rd_mols, highlightAtomLists=atoms)"
]
},
{
"cell_type": "markdown",
"id": "03c7cf3b",
"metadata": {},
"source": [
"We can see that BespokeFit has generated two bespoke SMARTS patterns, as determined by the\n",
"symmetry of the atoms around the central bond, covering all three torsion terms.\n",
"\n",
"## Setting up the BespokeExecutor\n",
"\n",
"In the last blog post we made use of the invaluable QCFractal infrastructure to compute our target data. Since then however,\n",
"we have been hard at work to offer a local execution pathway as an alternative reference generation method.\n",
"Here we still use QCEngine to execute the tasks, but spin up simple [Celery](https://docs.celeryproject.org/en/latest/index.html)\n",
"workers to perform the jobs. In fact the BespokeFit executor has been totally reworked to allow for a more flexible and scalable\n",
"fitting solution. Each of the major stages: fragmentation, reference generation and optimization are carried out\n",
"by dedicated workers which can be on the same local machine or distributed over multiple. Here we will be spinning up\n",
"the `BespokeExecutor` and three workers, one for each of the stages, in a single line! We will then submit our optimization\n",
"task via the RESTful API and let BespokeFit take care of the rest.\n",
"\n"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "c4f42b06",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n",
"Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n",
"Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n",
"Warning: importing 'simtk.openmm' is deprecated. Import 'openmm' instead.\n",
"Warning: Unable to load toolkit 'AmberTools'. \n",
"Warning: Unable to load toolkit 'AmberTools'. \n",
"Warning: Unable to load toolkit 'AmberTools'. \n",
"Warning: Unable to load toolkit 'AmberTools'. \n"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Output()"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
],
"text/plain": []
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"><span style=\"font-weight: bold\">[</span><span style=\"color: #008000; text-decoration-color: #008000\">✓</span><span style=\"font-weight: bold\">]</span> fragmentation successful\n",
"</pre>\n"
],
"text/plain": [
"\u001b[1m[\u001b[0m\u001b[32m✓\u001b[0m\u001b[1m]\u001b[0m fragmentation successful\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Output()"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
],
"text/plain": []
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"><span style=\"font-weight: bold\">[</span><span style=\"color: #008000; text-decoration-color: #008000\">✓</span><span style=\"font-weight: bold\">]</span> qc-generation successful\n",
"</pre>\n"
],
"text/plain": [
"\u001b[1m[\u001b[0m\u001b[32m✓\u001b[0m\u001b[1m]\u001b[0m qc-generation successful\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Output()"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"></pre>\n"
],
"text/plain": []
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<pre style=\"white-space:pre;overflow-x:auto;line-height:normal;font-family:Menlo,'DejaVu Sans Mono',consolas,'Courier New',monospace\"><span style=\"font-weight: bold\">[</span><span style=\"color: #008000; text-decoration-color: #008000\">✓</span><span style=\"font-weight: bold\">]</span> optimization successful\n",
"</pre>\n"
],
"text/plain": [
"\u001b[1m[\u001b[0m\u001b[32m✓\u001b[0m\u001b[1m]\u001b[0m optimization successful\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\n",
"worker: Warm shutdown (MainProcess)\n",
"\n",
"worker: Warm shutdown (MainProcess)\n",
"\n",
"worker: Warm shutdown (MainProcess)\n"
]
}
],
"source": [
"from openff.bespokefit.executor import BespokeExecutor, wait_until_complete\n",
"import os\n",
"\n",
"# set keep files to true so we can view the results\n",
"os.environ[\"BEFLOW_KEEP_FILES\"] = \"True\"\n",
"\n",
"# launch the executor\n",
"with BespokeExecutor(\n",
" n_fragmenter_workers=1,\n",
" n_qc_compute_workers=1,\n",
" n_optimizer_workers=1) as executor:\n",
" # grab the task id and wait for the task to finish\n",
" task = executor.submit(input_schema=schema)\n",
" result = wait_until_complete(optimization_id=task.id)"
]
},
{
"cell_type": "markdown",
"id": "84b22a7f",
"metadata": {},
"source": [
"As our optimization task works through each stage you should see them complete with a tick and once all tasks are complete\n",
"the executor will shut down for us and clean up the workers as well! We can then check the result object to make sure\n",
"each stage did finish with no errors."
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "435bea97",
"metadata": {
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"fragmentation success\n",
"qc-generation success\n",
"optimization success\n"
]
}
],
"source": [
"for stage in result.stages:\n",
" print(stage.type, stage.status)"
]
},
{
"cell_type": "markdown",
"id": "b45f0b79",
"metadata": {},
"source": [
"## Inspecting the Optimized Parameters\n",
"\n",
"Once the torsion optimization is complete a result schema is returned. The schema contains the final optimized\n",
"parameters that can be used with the `openff-toolkit` in workflows to parameterize molecules and set up systems in\n",
"OpenMM to run dynamics. The result schema also contains all provenance information which can help with reproducibility.\n",
"\n",
"Now to convince ourselves that the optimization was successful and has lead to an improvement in the parameters - lets load\n",
"the final output from ForceBalance."
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "68295664",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\n",
" <iframe\n",
" width=\"900\"\n",
" height=\"600\"\n",
" src=\"bespoke-executor/85c79a09-2897-43db-9abe-acba00146076/optimize.tmp/torsion-0/iter_0002/plot_torsion.pdf\"\n",
" frameborder=\"0\"\n",
" allowfullscreen\n",
" \n",
" ></iframe>\n",
" "
],
"text/plain": [
"<IPython.lib.display.IFrame at 0x1c99c17f0>"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from IPython.display import IFrame\n",
"\n",
"opt_id = f\"bespoke-executor/{result.results.input_schema.id}/optimize.tmp/torsion-0/iter_0002/plot_torsion.pdf\"\n",
"\n",
"IFrame(\n",
" opt_id, \n",
" width=900, \n",
" height=600)"
]
},
{
"cell_type": "markdown",
"id": "f9ab4fa2",
"metadata": {},
"source": [
"BespokeFit has been successful as there is a clear improvement in the PES around this torsion with respect to the reference data.\n",
"We can also plot how the parameters have changed during optimization as the result schema stores the initial and\n",
"final torsion parameters."
]
},
{
"cell_type": "code",
"execution_count": 19,
"id": "8bb9d4b5",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#1H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#35H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=2.977847229905e-05, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=-0.05907705944855, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.2124075556718, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=5.166921989902e-06, unit=kilocalorie/mole)},\n",
" ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#35H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#1H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=-2.134326655309, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=1.099572748702e-05, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.6985897230115, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=-4.932432892583e-06, unit=kilocalorie/mole)}}"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"result.results.refit_parameter_values"
]
},
{
"cell_type": "code",
"execution_count": 20,
"id": "f688cb80",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"{ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#1H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#35H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.6985935181583, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=0, unit=kilocalorie/mole)},\n",
" ProperTorsionSMIRKS(type='ProperTorsions', smirks='[#35H0X1x0!r+0A:1]-;!@[#6H2X4x0!r+0A:2](-;!@[#1H0X1x0!r+0A])(-;!@[#1H0X1x0!r+0A])-;!@[#8H1X2x0!r+0A:3]-;!@[#1H0X1x0!r+0A:4]', attributes={'k2', 'k4', 'k3', 'k1'}): {'k2': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k4': Quantity(value=0, unit=kilocalorie/mole),\n",
" 'k3': Quantity(value=0.6985935181583, unit=kilocalorie/mole),\n",
" 'k1': Quantity(value=0, unit=kilocalorie/mole)}}"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"result.results.input_schema.initial_parameter_values"
]
},
{
"cell_type": "code",
"execution_count": 21,
"id": "bf25a76c",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>parameter</th>\n",
" <th>before</th>\n",
" <th>after</th>\n",
" <th>change</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>0</th>\n",
" <td>smirks_0_k1</td>\n",
" <td>0.000000</td>\n",
" <td>0.000005</td>\n",
" <td>0.000005</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <td>smirks_0_k2</td>\n",
" <td>0.000000</td>\n",
" <td>0.000030</td>\n",
" <td>0.000030</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <td>smirks_0_k3</td>\n",
" <td>0.698594</td>\n",
" <td>0.212408</td>\n",
" <td>-0.486186</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <td>smirks_0_k4</td>\n",
" <td>0.000000</td>\n",
" <td>-0.059077</td>\n",
" <td>-0.059077</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <td>smirks_1_k1</td>\n",
" <td>0.000000</td>\n",
" <td>-0.000005</td>\n",
" <td>-0.000005</td>\n",
" </tr>\n",
" <tr>\n",
" <th>5</th>\n",
" <td>smirks_1_k2</td>\n",
" <td>0.000000</td>\n",
" <td>-2.134327</td>\n",
" <td>-2.134327</td>\n",
" </tr>\n",
" <tr>\n",
" <th>6</th>\n",
" <td>smirks_1_k3</td>\n",
" <td>0.698594</td>\n",
" <td>0.698590</td>\n",
" <td>-0.000004</td>\n",
" </tr>\n",
" <tr>\n",
" <th>7</th>\n",
" <td>smirks_1_k4</td>\n",
" <td>0.000000</td>\n",
" <td>0.000011</td>\n",
" <td>0.000011</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" parameter before after change\n",
"0 smirks_0_k1 0.000000 0.000005 0.000005\n",
"1 smirks_0_k2 0.000000 0.000030 0.000030\n",
"2 smirks_0_k3 0.698594 0.212408 -0.486186\n",
"3 smirks_0_k4 0.000000 -0.059077 -0.059077\n",
"4 smirks_1_k1 0.000000 -0.000005 -0.000005\n",
"5 smirks_1_k2 0.000000 -2.134327 -2.134327\n",
"6 smirks_1_k3 0.698594 0.698590 -0.000004\n",
"7 smirks_1_k4 0.000000 0.000011 0.000011"
]
},
"execution_count": 21,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"import pandas as pd\n",
"from simtk import unit\n",
"\n",
"parameter_data = []\n",
"for i, (parameter, initial_values) in enumerate(result.results.input_schema.initial_parameter_values.items()):\n",
" for final_parameter, final_values in result.results.refit_parameter_values.items():\n",
" if parameter.smirks == final_parameter.smirks:\n",
" for term in range(1, 5):\n",
" k_before = initial_values[f\"k{term}\"].value_in_unit(unit.kilocalorie_per_mole)\n",
" k_after = final_values[f\"k{term}\"].value_in_unit(unit.kilocalorie_per_mole)\n",
" parameter_data.append([f\"smirks_{i}_k{term}\", k_before, k_after, k_after - k_before])\n",
"\n",
"# make a pandas dataframe\n",
"df = pd.DataFrame(parameter_data, columns=[\"parameter\", \"before\", \"after\", \"change\"])\n",
"df"
]
},
{
"cell_type": "code",
"execution_count": 22,
"id": "ccc48b6c",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAaIAAAHqCAYAAABGAEpuAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAABNQUlEQVR4nO3de5hcVZnv8e9PEiRpSQcIaogOQTE6iDgOgjAqJA6CioIM7YhCBNEBBQQOXiAkowwkBGEIIjN44KAERBwghIuIyjVc1JHLKDiOoAKRkTQkgaSB7oTc3vPH2gWVSlf17kpV7a707/M89XTX3qvWequ6ut5aa6+9tiICMzOzoryq6ADMzGx4cyIyM7NCORGZmVmhnIjMzKxQTkRmZlYoJyIzMyvUiKIDsCHH8/nN6tDX1wfA6NGjC45kyFLVHT6PyCr4DWFmzVA1EXlozszMCuVEZGbWAH19fS8Pz9ng+BiRmVkD9PT0AD5GVA/3iMzMrFBORGZmVignIjMzK5QTkZmZFcqJyMzMCuVEZGZmhfLKClbJbwgzawavrGBmZkOTE5GZmRXKicjMrAGWLFnCkiVLig6jLXmJHzOzjdTT08Pll1/O4sWLmTRpEl1dXXR2dhYdVtvwZAWr5DeE2SDMmjWL2bNn09vb+/K2jo4Opk2bxvTp0wuMbMjx9YgsN78hzHKaNWsWM2bMqLp/5syZTkavcCKy3PyGMMuhp6eHCRMmrNcTqtTR0cGiRYsYM2ZMCyMbsobH9G1JcyUtzFn2CEkhaccmh4Wkt0u6RdKLkp6VdKmkrQdZR654JR0u6VpJf87Kz92o4M2sX/PmzauZhAB6e3uZN29eiyJqX5tUIgLOAA4qOohykrYDFgCjgC7gWGAf4CZJzXj9DwPeDNwKPN+E+s0M6O7ubmi54WyTmjUXEY8NVEbSSGBNC8Ip+SowEvhYRCzPYlgE3AV8HJjf4Pb2i4h1WTsfanDdZpYZP358Q8sNZy3vEUmaJOk6SYslrZT0pKRrJI2QNDkbTvq4pIskPSdpmaTzJG0maTdJ90rqlfQ7SftV1L3e0JykiVl9x0g6O0sALwFjq8S2q6RnJM2XtEW27QRJv5e0IovlAUmD6XUdAPy4lIQAIuJu4EngwEHUkyveUhIys+bq6uqio6OjZpmOjg66urpaFFH7KmJo7iZgAvBFYD/gFFJyKI/lW0Av8Eng34ATs22XA98D/gF4DpgvaVyONqcDk4CjSEN3KysLSNqXNIR2HfCJiFgp6VDgXOCHwEeAQ4F5QK7jO5JGATsA/93P7t8BO+Wpp0rdG8Rbb11mNnidnZ1MmzatZplp06Z5okIOLR2ay5LGW4ADI+LGsl1XZvtL9++IiJOy32+VtD9wHPD+iLg3K9sNPATsD1w2QNPPAAdF2RTBsrbIEs6lwFkR8fWyx+0JPBwRp5dtuznHUy3ZijRTZFk/+54D3jqIul5WI966SDqKlKSZM2cOU6dO3dgqzYaFo48+mt7eXs4//3z6+vpe3j569GhOOOEEjj76aJYuXVpghEPHuHHV+wytPkb0LPA4cJak1wELIuKP/ZT7ScX9R4BJpSRUtg3gjTnavT6qz1M/ETgaOD4ivlOx737gGEkXADcAv4iIPvIrZbv+2q46lXEAJ1I93rpExMXAxaW7jajTbLg488wzOeWUU7jkkkvWW1nBPaH8WpqIIiIkfRA4DZgNbCPpCeCcig/Vyh7EKmB5RV2rsl7NFjmarjVt5RDgKeDafvZdntX/OeAYYLWkm4GTImJhjnaXkT7Y+xvK24rUKxqsWvGaWQHGjBnDpz71KcCTE+rR8mNEEfF4RHwG2BZ4F3AHcKGkDzez2Rr7DiYdo1og6fXrPSi5KCJ2B8YBhwO7A1flajT1nhYCb+9n907A/+SpJ2+8ZmbtqLDziLIP+d8ApWNBOxcUylPAZNJrcaekfr/ORMSyiLgKuJrBxXojsL+kl1dAlPQ+YPtsX1PiNTNrF62erLALcD6pR/EnYDPgCNJ5PXcAW7YynpKI6JY0Gbid1NOYEhGLJF0MvAD8ElhMmnk3FbhlENWfQzrJ9EZJs4FO4GzgPtKMt4bFCyBpJ16ZjTcK2F5Saf7oXRHhderNbEhpdY/oadL5MyeRegM/BLYDPhoRD7Y4lvVExNPAFNLxqAWSJgA/B3YFLiStVDAduII0RJe33qfK6r0W+L/AncBHNuacnyrxAvwjcE1225rUeyrd72+I0MysUF701Cr5DWFWh9JSPp6sUJVX37bc/IYwq8Pq1asBGDlyZMGRDFlORI2mNHd8s1plIiLXmnaSBjpWt7bGeVCN5jeEmTXD8LgMRIsdDqwe4DYgSRNz1LN3Y0M3Mxs63COqk6RtSOvIVRURD+SoZ3NglwGKPRoRLwwivI3hN4RZHZYvXw7A2LFjC41jCPPQnOXmN4RZHTxZYUAemjMzs6HJicjMzArlRGRmZoVyIjIzs0I5EZmZWaFafWE8M7NN0ogR/jitl6dvWyW/IcysGTx928zMhiYnIjMzK5QTkZlZA3R3d7+8uoINjhORmZkVyonIzMwK5URkZmaFciIyM7NCORGZmVmhnIjMzKxQXpPCzKwBOjs7iw6hbXmJH6vkN4SZNYOX+DEzs6HJicjMrAH6+vro6+srOoy25GNEZmYN0NPTA8Do0aMLjqT9uEdkZmaFciIyM7NCORGZmVmhnIjMzKxQTkRmZlYoJyIzMyuUV1awSn5DmFkzeGUFMzMbmpyIzMysUE5EZmYNsGTJEpYsWVJ0GG3JS/yYmTXAmjVrig6hbblHZGZmhXIiMjOzQjkRmZlZoZyIzMysUJ6sYGa2kXp6erjyyitZvHgxkyZNoquri87OzqLDahubVI9I0lxJC3OWPUJSSNqxyWEh6e2SbpH0oqRnJV0qaetB1jFgvJLGS5ot6QFJPZKWSLpd0l4b/yzMrD+zZs1iwoQJfOUrX+Hss8/m85//PBMmTGDWrFlFh9Y2NqlEBJwBHFR0EOUkbQcsAEYBXcCxwD7ATZIa/frvCnwSuCFr6whgJbBA0kcb3JbZsDdr1ixmzJhBb2/vett7e3uZMWOGk1FOw26tOUkjgTXA4cClwFsi4k9NbO884LPAxIhYnm3bC7gLODgi5ues5wgGiFfSWODFiFhTtm0E8DvgmYjI0zMaXm8Iszr19PQwYcKEDZJQuY6ODhYtWsSYMWNaGNmQNXTWmpM0SdJ1khZLWinpSUnXSBohaXI2/PRxSRdJek7SMknnSdpM0m6S7pXUK+l3kvarqHu9oTlJE7P6jpF0tqRFwEvA2Cqx7SrpGUnzJW2RbTtB0u8lrchieUDSYHpdBwA/LiUhgIi4G3gSOHAQ9QwYb0QsL09CWVtrgN8AEzamLTNb37x582omIUg9o3nz5rUoovZVxGSFm4DlwBeBpaQPyI+wflL8FjCfNMy0FzCDFOs+wDnAU9m2+ZK2j4ilA7Q5HbgfOArYjDRctR5J+wLXAj8Ajo2ItZIOBc4FTgfuIQ2v7QLkOr4jaRSwA3BJP7t/B+yUp54qdW8Qb5VymwN7Ag/X25aZbai7u7uh5YazliYiSeOAtwAHRsSNZbuuzPaX7t8RESdlv98qaX/gOOD9EXFvVrYbeAjYH7hsgKafAQ6KsnHIsrbIEs6lwFkR8fWyx+0JPBwRp5dtuznHUy3ZitQdXdbPvueAtw6irpfViLc/pwFvAA6tpy0z69/48eMbWm44a3WP6FngceAsSa8DFkTEH/sp95OK+48Ak0pJqGwbwBtztHt9VD8YdiJwNHB8RHynYt/9wDGSLiBNAPhFRPTlaK+klO36a7vqeOkATqR6vOs3IH0aOAU4IyLuqVHuKFJvkTlz5jB16tQ6QzMbPqZMmcLo0aPp66v+kTB69GimTJnC0qUDDdps+saNG1d1X0sTUUSEpA+SvqXPBraR9ARwTsWHamUPYhVpOK+8rlVZr2aLHE3X6hsfQhrqu7affZdn9X8OOAZYLelm4KSIWJij3WWkJNTfUN5WpF7RYNWK92WSPgbMBb4bEd+oVTYiLgYuLt2tIyazYWfcuHGceuqpzJgxo2qZU089lR122KGFUbWnlk9WiIjHI+IzwLbAu4A7gAslfbiZzdbYdzBpAsMCSa9f70HJRRGxOzCONNNud+CqXI2m3tNC4O397N4J+J889eSNt0TS3wPXANeRek9m1gTTp09n5syZdHR0rLe9o6ODmTNnMn369IIiay+FnUeUfcj/BigdC9q5oFCeAiaTXos7JfU7oBsRyyLiKuBqBhfrjcD+kl4+zVrS+4Dts30NjVfSnqRhxNuBwyJiXR1tmFlO06dPZ9GiRZx77rmcfPLJfPe732XRokVOQoPQ6skKuwDnk3oUfyLNYDuCdF7PHcCWrYynJCK6JU0mfXgvkDQlIhZJuhh4AfglsBiYBEwFbhlE9ecAhwE3SpoNdAJnA/eReiyNjPdtwI9JsxHPAXYtn5QREf9ZT3tmVtuYMWP41Kc+BXhyQj1aPVnhadL5MyeRZnKtBH4LfDQiHsw+XAsREU9LmkLZhzvwc9LJqFNJCWQRcAVQ85hLRb1PZXXNIR3XWUXqsXx5Y3orVeLdg3TsaSvgzn4eVu8ECTMbQK2D8VbbsFtZwQbkN4SZNcPQWVnBzMysnC8DUSelgy+b1SpTudxOjboG+jusrXEelJkNAcuXLwdg7NixhcbRjtwjqt/hwOoBbgOSNDFHPXs3NnQza7QVK1awYsWKosNoS+4R1e9HwG4NqGdRjnoebUA7ZmZDkhNRnSLiWdKSRRtbzyrggY2PyMysPXlozszMCuVEZGZmhXIiMjOzQvkYkZlZA4wY4Y/TenllBavkN4SZNYNXVjAzs6HJicjMzArlRGRm1gDd3d10d9e6GLRV40RkZmaFciIyM7NCORGZmVmhnIjMzKxQTkRmZlYoJyIzMyuU16QwM2uAzs7OokNoW17ixyr5DWFmzeAlfszMbGhyIjIza4C+vj76+vqKDqMt+RiRmVkD9PT0ADB69OiCI2k/7hGZmVmhnIjMzKxQTkRmZlYoJyIzMyuUE5GZmRXKicjMzArllRWskt8QZtYMXlnBzMyGJiciMzMrlBORmVkDLFmyhCVLlhQdRlvyEj9mZg2wZs2aokNoW+4RmZlZoZyIzMysUE5EZmZWKCciMzMrlBORmdlQEAHXXw/TphUdSct5ZQWr5DeEWR2WL18OwNixYwf3wHXrUgL62tfg6adhxQpYu7bR4TVET08P8+bNo7u7m/Hjx9PV1UVnZ2feh1ddWcGJyCr5DWHWCuvWwfz5cPLJsHgxvPhi2v6qVw3JRDRr1ixmz55Nb2/vy9s6OjqYNm0a06dPz1OFE5Hl5jeEWTOtWwfXXpt6QEuXvpKASoZgIpo1axYzZsyoun/mzJl5ktHwWGtO0lxJC3OWPUJSSNqxyWEh6e2SbpH0oqRnJV0qaetB1pEr3qzu30t6PmvvIUlfkrTZxj0LM6tl9erVrF69unqBtWvhqqvgTW+CI4+EhQs3TEJDUE9PD7Nnz65ZZvbs2Tz//PN1t7GpraxwBnB+0UGUk7QdsAB4BOgCxgLnADdJel9ErGtwk6OAC4DHSL2b/UivyY7ACQ1uy8wyS5cuBWD8+PHr71i7Fq6+Og3BLVs2cPKJgAceaFKUgzfvhhvWG47rT29vL/PmzePII4+sq41NKhFFxGMDlZE0EmjlWhxfBUYCH4uI5VkMi4C7gI8D8xvZWEQcUrHpliwZHokTkVnrlHpAJ58My5fn7/10dMA++zQ1tMHoXrkyX7nu7rrbaHkikjQJ+CbwXmAMsBj4FfAp4H3AncBBwIeBT5DGFecCXwH+FjgPeBewEDgpIn5WVvdcYHJETMzuTwSeAI4FJgKHAa8HtqkS267AzcDPgU9HxEpJJwBfyB6/ktTTmBUR1+V8ygcAPy4lIYCIuFvSk8CBbEQi6i/eKkWfpbXJ12z4WrMG/uM/4JRToKdn8MNvQ2y4bvzARVK5yp7gIBTRI7oJWA58EVgKTAA+wvrHq75F+oD+JLAXMIMU6z6kYa2nsm3zJW0fEUsHaHM6cD9wFLAZKaGsR9K+wLXAD4BjI2KtpEOBc4HTgXtIw167ALmO70gaBewAXNLP7t8BO+Wpp0rdG8Rbtk+k5/ka4O+Bw4Gz623LzPIZ+atfwbHHwvPPD7mEUq8u0lBKrcG5jo4Ourq66m6jpYlI0jjgLcCBEXFj2a4rs/2l+3dExEnZ77dK2h84Dnh/RNyble0GHgL2By4boOlngIOibIpgWVtkCedS4KyI+HrZ4/YEHo6I08u23ZzjqZZsRerRLetn33PAWwdR18tqxFuyP/Cj7PfIyp1Ro76jSEmaOXPmMHXq1HrCMhv21rztbSy99NJ0Z12jD/8W5/gf/pDZl19eff/xx7Nq1aqXj5P1Z9y4cVX3tbpH9CzwOHCWpNcBCyLij/2U+0nF/UeASaUkVLYN4I052r0+qs9TPxE4Gjg+Ir5Tse9+4BhJFwA3AL+IiL4c7ZWUsl1/bVedyjiAE6keb8k9wG5AJ6lH9BVJERH9zq+MiIuBi2vEamYD6O7uJjo7Gbf77jBnTrqtW5dOUG1zZwIdwGzW7xkN8jyiqlqaiCIiJH0QOI30nLaR9ARwTsWHamUPYhVpOK+8rlVZr2aLHE3XOop2CGmo79p+9l2e1f854BhgtaSbScemFuZodxnpg72/obytSL2iwaoVLwAR0QOUpt3cLmkV8M+SLoyIp+po08zyGjsWTj8dvvxlOO88OPfcNHEhb0IaNQpGDL15ZNOBL0Uwb+VKuidNYvyXv0xXVxdjxozZ6Lpb/mwj4nHgM9lxjHeShtwuzM7/adZXh1rf8g8m9QYWSPpARDxdFmsAFwEXSdoK2Jd0zOgq4D0DNhrRlz2vt/ezeyfSzLnBqhpvDQ+QjsHtQEpiZtZgGww9dXbCaafBSSelhPSv/5p6SH0DDKqsXJkmOwxBY0jTb3nve2HChIbVW1jazT7kfyPpJFKPY2fSUFirPQVMJs3WuzP7cN+gBxURy4CrJL2HNDSW143A4ZI6s54Kkt4HbJ/ta0q8FfYmJePH62jPzHIYOXJk/zvGjIFvfCMlpPPPh29+s3ZCkuAf/7F5gQ5BLV1ZQdIuku6U9AVJ+0jaj9TjWAPc0cpYymUf5JOBdaSexnZZvBdLOldSl6S9JH0emArcMojqzwHWAjdK+pCkT5Jmut0H5J0Cnjfe/SXNk3S4pCmSDpD0HeAk4KKIWFRPe2bWAFtuCTNmQHc3TJ+e7o8eXXRUQ0Krl/h5GniS9MF4I/BDYDvgoxHxYItjWU82xDWFdDxqgaQJpPNzdgUuBG4lDZNeQZoOnbfep8rqvRb4v6TezEc2ZlWFKvE+RvqbzgR+Cvw/4B3AZ0jnUplZkyxfvvzlFbhres1r4NRTYdEi+Od/Tj2mYZ6QvOipVfIbwqwOpZUFBn1iZ28vXHghzJoFq1enY0RDbNHTBvHq25ab3xBmdag7EZX09aWE9MtfptW5Nz1ORI1WtnpBVRGRa1kdSQNNGllb4zyoRvMbwqwOG52INn31XwZC0maS3ilp28bG1PYOB1YPcBtQth7eQPXs3djQzcyGjjzTt4N0Hsr+DG622KbuR6TVCzbWohz1PNqAdszMhqQBE1FErJP0v6QVHiwTEc+Sliza2HpW8coqCGZmw07e6dsXASdK2ryZwZiZtasRI0YwYgguzdMO8r5qWwJvBh6X9FPS2m3lB7UjIr7R6ODMzNrFttv6MHq9cs2akzTQiZcRETVnkFnb8Kw5M2sGT9+23PyGMLNmqH/6tpmZDay7u/vlc4lscHInIiUHSPpXSZdK2j7bvndp0U0zM7PByjVZIbsWz82ka/A8T5q8cAHwZ+CfSBd4O75JMZqZ2SYsb4/oHNIlud8LjGP9sb7bSJejNjMzG7S807cPBL4SEb+UVDk77klSkjIzMxu0vD2i11D9EtNbUGM2hJmZWS15E9GjwL5V9u0N/LYx4ZiZ2XCTd2ju34F/l9QDXJltGyvps8BxwFHNCM7MrF10dnYWHULbyn1Cq6SzgK+QhuFEOvFxHXB2RExvWoTWaj6h1cyaoTErK2TnDn0QeC1p5elbI+LxjQ7PhhInIjNrho1LRJL2Av4rIl7sZ99rgL+NiLs3KkQbKpyIzOrQ19cHwOjRowuOZMja6CV+7gR2qrLvrdl+M7Nhq6enh56enqLDaEt5E1Gt6dmvBtY2IBYzMxuGqs6akzQReFPZpndnw3DlRgFHkk5qNTMzG7Ra07cPB75BOmYQpLXlyntGkd1fAxzbrADNzGzTVisRzQUWkJLNHaRk8z8VZV4C/hARzzUjODMz2/RVTUQR8WfS6tpImkKaNfdCqwIzM7PhIdfKChFxF4CkXYC9gG2AiyLiaUk7As84SZmZWT3ynkf0auAK4B94ZVWF3SLivyTNJw3PndLUSK1VfB6RmTXDRp9HNAvYB5gKvK6iwp8A+9UdmpmZDWt5Fz39FDAjIq7s53pETwATGxqVmZkNG3l7RNsAv69Rx6sbE46ZWXtasmQJS5YsKTqMtpQ3ET0B7Fll3+6k6xWZmQ1ba9asYc2aNUWH0ZbyJqLLgVMkHQpsnm2LbFr3/wG+14zgzMxs05d31txmwA+AfySdxPpqYAXpMuH/ERGHNjNIaynPmjOrQ3d3NwDjx48vOJIhq+qsubznEa0FDpH076QZcqXrEf20dI6RmZlZPfLOmgMgIu4B7mlSLGZmNgwNKhFJEjCeNCS3Hl+p1czM6pErEUnaBvh34KAaj6k8v8jMbFjo6elh3rx5PPPMM+ywww50dXXR2dlZdFhtI+9kheuBKcAlwCPAqsoyEXFZo4OzQniyglktEfDww7DLLiAxa9YsZs+eTW9v78tFOjo6mDZtGtOnTy8w0CGn6mSFvImoBzghIuY2MCgbmpyIzPoTATfdBF/9Kjz6KNxyC7Puu48ZM2ZUfcjMmTOdjF6x0YnoCeCYiPhJI6OyIcmJyKxcBNx4I3zta/DUU9DbC2PG0PO97zHh8MPX6wlV6ujoYNGiRYwZM6aFAQ9ZG73o6QXAF7LJCmZmm7516+C662DSJDjsMPjDH1ISApCYd++9NZMQQG9vL/PmzWtBsO0t73lEcyRtB/yPpNuAZRsWiW80PLpBkjQXmBwRE3OUPQK4FHhLRPypyXG9HTgP+DvSCcE3Al8ezJVt64lX0t8B95K+iYyMCK8/YjaQdevg+uvTENzixfDii/0W634u379v6URXqy7vrLmPkC4V/mrgrf0UCaDwRAScAZxfdBDlsgS+gDTJowsYC5wD3CTpfRGxrkntjgQuAp4BXt+MNsw2KevWwfz5aQhuyZKqCahk/NZb56rWKy0MLO/Q3BzgfuCdwKsj4lUVtyExdTsiHouIX9cqI2lki4cYvwqMBD4WET+NiP8ADiUtIvvxJrcrvA6gWW3r1sHVV8Ob3wyf/Sw88cSASYi1a+m67DI6Bvgo6ejooKurq4HBbpryJqK/AmZGxG8jYvXGNChpkqTrJC2WtFLSk5KukTRC0mRJIenjki6S9JykZZLOk7SZpN0k3SupV9LvJO1XUfdcSQvL7k/M6jtG0tmSFpGGxsZWiW1XSc9Imi9pi2zbCZJ+L2lFFssDkg4axFM+APhxRCwvbYiIu4EngQMHUU+ueLPtbwamA8cAG/X3MttkrV0LV10FO+wAn/scLFw4cAIqefFFOpctY9oAk72mTZvmiQo55F1Z4dfAdg1q8yZgOfBFYCkwAfgI6yfFbwHzgU8CewEzSLHuQxrWeirbNl/S9hGxdIA2p5N6dEeRTrxdWVlA0r7AtaTFXY+NiLXZauPnAqeTljYaBewC5OqTSxoF7EA6/6rS74Cd8tRTpe4N4i3b/R1gXkTcLekD9bZhtklauzb1gE4+GZYty598+lGamD0bKJ+24POIBidvIjoeuEzSHyPi5/U2Jmkc8BbgwIi4sWzXldn+0v07IuKk7PdbJe0PHAe8PyLuzcp2Aw8B+wMDnUz7DHBQlM1VLx+dyxLOpcBZEfH1ssftCTwcEaeXbbs5x1Mt2Yo0PFY5uQPgOfo/3jagGvEi6TDg3cDbBlHfUaQkzZw5c5g6dWo9YZm1hxUrYKut4OKLG1Ld0cCne3u58Re/oDuC8dtvzwEHHMCWW27J0qUDfUcePsaNG1d1X95EdD0wBrhbUi+pR1MuImL7HPU8CzwOnCXpdcCCiPhjP+Uqz1d6BJhUSkJl2wDemKPd68uTUIUTSe+l4yPiOxX77geOkXQBcAPwi4joy9FeSSnb9dd2vcepTqRKvJK2JvXgTo2IxXkrjIiLgdJ/pc8jsk3fU0+lSQkPPpgSU47zKWsZBxzf2cnaH/2IEe9/f2NiHEbyJqLbacAHVESEpA8Cp5F6s9tkJ8ueU/GhWtmDWEVF8ouIVVmvZoMFWPtRa/7kIaShvmv72Xd5Vv/nyI63SLoZOCkiFuZodxnpdetvKG8rUq9osGrFO5PU+7ta0thsW+n16ZS0MiJqn/hgNhzssQfcfTfcd19KSPffv9EJSatXM2LEoNaRtkze84iOaFSD2Srdn8lmrr2TNOR2YTbJYEWj2qlstsa+g0m9gQWSPhART5fFGqQp0BdJ2grYl9TjuAp4z4CNRvRlz+vt/ezeCajnWk5V483qfAep51lpKalX9/E62jTbNO2+OyxYkBLRySfDr36VLyFtuSW88Y1QPoX7Na+Bt9Y12j7s5Z0113CR/AYoHQvauaBQngImk16LOyX1O+k/IpZFxFXA1Qwu1huB/SW9vBSvpPcB22f7GhnviaTFactvpeNn+5AmeJhZpd12gzvuSL2kD3wARo2CWlOzX/UqOPdcuOeel2/Lf/hDlr+qsI/UtjbY6xG9k3SAvb/rEV2e4/G7kE44vQr4E2kG2xHAGuAOYMvBxNMoEdEtaTJpCHKBpCkRsUjSxcALwC+BxcAkYCpwyyCqPwc4DLhR0mygEzgbuA+4rpHxZol9PVk5gLu8soLZAHbdFW67DX7969RDuvdeeOmldK7RAFasSAM6Y8eObXKQm568KyuMBX4M7FHalP0s778OmIiAp0nnz5wEvIE0jfq3wEcj4sGyD82Wi4inJU2h7MMd+DnwWVLy6QQWAVcwiFUkIuKprK45pOM6q0hDZF/emFUV+os3Ip6qtz4zK/Oud8Ett8BvfgOnnJJ6SjkTkg1e3tW3LwQ+QDpofw/pAnk9wJGkKc6HRMSDTYzTWsez5swqPfxwSkgLFqSE1NGRzkX60IdeLlJaU85L+lS10atv7wecCfxndv8vEbEgIj4D3AacsHHxmZkNYbvsAjffnCYzfPjD6STY17ym6Kg2GXmPEY0HHs9WG1jJ+sdy5gP/0fDIhrhs1l/NNfbyHpORNNDfYW2N86DMrFXe8Y50cbxly9JJsdYQeXtET/PK+mx/Jg3HlezYyIDayOGkddxq3QYkaWKOevZubOhmtlGchBoqb4/oXlLyuQn4PvCN7AN0DekDuZ5pyO3uR8BuDahnUY56Hm1AO2bWRD6ZtX55Jyu8GdguIu7JrnNzFmlB0tHAT4EvRUR/J1Fa+/EQoJk1Q9XJCrkSkQ0rfkOYWTPUP2tO0ubZdYEOaGxMZmZmORJRRKwiHQva4Bo+ZmaWdHd3v3wukQ1O3llz1wO+3q2ZmTVc3mkePwG+LWkeKSl1U3EsISLuaGxoZmY2HOSdNVdtgaUgHYCKiKh5cqe1DU9WMKuDl/gZUNXJCnl7RFMaFIiZmdl68l4Yr54LuJmZmQ3IV3EyM7NC5V6TQtLOpMtA9HdhvIiIv29kYGZm7aSzs3PgQtavvBfGew9wF7AQeAvwMLAV8FfAX0hXWzUzG7ZGjx5ddAhtK+/Q3Jmkyz28nTTz4XMRMRHYh3QphJlNic7MzDZ5eRPRLqRLZJem9m4GL587NBOY3fjQzMzaR19fH319fUWH0ZbyHiMaCfRGxDpJz5EulFfyKLBzwyMzM2sjPT09gIfo6pG3R/QYMCH7/WHgSEmvkvQq4LOkC+eZmZkNWt4e0U3AZOBK0vGiHwPPA2uB1wDHNyM4MzPb9NV1PSJJ7wIOJrswXkTc0ujArDBe4sesDl7iZ0D1LfEjaRxwGLAjsAy4NiJ+ExG/Bn7d0BDNzGxYqpqIJL0VuBvYtmzzKZK6IuKGpkdmZmbDQq3JCjNJF8ObDHQA7wDuA+Y0PywzMxsuqh4jkvQk8M8RcVnZtp2Bh4DXR8SS1oRoLeZjRGbWDFWPEdXqEU0gnSNU7tGssu0aEJSZmVnNRCTS9OxypQvkedVuMzNriFpDc+uAnwJLyzcDhwI3A8+VbY+IOLxZQVpLeWjOrA5LlqSjFdtuu+0AJYetuqZvPwn8dT/b/0xa/LScP7zMbFhbs2ZN0SG0raqJKFtd28zMrKl8rMfMzArlRGRmZoVyIjIzs0I5EZmZWaHyXgbCzMxqGDVqVNEhtK26LgNhmzS/IcysGepa4sfMzKzpciUiSZtL+oakRyT1SVpbcfOZXGY2rK1evZrVq1cXHUZbynuM6BzgWOAnwHzgpaZFZGbWhpYuTauh+Qqtg5c3EXUB34iIWc0MxszMhp+8x4heA/yymYGYmdnwlDcR/QjYq5mBNIKkuZIW5ix7hKSQtGOTY9pZ0kWSHpS0SlJds9LyxivpcEnXSvpzVn5uXYGbmbVI3qG5C4DLs0tDVF4CAoCIeLyRgdXpDOD8ooOosCvwEeAB0rG1PZvc3mHAtsCtwCea3JaZ2UbLm4hKw3KnAd+oUmazjY5mI0XEYwOVkTQSaOUsv++XLrcuaSbNT0T7RcS6rL0PNbktM7ONljcRHUmDTnSUNAn4JvBeYAywGPgV8CngfcCdwEHAh0nf6AXMBb4C/C1wHvAuYCFwUkT8rKzuucDk0iUsJE0EniDN+JtI6i28HtimSmy7knp8Pwc+HRErJZ0AfCF7/ErgMWBWRFyX5/mWkkIz9BdvM9szs/719PRw5ZVXsnjxYiZNmkRXVxednZ1Fh9U2ciWiiJjbwDZvApYDXyRd/XUCaeiq/HjVt0jTxD9JOjY1gxTrPqSp5E9l2+ZL2j4iyq8i25/pwP3AUaSe28rKApL2Ba4FfgAcGxFrJR0KnAucDtwDjAJ2AbYe5HNuuP7iLTgks2Fp1qxZzJ49m97e3pe3nXDCCUybNo3p06cXGFn7GNRac5IE7ET6IH4W+H0MYo0gSeOAtwAHRsSNZbuuzPaX7t8RESdlv98qaX/gOOD9EXFvVrYbeAjYH7hsgKafAQ4qj7WsLbKEcylwVkR8vexxewIPR8TpZdtuzvFUm6pGvGbWQrNmzWLGjBkbbO/t7X15u5PRwHInIkmfB2aSDoSXLJY0IyK+m7OaZ4HHgbMkvQ5YEBF/7KfcTyruPwJMKiWhsm0Ab8zR7vU1EuaJwNHA8RHxnYp99wPHSLoAuAH4RUT05WivmU6kerx1kXQUqbfInDlzmDp1aiOqNdukPf/885x55pk1y5x55pl8+tOfZsstt2xRVEPXuHHjqu7LlYiyb+AXA7cDVwBPk461HApcLKkvIn44UD0REZI+SJr0MBvYRtITwDkVH6rLKh66ijScV17XqqxXs0WOp9BdY98hpKG+a/vZd3lW/+eAY4DVkm4mHZtamKPdZqgVb10i4mLS3xe86KlZLjfccAN9fbW/l/b19XHnnXdy5JFHtiiq9pT3PKKvAT+IiA9GxGUR8bPs576kYbWT8zYYEY9HxGdIPat3AXcAF0r68GCDH4RaH64Hk6ZVL5D0+vUelFwUEbsD44DDgd2Bq5oW6cCqxmtmrdPdXev77eDLDWd5E9FbST2h/lyR7R+U7EP+N0DpWNDOg62jQZ4CJpNeizsl9btQVEQsi4irgKspLlbIGa+ZNVfeNeW89tzA8h4jegF4Q5V9b8j2D0jSLqQTTq8C/kSawXYE6byeO4BCBlIjolvSZNLQ4wJJUyJikaSLSc/tl6Rp5pOAqcAteeuWNJo0KxDgbdm2ruz+woh4oFHxZnXvRJpQAmmW3/Zl7d0VEUsG256Zbairq4sTTjhhvdlylTo6Oujq6qq635K8PaKfAGdKen/5Rkl7kiYwVE4uqOZp4ElSL+hG4IfAdsBHI+LBnHU0RUQ8DUwhHY9aIGkC6fycXYELSSsVTCf1AA8fRNWvBa7Jbgdn20r3j2twvAD/WFb/1qTeU+n+2+ttz8zW19nZybRp02qWmTZtGmPGjGlRRO0r1xVas2MRdwNvJg0NdZMmK7yB1LPZKyKeaWKc1jqerGA2CP2dR9TR0eHziDZU9QqtuS8Vng0xHQm8n/RN+zngLmDuEJjSbI3jRGQ2SM8//zyXXHLJeisruCe0gY1PRLa+7OTemuvrRUSuNe0kDXSsbu1gThzeSH5DmNVhyZJ0+HXbbbcdoOSwVTUR5T1GZBs6HFg9wG1A2Xp4A9Wzd2NDN7NG23bbbZ2E6lS1RyTpcdKyOA9lJ53W+qYcEfHmZgQ4VEnaBtihVpk8M+IkbU5av66WRyMi18zEBnCPyMyaoWqPqNaQ0F3A82W/+wOqTEQ8S1qyaGPrWUW6VpGZ2bDkY0RWyW8IszqUVlDwCaxVNecYUTY8ZWZmVrdciUjSP0n6atn9d0j6C2n17Qe85pmZmdUrb4/oS8CKsvtzSKthnwh0ki4cZ2ZmNmh515r7K7Lr/0jqJE0n/nhE3CzpWdIlHczMzAYtb49oM2Bd9vv7SAe0F2T3/5e0npqZmdmg5U1EfyRdkhvShdnKr1S6HWm5HzMzs0HLOzT3r8D3JR0ObAV8omzfFODhRgdmZtZOOjs7iw6hbeVKRBFxpaQ/A3sA90fE3WW7nyFd0sHMbNgaPXp00SG0rQFPaM2WoPkmcGVE3N+SqKxIPqHVzJqh/hNasyVojiZd7dPMzPrR19dHX5+viFOPvJMVfg28o5mBmJm1s56eHnp6eooOoy3lTURfBr4i6aPZdXjMzMwaIu+suWtIKyjcAKyRtJj1jyVERGzf6ODMzGzTlzcR3Y4PYpuZWRPknb59RJPjMDOzYcqXCjczs0LlTkTZpR/mSVoiaY2kxZKuluTZdGZmVrdcV2iVtBvpcuErSKsoPA28HvgY6fyivSLiwSbGaa3jY4Fm1gxVZ1znTUS3AWOAv4+IF8q2bwncBvRExL4NCNSK50RkZs2w0ZcK3wOYXZ6EALL73wT2rD82MzMbzvImooG+JftbtJkNa0uWLGHJkiVFh9GW8iaiXwGnZkNxL5PUAZwM/GejAzMzaydr1qxhzZo1RYfRlqqeRyTpceCgiHgIOJV0RdY/S7oJ6CZNVtifNFlhctMjNTOzTVKtE1onAq8GiIj7JO0BfB3YD9iadFXWO4AzIuK3TY7TzMw2UXmX+CEiHga6mhiLmZkNQwMdI/IkBDMza6qBekT/ImlpjnoiIg5vREBmZja8DJSI/gZ4KUc97jmZ2bA2apQvYl2vqisrSFoH7BER97U2JCuYv1SYWTNs9MoKZmZmTeFEZGbWAKtXr2b16tVFh9GWnIjMzBpg6dKlLF2aZ26XVao6WSEinKTMzKzpnGysOOvWwTXXwM9+VnQkZlYgJyJrvbVr4aqr4E1vgkMPhRkzio7IzAqUe4kfs422di1cfTWcfDIsWwYvvlh0RGY2BDgRWfOVekAnnwzLlzsBmdl6NqmhOUlzJS3MWfYISSFpxybHtLOkiyQ9KGmVpLpOGM0Tr6TxkmZLekBSj6Qlkm6XtFf9z2AjrFkDV1wB228PRx8Nf/mLk5CZbWBT6xGdAZxfdBAVdgU+AjxAWi6pmZdV3xX4JHAp6WKFmwPHAAskHRARNzWx7VesWQNXXgnTpsHzzw+cfFavhu7uloRm1izjAMaNKzqMtrRJJaKIeGygMpJGAq28jOL3I+KyrO2ZNDcR3QtMioiXn5+knwG/A74GNDcRrVkD3/8+TJ8OL7yQv/fzhz/Am9/c1NDMmm3kmjXw0EPw139ddChtp+WJSNIk4JvAe4ExwGLSpcg/BbwPuBM4CPgw8AnS+kRzga8AfwucB7wLWAicFBE/K6t7LjA5IiZm9ycCTwDHki70dxjpyrLbVIltV+Bm4OfApyNipaQTgC9kj18JPAbMiojr8jzfiFiXp1w9+ol3eT/tr5H0G+DdzYqD1avh8stTAurtHfzw24oVzYnLrJU6O9MxUBu0InpENwHLgS8CS4EJpKGr8uNV3wLmk4aZ9gJmkGLdBzgHeCrbNl/S9hEx0OnM04H7gaOAzUgJZT2S9gWuBX4AHBsRayUdCpwLnA7cQ7os+i6kK9QWqr94q5TbnNQLe7gpgVxzDRx3HPT1+fiPDWvrIuh94QW2LDqQNtTSRCRpHPAW4MCIuLFs15XZ/tL9OyLipOz3WyXtDxwHvD8i7s3KdgMPAfsDlw3Q9DPAQVG21HhZW2QJ51LgrIj4etnj9gQejojTy7bdnOOpNlWNePtzGvAG4NAa9R1FStLMmTOHqVOn5g9m4sTUG6qyirvZcBGbbcbaHXbgJS/z069xNY6ftbpH9CzwOHCWpNcBCyLij/2U+0nF/UdIxz7urdgG8MYc7V4f1a53AScCRwPHR8R3KvbdDxwj6QLgBuAXEdGXo71mOpHq8a5H0qeBU4AzIuKeauUi4mLg4tLdQUUzbhz88pfw1a/Cr3+dhtmclGyY6QGufvWrebyrix2nTKGrq4vOzs6iw2obLU1EERGSPkj6lj4b2EbSE8A5FR+qyyoeuoo0nFde16qsV7NFjqZrTck6hDTUd20/+y7P6v8cafbZakk3k45NLczRbjPUivdlkj5GOrb23Yj4RlMj2nNPuPde+M//hK99DR58cHAJaeRIUNVLlZgNabPWrmX22rX0vvQS/OAH8IMfcMIJJzBt2jSmT59edHhtoeXHiCLiceAzSlnknaQhtwuz83+addS61ifiwaTewAJJH4iIp8tiDeAi4CJJWwH7ko4ZXQW8p0mxDqRqvCWS/h64BriO1HtqjT32gLvvhvvuSwnp/vvzJaTXvQ6+/OXWxGjWQLNuvZUZN284Wt/b28uMbOkqJ6McIqLQG2nmXABfBSZnv+9TUWYu8Jd+HhvAzIpyC8vuT8zKfL6fxx6R7dsRGE8a6vs9MH6AeOcAvXU+15lk+a2Ox+aKl3Rc60Xgx8DIOtpqnPvvj5gyJWL06AgpIqWkDW/vfndDmzVrheXLl0dHR0dk/5f93jo6OqKnp6foUIeKqp87LV1ZQdIuku6U9AVJ+0jaj9TjWAPc0cpYykVENykJriP1NLbL4r1Y0rmSuiTtJenzwFTglrx1SxqdPb4LeFu2rSu71TWluka8byMloKWk2YW7StqjdKunrY3y7nfDHXekXtIHPgCjRnkIzjYZ8+bNo7e3t2aZ3t5e5s2b16KI2lerh+aeBp4ETiLN5FoJ/Bb4aEQ8KGlyi+N5WUQ8LWkKcDvpw30K6fycz5KSTyewCLgCGMwxl9eShsnKle5fRurpNCrePYCtstud/TysmCyw665w221pMsPJJ6fjSS+9lC4DYdamunOuBpK33HCm8AwnW1/z3xAPPZQS0t13w8qVKVHdf3/TmzVrpO9+97t8/vOfz1XuyCOPbEFEQ17VL8JORFapdW+Ihx+GU06BSZPgW99qWbNmjdDT08OECRNqDs91dHSwaNEixowZ08LIhqyqiWiTWn27lZSMqHUbRF0165E20QMru+wCN9/sJGRtqbOzk2nTptUsM23aNCehHJyI6nc4sHqA24Cy9fAGqmfvxoZuZo0wffp0Zs6cSUdHx3rbOzo6mDlzpqdu5+ShuTpJ2gbYoVaZiHggRz2bk9avq+XRiHhhEOFtDL8hzAbp+eef55JLLmHx4sVMmjSJrq4u94Q25GNElpvfEGZ1KM2OGz9+fMGRDFk+RmRmZkOTE5GZmRXKicjMzArlRGRmZoUq4gqtZmabHF9/qH6eNWeV/IYws2bwrDkzMxuanIjMzBqgr6+Pvr6+osNoSz5GZGbWAD09PQCMHj264Ejaj3tEZmZWKCciMzMrlBORmZkVyonIzMwK5URkZmaFciIyM7NCeWUFq+Q3hJk1g1dWMDOzocmJyMzMCuVEZGbWAEuWLGHJkiVFh9GWvMSPmVkDrFmzpugQ2pZ7RGZmVignIjMzK5QTkZmZFcqJyMzMCuVEZGZmhfKsOTOzBhg1alTRIbQtL/FjlfyGMLNm8BI/ZmY2NDkRmZk1wOrVq1m9enXRYbQlJyIzswZYunQpS5cuLTqMtuREZGZmhXIiMjOzQjkRmZlZoZyIzMysUD6h1czMcunp6WHevHl0d3czfvx4urq66Ozs3Oh6fUKrVfIbwqwO3d3dAIwfP77gSJpj1qxZzJ49m97e3pe3dXR0MG3aNKZPn56niqontDoRWSW/IczqUDqHaOTIkQVH0nizZs1ixowZVffPnDkzTzIaHolI0lxgckRMzFH2COBS4C0R8acmxrQz8CXg3cA7gJERUfUPUqOeI8gRr6RLgT2ACaRjgI8BlwAXRsTaHE1tOm8IM9toPT09TJgwYb2eUKWOjg4WLVrEmDFjalU1bJb4OQM4qOggKuwKfAR4EnigBe2NAi4APgH8A3AbcD4wpwVtm9kmZt68eTWTEEBvby/z5s2ru41NarJCRDw2UBlJI4FWXlz++xFxWdb2TGDPZjYWEYdUbLpF0nbAkcAJzWzbbDhbvnw5AGPHji00jkYrHftqVLn+tLxHJGmSpOskLZa0UtKTkq6RNELSZEkh6eOSLpL0nKRlks6TtJmk3STdK6lX0u8k7VdR91xJC8vuT8zqO0bS2ZIWAS8BY6vEtqukZyTNl7RFtu0ESb+XtCKL5QFJuXtdEbGuntcpj/7ireJZWpt8zYadFStWsGLFiqLDaLi8ky82ZpJGEUNzN5GOX3wR2A84hZQcymP5FtALfBL4N+DEbNvlwPdIQ07PAfMljcvR5nRgEnAUaehuZWUBSfsCC4DrgE9ExEpJhwLnAj8kDa8dCswDts79bJukv3jL9ilL7GMlHQwcjofmzKwOXV1ddHR01CzT0dFBV1dX3W20dGguSxpvAQ6MiBvLdl2Z7S/dvyMiTsp+v1XS/sBxwPsj4t6sbDfwELA/cNkATT8DHBRlMzPK2iJLOJcCZ0XE18setyfwcEScXrbt5hxPtalqxFuyP/Cj7PfIyp3RqvjMbNPR2dnJtGnTas6amzZt2kATFWpq9TGiZ4HHgbMkvQ5YEBF/7KfcTyruPwJMKiWhsm0Ab8zR7vVRfXrgicDRwPER8Z2KffcDx0i6ALgB+EVE9OVor5lOpHq8JfcAuwGdwN8DX5EUEdHv/EpJR5F6i8yZM4epU6c2PGiz4WJTXIH76KOPpre3l/PPP5++vlc+AkePHs0JJ5zA0UcfPeDzHjeu+uBVSxNRRISkDwKnAbOBbSQ9AZxT8aG6rOKhq4DlFXWtyno1tY6NlNQ6inYI8BRwbT/7Ls/q/xxwDLBa0s3ASRGxMEe7zVArXgAioodXZujdLmkV8M+SLoyIp/opfzFwcelug+M1GxZKB+trfeC2szPPPJNTTjllg5UVNqYnVNLyWXMR8TjwGaUs8k7SkNuF2SSDZh3pq/XhejDpQ3iBpA9ExNNlsQZwEXCRpK2AfUnHjK4C3tOkWAdSNd4aHiAdg9uBlMTMzAZtzJgxHHnkkQ2vt7DziCL5DVA6FrRzQaE8BUwmvRZ3Sup36kdELIuIq4CrKS5WyBlvhb1JyfjxJsZlNqyNGDGCESM2qTNiWqbVkxV2IZ1ceRXwJ2Az4AjS1OI7gC1bGU9JRHRLmgzcTuppTImIRZIuBl4AfgksJs28mwrckrduSaNJM+4A3pZtK00vWRgRgz7JtUa8+wOfJU1UeJL0en6YdPznoohYNNi2zCyfbbfdtugQ2lar0/fTpA/Ik4A3kKZR/xb4aEQ8mH24FiIinpY0hbIPd+DnpA/2qaQD/4uAK4BvDKLq1wLXVGwr3b+MlIgbFe9jpJ7SzKzd5cAfgc+QpqCbmQ05m9Rac9YQfkOYWTMMm7XmzMwK0d3dvVHL3AxnPrJWp2zW32a1ykRErmV1JA30d1hb4zwoM7O25h5R/Q4HVg9wG5CkiTnq2buxoZuZDR3uEdXvR6TVCzbWohz1PNqAdszMhiQnojpFxLOkJYs2tp5VtOY6RWZmQ5KH5szMrFBORGZmVigPzZmZNUBnZ2fRIbQtn9BqlfyGMLNm8AmtZmY2NDkRmZk1QF9f33oXjbP8fIzIzKwBenp6gHTVUhsc94jMzKxQTkRmZlYoJyIzMyuUE5GZmRXKicjMzArlRGRmZoXyygpWyW8IM2sGr6xgZmZDkxORmZkVyonIzKwBlixZwpIlS4oOoy15iR8zswZYs2ZN0SG0LfeIzMysUE5EZmZWKCciMzMrlBORmZkVyonIzMwK5VlzZmYNMGrUqKJDaFte4scq+Q1hZs3gJX7MzGxociIyM2uA1atXs3r16qLDaEtORGZmDbB06VKWLl1adBhtyYnIzMwK5URkZmaFciIyM7NCORGZmVmhnIjMzKxQTkRmZlYor6xglfyGMKtD6RyikSNHFhzJkFV1ZQUnIqvkN4SZNYOX+DEzs6Fpk0pEkuZKWpiz7BGSQtKOTY5pZ0kXSXpQ0ipJdfU46olX0t9JWpc9ziutmzXR8uXLWb58edFhtKVNKhEBZwAHFR1EhV2BjwBPAg+0qlFJI4GLgGda1abZcLZixQpWrFhRdBhtaZNKRBHxWET8ulYZSSMlVR2rbILvR8QbI+Ig4I4WtvtV0pjs91rYppnZoLU8EUmaJOk6SYslrZT0pKRrJI2QNDkbRvp4Npz1nKRlks6TtJmk3STdK6lX0u8k7VdR93pDc5ImZvUdI+lsSYuAl4CxVWLbVdIzkuZL2iLbdoKk30takcXygKTcva6IWFfP65RHf/Fm298MTAeOAbwcsJkNaUUcN7gJWA58EVgKTCANXZUnxW8B84FPAnsBM0ix7gOcAzyVbZsvafuIGGjJ2+nA/cBRwGbAysoCkvYFrgV+ABwbEWslHQqcC5wO3AOMAnYBth7kc264/uIt2/0dYF5E3C3pA4UEaDaM9PT0cOWVV7J48WImTZpEV1cXnZ2dRYfVPiKiZTdgHGl68AFV9k/O9n+vYvt/ZdvfV7Ztl2zb4WXb5gILy+5PzMr8F9lU9bJ9R2T7dgQOBVYBp1eU+Tfgvxr4/Geml7yuxw4Yb1buMOA54LXZ/dOyx43I2ZaZDcLMmTOjo6Mjsv+zAKKjoyNmzpxZdGhDTdXPnVb3iJ4FHgfOkvQ6YEFE/LGfcj+puP8IMCki7q3YBvDGHO1eH1H1hKkTgaOB4yPiOxX77geOkXQBcAPwi4joy9FeM51IlXglbU3qwZ0aEYvzVijpKFJvkTlz5jB16tTGRWu2CZszZw6zZ8/eYHtvby8zZsygt7eXk046qYDIhp5x48ZV3dfSRBQRIemDpG/ps4FtJD0BnFPxobqs4qGrSMN55XWtyuYcbMHAumvsO4Q01HdtP/suz+r/HNnxFkk3AydFxMIc7TZDrXhnkmbJXS1pbLat9Pp0SloZEb2VD4qIi4GLS3cbG67Zpqmnp4dvf/vbNct8+9vf5pRTTmHMmDEtiqo9tXyyQkQ8HhGfAbYF3kWaSXahpA83s9ka+w4mTWBYIOn16z0ouSgidicNKx4O7A5c1bRIB1Y1XmAn4B2knuey7HZytm8p6XiSmTXAvHnz6O3d4Hvdenp7e5k3b16LImpfhU3fzj7kfwOU+q07FxTKU6RjU68C7pQ0vr9CEbEsIq4Crqa4WKF2vCcCUypul2X79iFN8DCzBujurjXQMvhyw1lLh+Yk7QKcT+pR/Ik0g+0IYA2pZ7RlK+MpiYhuSZOB20k9jSkRsUjSxcALwC+BxcAkYCpwS966JY0mzQoEeFu2rSu7vzAiBn2Sa7V4s8Re2f7k7Ne7ImLNYNsys/6NH9/vd9a6yw1nre4RPU1aYeAk4Ebgh8B2wEcj4sEWx7KeiHia1INYRfpwnwD8nLQywoXAraRp4FeQhujyei1wTXY7ONtWun9cg+M1sxbp6uqio6OjZpmOjg66urpqljGvvm0b8hvCLKdZs2YxY0b1Ee+ZM2cyffr0FkY0pFVd0cYLYZqZ1amUZGbPnr3exIWOjg6mTZvmJJSTe0R1ytar26xWmbzHZHKsjL22xnlQjeY3hNkgPf/881xyySXrrazgKdsb8IXxGk3SEcCltcpExICLq0qaCDwxQLEpEbEgb2wbyW8IszqUZsd5ckJVHpprgh8BuzWgnkU56nm0Ae2YmQ1JTkR1iohnSSeObmw9q2jhdYrMzIYaJyIzswbwatv18zEiq+Q3hJk1Q9VjRJvUFVrNzKz9OBGZmTVAX18ffX1FXyWmPfkYkZlZA/T09AAwevTogiNpP+4RmZlZoZyIzMysUE5EZmZWKCciMzMrlBORmZkVyonIzMwK5ZUVrK1IOioiLi46jk2FX8/G8utZH/eIrN0cVXQAmxi/no3l17MOTkRmZlYoJyIzMyuUE5G1G4+/N5Zfz8by61kHT1YwM7NCuUdkZmaFciKytiRpkqTzJT0s6UVJ3ZJulPTOomNrV5JOkvSj7LUMSacVHVM7kPRGSfMk9Uh6XtJ8SX9VdFztxInI2tW+wBTgMuBjwDHAtsCvJO1aZGBt7J+A1wLXFxxH25A0GrgDeBtwODAVeAtwp6SOImNrJz5GZG1J0jjg2Sh7A0vqBBYCP4qIzxQVW7uS9KqIWCdpBLAa+JeIOK3gsIY0SScAc4C3RsSfsm07AH8EvhYRc4qMr124R2RtKSKWRsW3qIjoAf4ATCgmqvYWEeuKjqENHQD8ZykJAUTEE8DPgQMLi6rNOBHZJkPS1sDOwO+LjsWGjbcD/93P9t8BO7U4lrblRGSbkgsAAd8qOA4bPrYGlvWz/TlgqxbH0raciGxIkLRPNlNroNuCKo+fBnwaOK58mGS42tjX0walvwPtankUbWxE0QGYZX4B/HWOcn2VGyR9ATgTmBER32t0YG2q7tfTBmUZqVdUaSv67ylZP5yIbEiIiD7gkcE+TtJU4ELg3IiY1fDA2lS9r6cN2u9Ix4kq7QT8T4tjaVsemrO2Jekg4FLgkoj4StHx2LB0I7CHpDeVNkiaCLw322c5+Dwia0uS9gJuIX3rPA4on3r8UkT8upDA2pikdwMTSV9QrwKuAa7Odt+c9bKsTHbS6kPACmAG6XjRGcCWwC4R8WKB4bUNJyJrS9nyM9+osvvPETGxddFsGiTNJa0O0J8dImJh66JpH9lyPucBHyRNUrgdONGvV35ORGZmVigfIzIzs0I5EZmZWaGciMzMrFBORGZmVignIjMzK5QTkZmZFcqJyMzMCuVEZEOepCNqrB69T9Hx5SWpQ9I0Sf8l6QVJKyU9KunfJO1YYFxHSDqyBe38jaTTsutG5X1Mh6RuSQeXbZsr6S/NibLfGOZKWlhl388k/b8WxlL6X5iY3ZekX0v6aqtiaAYvemrt5BNA5QdQWywsKWk8cBuwHfBvwL3AKtLimEeS1iZ7V0HhHUH6LGj2yuV/Q1oN4wrS9Xry+DKwFJjfpJjqJmkMMBn4h6JiiIiQdDrwPUnfjYi8r+uQ4kRk7eQ3jb7WkKSRwJrKy443wfeB8cDuEfHHsu13SroQX1Z6A5I2B74EnNaCv0899id9mbit4DhuBFYCnwfOLjiWunhozjYJkkZKmilpoaRV2c+ZWaIplZmYDWscI+lsSYuAl4Cx2f6DJP1c0ouSnpd0n6QDyh4/Ihtae0TSS5IWSTpX0hYDxLY78PfAmRVJCEjfaiPi+jqfy9GSTs+Gr5ZL+pGkN1S0/+ls+OZFST2Sfivp6GzfAmBv4L2VF8uTtK2kiyT9QVKfpP+VdKWkCRX1n5Y97i2Sfpy182dJX5f0qqzMEaSV0gH+WNbWxBov3UGka/1cVev1zer/rKTVkk4p27atpAuzuF/Kfn5f0quz/Ttm95+QtELS45K+IynvlVU/Dvw0Il7K6lsg6V5JH5L0m6zOX0t6T/beOTP7Oz2XDfd1VDyH8ZIul7Q0i/dhSYcNFERErCUtUPv5nHEPOe4RWTvZTFL5ezayf0KAy4B/JF0g715gT9JqyG8iXbm13HTgfuAoYDNgpaQvAd8Grict/Pki8Lek1ahLrgA+BnyTVy48d0ZW5mCqKx3HyntZgME8l2lZLEcCrwXOBX5ASi5Iel8W97eBr5K+fL6NLPkCx2T7NwOOzrY9n/3cmvRNexqwhDSs+GXg55LeFhErK2K5jpRsziO9Tv8C/G+27cfAzOx5lA+xdtd4HT4E/D4iltYoU7o6778A/xQRc7NtW2Wvy9ZZuw9nr8+BwOakLyDbZXGcSLqI3ZuAU4GbSa95rTY3z+I7pmLXjsA5wCzSe+hs0t/9RtLn7RGk9805wGLga1l9HcBdpAvqnUp63Q4Dvi9pdERcXCse4G7gS5LeFBGPD1B26IkI33wb0jfSP2/0c7s3279zdv+0iseVluXfJbs/Mbv/X2QL/mbbxwAvAPNrxPD+7LGfqdh+aLb9b2o89jtZmVfneK6DfS53VZT7SrZ9u7L7zw3Q5oLSazlAuc2AN2b1H1S2/bRs22cryv8WuKWfv+OOOf/uvwd+0M/2uaQE8irgAqAX2L+izOnAWuBdg3ifjQDel8X4ror2FlaU/TBpWG5sxeu4GnhT2bYDsvpuq3j8fOCJsvvHZeUmV5S7jZSwNqt4DSdWlHtztv3Tjf7/a8XNQ3PWTg4Cdiu7fS7bvlf284qK8qX7e1dsvz6y/97M3wGvAWp96/wQ6YPn2myYZUTWO7ulIoaNNdjn8uOK+7/Nfv5V9vN+YCtJV0j6qKSxgwlG0hclPSTpRWAN8GS26639FK+M5b/L4qjHdqSeWH9GAP9B6iHuExGVbe8L3B81rkslaXNJp2ZDrStISeSebHd/z6/cx0lfApZXbP9DrN8jKV0l92cV5R4B3iBJ2f29gKciYkFFuSuAbUmTWmopvU7bDVBuSHIisnby3xHxQNnt0Wx7aTpw5TDP0xX7qVJum+xnrSnBryUN6bxI+sAq3RZX1NGf/81+bl+jTMlgn0vlLKmXsp9bAETEXaShsDeShs6WSLpN0i4DBZINV15I+lb+D8DuwB7l9eeIpebxswFswSvPp9IY0mSBXwD39bN/G2r/PQFmk3pzV2R17c4rM+Cqxp0lj4+RhnErLau4v6rG9hGkXiakv2t/w5TV/u6VVmQ/Rw1QbkhyIrJNQekD8PUV20v3n63YXjkDq3QMYgLVPUs6XrJbldtFNR5bmlX1sRplSgb7XAYUEfMiYm/S8YeDSLP3flqaSFDDIcDtEfHliLglIu7nlcTbCs+SYu7Pc6TkMQX4YcWxQ0h/01p/T0jP7/KImBkRd2TPb3mOuPYgvYY35Cib13Ns+DeH/H/3UqKqeTxtqHIisk3BXdnPQyq2H5r9vHuAx/+C1NM5qkaZn5K+JXdW9MpKt0XVHhgR95Gu2nmqqpy4Kqk0fXtjn0tVEfFiRNxESprjeaUX9xL9f5MeTer1lftsve3zSu8m77f2R0gTCPqVDWN9OLv9R0UyugXYXdI7a9Rf7/P7OPBARDTypNq7SEN1763Y/mlS8v/9AI/fIfv5aM1SQ5RnzVnbi4jfSfohcFr2YfQL0qynfwZ+GBEPD/D4F7KZVxdIupY06+wF0gmYKyPigohYkLUxT9Ic0nDQOtKkgY8AJ0fEH2o0M5XUM7pf0gW8ckLr20gz3kYCN2zsc6mkdLLj64A7gUXAG4DjSedklY4r/A9wjKRPAo8BL2TDnj8FTpZ0avZ8PwB0Dab9CqWTj4+VdBkpCTwcEauqlL8bOFHSqyJiXX8FIuIeSR8CfgJcJemQiFhNmrn3aeA2STNJx87GkWbNfSEiXsie3+GSfgv8iTQs93c5nseBpPPCGmkucAIwX9J00rDioaTLjx8dr8wOreY9pNfzPxscV2sUPVvCN98GupFjthXpg3wm8GfSP+Sfs/sjy8pMzOr5fJU6uoBfkcbbn89+/2jZ/leRPiweIg3T9WS/n03qKQ30PF5Dmpr7a9JMr5dI32DPZ/2ZVnU/F9KZ/i/PviINX/2MdPzhJdLxqu+SzarLyryeNGX5heyxC7Lto0gz/pZk+24iffNeb1Yfr8yaG1ERy1w2nG32DeAp0oy2DWZ/VZT966zM3v3U+5eKbXtmf4/rgc2zba8lTUDpJiX9/yVNjX91tn8cacLDsuz2A9IwawBH9Pc8SF8cAnh7P/EuoGL2YY2/0wavGamX+n3S8NpLpCnnh1X5X5hYsf1WYF7R/6v13pQ9CTOzISc7ufZPETEkTtbMTpj9XES8pehYSiRtR5rNuF9E3F50PPVwIjKzISs7ZnIbqTf8VNHxDEWSzgPeGREfKDqWenmygpkNWRHxc+D/kG/q+3DVDRxbdBAbwz0iMzMrlHtEZmZWKCciMzMrlBORmZkVyonIzMwK5URkZmaF+v8qsYGYGswoowAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 360x576 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import matplotlib.pyplot as plt\n",
"import seaborn as sns\n",
"%matplotlib inline\n",
"\n",
"plt.rc('font', size=16)\n",
"\n",
"\n",
"ax = plt.figure(figsize=(5,8))\n",
"#add start points\n",
"ax = sns.stripplot(data=df, \n",
" x='before', \n",
" y='parameter', \n",
" orient='h', \n",
" order=df['parameter'], \n",
" size=10, \n",
" color='black')\n",
"ax.grid(axis='y', color='0.9') \n",
"#add arrows to plot only if the parameter changed by more than 1e-3 kcal/mol\n",
"for i in range(len(df.index)):\n",
" term = df.iloc[i]\n",
" if abs(term[\"change\"]) > 1e-3:\n",
" if term[\"after\"] > term[\"before\"]:\n",
" arrow_color = '#347768'\n",
" elif term[\"after\"] < term[\"before\"]:\n",
" arrow_color = 'red'\n",
" else:\n",
" arrow_color = 'black'\n",
" ax.arrow(term[\"before\"], \n",
" i, \n",
" term[\"change\"], \n",
" 0, \n",
" head_width=0.3, \n",
" head_length=0.2, \n",
" width=0.1, \n",
" fc=arrow_color, \n",
" ec=arrow_color) \n",
"ax.axvline(x=0, color='0.9', ls='--', lw=2, zorder=0)\n",
"ax.set_xlabel('Force Constant (kcal/mol)') \n",
"ax.set_ylabel('Torsion Parameter') \n",
"sns.despine(left=True, bottom=True)"
]
},
{
"cell_type": "markdown",
"id": "12baaceb",
"metadata": {},
"source": [
"## Conclusion\n",
"\n",
"To recap in this blog post we have:\n",
"- Created and configured a general BespokeFit workflow\n",
"- Built a molecule specific optimization schema\n",
"- Spun up the BespokeFitExecutor and workers\n",
"- Fragmented the molecule\n",
"- Automatically generated torsiondrive reference data using QCEngine\n",
"- Automatically generated bespoke SMARTS patterns for the molecule\n",
"- Optimized the torsion parameters using ForceBalance\n",
"\n",
"Hopefully this demonstration has shown just how easy it is to set up a BespokeFit general fitting pipeline with a custom\n",
"configuration and train bespoke parameters to reference data generated on the fly from just a few imports. Better yet\n",
"this whole process can now be routinely applied to new molecules from the CLI using our serialized workflow.\n",
"\n",
"`openff-bespoke executor run --input BrCO.sdf --spec-file workflow.json`"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "faac0051",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.9.7"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment