Chack the flux only on Neumann sides

rafavzqz · Jul 16, 2024 · cbb8296 · cbb8296
1 parent 455f6d9
commit cbb8296
Showing 1 changed file with 3 additions and 4 deletions.
diff --git a/geopdes/inst/solve/solve_cahn_hilliard.m b/geopdes/inst/solve/solve_cahn_hilliard.m
@@ -156,7 +156,7 @@
   udot_n = zeros(space.ndof, 1);
 end
 
-flux_initial = check_flux_phase_field(space, msh, u_n, zeros(size(u_n)));
+flux_initial = check_flux_phase_field(space, msh, u_n, zeros(size(u_n)), nmnn_sides);
 disp(strcat('initial flux =',num2str(flux_initial)))
 
 %%-------------------------------------------------------------------------
@@ -182,7 +182,7 @@
                     mass_mat, lapl_mat, bnd_mat, Pen, pen_rhs, space, msh);
 
   % check flux through the boundary
-  flux = check_flux_phase_field(space, msh, u_n1, zeros(size(u_n1)));
+  flux = check_flux_phase_field(space, msh, u_n1, zeros(size(u_n1)), nmnn_sides);
   disp(strcat('flux norm =',num2str(flux)))
 
   % Time step update
@@ -374,9 +374,8 @@
 % Check flux through the boundaries
 %--------------------------------------------------------------------------
 
-function flux = check_flux_phase_field(space, msh, uhat, uhat0)
+function flux = check_flux_phase_field(space, msh, uhat, uhat0, sides)
 
-sides = [1,2,3,4]; % all the boundaries
 flux = 0;
 
 for iside = 1:numel(sides)