Mercurial > repos > public > sbplib

comparison +multiblock/DiffOp.m @ 592:4422c4476650 feature/utux2D

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge with feature/grids
author Martin Almquist <martin.almquist@it.uu.se>
date Mon, 11 Sep 2017 14:17:15 +0200
parents b0386d2c180d
children c360bbecf260
comparison
equal deleted inserted replaced
591:39554f2de783 592:4422c4476650
118 % Splits a matrix operator into a cell-matrix of matrix operators for 118 % Splits a matrix operator into a cell-matrix of matrix operators for
119 % each grid. 119 % each grid.
120 ops = sparse2cell(op, obj.NNN); 120 ops = sparse2cell(op, obj.NNN);
121 end 121 end
122 122
123 function op = getBoundaryOperator(obj, op, boundary) 123 % Get a boundary operator specified by opName for the given boundary/BoundaryGroup
124 if iscell(boundary) 124 function op = getBoundaryOperator(obj, opName, boundary)
125 localOpName = [op '_' boundary{2}]; 125 switch class(boundary)
126 blockId = boundary{1}; 126 case 'cell'
127 localOp = obj.diffOps{blockId}.(localOpName); 127 localOpName = [opName '_' boundary{2}];
128 128 blockId = boundary{1};
129 div = {obj.blockmatrixDiv{1}, size(localOp,2)}; 129 localOp = obj.diffOps{blockId}.(localOpName);
130 blockOp = blockmatrix.zero(div); 130
131 blockOp{blockId,1} = localOp; 131 div = {obj.blockmatrixDiv{1}, size(localOp,2)};
132 op = blockmatrix.toMatrix(blockOp); 132 blockOp = blockmatrix.zero(div);
133 return 133 blockOp{blockId,1} = localOp;
134 else 134 op = blockmatrix.toMatrix(blockOp);
135 % Boundary är en sträng med en boundary group i. 135 return
136 case 'multiblock.BoundaryGroup'
137 op = sparse(size(obj.D,1),0);
138 for i = 1:length(boundary)
139 op = [op, obj.getBoundaryOperator(opName, boundary{i})];
140 end
141 otherwise
142 error('Unknown boundary indentifier')
136 end 143 end
137 end 144 end
138 145
139 % Creates the closure and penalty matrix for a given boundary condition, 146 % Creates the closure and penalty matrix for a given boundary condition,
140 % boundary -- the name of the boundary on the form {id,name} where 147 % boundary -- the name of the boundary on the form {id,name} where
141 % id is the number of a block and name is the name of a 148 % id is the number of a block and name is the name of a
142 % boundary of that block example: {1,'s'} or {3,'w'} 149 % boundary of that block example: {1,'s'} or {3,'w'}. It
150 % can also be a boundary group
143 function [closure, penalty] = boundary_condition(obj, boundary, type) 151 function [closure, penalty] = boundary_condition(obj, boundary, type)
152 switch class(boundary)
153 case 'cell'
154 [closure, penalty] = obj.singleBoundaryCondition(boundary, type);
155 case 'multiblock.BoundaryGroup'
156 [n,m] = size(obj.D);
157 closure = sparse(n,m);
158 penalty = sparse(n,0);
159 for i = 1:length(boundary)
160 [closurePart, penaltyPart] = obj.boundary_condition(boundary{i}, type);
161 closure = closure + closurePart;
162 penalty = [penalty, penaltyPart];
163 end
164 otherwise
165 error('Unknown boundary indentifier')
166 end
167
168 end
169
170 function [closure, penalty] = singleBoundaryCondition(obj, boundary, type)
144 I = boundary{1}; 171 I = boundary{1};
145 name = boundary{2}; 172 name = boundary{2};
146 173
147 % Get the closure and penaly matrices 174 % Get the closure and penaly matrices
148 [blockClosure, blockPenalty] = obj.diffOps{I}.boundary_condition(name, type); 175 [blockClosure, blockPenalty] = obj.diffOps{I}.boundary_condition(name, type);
175 end 202 end
176 end 203 end
177 end 204 end
178 205
179 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) 206 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)
180 207 error('not implemented')
181 end 208 end
182 209
183 % Size returns the number of degrees of freedom 210 % Size returns the number of degrees of freedom
184 function N = size(obj) 211 function N = size(obj)
185 N = 0; 212 N = 0;