NMPDE/btcs.m at master · shardoolk/NMPDE · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
clear all
close all
clc
N = 10;

format long;
for p = 1:5
    x0 = 0;
    xN = 1;
     t0 = 0;
    tf = 1;
     delt = 0.0001;
M = fix((tf-t0)/delt);
   h(p) = (xN - x0)/N;

lambda = delt/(h(p)^2);

    for i = 1:N-1
    x(i) = i*h(p);
    U0(i) = cos(pi*(x(i) - 0.5));
end
A = zeros(N-1,N-1);
t = t0;
b = zeros(N-1,1);
for j = 1:N-1
    A(j,j) = 1 + 2*lambda;
end
for i = 1:N-1
    b(i) = U0(i);
end
for j = 1:N-2
    A(j+1,j) = -1*lambda;
    A(j,j+1) = -1*lambda;
end

for i = 1:M
    t = t + delt;
    for k = 1:N-1
        b(k) = U0(k) + delt*fcts(t,x(k));
    end

    Asol = A\b;
    U0 = Asol;
end


for i = 1:N-1
    exact(i) = exp(-t)*cos(pi*(x(i) - 0.5));
end


[Asol , exact']

error(p) = max(abs(exact' - Asol));

N = N*2;

end

for i = 1:p-1
    order(i) = log(error(i)/error(i+1))/log(2);
end
order

plot(x,Asol, '*', x,exact)

plot(log(h),log(error));