% Generate efficiency plot.

clear;

%% DATA
log_instantaneous_acceptance = -76.5;
nsteps_list = [ 1        2        4        8       16       32       64      128      256      512     1024     2048 ];
log_mean_acceptance = [-76.4    -78.1    -80.3    -59.9    -41.1    -28.9    -15.3    -10.6     -6.8     -5.0     -3.4     -2.4];
log_relative_efficiency = [ 0.1     -2.3     -5.2     14.5     32.6     44.2     57.0     61.1     64.1     65.2     66.2     66.4 ];

log_instantaneous_acceptance = -25.7;
nsteps_list = [1        2        4        8       16       32       64      128      256      512     1024     2048];
log_mean_acceptance = [-25.7    -26.4    -27.6    -30.4    -33.7    -36.1    -18.0    -14.7     -6.1     -4.7     -3.1     -2.2];
log_mean_acceptance_error = [1.000    1.000    1.000    1.000    1.000    0.721    0.783    0.302    0.370    0.156    0.073    0.042];
log_relative_efficiency = [-0.7     -1.7     -3.5     -6.8    -10.8    -13.8      3.5      6.2     14.0     14.8     15.7     15.9 ];

% 1770 iterations
log_instantaneous_acceptance = -70.1;
log_mean_acceptance = [-70.1    -70.7    -72.9    -74.5    -60.5    -33.0    -16.9    -12.8     -8.2     -4.5     -3.1     -2.0];
log_mean_acceptance_error = [0.000    0.000    0.000    0.000    0.000    0.989    0.990    0.724    0.578    0.184    0.091    0.051];
log_relative_efficiency = [-0.6     -1.7     -4.4     -6.6      6.8     33.6     49.1     52.5     56.4     59.4     60.1     60.5];

% 4395 iterations
log_instantaneous_acceptance = -41.1;
log_mean_acceptance = [ -41.2    -42.7    -45.8    -53.2    -45.5    -24.7    -17.6    -10.1     -7.2     -4.5     -3.1     -2.1];
log_mean_acceptance_error = [1.000    1.000    1.000    1.000    1.000    1.000    0.829    0.775    0.382    0.116    0.058    0.033];
log_relative_efficiency = [ -0.7     -2.7     -6.3    -14.3     -7.2     12.9     19.3     26.1     28.4     30.4     31.1     31.4];

% 5190 iterations
log_instantaneous_acceptance = -41.3;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048 ];
log_mean_acceptance = [ -41.3 -42.8 -46.0 -53.4 -45.7 -24.9 -17.8 -10.3  -7.1  -4.4  -3.1  -2.1 ];
log_mean_acceptance_lower = [ -81.3 -83.0 -86.6 -80.5 -61.5 -38.7 -21.1 -12.9  -8.0  -4.6  -3.2  -2.2 ];
log_mean_acceptance_upper = [ -40.2 -41.7 -44.9 -52.3 -44.6 -23.8 -16.8  -9.3  -6.5  -4.2  -3.0  -2.0 ];
log_relative_efficiency = [  -0.7  -2.7  -6.3 -14.3  -7.2  12.9  19.3  26.1  28.7  30.6  31.2  31.6 ];
log_relative_efficiency_lower = [  -0.7  -2.7  -6.3 -14.3 -23.1  -0.8  15.8  23.3  27.2  29.4  30.0  30.4 ];
log_relative_efficiency_upper = [  -0.6  -1.7  -4.3   5.7  32.0  48.8  59.3  65.7  68.7  70.7  71.3  71.6 ];

% 6085 iterations
log_instantaneous_acceptance = -41.4;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048 ];
log_mean_acceptance = [ -41.5 -43.0 -46.2 -53.6 -45.8 -25.1 -17.8 -10.4  -7.0  -4.5  -3.1  -2.1 ];
log_mean_acceptance_lower = [ -73.4 -74.4 -77.7 -75.7 -61.7 -38.7 -20.4 -13.0  -7.8  -4.7  -3.2  -2.1 ];
log_mean_acceptance_upper = [ -40.4 -41.9 -45.1 -52.5 -44.7 -24.0 -16.8  -9.5  -6.5  -4.3  -3.1  -2.0 ];
log_relative_efficiency = [  -0.7  -2.7  -6.3 -14.3  -7.2  12.9  19.4  26.1  28.9  30.7  31.4  31.7 ];
log_relative_efficiency_lower = [  -0.7  -2.7  -6.3 -14.3 -23.1  -0.5  16.6  23.3  27.4  29.5  30.2  30.6 ];
log_relative_efficiency_upper = [  -0.5  -0.1   0.5   0.2  24.4  44.9  51.4  58.0  60.7  62.7  63.4  63.8 ];

% 6956 iterations
log_instantaneous_acceptance = -41.6;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048 ];
log_mean_acceptance = [ -41.6 -43.1 -46.3 -53.7 -46.0 -25.2 -18.0 -10.5  -6.7  -4.5  -3.1  -2.1 ];
log_mean_acceptance_lower = [ -73.5 -74.2 -74.2 -71.2 -60.2 -38.3 -20.5 -12.7  -7.4  -4.7  -3.2  -2.1 ];
log_mean_acceptance_upper = [ -40.5 -42.0 -45.2 -52.6 -44.9 -24.1 -16.9  -9.6  -6.3  -4.3  -3.0  -2.0 ];
log_relative_efficiency = [  -0.7  -2.7  -6.3 -14.3  -7.2  12.9  19.4  26.2  29.3  30.8  31.5  31.9 ];
log_relative_efficiency_lower = [  -0.7  -2.7  -6.3 -14.3 -17.4   0.2  16.6  23.5  27.9  29.6  30.3  30.7 ];
log_relative_efficiency_upper = [  -0.5  -0.1   1.2   8.6  24.8  44.9  51.5  58.3  61.4  62.9  63.6  63.9 ];

% 8108 iterations
log_instantaneous_acceptance = -41.7;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048 ];
log_mean_acceptance = [ -41.8 -43.3 -46.5 -53.8 -46.1 -25.4 -18.1 -10.7  -6.8  -4.5  -3.1  -2.1 ];
log_mean_acceptance_lower = [ -73.7 -74.3 -76.9 -71.4 -59.6 -38.4 -20.9 -12.8  -7.3  -4.7  -3.2  -2.1 ];
log_mean_acceptance_upper = [ -40.7 -42.2 -45.4 -52.7 -45.0 -24.3 -17.1  -9.7  -6.4  -4.4  -3.0  -2.0 ];
log_relative_efficiency = [  -0.7  -2.7  -6.3 -14.3  -7.2  12.9  19.4  26.2  29.4  30.9  31.7  32.0 ];
log_relative_efficiency_lower = [  -0.7  -2.7  -6.3 -14.3 -20.5   2.0  16.7  23.7  28.1  29.8  30.5  30.9 ];
log_relative_efficiency_upper = [  -0.5  -0.1   0.5   7.8  24.3  44.5  51.3  58.1  61.5  63.0  63.7  64.1 ];

% 10000 iterations
log_instantaneous_acceptance = -41.93;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048 ];
log_mean_acceptance = [ -41.99 -43.50 -46.66 -54.06 -46.32 -25.56 -18.32 -10.76  -6.73  -4.55  -3.13  -2.07 ];
log_mean_acceptance_lower = [ -73.90 -74.50 -74.66 -71.59 -55.44 -38.46 -21.24 -12.41  -7.24  -4.71  -3.21  -2.11 ];
log_mean_acceptance_upper = [ -40.89 -42.40 -45.56 -52.96 -45.23 -24.46 -17.23  -9.91  -6.34  -4.40  -3.06  -2.03 ];
log_relative_efficiency = [  -0.75  -2.66  -6.34 -14.32  -7.22  12.87  19.44  26.31  29.65  31.14  31.87  32.23 ];
log_relative_efficiency_lower = [  -0.75  -2.66  -6.34 -14.32 -16.98   0.67  16.59  24.21  28.34  29.91  30.69  31.11 ];
log_relative_efficiency_upper = [  -0.54  -0.07   0.52   8.64  24.21  44.31  51.34  58.41  61.62  63.18  63.92  64.28 ];

% 79 iterations - up to 8192 steps
%log_instantaneous_acceptance = -243.61;
%nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
%log_mean_acceptance = [ -237.38 -215.51 -180.71 -128.30 -88.61 -51.04 -32.05 -18.56 -10.58  -3.96  -3.19  -2.30  -1.35  -0.84 ];
%log_mean_acceptance_lower = [ -401.84 -373.48 -321.03 -254.58 -123.04 -70.63 -39.64 -22.86 -13.22  -6.12  -4.42  -2.80  -1.71  -1.04 ];
%log_mean_acceptance_upper = [ -236.28 -214.41 -179.61 -127.20 -87.51 -49.94 -30.96 -17.49  -9.70  -3.15  -2.50  -1.91  -1.10  -0.68 ];
%log_relative_efficiency = [   5.54  27.00  61.29 113.12 152.17 189.08 207.39 220.19 227.49 233.41 233.49 233.69 233.95 233.76 ];
%log_relative_efficiency_lower = [   5.54  14.16  22.26  60.81 152.17 170.37 199.80 215.73 224.78 231.17 231.80 232.35 232.83 232.57 ];
%log_relative_efficiency_upper = [  58.20  86.15 148.74 249.90 336.54 392.65 420.82 436.60 444.26 450.24 450.63 450.66 450.81 450.84 ];

% 1159 iterations - up to 8192 steps
log_instantaneous_acceptance = -105.73;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -105.44 -106.59 -106.87 -101.26 -71.32 -36.42 -22.56 -11.58  -6.44  -4.34  -3.02  -2.04  -1.42  -0.93 ];
log_mean_acceptance_lower = [ -168.83 -168.21 -163.76 -130.98 -82.96 -47.07 -27.71 -16.30  -8.98  -4.82  -3.25  -2.16  -1.50  -0.99 ];
log_mean_acceptance_upper = [ -104.34 -105.49 -105.77 -100.16 -70.24 -35.32 -21.48 -10.68  -5.60  -4.01  -2.82  -1.92  -1.34  -0.88 ];
log_relative_efficiency = [  -0.40  -1.95  -2.74   2.28  31.58  65.82  79.00  89.30  93.75  95.15  95.78  96.07  96.00  95.79 ];
log_relative_efficiency_lower = [  -0.40  -3.30  -2.74   2.28  24.61  55.11  73.75  84.57  91.29  93.79  94.54  94.89  94.84  94.65 ];
log_relative_efficiency_upper = [   2.64   6.61  19.80  50.46  93.01 130.18 140.25 153.77 158.09 159.80 160.31 160.68 160.64 160.36 ];

% 2268 iterations - up to 8192 steps
log_instantaneous_acceptance = -79.44;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -79.11 -77.63 -77.33 -79.60 -71.81 -37.08 -23.06 -12.16  -6.55  -4.58  -3.08  -2.08  -1.38  -0.94 ];
log_mean_acceptance_lower = [ -128.81 -129.90 -125.27 -102.36 -74.59 -46.69 -25.79 -14.72  -7.67  -4.98  -3.26  -2.17  -1.43  -0.98 ];
log_mean_acceptance_upper = [ -77.73 -76.25 -75.94 -78.52 -70.81 -35.99 -22.04 -11.31  -5.94  -4.28  -2.94  -2.00  -1.32  -0.91 ];
log_relative_efficiency = [  -0.36   0.71   0.51  -2.35   4.80  38.87  52.21  62.43  67.35  68.62  69.43  69.74  69.75  69.49 ];
log_relative_efficiency_lower = [  -0.65  -2.18  -2.74  -4.99   1.84  31.92  49.20  59.59  65.74  67.31  68.12  68.38  68.40  68.12 ];
log_relative_efficiency_upper = [   0.54   5.23  13.15  24.27  55.16  88.24 102.57 112.38 117.26 118.34 118.97 119.22 119.14 118.90 ];

% 3326 iterations - up to 8192 steps
log_instantaneous_acceptance = -79.83;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -79.50 -78.01 -77.71 -79.97 -65.12 -37.43 -18.37 -12.49  -6.68  -4.63  -3.10  -2.11  -1.40  -0.94 ];
log_mean_acceptance_lower = [ -120.90 -118.28 -108.90 -101.60 -72.18 -42.29 -25.35 -14.74  -7.51  -4.95  -3.23  -2.18  -1.44  -0.98 ];
log_mean_acceptance_upper = [ -78.40 -76.92 -76.61 -78.90 -64.01 -36.35 -16.99 -11.67  -6.11  -4.39  -2.97  -2.04  -1.35  -0.91 ];
log_relative_efficiency = [  -0.36   0.71   0.51  -2.34  11.88  38.90  57.28  62.47  67.60  68.95  69.80  70.09  70.11  69.87 ];
log_relative_efficiency_lower = [  -0.40  -1.95  -2.74  -4.99   4.82  33.88  50.16  59.92  66.18  67.73  68.64  68.96  69.00  68.75 ];
log_relative_efficiency_upper = [   0.54   5.23  13.15  28.70  49.47  76.49  94.60 102.86 109.60 110.74 111.83 112.18 112.16 111.90 ];

% 4466 iterations - up to 8192 steps
log_instantaneous_acceptance = -68.84;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -68.31 -63.65 -56.19 -46.83 -38.40 -35.42 -18.66 -10.74  -6.49  -4.59  -3.07  -2.11  -1.39  -0.95 ];
log_mean_acceptance_lower = [ -105.80 -99.59 -91.17 -82.91 -67.67 -40.95 -24.95 -13.36  -7.18  -4.85  -3.19  -2.18  -1.43  -0.98 ];
log_mean_acceptance_upper = [ -67.21 -62.55 -55.09 -45.73 -37.32 -34.36 -17.57  -9.76  -6.00  -4.36  -2.97  -2.05  -1.35  -0.92 ];
log_relative_efficiency = [  -0.17   4.09  11.04  19.81  27.60  29.91  46.00  53.24  56.80  58.01  58.83  59.10  59.13  58.88 ];
log_relative_efficiency_lower = [  -0.36   0.71   0.51  -4.91   9.59  24.27  39.32  50.39  55.27  56.73  57.59  57.92  57.98  57.74 ];
log_relative_efficiency_upper = [   0.56  10.13  26.25  47.00  65.16  67.72  83.65  91.32  94.65  95.70  96.57  96.85  96.80  96.56 ];

% 9247 iterations - up to 8192 steps
log_instantaneous_acceptance = -62.02;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -62.06 -63.78 -56.92 -47.55 -39.13 -21.63 -18.47  -8.79  -6.46  -4.50  -3.07  -2.10  -1.38  -0.96 ];
log_mean_acceptance_lower = [ -69.04 -70.91 -78.29 -78.28 -52.17 -33.94 -20.12 -12.34  -6.97  -4.68  -3.15  -2.15  -1.41  -0.97 ];
log_mean_acceptance_upper = [ -60.93 -62.81 -55.82 -46.45 -38.05 -20.54 -17.70  -7.74  -6.10  -4.36  -3.00  -2.06  -1.36  -0.93 ];
log_relative_efficiency = [  -0.74  -2.86   3.49  12.27  20.05  36.89  39.37  48.37  50.01  51.27  52.02  52.29  52.32  52.05 ];
log_relative_efficiency_lower = [  -0.80  -5.45 -15.20 -16.40   6.59  24.61  37.49  44.69  48.79  50.14  50.90  51.15  51.20  50.94 ];
log_relative_efficiency_upper = [  -0.17   4.09  11.04  19.81  27.60  45.12  47.23  56.33  57.72  58.91  59.61  59.84  59.88  59.61 ];

% 10000 iterations - up to 8192 steps
log_instantaneous_acceptance = -62.10;
nsteps_list = [      1      2      4      8     16     32     64    128    256    512   1024   2048   4096   8192 ];
log_mean_acceptance = [ -62.14 -63.86 -56.99 -47.63 -39.21 -21.71 -18.55  -8.87  -6.52  -4.51  -3.07  -2.10  -1.38  -0.96 ];
log_mean_acceptance_lower = [ -69.12 -70.99 -78.02 -74.79 -52.25 -34.23 -20.13 -12.40  -7.00  -4.69  -3.15  -2.15  -1.41  -0.98 ];
log_mean_acceptance_upper = [ -61.01 -62.86 -55.90 -46.53 -38.11 -20.61 -17.84  -7.86  -6.14  -4.37  -3.00  -2.06  -1.36  -0.94 ];
log_relative_efficiency = [  -0.74  -2.86   3.49  12.27  20.05  36.89  39.37  48.37  50.03  51.35  52.09  52.37  52.40  52.13 ];
log_relative_efficiency_lower = [  -0.80  -5.45 -15.20 -12.92   6.09  24.45  37.57  44.78  48.82  50.24  50.97  51.26  51.29  51.02 ];
log_relative_efficiency_upper = [  -0.17   4.09  11.04  19.81  27.60  44.43  46.53  55.49  57.47  58.84  59.63  59.91  59.93  59.68 ];

% Compute modelled correlation times.
tau_md = 299.8; % correlation time of MD alone
tau_ncmc = (2 * exp(log_mean_acceptance)).^(-1);
%tau_eff = (tau_md.^(-1) + tau_ncmc.^(-1) + tau_md.^(-1) .* tau_ncmc.^(-1)).^(-1);
tau_eff = tau_md .* (tau_ncmc ./ (tau_md + tau_ncmc));
g_md = 1 + 2*tau_md;
g_ncmc = 1 + 2*tau_ncmc;
g_eff = 1 + 2*tau_eff;
nsteps_md = 500;
efficiency_tau = (g_md ./ g_eff) .* nsteps_md ./ (nsteps_md + nsteps_list);
efficiency_ncmc_alone = (g_md ./ g_ncmc) .* nsteps_md ./ (nsteps_md + nsteps_list);
efficiency_bound = (g_md / 1) .* nsteps_md ./ (nsteps_md + nsteps_list);
force_evals_per_uncorrelated_sample = g_eff .* (nsteps_md + nsteps_list);
force_evals_per_uncorrelated_sample_ncmc = g_ncmc .* nsteps_list;
force_evals_per_uncorrelated_sample_md = g_md .* nsteps_md;

%% Generate plot.
fontsize = 9;

figure(1);
clf;
scale = 1.4; % original Fig. 3
%scale = 2.0;

% Plot acceptance probabilities.
set(gca, 'fontsize', fontsize);

xmin = 0.5;
xmax = nsteps_list(end) * 2;
ymin = log_instantaneous_acceptance - 1;
ymax = 0.0;

%ymin = exp(min(log_mean_acceptance) - 15);
%ymax = 1.0;

ymin = 1e-35;
ymax = 1e-0;

%plot(nsteps_list, log_mean_acceptance, 'k.');
loglog(nsteps_list, exp(log_mean_acceptance), 'k.', 'markersize', 12);

% error bars
hold on;
%loglog([nsteps_list;nsteps_list], [exp(log_mean_acceptance+2*log_mean_acceptance_error);exp(log_mean_acceptance-2*log_mean_acceptance_error)], 'k-');
loglog([nsteps_list;nsteps_list], [exp(log_mean_acceptance_upper);exp(log_mean_acceptance_lower)], 'k-');

axis([xmin xmax ymin ymax]);
hold on;
loglog([xmin xmax], exp(log_instantaneous_acceptance) * [1 1], 'r-');
xlabel('NCMC switching steps');
ylabel('acceptance probability');
set(gca, 'XTick', nsteps_list)
set(gca, 'YTick', [1e-30, 1e-20, 1e-10, 1e-0]);

% INSET
axes('pos', [0.6 0.19 0.305 0.50]);
i = 10;
semilogx(nsteps_list(i:end), exp(log_mean_acceptance(i:end)), 'k.', 'markersize', 12);
hold on;
semilogx([nsteps_list(i:end);nsteps_list(i:end)], [exp(log_mean_acceptance_upper(i:end));exp(log_mean_acceptance_lower(i:end))], 'k-');
set(gca, 'XTick', nsteps_list(i:end));
set(gca, 'XTickLabel', '');
set(gca, 'YTick', [0 0.1 0.2 0.3 0.4 0.5]);
axis([nsteps_list(i)/2 nsteps_list(end)*2 0 0.5]);
% END INSET

% Generate PDF.
filename = 'wca-dimer-acceptance.eps';
%print('-depsc', filename);
exportfig(gcf, filename, 'color', 'cmyk', 'width', scale*4.5, 'height', scale*3.5/2.0, 'renderer', 'painters');
system(sprintf('epstopdf %s', filename));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

figure(2);
clf;

% Plot modelled correlation times.
subplot(2,1,1);
set(gca, 'fontsize', fontsize);
%plot(nsteps_list, tau_eff, 'k.');
semilogx(nsteps_list, tau_eff, 'k.', 'markersize', 12);
%loglog(nsteps_list, tau_eff, 'k.', 'markersize', 12);

% reference data from simulations
%hold on;
%semilogx(2048, 4.0, 'r.', 'markersize', 12);
%loglog(2048, 4.0, 'r.', 'markersize', 12);

% Show NCMC alone (DEBUG).
%hold on;
%semilogx(nsteps_list, tau_ncmc, 'r.', 'markersize', 12);
%legend('MD+NCMC', 'NCMC alone');

% error bars
%hold on;
%loglog([nsteps_list;nsteps_list], [exp(log_relative_efficiency_upper);exp(log_relative_efficiency_lower)], 'k-');

xmin = 1;
xmin = 0.5;
xmax = nsteps_list(end) * 2;
%ymin = exp(min(log_relative_efficiency) - 1);
%ymax = exp(max(log_relative_efficiency) + 10);
%ymin = 1e-10;
%ymax = 1e+30;
ymin = 0;
ymax = 310;

axis([xmin xmax ymin ymax]);

% Show MD correlation time.
hold on;
semilogx([xmin xmax], tau_md*[1 1], 'r-', 'markersize', 12);

xlabel('NCMC switching steps');
ylabel('correlation time \tau_{eff}');
set(gca, 'XTick', nsteps_list)
set(gca, 'YTick', [0 100 200 300 400]);

% Plot actual efficiency.
subplot(2,1,2);
set(gca, 'fontsize', fontsize);
%plot(nsteps_list, tau_eff, 'k.');
%semilogx(nsteps_list, efficiency_tau, 'k.', 'markersize', 12);
semilogx(nsteps_list, force_evals_per_uncorrelated_sample_md ./ force_evals_per_uncorrelated_sample, 'k.', 'markersize', 12);
%loglog(nsteps_list, force_evals_per_uncorrelated_sample_md ./ force_evals_per_uncorrelated_sample, 'k.', 'markersize', 12);

% Show eNCMC alone (DEBUG).
%hold on;
%semilogx(nsteps_list, force_evals_per_uncorrelated_sample_md ./ force_evals_per_uncorrelated_sample_ncmc, 'r.', 'markersize', 12);
%legend('MD+NCMC', 'NCMC alone');

% error bars
%hold on;
%loglog([nsteps_list;nsteps_list], [exp(log_relative_efficiency_upper);exp(log_relative_efficiency_lower)], 'k-');

xmin = 1;
xmin = 0.5;
xmax = nsteps_list(end) * 2;
%ymin = exp(min(log_relative_efficiency) - 1);
%ymax = exp(max(log_relative_efficiency) + 10);
%ymin = 1e-10;
%ymax = 1e+30;
ymin = 0;
%ymax = max(force_evals_per_uncorrelated_sample);
ymax = 16;

% plot bound
hold on;
%semilogx(nsteps_list, efficiency_bound, 'r-', 'markersize', 12);
%semilogx([nsteps_list(1) nsteps_list(end)], force_evals_per_uncorrelated_sample_md * [1 1], 'r-', 'markersize', 12);
semilogx([xmin xmax], [1 1], 'r-', 'markersize', 12);
%loglog([nsteps_list(1) nsteps_list(end)], [1 1], 'r-', 'markersize', 12);

axis([xmin xmax ymin ymax]);

xlabel('NCMC switching steps');
ylabel('relative efficiency');
set(gca, 'XTick', nsteps_list)
%set(gca, 'YTick', [0 5 10 15 20]);
set(gca, 'YTick', [0 5 10 15 20]);
%set(gca, 'YTick', [0.5 1 2 4 8 16]);


% INSET
axes('pos', [0.172 0.20 0.452 0.22]);
i = 9;
semilogx(nsteps_list, force_evals_per_uncorrelated_sample_md ./ force_evals_per_uncorrelated_sample, 'k.', 'markersize', 12);
hold on;
semilogx([xmin xmax], [1 1], 'r-', 'markersize', 12);
set(gca, 'XTick', nsteps_list(1:i));
set(gca, 'YTick', [0 0.5 1 1.5 2]);
set(gca, 'XTickLabel', '');
set(gca, 'YTickLabel', {'0', '', '1', '', '2'});
%set(gca, 'YTick', [0 0.1 0.2 0.3 0.4 0.5]);
axis([nsteps_list(1)-0.1 nsteps_list(i)*1.5 0 2]);
% END INSET

% Generate PDF.
filename = 'wca-dimer-efficiency-tau.eps';
%print('-depsc', filename);
exportfig(gcf, filename, 'color', 'cmyk', 'width', scale*4.5, 'height', scale*3.5, 'renderer', 'painters');
system(sprintf('epstopdf %s', filename));

