<jupyter_start><jupyter_text>EXC RB simulation<NAME>, <NAME><jupyter_code># ~~~
# This file is part of the paper:
#
#  "A NON-CONFORMING DUAL APPROACH FOR ADAPTIVE TRUST-REGION REDUCED BASIS
#           APPROXIMATION OF PDE-CONSTRAINED OPTIMIZATION"
#
#   https://github.com/TiKeil/NCD-corrected-TR-RB-approach-for-pde-opt
#
# Copyright 2019-2020 all developers. All rights reserved.
# License: Licensed as BSD 2-Clause License (http://opensource.org/licenses/BSD-2-Clause)
# ~~~<jupyter_output><empty_output><jupyter_text>Preparations  details<jupyter_code>import sys
path = '../../'
sys.path.append(path)
import numpy as np

from matplotlib import pyplot as plt

from pymor.basic import *
set_log_levels({'pymor': 'WARN'})
from pymor.core.logger import set_log_levels, getLogger
set_log_levels({'pymor': 'ERROR',
                'distributed_adaptive_discretizations': 'DEBUG',
                'notebook': 'INFO'})
logger = getLogger('notebook.notebook')
import matplotlib as mpl
mpl.rcParams['figure.figsize'] = (12.0, 8.0)
mpl.rcParams['font.size'] = 12
mpl.rcParams['savefig.dpi'] = 300

data_path = '../../../EXC_data'
# domain of interest
bounding_box = [[0,0],[2,1]]
domain_of_interest = BitmapFunction('{}/Domain_of_interest.png'.format(data_path), range=[1,0], bounding_box=bounding_box)<jupyter_output><empty_output><jupyter_text>problem definition<jupyter_code>from pdeopt.problems import EXC_problem, set_input_dict
from pdeopt.discretizer import discretize_quadratic_pdeopt_stationary_cg

parametric_quantities = {'walls': [1,4,9], 'windows': [], 'doors': [], 'heaters': [1,3,5,6,7,8,9]}
inactive_quantities = {'removed_walls': [], 'open_windows': [], 'open_doors': [1,2,3,4,5,6,7,10], 'active_heaters': []}
summed_quantities = {'walls': [[1,2,3,8],[4,5,6,7]], 'windows': [], 'doors': [], 'heaters': [[1,2],[3,4],[9,10,11,12]]}

coefficient_expressions = None

parameters_in_q = True
input_dict = set_input_dict(parametric_quantities, inactive_quantities, coefficient_expressions, summed_quantities, parameters_in_q,
                            ac=0.5, owc=[0.025,0.1], iwc= [0.025,0.1], idc=[0.005], wc=[0.0005], ht=[0,100],
                                    owc_c=0.001,  iwc_c= 0.025,     idc_c=0.01,  wc_c=0.025,  ht_c=80)


parameter_scaling = False
u_out = 5

problem, parameter_scales = EXC_problem(input_dict, summed_quantities, outside_temperature=u_out, #, q_inverse=0.0001
                                        data_path = data_path,parameters_in_q=parameters_in_q, 
                                        parameter_scaling=parameter_scaling,
                                        coefficient_expressions=coefficient_expressions)

u_d = 18  
mu_d = None 
print('desired_parameter: ', mu_d)

sigma_d = 100
weights = {'walls': [0.5,0.25,0.05], 'doors': 1, 'heaters': [0.002,0.002,0.001,0.001,0.001,0.001,0.004], 'windows': 1, 'state': sigma_d}

diameter = np.sqrt(2)/200.
opt_fom, data, mu_bar = discretize_quadratic_pdeopt_stationary_cg(problem, diameter, weights, parameter_scales, 
                                                          domain_of_interest, desired_temperature=u_d, 
                                                          mu_for_u_d=mu_d, mu_for_tikhonov=mu_d,
                                                          parameters_in_q=parameters_in_q, product='fixed_energy',
                                                          use_corrected_gradient= True)
print('information on the grid:')
print(data['grid'])

N = 16
validation_set_size = 100
tau_global_RB_J = 5e-4
tau_global_RB_DJ = 1e-1<jupyter_output>information on the grid:
Rect-Grid on domain [0,2] x [0,1]
x0-intervals: 400, x1-intervals: 200
faces: 80000, edges: 160600, vertices: 80601<jupyter_text>Classical Model Order Reduction  BFGS-Greedy for non-corrected functional We now construct a simple RB basis for primal, dual and all sensitivities. For this, we start with an empty basis<jupyter_code>params = []
from pdeopt.model import build_initial_basis
RBbasis, dual_RBbasis = build_initial_basis(opt_fom, params, build_sensitivities=False)

from pdeopt.reductor import QuadraticPdeoptStationaryCoerciveReductor

from pymor.parameters.functionals import MinThetaParameterFunctional

ce = MinThetaParameterFunctional(opt_fom.primal_model.operator.coefficients, mu_bar)

opt_fom = opt_fom.with_(use_corrected_functional=False)
opt_fom = opt_fom.with_(use_corrected_gradient=False)
opt_fom = opt_fom.with_(adjoint_approach=False)

pdeopt_reductor = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                                RBbasis, dual_RBbasis, 
                                                                opt_product=opt_fom.opt_product,
                                                                coercivity_estimator=ce, 
                                                                prepare_for_gradient_estimate=True,
                                                                mu_bar=mu_bar)<jupyter_output>Starting with two bases. Primal and dual have length 0 and 0<jupyter_text>We start a greedy for the whole domain<jupyter_code>set_log_levels({'pymor': 'WARN'})  # <-- set this to 'INFO' if you want to have further details 
from pdeopt.greedy import pdeopt_adaptive_greedy

result_J, result_DJ = pdeopt_adaptive_greedy(opt_fom, pdeopt_reductor, opt_fom.parameter_space, validation_mus=-100, 
                                max_extensions=N, J_atol=tau_global_RB_J, DJ_atol=tau_global_RB_DJ)
opt_rom = result_DJ['rom']
tictoc = result_DJ['time'] + result_J['time']
print('Greedy took {}'.format(tictoc))
picked_mus = result_J['max_err_mus'][:-1]
picked_mus.extend(result_DJ['max_err_mus'])
print('training set sizes: ', result_J['training_set_sizes'], result_DJ['training_set_sizes'])
print('Before the {}th extension of J goal, the max error was {}'.format(result_J['extensions'],result_J['max_errs'][-2]))
print('Before the {}th extension of DJ goal, the max error was {}'.format(result_DJ['extensions'],result_DJ['max_errs'][-1]))
params = picked_mus
RBbasis, dual_RBbasis, RBPrimalSens, RBDualSens = build_initial_basis(opt_fom, params, build_sensitivities=True)<jupyter_output><empty_output><jupyter_text>Model for 1a<jupyter_code>ce = MinThetaParameterFunctional(opt_fom.primal_model.operator.coefficients, mu_bar)

opt_fom = opt_fom.with_(use_corrected_functional=False)
opt_fom = opt_fom.with_(use_corrected_gradient=False)
opt_fom = opt_fom.with_(adjoint_approach=False)

pdeopt_reductor = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                                RBbasis.copy(), dual_RBbasis.copy(), 
                                                                opt_product=opt_fom.opt_product,
                                                                coercivity_estimator=ce, 
                                                                prepare_for_gradient_estimate=True,
                                                                mu_bar=mu_bar)

opt_rom_1a = pdeopt_reductor.reduce()<jupyter_output>Starting with two bases. Primal and dual have length 15 and 15
GRAD J ESTIMATOR: non corrected estimator<jupyter_text>Model for 2a<jupyter_code>opt_fom = opt_fom.with_(use_corrected_functional=True)
opt_fom = opt_fom.with_(use_corrected_gradient=False)
opt_fom = opt_fom.with_(adjoint_approach=False)

pdeopt_reductor_2a = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                            RBbasis.copy(), dual_RBbasis.copy(), 
                                                            opt_product=opt_fom.opt_product,
                                                            coercivity_estimator=ce, mu_bar=mu_bar,
                                                            prepare_for_gradient_estimate=True,
                                                            true_lagrange=True)

opt_rom_2a = pdeopt_reductor_2a.reduce()<jupyter_output>Starting with two bases. Primal and dual have length 15 and 15
GRAD J ESTIMATOR: non corrected estimator<jupyter_text>Model for 3a<jupyter_code>opt_fom = opt_fom.with_(use_corrected_functional=True)
opt_fom = opt_fom.with_(use_corrected_gradient=False)
opt_fom = opt_fom.with_(adjoint_approach=True)

pdeopt_reductor_3a = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                            RBbasis.copy(), dual_RBbasis.copy(), 
                                                            opt_product=opt_fom.opt_product,
                                                            coercivity_estimator=ce, mu_bar=mu_bar,
                                                            adjoint_estimate=True,
                                                            prepare_for_gradient_estimate=True,
                                                            true_lagrange=True)

opt_rom_3a = pdeopt_reductor_3a.reduce()<jupyter_output>Starting with two bases. Primal and dual have length 15 and 15
GRAD J ESTIMATOR: adjoint approach - adjoint estimate (no sensitivities)<jupyter_text>Model for 4a<jupyter_code>opt_fom = opt_fom.with_(use_corrected_functional=True)
opt_fom = opt_fom.with_(use_corrected_gradient=False)
opt_fom = opt_fom.with_(adjoint_approach=True)

pdeopt_reductor_4a = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                            RBbasis.copy(), dual_RBbasis.copy(), 
                                                            opt_product=opt_fom.opt_product,
                                                            coercivity_estimator=ce, mu_bar=mu_bar,
                                                            prepare_for_gradient_estimate=True,
                                                            prepare_for_sensitivity_estimate=True,
                                                            true_lagrange=True)

opt_rom_4a = pdeopt_reductor_4a.reduce()<jupyter_output>Starting with two bases. Primal and dual have length 15 and 15
building a SINGLE sensitivtiy model for any direction...
GRAD J ESTIMATOR: adjoint approach - sensitivity estimate (d_mu) - Also valid for unique<jupyter_text>Model for 5a<jupyter_code>opt_fom = opt_fom.with_(use_corrected_functional=True)
opt_fom = opt_fom.with_(use_corrected_gradient=True)
opt_fom = opt_fom.with_(adjoint_approach=False)

pdeopt_reductor_5a = QuadraticPdeoptStationaryCoerciveReductor(opt_fom, 
                                                            RBbasis.copy(), dual_RBbasis.copy(), 
                                                            RBPrimalSens.copy(), RBDualSens.copy(),
                                                            opt_product=opt_fom.opt_product,
                                                            coercivity_estimator=ce, mu_bar=mu_bar,
                                                            prepare_for_gradient_estimate=True,
                                                            prepare_for_sensitivity_estimate=True,
                                                            true_lagrange=True)

opt_rom_5a = pdeopt_reductor_5a.reduce()<jupyter_output>Starting with separated sensitivity bases. Primal and dual have length 15 and 15
building MULTIPLE sensitivity models for 10 directions...
GRAD J ESTIMATOR: standard (d_mu-tilde sensitivity correction optional)<jupyter_text>compute true errors and estimators<jupyter_code>validation_set = opt_fom.parameter_space.sample_randomly(validation_set_size, seed=0)
from pdeopt.tools import compute_all_errors_and_estimators_for_all_ROMS

J_errors_1a, DJ_errors_1a, rel_J_errors_1a, rel_DJ_errors_1a, J_estimators_1a, DJ_estimators_1a, effectivities_J_1a, effectivities_DJ_1a, \
J_errors_2a, DJ_errors_2a, rel_J_errors_2a, rel_DJ_errors_2a, J_estimators_2a, DJ_estimators_2a, effectivities_J_2a, effectivities_DJ_2a, \
J_errors_3a, DJ_errors_3a, rel_J_errors_3a, rel_DJ_errors_3a, J_estimators_3a, DJ_estimators_3a, effectivities_J_3a, effectivities_DJ_3a, \
J_errors_4a, DJ_errors_4a, rel_J_errors_4a, rel_DJ_errors_4a, J_estimators_4a, DJ_estimators_4a, effectivities_J_4a, effectivities_DJ_4a, \
J_errors_5a, DJ_errors_5a, rel_J_errors_5a, rel_DJ_errors_5a, J_estimators_5a, DJ_estimators_5a, effectivities_J_5a, effectivities_DJ_5a, \
                                    J, DJ, \
u_mu_errors_4a, rel_u_mu_errors_4a, u_mu_estimators_4a, effectivities_u_mu_4a, \
u_mu_errors_5a, rel_u_mu_errors_5a, u_mu_estimators_5a, effectivities_u_mu_5a, \
p_mu_errors_4a, rel_p_mu_errors_4a, p_mu_estimators_4a, effectivities_p_mu_4a, \
p_mu_errors_5a, rel_p_mu_errors_5a, p_mu_estimators_5a, effectivities_p_mu_5a \
        = compute_all_errors_and_estimators_for_all_ROMS(
            validation_set, opt_fom, opt_rom_1a, opt_rom_2a, opt_rom_3a, opt_rom_4a, opt_rom_5a, 
            pdeopt_reductor_4a, pdeopt_reductor_5a)
#J
max_J_error_1a = max(J_errors_1a)
min_J_error_1a = min(J_errors_1a)

max_J_error_2a = max(J_errors_2a)
min_J_error_2a = min(J_errors_2a)

max_J_error_3a = max(J_errors_3a)
min_J_error_3a = min(J_errors_3a)

max_J_error_4a = max(J_errors_4a)
min_J_error_4a = min(J_errors_4a)

max_J_error_5a = max(J_errors_5a)
min_J_error_5a = min(J_errors_5a)

#DJ
max_DJ_error_1a = max(DJ_errors_1a)
min_DJ_error_1a = min(DJ_errors_1a)

max_DJ_error_2a = max(DJ_errors_2a)
min_DJ_error_2a = min(DJ_errors_2a)

max_DJ_error_3a = max(DJ_errors_3a)
min_DJ_error_3a = min(DJ_errors_3a)

max_DJ_error_4a = max(DJ_errors_4a)
min_DJ_error_4a = min(DJ_errors_4a)

max_DJ_error_5a = max(DJ_errors_5a)
min_DJ_error_5a = min(DJ_errors_5a)

#J estimator
max_J_estimators_1a = max(J_estimators_1a)
min_J_estimators_1a = min(J_estimators_1a)

max_J_estimators_2a = max(J_estimators_2a)
min_J_estimators_2a = min(J_estimators_2a)

max_J_estimators_3a = max(J_estimators_3a)
min_J_estimators_3a = min(J_estimators_3a)

max_J_estimators_4a = max(J_estimators_4a)
min_J_estimators_4a = min(J_estimators_4a)

max_J_estimators_5a = max(J_estimators_5a)
min_J_estimators_5a = min(J_estimators_5a)

#DJ estimator
max_DJ_estimators_1a = max(DJ_estimators_1a)
min_DJ_estimators_1a = min(DJ_estimators_1a)

max_DJ_estimators_2a = max(DJ_estimators_2a)
min_DJ_estimators_2a = min(DJ_estimators_2a)

max_DJ_estimators_3a = max(DJ_estimators_3a)
min_DJ_estimators_3a = min(DJ_estimators_3a)

max_DJ_estimators_4a = max(DJ_estimators_4a)
min_DJ_estimators_4a = min(DJ_estimators_4a)

max_DJ_estimators_5a = max(DJ_estimators_5a)
min_DJ_estimators_5a = min(DJ_estimators_5a)


median_effectivities_J_1a = np.sum(effectivities_J_1a)/len(effectivities_J_1a)
median_effectivities_J_2a = np.sum(effectivities_J_2a)/len(effectivities_J_1a)
median_effectivities_J_3a = np.sum(effectivities_J_3a)/len(effectivities_J_1a)
median_effectivities_J_4a = np.sum(effectivities_J_4a)/len(effectivities_J_1a)
median_effectivities_J_5a = np.sum(effectivities_J_5a)/len(effectivities_J_1a)

median_effectivities_DJ_1a = np.sum(effectivities_DJ_1a)/len(effectivities_J_1a)
median_effectivities_DJ_2a = np.sum(effectivities_DJ_2a)/len(effectivities_J_1a)
median_effectivities_DJ_3a = np.sum(effectivities_DJ_3a)/len(effectivities_J_1a)
median_effectivities_DJ_4a = np.sum(effectivities_DJ_4a)/len(effectivities_J_1a)
median_effectivities_DJ_5a = np.sum(effectivities_DJ_5a)/len(effectivities_J_1a)

median_errors_J_1a = np.sum(J_errors_1a)/len(effectivities_J_1a)
median_errors_J_2a = np.sum(J_errors_2a)/len(effectivities_J_1a)
median_errors_J_3a = np.sum(J_errors_3a)/len(effectivities_J_1a)
median_errors_J_4a = np.sum(J_errors_4a)/len(effectivities_J_1a)
median_errors_J_5a = np.sum(J_errors_5a)/len(effectivities_J_1a)

median_estimators_J_1a = np.sum(J_estimators_1a)/len(effectivities_J_1a)
median_estimators_J_2a = np.sum(J_estimators_2a)/len(effectivities_J_1a)
median_estimators_J_3a = np.sum(J_estimators_3a)/len(effectivities_J_1a)
median_estimators_J_4a = np.sum(J_estimators_4a)/len(effectivities_J_1a)
median_estimators_J_5a = np.sum(J_estimators_5a)/len(effectivities_J_1a)

median_errors_DJ_1a = np.sum(DJ_errors_1a)/len(effectivities_J_1a)
median_errors_DJ_2a = np.sum(DJ_errors_2a)/len(effectivities_J_1a)
median_errors_DJ_3a = np.sum(DJ_errors_3a)/len(effectivities_J_1a)
median_errors_DJ_4a = np.sum(DJ_errors_4a)/len(effectivities_J_1a)
median_errors_DJ_5a = np.sum(DJ_errors_5a)/len(effectivities_J_1a)

median_estimators_DJ_1a = np.sum(DJ_estimators_1a)/len(effectivities_J_1a)
median_estimators_DJ_2a = np.sum(DJ_estimators_2a)/len(effectivities_J_1a)
median_estimators_DJ_3a = np.sum(DJ_estimators_3a)/len(effectivities_J_1a)
median_estimators_DJ_4a = np.sum(DJ_estimators_4a)/len(effectivities_J_1a)
median_estimators_DJ_5a = np.sum(DJ_estimators_5a)/len(effectivities_J_1a)<jupyter_output><empty_output><jupyter_text>print tables<jupyter_code>from tabulate import tabulate

# tabulate for the output functional
print('J estimator comparison')
print()
headers = ['Method', 'max J error', 'min J error', 'max J estimate', 'min J estimate', 'average effectivity']
table = [
    ['1a not-corr ', max_J_error_1a, min_J_error_1a, max_J_estimators_1a, min_J_estimators_1a, median_effectivities_J_1a],
    ['2a semi-corr', max_J_error_2a, min_J_error_2a, max_J_estimators_2a, min_J_estimators_2a, median_effectivities_J_2a],
    ['3a AA A-Est ', max_J_error_3a, min_J_error_3a, max_J_estimators_3a, min_J_estimators_3a, median_effectivities_J_3a],
    ['4a AA S-Est ', max_J_error_4a, min_J_error_4a, max_J_estimators_4a, min_J_estimators_4a, median_effectivities_J_4a],
    ['5a Corr SA  ', max_J_error_5a, min_J_error_5a, max_J_estimators_5a, min_J_estimators_5a, median_effectivities_J_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()
# tabulate for the gradient of the output functional
print('DJ estimator comparison')
print()
headers = ['Method', 'max DJ error', 'min DJ error', 'max DJ estimate', 'min DJ estimate', 'average effectivity']
table = [
    ['1a not-corr ', max_DJ_error_1a, min_DJ_error_1a, max_DJ_estimators_1a, min_DJ_estimators_1a, median_effectivities_DJ_1a],
    ['2a semi-corr', max_DJ_error_2a, min_DJ_error_2a, max_DJ_estimators_2a, min_DJ_estimators_2a, median_effectivities_DJ_2a],
    ['3a AA A-Est ', max_DJ_error_3a, min_DJ_error_3a, max_DJ_estimators_3a, min_DJ_estimators_3a, median_effectivities_DJ_3a],
    ['4a AA S-Est ', max_DJ_error_4a, min_DJ_error_4a, max_DJ_estimators_4a, min_DJ_estimators_4a, median_effectivities_DJ_4a],
    ['5a Corr SA  ', max_DJ_error_5a, min_DJ_error_5a, max_DJ_estimators_5a, min_DJ_estimators_5a, median_effectivities_DJ_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f'))
from tabulate import tabulate

# tabulate for the output functional
print('J estimator comparison')
print()
headers = ['Method', 'median J error', 'median J estimate', 'average effectivity']
table = [
    ['1a not-corr ', median_errors_J_1a, median_estimators_J_1a, median_effectivities_J_1a],
    ['2a semi-corr', median_errors_J_2a, median_estimators_J_2a, median_effectivities_J_2a],
    ['3a AA A-Est ', median_errors_J_3a, median_estimators_J_3a, median_effectivities_J_3a],
    ['4a AA S-Est ', median_errors_J_4a, median_estimators_J_4a, median_effectivities_J_4a],
    ['5a Corr SA  ', median_errors_J_5a, median_estimators_J_5a, median_effectivities_J_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()
# tabulate for the gradient of the output functional
print('DJ estimator comparison')
print()
headers = ['Method', 'median DJ error', 'median DJ estimate', 'average effectivity']
table = [
    ['1a not-corr ', median_errors_DJ_1a, median_estimators_DJ_1a, median_effectivities_DJ_1a],
    ['2a semi-corr', median_errors_DJ_2a, median_estimators_DJ_2a, median_effectivities_DJ_2a],
    ['3a AA A-Est ', median_errors_DJ_3a, median_estimators_DJ_3a, median_effectivities_DJ_3a],
    ['4a AA S-Est ', median_errors_DJ_4a, median_estimators_DJ_4a, median_effectivities_DJ_4a],
    ['5a Corr SA  ', median_errors_DJ_5a, median_estimators_DJ_5a, median_effectivities_DJ_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f'))<jupyter_output>DJ estimator comparison

| Method       |   median DJ error |   median DJ estimate |   average effectivity |
|--------------|-------------------|----------------------|-----------------------|
| 1a not-corr  |         8.0175611 |         7445.1147871 |          2436.7337869 |
| 2a semi-corr |         8.0175611 |         7445.1147871 |          2436.7337869 |
| 3a AA A-Est  |         0.8402538 |        18199.1093634 |         58503.3881541 |
| 4a AA S-Est  |         0.8402538 |           92.3805535 |           358.3327077 |
| 5a Corr SA   |         0.2032704 |           69.3861443 |          2531.2768237 |<jupyter_text>primal sensitivities<jupyter_code>max_u_error_4a = max(u_mu_errors_4a)
min_u_error_4a = min(u_mu_errors_4a)

max_u_error_5a = max(u_mu_errors_5a)
min_u_error_5a = min(u_mu_errors_5a)

max_u_estimators_4a = max(u_mu_estimators_4a)
min_u_estimators_4a = min(u_mu_estimators_4a)

max_u_estimators_5a = max(u_mu_estimators_5a)
min_u_estimators_5a = min(u_mu_estimators_5a)

median_effectivities_u_4a = np.sum(effectivities_u_mu_4a)/len(effectivities_J_1a)
median_effectivities_u_5a = np.sum(effectivities_u_mu_5a)/len(effectivities_J_1a)

median_errors_u_4a = np.sum(u_mu_errors_4a)/len(effectivities_J_1a)
median_errors_u_5a = np.sum(u_mu_errors_5a)/len(effectivities_J_1a)

median_estimators_u_4a = np.sum(u_mu_estimators_4a)/len(effectivities_J_1a)
median_estimators_u_5a = np.sum(u_mu_estimators_5a)/len(effectivities_J_1a)
# tabulate for the output functional
print('u estimator comparison')
print()
headers = ['Method', 'max error', 'min error', 'max estimator', 'min estimator']
table = [
    ['4a AA S-Est ', max_u_error_4a, min_u_error_4a, max_u_estimators_4a, min_u_estimators_4a],
    ['5a Corr SA  ', max_u_error_5a, min_u_error_5a, max_u_estimators_5a, min_u_estimators_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()
headers = ['Method', 'average error', 'average estimate', 'average effectivity']
table = [
    ['4a AA S-Est ', median_errors_u_4a, median_estimators_u_4a, median_effectivities_u_4a],
    ['5a Corr SA  ', median_errors_u_5a, median_estimators_u_5a, median_effectivities_u_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()<jupyter_output>| Method      |   average error |   average estimate |   average effectivity |
|-------------|-----------------|--------------------|-----------------------|
| 4a AA S-Est |       0.6694017 |          2.8630712 |            18.5036803 |
| 5a Corr SA  |       0.2373445 |          2.1682993 |           141.1353100 |<jupyter_text>dual sensitivities<jupyter_code>max_p_error_4a = max(p_mu_errors_4a)
min_p_error_4a = min(p_mu_errors_4a)

max_p_error_5a = max(p_mu_errors_5a)
min_p_error_5a = min(p_mu_errors_5a)

#J estimator
max_p_estimators_4a = max(p_mu_estimators_4a)
min_p_estimators_4a = min(p_mu_estimators_4a)

max_p_estimators_5a = max(p_mu_estimators_5a)
min_p_estimators_5a = min(p_mu_estimators_5a)

median_effectivities_p_4a = np.sum(effectivities_p_mu_4a)/len(effectivities_J_1a)
median_effectivities_p_5a = np.sum(effectivities_p_mu_5a)/len(effectivities_J_1a)

median_errors_p_4a = np.sum(p_mu_errors_4a)/len(effectivities_J_1a)
median_errors_p_5a = np.sum(p_mu_errors_5a)/len(effectivities_J_1a)

median_estimators_p_4a = np.sum(p_mu_estimators_4a)/len(effectivities_J_1a)
median_estimators_p_5a = np.sum(p_mu_estimators_5a)/len(effectivities_J_1a)
# tabulate for the output functional
print('p estimator comparison')
print()
headers = ['Method', 'max error', 'min error', 'max estimator', 'min estimator']
table = [
    ['4a AA S-Est ', max_p_error_4a, min_p_error_4a, max_p_estimators_4a, min_p_estimators_4a],
    ['5a Corr SA  ', max_p_error_5a, min_p_error_5a, max_p_estimators_5a, min_p_estimators_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()
headers = ['Method', 'average error', 'average estimate', 'average effectivity']
table = [
    ['4a AA S-Est ', median_errors_p_4a, median_estimators_p_4a, median_effectivities_p_4a],
    ['5a Corr SA  ', median_errors_p_5a, median_estimators_p_5a, median_effectivities_p_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f')) 
print()
max_eff_J_1a = max(effectivities_J_1a)
max_eff_J_2a = max(effectivities_J_2a)
max_eff_J_3a = max(effectivities_J_3a)
max_eff_J_4a = max(effectivities_J_4a)
max_eff_J_5a = max(effectivities_J_5a)

max_eff_DJ_1a = max(effectivities_DJ_1a)
max_eff_DJ_2a = max(effectivities_DJ_2a)
max_eff_DJ_3a = max(effectivities_DJ_3a)
max_eff_DJ_4a = max(effectivities_DJ_4a)
max_eff_DJ_5a = max(effectivities_DJ_5a)

max_eff_u_4a = max(effectivities_u_mu_4a)
max_eff_u_5a = max(effectivities_u_mu_5a)

max_eff_p_4a = max(effectivities_p_mu_4a)
max_eff_p_5a = max(effectivities_p_mu_5a)
print('max effectivities')

headers = ['Method', 'max eff for J', 'max eff for DJ', 'max eff for u_mu', 'max eff for p_mu']
table = [
    ['1a not-corr ', max_eff_J_1a, max_eff_DJ_1a, '-', '-'],
    ['2a semi-corr', max_eff_J_2a, max_eff_DJ_2a, '-', '-'],
    ['3a AA A-Est ', max_eff_J_3a, max_eff_DJ_3a, '-', '-'],
    ['4a AA S-Est ', max_eff_J_4a, max_eff_DJ_4a, max_eff_u_4a, max_eff_p_4a],
    ['5a Corr SA  ', max_eff_J_5a, max_eff_DJ_5a, max_eff_u_5a, max_eff_p_5a]]

print(tabulate(table, headers=headers, tablefmt='github', floatfmt='.7f'))<jupyter_output>max effectivities
| Method       |   max eff for J |   max eff for DJ | max eff for u_mu   | max eff for p_mu   |
|--------------|-----------------|------------------|--------------------|--------------------|
| 1a not-corr  |    8744.1808226 |    43287.5616865 | -                  | -                  |
| 2a semi-corr |    8930.1670963 |    43287.5616866 | -                  | -                  |
| 3a AA A-Est  |    8930.1670963 |   684164.9871887 | -                  | -                  |
| 4a AA S-Est  |    8930.1670963 |     6818.7796795 | 55.455289395002836 | 4427.108668360115  |
| 5a Corr SA   |    8930.1670963 |    94919.3914497 | 1890.945110056238  | 17523.007090085026 |
{'fraction_non_alphanumeric': 0.10333917816573912, 'fraction_numerical': 0.04311199207135778, 'mean_word_length': 3.6147757255936677, 'pattern_counts': {'":': 0, '<': 50, '<?xml version=': 0, '>': 49, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 85, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '14', 'dir': 'jupyter-structured-clean-dedup', 'id': '73989', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 10654, 'n_tokens_neox': 9737, 'n_words': 1444, 'path': 'notebooks/Paper1_simulations/Model_Problem_2_Estimator_study/EXC10-estimator_study_with_basis_size_16.ipynb', 'repo_name': 'TiKeil/NCD-corrected-TR-RB-approach-for-pde-opt', 'size': '26243'}