Perovskite light-intensity dependant JV fits with SIMsalabim (fake data)
This notebook is a demonstration of how to fit light-intensity dependent JV curves with drift-diffusion models using the SIMsalabim package.
[1]:
# Import necessary libraries
import warnings, os, sys, shutil
# remove warnings from the output
os.environ["PYTHONWARNINGS"] = "ignore"
warnings.filterwarnings(action='ignore', category=FutureWarning)
warnings.filterwarnings(action='ignore', category=UserWarning)
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
from numpy.random import default_rng
import torch, copy, uuid
import pySIMsalabim as sim
from pySIMsalabim.experiments.JV_steady_state import *
import ax, logging
from ax.utils.notebook.plotting import init_notebook_plotting, render
init_notebook_plotting() # for Jupyter notebooks
try:
from optimpv import *
from optimpv.axBOtorch.axUtils import *
except Exception as e:
sys.path.append('../') # add the path to the optimpv module
from optimpv import *
from optimpv.axBOtorch.axUtils import *
[INFO 08-13 17:56:35] ax.utils.notebook.plotting: Injecting Plotly library into cell. Do not overwrite or delete cell.
[INFO 08-13 17:56:35] ax.utils.notebook.plotting: Please see
(https://ax.dev/tutorials/visualizations.html#Fix-for-plots-that-are-not-rendering)
if visualizations are not rendering.
Define the parameters for the simulation
[2]:
params = [] # list of parameters to be optimized
mun = FitParam(name = 'l2.mu_n', value = 6e-4, bounds = [1e-5,1e-3], values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$\mu_n$', unit='m$^2$ V$^{-1}$s$^{-1}$', axis_type = 'log', std = 0,encoding = None,force_log = True)
params.append(mun)
mup = FitParam(name = 'l2.mu_p', value = 4e-4, bounds = [1e-5,1e-3], values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$\mu_p$', unit=r'm$^2$ V$^{-1}$s$^{-1}$', axis_type = 'log', std = 0,encoding = None,force_log = True)
params.append(mup)
bulk_tr = FitParam(name = 'l2.N_t_bulk', value = 1e20, bounds = [1e19,1e21], values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$N_{T}$', unit=r'm$^{-3}$', axis_type = 'log', std = 0,encoding = None,force_log = False)
params.append(bulk_tr)
HTL_int_trap = FitParam(name = 'l1.N_t_int', value = 5e11, bounds = [1e11,1e13], values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$N_{T,int}^{HTL}$', unit='m$^{-2}$', axis_type = 'log', std = 0,encoding = None,force_log = False)
params.append(HTL_int_trap)
ETL_int_trap = FitParam(name = 'l2.N_t_int', value = 4e12, bounds = [1e11,1e13], values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$N_{T,int}^{ETL}$', unit='m$^{-2}$', axis_type = 'log', std = 0,encoding = None,force_log = False)
params.append(ETL_int_trap)
Nions = FitParam(name = 'l2.N_ions', value = 1e22, bounds = [1e20,5e22], type='range', values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$C_{ions}$', unit='m$^{-3}$', axis_type = 'log', std = 0,encoding = None,force_log = False)
params.append(Nions)
R_series = FitParam(name = 'R_series', value = 1e-4, bounds = [1e-5,1e-3], type='range', values = None, start_value = None, log_scale = True, value_type = 'float', fscale = None, rescale = False, stepsize = None, display_name=r'$R_{series}$', unit=r'$\Omega$ m$^2$', axis_type = 'log', std = 0,encoding = None,force_log = False)
params.append(R_series)
# save the original parameters for later
params_orig = copy.deepcopy(params)
Generate some fake data
Here we generate some fake data to fit. The data is generated using the same model as the one used for the fitting, so it is a good test of the fitting procedure. For more information on how to run SIMsalabim from python see the pySIMsalabim package.
[3]:
# Set the session path for the simulation and the input files
session_path = os.path.join(os.path.join(os.path.abspath('../'),'SIMsalabim','SimSS'))
input_path = os.path.join(os.path.join(os.path.join(os.path.abspath('../'),'Data','simsalabim_test_inputs','fakePerovskite')))
simulation_setup_filename = 'simulation_setup_fakePerovskite.txt'
simulation_setup = os.path.join(session_path, simulation_setup_filename)
optical_files = ['nk_glass.txt','nk_ITO.txt','nk_PTAA.txt','nk_FACsPbIBr.txt','nk_C60_1.txt','nk_Au.txt']
# path to the layer files defined in the simulation_setup file
l1 = 'PTAA.txt'
l2 = 'fakePerovskite.txt'
l3 = 'C60.txt'
l1 = os.path.join(input_path, l1)
l2 = os.path.join(input_path, l2)
l3 = os.path.join(input_path, l3)
# copy this files to session_path
force_copy = True
if not os.path.exists(session_path):
os.makedirs(session_path)
for file in [l1,l2,l3,simulation_setup_filename]:
file = os.path.join(input_path, os.path.basename(file))
if force_copy or not os.path.exists(os.path.join(session_path, os.path.basename(file))):
shutil.copyfile(file, os.path.join(session_path, os.path.basename(file)))
else:
print('File already exists: ',file)
# copy the optical files to the session path
for file in optical_files:
file = os.path.join(input_path, os.path.basename(file))
if force_copy or not os.path.exists(os.path.join(session_path, os.path.basename(file))):
shutil.copyfile(file, os.path.join(session_path, os.path.basename(file)))
# Show the device structure
fig = sim.plot_band_diagram(simulation_setup, session_path)
# reset simss
# Set the JV parameters
Gfracs = [0.1,0.5,1.0] # Fractions of the generation rate to simulate (None if you want only one light intensity as define in the simulation_setup file)
UUID = str(uuid.uuid4()) # random UUID to avoid overwriting files
cmd_pars = [] # see pySIMsalabim documentation for the command line parameters
# Add the parameters to the command line arguments
for param in params:
if param.name == 'l2.N_ions':
cmd_pars.append({'par':'l2.N_cation', 'val':str(param.value)})
cmd_pars.append({'par':'l2.N_anion', 'val':str(param.value)})
else:
cmd_pars.append({'par':param.name, 'val':str(param.value)})
# Add the layer files to the command line arguments
cmd_pars.append({'par':'l1', 'val':'PTAA.txt'})
cmd_pars.append({'par':'l2', 'val':'fakePerovskite.txt'})
cmd_pars.append({'par':'l3', 'val':'C60.txt'})
# Run the JV simulation
ret, mess = run_SS_JV(simulation_setup, session_path, JV_file_name = 'JV.dat', G_fracs = Gfracs, parallel = True, max_jobs = 3, UUID=UUID, cmd_pars=cmd_pars)
# save data for fitting
X,y = [],[]
X_orig,y_orig = [],[]
if Gfracs is None:
data = pd.read_csv(os.path.join(session_path, 'JV_'+UUID+'.dat'), sep=r'\s+') # Load the data
Vext = np.asarray(data['Vext'].values)
Jext = np.asarray(data['Jext'].values)
G = np.ones_like(Vext)
rng = default_rng()#
noise = rng.standard_normal(Jext.shape) * 0.01 * Jext
Jext = Jext + noise
X = Vext
y = Jext
plt.figure()
plt.plot(X,y)
plt.show()
else:
for Gfrac in Gfracs:
data = pd.read_csv(os.path.join(session_path, 'JV_Gfrac_'+str(Gfrac)+'_'+UUID+'.dat'), sep=r'\s+') # Load the data
Vext = np.asarray(data['Vext'].values)
Jext = np.asarray(data['Jext'].values)
G = np.ones_like(Vext)*Gfrac
rng = default_rng()#
noise = rng.standard_normal(Jext.shape) * 0.005 * Jext
if len(X) == 0:
X = np.vstack((Vext,G)).T
y = Jext + noise
y_orig = Jext
else:
X = np.vstack((X,np.vstack((Vext,G)).T))
y = np.hstack((y,Jext+ noise))
y_orig = np.hstack((y_orig,Jext))
# remove all the current where Jext is higher than a given value
X = X[y<200]
X_orig = copy.deepcopy(X)
y_orig = y_orig[y<200]
y = y[y<200]
plt.figure()
for Gfrac in Gfracs:
plt.plot(X[X[:,1]==Gfrac,0],y[X[:,1]==Gfrac],label='Gfrac = '+str(Gfrac))
plt.xlabel('Voltage [V]')
plt.ylabel('Current density [A/m$^2$]')
plt.legend()
plt.show()
Run the optimization
[4]:
# Define the Agent and the target metric/loss function
from optimpv.DDfits.JVAgent import JVAgent
metric = 'nrmse' # can be 'nrmse', 'mse', 'mae'
loss = 'linear' # can be 'linear', 'huber', 'soft_l1'
jv = JVAgent(params, X, y, session_path, simulation_setup, parallel = True, max_jobs = 3, metric = metric, loss = loss)
# Calulate the target metric for the original parameters
best_fit_possible = loss_function(calc_metric(y,y_orig, metric_name = metric),loss)
print('Best fit: ',best_fit_possible)
Best fit: 0.0014083591744227725
[5]:
from optimpv.axBOtorch.axBOtorchOptimizer import axBOtorchOptimizer
from botorch.acquisition.logei import qLogNoisyExpectedImprovement
from ax.adapter.transforms.standardize_y import StandardizeY
from ax.adapter.transforms.unit_x import UnitX
from ax.adapter.transforms.remove_fixed import RemoveFixed
from ax.adapter.transforms.log import Log
from ax.generators.torch.botorch_modular.utils import ModelConfig
from ax.generators.torch.botorch_modular.surrogate import SurrogateSpec
from gpytorch.kernels import MaternKernel
from gpytorch.kernels import ScaleKernel
from botorch.models import SingleTaskGP
model_gen_kwargs_list = None
parameter_constraints = None
model_kwargs_list = [{},{"torch_device":torch.device("cuda" if torch.cuda.is_available() else "cpu"),'botorch_acqf_class':qLogNoisyExpectedImprovement,'transforms':[RemoveFixed, Log,UnitX, StandardizeY],'surrogate_spec':SurrogateSpec(model_configs=[ModelConfig(botorch_model_class=SingleTaskGP,covar_module_class=ScaleKernel, covar_module_options={'base_kernel':MaternKernel(nu=2.5, ard_num_dims=len(params))})])}]
# Define the optimizer
optimizer = axBOtorchOptimizer(params = params, agents = jv, models = ['SOBOL','BOTORCH_MODULAR'],n_batches = [1,45], batch_size = [10,2], ax_client = None, max_parallelism = -1, model_kwargs_list = model_kwargs_list, model_gen_kwargs_list = None, name = 'ax_opti',parallel=True)
[6]:
# optimizer.optimize() # run the optimization with ax
optimizer.optimize_turbo() # run the optimization with turbo
[INFO 08-13 17:56:37] optimpv.axBOtorchOptimizer: Starting optimization with 46 batches and a total of 100 trials
[INFO 08-13 17:56:37] optimpv.axBOtorchOptimizer: Starting Sobol batch 1 with 10 trials
[INFO 08-13 17:56:38] optimpv.axBOtorchOptimizer: Finished Sobol
[INFO 08-13 17:56:40] optimpv.axBOtorchOptimizer: Finished Turbo batch 2 with 2 trials with current best value: 6.02e-02, TR length: 8.00e-01
[INFO 08-13 17:56:41] optimpv.axBOtorchOptimizer: Finished Turbo batch 3 with 2 trials with current best value: 5.34e-02, TR length: 8.00e-01
[INFO 08-13 17:56:43] optimpv.axBOtorchOptimizer: Finished Turbo batch 4 with 2 trials with current best value: 4.92e-02, TR length: 8.00e-01
[INFO 08-13 17:56:44] optimpv.axBOtorchOptimizer: Finished Turbo batch 5 with 2 trials with current best value: 4.92e-02, TR length: 8.00e-01
[INFO 08-13 17:56:45] optimpv.axBOtorchOptimizer: Finished Turbo batch 6 with 2 trials with current best value: 4.24e-02, TR length: 8.00e-01
[INFO 08-13 17:56:46] optimpv.axBOtorchOptimizer: Finished Turbo batch 7 with 2 trials with current best value: 4.24e-02, TR length: 8.00e-01
[INFO 08-13 17:56:47] optimpv.axBOtorchOptimizer: Finished Turbo batch 8 with 2 trials with current best value: 2.89e-02, TR length: 8.00e-01
[INFO 08-13 17:56:48] optimpv.axBOtorchOptimizer: Finished Turbo batch 9 with 2 trials with current best value: 2.50e-02, TR length: 8.00e-01
[INFO 08-13 17:56:49] optimpv.axBOtorchOptimizer: Finished Turbo batch 10 with 2 trials with current best value: 2.26e-02, TR length: 1.60e+00
[INFO 08-13 17:56:50] optimpv.axBOtorchOptimizer: Finished Turbo batch 11 with 2 trials with current best value: 2.26e-02, TR length: 1.60e+00
[INFO 08-13 17:56:51] optimpv.axBOtorchOptimizer: Finished Turbo batch 12 with 2 trials with current best value: 2.26e-02, TR length: 1.60e+00
[INFO 08-13 17:56:53] optimpv.axBOtorchOptimizer: Finished Turbo batch 13 with 2 trials with current best value: 1.79e-02, TR length: 1.60e+00
[INFO 08-13 17:56:54] optimpv.axBOtorchOptimizer: Finished Turbo batch 14 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:56:55] optimpv.axBOtorchOptimizer: Finished Turbo batch 15 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:56:56] optimpv.axBOtorchOptimizer: Finished Turbo batch 16 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:56:57] optimpv.axBOtorchOptimizer: Finished Turbo batch 17 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:56:59] optimpv.axBOtorchOptimizer: Finished Turbo batch 18 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:00] optimpv.axBOtorchOptimizer: Finished Turbo batch 19 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:01] optimpv.axBOtorchOptimizer: Finished Turbo batch 20 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:02] optimpv.axBOtorchOptimizer: Finished Turbo batch 21 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:03] optimpv.axBOtorchOptimizer: Finished Turbo batch 22 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:04] optimpv.axBOtorchOptimizer: Finished Turbo batch 23 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:06] optimpv.axBOtorchOptimizer: Finished Turbo batch 24 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:07] optimpv.axBOtorchOptimizer: Finished Turbo batch 25 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:08] optimpv.axBOtorchOptimizer: Finished Turbo batch 26 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:09] optimpv.axBOtorchOptimizer: Finished Turbo batch 27 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:10] optimpv.axBOtorchOptimizer: Finished Turbo batch 28 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:11] optimpv.axBOtorchOptimizer: Finished Turbo batch 29 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:12] optimpv.axBOtorchOptimizer: Finished Turbo batch 30 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:14] optimpv.axBOtorchOptimizer: Finished Turbo batch 31 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:15] optimpv.axBOtorchOptimizer: Finished Turbo batch 32 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:16] optimpv.axBOtorchOptimizer: Finished Turbo batch 33 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:17] optimpv.axBOtorchOptimizer: Finished Turbo batch 34 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:18] optimpv.axBOtorchOptimizer: Finished Turbo batch 35 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:19] optimpv.axBOtorchOptimizer: Finished Turbo batch 36 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:21] optimpv.axBOtorchOptimizer: Finished Turbo batch 37 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:22] optimpv.axBOtorchOptimizer: Finished Turbo batch 38 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:23] optimpv.axBOtorchOptimizer: Finished Turbo batch 39 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:24] optimpv.axBOtorchOptimizer: Finished Turbo batch 40 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:25] optimpv.axBOtorchOptimizer: Finished Turbo batch 41 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:27] optimpv.axBOtorchOptimizer: Finished Turbo batch 42 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:28] optimpv.axBOtorchOptimizer: Finished Turbo batch 43 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:29] optimpv.axBOtorchOptimizer: Finished Turbo batch 44 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:30] optimpv.axBOtorchOptimizer: Finished Turbo batch 45 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:31] optimpv.axBOtorchOptimizer: Finished Turbo batch 46 with 2 trials with current best value: 3.63e-03, TR length: 1.60e+00
[INFO 08-13 17:57:33] optimpv.axBOtorchOptimizer: Finished Turbo batch 47 with 2 trials with current best value: 3.47e-03, TR length: 1.60e+00
[INFO 08-13 17:57:33] ax.api.client: Trial 0 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 1 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 2 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 3 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 4 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 5 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 6 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 7 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 8 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 9 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 10 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 11 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 12 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 13 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 14 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 15 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 16 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 17 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 18 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 19 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 20 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 21 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 22 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 23 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 24 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 25 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 26 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 27 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 28 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 29 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 30 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 31 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 32 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 33 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 34 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 35 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 36 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 37 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 38 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 39 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 40 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 41 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 42 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 43 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 44 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 45 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 46 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 47 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 48 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 49 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 50 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 51 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 52 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 53 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 54 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 55 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 56 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 57 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 58 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 59 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 60 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 61 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 62 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 63 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 64 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 65 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 66 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 67 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 68 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 69 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 70 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 71 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 72 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 73 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 74 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 75 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 76 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 77 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 78 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 79 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 80 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 81 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 82 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 83 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 84 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 85 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 86 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 87 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 88 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 89 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 90 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 91 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 92 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 93 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 94 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 95 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 96 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 97 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 98 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 99 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 100 marked COMPLETED.
[INFO 08-13 17:57:33] ax.api.client: Trial 101 marked COMPLETED.
[INFO 08-13 17:57:33] optimpv.axBOtorchOptimizer: Turbo is terminated as the max number (100) of trials is reached
[7]:
# get the best parameters and update the params list in the optimizer and the agent
ax_client = optimizer.ax_client # get the ax client
optimizer.update_params_with_best_balance() # update the params list in the optimizer with the best parameters
jv.params = optimizer.params # update the params list in the agent with the best parameters
# print the best parameters
print('Best parameters:')
for p,po in zip(optimizer.params, params_orig):
print(p.name, 'fitted value:', p.value, 'original value:', po.value)
print('\nSimSS command line:')
print(jv.get_SIMsalabim_clean_cmd(jv.params)) # print the simss command line with the best parameters
Best parameters:
l2.mu_n fitted value: 8.234637143291717e-05 original value: 0.0006
l2.mu_p fitted value: 0.0006940768150397674 original value: 0.0004
l2.N_t_bulk fitted value: 6.184054359793663e+19 original value: 1e+20
l1.N_t_int fitted value: 1039488674349.9207 original value: 500000000000.0
l2.N_t_int fitted value: 4027796411514.282 original value: 4000000000000.0
l2.N_ions fitted value: 4.521237521869039e+22 original value: 1e+22
R_series fitted value: 0.00014029361686480063 original value: 0.0001
SimSS command line:
./simss -l2.mu_n 8.234637143291717e-05 -l2.mu_p 0.0006940768150397674 -l2.N_t_bulk 6.184054359793663e+19 -l1.N_t_int 1039488674349.9207 -l2.N_t_int 4027796411514.282 -l2.N_anion 4.521237521869039e+22 -l2.N_cation 4.521237521869039e+22 -R_series 0.00014029361686480063
[8]:
# Plot optimization results
data = ax_client.summarize()
all_metrics = optimizer.all_metrics
plt.figure()
plt.plot(np.minimum.accumulate(data[all_metrics]), label="Best value seen so far")
plt.axhline(y=best_fit_possible, color='red', linestyle='--', label="Best fit possible")
plt.yscale("log")
plt.xlabel("Iteration")
plt.ylabel("Target metric: " + metric + " with " + loss + " loss")
plt.legend()
plt.title("Best value seen so far")
print("Best value seen so far is ", min(data[all_metrics]), "at iteration ", int(data[all_metrics].idxmin()))
print("Best value possible is ", best_fit_possible)
plt.show()
Best value seen so far is JV_JV_nrmse_linear at iteration 101
Best value possible is 0.0014083591744227725
[10]:
# Plot the density of the exploration of the parameters
# this gives a nice visualization of where the optimizer focused its exploration and may show some correlation between the parameters
plot_dens = True
if plot_dens:
from optimpv.posterior.posterior import *
params_orig_dict, best_parameters = {}, {}
for p in params_orig:
params_orig_dict[p.name] = p.value
for p in optimizer.params:
best_parameters[p.name] = p.value
fig_dens, ax_dens = plot_density_exploration(params, optimizer = optimizer, best_parameters = best_parameters, params_orig = params_orig_dict, optimizer_type = 'ax')
[11]:
# rerun the simulation with the best parameters
yfit = jv.run(parameters={}) # run the simulation with the best parameters
viridis = plt.get_cmap('viridis', len(Gfracs))
plt.figure(figsize=(10,10))
linewidth = 2
for idx, Gfrac in enumerate(Gfracs[::-1]):
plt.plot(X[X[:,1]==Gfrac,0],y[X[:,1]==Gfrac],label='Gfrac = '+str(Gfrac),color=viridis(idx),alpha=0.5,linewidth=linewidth)
plt.plot(X[X[:,1]==Gfrac,0],yfit[X[:,1]==Gfrac],label='Gfrac = '+str(Gfrac)+' fit',linestyle='--',color=viridis(idx),linewidth=linewidth)
plt.xlabel('Voltage [V]')
plt.ylabel('Current density [A m$^{-2}$]')
plt.legend()
plt.show()
[12]:
# Clean up the output files (comment out if you want to keep the output files)
sim.clean_all_output(session_path)
sim.delete_folders('tmp',session_path)
# uncomment the following lines to delete specific files
sim.clean_up_output('PTAA',session_path)
sim.clean_up_output('fakePerovskite',session_path)
sim.clean_up_output('C60',session_path)
sim.clean_up_output('simulation_setup_fakePerovskite',session_path)
for file in optical_files:
sim.clean_up_output(file,session_path)