Mercurial > repos > public > sbplib

annotate +multiblock/stitchSchemes.m @ 714:bca7a52550d1 feature/optim

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
closing feature/optim
author Ylva Rydin <ylva.rydin@telia.com>
date Thu, 08 Mar 2018 09:42:49 +0100
parents a5d4e3ced9a5
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
1 % Stitch schemes together given connection matrix and BC vector.
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
2 % schmHand - function_handle to a Scheme constructor
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
3 % order - order of accuracy
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
4 % schmParam - cell array of extra parameters sent to each Scheme stored as cell arrays
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
5 % blocks - block definitions, On whatever form the scheme expects.
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
6 % ms - grid points in each direction for each block. Ex {[10,10], [10, 20]}
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
7 % conn - connection matrix
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
8 % bound - boundary condition vector, array of structs with fields w,e,s,n
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
9 % each field with a parameter array that is sent to schm.boundary_condition
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
10 %
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
11 % Output parameters are cell arrays and cell matrices.
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
12 %
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
13 % Ex: [schms, D, H] = stitchSchemes(schmHand, order, schmParam, blocks, ms, conn, bound)
692
dd3e9a0f5032 added jacobian to stichshemes
Ylva Rydin <ylva.rydin@telia.com>
parents: 498
diff changeset
14 function [schms, D, H,penalty,J] = stitchSchemes(schmHand, order, schmParam, grids, conn, bound,timeDep)
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
15 default_arg('schmParam',[]);
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
16 default_arg('timeDep','N');
181
419ec303e97d Made some local changes to stichSchemes.
Jonatan Werpers <jonatan@werpers.com>
parents: 120
diff changeset
17 n_blocks = numel(grids);
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
18
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
19 % Creating Schemes
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
20 for i = 1:n_blocks
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
21 if isempty(schmParam);
181
419ec303e97d Made some local changes to stichSchemes.
Jonatan Werpers <jonatan@werpers.com>
parents: 120
diff changeset
22 schms{i} = schmHand(grids{i},order,[]);
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
23 elseif ~iscell(schmParam)
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
24 param = schmParam(i);
181
419ec303e97d Made some local changes to stichSchemes.
Jonatan Werpers <jonatan@werpers.com>
parents: 120
diff changeset
25 schms{i} = schmHand(grids{i},order,param);
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
26 else
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
27 param = schmParam{i};
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
28 if iscell(param)
181
419ec303e97d Made some local changes to stichSchemes.
Jonatan Werpers <jonatan@werpers.com>
parents: 120
diff changeset
29 schms{i} = schmHand(grids{i},order,param{:});
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
30 else
181
419ec303e97d Made some local changes to stichSchemes.
Jonatan Werpers <jonatan@werpers.com>
parents: 120
diff changeset
31 schms{i} = schmHand(grids{i},order,param);
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
32 end
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
33 end
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
34
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
35 % class(schmParam)
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
36 % class(ms)
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
37 % class(blocks)
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
38 % class(schmParam{i})
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
39 % class(ms)
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
40
76
5c569cbef49e Added more input parameter handling and fixed some bugs.
Jonatan Werpers <jonatan@werpers.com>
parents: 0
diff changeset
41
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
42 end
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
43
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
44
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
45 % Total norm
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
46 H = cell(n_blocks,n_blocks);
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
47 for i = 1:n_blocks
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
48 H{i,i} = schms{i}.H;
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
49 end
693
a5d4e3ced9a5 changed my mond about the jacobian
Ylva Rydin <ylva.rydin@telia.com>
parents: 692
diff changeset
50
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
51 %% Total system matrix
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
52
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
53 % Differentiation terms
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
54 D = cell(n_blocks,n_blocks);
498
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
55 penalty = cell(n_blocks,n_blocks);
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
56
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
57
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
58 for i = 1:n_blocks
498
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
59 penalty{i,i}= @(t)0;
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
60 D{i,i} = schms{i}.D;
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
61 end
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
62
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
63 % Boundary penalty terms
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
64 switch timeDep
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
65 case {'n','no','N','No'}
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
66 for i = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
67 if ~isstruct(bound{i})
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
68 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
69 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
70
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
71 fn = fieldnames(bound{i});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
72 for j = 1:length(fn)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
73 bc = bound{i}.(fn{j});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
74 if isempty(bc)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
75 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
76 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
77
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
78 [closure, ~] = schms{i}.boundary_condition(fn{j},bc{:});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
79 D{i,i} = D{i,i}+closure;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
80 end
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
81 end
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
82
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
83 % Interface penalty terms
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
84 for i = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
85 for j = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
86 intf = conn{i,j};
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
87 if isempty(intf)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
88 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
89 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
90
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
91 [uu,uv,vv,vu] = schms{i}.interface_coupling(schms{i},intf{1},schms{j},intf{2});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
92 D{i,i} = D{i,i} + uu;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
93 D{i,j} = uv;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
94 D{j,j} = D{j,j} + vv;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
95 D{j,i} = vu;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
96 end
0
48b6fb693025 Initial commit.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
97 end
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
98 case {'y','yes','Y','Yes'}
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
99 for i = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
100 if ~isstruct(bound{i})
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
101 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
102 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
103
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
104 fn = fieldnames(bound{i});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
105 for j = 1:length(fn)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
106 bc = bound{i}.(fn{j});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
107 if isempty(bc)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
108 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
109 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
110
498
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
111 [closure, penaltyi] = schms{i}.boundary_condition(fn{j},bc{:});
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
112 D{i,i} =@(t) D{i,i}(t) + closure(t);
498
324c927d8b1d chnaged sbp interfacein 1d among many things
Ylva Rydin <ylva.rydin@telia.com>
parents: 496
diff changeset
113 penalty{i,i} = @(t) penalty{i,i}(t) + penaltyi(t);
496
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
114 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
115 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
116
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
117 % Interface penalty terms
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
118 for i = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
119 for j = 1:n_blocks
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
120 intf = conn{i,j};
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
121 if isempty(intf)
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
122 continue
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
123 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
124
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
125 [uu,uv,vv,vu] = schms{i}.interface_coupling(schms{i},intf{1},schms{j},intf{2});
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
126 D{i,i} = @(t) D{i,i}(t) + uu(t);
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
127 D{i,j} = uv;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
128 D{j,j} = @(t)D{j,j}(t) + vv(t);
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
129 D{j,i} = vu;
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
130 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
131 end
437fad4a47e1 Add 1d interface for shrodinger in 1d
Ylva Rydin <ylva.rydin@telia.com>
parents: 181
diff changeset
132 end