Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
OpenFOAM
matrices
lduMatrix
lduMatrix
lduMatrix.H
Go to the documentation of this file.
1
/*---------------------------------------------------------------------------*\
2
========= |
3
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4
\\ / O peration |
5
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6
\\/ M anipulation |
7
-------------------------------------------------------------------------------
8
License
9
This file is part of OpenFOAM.
10
11
OpenFOAM is free software: you can redistribute it and/or modify it
12
under the terms of the GNU General Public License as published by
13
the Free Software Foundation, either version 3 of the License, or
14
(at your option) any later version.
15
16
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19
for more details.
20
21
You should have received a copy of the GNU General Public License
22
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23
24
Class
25
Foam::lduMatrix
26
27
Description
28
lduMatrix is a general matrix class in which the coefficients are
29
stored as three arrays, one for the upper triangle, one for the
30
lower triangle and a third for the diagonal.
31
32
Addressing arrays must be supplied for the upper and lower triangles.
33
34
It might be better if this class were organised as a hierachy starting
35
from an empty matrix, then deriving diagonal, symmetric and asymmetric
36
matrices.
37
38
SourceFiles
39
lduMatrixATmul.C
40
lduMatrix.C
41
lduMatrixTemplates.C
42
lduMatrixOperations.C
43
lduMatrixSolver.C
44
lduMatrixPreconditioner.C
45
lduMatrixTests.C
46
lduMatrixUpdateMatrixInterfaces.C
47
48
\*---------------------------------------------------------------------------*/
49
50
#ifndef lduMatrix_H
51
#define lduMatrix_H
52
53
#include <
OpenFOAM/lduMesh.H
>
54
#include <
OpenFOAM/primitiveFieldsFwd.H
>
55
#include <
OpenFOAM/FieldField.H
>
56
#include <
OpenFOAM/lduInterfaceFieldPtrsList.H
>
57
#include <
OpenFOAM/typeInfo.H
>
58
#include <
OpenFOAM/autoPtr.H
>
59
#include <
OpenFOAM/runTimeSelectionTables.H
>
60
61
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62
63
namespace
Foam
64
{
65
66
// Forward declaration of friend functions and operators
67
68
class
lduMatrix;
69
Ostream&
operator<<
(Ostream&,
const
lduMatrix&);
70
71
72
/*---------------------------------------------------------------------------*\
73
Class lduMatrix Declaration
74
\*---------------------------------------------------------------------------*/
75
76
class
lduMatrix
77
{
78
// private data
79
80
//- LDU mesh reference
81
const
lduMesh
& lduMesh_;
82
83
//- Coefficients (not including interfaces)
84
scalarField
*lowerPtr_, *diagPtr_, *upperPtr_;
85
86
87
public
:
88
89
//- Class returned by the solver, containing performance statistics
90
class
solverPerformance
91
{
92
word
solverName_;
93
word
fieldName_;
94
scalar initialResidual_;
95
scalar finalResidual_;
96
label noIterations_;
97
bool
converged_;
98
bool
singular_;
99
100
101
public
:
102
103
// Constructors
104
105
solverPerformance
()
106
:
107
initialResidual_(0),
108
finalResidual_(0),
109
noIterations_(0),
110
converged_(false),
111
singular_(false)
112
{}
113
114
115
solverPerformance
116
(
117
const
word
&
solverName
,
118
const
word
& fieldName,
119
const
scalar iRes = 0,
120
const
scalar fRes = 0,
121
const
label nIter = 0,
122
const
bool
converged
=
false
,
123
const
bool
singular
=
false
124
)
125
:
126
solverName_(solverName),
127
fieldName_(fieldName),
128
initialResidual_(iRes),
129
finalResidual_(fRes),
130
noIterations_(nIter),
131
converged_(
converged
),
132
singular_(
singular
)
133
{}
134
135
136
// Member functions
137
138
//- Return solver name
139
const
word
&
solverName
()
const
140
{
141
return
solverName_;
142
}
143
144
//- Return initial residual
145
scalar
initialResidual
()
const
146
{
147
return
initialResidual_;
148
}
149
150
//- Return initial residual
151
scalar&
initialResidual
()
152
{
153
return
initialResidual_;
154
}
155
156
157
//- Return final residual
158
scalar
finalResidual
()
const
159
{
160
return
finalResidual_;
161
}
162
163
//- Return final residual
164
scalar&
finalResidual
()
165
{
166
return
finalResidual_;
167
}
168
169
170
//- Return number of iterations
171
label
nIterations
()
const
172
{
173
return
noIterations_;
174
}
175
176
//- Return number of iterations
177
label&
nIterations
()
178
{
179
return
noIterations_;
180
}
181
182
183
//- Has the solver converged?
184
bool
converged
()
const
185
{
186
return
converged_;
187
}
188
189
//- Is the matrix singular?
190
bool
singular
()
const
191
{
192
return
singular_;
193
}
194
195
//- Convergence test
196
bool
checkConvergence
197
(
198
const
scalar tolerance,
199
const
scalar relTolerance
200
);
201
202
//- Singularity test
203
bool
checkSingularity
(
const
scalar
residual
);
204
205
//- Print summary of solver performance
206
void
print
()
const
;
207
};
208
209
210
//- Abstract base-class for lduMatrix solvers
211
class
solver
212
{
213
protected
:
214
215
// Protected data
216
217
word
fieldName_
;
218
const
lduMatrix
&
matrix_
;
219
const
FieldField<Field, scalar>
&
interfaceBouCoeffs_
;
220
const
FieldField<Field, scalar>
&
interfaceIntCoeffs_
;
221
lduInterfaceFieldPtrsList
interfaces_
;
222
223
//- dictionary of controls
224
dictionary
controlDict_
;
225
226
//- Maximum number of iterations in the solver
227
label
maxIter_
;
228
229
//- Final convergence tolerance
230
scalar
tolerance_
;
231
232
//- Convergence tolerance relative to the initial
233
scalar
relTol_
;
234
235
236
// Protected Member Functions
237
238
//- Read the control parameters from the controlDict_
239
virtual
void
readControls
();
240
241
242
public
:
243
244
//- Runtime type information
245
virtual
const
word
&
type
()
const
= 0;
246
247
248
// Declare run-time constructor selection tables
249
250
declareRunTimeSelectionTable
251
(
252
autoPtr
,
253
solver
,
254
symMatrix,
255
(
256
const
word
&
fieldName
,
257
const
lduMatrix
&
matrix
,
258
const
FieldField<Field, scalar>
&
interfaceBouCoeffs
,
259
const
FieldField<Field, scalar>
&
interfaceIntCoeffs
,
260
const
lduInterfaceFieldPtrsList
&
interfaces
,
261
const
dictionary
& solverControls
262
),
263
(
264
fieldName,
265
matrix,
266
interfaceBouCoeffs,
267
interfaceIntCoeffs,
268
interfaces,
269
solverControls
270
)
271
);
272
273
declareRunTimeSelectionTable
274
(
275
autoPtr
,
276
solver
,
277
asymMatrix,
278
(
279
const
word
& fieldName,
280
const
lduMatrix
& matrix,
281
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
282
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
283
const
lduInterfaceFieldPtrsList
& interfaces,
284
const
dictionary
& solverControls
285
),
286
(
287
fieldName,
288
matrix,
289
interfaceBouCoeffs,
290
interfaceIntCoeffs,
291
interfaces,
292
solverControls
293
)
294
);
295
296
297
// Constructors
298
299
solver
300
(
301
const
word
& fieldName,
302
const
lduMatrix
& matrix,
303
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
304
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
305
const
lduInterfaceFieldPtrsList
& interfaces,
306
const
dictionary
& solverControls
307
);
308
309
// Selectors
310
311
//- Return a new solver
312
static
autoPtr<solver>
New
313
(
314
const
word
& fieldName,
315
const
lduMatrix
& matrix,
316
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
317
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
318
const
lduInterfaceFieldPtrsList
& interfaces,
319
const
dictionary
& solverControls
320
);
321
322
323
324
// Destructor
325
326
virtual
~solver
()
327
{}
328
329
330
// Member functions
331
332
// Access
333
334
const
word
&
fieldName
()
const
335
{
336
return
fieldName_
;
337
}
338
339
const
lduMatrix
&
matrix
()
const
340
{
341
return
matrix_
;
342
}
343
344
const
FieldField<Field, scalar>
&
interfaceBouCoeffs
()
const
345
{
346
return
interfaceBouCoeffs_
;
347
}
348
349
const
FieldField<Field, scalar>
&
interfaceIntCoeffs
()
const
350
{
351
return
interfaceIntCoeffs_
;
352
}
353
354
const
lduInterfaceFieldPtrsList
&
interfaces
()
const
355
{
356
return
interfaces_
;
357
}
358
359
360
//- Read and reset the solver parameters from the given stream
361
virtual
void
read
(
const
dictionary
&);
362
363
virtual
solverPerformance
solve
364
(
365
scalarField
&
psi
,
366
const
scalarField
& source,
367
const
direction
cmpt=0
368
)
const
= 0;
369
370
//- Return the matrix norm used to normalise the residual for the
371
// stopping criterion
372
scalar
normFactor
373
(
374
const
scalarField
&
psi
,
375
const
scalarField
& source,
376
const
scalarField
& Apsi,
377
scalarField
& tmpField
378
)
const
;
379
};
380
381
382
//- Abstract base-class for lduMatrix smoothers
383
class
smoother
384
{
385
protected
:
386
387
// Protected data
388
389
word
fieldName_
;
390
const
lduMatrix
&
matrix_
;
391
const
FieldField<Field, scalar>
&
interfaceBouCoeffs_
;
392
const
FieldField<Field, scalar>
&
interfaceIntCoeffs_
;
393
const
lduInterfaceFieldPtrsList
&
interfaces_
;
394
395
396
public
:
397
398
//- Find the smoother name (directly or from a sub-dictionary)
399
static
word
getName
(
const
dictionary
&);
400
401
//- Runtime type information
402
virtual
const
word
&
type
()
const
= 0;
403
404
405
// Declare run-time constructor selection tables
406
407
declareRunTimeSelectionTable
408
(
409
autoPtr
,
410
smoother
,
411
symMatrix,
412
(
413
const
word
& fieldName,
414
const
lduMatrix
& matrix,
415
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
416
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
417
const
lduInterfaceFieldPtrsList
& interfaces
418
),
419
(
420
fieldName,
421
matrix,
422
interfaceBouCoeffs,
423
interfaceIntCoeffs,
424
interfaces
425
)
426
);
427
428
declareRunTimeSelectionTable
429
(
430
autoPtr
,
431
smoother
,
432
asymMatrix,
433
(
434
const
word
& fieldName,
435
const
lduMatrix
& matrix,
436
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
437
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
438
const
lduInterfaceFieldPtrsList
& interfaces
439
),
440
(
441
fieldName,
442
matrix,
443
interfaceBouCoeffs,
444
interfaceIntCoeffs,
445
interfaces
446
)
447
);
448
449
450
// Constructors
451
452
smoother
453
(
454
const
word
& fieldName,
455
const
lduMatrix
& matrix,
456
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
457
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
458
const
lduInterfaceFieldPtrsList
& interfaces
459
);
460
461
462
// Selectors
463
464
//- Return a new smoother
465
static
autoPtr<smoother>
New
466
(
467
const
word
& fieldName,
468
const
lduMatrix
& matrix,
469
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
470
const
FieldField<Field, scalar>
& interfaceIntCoeffs,
471
const
lduInterfaceFieldPtrsList
& interfaces,
472
const
dictionary
& solverControls
473
);
474
475
476
// Destructor
477
478
virtual
~smoother
()
479
{}
480
481
482
// Member functions
483
484
// Access
485
486
const
word
&
fieldName
()
const
487
{
488
return
fieldName_
;
489
}
490
491
const
lduMatrix
&
matrix
()
const
492
{
493
return
matrix_
;
494
}
495
496
const
FieldField<Field, scalar>
&
interfaceBouCoeffs
()
const
497
{
498
return
interfaceBouCoeffs_
;
499
}
500
501
const
FieldField<Field, scalar>
&
interfaceIntCoeffs
()
const
502
{
503
return
interfaceIntCoeffs_
;
504
}
505
506
const
lduInterfaceFieldPtrsList
&
interfaces
()
const
507
{
508
return
interfaces_
;
509
}
510
511
512
//- Smooth the solution for a given number of sweeps
513
virtual
void
smooth
514
(
515
scalarField
&
psi
,
516
const
scalarField
& source,
517
const
direction
cmpt,
518
const
label nSweeps
519
)
const
= 0;
520
};
521
522
523
//- Abstract base-class for lduMatrix preconditioners
524
class
preconditioner
525
{
526
protected
:
527
528
// Protected data
529
530
//- Reference to the base-solver this preconditioner is used with
531
const
solver
&
solver_
;
532
533
534
public
:
535
536
//- Find the preconditioner name (directly or from a sub-dictionary)
537
static
word
getName
(
const
dictionary
&);
538
539
//- Runtime type information
540
virtual
const
word
&
type
()
const
= 0;
541
542
543
// Declare run-time constructor selection tables
544
545
declareRunTimeSelectionTable
546
(
547
autoPtr
,
548
preconditioner
,
549
symMatrix,
550
(
551
const
solver
& sol,
552
const
dictionary
& solverControls
553
),
554
(sol, solverControls)
555
);
556
557
declareRunTimeSelectionTable
558
(
559
autoPtr
,
560
preconditioner
,
561
asymMatrix,
562
(
563
const
solver
& sol,
564
const
dictionary
& solverControls
565
),
566
(sol, solverControls)
567
);
568
569
570
// Constructors
571
572
preconditioner
573
(
574
const
solver
& sol
575
)
576
:
577
solver_
(sol)
578
{}
579
580
581
// Selectors
582
583
//- Return a new preconditioner
584
static
autoPtr<preconditioner>
New
585
(
586
const
solver
& sol,
587
const
dictionary
& solverControls
588
);
589
590
591
// Destructor
592
593
virtual
~preconditioner
()
594
{}
595
596
597
// Member functions
598
599
//- Read and reset the preconditioner parameters
600
// from the given stream
601
virtual
void
read
(
const
dictionary
&)
602
{}
603
604
//- Return wA the preconditioned form of residual rA
605
virtual
void
precondition
606
(
607
scalarField
& wA,
608
const
scalarField
& rA,
609
const
direction
cmpt=0
610
)
const
= 0;
611
612
//- Return wT the transpose-matrix preconditioned form of
613
// residual rT.
614
// This is only required for preconditioning asymmetric matrices.
615
virtual
void
preconditionT
616
(
617
scalarField
& wT,
618
const
scalarField
& rT,
619
const
direction
cmpt=0
620
)
const
621
{
622
notImplemented
623
(
624
type
() +
"::preconditionT"
625
"(scalarField& wT, const scalarField& rT, "
626
"const direction cmpt)"
627
);
628
}
629
};
630
631
632
// Static data
633
634
// Declare name of the class and its debug switch
635
ClassName
(
"lduMatrix"
);
636
637
//- Large scalar for the use in solvers
638
static
const
scalar
great_
;
639
640
//- Small scalar for the use in solvers
641
static
const
scalar
small_
;
642
643
644
// Constructors
645
646
//- Construct given an LDU addressed mesh.
647
// The coefficients are initially empty for subsequent setting.
648
lduMatrix
(
const
lduMesh
&);
649
650
//- Construct as copy
651
lduMatrix
(
const
lduMatrix
&);
652
653
//- Construct as copy or re-use as specified.
654
lduMatrix
(
lduMatrix
&,
bool
reUse);
655
656
//- Construct given an LDU addressed mesh and an Istream
657
// from which the coefficients are read
658
lduMatrix
(
const
lduMesh
&,
Istream
&);
659
660
661
// Destructor
662
663
~lduMatrix
();
664
665
666
// Member functions
667
668
// Access to addressing
669
670
//- Return the LDU mesh from which the addressing is obtained
671
const
lduMesh
&
mesh
()
const
672
{
673
return
lduMesh_;
674
}
675
676
//- Return the LDU addressing
677
const
lduAddressing
&
lduAddr
()
const
678
{
679
return
lduMesh_.
lduAddr
();
680
}
681
682
//- Return the patch evaluation schedule
683
const
lduSchedule
&
patchSchedule
()
const
684
{
685
return
lduAddr
().
patchSchedule
();
686
}
687
688
689
// Access to coefficients
690
691
scalarField
&
lower
();
692
scalarField
&
diag
();
693
scalarField
&
upper
();
694
695
const
scalarField
&
lower
()
const
;
696
const
scalarField
&
diag
()
const
;
697
const
scalarField
&
upper
()
const
;
698
699
bool
hasDiag
()
const
700
{
701
return
(diagPtr_);
702
}
703
704
bool
hasUpper
()
const
705
{
706
return
(upperPtr_);
707
}
708
709
bool
hasLower
()
const
710
{
711
return
(lowerPtr_);
712
}
713
714
bool
diagonal
()
const
715
{
716
return
(diagPtr_ && !lowerPtr_ && !upperPtr_);
717
}
718
719
bool
symmetric
()
const
720
{
721
return
(diagPtr_ && (!lowerPtr_ && upperPtr_));
722
}
723
724
bool
asymmetric
()
const
725
{
726
return
(diagPtr_ && lowerPtr_ && upperPtr_);
727
}
728
729
730
// operations
731
732
void
sumDiag
();
733
void
negSumDiag
();
734
735
void
sumMagOffDiag
(
scalarField
& sumOff)
const
;
736
737
//- Matrix multiplication with updated interfaces.
738
void
Amul
739
(
740
scalarField
&,
741
const
tmp<scalarField>
&,
742
const
FieldField<Field, scalar>
&,
743
const
lduInterfaceFieldPtrsList
&,
744
const
direction
cmpt
745
)
const
;
746
747
//- Matrix transpose multiplication with updated interfaces.
748
void
Tmul
749
(
750
scalarField
&,
751
const
tmp<scalarField>
&,
752
const
FieldField<Field, scalar>
&,
753
const
lduInterfaceFieldPtrsList
&,
754
const
direction
cmpt
755
)
const
;
756
757
758
//- Sum the coefficients on each row of the matrix
759
void
sumA
760
(
761
scalarField
&,
762
const
FieldField<Field, scalar>
&,
763
const
lduInterfaceFieldPtrsList
&
764
)
const
;
765
766
767
void
residual
768
(
769
scalarField
& rA,
770
const
scalarField
&
psi
,
771
const
scalarField
& source,
772
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
773
const
lduInterfaceFieldPtrsList
& interfaces,
774
const
direction
cmpt
775
)
const
;
776
777
tmp<scalarField>
residual
778
(
779
const
scalarField
&
psi
,
780
const
scalarField
& source,
781
const
FieldField<Field, scalar>
& interfaceBouCoeffs,
782
const
lduInterfaceFieldPtrsList
& interfaces,
783
const
direction
cmpt
784
)
const
;
785
786
787
//- Initialise the update of interfaced interfaces
788
// for matrix operations
789
void
initMatrixInterfaces
790
(
791
const
FieldField<Field, scalar>
& interfaceCoeffs,
792
const
lduInterfaceFieldPtrsList
& interfaces,
793
const
scalarField
& psiif,
794
scalarField
& result,
795
const
direction
cmpt
796
)
const
;
797
798
//- Update interfaced interfaces for matrix operations
799
void
updateMatrixInterfaces
800
(
801
const
FieldField<Field, scalar>
& interfaceCoeffs,
802
const
lduInterfaceFieldPtrsList
& interfaces,
803
const
scalarField
& psiif,
804
scalarField
& result,
805
const
direction
cmpt
806
)
const
;
807
808
809
template
<
class
Type>
810
tmp<Field<Type>
>
H
(
const
Field<Type>
&)
const
;
811
812
template
<
class
Type>
813
tmp<Field<Type>
>
H
(
const
tmp
<
Field<Type>
>&)
const
;
814
815
tmp<scalarField>
H1
()
const
;
816
817
template
<
class
Type>
818
tmp<Field<Type>
>
faceH
(
const
Field<Type>
&)
const
;
819
820
template
<
class
Type>
821
tmp<Field<Type>
>
faceH
(
const
tmp
<
Field<Type>
>&)
const
;
822
823
824
// Member operators
825
826
void
operator=
(
const
lduMatrix
&);
827
828
void
negate
();
829
830
void
operator+=
(
const
lduMatrix
&);
831
void
operator-=
(
const
lduMatrix
&);
832
833
void
operator*=
(
const
scalarField
&);
834
void
operator*=
(scalar);
835
836
837
// Ostream operator
838
839
friend
Ostream
&
operator<<
(
Ostream
&,
const
lduMatrix
&);
840
};
841
842
843
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
844
845
}
// End namespace Foam
846
847
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
848
849
#ifdef NoRepository
850
# include <
OpenFOAM/lduMatrixTemplates.C
>
851
#endif
852
853
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
854
855
#endif
856
857
// ************************ vim: set sw=4 sts=4 et: ************************ //