FEAT 3
Finite Element Analysis Toolbox
Loading...
Searching...
No Matches
stokes_blocked.hpp
1// FEAT3: Finite Element Analysis Toolbox, Version 3
2// Copyright (C) 2010 by Stefan Turek & the FEAT group
3// FEAT3 is released under the GNU General Public License version 3,
4// see the file 'copyright.txt' in the top level directory for details.
5
6#pragma once
7
9#include <kernel/lafem/dense_vector.hpp>
10#include <kernel/lafem/dense_vector_blocked.hpp>
11#include <kernel/lafem/filter_chain.hpp>
12#include <kernel/lafem/filter_sequence.hpp>
13#include <kernel/lafem/sparse_matrix_csr.hpp>
14#include <kernel/lafem/sparse_matrix_bcsr.hpp>
15#include <kernel/lafem/sparse_matrix_bwrappedcsr.hpp>
16#include <kernel/lafem/saddle_point_matrix.hpp>
17#include <kernel/lafem/unit_filter.hpp>
18#include <kernel/lafem/unit_filter_blocked.hpp>
19#include <kernel/lafem/mean_filter.hpp>
20#include <kernel/lafem/none_filter.hpp>
21#include <kernel/lafem/slip_filter.hpp>
22#include <kernel/lafem/tuple_filter.hpp>
23#include <kernel/lafem/tuple_mirror.hpp>
24#include <kernel/lafem/tuple_diag_matrix.hpp>
25#include <kernel/lafem/transfer.hpp>
26#include <kernel/assembly/mirror_assembler.hpp>
27#include <kernel/assembly/symbolic_assembler.hpp>
28#include <kernel/assembly/bilinear_operator_assembler.hpp>
29#include <kernel/assembly/domain_assembler_basic_jobs.hpp>
30#include <kernel/assembly/common_operators.hpp>
31#include <kernel/assembly/gpdv_assembler.hpp>
32#include <kernel/assembly/grid_transfer.hpp>
33#include <kernel/assembly/mean_filter_assembler.hpp>
34#include <kernel/assembly/unit_filter_assembler.hpp>
35#include <kernel/global/gate.hpp>
36#include <kernel/global/muxer.hpp>
37#include <kernel/global/vector.hpp>
38#include <kernel/global/matrix.hpp>
39#include <kernel/global/filter.hpp>
40#include <kernel/global/mean_filter.hpp>
41#include <kernel/global/transfer.hpp>
42#include <kernel/global/splitter.hpp>
43
44#include <control/domain/domain_control.hpp>
45#include <control/asm/gate_asm.hpp>
46#include <control/asm/muxer_asm.hpp>
47#include <control/asm/splitter_asm.hpp>
48#include <control/asm/transfer_asm.hpp>
49#include <control/asm/transfer_voxel_asm.hpp>
50#include <control/asm/unit_filter_asm.hpp>
51#include <control/asm/mean_filter_asm.hpp>
52#include <control/asm/slip_filter_asm.hpp>
53
54
55namespace FEAT
56{
57 namespace Control
58 {
59 template
60 <
61 int dim_,
62 typename DataType_ = Real,
63 typename IndexType_ = Index,
64 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
65 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
66 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
67 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
68 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
69 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
70 typename TransferMatrixS_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
71 >
73 {
74 public:
75 // basic types
76 typedef DataType_ DataType;
77 typedef IndexType_ IndexType;
78 static constexpr int dim = dim_;
79
80 // scalar types
81 typedef ScalarMatrix_ LocalScalarMatrix;
82 typedef typename LocalScalarMatrix::VectorTypeL LocalScalarVector;
83
84 // define local blocked matrix type
85 typedef MatrixBlockA_ LocalMatrixBlockA;
86 typedef MatrixBlockB_ LocalMatrixBlockB;
87 typedef MatrixBlockD_ LocalMatrixBlockD;
88 typedef LocalScalarMatrix LocalSchurMatrix;
90
91 // define local vector types
92 typedef typename LocalMatrixBlockB::VectorTypeL LocalVeloVector;
93 typedef typename LocalMatrixBlockD::VectorTypeL LocalPresVector;
95
96 // define local transfer matrix types
97 typedef TransferMatrixV_ LocalVeloTransferMatrix;
98 typedef TransferMatrixP_ LocalPresTransferMatrix;
99 typedef TransferMatrixS_ LocalScalarTransferMatrix;
101
102 // define local transfer operators
107
108 // define mirror types
110 typedef ScalarMirror VeloMirror;
111 typedef ScalarMirror PresMirror;
113
114 // define gates
119
120 // define muxers
125
126 // define splitters
131
132 // define global vector types
136
137 // define global matrix types
143
144 // define global transfer types
149
150 /* ***************************************************************************************** */
151
154 PresGate gate_pres;
155 SystemGate gate_sys;
156 ScalarGate gate_scalar_velo;
157
160 PresMuxer coarse_muxer_pres;
161 SystemMuxer coarse_muxer_sys;
162 ScalarMuxer coarse_muxer_scalar_velo;
163
166 PresSplitter base_splitter_pres;
167 SystemSplitter base_splitter_sys;
168 ScalarSplitter base_splitter_scalar_velo;
169
172 GlobalMatrixBlockA matrix_a;
173 GlobalMatrixBlockB matrix_b;
174 GlobalMatrixBlockD matrix_d;
175 GlobalSchurMatrix matrix_s;
176
179 GlobalPresTransfer transfer_pres;
180 GlobalSystemTransfer transfer_sys;
181 GlobalScalarTransfer transfer_scalar_velo;
182
185
188 matrix_sys(&gate_sys, &gate_sys),
189 matrix_a(&gate_velo, &gate_velo),
190 matrix_b(&gate_velo, &gate_pres),
191 matrix_d(&gate_pres, &gate_velo),
192 matrix_s(&gate_pres, &gate_pres),
194 transfer_pres(&coarse_muxer_pres),
195 transfer_sys(&coarse_muxer_sys),
196 transfer_scalar_velo(&coarse_muxer_scalar_velo)
197 {
198 }
199
200 // no copies, no problems
202 StokesBlockedSystemLevel& operator=(const StokesBlockedSystemLevel&) = delete;
203
205 {
206 }
207
209 std::size_t bytes() const
210 {
211 return this->matrix_sys.local().bytes () + this->matrix_s.local().bytes() + transfer_sys.bytes();
212 }
213
214 void compile_system_matrix()
215 {
219 }
220
221 void compile_system_transfer()
222 {
223 transfer_sys.get_mat_prol().template at<0,0>() = transfer_velo.get_mat_prol().clone(LAFEM::CloneMode::Shallow);
224 transfer_sys.get_mat_rest().template at<0,0>() = transfer_velo.get_mat_rest().clone(LAFEM::CloneMode::Shallow);
225 transfer_sys.get_mat_prol().template at<1,1>() = transfer_pres.get_mat_prol().clone(LAFEM::CloneMode::Shallow);
226 transfer_sys.get_mat_rest().template at<1,1>() = transfer_pres.get_mat_rest().clone(LAFEM::CloneMode::Shallow);
227 transfer_sys.get_mat_trunc().template at<0,0>() = transfer_velo.get_mat_trunc().clone(LAFEM::CloneMode::Shallow);
228 transfer_sys.get_mat_trunc().template at<1,1>() = transfer_pres.get_mat_trunc().clone(LAFEM::CloneMode::Shallow);
229 transfer_sys.compile();
230 }
231
232 void compile_scalar_transfer()
233 {
234 transfer_scalar_velo.get_mat_prol().clone(this->transfer_velo.get_mat_prol().unwrap(), LAFEM::CloneMode::Shallow);
235 transfer_scalar_velo.get_mat_rest().clone(this->transfer_velo.get_mat_rest().unwrap(), LAFEM::CloneMode::Shallow);
236 transfer_scalar_velo.get_mat_trunc().clone(this->transfer_velo.get_mat_trunc().unwrap(), LAFEM::CloneMode::Shallow);
237 transfer_scalar_velo.compile();
238 }
239
240 template<typename D_, typename I_, typename SMA_, typename SMB_, typename SMD_, typename SM_, typename TV_, typename TP_>
241 void convert(const StokesBlockedSystemLevel<dim_, D_, I_, SMA_, SMB_, SMD_, SM_, TV_, TP_> & other)
242 {
243 gate_velo.convert(other.gate_velo);
244 gate_pres.convert(other.gate_pres);
245 Asm::build_gate_tuple(this->gate_sys, this->gate_velo, this->gate_pres);
246 this->gate_scalar_velo.convert(this->gate_velo, LocalScalarVector(this->gate_velo.get_freqs().template size<LAFEM::Perspective::native>()));
247
248 coarse_muxer_velo.convert(other.coarse_muxer_velo);
249 coarse_muxer_pres.convert(other.coarse_muxer_pres);
250 Asm::build_muxer_tuple(this->coarse_muxer_sys, this->gate_sys.get_freqs(), this->coarse_muxer_velo, this->coarse_muxer_pres);
251 this->coarse_muxer_scalar_velo.convert(this->coarse_muxer_velo, this->gate_scalar_velo.get_freqs().clone(LAFEM::CloneMode::Shallow));
252
253 base_splitter_velo.convert(other.base_splitter_velo);
254 base_splitter_pres.convert(other.base_splitter_pres);
255 Asm::build_splitter_tuple(this->base_splitter_sys, this->gate_sys.get_freqs(), this->base_splitter_velo, this->base_splitter_pres);
256 this->base_splitter_scalar_velo.convert(this->base_splitter_velo, this->gate_scalar_velo.get_freqs().clone(LAFEM::CloneMode::Shallow));
257
259 transfer_pres.convert(&coarse_muxer_pres, other.transfer_pres);
260 this->compile_system_transfer();
261 this->compile_scalar_transfer();
262
263 matrix_a.convert(&gate_velo, &gate_velo, other.matrix_a);
264 matrix_b.convert(&gate_velo, &gate_pres, other.matrix_b);
265 matrix_d.convert(&gate_pres, &gate_velo, other.matrix_d);
266 matrix_s.convert(&gate_pres, &gate_pres, other.matrix_s);
267 this->compile_system_matrix();
268 }
269
270 GlobalSystemVector create_global_vector_sys() const
271 {
272 return GlobalSystemVector(&this->gate_sys, this->gate_sys.get_freqs().clone(LAFEM::CloneMode::Layout));
273 }
274
275 GlobalVeloVector create_global_vector_velo() const
276 {
277 return GlobalVeloVector(&this->gate_velo, this->gate_velo.get_freqs().clone(LAFEM::CloneMode::Layout));
278 }
279
280 GlobalPresVector create_global_vector_pres() const
281 {
282 return GlobalPresVector(&this->gate_pres, this->gate_pres.get_freqs().clone(LAFEM::CloneMode::Layout));
283 }
284
285 template<typename DomainLevel_>
286 void assemble_gates(const Domain::VirtualLevel<DomainLevel_>& virt_dom_lvl)
287 {
288 Asm::asm_gate(virt_dom_lvl, virt_dom_lvl->space_velo, this->gate_velo, true);
289 Asm::asm_gate(virt_dom_lvl, virt_dom_lvl->space_pres, this->gate_pres, true);
290 Asm::build_gate_tuple(this->gate_sys, this->gate_velo, this->gate_pres);
291 this->gate_scalar_velo.convert(this->gate_velo, LocalScalarVector(virt_dom_lvl->space_velo.get_num_dofs()));
292 }
293
294 template<typename DomainLevel_>
295 void assemble_coarse_muxers(const Domain::VirtualLevel<DomainLevel_>& virt_lvl_coarse)
296 {
297 Asm::asm_muxer(virt_lvl_coarse, [](const DomainLevel_& dl){return &dl.space_velo;}, this->coarse_muxer_velo);
298 Asm::asm_muxer(virt_lvl_coarse, [](const DomainLevel_& dl){return &dl.space_pres;}, this->coarse_muxer_pres);
299 Asm::build_muxer_tuple(this->coarse_muxer_sys, this->gate_sys.get_freqs(), this->coarse_muxer_velo, this->coarse_muxer_pres);
300 this->coarse_muxer_scalar_velo.convert(this->coarse_muxer_velo, this->gate_scalar_velo.get_freqs().clone(LAFEM::CloneMode::Shallow));
301 }
302
303 template<typename DomainLevel_>
304 void assemble_base_splitters(const Domain::VirtualLevel<DomainLevel_>& virt_lvl)
305 {
306 Asm::asm_splitter(virt_lvl, [](const DomainLevel_& dl){return &dl.space_velo;}, this->base_splitter_velo);
307 Asm::asm_splitter(virt_lvl, [](const DomainLevel_& dl){return &dl.space_pres;}, this->base_splitter_pres);
308 Asm::build_splitter_tuple(this->base_splitter_sys, this->gate_sys.get_freqs(), this->base_splitter_velo, this->base_splitter_pres);
309 this->base_splitter_scalar_velo.convert(this->base_splitter_velo, this->gate_scalar_velo.get_freqs().clone(LAFEM::CloneMode::Shallow));
310 }
311
312 template<typename DomainLevel_>
313 void assemble_transfers(
314 const StokesBlockedSystemLevel& sys_lvl_coarse,
315 const Domain::VirtualLevel<DomainLevel_>& virt_lvl_fine,
316 const Domain::VirtualLevel<DomainLevel_>& virt_lvl_coarse,
317 const String& cubature, bool trunc_v = false, bool trunc_p = false, bool shrink = true)
318 {
319 Asm::asm_transfer_blocked(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_v, shrink,
320 [](const DomainLevel_& dl) {return &dl.space_velo;},
321 this->transfer_velo.local(), this->coarse_muxer_velo, this->gate_velo, sys_lvl_coarse.gate_velo);
322 Asm::asm_transfer_scalar(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_p, shrink,
323 [](const DomainLevel_& dl) {return &dl.space_pres;},
324 this->transfer_pres.local(), this->coarse_muxer_pres, this->gate_pres, sys_lvl_coarse.gate_pres);
325
326 this->transfer_velo.compile();
327 this->transfer_pres.compile();
328 this->compile_system_transfer();
329 this->compile_scalar_transfer();
330 }
331
332 template<typename DomainLevel_>
333 void assemble_transfers(
334 const Domain::VirtualLevel<DomainLevel_>& virt_lvl_fine,
335 const Domain::VirtualLevel<DomainLevel_>& virt_lvl_coarse,
336 const String& cubature, bool trunc_v = false, bool trunc_p = false, bool shrink = true)
337 {
338 // if the coarse level is a parent, then we need the coarse system level
339 XASSERT(!virt_lvl_coarse.is_parent());
340
341 Asm::asm_transfer_blocked(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_v, shrink,
342 [](const DomainLevel_& dl) {return &dl.space_velo;},
343 this->transfer_velo.local(), this->coarse_muxer_velo, this->gate_velo, this->gate_velo);
344 Asm::asm_transfer_scalar(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_p, shrink,
345 [](const DomainLevel_& dl) {return &dl.space_pres;},
346 this->transfer_pres.local(), this->coarse_muxer_pres, this->gate_pres, this->gate_pres);
347
348 this->transfer_velo.compile();
349 this->transfer_pres.compile();
350 this->compile_system_transfer();
351 this->compile_scalar_transfer();
352 }
353
354 template<typename DomainLevel_, template<typename> class SlagDomainLevelWrapper_>
355 void assemble_transfers_voxel(
356 const StokesBlockedSystemLevel& sys_lvl_coarse,
357 const Domain::VirtualLevel<SlagDomainLevelWrapper_<DomainLevel_>>& virt_lvl_fine,
358 const Domain::VirtualLevel<SlagDomainLevelWrapper_<DomainLevel_>>& virt_lvl_coarse,
359 const String& cubature, bool trunc_v = false, bool trunc_p = false, bool shrink = true)
360 {
361 Asm::asm_transfer_voxel_blocked(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_v, shrink,
362 [](const DomainLevel_& dl) {return &dl.space_velo;},
363 this->transfer_velo.local(), this->coarse_muxer_velo, this->gate_velo, sys_lvl_coarse.gate_velo);
364 Asm::asm_transfer_voxel_scalar(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_p, shrink,
365 [](const DomainLevel_& dl) {return &dl.space_pres;},
366 this->transfer_pres.local(), this->coarse_muxer_pres, this->gate_pres, sys_lvl_coarse.gate_pres);
367
368 this->transfer_velo.compile();
369 this->transfer_pres.compile();
370 this->compile_system_transfer();
371 this->compile_scalar_transfer();
372 }
373
374 template<typename DomainLevel_, template<typename> class SlagDomainLevelWrapper_>
375 void assemble_transfers_voxel(
376 const Domain::VirtualLevel<SlagDomainLevelWrapper_<DomainLevel_>>& virt_lvl_fine,
377 const Domain::VirtualLevel<SlagDomainLevelWrapper_<DomainLevel_>>& virt_lvl_coarse,
378 const String& cubature, bool trunc_v = false, bool trunc_p = false, bool shrink = true)
379 {
380 // if the coarse level is a parent, then we need the coarse system level
381 XASSERT(!virt_lvl_coarse.is_parent());
382
383 Asm::asm_transfer_voxel_blocked(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_v, shrink,
384 [](const DomainLevel_& dl) {return &dl.space_velo;},
385 this->transfer_velo.local(), this->coarse_muxer_velo, this->gate_velo, this->gate_velo);
386 Asm::asm_transfer_voxel_scalar(virt_lvl_fine, virt_lvl_coarse, cubature, trunc_p, shrink,
387 [](const DomainLevel_& dl) {return &dl.space_pres;},
388 this->transfer_pres.local(), this->coarse_muxer_pres, this->gate_pres, this->gate_pres);
389
390 this->transfer_velo.compile();
391 this->transfer_pres.compile();
392 this->compile_system_transfer();
393 this->compile_scalar_transfer();
394 }
395
396 void clear_transfers()
397 {
398 this->transfer_velo.local() = LocalVeloTransfer();
399 this->transfer_pres.local() = LocalPresTransfer();
400 this->transfer_scalar_velo.local() = LocalScalarTransfer();
401 this->compile_system_transfer();
402 this->compile_scalar_transfer();
403 }
404
405 template<typename SpaceVelo_, typename SpacePres_, typename Cubature_>
406 void assemble_grad_div_matrices(const SpaceVelo_& space_velo, const SpacePres_& space_pres, const Cubature_& cubature)
407 {
408 Assembly::GradPresDivVeloAssembler::assemble(this->matrix_b.local(), this->matrix_d.local(), space_velo, space_pres, cubature);
409 }
410
411 template<typename Trafo_, typename SpaceVelo_, typename SpacePres_>
412 void assemble_grad_div_matrices(Assembly::DomainAssembler<Trafo_>& dom_asm,
413 const SpaceVelo_& space_velo, const SpacePres_& space_pres, const String& cubature)
414 {
415 // assemble matrix structure of B
416 if(this->matrix_b.local().empty())
417 Assembly::SymbolicAssembler::assemble_matrix_std2(this->matrix_b.local(), space_velo, space_pres);
418
419 // assemble matrix B
420 this->matrix_b.local().format();
421 Assembly::Common::GradientTestOperatorBlocked<dim_> grad_op;
423 dom_asm, this->matrix_b.local(), grad_op, space_velo, space_pres, cubature, -DataType(1));
424
425 // transpose to obtain matrix D
426 this->matrix_d.local().transpose(this->matrix_b.local());
427 }
428
429 template<typename Trafo_, typename SpaceVelo_, typename SpacePres_, typename AsmDT_>
430 void assemble_grad_div_matrices_high_prec(Assembly::DomainAssembler<Trafo_>& dom_asm,
431 const SpaceVelo_& space_velo, const SpacePres_& space_pres, const String& cubature, const AsmDT_ asm_weight)
432 {
433 // assemble matrix structure of B
434 if(this->matrix_b.local().empty())
435 Assembly::SymbolicAssembler::assemble_matrix_std2(this->matrix_b.local(), space_velo, space_pres);
436
437 // assemble matrix B
438 this->matrix_b.local().format();
439 {
440 // create local matrix with asm datatype
441 typename LocalMatrixBlockB::template ContainerTypeByDI<AsmDT_, IndexType_> tmp_mat_b;
442 tmp_mat_b.convert(this->matrix_b.local());
443 tmp_mat_b.format();
444
445 Assembly::Common::GradientTestOperatorBlocked<dim_> grad_op;
447 dom_asm, tmp_mat_b, grad_op, space_velo, space_pres, cubature, -asm_weight);
448
449 // convert matrix
450 this->matrix_b.local().convert(tmp_mat_b);
451 }
452
453 // transpose to obtain matrix D
454 this->matrix_d.local().transpose(this->matrix_b.local());
455 }
456
457 template<typename SpaceVelo_>
458 void assemble_velo_struct(const SpaceVelo_& space_velo)
459 {
460 // assemble matrix structure
461 Assembly::SymbolicAssembler::assemble_matrix_std1(this->matrix_a.local(), space_velo);
462 this->matrix_a.local().format();
463 }
464
465 template<typename SpacePres_>
466 void assemble_pres_struct(const SpacePres_& space_pres)
467 {
468 // assemble matrix structure
469 Assembly::SymbolicAssembler::assemble_matrix_std1(this->matrix_s.local(), space_pres);
470 this->matrix_s.local().format();
471 }
472
473 void compile_local_matrix_sys_type1()
474 {
475 this->local_matrix_sys_type1.block_a() = matrix_a.convert_to_1();
476 this->local_matrix_sys_type1.block_b() = matrix_b.local().clone(LAFEM::CloneMode::Weak);
477 this->local_matrix_sys_type1.block_d() = matrix_d.local().clone(LAFEM::CloneMode::Weak);
478 }
479 }; // class StokesBlockedSystemLevel<...>
480
481 template
482 <
483 int dim_,
484 typename DataType_ = Real,
485 typename IndexType_ = Index,
486 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
487 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
488 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
489 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
490 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
491 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
492 >
494 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
495 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
496 {
497 public:
498 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
499 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
500
501 // define local filter types
505
506 // define global filter types
510
511 // (global) filters
512 GlobalSystemFilter filter_sys;
513 GlobalVeloFilter filter_velo;
514 GlobalPresFilter filter_pres;
515
517 std::size_t bytes() const
518 {
519 return this->filter_sys.bytes() + BaseClass::bytes();
520 }
521
522 void compile_system_filter()
523 {
524 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
525 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
526 }
527
528 void compile_local_matrix_sys_type1()
529 {
530 BaseClass::compile_local_matrix_sys_type1();
531
532 // apply filter to A and B
533 this->filter_velo.local().filter_mat(this->local_matrix_sys_type1.block_a());
534 this->filter_velo.local().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
535 }
536 }; // class StokesBlockedUnitVeloNonePresSystemLevel
537
538 template
539 <
540 int dim_,
541 typename DataType_ = Real,
542 typename IndexType_ = Index,
543 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
544 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
545 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
546 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
547 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
548 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
549 >
551 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
552 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
553 {
554 public:
555 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
556 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
557
558 // define local filter types
564
565 // define global filter types
569
570 // (global) filters
571 GlobalSystemFilter filter_sys;
572 GlobalVeloFilter filter_velo;
573 GlobalPresFilter filter_pres;
574
576 std::size_t bytes() const
577 {
578 return this->filter_sys.bytes() + BaseClass::bytes();
579 }
580
581 void compile_system_filter()
582 {
583 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
584 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
585 }
586
587 void compile_local_matrix_sys_type1()
588 {
589 BaseClass::compile_local_matrix_sys_type1();
590
591 // apply filter to A and B
592 this->filter_velo.local().template at<1>().filter_mat(this->local_matrix_sys_type1.block_a());
593 this->filter_velo.local().template at<1>().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
594 }
595 }; // class StokesBlockedSlipUnitVeloNonePresSystemLevel
596
602 template
603 <
604 int dim_,
605 typename DataType_ = Real,
606 typename IndexType_ = Index,
607 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
608 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
609 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
610 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
611 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
612 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
613 >
615 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
616 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
617 {
618 public:
619 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
620 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
621
622 // define local filter types
626
627 // define global filter types
631
632 // (global) filters
633 GlobalSystemFilter filter_sys;
634 GlobalVeloFilter filter_velo;
635 GlobalPresFilter filter_pres;
636
638 std::size_t bytes() const
639 {
640 return this->matrix_sys.local().bytes () + this->matrix_s.local().bytes() + filter_sys.local().bytes();
641 }
642
643 void compile_system_filter()
644 {
645 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
646 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
647 }
648
649 template<typename D_, typename I_, typename SM_>
650 void convert(const StokesBlockedUnitVeloMeanPresSystemLevel<dim_, D_, I_, SM_> & other)
651 {
652 BaseClass::convert(other);
653 filter_velo.convert(other.filter_velo);
654 filter_pres.convert(other.filter_pres);
655
656 compile_system_filter();
657 }
658
659 template<typename SpacePres_, typename Cubature_>
660 void assemble_pressure_mean_filter(const SpacePres_& space_pres, const Cubature_& cubature)
661 {
662 // get our local pressure filter
663 LocalPresFilter& fil_loc_p = this->filter_pres.local();
664
665 // create two global vectors
666 typename BaseClass::GlobalPresVector vec_glob_v(&this->gate_pres), vec_glob_w(&this->gate_pres);
667
668 // fetch the local vectors
669 typename BaseClass::LocalPresVector& vec_loc_v = vec_glob_v.local();
670 typename BaseClass::LocalPresVector& vec_loc_w = vec_glob_w.local();
671
672 // fetch the frequency vector of the pressure gate
673 typename BaseClass::LocalPresVector& vec_loc_f = this->gate_pres._freqs;
674
675 // assemble the mean filter
676 Assembly::MeanFilterAssembler::assemble(vec_loc_v, vec_loc_w, space_pres, cubature);
677
678 // synchronize the vectors
679 vec_glob_v.sync_1();
680 vec_glob_w.sync_0();
681
682 // build the mean filter
683 fil_loc_p = LocalPresFilter(vec_loc_v.clone(), vec_loc_w.clone(), vec_loc_f.clone(), this->gate_pres.get_comm());
684 }
685
686 void compile_local_matrix_sys_type1()
687 {
688 BaseClass::compile_local_matrix_sys_type1();
689
690 // apply filter to A and B
691 this->filter_velo.local().filter_mat(this->local_matrix_sys_type1.block_a());
692 this->filter_velo.local().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
693 }
694 }; // struct StokesBlockedUnitVeloMeanPresSystemLevel<...>
695
696 template
697 <
698 int dim_,
699 typename DataType_ = Real,
700 typename IndexType_ = Index,
701 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
702 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
703 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
704 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
705 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
706 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
707 >
709 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
710 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
711 {
712 public:
713 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
714 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
715
716 // define local filter types
722
723 // define global filter types
727
728 // (global) filters
729 GlobalSystemFilter filter_sys;
730 GlobalVeloFilter filter_velo;
731 GlobalPresFilter filter_pres;
732
734 std::size_t bytes() const
735 {
736 return this->filter_sys.bytes() + BaseClass::bytes();
737 }
738
739 void compile_system_filter()
740 {
741 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
742 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
743 }
744
745 template<typename D_, typename I_, typename SM_>
746 void convert(const StokesBlockedUnitVeloMeanPresSystemLevel<dim_, D_, I_, SM_> & other)
747 {
748 BaseClass::convert(other);
749 filter_velo.convert(other.filter_velo);
750 filter_pres.convert(other.filter_pres);
751
752 compile_system_filter();
753 }
754
755 template<typename SpacePres_, typename Cubature_>
756 void assemble_global_filters(const SpacePres_& space_pres, const Cubature_& cubature)
757 {
758 // get our local pressure filter
759 LocalPresFilter& fil_loc_p = this->filter_pres.local();
760
761 // create two global vectors
762 typename BaseClass::GlobalPresVector vec_glob_v(&this->gate_pres), vec_glob_w(&this->gate_pres);
763
764 // fetch the local vectors
765 typename BaseClass::LocalPresVector& vec_loc_v = vec_glob_v.local();
766 typename BaseClass::LocalPresVector& vec_loc_w = vec_glob_w.local();
767
768 // fetch the frequency vector of the pressure gate
769 typename BaseClass::LocalPresVector& vec_loc_f = this->gate_pres._freqs;
770
771 // assemble the mean filter
772 Assembly::MeanFilterAssembler::assemble(vec_loc_v, vec_loc_w, space_pres, cubature);
773
774 // synchronize the vectors
775 vec_glob_v.sync_1();
776 vec_glob_w.sync_0();
777
778 // build the mean filter
779 fil_loc_p = LocalPresFilter(vec_loc_v.clone(), vec_loc_w.clone(), vec_loc_f.clone(), this->gate_pres.get_comm());
780
781 // synchronize the slip filter
782 Asm::sync_slip_filter(this->gate_velo, filter_velo.local().template at<0>());
783 }
784
785 void compile_local_matrix_sys_type1()
786 {
787 BaseClass::compile_local_matrix_sys_type1();
788
789 // apply filter to A and B
790 this->filter_velo.local().template at<1>().filter_mat(this->local_matrix_sys_type1.block_a());
791 this->filter_velo.local().template at<1>().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
792 }
793 }; // class StokesBlockedSlipUnitVeloMeanPresSystemLevel
794
800 template
801 <
802 int dim_,
803 typename DataType_ = Real,
804 typename IndexType_ = Index,
805 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
806 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
807 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
808 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
809 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
810 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
811 >
813 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
814 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
815 {
816 public:
817 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
818 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
819
820 // define local velocity filter chain
826
827 // define local pressure filter chain
831
832 // define local system filter
834
835 // define global filter types
839
840 // (global) filters
841 GlobalSystemFilter filter_sys;
842 GlobalVeloFilter filter_velo;
843 GlobalPresFilter filter_pres;
844
846 std::size_t bytes() const
847 {
848 return this->filter_sys.bytes() + BaseClass::bytes();
849 }
850
851 void compile_system_filter()
852 {
853 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
854 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
855 }
856
857 template<typename D_, typename I_, typename SM_>
858 void convert(const StokesBlockedCombinedSystemLevel<dim_, D_, I_, SM_> & other)
859 {
860 BaseClass::convert(other);
861 filter_velo.convert(other.filter_velo);
862 filter_pres.convert(other.filter_pres);
863
864 compile_system_filter();
865 }
866
867 LocalVeloSlipFilterSeq& get_local_velo_slip_filter_seq()
868 {
869 return this->filter_velo.local().template at<0>();
870 }
871
872 LocalVeloUnitFilterSeq& get_local_velo_unit_filter_seq()
873 {
874 return this->filter_velo.local().template at<1>();
875 }
876
877 LocalPresMeanFilter& get_local_pres_mean_filter()
878 {
879 return this->filter_pres.local().template at<0>();
880 }
881
882 LocalPresUnitFilter& get_local_pres_unit_filter()
883 {
884 return this->filter_pres.local().template at<1>();
885 }
886
887 template<typename DomainLevel_, typename Space_>
888 void assemble_velocity_unit_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts)
889 {
890 auto& loc_filter = get_local_velo_unit_filter_seq().find_or_add(name);
891 loc_filter.clear();
892 Asm::asm_unit_filter_blocked_homogeneous(loc_filter, dom_level, space, mesh_parts);
893 }
894
895 template<typename DomainLevel_, typename Space_, typename Function_>
896 void assemble_velocity_unit_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts, const Function_& function)
897 {
898 auto& loc_filter = get_local_velo_unit_filter_seq().find_or_add(name);
899 loc_filter.clear();
900 Asm::asm_unit_filter_blocked(loc_filter, dom_level, space, mesh_parts, function);
901 }
902
903 template<typename DomainLevel_, typename Space_>
904 void assemble_velocity_slip_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts)
905 {
906 auto& loc_filter = get_local_velo_slip_filter_seq().find_or_add(name);
907 loc_filter.clear();
908 Asm::asm_slip_filter(loc_filter, dom_level, space, mesh_parts);
909 Asm::sync_slip_filter(this->gate_velo, loc_filter);
910 }
911
912 template<typename SpacePres_>
913 void assemble_pressure_mean_filter(const SpacePres_& space_pres, const String& cubature_name)
914 {
915 this->get_local_pres_mean_filter() = Asm::asm_mean_filter(this->gate_pres, space_pres, cubature_name);
916 }
917
918 void sync_velocity_slip_filters()
919 {
920 for(auto& it : get_local_velo_slip_filter_seq())
921 Asm::sync_slip_filter(this->gate_velo, it.second);
922 }
923
924 void compile_local_matrix_sys_type1()
925 {
926 BaseClass::compile_local_matrix_sys_type1();
927 for(const auto& loc_fil : get_local_velo_unit_filter_seq())
928 {
929 loc_fil.second.filter_mat(this->local_matrix_sys_type1.block_a());
930 loc_fil.second.filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
931 }
932 }
933 }; // class StokesBlockedCombinedSystemLevel<...>
934 } // namespace Control
935} // namespace FEAT
#define XASSERT(expr)
Assertion macro definition.
Definition: assertion.hpp:262
FEAT Kernel base header.
static void assemble(LAFEM::SparseMatrixBCSR< DataType_, IndexType_, dim_, 1 > &matrix_b, LAFEM::SparseMatrixBCSR< DataType_, IndexType_, 1, dim_ > &matrix_d, const SpaceVelo_ &space_velo, const SpacePres_ &space_pres, const String &cubature_name, const DataType_ scale_b=-DataType_(1), const DataType_ scale_d=-DataType_(1))
Assembles the B and D matrices.
static void assemble(LAFEM::DenseVector< DataType_, IndexType_ > &vec_prim, LAFEM::DenseVector< DataType_, IndexType_ > &vec_dual, const Space_ &space, const String &cubature_name)
Assembles an integral mean filter.
static void assemble_matrix_std1(MatrixType_ &matrix, const Space_ &space)
Assembles a standard matrix structure from a single space.
static void assemble_matrix_std2(MatrixType_ &matrix, const TestSpace_ &test_space, const TrialSpace_ &trial_space)
Assembles a standard matrix structure from a test-trial-space pair.
Stokes blocked System level combining all supported types of boundary conditions.
std::size_t bytes() const
Returns the total amount of bytes allocated.
std::size_t bytes() const
Returns the total amount of bytes allocated.
std::size_t bytes() const
Returns the total amount of bytes allocated.
GlobalSystemMatrix matrix_sys
our global system matrix
VeloSplitter base_splitter_velo
our base-mesh multiplexer
LocalSystemMatrix local_matrix_sys_type1
local type-1 system matrix
std::size_t bytes() const
Returns the total amount of bytes allocated.
GlobalVeloTransfer transfer_velo
our global transfer operator
VeloMuxer coarse_muxer_velo
our coarse-level system muxer
std::size_t bytes() const
Returns the total amount of bytes allocated.
Global Filter wrapper class template.
Definition: filter.hpp:21
std::size_t bytes() const
Returns the total amount of bytes allocated.
Definition: filter.hpp:77
Global gate implementation.
Definition: gate.hpp:51
void convert(const Gate< LVT2_, MT2_ > &other)
Conversion function for same vector container type but with different MDI-Type.
Definition: gate.hpp:191
const LocalVector_ & get_freqs() const
Returns a const reference to the frequencies vector.
Definition: gate.hpp:179
LocalVector_ _freqs
frequency vector
Definition: gate.hpp:75
Global Matrix wrapper class template.
Definition: matrix.hpp:40
LocalMatrix_ convert_to_1() const
Computes and returns the type-1 conversion of this matrix as a local matrix.
Definition: matrix.hpp:617
LocalMatrix_ & local()
Returns a reference to the internal local LAFEM matrix object.
Definition: matrix.hpp:126
Mean Filter class template.
Definition: mean_filter.hpp:23
MeanFilter clone(LAFEM::CloneMode clone_mode=LAFEM::CloneMode::Deep) const
Creates a clone of itself.
Global multiplexer/demultiplexer implementation.
Definition: muxer.hpp:68
void convert(const Muxer< LVT2_, MT2_ > &other)
Conversion function for same vector container type but with different MDI-Type.
Definition: muxer.hpp:155
Global base-mesh vector splitter (and joiner) implementation.
Definition: splitter.hpp:59
void convert(const Splitter< LVT2_, MT2_ > &other)
Conversion function for same vector container type but with different MDI-Type.
Definition: splitter.hpp:124
Global grid-transfer operator class template.
Definition: transfer.hpp:23
void convert(MuxerType *coarse_muxer, const Transfer< LocalTransfer2_, Mirror2_ > &other)
container conversion function
Definition: transfer.hpp:107
Transfer clone(LAFEM::CloneMode clone_mode=LAFEM::CloneMode::Weak) const
Creates a clone of this object.
Definition: transfer.hpp:125
std::size_t bytes() const
Definition: transfer.hpp:143
LocalTransfer_ & local()
Definition: transfer.hpp:131
Global vector wrapper class template.
Definition: vector.hpp:68
Filter Chainclass template.
FilterChain clone(CloneMode clone_mode=CloneMode::Deep) const
Creates a clone of itself.
void filter_mat(Matrix_ &matrix) const
Applies the filter onto a system matrix.
Sequence of filters of the same type.
None Filter class template.
Definition: none_filter.hpp:29
NoneFilter clone(CloneMode=CloneMode::Deep) const
Creates a (empty) clone of itself.
Definition: none_filter.hpp:51
Saddle-Point matrix meta class template.
MatrixTypeA & block_a()
Returns the sub-matrix block A.
std::size_t bytes() const
Returns the total amount of bytes allocated.
MatrixTypeB & block_b()
Returns the sub-matrix block B.
MatrixTypeD & block_d()
Returns the sub-matrix block D.
Slip Filter class template.
Definition: slip_filter.hpp:67
Grid-Transfer operator class template.
Definition: transfer.hpp:27
Tuple-Diag-Matrix meta class template.
TupleVector meta-filter class template.
std::size_t bytes() const
Returns the total amount of bytes allocated.
TupleVector meta-mirror class template.
Variadic TupleVector class template.
TupleVector clone(LAFEM::CloneMode mode=LAFEM::CloneMode::Weak) const
Creates and returns a copy of this vector.
Unit Filter Blocked class template.
void filter_offdiag_row_mat(MatrixType &matrix) const
Filter the non-diagonal row entries.
UnitFilterBlocked clone(CloneMode clone_mode=CloneMode::Deep) const
Creates a clone of itself.
void filter_mat(MatrixType &matrix) const
Applies the filter onto a system matrix.
Unit Filter class template.
Definition: unit_filter.hpp:29
Handles vector prolongation, restriction and serialization.
void assemble_bilinear_operator_matrix_2(DomainAssembler< Trafo_ > &dom_asm, Matrix_ &matrix, const BilOp_ &bilinear_operator, const TestSpace_ &test_space, const TrialSpace_ &trial_space, const String &cubature, const typename Matrix_::DataType alpha=typename Matrix_::DataType(1))
Assembles a bilinear operator into a matrix with different test- and trial-spaces.
@ other
generic/other permutation strategy
FEAT namespace.
Definition: adjactor.hpp:12
double Real
Real data type.
std::uint64_t Index
Index data type.
System level using a MeanFilter for the pressure.
std::size_t bytes() const
Returns the total amount of bytes allocated.