% b = @(chi, nu) nu^(1/3)*sqrt(chi/(chi-2)); % e = @(chi) sqrt(3/8*(chi-2)); % t = linspace(0, 10, 100); % T = Ternary_model(0, 'FRAP', {-0.1, b(7.9/3, 10^-12), 0.5, e(7.9/3),... % 0, 1, 30, 7, 0, 'Constituent', 0},... % t, 2); % T.solve_tern_frap() load ~/Desktop/f_eos_parts_local.mat fit_seed_part = fit_seed; load ~/Desktop/f_eos.mat %% Select individual simulation jj = find(n_f==25 & fix_params(3, :)'==2); jj = jj(20); D = abs(f_temp{n_f(jj)}.conv_factor_D*f_temp{n_f(jj)}.x_seed); [cost, T] = to_min(fit_seed(jj), fit, [D; fix_params(2:3, jj)], fixed,... f_temp{n_f(jj)}); %% Analytical solution for steady state (see note from 03/11/20 in FRAP D_out = (1-T.phi_t(end))*T.ga0(end); L = sqrt(D_out)*80; R = -T.a+10*T.b; t_ind = 11; ind_R = find(T.x>R, 1, 'first'); x = T.x(ind_R:end); j = -D_out*(T.sol(t_ind, ind_R+1)-T.sol(t_ind, ind_R))/... (T.x(ind_R+1)-T.x(ind_R)); B = j*R^2/D_out; % infinite boundary A = T.phi_t(end); % finite boundary % ind_L = find(T.x>L-T.a, 1); % A = T.sol(t_ind, ind_L)-B/T.x(ind_L); c = A+B./x; T.plot_sim('plot', 10, 'green', 1); plot(x, c, 'r', 'LineWidth', 3); diff = sum(T.sol(t_ind, ind_R:end)-c)/length(c); axis([-inf, inf, -inf, inf]); j_qs = -D_out*(c(2)-c(1))/(x(2)-x(1));