{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Run the fusion EasyVVUQ campaign using SC\n",
"\n",
"Run an EasyVVUQ campaign to analyze the sensitivity of the temperature\n",
"profile predicted by a simplified model of heat conduction in a\n",
"tokamak plasma.\n",
"\n",
"This is done with SC."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:31.713466Z",
"start_time": "2021-06-04T20:24:40.688493Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:33.409694Z",
"iopub.status.busy": "2025-07-19T22:23:33.409620Z",
"iopub.status.idle": "2025-07-19T22:23:35.128951Z",
"shell.execute_reply": "2025-07-19T22:23:35.128693Z",
"shell.execute_reply.started": "2025-07-19T22:23:33.409687Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Volumes/UserData/dpc/GIT/EasyVVUQ/env_3.12/lib/python3.12/site-packages/chaospy/__init__.py:9: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81.\n",
" import pkg_resources\n"
]
}
],
"source": [
"# import packages that we will use\n",
"%matplotlib inline\n",
"import os\n",
"import easyvvuq as uq\n",
"import chaospy as cp\n",
"import pickle\n",
"import time\n",
"import numpy as np\n",
"import pandas as pd\n",
"import matplotlib\n",
"if not os.getenv(\"DISPLAY\"): matplotlib.use('Agg')\n",
"import matplotlib.pylab as plt\n",
"from IPython.display import display\n",
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.272217Z",
"start_time": "2021-06-04T20:25:31.716919Z"
},
"code_folding": [
0,
2,
4
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.129386Z",
"iopub.status.busy": "2025-07-19T22:23:35.129227Z",
"iopub.status.idle": "2025-07-19T22:23:35.462597Z",
"shell.execute_reply": "2025-07-19T22:23:35.462303Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.129379Z"
}
},
"outputs": [],
"source": [
"# we need fipy -- install if not already available\n",
"\n",
"try:\n",
" import fipy\n",
"except ModuleNotFoundError:\n",
" ! pip install future\n",
" ! pip install fipy\n",
" import fipy"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.279881Z",
"start_time": "2021-06-04T20:25:34.275176Z"
},
"code_folding": [
0,
2
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.463084Z",
"iopub.status.busy": "2025-07-19T22:23:35.463019Z",
"iopub.status.idle": "2025-07-19T22:23:35.465012Z",
"shell.execute_reply": "2025-07-19T22:23:35.464783Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.463076Z"
}
},
"outputs": [],
"source": [
"# routine to write out (if needed) the fusion .template file\n",
"\n",
"def write_template(params):\n",
" str = \"\"\n",
" first = True\n",
" for k in params.keys():\n",
" if first:\n",
" str += '{\"%s\": \"$%s\"' % (k,k) ; first = False\n",
" else:\n",
" str += ', \"%s\": \"$%s\"' % (k,k)\n",
" str += '}'\n",
" print(str, file=open('fusion.template','w'))"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.304016Z",
"start_time": "2021-06-04T20:25:34.282193Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.465428Z",
"iopub.status.busy": "2025-07-19T22:23:35.465351Z",
"iopub.status.idle": "2025-07-19T22:23:35.467656Z",
"shell.execute_reply": "2025-07-19T22:23:35.467474Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.465423Z"
}
},
"outputs": [],
"source": [
"# define parameters of the fusion model\n",
"def define_params():\n",
" return {\n",
" \"Qe_tot\": {\"type\": \"float\", \"min\": 1.0e6, \"max\": 50.0e6, \"default\": 2e6},\n",
" \"H0\": {\"type\": \"float\", \"min\": 0.00, \"max\": 1.0, \"default\": 0},\n",
" \"Hw\": {\"type\": \"float\", \"min\": 0.01, \"max\": 100.0, \"default\": 0.1},\n",
" \"Te_bc\": {\"type\": \"float\", \"min\": 10.0, \"max\": 1000.0, \"default\": 100},\n",
" \"chi\": {\"type\": \"float\", \"min\": 0.01, \"max\": 100.0, \"default\": 1},\n",
" \"a0\": {\"type\": \"float\", \"min\": 0.2, \"max\": 10.0, \"default\": 1},\n",
" \"R0\": {\"type\": \"float\", \"min\": 0.5, \"max\": 20.0, \"default\": 3},\n",
" \"E0\": {\"type\": \"float\", \"min\": 1.0, \"max\": 10.0, \"default\": 1.5},\n",
" \"b_pos\": {\"type\": \"float\", \"min\": 0.95, \"max\": 0.99, \"default\": 0.98},\n",
" \"b_height\": {\"type\": \"float\", \"min\": 3e19, \"max\": 10e19, \"default\": 6e19},\n",
" \"b_sol\": {\"type\": \"float\", \"min\": 2e18, \"max\": 3e19, \"default\": 2e19},\n",
" \"b_width\": {\"type\": \"float\", \"min\": 0.005, \"max\": 0.025, \"default\": 0.01},\n",
" \"b_slope\": {\"type\": \"float\", \"min\": 0.0, \"max\": 0.05, \"default\": 0.01},\n",
" \"nr\": {\"type\": \"integer\", \"min\": 10, \"max\": 1000, \"default\": 100},\n",
" \"dt\": {\"type\": \"float\", \"min\": 1e-3, \"max\": 1e3, \"default\": 100},\n",
" \"out_file\": {\"type\": \"string\", \"default\": \"output.csv\"}\n",
" }"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.320552Z",
"start_time": "2021-06-04T20:25:34.305858Z"
},
"code_folding": [
0,
2,
17,
21,
28
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.468065Z",
"iopub.status.busy": "2025-07-19T22:23:35.467981Z",
"iopub.status.idle": "2025-07-19T22:23:35.470305Z",
"shell.execute_reply": "2025-07-19T22:23:35.470129Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.468059Z"
}
},
"outputs": [],
"source": [
"# define varying quantities\n",
"def define_vary():\n",
" vary_all = {\n",
" \"Qe_tot\": cp.Uniform(1.8e6, 2.2e6),\n",
" \"H0\": cp.Uniform(0.0, 0.2),\n",
" \"Hw\": cp.Uniform(0.1, 0.5),\n",
" \"chi\": cp.Uniform(0.8, 1.2),\n",
" \"Te_bc\": cp.Uniform(80.0, 120.0),\n",
" \"a0\": cp.Uniform(0.9, 1.1),\n",
" \"R0\": cp.Uniform(2.7, 3.3),\n",
" \"E0\": cp.Uniform(1.4, 1.6),\n",
" \"b_pos\": cp.Uniform(0.95, 0.99),\n",
" \"b_height\": cp.Uniform(5e19, 7e19),\n",
" \"b_sol\": cp.Uniform(1e19, 3e19),\n",
" \"b_width\": cp.Uniform(0.015, 0.025),\n",
" \"b_slope\": cp.Uniform(0.005, 0.020)\n",
" }\n",
" vary_2 = {\n",
" \"Qe_tot\": cp.Uniform(1.8e6, 2.2e6),\n",
" \"Te_bc\": cp.Uniform(80.0, 120.0)\n",
" }\n",
" vary_5 = {\n",
" \"Qe_tot\": cp.Uniform(1.8e6, 2.2e6),\n",
" \"H0\": cp.Uniform(0.0, 0.2),\n",
" \"Hw\": cp.Uniform(0.1, 0.5),\n",
" \"chi\": cp.Uniform(0.8, 1.2),\n",
" \"Te_bc\": cp.Uniform(80.0, 120.0)\n",
" }\n",
" vary_10 = {\n",
" \"Qe_tot\": cp.Uniform(1.8e6, 2.2e6),\n",
" \"H0\": cp.Uniform(0.0, 0.2),\n",
" \"Hw\": cp.Uniform(0.1, 0.5),\n",
" \"chi\": cp.Uniform(0.8, 1.2),\n",
" \"Te_bc\": cp.Uniform(80.0, 120.0),\n",
" \"b_pos\": cp.Uniform(0.95, 0.99),\n",
" \"b_height\": cp.Uniform(5e19, 7e19),\n",
" \"b_sol\": cp.Uniform(1e19, 3e19),\n",
" \"b_width\": cp.Uniform(0.015, 0.025),\n",
" \"b_slope\": cp.Uniform(0.005, 0.020)\n",
" }\n",
" return vary_5"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.335304Z",
"start_time": "2021-06-04T20:25:34.322416Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.471544Z",
"iopub.status.busy": "2025-07-19T22:23:35.471478Z",
"iopub.status.idle": "2025-07-19T22:23:35.474090Z",
"shell.execute_reply": "2025-07-19T22:23:35.473798Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.471539Z"
}
},
"outputs": [],
"source": [
"# define a model to run the fusion code directly from python, expecting a dictionary and returning a dictionary\n",
"def run_fusion_model(input):\n",
" import json\n",
" import fusion\n",
" qois = [\"te\", \"ne\", \"rho\", \"rho_norm\"]\n",
" del input['out_file']\n",
" return {q: v for q,v in zip(qois, [t.tolist() for t in fusion.solve_Te(**input, plots=False, output=False)])}"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.351660Z",
"start_time": "2021-06-04T20:25:34.337059Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.474632Z",
"iopub.status.busy": "2025-07-19T22:23:35.474536Z",
"iopub.status.idle": "2025-07-19T22:23:35.479503Z",
"shell.execute_reply": "2025-07-19T22:23:35.479220Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.474626Z"
}
},
"outputs": [],
"source": [
"# routine to run a SC campaign\n",
"\n",
"def run_sc_case(sc_order=2, local=True, dask=True, batch_size=os.cpu_count(), use_files=True):\n",
" \"\"\"\n",
" Inputs:\n",
" sc_order: order of the sc expansion\n",
" local: if using Dask, whether to use the local option (True)\n",
" dask: whether to use dask (True)\n",
" batch_size: for the non Dask option, number of cases to run in parallel (16)\n",
" Outputs:\n",
" results_df: Pandas dataFrame containing inputs to and output from the model\n",
" results: Results of the sc analysis\n",
" times: Information about the elapsed time for the various phases of the calculation\n",
" sc_order: sc_order \n",
" count: number of sc samples\n",
" \"\"\"\n",
" \n",
" if dask:\n",
" if local:\n",
" print('Running locally')\n",
" import multiprocessing.popen_spawn_posix\n",
" from dask.distributed import Client, LocalCluster\n",
" cluster = LocalCluster(threads_per_worker=1)\n",
" client = Client(cluster) # processes=True, threads_per_worker=1)\n",
" else:\n",
" print('Running using SLURM')\n",
" from dask.distributed import Client\n",
" from dask_jobqueue import SLURMCluster\n",
" cluster = SLURMCluster(\n",
" job_extra=['--qos=p.tok.openmp.2h', '--mail-type=end', '--mail-user=dpc@rzg.mpg.de', '-t 2:00:00'], \n",
" queue='p.tok.openmp', \n",
" cores=8, \n",
" memory='8 GB',\n",
" processes=8)\n",
" cluster.scale(32)\n",
" print(cluster)\n",
" print(cluster.job_script())\n",
" client = Client(cluster)\n",
" print(client)\n",
"\n",
" else:\n",
" import concurrent.futures\n",
"# client = concurrent.futures.ProcessPoolExecutor(max_workers=batch_size)\n",
" client = concurrent.futures.ThreadPoolExecutor(max_workers=batch_size)\n",
"# client = None\n",
" \n",
" times = np.zeros(7)\n",
"\n",
" time_start = time.time()\n",
" time_start_whole = time_start\n",
" # Set up a fresh campaign called \"fusion_sc.\"\n",
" my_campaign = uq.Campaign(name='fusion_sc.') \n",
"\n",
" # Define parameter space\n",
" params = define_params()\n",
"\n",
" # Create an encoder and decoder for sc test app\n",
" if use_files:\n",
" encoder = uq.encoders.GenericEncoder(template_fname='fusion.template',\n",
" delimiter='$',\n",
" target_filename='fusion_in.json')\n",
"\n",
"\n",
" decoder = uq.decoders.SimpleCSV(target_filename=\"output.csv\",\n",
" output_columns=[\"te\", \"ne\", \"rho\", \"rho_norm\"])\n",
"\n",
" execute = uq.actions.ExecuteLocal('python3 %s/fusion_model.py fusion_in.json' % (os.getcwd()))\n",
"\n",
" actions = uq.actions.Actions(uq.actions.CreateRunDirectory('/tmp'), \n",
" uq.actions.Encode(encoder), execute, uq.actions.Decode(decoder))\n",
" else:\n",
" actions = uq.actions.Actions(uq.actions.ExecutePython(run_fusion_model))\n",
"\n",
"\n",
" # Add the app (automatically set as current app)\n",
" my_campaign.add_app(name=\"fusion\", params=params, actions=actions)\n",
"\n",
" time_end = time.time()\n",
" times[1] = time_end-time_start\n",
" print('Time for phase 1 = %.3f' % (times[1]))\n",
"\n",
" time_start = time.time()\n",
" # Associate a sampler with the campaign\n",
" my_campaign.set_sampler(uq.sampling.SCSampler(vary=define_vary(), polynomial_order=sc_order))\n",
" my_campaign.draw_samples()\n",
" print('Number of samples = %s' % my_campaign.get_active_sampler().count)\n",
"\n",
" time_end = time.time()\n",
" times[2] = time_end-time_start\n",
" print('Time for phase 2 = %.3f' % (times[2]))\n",
"\n",
" time_start = time.time()\n",
" # Perform the actions\n",
" my_campaign.execute(pool=client).collate(progress_bar=True)\n",
"\n",
" if dask:\n",
" client.close()\n",
" client.shutdown()\n",
"\n",
" time_end = time.time()\n",
" times[3] = time_end-time_start\n",
" print('Time for phase 3 = %.3f' % (times[3]))\n",
"\n",
" time_start = time.time()\n",
" # Collate the results\n",
" results_df = my_campaign.get_collation_result()\n",
"\n",
" time_end = time.time()\n",
" times[4] = time_end-time_start\n",
" print('Time for phase 4 = %.3f' % (times[4]))\n",
"\n",
" time_start = time.time()\n",
" # Post-processing analysis\n",
" results = my_campaign.analyse(qoi_cols=[\"te\", \"ne\", \"rho\", \"rho_norm\"])\n",
"\n",
" time_end = time.time()\n",
" times[5] = time_end-time_start\n",
" print('Time for phase 5 = %.3f' % (times[5]))\n",
"\n",
" time_start = time.time()\n",
" # Save the results\n",
" pickle.dump(results, open('fusion_results.pickle','bw'))\n",
" time_end = time.time()\n",
" times[6] = time_end-time_start\n",
" print('Time for phase 6 = %.3f' % (times[6]))\n",
"\n",
" times[0] = time_end - time_start_whole\n",
"\n",
" return results_df, results, times, sc_order, my_campaign.get_active_sampler().count"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-04T20:25:34.378499Z",
"start_time": "2021-06-04T20:25:34.354265Z"
},
"code_folding": [
0,
2,
24,
46,
58,
72,
84
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.480057Z",
"iopub.status.busy": "2025-07-19T22:23:35.479950Z",
"iopub.status.idle": "2025-07-19T22:23:35.486738Z",
"shell.execute_reply": "2025-07-19T22:23:35.486487Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.480051Z"
}
},
"outputs": [],
"source": [
"# routines for plotting the results\n",
"\n",
"def plot_Te(results, title=None):\n",
" # plot the calculated Te: mean, with std deviation, 1, 10, 90 and 99%\n",
" plt.figure()\n",
" rho = results.describe('rho', 'mean')\n",
" plt.plot(rho, results.describe('te', 'mean'), 'b-', label='Mean')\n",
" plt.plot(rho, results.describe('te', 'mean')-results.describe('te', 'std'), 'b--', label='+1 std deviation')\n",
" plt.plot(rho, results.describe('te', 'mean')+results.describe('te', 'std'), 'b--')\n",
" plt.fill_between(rho, results.describe('te', 'mean')-results.describe('te', 'std'), results.describe('te', 'mean')+results.describe('te', 'std'), color='b', alpha=0.2)\n",
" try:\n",
" plt.plot(rho, results.describe('te', '10%'), 'b:', label='10 and 90 percentiles')\n",
" plt.plot(rho, results.describe('te', '90%'), 'b:')\n",
" plt.fill_between(rho, results.describe('te', '10%'), results.describe('te', '90%'), color='b', alpha=0.1)\n",
" plt.fill_between(rho, results.describe('te', '1%'), results.describe('te', '99%'), color='b', alpha=0.05)\n",
" except:\n",
" print('Problem with some of the percentiles')\n",
" plt.legend(loc=0)\n",
" plt.xlabel('rho [$m$]')\n",
" plt.ylabel('Te [$eV$]')\n",
" if not title is None: plt.title(title)\n",
" plt.savefig('Te.png')\n",
" plt.savefig('Te.pdf')\n",
"\n",
"def plot_ne(results, title=None):\n",
" # plot the calculated ne: mean, with std deviation, 1, 10, 90 and 99%\n",
" plt.figure()\n",
" rho = results.describe('rho', 'mean')\n",
" plt.plot(rho, results.describe('ne', 'mean'), 'b-', label='Mean')\n",
" plt.plot(rho, results.describe('ne', 'mean')-results.describe('ne', 'std'), 'b--', label='+1 std deviation')\n",
" plt.plot(rho, results.describe('ne', 'mean')+results.describe('ne', 'std'), 'b--')\n",
" plt.fill_between(rho, results.describe('ne', 'mean')-results.describe('ne', 'std'), results.describe('ne', 'mean')+results.describe('ne', 'std'), color='b', alpha=0.2)\n",
" try:\n",
" plt.plot(rho, results.describe('ne', '10%'), 'b:', label='10 and 90 percentiles')\n",
" plt.plot(rho, results.describe('ne', '90%'), 'b:')\n",
" plt.fill_between(rho, results.describe('ne', '10%'), results.describe('ne', '90%'), color='b', alpha=0.1)\n",
" plt.fill_between(rho, results.describe('ne', '1%'), results.describe('ne', '99%'), color='b', alpha=0.05)\n",
" except:\n",
" print('Problem with some of the percentiles')\n",
" plt.legend(loc=0)\n",
" plt.xlabel('rho [$m$]')\n",
" plt.ylabel('ne [$m^{-3}$]')\n",
" if not title is None: plt.title(title)\n",
" plt.savefig('ne.png')\n",
" plt.savefig('ne.pdf')\n",
"\n",
"def plot_sobols_first(results, title=None, field='te'):\n",
" # plot the first Sobol results\n",
" plt.figure()\n",
" rho = results.describe('rho', 'mean')\n",
" for k in results.sobols_first()[field].keys(): plt.plot(rho, results.sobols_first()[field][k], label=k)\n",
" plt.legend(loc=0)\n",
" plt.xlabel('rho [$m$]')\n",
" plt.ylabel('sobols_first')\n",
" if not title is None: plt.title(field + ': ' + title)\n",
" plt.savefig('sobols_first_%s.png' % (field))\n",
" plt.savefig('sobols_first_%s.pdf' % (field))\n",
"\n",
"def plot_sobols_second(results, title=None, field='te'):\n",
" # plot the second Sobol results\n",
" plt.figure()\n",
" rho = results.describe('rho', 'mean')\n",
" for k1 in results.sobols_second()[field].keys():\n",
" for k2 in results.sobols_second()[field][k1].keys():\n",
" plt.plot(rho, results.sobols_second()[field][k1][k2], label=k1+'/'+k2)\n",
" plt.legend(loc=0, ncol=2)\n",
" plt.xlabel('rho [$m$]')\n",
" plt.ylabel('sobols_second')\n",
" if not title is None: plt.title(field + ': ' + title)\n",
" plt.savefig('sobols_second_%s.png' % (field))\n",
" plt.savefig('sobols_second_%s.pdf' % (field))\n",
"\n",
"def plot_sobols_total(results, title=None, field='te'):\n",
" # plot the total Sobol results\n",
" plt.figure()\n",
" rho = results.describe('rho', 'mean')\n",
" for k in results.sobols_total()[field].keys(): plt.plot(rho, results.sobols_total()[field][k], label=k)\n",
" plt.legend(loc=0)\n",
" plt.xlabel('rho [$m$]')\n",
" plt.ylabel('sobols_total')\n",
" if not title is None: plt.title(field + ': ' + title)\n",
" plt.savefig('sobols_total_%s.png' % (field))\n",
" plt.savefig('sobols_total_%s.pdf' % (field))\n",
"\n",
"def plot_distribution(results, results_df, title=None):\n",
" te_dist = results.raw_data['output_distributions']['te']\n",
" rho_norm = results.describe('rho_norm', 'mean')\n",
" for i in [np.maximum(0, int(i-1)) \n",
" for i in np.linspace(0,1,5) * rho_norm.shape]:\n",
" plt.figure()\n",
" pdf_raw_samples = cp.GaussianKDE(results_df.te[i])\n",
" pdf_kde_samples = cp.GaussianKDE(te_dist.samples[i])\n",
" plt.hist(results_df.te[i], density=True, bins=50, label='histogram of raw samples', alpha=0.25)\n",
" if hasattr(te_dist, 'samples'):\n",
" plt.hist(te_dist.samples[i], density=True, bins=50, label='histogram of kde samples', alpha=0.25)\n",
"\n",
" plt.plot(np.linspace(pdf_raw_samples.lower, pdf_raw_samples.upper), pdf_raw_samples.pdf(np.linspace(pdf_raw_samples.lower, pdf_raw_samples.upper)), label='PDF (raw samples)')\n",
" plt.plot(np.linspace(pdf_kde_samples.lower, pdf_kde_samples.upper), pdf_kde_samples.pdf(np.linspace(pdf_kde_samples.lower, pdf_kde_samples.upper)), label='PDF (kde samples)')\n",
"\n",
" plt.legend(loc=0)\n",
" plt.xlabel('Te [$eV$]')\n",
" if title is None:\n",
" plt.title('Distributions for rho_norm = %0.4f' % (rho_norm[i]))\n",
" else:\n",
" plt.title('%s\\nDistributions for rho_norm = %0.4f' % (title, rho_norm[i]))\n",
" plt.savefig('distribution_function_rho_norm=%0.4f.png' % (rho_norm[i]))\n",
" plt.savefig('distribution_function_rho_norm=%0.4f.pdf' % (rho_norm[i]))"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"ExecuteTime": {
"end_time": "2021-07-16T11:18:30.251249Z",
"start_time": "2021-06-04T20:25:34.380415Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:23:35.487185Z",
"iopub.status.busy": "2025-07-19T22:23:35.487072Z",
"iopub.status.idle": "2025-07-19T22:32:36.791508Z",
"shell.execute_reply": "2025-07-19T22:32:36.791242Z",
"shell.execute_reply.started": "2025-07-19T22:23:35.487177Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Volumes/UserData/dpc/GIT/EasyVVUQ/env_3.12/lib/python3.12/site-packages/cerberus/validator.py:618: UserWarning: These types are defined both with a method and in the'types_mapping' property of this validator: {'integer'}\n",
" warn(\n",
"/Volumes/UserData/dpc/GIT/EasyVVUQ/env_3.12/lib/python3.12/site-packages/cerberus/validator.py:618: UserWarning: These types are defined both with a method and in the'types_mapping' property of this validator: {'integer'}\n",
" warn(\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 1 = 0.017\n",
"Number of samples = 32\n",
"Time for phase 2 = 0.035\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 32/32 [00:00<00:00, 74.02it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 0.448\n",
"Time for phase 4 = 0.006\n",
"Time for phase 5 = 0.030\n",
"Time for phase 6 = 0.001\n",
"Time for phase 1 = 0.005\n",
"Number of samples = 243\n",
"Time for phase 2 = 0.115\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 243/243 [00:02<00:00, 88.70it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 3.087\n",
"Time for phase 4 = 0.024\n",
"Time for phase 5 = 0.214\n",
"Time for phase 6 = 0.002\n",
"Time for phase 1 = 0.004\n",
"Number of samples = 1024\n",
"Time for phase 2 = 0.381\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 1024/1024 [00:12<00:00, 80.66it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 13.070\n",
"Time for phase 4 = 0.097\n",
"Time for phase 5 = 1.204\n",
"Time for phase 6 = 0.010\n",
"Time for phase 1 = 0.006\n",
"Number of samples = 3125\n",
"Time for phase 2 = 1.143\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 3125/3125 [00:39<00:00, 78.79it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 40.277\n",
"Time for phase 4 = 0.287\n",
"Time for phase 5 = 6.339\n",
"Time for phase 6 = 0.023\n",
"Time for phase 1 = 0.005\n",
"Number of samples = 7776\n",
"Time for phase 2 = 2.637\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 7776/7776 [01:39<00:00, 77.86it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 100.147\n",
"Time for phase 4 = 0.747\n",
"Time for phase 5 = 28.585\n",
"Time for phase 6 = 0.057\n",
"Time for phase 1 = 0.004\n",
"Number of samples = 16807\n",
"Time for phase 2 = 5.628\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"100%|████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 16807/16807 [03:33<00:00, 78.86it/s]\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Time for phase 3 = 213.790\n",
"Time for phase 4 = 1.743\n",
"Time for phase 5 = 121.002\n",
"Time for phase 6 = 0.131\n"
]
}
],
"source": [
"# Calculate the stochastic collocation expansion for a range of orders\n",
"\n",
"if __name__ == '__main__':\n",
" local = False # if True, use local cores; if False, use SLURM\n",
" dask = False # if True, use DASK; if False, use a fall-back non-DASK option\n",
"\n",
" R = {}\n",
" for sc_order in range(1, 7):\n",
" R[sc_order] = {}\n",
" (R[sc_order]['results_df'], \n",
" R[sc_order]['results'], \n",
" R[sc_order]['times'], \n",
" R[sc_order]['order'], \n",
" R[sc_order]['number_of_samples']) = run_sc_case(sc_order=sc_order, \n",
" local=local, dask=dask, \n",
" batch_size=7, use_files=False)"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"ExecuteTime": {
"end_time": "2021-07-16T11:20:48.322921Z",
"start_time": "2021-07-16T11:18:30.253551Z"
},
"code_folding": [],
"execution": {
"iopub.execute_input": "2025-07-19T22:32:36.791968Z",
"iopub.status.busy": "2025-07-19T22:32:36.791894Z",
"iopub.status.idle": "2025-07-19T22:32:37.078940Z",
"shell.execute_reply": "2025-07-19T22:32:37.078592Z",
"shell.execute_reply.started": "2025-07-19T22:32:36.791960Z"
}
},
"outputs": [],
"source": [
"# save the results\n",
"\n",
"if __name__ == '__main__':\n",
"\n",
" pickle.dump(R, open('collected_results.pickle','bw'))"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"ExecuteTime": {
"end_time": "2021-07-16T11:20:52.573368Z",
"start_time": "2021-07-16T11:20:48.326442Z"
},
"code_folding": [
0
],
"execution": {
"iopub.execute_input": "2025-07-19T22:32:37.079391Z",
"iopub.status.busy": "2025-07-19T22:32:37.079312Z",
"iopub.status.idle": "2025-07-19T22:32:37.089674Z",
"shell.execute_reply": "2025-07-19T22:32:37.089473Z",
"shell.execute_reply.started": "2025-07-19T22:32:37.079384Z"
}
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" Total | \n",
" Phase 1 | \n",
" Phase 2 | \n",
" Phase 3 | \n",
" Phase 4 | \n",
" Phase 5 | \n",
" Phase 6 | \n",
"
\n",
" \n",
" \n",
" \n",
" | 1 | \n",
" 0.536658 | \n",
" 0.017182 | \n",
" 0.034876 | \n",
" 0.448066 | \n",
" 0.006149 | \n",
" 0.029729 | \n",
" 0.000557 | \n",
"
\n",
" \n",
" | 2 | \n",
" 3.446937 | \n",
" 0.004853 | \n",
" 0.114904 | \n",
" 3.087407 | \n",
" 0.024183 | \n",
" 0.213731 | \n",
" 0.001767 | \n",
"
\n",
" \n",
" | 3 | \n",
" 14.766615 | \n",
" 0.004219 | \n",
" 0.381222 | \n",
" 13.070184 | \n",
" 0.096659 | \n",
" 1.204292 | \n",
" 0.009901 | \n",
"
\n",
" \n",
" | 4 | \n",
" 48.073970 | \n",
" 0.005606 | \n",
" 1.143053 | \n",
" 40.276773 | \n",
" 0.286804 | \n",
" 6.338882 | \n",
" 0.022661 | \n",
"
\n",
" \n",
" | 5 | \n",
" 132.179250 | \n",
" 0.004542 | \n",
" 2.637365 | \n",
" 100.147288 | \n",
" 0.747055 | \n",
" 28.585327 | \n",
" 0.057435 | \n",
"
\n",
" \n",
" | 6 | \n",
" 342.298249 | \n",
" 0.004455 | \n",
" 5.628100 | \n",
" 213.790076 | \n",
" 1.743332 | \n",
" 121.001515 | \n",
" 0.130513 | \n",
"
\n",
" \n",
"
\n",
"