Mercurial > repos > public > sbplib

comparison +time/SBPInTimeSecondOrderFormImplicit.m @ 797:5cf9fdf4c98f feature/poroelastic

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge with feature/grids and bugfix bcSetup
author Martin Almquist <malmquist@stanford.edu>
date Thu, 26 Jul 2018 10:53:05 -0700
parents 66eb4a2bbb72
children 8894e9c49e40
comparison
equal deleted inserted replaced
796:aa1ed37a1b56 797:5cf9fdf4c98f
12 12
13 methods 13 methods
14 % Solves A*u_tt + B*u_t + C*u = f(t) 14 % Solves A*u_tt + B*u_t + C*u = f(t)
15 % A, B can either both be constants or both be function handles, 15 % A, B can either both be constants or both be function handles,
16 % They can also be omitted by setting them equal to the empty matrix. 16 % They can also be omitted by setting them equal to the empty matrix.
17 function obj = SBPInTimeSecondOrderFormImplicit(A, B, C, f, k, t0, v0, v0t, TYPE, order, blockSize) 17 function obj = SBPInTimeSecondOrderFormImplicit(A, B, C, f, k, t0, v0, v0t, do_scaling, TYPE, order, blockSize)
18 default_arg('f', []); 18 default_arg('f', []);
19 default_arg('TYPE', []); 19 default_arg('TYPE', []);
20 default_arg('order', []); 20 default_arg('order', []);
21 default_arg('blockSize',[]); 21 default_arg('blockSize',[]);
22 default_arg('do_scaling', false);
22 23
23 m = length(v0); 24 m = length(v0);
24 25
25 default_arg('A', sparse(m, m)); 26 default_arg('A', speye(m, m));
26 default_arg('B', sparse(m, m)); 27 default_arg('B', sparse(m, m));
27 default_arg('C', sparse(m, m)); 28 default_arg('C', sparse(m, m));
28 29
29 I = speye(m); 30 I = speye(m);
30 O = sparse(m,m); 31 O = sparse(m,m);
54 55
55 obj.k = k; 56 obj.k = k;
56 obj.t = t0; 57 obj.t = t0;
57 obj.n = 0; 58 obj.n = 0;
58 59
59 obj.firstOrderTimeStepper = time.SBPInTimeImplicitFormulation(obj.AA, obj.BB, obj.ff, obj.k, obj.t, w0, TYPE, order, blockSize); 60 if do_scaling
61 scaling = [ones(m,1); sqrt(diag(C))];
62 obj.firstOrderTimeStepper = time.SBPInTimeScaled(obj.AA, obj.BB, obj.ff, obj.k, obj.t, w0, scaling, TYPE, order, blockSize);
63 else
64 obj.firstOrderTimeStepper = time.SBPInTimeImplicitFormulation(obj.AA, obj.BB, obj.ff, obj.k, obj.t, w0, TYPE, order, blockSize);
65 end
60 end 66 end
61 67
62 function [v,t] = getV(obj) 68 function [v,t] = getV(obj)
63 w = obj.firstOrderTimeStepper.getV(); 69 w = obj.firstOrderTimeStepper.getV();
64 v = w(1:end/2); 70 v = w(1:end/2);