From 1c550bd2e5834874e18ec616168877a806d61a57 Mon Sep 17 00:00:00 2001
From: "lisa.pizzo" <lisa.pizzo@noreply.imsc.uni-graz.at>
Date: Sun, 25 Jan 2026 11:49:50 +0100
Subject: [PATCH] Additional Code

---
 Project/AdditionalPlotCodes.txt | 106 ++++++++++++++++++++++++++++++++
 1 file changed, 106 insertions(+)
 create mode 100644 Project/AdditionalPlotCodes.txt

diff --git a/Project/AdditionalPlotCodes.txt b/Project/AdditionalPlotCodes.txt
new file mode 100644
index 0000000..1215693
--- /dev/null
+++ b/Project/AdditionalPlotCodes.txt
@@ -0,0 +1,106 @@
+%% TASK 5: 3D plot and video 
+% Nodes and solution from Task 5
+r = nodes(1,:)';  % r-coordinate
+z = nodes(2,:)';  % z-coordinate
+u_sol = u;        % solution vector
+
+% Create a fine phi grid for rotation
+phi = linspace(0, 2*pi, 50);
+
+% Create 3D coordinates
+[Phi, R] = meshgrid(phi, r);
+Z = repmat(z, 1, length(phi));
+X = R .* cos(Phi);
+Y = R .* sin(Phi);
+
+% Interpolate u onto this grid (simple nearest neighbor)
+U = repmat(u_sol, 1, length(phi));
+
+% Plot
+figure
+surf(X, Y, Z, U, 'EdgeColor','none');
+axis equal
+xlabel('x'); ylabel('y'); zlabel('z');
+title('Axisymmetric 3D solution');
+colorbar
+view(30,30)
+
+% Assume you already have the 3D plot (X,Y,Z,U)
+figure
+h = surf(X,Y,Z,U, 'EdgeColor','none');
+axis equal
+colorbar
+xlabel('x'); ylabel('y'); zlabel('z');
+title('Axisymmetric 3D solution');
+
+v = VideoWriter('mug_rotation.avi'); % create video file
+open(v);
+
+for angle = 0:2:360
+    view(angle, 30); % rotate camera around z-axis
+    drawnow;
+    frame = getframe(gcf); % capture current frame
+    writeVideo(v, frame);
+end
+close(v);
+
+
+%% TASK 8: 9 snapshot of the solution over time: 
+
+% --- Snapshot times (seconds)
+snapshot_times = [0 50 100 150 200 250 300 350 400];
+snapshot_steps = round(snapshot_times / tau) + 1;
+
+Nsnap = length(snapshot_steps);
+snapshots = zeros(Nnodes, Nsnap);
+
+% Initial condition (t = 0)
+snapshots(:,1) = u0;
+snapshot_counter = 2;
+
+% --- Time stepping parameters
+T_end = snapshot_times(end);
+Nt = ceil(T_end / tau);
+
+% --- System matrix
+A = (1/tau)*M + K;
+
+% LU factorization (constant matrix)
+[L,U,P,Q] = lu(A);
+
+% --- Initialize solution
+u = u0;
+
+% --- Time stepping loop
+for k = 1:Nt
+    % Right-hand side
+    b = (1/tau)*M*u + F;
+
+    % Solve system
+    u_next = Q * (U \ (L \ (P * b)));
+    u = u_next;
+
+    % Store snapshots at requested times
+    if snapshot_counter <= Nsnap
+        if k+1 == snapshot_steps(snapshot_counter)
+            snapshots(:,snapshot_counter) = u;
+            snapshot_counter = snapshot_counter + 1;
+        end
+    end
+end
+
+%% --- 3x3 grid plot (overview figure)
+figure(20)
+for i = 1:Nsnap
+    subplot(3,3,i)
+    pdeplot(model, 'XYData', snapshots(:,i), 'Mesh','off');
+    axis equal
+    title(['T = ', num2str(snapshot_times(i)), ' s'])
+end
+colorbar
+sgtitle('Transient temperature evolution')
+
+% Save combined figure
+saveas(gcf, 'Task8_Transient_3x3.png');
+
+