comparison +scheme/bcSetup.m @ 1033:037f203b9bf5 feature/burgers1d

Merge with branch feature/advectioRV to utilize the +rv package
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 17 Jan 2019 10:44:12 +0100
parents b45a6dcb61ac
children 386ef449df51
comparison
equal deleted inserted replaced
854:18162a0a5bb5 1033:037f203b9bf5
1 % function [closure, S] = bcSetup(diffOp, bc)
2 % Takes a diffOp and a cell array of boundary condition definitions. 1 % Takes a diffOp and a cell array of boundary condition definitions.
3 % Each bc is a struct with the fields 2 % Each bc is a struct with the fields
4 % * type -- Type of boundary condition 3 % * type -- Type of boundary condition
5 % * boundary -- Boundary identifier 4 % * boundary -- Boundary identifier
6 % * data -- A function_handle for a function which provides boundary data.(see below) 5 % * data -- A function_handle for a function which provides boundary data.(see below)
7 % Also takes S_sign which modifies the sign of S, [-1,1] 6 % Also takes S_sign which modifies the sign of the penalty function, [-1,1]
8 % Returns a closure matrix and a forcing function S. 7 % Returns a closure matrix and a forcing function S.
9 % 8 %
10 % The boundary data function can either be a function of time or a function of time and space coordinates. 9 % The boundary data function can either be a function of time or a function of time and space coordinates.
11 % In the case where it only depends on time it should return the data as grid function for the boundary. 10 % In the case where it only depends on time it should return the data as grid function for the boundary.
12 % In the case where it also takes space coordinates the number of space coordinates should match the number of dimensions of the problem domain. 11 % In the case where it also takes space coordinates the number of space coordinates should match the number of dimensions of the problem domain.
14 function [closure, S] = bcSetup(diffOp, bcs, S_sign) 13 function [closure, S] = bcSetup(diffOp, bcs, S_sign)
15 default_arg('S_sign', 1); 14 default_arg('S_sign', 1);
16 assertType(bcs, 'cell'); 15 assertType(bcs, 'cell');
17 assert(S_sign == 1 || S_sign == -1, 'S_sign must be either 1 or -1'); 16 assert(S_sign == 1 || S_sign == -1, 'S_sign must be either 1 or -1');
18 17
19 verifyBcFormat(bcs, diffOp); 18 [closure, penalties] = scheme.bc.closureSetup(diffOp, bcs);
20 19 S = scheme.bc.forcingSetup(diffOp, penalties, bcs, S_sign);
21 % Setup storage arrays
22 closure = spzeros(size(diffOp));
23 gridData = {};
24 symbolicData = {};
25
26 % Collect closures, penalties and data
27 for i = 1:length(bcs)
28 [localClosure, penalty] = diffOp.boundary_condition(bcs{i}.boundary, bcs{i}.type);
29 closure = closure + localClosure;
30
31 [ok, isSymbolic, data] = parseData(bcs{i}, penalty, diffOp.grid);
32
33 if ~ok
34 % There was no data
35 continue
36 end
37
38 if isSymbolic
39 symbolicData{end+1} = data;
40 else
41 gridData{end+1} = data;
42 end
43 end
44
45 % Setup penalty function
46 O = spzeros(size(diffOp),1);
47 function v = S_fun(t)
48 v = O;
49 for i = 1:length(gridData)
50 v = v + gridData{i}.penalty*gridData{i}.func(t);
51 end
52
53 for i = 1:length(symbolicData)
54 v = v + symbolicData{i}.penalty*symbolicData{i}.func(t, symbolicData{i}.coords{:});
55 end
56
57 v = S_sign * v;
58 end
59 S = @S_fun;
60 end 20 end
61
62 function verifyBcFormat(bcs, diffOp)
63 for i = 1:length(bcs)
64 assertType(bcs{i}, 'struct');
65 assertStructFields(bcs{i}, {'type', 'boundary'});
66
67 if ~isfield(bcs{i}, 'data') || isempty(bcs{i}.data)
68 continue
69 end
70
71 if ~isa(bcs{i}.data, 'function_handle')
72 error('bcs{%d}.data should be a function of time or a function of time and space',i);
73 end
74
75 b = diffOp.grid.getBoundary(bcs{i}.boundary);
76
77 dim = size(b,2);
78
79 if nargin(bcs{i}.data) == 1
80 % Grid data (only function of time)
81 assertSize(bcs{i}.data(0), 1, size(b));
82 elseif nargin(bcs{i}.data) ~= 1+dim
83 error('sbplib:scheme:bcSetup:DataWrongNumberOfArguments', 'bcs{%d}.data has the wrong number of input arguments. Must be either only time or time and space.', i);
84 end
85 end
86 end
87
88 function [ok, isSymbolic, dataStruct] = parseData(bc, penalty, grid)
89 if ~isfield(bc,'data') || isempty(bc.data)
90 isSymbolic = [];
91 dataStruct = struct();
92 ok = false;
93 return
94 end
95 ok = true;
96
97 nArg = nargin(bc.data);
98
99 if nArg > 1
100 % Symbolic data
101 isSymbolic = true;
102 coord = grid.getBoundary(bc.boundary);
103 dataStruct.penalty = penalty;
104 dataStruct.func = bc.data;
105 dataStruct.coords = num2cell(coord, 1);
106 else
107 % Grid data
108 isSymbolic = false;
109 dataStruct.penalty = penalty;
110 dataStruct.func = bcs{i}.data;
111 end
112 end