Metodi numerici Cheat Sheet (DRAFT) by Unobia

n = size(A,1);
D = spdiags(diag(A),0,n,n);
N = A-D;
for k = 1:kmax
    x = D\(b-N*x);
    res_rel = norm(b-A*x)/norm(b);
    if res_rel <= toll
        break
    end
end

M = tril(A);
N = A-M;
for k = 1:kmax
    x = M\(b-N*x);
    res_rel = norm(b-A*x)/norm(b);
    if res_rel <= toll
        break
    end
end

r = b-A*x;
for k = 1:kmax
    z = A*r;
    alpha = (r'r)/(r'z);
    x = x + alpha*r;
    r = r - alpha*z;
    res_rel = norm(r)/norm(b);
    if res_rel <= tol
        break
    end
end

x = linspace(a,b,n+1);
y = f(x);
h = (b-a)/n;
Itrap = (h/2)(y(1) + 2sum(y(2:end-1)) + y(end));

x = linspace(a,b,2*n+1);
y = f(x);
h = (b-a)/n;
Isimpson = (h/6)(y(1) + 4sum(y(2:2:end-1)) + 2*sum(y(3:2:end-2)) + y(end));

N = length(x) - 1;
h = x(2) - x(1);
xM = (x(1:end-1)+x(2:end))/2;
muVec = mu(xM);

d0 = 1/h^2*(muVec(1:end-1) + muVec(2:end));
d1 = -1/h^2*muVec(2:end);
d_1 = -1/h^2*muVec(1:end-1);

A = spdiags([d_1, d0, d1], [-1, 0, 1], N-1, N-1);
b = f(x(2:end-1));
b(1) = b(1) + (1/h^2)*muVec(1) *BC(1);
b(end) = b(end) + (1/h^2)*muVec(end)*BC(2);

u = A\b;
u = [BC(1); u; BC(2)];

N = length(x) - 1;
h = x(2) - x(1);
xM = (x(1:end-1)+x(2:end))/2;

d0 = mu/h^2*([1; 2*ones(N-1,1)]);
d1 = -mu/h^2*ones(N,1);
d_1 = -mu/h^2*ones(N,1);

A = spdiags([d_1, d0, d1], [-1, 0, 1], N, N);
b = f(x(1:end-1));
b(1) = 1/2*b(1)-BC(1)/h;
b(end) = b(end) + 1/h^2*mu*BC(2);

u = A\b;
u = [u; BC(2)];

N = length(x) - 1;
h = diff(x);
xM = x(1:N) + 0.5*h;
muVec = mu(xM)./h;

d = muVec(1:N-1) + muVec(2:N);
d_1 = -muVec(2:N);
d1 = -muVec(1:N-1);

A = spdiags([d_1 d d1], [-1 0 1], N-1, N-1);

b = f(x(2:N)).*(h(1:N-1) + h(2:N))*0.5;
b(1) = b(1) + muVec(1)*BC(1);
b(N-1) = b(N-1) + muVec(N-1)*BC(2);

uDof = A\b;
u = [BC(1); uDof; BC(2)];

N = length(x) - 1;
h = diff(x);
xM = x(1:N) + 0.5*h;
muVec = mu(xM)./h;

d = [muVec(1); muVec(1:N-1) + muVec(2:N)];
d_1 = -[muVec(1); muVec(2:N)];
d1 = -[muVec(1); muVec(1:N-1)];

A = spdiags([d_1 d d1], [-1 0 1], N, N);

b = f(x(2:N)).*(h(1:N-1) + h(2:N))*0.5;
b = [0.5*f(x(1))*h(1) - BC(1); b];
b(end) = b(end) + muVec(end)*BC(2);

uDof = A\b;
u = [uDof; BC(2)];

for k=1:kmax
    x = x0 - m*f(x0)/df(x0);
    errx(k) = abs(x-x0)/abs(x);
    val_f = f(x);
    if val_f<=tolf
        break
    end
    if errx(k)<=tolr
        break
    end
    x0 = x;
end

x_vec(1) = x0; 
for k=1:kmax
    x = phi(x0);
    x_vec(k+1) = x;
    errx(k) = abs(x-x0)/x;
    if errx(k)<=tolr
        break
    end
    val_f = abs(phi(x)-phi(x0));
    if val_f<=tolf
        break
    end
    x0 = x;
end

t = t0:dt:T;
u = zeros(length(t),1);
u(1) = u0;
for k = 1:length(t) - 1
    u(k+1) = u(k) + dt*f(t(k),u(k));
end
t = t';

t = t0:dt:T;
u = zeros(length(t),1);
u(1) = u0;
for k = length(t)-1
    u(k+1) = 1/(1 - dtlambda)u(k);
end
t = t';

t = t0:dt:T;
dim_u = size(u0,1);
u = zeros(length(t),dim_u);
u(1,:) = u0;
for k = 1:length(t) - 1
    u(k+1,:) = u(k,:) + dt*f(t(k),u(k,:))';
end
t = t';