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_>
355 void assemble_transfers_voxel(
356 const StokesBlockedSystemLevel& sys_lvl_coarse,
357 const Domain::VirtualLevel<Domain::VoxelDomainLevelWrapper<DomainLevel_>>& virt_lvl_fine,
358 const Domain::VirtualLevel<Domain::VoxelDomainLevelWrapper<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_>
375 void assemble_transfers_voxel(
376 const Domain::VirtualLevel<Domain::VoxelDomainLevelWrapper<DomainLevel_>>& virt_lvl_fine,
377 const Domain::VirtualLevel<Domain::VoxelDomainLevelWrapper<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 template<typename SpaceVelo_, typename SpacePres_, typename Cubature_>
397 void assemble_grad_div_matrices(const SpaceVelo_& space_velo, const SpacePres_& space_pres, const Cubature_& cubature)
398 {
399 Assembly::GradPresDivVeloAssembler::assemble(this->matrix_b.local(), this->matrix_d.local(), space_velo, space_pres, cubature);
400 }
401
402 template<typename Trafo_, typename SpaceVelo_, typename SpacePres_>
403 void assemble_grad_div_matrices(Assembly::DomainAssembler<Trafo_>& dom_asm,
404 const SpaceVelo_& space_velo, const SpacePres_& space_pres, const String& cubature)
405 {
406 // assemble matrix structure of B
407 if(this->matrix_b.local().empty())
408 Assembly::SymbolicAssembler::assemble_matrix_std2(this->matrix_b.local(), space_velo, space_pres);
409
410 // assemble matrix B
411 this->matrix_b.local().format();
412 Assembly::Common::GradientTestOperatorBlocked<dim_> grad_op;
414 dom_asm, this->matrix_b.local(), grad_op, space_velo, space_pres, cubature, -DataType(1));
415
416 // transpose to obtain matrix D
417 this->matrix_d.local().transpose(this->matrix_b.local());
418 }
419
420 template<typename Trafo_, typename SpaceVelo_, typename SpacePres_, typename AsmDT_>
421 void assemble_grad_div_matrices_high_prec(Assembly::DomainAssembler<Trafo_>& dom_asm,
422 const SpaceVelo_& space_velo, const SpacePres_& space_pres, const String& cubature, const AsmDT_ asm_weight)
423 {
424 // assemble matrix structure of B
425 if(this->matrix_b.local().empty())
426 Assembly::SymbolicAssembler::assemble_matrix_std2(this->matrix_b.local(), space_velo, space_pres);
427
428 // assemble matrix B
429 this->matrix_b.local().format();
430 {
431 // create local matrix with asm datatype
432 typename LocalMatrixBlockB::template ContainerTypeByDI<AsmDT_, IndexType_> tmp_mat_b;
433 tmp_mat_b.convert(this->matrix_b.local());
434 tmp_mat_b.format();
435
436 Assembly::Common::GradientTestOperatorBlocked<dim_> grad_op;
438 dom_asm, tmp_mat_b, grad_op, space_velo, space_pres, cubature, -asm_weight);
439
440 // convert matrix
441 this->matrix_b.local().convert(tmp_mat_b);
442 }
443
444 // transpose to obtain matrix D
445 this->matrix_d.local().transpose(this->matrix_b.local());
446 }
447
448 template<typename SpaceVelo_>
449 void assemble_velo_struct(const SpaceVelo_& space_velo)
450 {
451 // assemble matrix structure
452 Assembly::SymbolicAssembler::assemble_matrix_std1(this->matrix_a.local(), space_velo);
453 this->matrix_a.local().format();
454 }
455
456 template<typename SpacePres_>
457 void assemble_pres_struct(const SpacePres_& space_pres)
458 {
459 // assemble matrix structure
460 Assembly::SymbolicAssembler::assemble_matrix_std1(this->matrix_s.local(), space_pres);
461 this->matrix_s.local().format();
462 }
463
464 void compile_local_matrix_sys_type1()
465 {
466 this->local_matrix_sys_type1.block_a() = matrix_a.convert_to_1();
467 this->local_matrix_sys_type1.block_b() = matrix_b.local().clone(LAFEM::CloneMode::Weak);
468 this->local_matrix_sys_type1.block_d() = matrix_d.local().clone(LAFEM::CloneMode::Weak);
469 }
470 }; // class StokesBlockedSystemLevel<...>
471
472 template
473 <
474 int dim_,
475 typename DataType_ = Real,
476 typename IndexType_ = Index,
477 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
478 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
479 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
480 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
481 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
482 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
483 >
485 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
486 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
487 {
488 public:
489 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
490 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
491
492 // define local filter types
496
497 // define global filter types
501
502 // (global) filters
503 GlobalSystemFilter filter_sys;
504 GlobalVeloFilter filter_velo;
505 GlobalPresFilter filter_pres;
506
508 std::size_t bytes() const
509 {
510 return this->filter_sys.bytes() + BaseClass::bytes();
511 }
512
513 void compile_system_filter()
514 {
515 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
516 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
517 }
518
519 void compile_local_matrix_sys_type1()
520 {
521 BaseClass::compile_local_matrix_sys_type1();
522
523 // apply filter to A and B
524 this->filter_velo.local().filter_mat(this->local_matrix_sys_type1.block_a());
525 this->filter_velo.local().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
526 }
527 }; // class StokesBlockedUnitVeloNonePresSystemLevel
528
529 template
530 <
531 int dim_,
532 typename DataType_ = Real,
533 typename IndexType_ = Index,
534 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
535 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
536 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
537 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
538 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
539 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
540 >
542 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
543 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
544 {
545 public:
546 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
547 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
548
549 // define local filter types
555
556 // define global filter types
560
561 // (global) filters
562 GlobalSystemFilter filter_sys;
563 GlobalVeloFilter filter_velo;
564 GlobalPresFilter filter_pres;
565
567 std::size_t bytes() const
568 {
569 return this->filter_sys.bytes() + BaseClass::bytes();
570 }
571
572 void compile_system_filter()
573 {
574 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
575 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
576 }
577
578 void compile_local_matrix_sys_type1()
579 {
580 BaseClass::compile_local_matrix_sys_type1();
581
582 // apply filter to A and B
583 this->filter_velo.local().template at<1>().filter_mat(this->local_matrix_sys_type1.block_a());
584 this->filter_velo.local().template at<1>().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
585 }
586 }; // class StokesBlockedSlipUnitVeloNonePresSystemLevel
587
593 template
594 <
595 int dim_,
596 typename DataType_ = Real,
597 typename IndexType_ = Index,
598 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
599 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
600 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
601 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
602 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
603 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
604 >
606 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
607 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
608 {
609 public:
610 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
611 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
612
613 // define local filter types
617
618 // define global filter types
622
623 // (global) filters
624 GlobalSystemFilter filter_sys;
625 GlobalVeloFilter filter_velo;
626 GlobalPresFilter filter_pres;
627
629 std::size_t bytes() const
630 {
631 return this->matrix_sys.local().bytes () + this->matrix_s.local().bytes() + filter_sys.local().bytes();
632 }
633
634 void compile_system_filter()
635 {
636 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
637 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
638 }
639
640 template<typename D_, typename I_, typename SM_>
641 void convert(const StokesBlockedUnitVeloMeanPresSystemLevel<dim_, D_, I_, SM_> & other)
642 {
643 BaseClass::convert(other);
644 filter_velo.convert(other.filter_velo);
645 filter_pres.convert(other.filter_pres);
646
647 compile_system_filter();
648 }
649
650 template<typename SpacePres_, typename Cubature_>
651 void assemble_pressure_mean_filter(const SpacePres_& space_pres, const Cubature_& cubature)
652 {
653 // get our local pressure filter
654 LocalPresFilter& fil_loc_p = this->filter_pres.local();
655
656 // create two global vectors
657 typename BaseClass::GlobalPresVector vec_glob_v(&this->gate_pres), vec_glob_w(&this->gate_pres);
658
659 // fetch the local vectors
660 typename BaseClass::LocalPresVector& vec_loc_v = vec_glob_v.local();
661 typename BaseClass::LocalPresVector& vec_loc_w = vec_glob_w.local();
662
663 // fetch the frequency vector of the pressure gate
664 typename BaseClass::LocalPresVector& vec_loc_f = this->gate_pres._freqs;
665
666 // assemble the mean filter
667 Assembly::MeanFilterAssembler::assemble(vec_loc_v, vec_loc_w, space_pres, cubature);
668
669 // synchronize the vectors
670 vec_glob_v.sync_1();
671 vec_glob_w.sync_0();
672
673 // build the mean filter
674 fil_loc_p = LocalPresFilter(vec_loc_v.clone(), vec_loc_w.clone(), vec_loc_f.clone(), this->gate_pres.get_comm());
675 }
676
677 void compile_local_matrix_sys_type1()
678 {
679 BaseClass::compile_local_matrix_sys_type1();
680
681 // apply filter to A and B
682 this->filter_velo.local().filter_mat(this->local_matrix_sys_type1.block_a());
683 this->filter_velo.local().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
684 }
685 }; // struct StokesBlockedUnitVeloMeanPresSystemLevel<...>
686
687 template
688 <
689 int dim_,
690 typename DataType_ = Real,
691 typename IndexType_ = Index,
692 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
693 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
694 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
695 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
696 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
697 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
698 >
700 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
701 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
702 {
703 public:
704 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
705 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
706
707 // define local filter types
713
714 // define global filter types
718
719 // (global) filters
720 GlobalSystemFilter filter_sys;
721 GlobalVeloFilter filter_velo;
722 GlobalPresFilter filter_pres;
723
725 std::size_t bytes() const
726 {
727 return this->filter_sys.bytes() + BaseClass::bytes();
728 }
729
730 void compile_system_filter()
731 {
732 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
733 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
734 }
735
736 template<typename D_, typename I_, typename SM_>
737 void convert(const StokesBlockedUnitVeloMeanPresSystemLevel<dim_, D_, I_, SM_> & other)
738 {
739 BaseClass::convert(other);
740 filter_velo.convert(other.filter_velo);
741 filter_pres.convert(other.filter_pres);
742
743 compile_system_filter();
744 }
745
746 template<typename SpacePres_, typename Cubature_>
747 void assemble_global_filters(const SpacePres_& space_pres, const Cubature_& cubature)
748 {
749 // get our local pressure filter
750 LocalPresFilter& fil_loc_p = this->filter_pres.local();
751
752 // create two global vectors
753 typename BaseClass::GlobalPresVector vec_glob_v(&this->gate_pres), vec_glob_w(&this->gate_pres);
754
755 // fetch the local vectors
756 typename BaseClass::LocalPresVector& vec_loc_v = vec_glob_v.local();
757 typename BaseClass::LocalPresVector& vec_loc_w = vec_glob_w.local();
758
759 // fetch the frequency vector of the pressure gate
760 typename BaseClass::LocalPresVector& vec_loc_f = this->gate_pres._freqs;
761
762 // assemble the mean filter
763 Assembly::MeanFilterAssembler::assemble(vec_loc_v, vec_loc_w, space_pres, cubature);
764
765 // synchronize the vectors
766 vec_glob_v.sync_1();
767 vec_glob_w.sync_0();
768
769 // build the mean filter
770 fil_loc_p = LocalPresFilter(vec_loc_v.clone(), vec_loc_w.clone(), vec_loc_f.clone(), this->gate_pres.get_comm());
771
772 // synchronize the slip filter
773 Asm::sync_slip_filter(this->gate_velo, filter_velo.local().template at<0>());
774 }
775
776 void compile_local_matrix_sys_type1()
777 {
778 BaseClass::compile_local_matrix_sys_type1();
779
780 // apply filter to A and B
781 this->filter_velo.local().template at<1>().filter_mat(this->local_matrix_sys_type1.block_a());
782 this->filter_velo.local().template at<1>().filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
783 }
784 }; // class StokesBlockedSlipUnitVeloMeanPresSystemLevel
785
791 template
792 <
793 int dim_,
794 typename DataType_ = Real,
795 typename IndexType_ = Index,
796 typename MatrixBlockA_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, dim_>,
797 typename MatrixBlockB_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, dim_, 1>,
798 typename MatrixBlockD_ = LAFEM::SparseMatrixBCSR<DataType_, IndexType_, 1, dim_>,
799 typename ScalarMatrix_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>,
800 typename TransferMatrixV_ = LAFEM::SparseMatrixBWrappedCSR<DataType_, IndexType_, dim_>,
801 typename TransferMatrixP_ = LAFEM::SparseMatrixCSR<DataType_, IndexType_>
802 >
804 public StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
805 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_>
806 {
807 public:
808 typedef StokesBlockedSystemLevel<dim_, DataType_, IndexType_,
809 MatrixBlockA_, MatrixBlockB_, MatrixBlockD_, ScalarMatrix_, TransferMatrixV_, TransferMatrixP_> BaseClass;
810
811 // define local velocity filter chain
817
818 // define local pressure filter chain
822
823 // define local system filter
825
826 // define global filter types
830
831 // (global) filters
832 GlobalSystemFilter filter_sys;
833 GlobalVeloFilter filter_velo;
834 GlobalPresFilter filter_pres;
835
837 std::size_t bytes() const
838 {
839 return this->filter_sys.bytes() + BaseClass::bytes();
840 }
841
842 void compile_system_filter()
843 {
844 filter_sys.local().template at<0>() = filter_velo.local().clone(LAFEM::CloneMode::Shallow);
845 filter_sys.local().template at<1>() = filter_pres.local().clone(LAFEM::CloneMode::Shallow);
846 }
847
848 template<typename D_, typename I_, typename SM_>
849 void convert(const StokesBlockedCombinedSystemLevel<dim_, D_, I_, SM_> & other)
850 {
851 BaseClass::convert(other);
852 filter_velo.convert(other.filter_velo);
853 filter_pres.convert(other.filter_pres);
854
855 compile_system_filter();
856 }
857
858 LocalVeloSlipFilterSeq& get_local_velo_slip_filter_seq()
859 {
860 return this->filter_velo.local().template at<0>();
861 }
862
863 LocalVeloUnitFilterSeq& get_local_velo_unit_filter_seq()
864 {
865 return this->filter_velo.local().template at<1>();
866 }
867
868 LocalPresMeanFilter& get_local_pres_mean_filter()
869 {
870 return this->filter_pres.local().template at<0>();
871 }
872
873 LocalPresUnitFilter& get_local_pres_unit_filter()
874 {
875 return this->filter_pres.local().template at<1>();
876 }
877
878 template<typename DomainLevel_, typename Space_>
879 void assemble_velocity_unit_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts)
880 {
881 auto& loc_filter = get_local_velo_unit_filter_seq().find_or_add(name);
882 loc_filter.clear();
883 Asm::asm_unit_filter_blocked_homogeneous(loc_filter, dom_level, space, mesh_parts);
884 }
885
886 template<typename DomainLevel_, typename Space_, typename Function_>
887 void assemble_velocity_unit_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts, const Function_& function)
888 {
889 auto& loc_filter = get_local_velo_unit_filter_seq().find_or_add(name);
890 loc_filter.clear();
891 Asm::asm_unit_filter_blocked(loc_filter, dom_level, space, mesh_parts, function);
892 }
893
894 template<typename DomainLevel_, typename Space_>
895 void assemble_velocity_slip_filter(const DomainLevel_& dom_level, const Space_& space, const String& name, const String& mesh_parts)
896 {
897 auto& loc_filter = get_local_velo_slip_filter_seq().find_or_add(name);
898 loc_filter.clear();
899 Asm::asm_slip_filter(loc_filter, dom_level, space, mesh_parts);
900 Asm::sync_slip_filter(this->gate_velo, loc_filter);
901 }
902
903 template<typename SpacePres_>
904 void assemble_pressure_mean_filter(const SpacePres_& space_pres, const String& cubature_name)
905 {
906 this->get_local_pres_mean_filter() = Asm::asm_mean_filter(this->gate_pres, space_pres, cubature_name);
907 }
908
909 void sync_velocity_slip_filters()
910 {
911 for(auto& it : get_local_velo_slip_filter_seq())
912 Asm::sync_slip_filter(this->gate_velo, it.second);
913 }
914
915 void compile_local_matrix_sys_type1()
916 {
917 BaseClass::compile_local_matrix_sys_type1();
918 for(const auto& loc_fil : get_local_velo_unit_filter_seq())
919 {
920 loc_fil.second.filter_mat(this->local_matrix_sys_type1.block_a());
921 loc_fil.second.filter_offdiag_row_mat(this->local_matrix_sys_type1.block_b());
922 }
923 }
924 }; // class StokesBlockedCombinedSystemLevel<...>
925 } // namespace Control
926} // 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.