.

Sheet6/flux_jump.m

@@ -0,0 +1,10 @@
+function eta = flux_jump(x,u,lambda)
+N = length(x);
+eta = zeros(N,1);
+
+for j = 2:N-1
+    ul = (u(j)-u(j-1))/(x(j)-x(j-1)); %left derivative
+    ur = (u(j+1)-u(j))/(x(j+1)-x(j)); %right derivative 
+    eta(j) = abs(lambda(x(j))*ur - lambda(x(j))*ul); %error indicator
+end
+end

Sheet6/h_adapt.m

@@ -0,0 +1,10 @@
+function xnew = h_adapt(x,eta,theta) %eta is local error
+xnew = x(1);
+for i = 1:length(x)-1
+    if eta(i) > theta*max(eta) %if local error is large, inster a midpoint
+        xnew = [xnew, (x(i)+x(i+1))/2];
+    end
+    xnew = [xnew, x(i+1)]; %if local error is small, keep the element unchanged 
+end
+xnew = unique(xnew);
+end

Sheet6/r_adapt.m

@@ -0,0 +1,11 @@
+function xnew = r_adapt(x,eta_elem)
+N = length(x);
+w = abs(eta_elem) + 1e-10; %error based monitor function
+s = zeros(N,1); 
+for i = 2:N
+    s(i) = s(i-1) + w(i-1); %cumulative sum
+end
+s = s / s(end);   % normalize to [0,1]
+s_new = linspace(0,1,N)';
+xnew = interp1(s, x, s_new, 'linear');
+end