Skip to content

Commit

Permalink
Changed the way to store results
Browse files Browse the repository at this point in the history
  • Loading branch information
@rafavzqz
rafavzqz committed Jul 23, 2024
1 parent b16c6a8 commit 2a9a89a
Show file tree
Hide file tree
Showing 2 changed files with 52 additions and 44 deletions.
48 changes: 26 additions & 22 deletions geopdes/inst/multipatch/mp_solve_cahn_hilliard_C1.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,8 @@
% - dt: time step size for generalized-alpha method
% - rho_inf_gen_alpha: parameter in [0,1], which governs numerical damping of the generalized alpha method
%
% save_info: time values at which the solution should be saved (see "results" in the output)
%
% OUTPUT:
%
% geometry: geometry structure (see mp_geo_load)
Expand Down Expand Up @@ -131,7 +133,7 @@

%%-------------------------------------------------------------------------
% Initial conditions
time = 0;
time = problem_data.initial_time;
if (exist('fun_u', 'var') && ~isempty(fun_u))
rhs = op_f_v_mp (space, msh, fun_u);
u_n = (mass_mat + Cpen_projection/Cpen_nitsche * Pen)\rhs;
Expand All @@ -148,16 +150,22 @@

%%-------------------------------------------------------------------------
% Initialize structure to store the results
save_id = 1;
results.u = zeros(length(u_n), length(save_info)+1);
results.udot = zeros(length(u_n), length(save_info)+1);
results.time = zeros(length(save_info)+1,1);
flag_stop_save = false;
save_info = save_info(save_info>=problem_data.initial_time & save_info<=problem_data.Time_max);

results.u = zeros(length(u_n), length(save_info));
results.udot = zeros(length(u_n), length(save_info));
results.time = zeros(length(save_info), 1);
save_info(end+1) = problem_data.Time_max + 1e5;

% Save initial conditions
results.u(:,1) = u_n;
results.udot(:,1) = udot_n;
results.time(1) = time;
save_id = 1;
if (time >= save_info(1))
results.u(:,1) = u_n;
results.udot(:,1) = udot_n;
results.time(1) = time;
save_id = 2;
save_info = save_info(save_info > time);
end

%%-------------------------------------------------------------------------
% Loop over time steps
Expand All @@ -179,16 +187,12 @@
end

% Store results
if (flag_stop_save == false)
if (time >= save_info(save_id))
save_id = save_id + 1;
results.u(:,save_id) = u_n;
results.udot(:,save_id) = udot_n;
results.time(save_id) = time;
if (save_id > length(save_info))
flag_stop_save = true;
end
end
if (time >= save_info(1))
results.u(:,save_id) = u_n;
results.udot(:,save_id) = udot_n;
results.time(save_id) = time;
save_id = save_id + 1;
save_info = save_info(save_info > time);
end
end

Expand All @@ -197,9 +201,9 @@
disp('----------------------------------------------------------------')

% Crop results
results.u = results.u(:,1:save_id);
results.udot = results.udot(:,1:save_id);
results.time = results.time(1:save_id);
results.u = results.u(:,1:save_id-1);
results.udot = results.udot(:,1:save_id-1);
results.time = results.time(1:save_id-1);

end

Expand Down
48 changes: 26 additions & 22 deletions geopdes/inst/solve/solve_cahn_hilliard.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,8 @@
% - Cpen_nitsche: penalization parameter for Nitsche's method, to impose Neumann conditions
% - Cpen_projection: penalization parameter to impose zero flux for the initial condition
%
% save_info: time values at which the solution should be saved (see "results" in the output)
%
% OUTPUT:
%
% geometry: geometry structure (see geo_load)
Expand Down Expand Up @@ -141,7 +143,7 @@

%%-------------------------------------------------------------------------
% Initial conditions
time = 0;
time = problem_data.initial_time;
if (exist('fun_u', 'var') && ~isempty(fun_u))
rhs = op_f_v_tp (space, msh, fun_u);
u_n = (mass_mat + Cpen_projection/Cpen_nitsche * Pen)\rhs;
Expand All @@ -158,16 +160,22 @@

%%-------------------------------------------------------------------------
% Initialize structure to store the results
save_id = 1;
results.u = zeros(length(u_n), length(save_info)+1);
results.udot = zeros(length(u_n), length(save_info)+1);
results.time = zeros(length(save_info)+1,1);
flag_stop_save = false;
save_info = save_info(save_info>=problem_data.initial_time & save_info<=problem_data.Time_max);

results.u = zeros(length(u_n), length(save_info));
results.udot = zeros(length(u_n), length(save_info));
results.time = zeros(length(save_info), 1);
save_info(end+1) = problem_data.Time_max + 1e5;

% Save initial conditions
results.u(:,1) = u_n;
results.udot(:,1) = udot_n;
results.time(1) = time;
save_id = 1;
if (time >= save_info(1))
results.u(:,1) = u_n;
results.udot(:,1) = udot_n;
results.time(1) = time;
save_id = 2;
save_info = save_info(save_info > time);
end

%%-------------------------------------------------------------------------
% Loop over time steps
Expand All @@ -189,16 +197,12 @@
end

% Store results
if (flag_stop_save == false)
if (time >= save_info(save_id))
save_id = save_id + 1;
results.u(:,save_id) = u_n;
results.udot(:,save_id) = udot_n;
results.time(save_id) = time;
if (save_id > length(save_info))
flag_stop_save = true;
end
end
if (time >= save_info(1))
results.u(:,save_id) = u_n;
results.udot(:,save_id) = udot_n;
results.time(save_id) = time;
save_id = save_id + 1;
save_info = save_info(save_info > time);
end
end

Expand All @@ -207,9 +211,9 @@
disp('----------------------------------------------------------------')

% Crop results
results.u = results.u(:,1:save_id);
results.udot = results.udot(:,1:save_id);
results.time = results.time(1:save_id);
results.u = results.u(:,1:save_id-1);
results.udot = results.udot(:,1:save_id-1);
results.time = results.time(1:save_id-1);

end

Expand Down

0 comments on commit 2a9a89a

Please sign in to comment.