Mercurial > repos > public > sbplib
diff +sbp/+implementations/d4_variable_6_min_boundary_points.m @ 316:203afa156f59 feature/beams
Collected boundary operators.
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Fri, 23 Sep 2016 23:10:44 +0200 |
parents | 9230c056a574 |
children | 99005a80b4c2 |
line wrap: on
line diff
--- a/+sbp/+implementations/d4_variable_6_min_boundary_points.m Fri Sep 23 22:58:45 2016 +0200 +++ b/+sbp/+implementations/d4_variable_6_min_boundary_points.m Fri Sep 23 23:10:44 2016 +0200 @@ -22,6 +22,27 @@ H(1:6,1:6) = diag([13649/43200,12013/8640,2711/4320,5359/4320,7877/8640, 43801/43200]); H(m-5:m,m-5:m) = fliplr(flipud(diag([13649/43200,12013/8640, 2711/4320,5359/4320,7877/8640,43801/43200]))); + e_1 = zeros(m,1);e_1(1)=1; + e_m = zeros(m,1);e_m(m)=1; + + S_U = [-25/12, 4, -3, 4/3, -1/4]/h; + S_1 = zeros(1,m); + S_1(1:5) = S_U; + S_m = zeros(1,m); + S_m(m-4:m) = fliplr(-S_U); + + S2_U = [0.35e2/0.12e2 -0.26e2/0.3e1 0.19e2/0.2e1 -0.14e2/0.3e1 0.11e2/0.12e2;]/h^2; + S2_1 = zeros(1,m); + S2_1(1:5) = S2_U; + S2_m = zeros(1,m); + S2_m(m-4:m) = fliplr(S2_U); + + S3_U = [-0.5e1/0.2e1 9 -12 7 -0.3e1/0.2e1;]/h^3; + S3_1 = zeros(1,m); + S3_1(1:5) = S3_U; + S3_m = zeros(1,m); + S3_m(m-4:m) = fliplr(-S3_U); + x1=0.70127127127127; @@ -39,21 +60,6 @@ D1(m-5:m,m-8:m) = flipud( fliplr(-D1(1:6,1:9))); D1 = D1/h; - e_1 = zeros(m,1);e_1(1)=1; - e_m = zeros(m,1);e_m(m)=1; - - S_U = [-25/12, 4, -3, 4/3, -1/4]/h; - S_1 = zeros(1,m); - S_1(1:5) = S_U; - S_m = zeros(1,m); - S_m(m-4:m) = fliplr(-S_U); - - - %DS=zeros(m,m); - %DS(1,1:5)=-[-25/12, 4, -3, 4/3, -1/4]; - %DS(m,m-4:m)=fliplr(-[-25/12, 4, -3, 4/3, -1/4]); - %DS=diag(c)*DS/h; - H = h*H; HI = inv(H); @@ -97,15 +103,6 @@ D2 = HI*(-M-diag(c)*e_1*S_1+diag(c)*e_m*S_m); - S2_U = [0.35e2/0.12e2 -0.26e2/0.3e1 0.19e2/0.2e1 -0.14e2/0.3e1 0.11e2/0.12e2;]/h^2; - S2_1 = zeros(1,m); - S2_1(1:5) = S2_U; - S2_m = zeros(1,m); - S2_m(m-4:m) = fliplr(S2_U); - - - - % Fourth derivative, 1th order accurate at first 8 boundary points (still % yield 5th order convergence if stable: for example u_tt=-u_xxxx @@ -134,12 +131,6 @@ M4(m-5:m,m-5:m) = flipud( fliplr( M4_U ) ); M4 = M4/h^3; - S3_U = [-0.5e1/0.2e1 9 -12 7 -0.3e1/0.2e1;]/h^3; - S3_1 = zeros(1,m); - S3_1(1:5) = S3_U; - S3_m = zeros(1,m); - S3_m(m-4:m) = fliplr(-S3_U); - D4 = HI*(M4-e_1*S3_1+e_m*S3_m + S_1'*S2_1-S_m'*S2_m); end