Mercurial > repos > public > sbplib

comparison +sbp/+implementations/d2_variable_hollow_2.m @ 1305:b5907140c069 feature/poroelastic

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Clean up hollow D2 generation by using precomputed coefficients stored in mat files
author Martin Almquist <malmquist@stanford.edu>
date Mon, 13 Jul 2020 14:21:41 -0700
parents 49e3870335ef
children
comparison
equal deleted inserted replaced
1304:a38e80fdbf60 1305:b5907140c069
32 D1(m,m-1)=-1;D1(m,m)=1; 32 D1(m,m-1)=-1;D1(m,m)=1;
33 D1=D1/h; 33 D1=D1/h;
34 %Q=H*D1 + 1/2*(e_1*e_1') - 1/2*(e_m*e_m'); 34 %Q=H*D1 + 1/2*(e_1*e_1') - 1/2*(e_m*e_m');
35 35
36 36
37 M=sparse(m,m);
38 nBP = 2; 37 nBP = 2;
39 38 M = sparse(m,m);
40 scheme_width = 3; 39 coeffs = load('sbplib/+sbp/+implementations/coeffs_d2_variable_2.mat');
41 scheme_radius = (scheme_width-1)/2;
42 40
43 function D2 = D2_fun(c) 41 function D2 = D2_fun(c)
42 M_l = zeros(nBP, coeffs.nBPC);
43 M_r = zeros(nBP, coeffs.nBPC);
44 44
45 for r = 1+scheme_radius:nBP 45 for i=1:coeffs.nBPC
46 Mm1 = -c(r-1)/2 - c(r)/2; 46 M_l = M_l + coeffs.C_l{i}*c(i);
47 M0 = c(r-1)/2 + c(r) + c(r+1)/2; 47 M_r = M_r + coeffs.C_r{i}*c(m-coeffs.nBPC+i);
48 Mp1 = -c(r)/2 - c(r+1)/2;
49
50 M(r, r-1) = Mm1;
51 M(r, r) = M0;
52 M(r, r+1) = Mp1;
53 end 48 end
54 49
55 for r = m-nBP+1:m-scheme_radius 50 M(1:nBP, 1:coeffs.nBPC) = M_l;
56 Mm1 = -c(r-1)/2 - c(r)/2; 51 M(m-nBP+1:m, m-coeffs.nBPC+1:m) = M_r;
57 M0 = c(r-1)/2 + c(r) + c(r+1)/2;
58 Mp1 = -c(r)/2 - c(r+1)/2;
59 52
60 M(r, r-1) = Mm1; 53 D2 = M/h^2;
61 M(r, r) = M0;
62 M(r, r+1) = Mp1;
63 end
64
65 M(1:2,1:2)=[c(1)/2 + c(2)/2 -c(1)/2 - c(2)/2; -c(1)/2 - c(2)/2 c(1)/2 + c(2) + c(3)/2;];
66 M(m-1:m,m-1:m)=[c(m-2)/2 + c(m-1) + c(m)/2 -c(m-1)/2 - c(m)/2; -c(m-1)/2 - c(m)/2 c(m-1)/2 + c(m)/2;];
67 M=M/h;
68
69 D2=HI*(-M-c(1)*e_l*d1_l'+c(m)*e_r*d1_r');
70 end 54 end
71 D2 = @D2_fun; 55 D2 = @D2_fun;
72 end 56 end