Works also for precision=2.

prob_laplace.m

@@ -98,12 +98,12 @@ b = @(chi, nu) nu^(1/3)*sqrt(chi/(chi-2));
 e = @(chi) sqrt(3/8*(chi-2));
 t = [0, 0.05, 0.1];
 T1 = Ternary_model(0, 'Gauss', [-5, b(chi, nu), 0.5, e(chi), ...
-                     0.16, 0.2, 10, a], t, 3);
-x0 = T1.x(find(T1.x>a, 1):end);
-x0 = [x0(1:20:8000), x0(8001:2:end)];
+                     0.16, 0.2, 10, a], t, 2);
+% x0 = T1.x(find(T1.x>a, 1):end);
+% x0 = [x0(1:20:8000), x0(8001:2:end)];
 x0 = x0(x0<7.1);
 % x0 = x0(1:10:end);
-% x0 = 0.0:0.001:0.1;
+x0 = 5.0:0.001:6.01;
 %%
 T = {};
 tic
@@ -111,11 +111,11 @@ parfor i = 1:length(x0)
 b = @(chi, nu) nu^(1/3)*sqrt(chi/(chi-2));
 e = @(chi) sqrt(3/8*(chi-2));
 T{i} = Ternary_model(0, 'Gauss', [-a, b(chi, nu), 0.5, e(chi), ...
-                     0.16, 0.2, 10, x0(i)], t, 3);
+                     0.16, 0.2, 10, x0(i)], t, 2);
 T{i}.solve_tern_frap();
 end
 toc
-save prob_laplace_X_7_7
+save prob_laplace_X_7_7_short
 %%
 ls = 0.0005:0.005:1;
 p = nan(1, length(ls));
@@ -124,7 +124,8 @@ p(i) = int_prob(ls(i), T, x0);
 end
 %% Do we need to look at the left side as well?
 figure; hold on;
-plot(ls, p);
+plot(ls, N*p);
+% plot(lx, xt);
 % plot(ls, q/50000);
 %%
 int_prob(0.01, T, x0+a)
@@ -136,11 +137,11 @@ for i = 1:2:length(T)
 end
 %% Normalization factor
 tic
-N = normalization(T, x0+a);
+N = normalization(T, x0);
 toc
 %% Write to file
-csvwrite('jump_length.csv', ls)
-csvwrite('prob.csv', q);
+csvwrite('jump_length_7_7_3.csv', ls)
+csvwrite('prob_7_7_3.csv', p);
 %% distribution for x0 can be taken from phi_tot (steady state)
 function p = int_prob(l, T, x0)
 delta_x0 = diff(x0);