template with coefficient matrix and rhs vector using pointer for a 2d unsteady heat conduction equation. implicit euler for time and central differencing for space. use a jacobi solver. in matlab

Here's an example of how you can template a solution for the 2D Unsteady Heat Conduction equation. Below I will represent the Heat equation discretized with implicit Euler in time and central difference in space.

Where the coefficient matrix can be represented as

Here's an example code:


main.m
%Parameters
Lx = 1;
Ly = 1;
Nx = 30;
Ny = Nx;
dx = Lx/(Nx+1);
dy = Ly/(Ny+1);
dt = 0.01;
t_end = 1;
alpha = 1;

%Initial Conditions
T = zeros(Nx,Ny);
T_old = T;

%Coefficient Matrix calculation
p = -alpha*dt/dx^2;
q = -alpha*dt/dy^2;
r = 1+2*(alpha*dt/dx^2) + 2*(alpha*dt/dy^2);
C = zeros(Nx*Ny,Nx*Ny);
RHS = zeros(Nx*Ny,1);
for i = 1:Nx
    for j = 1:Ny
        n = j+(i-1)*Ny;
        if i == 1 %South Boundary Condition
            C(n,n) = 1;
            RHS(n) = 0;
        elseif i == Nx %North Boundary Condition
            C(n,n) = 1;
            RHS(n) = 0;
        elseif j == 1 %West Boundary Condition
            C(n,n) = 1;
            RHS(n) = 0;
        elseif j == Ny %East Boundary Condition
            C(n,n) = 1;
            RHS(n) = 0;
        else %Interior Nodes
            %Central Differencing
            C(n,n) = r;
            C(n,n-1) = p;
            C(n,n+1) = p;
            C(n,n-Ny) = q;
            C(n,n+Ny) = q;
        end
    end    
end

%Time Integration
time = 0;
while time<t_end %March in time
    T_old = T;
    
    %Jacobi iteration
    for i = 1:50 %Number of iterations
        T = C\RHS;
        RHS = zeros(Nx*Ny,1); %Reset RHS vector
        for i = 1:Nx
            for j = 1:Ny
                n = j+(i-1)*Ny;
                if i == 1 %South Boundary Condition
                    RHS(n) = T_old(i,j);
                elseif i == Nx %North Boundary Condition
                    RHS(n) = T_old(i,j);
                elseif j == 1 %West Boundary Condition
                    RHS(n) = T_old(i,j);
                elseif j == Ny %East Boundary Condition
                    RHS(n) = T_old(i,j);
                else %Interior Nodes
                    RHS(n) = T_old(i,j);
                end
            end    
        end
    end
    
    time = time+dt;
end
1818 chars
77 lines

Here we have used a Jacobi Solver to solve the system of equations. Jacobi iteration method is one of the simplest methods to solve a system of linear equations. It starts with an initial guess and performs iterations to approach the solution.

This code is a good starting point but you may need to modify it depending on your boundary conditions, domain size, and other requirements.