From a9c0d6b59de1e3fd22667ee894d44d7072c74b7a Mon Sep 17 00:00:00 2001
From: Lars Hubatsch <hubatsch@pks.mpg.de>
Date: Tue, 29 Sep 2020 15:12:50 +0200
Subject: [PATCH] Comparison ML Python works. Spherical still off. Slab fine.
---
@Ternary_model/create_mesh.m | 2 +-
@Ternary_model/plot_sim.m | 14 +++++------
@Ternary_model/solve_tern_frap.m | 2 ++
FloryHugg_DiffUnbleached.ipynb | 43 ++++++++++++++++++++++++++------
prob_laplace.m | 19 +++++++++++++-
5 files changed, 63 insertions(+), 17 deletions(-)
diff --git a/@Ternary_model/create_mesh.m b/@Ternary_model/create_mesh.m
index 74c139c..5282414 100644
--- a/@Ternary_model/create_mesh.m
+++ b/@Ternary_model/create_mesh.m
@@ -4,7 +4,7 @@ function create_mesh(T)
if strcmp(T.mode, 'Constituent') | T.v == 0
% Start with left side of mesh, at r=0/x=0, with a given step size
- x = linspace(0, 0.75, 40);
+ x = linspace(0, -0.75*T.a, 40);
x = [x, x(end)+0.001/T.precision];
while x(end)<-T.a
x_temp = x(end)-x(end-1);
diff --git a/@Ternary_model/plot_sim.m b/@Ternary_model/plot_sim.m
index 72f7768..d9de961 100644
--- a/@Ternary_model/plot_sim.m
+++ b/@Ternary_model/plot_sim.m
@@ -12,26 +12,26 @@ if strcmp(mode, 'mov')
for i = 1:nth:length(T.t)
cla;
axis([0, T.system_size, 0, max(T.sol(:))]);
- ax = gca;
- ax.FontSize = 12;
+ ax = gca; ax.FontSize = 12;
xlabel('position'); ylabel('probability [not normalized]');
plot(T.x, T.phi_tot(T.x, T.a, T.b, T.e, T.u0, T.t(i)*T.v), 'LineWidth', 2, ...
'LineStyle', '--');
plot(T.x, T.sol(i, :), 'LineWidth', 2);
text(abs(-T.a+2), 0.7, num2str(T.t(i)));
- shg
- pause();
+ shg; pause();
% print([num2str(i),'.png'],'-dpng')
end
elseif strcmp(mode, 'plot')
figure(20);
- cla;
+% cla;
hold on;
-% xlim([-inf, 6.0]); ylim([-inf, max(T.sol(:))]);
+ xlim([-inf, -2*T.a]); %ylim([-inf, max(T.sol(:))]);
ax = gca;
ax.FontSize = 12;
xlabel('x [\mum]'); ylabel('volume fraction');
- plot(T.x, T.sol([1, 2:nth:si(1)], :), 'LineWidth', 1, 'Color', color);
+% N = max(max(T.sol([1, 2:nth:si(1)], :)));
+ N = 1;
+ plot(T.x, T.sol([1, 1:nth:si(1)], :)/N, 'LineWidth', 1, 'Color', color);
% plot(T.x, T.sol(nth, :), 'LineWidth', 1, 'Color', color);
plot(T.x, Ternary_model.phi_tot(T.x, T.a, T.b, T.e, T.u0, 0),...
'LineWidth', 1, 'LineStyle', '--', 'Color', [0.97, 0.55, 0.03]);
diff --git a/@Ternary_model/solve_tern_frap.m b/@Ternary_model/solve_tern_frap.m
index 837fdbc..a1fde6a 100644
--- a/@Ternary_model/solve_tern_frap.m
+++ b/@Ternary_model/solve_tern_frap.m
@@ -1,12 +1,14 @@
function solve_tern_frap(T)
%SOLVE_TERN_FRAP solves Fokker Planck equation for ternary FRAP
% Assumes fixed profile of phi_tot
+tic
fh_ic = @(x) flory_ic(x, T.a, T.b, T.u0, T.e, T.ic, T.x0);
fh_bc = @(xl, ul, xr, ur, t) flory_bc(xl, ul, xr, ur, t, T.u0, T.e);
fh_pde = @(x, t, u, dudx) flory_hugg_pde(x, t, u, dudx, T.a, T.b, T.e,...
T.u0, T.e_g0, T.u_g0, T.v,...
T.mode);
T.sol = pdepe(T.geometry, fh_pde, fh_ic, fh_bc, T.x, T.t);
+toc
end
%% Function definitions for pde solver
diff --git a/FloryHugg_DiffUnbleached.ipynb b/FloryHugg_DiffUnbleached.ipynb
index c872512..63c1eca 100644
--- a/FloryHugg_DiffUnbleached.ipynb
+++ b/FloryHugg_DiffUnbleached.ipynb
@@ -94,8 +94,9 @@
"outputs": [],
"source": [
"bp = 0.1 # Boundary position at IC\n",
- "sym = 2 # Symmetry of the problem\n",
- "dt = 0.0001\n",
+ "sym = 0 # Symmetry of the problem\n",
+ "dt = 0.001\n",
+ "nt = 10\n",
"\n",
"c_tot_1 = create_func(F, p_tot(0.9, 0.9), 1)\n",
"c0_1 = create_func(F, 'x[0]<'+str(bp)+'? 0 :' + p_tot(0.9, 0.9), 1)\n",
@@ -121,12 +122,12 @@
"c0_9 = create_func(F, 'x[0]<'+str(bp)+'? 0 :' +p_tot(0.9, 0.1), 1)\n",
"Ga0_9 = create_func(F,p_tot(1, 1), 1)\n",
"\n",
- "c0_1 = calc_sim(c0_1, c_tot_1, Ga0_1, dt, 10, sym)\n",
- "c0_2 = calc_sim(c0_2, c_tot_2, Ga0_2, dt, 10, sym)\n",
- "c0_3 = calc_sim(c0_3, c_tot_3, Ga0_3, dt, 10, sym)\n",
- "c0_5 = calc_sim(c0_5, c_tot_5, Ga0_5, dt, 10, sym)\n",
- "c0_8 = calc_sim(c0_8, c_tot_8, Ga0_8, dt, 10, sym)\n",
- "c0_9 = calc_sim(c0_9, c_tot_9, Ga0_9, dt, 10, sym)"
+ "c0_1 = calc_sim(c0_1, c_tot_1, Ga0_1, dt, nt, sym)\n",
+ "# c0_2 = calc_sim(c0_2, c_tot_2, Ga0_2, dt, nt, sym)\n",
+ "# c0_3 = calc_sim(c0_3, c_tot_3, Ga0_3, dt, nt, sym)\n",
+ "# c0_5 = calc_sim(c0_5, c_tot_5, Ga0_5, dt, nt, sym)\n",
+ "# c0_8 = calc_sim(c0_8, c_tot_8, Ga0_8, dt, nt, sym)\n",
+ "c0_9 = calc_sim(c0_9, c_tot_9, Ga0_9, dt, nt, sym)"
]
},
{
@@ -270,6 +271,32 @@
"plt.ylim(0, 1.1)"
]
},
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "### Comparison with Matlab"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "ML_9 = np.genfromtxt('Matlab_workspaces/ML_9_1.csv', delimiter=',')\n",
+ "ML_1 = np.genfromtxt('Matlab_workspaces/ML_1_9.csv', delimiter=',')\n",
+ "plt.plot(ML_9[0, :], ML_9[-1, :])\n",
+ "plt.plot(x, eval_func(c0_9, x))\n",
+ "plt.xlim(0, 0.5)\n",
+ "# plt.show()\n",
+ "\n",
+ "plt.plot(ML_1[0, :], ML_1[-1, :])\n",
+ "plt.plot(x, eval_func(c0_1, x))\n",
+ "plt.xlim(0, 0.5)\n",
+ "plt.show()"
+ ]
+ },
{
"cell_type": "markdown",
"metadata": {},
diff --git a/prob_laplace.m b/prob_laplace.m
index be32369..e288f90 100644
--- a/prob_laplace.m
+++ b/prob_laplace.m
@@ -2,7 +2,24 @@
% Define parameters for tanh, according to Weber, Zwicker, Julicher, Lee.
b = @(chi, nu) nu^(1/3)*sqrt(chi/(chi-2));
e = @(chi) sqrt(3/8*(chi-2));
-t = [0, 0.0001, 0.02, 0.05, 0.1, 1, 50];
+% t = [0, 0.0001, 0.02, 0.05, 0.1, 1, 50];
+t = linspace(0.0, 0.01, 11);
+% t=0.01;
+% t = logspace(-4, -2, 20);
+%% Comparison with Python
+T1 = Ternary_model(0, 'FRAP', {-0.1, b(7/3, 10^-12), 0.5, 0.4,...
+ 0, 1, 300, 7, 0, 'Constituent'},...
+ t, 1);
+T1.solve_tern_frap();
+csvwrite('ML_9_1.csv', [T1.x; T1.sol]);
+
+
+T2 = Ternary_model(0, 'FRAP', {-0.1, b(7/3, 10^-12), 0.9, 0,...
+ -0.895, 1, 300, 7, 0, 'Constituent'},...
+ t, 1);
+T2.x = T1.x;
+T2.solve_tern_frap();
+csvwrite('ML_1_9.csv', [T2.x; T2.sol]);
%% Try out different precisions
% prec = [0.5, 0.5, 1, 2, 3.5, 5];
fac = [1, 5, 10];%, 10, 100];
--
GitLab