src/finiteArea/faMatrices/faMatrix/faMatrix.H

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     |
   5     \\  /    A nd           | For copyright notice see file Copyright
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by the
  13     Free Software Foundation, either version 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Class
  25     faMatrix<Type>
  26
  27 Description
  28     Finite-Area matrix.
  29
  30 SourceFiles
  31     faMatrix.C
  32     faMatrixSolve.C
  33
  34 \*---------------------------------------------------------------------------*/
  35
  36 #ifndef faMatrix_H
  37 #define faMatrix_H
  38
  39 #include "areaFields.H"
  40 #include "edgeFields.H"
  41 #include "lduMatrix.H"
  42 #include "tmp.H"
  43 #include "autoPtr.H"
  44 #include "dimensionedTypes.H"
  45 #include "className.H"
  46
  47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  48
  49 namespace Foam
  50 {
  51
  52 // * * * * * * Forward declaration of template friend fuctions * * * * * * * //
  53
  54 template<class Type>
  55 class faMatrix;
  56
  57 template<class Type>
  58 Ostream& operator<<(Ostream&, const faMatrix<Type>&);
  59
  60
  61 /*---------------------------------------------------------------------------*\
  62                            Class faMatrix Declaration
  63 \*---------------------------------------------------------------------------*/
  64
  65 template<class Type>
  66 class faMatrix
  67 :
  68     public refCount,
  69     public lduMatrix
  70 {
  71     // Private data
  72
  73         // Reference to GeometricField<Type, faPatchField, areaMesh>
  74         GeometricField<Type, faPatchField, areaMesh>& psi_;
  75
  76         //- Dimension set
  77         dimensionSet dimensions_;
  78
  79         //- Source term
  80         Field<Type> source_;
  81
  82         //- Boundary scalar field containing pseudo-matrix coeffs
  83         //  for internal faces
  84         FieldField<Field, Type> internalCoeffs_;
  85
  86         //- Boundary scalar field containing pseudo-matrix coeffs
  87         //  for boundary faces
  88         FieldField<Field, Type> boundaryCoeffs_;
  89
  90
  91         //- Face flux field for non-orthogonal correction
  92         mutable GeometricField<Type, faePatchField, edgeMesh>
  93             *faceFluxCorrectionPtr_;
  94
  95         direction solvingComponent;
  96
  97
  98     // Private member functions
  99
 100         //- Add patch contribution to internal field
 101         template<class Type2>
 102         void addToInternalField
 103         (
 104             const unallocLabelList& addr,
 105             const Field<Type2>& pf,
 106             Field<Type2>& intf
 107         ) const;
 108
 109         template<class Type2>
 110         void addToInternalField
 111         (
 112             const unallocLabelList& addr,
 113             const tmp<Field<Type2> >& tpf,
 114             Field<Type2>& intf
 115         ) const;
 116
 117         //- Subtract patch contribution from internal field
 118         template<class Type2>
 119         void subtractFromInternalField
 120         (
 121             const unallocLabelList& addr,
 122             const Field<Type2>& pf,
 123             Field<Type2>& intf
 124         ) const;
 125
 126         template<class Type2>
 127         void subtractFromInternalField
 128         (
 129             const unallocLabelList& addr,
 130             const tmp<Field<Type2> >& tpf,
 131             Field<Type2>& intf
 132         ) const;
 133
 134
 135         // Matrix completion functionality
 136
 137             void addBoundaryDiag
 138             (
 139                 scalarField& diag,
 140                 const direction cmpt
 141             ) const;
 142
 143             void addCmptAvBoundaryDiag(scalarField& diag) const;
 144
 145             void addBoundarySource
 146             (
 147                 Field<Type>& source,
 148                 const bool couples = true
 149             ) const;
 150
 151
 152 public:
 153
 154     //- Solver class returned by the solver function
 155     class faSolver
 156     {
 157         faMatrix<Type>& faMat_;
 158
 159         autoPtr<lduSolver> solver_;
 160
 161     public:
 162
 163         // Constructors
 164
 165             faSolver(faMatrix<Type>& faMat, autoPtr<lduSolver> sol)
 166             :
 167                 faMat_(faMat),
 168                 solver_(sol)
 169             {}
 170
 171
 172         // Member functions
 173
 174             //- Solve returning the solution statistics.
 175             //  Solver controls read from dictionary
 176             lduSolverPerformance solve(const dictionary&);
 177
 178             //- Solve returning the solution statistics.
 179             //  Solver controls read from faSolution
 180             lduSolverPerformance solve();
 181     };
 182
 183
 184     ClassName("faMatrix");
 185
 186
 187     // Constructors
 188
 189         //- Construct given a field to solve for
 190         faMatrix
 191         (
 192             GeometricField<Type, faPatchField, areaMesh>&,
 193             const dimensionSet&
 194         );
 195
 196         //- Construct as copy
 197         faMatrix(const faMatrix<Type>&);
 198
 199         //- Construct from Istream given field to solve for
 200         faMatrix(GeometricField<Type, faPatchField, areaMesh>&, Istream&);
 201
 202
 203     // Destructor
 204
 205         virtual ~faMatrix();
 206
 207
 208     // Member functions
 209
 210         // Access
 211
 212             const GeometricField<Type, faPatchField, areaMesh>& psi() const
 213             {
 214                 return psi_;
 215             }
 216
 217             GeometricField<Type, faPatchField, areaMesh>& psi()
 218             {
 219                 return psi_;
 220             }
 221
 222             const dimensionSet& dimensions() const
 223             {
 224                 return dimensions_;
 225             }
 226
 227             Field<Type>& source()
 228             {
 229                 return source_;
 230             }
 231
 232             const Field<Type>& source() const
 233             {
 234                 return source_;
 235             }
 236
 237             //- faBoundary scalar field containing pseudo-matrix coeffs
 238             //  for internal cells
 239             FieldField<Field, Type>& internalCoeffs()
 240             {
 241                 return internalCoeffs_;
 242             }
 243
 244             //- faBoundary scalar field containing pseudo-matrix coeffs
 245             //  for boundary cells
 246             FieldField<Field, Type>& boundaryCoeffs()
 247             {
 248                 return boundaryCoeffs_;
 249             }
 250
 251
 252             //- Declare return type of the faceFluxCorrectionPtr() function
 253             typedef GeometricField<Type, faePatchField, edgeMesh>
 254                 *edgeTypeFieldPtr;
 255
 256             //- Return pointer to face-flux non-orthogonal correction field
 257             edgeTypeFieldPtr& faceFluxCorrectionPtr()
 258             {
 259                 return faceFluxCorrectionPtr_;
 260             }
 261
 262
 263         // Operations
 264
 265             //- Set solution in given cells and eliminate corresponding
 266             //  equations from the matrix
 267             void setValues
 268             (
 269                 const labelList& faces,
 270                 const Field<Type>& values
 271             );
 272
 273             //- Set reference level for solution
 274             void setReference
 275             (
 276                 const label facei,
 277                 const Type& value
 278             );
 279
 280             //- Set reference level for a component of the solution
 281             //  on a given patch face
 282             void setComponentReference
 283             (
 284                 const label patchi,
 285                 const label facei,
 286                 const direction cmpt,
 287                 const scalar value
 288             );
 289
 290             //- Relax matrix (for steady-state solution).
 291             //  alpha = 1 : diagonally equal
 292             //  alpha < 1 :    ,,      dominant
 293             //  alpha = 0 : do nothing
 294             void relax(const scalar alpha);
 295
 296             //- Relax matrix (for steadty-state solution).
 297             //  alpha is read from controlDict
 298             void relax();
 299
 300             //- Construct and return the solver
 301             //  Solver controls read from Istream
 302             autoPtr<faSolver> solver(const dictionary&);
 303
 304             //- Construct and return the solver
 305             //  Solver controls read from faSolution
 306             autoPtr<faSolver> solver();
 307
 308             //- Solve returning the solution statistics.
 309             //  Solver controls read Istream
 310             lduSolverPerformance solve(const dictionary&);
 311
 312             //- Solve returning the solution statistics.
 313             //  Solver controls read from faSolution
 314             lduSolverPerformance solve();
 315
 316             //- Return the matrix residual
 317             tmp<Field<Type> > residual() const;
 318
 319             //- Return the matrix diagonal
 320             tmp<scalarField> D() const;
 321
 322             //- Return the central coefficient
 323             tmp<areaScalarField> A() const;
 324
 325             //- Return the H operation source
 326             tmp<GeometricField<Type, faPatchField, areaMesh> > H() const;
 327
 328             //- Return the face-flux field from the matrix
 329             tmp<GeometricField<Type, faePatchField, edgeMesh> > flux() const;
 330
 331
 332     // Member operators
 333
 334         void operator=(const faMatrix<Type>&);
 335         void operator=(const tmp<faMatrix<Type> >&);
 336
 337         void negate();
 338
 339         void operator+=(const faMatrix<Type>&);
 340         void operator+=(const tmp<faMatrix<Type> >&);
 341
 342         void operator-=(const faMatrix<Type>&);
 343         void operator-=(const tmp<faMatrix<Type> >&);
 344
 345         void operator+=(const GeometricField<Type,faPatchField,areaMesh>&);
 346         void operator+=(const tmp<GeometricField<Type,faPatchField,areaMesh> >&);
 347
 348         void operator-=(const GeometricField<Type,faPatchField,areaMesh>&);
 349         void operator-=(const tmp<GeometricField<Type,faPatchField,areaMesh> >&);
 350
 351         void operator+=(const dimensioned<Type>&);
 352         void operator-=(const dimensioned<Type>&);
 353
 354         void operator*=(const areaScalarField&);
 355         void operator*=(const tmp<areaScalarField>&);
 356
 357         void operator*=(const dimensioned<scalar>&);
 358
 359
 360     // Ostream operator
 361
 362         friend Ostream& operator<< <Type>
 363         (
 364             Ostream&,
 365             const faMatrix<Type>&
 366         );
 367 };
 368
 369
 370 // * * * * * * * * * * * * * * * Global functions  * * * * * * * * * * * * * //
 371
 372 template<class Type>
 373 void checkMethod
 374 (
 375     const faMatrix<Type>&,
 376     const faMatrix<Type>&,
 377     const char*
 378 );
 379
 380 template<class Type>
 381 void checkMethod
 382 (
 383     const faMatrix<Type>&,
 384     const GeometricField<Type, faPatchField, areaMesh>&,
 385     const char*
 386 );
 387
 388 template<class Type>
 389 void checkMethod
 390 (
 391     const faMatrix<Type>&,
 392     const dimensioned<Type>&,
 393     const char*
 394 );
 395
 396
 397 //- Solve returning the solution statistics given convergence tolerance
 398 //  Solver controls read Istream
 399 template<class Type>
 400 lduSolverPerformance solve(faMatrix<Type>&, Istream&);
 401
 402
 403 //- Solve returning the solution statistics given convergence tolerance,
 404 //  deleting temporary matrix after solution.
 405 //  Solver controls read Istream
 406 template<class Type>
 407 lduSolverPerformance solve(const tmp<faMatrix<Type> >&, Istream&);
 408
 409
 410 //- Solve returning the solution statistics given convergence tolerance
 411 //  Solver controls read faSolution
 412 template<class Type>
 413 lduSolverPerformance solve(faMatrix<Type>&);
 414
 415
 416 //- Solve returning the solution statistics given convergence tolerance,
 417 //  deleting temporary matrix after solution.
 418 //  Solver controls read faSolution
 419 template<class Type>
 420 lduSolverPerformance solve(const tmp<faMatrix<Type> >&);
 421
 422
 423 // * * * * * * * * * * * * * * * Global operators  * * * * * * * * * * * * * //
 424
 425 template<class Type>
 426 tmp<faMatrix<Type> > operator-
 427 (
 428     const faMatrix<Type>&
 429 );
 430
 431 template<class Type>
 432 tmp<faMatrix<Type> > operator-
 433 (
 434     const tmp<faMatrix<Type> >&
 435 );
 436
 437 template<class Type>
 438 tmp<faMatrix<Type> > operator+
 439 (
 440     const faMatrix<Type>&,
 441     const faMatrix<Type>&
 442 );
 443
 444 template<class Type>
 445 tmp<faMatrix<Type> > operator+
 446 (
 447     const tmp<faMatrix<Type> >&,
 448     const faMatrix<Type>&
 449 );
 450
 451 template<class Type>
 452 tmp<faMatrix<Type> > operator+
 453 (
 454     const faMatrix<Type>&,
 455     const tmp<faMatrix<Type> >&
 456 );
 457
 458 template<class Type>
 459 tmp<faMatrix<Type> > operator+
 460 (
 461     const tmp<faMatrix<Type> >&,
 462     const tmp<faMatrix<Type> >&
 463 );
 464
 465 template<class Type>
 466 tmp<faMatrix<Type> > operator-
 467 (
 468     const faMatrix<Type>&,
 469     const faMatrix<Type>&
 470 );
 471
 472 template<class Type>
 473 tmp<faMatrix<Type> > operator-
 474 (
 475     const tmp<faMatrix<Type> >&,
 476     const faMatrix<Type>&
 477 );
 478
 479 template<class Type>
 480 tmp<faMatrix<Type> > operator-
 481 (
 482     const faMatrix<Type>&,
 483     const tmp<faMatrix<Type> >&
 484 );
 485
 486 template<class Type>
 487 tmp<faMatrix<Type> > operator-
 488 (
 489     const tmp<faMatrix<Type> >&,
 490     const tmp<faMatrix<Type> >&
 491 );
 492
 493 template<class Type>
 494 tmp<faMatrix<Type> > operator==
 495 (
 496     const faMatrix<Type>&,
 497     const faMatrix<Type>&
 498 );
 499
 500 template<class Type>
 501 tmp<faMatrix<Type> > operator==
 502 (
 503     const tmp<faMatrix<Type> >&,
 504     const faMatrix<Type>&
 505 );
 506
 507 template<class Type>
 508 tmp<faMatrix<Type> > operator==
 509 (
 510     const faMatrix<Type>&,
 511     const tmp<faMatrix<Type> >&
 512 );
 513
 514 template<class Type>
 515 tmp<faMatrix<Type> > operator==
 516 (
 517     const tmp<faMatrix<Type> >&,
 518     const tmp<faMatrix<Type> >&
 519 );
 520
 521 template<class Type>
 522 tmp<faMatrix<Type> > operator+
 523 (
 524     const faMatrix<Type>&,
 525     const GeometricField<Type, faPatchField, areaMesh>&
 526 );
 527
 528 template<class Type>
 529 tmp<faMatrix<Type> > operator+
 530 (
 531     const tmp<faMatrix<Type> >&,
 532     const GeometricField<Type, faPatchField, areaMesh>&
 533 );
 534
 535 template<class Type>
 536 tmp<faMatrix<Type> > operator+
 537 (
 538     const faMatrix<Type>&,
 539     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 540 );
 541
 542 template<class Type>
 543 tmp<faMatrix<Type> > operator+
 544 (
 545     const tmp<faMatrix<Type> >&,
 546     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 547 );
 548
 549 template<class Type>
 550 tmp<faMatrix<Type> > operator+
 551 (
 552     const GeometricField<Type, faPatchField, areaMesh>&,
 553     const faMatrix<Type>&
 554 );
 555
 556 template<class Type>
 557 tmp<faMatrix<Type> > operator+
 558 (
 559     const GeometricField<Type, faPatchField, areaMesh>&,
 560     const tmp<faMatrix<Type> >&
 561 );
 562
 563 template<class Type>
 564 tmp<faMatrix<Type> > operator+
 565 (
 566     const tmp<GeometricField<Type, faPatchField, areaMesh> >&,
 567     const faMatrix<Type>&
 568 );
 569
 570 template<class Type>
 571 tmp<faMatrix<Type> > operator+
 572 (
 573     const tmp<GeometricField<Type, faPatchField, areaMesh> >&,
 574     const tmp<faMatrix<Type> >&
 575 );
 576
 577 template<class Type>
 578 tmp<faMatrix<Type> > operator-
 579 (
 580     const faMatrix<Type>&,
 581     const GeometricField<Type, faPatchField, areaMesh>&
 582 );
 583
 584 template<class Type>
 585 tmp<faMatrix<Type> > operator-
 586 (
 587     const tmp<faMatrix<Type> >&,
 588     const GeometricField<Type, faPatchField, areaMesh>&
 589 );
 590
 591 template<class Type>
 592 tmp<faMatrix<Type> > operator-
 593 (
 594     const faMatrix<Type>&,
 595     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 596 );
 597
 598 template<class Type>
 599 tmp<faMatrix<Type> > operator-
 600 (
 601     const tmp<faMatrix<Type> >&,
 602     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 603 );
 604
 605 template<class Type>
 606 tmp<faMatrix<Type> > operator-
 607 (
 608     const GeometricField<Type, faPatchField, areaMesh>&,
 609     const faMatrix<Type>&
 610 );
 611
 612 template<class Type>
 613 tmp<faMatrix<Type> > operator-
 614 (
 615     const GeometricField<Type, faPatchField, areaMesh>&,
 616     const tmp<faMatrix<Type> >&
 617 );
 618
 619 template<class Type>
 620 tmp<faMatrix<Type> > operator-
 621 (
 622     const tmp<GeometricField<Type, faPatchField, areaMesh> >&,
 623     const faMatrix<Type>&
 624 );
 625
 626 template<class Type>
 627 tmp<faMatrix<Type> > operator-
 628 (
 629     const tmp<GeometricField<Type, faPatchField, areaMesh> >&,
 630     const tmp<faMatrix<Type> >&
 631 );
 632
 633 template<class Type>
 634 tmp<faMatrix<Type> > operator+
 635 (
 636     const tmp<faMatrix<Type> >&,
 637     const dimensioned<Type>&
 638 );
 639
 640 template<class Type>
 641 tmp<faMatrix<Type> > operator+
 642 (
 643     const dimensioned<Type>&,
 644     const tmp<faMatrix<Type> >&
 645 );
 646
 647 template<class Type>
 648 tmp<faMatrix<Type> > operator-
 649 (
 650     const tmp<faMatrix<Type> >&,
 651     const dimensioned<Type>&
 652 );
 653
 654 template<class Type>
 655 tmp<faMatrix<Type> > operator-
 656 (
 657     const dimensioned<Type>&,
 658     const tmp<faMatrix<Type> >&
 659 );
 660
 661 template<class Type>
 662 tmp<faMatrix<Type> > operator==
 663 (
 664     const faMatrix<Type>&,
 665     const GeometricField<Type, faPatchField, areaMesh>&
 666 );
 667
 668 template<class Type>
 669 tmp<faMatrix<Type> > operator==
 670 (
 671     const tmp<faMatrix<Type> >&,
 672     const GeometricField<Type, faPatchField, areaMesh>&
 673 );
 674
 675 template<class Type>
 676 tmp<faMatrix<Type> > operator==
 677 (
 678     const faMatrix<Type>&,
 679     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 680 );
 681
 682 template<class Type>
 683 tmp<faMatrix<Type> > operator==
 684 (
 685     const tmp<faMatrix<Type> >&,
 686     const tmp<GeometricField<Type, faPatchField, areaMesh> >&
 687 );
 688
 689 template<class Type>
 690 tmp<faMatrix<Type> > operator==
 691 (
 692     const faMatrix<Type>&,
 693     const dimensioned<Type>&
 694 );
 695
 696 template<class Type>
 697 tmp<faMatrix<Type> > operator==
 698 (
 699     const tmp<faMatrix<Type> >&,
 700     const dimensioned<Type>&
 701 );
 702
 703
 704 template<class Type>
 705 tmp<faMatrix<Type> > operator*
 706 (
 707     const areaScalarField&,
 708     const faMatrix<Type>&
 709 );
 710
 711 template<class Type>
 712 tmp<faMatrix<Type> > operator*
 713 (
 714     const areaScalarField&,
 715     const tmp<faMatrix<Type> >&
 716 );
 717
 718 template<class Type>
 719 tmp<faMatrix<Type> > operator*
 720 (
 721     const tmp<areaScalarField>&,
 722     const faMatrix<Type>&
 723 );
 724
 725 template<class Type>
 726 tmp<faMatrix<Type> > operator*
 727 (
 728     const tmp<areaScalarField>&,
 729     const tmp<faMatrix<Type> >&
 730 );
 731
 732
 733 template<class Type>
 734 tmp<faMatrix<Type> > operator*
 735 (
 736     const dimensioned<scalar>&,
 737     const faMatrix<Type>&
 738 );
 739
 740 template<class Type>
 741 tmp<faMatrix<Type> > operator*
 742 (
 743     const dimensioned<scalar>&,
 744     const tmp<faMatrix<Type> >&
 745 );
 746
 747
 748 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 749
 750 } // End namespace Foam
 751
 752 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 753
 754 #ifdef NoRepository
 755 #   include "faMatrix.C"
 756 #endif
 757
 758 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 759
 760 #endif
 761
 762 // ************************************************************************* //