fvMatrix.C
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) 2011-2017 OpenFOAM Foundation
6  \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd.
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 \*---------------------------------------------------------------------------*/
25 
26 #include "volFields.H"
27 #include "surfaceFields.H"
30 #include "coupledFvPatchFields.H"
31 #include "UIndirectList.H"
32 
33 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
34 
35 template<class Type>
36 template<class Type2>
38 (
39  const labelUList& addr,
40  const Field<Type2>& pf,
41  Field<Type2>& intf
42 ) const
43 {
44  if (addr.size() != pf.size())
45  {
47  << "sizes of addressing and field are different"
48  << abort(FatalError);
49  }
50 
51  forAll(addr, facei)
52  {
53  intf[addr[facei]] += pf[facei];
54  }
55 }
56 
57 
58 template<class Type>
59 template<class Type2>
61 (
62  const labelUList& addr,
63  const tmp<Field<Type2>>& tpf,
64  Field<Type2>& intf
65 ) const
66 {
67  addToInternalField(addr, tpf(), intf);
68  tpf.clear();
69 }
70 
71 
72 template<class Type>
73 template<class Type2>
75 (
76  const labelUList& addr,
77  const Field<Type2>& pf,
78  Field<Type2>& intf
79 ) const
80 {
81  if (addr.size() != pf.size())
82  {
84  << "sizes of addressing and field are different"
85  << abort(FatalError);
86  }
87 
88  forAll(addr, facei)
89  {
90  intf[addr[facei]] -= pf[facei];
91  }
92 }
93 
94 
95 template<class Type>
96 template<class Type2>
98 (
99  const labelUList& addr,
100  const tmp<Field<Type2>>& tpf,
101  Field<Type2>& intf
102 ) const
103 {
104  subtractFromInternalField(addr, tpf(), intf);
105  tpf.clear();
106 }
107 
108 
109 template<class Type>
111 (
112  scalarField& diag,
113  const direction solveCmpt
114 ) const
115 {
116  forAll(internalCoeffs_, patchi)
117  {
118  addToInternalField
119  (
120  lduAddr().patchAddr(patchi),
121  internalCoeffs_[patchi].component(solveCmpt),
122  diag
123  );
124  }
125 }
126 
127 
128 template<class Type>
130 {
131  forAll(internalCoeffs_, patchi)
132  {
133  addToInternalField
134  (
135  lduAddr().patchAddr(patchi),
136  cmptAv(internalCoeffs_[patchi]),
137  diag
138  );
139  }
140 }
141 
142 
143 template<class Type>
145 (
146  Field<Type>& source,
147  const bool couples
148 ) const
149 {
150  forAll(psi_.boundaryField(), patchi)
151  {
152  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
153  const Field<Type>& pbc = boundaryCoeffs_[patchi];
154 
155  if (!ptf.coupled())
156  {
157  addToInternalField(lduAddr().patchAddr(patchi), pbc, source);
158  }
159  else if (couples)
160  {
161  const tmp<Field<Type>> tpnf = ptf.patchNeighbourField();
162  const Field<Type>& pnf = tpnf();
163 
164  const labelUList& addr = lduAddr().patchAddr(patchi);
165 
166  forAll(addr, facei)
167  {
168  source[addr[facei]] += cmptMultiply(pbc[facei], pnf[facei]);
169  }
170  }
171  }
172 }
173 
174 
175 template<class Type>
176 template<template<class> class ListType>
178 (
179  const labelUList& cellLabels,
180  const ListType<Type>& values
181 )
182 {
183  const fvMesh& mesh = psi_.mesh();
184 
185  const cellList& cells = mesh.cells();
186  const labelUList& own = mesh.owner();
187  const labelUList& nei = mesh.neighbour();
188 
189  scalarField& Diag = diag();
190  Field<Type>& psi =
191  const_cast
192  <
194  >(psi_).primitiveFieldRef();
195 
196  forAll(cellLabels, i)
197  {
198  const label celli = cellLabels[i];
199  const Type& value = values[i];
200 
201  psi[celli] = value;
202  source_[celli] = value*Diag[celli];
203 
204  if (symmetric() || asymmetric())
205  {
206  const cell& c = cells[celli];
207 
208  forAll(c, j)
209  {
210  const label facei = c[j];
211 
212  if (mesh.isInternalFace(facei))
213  {
214  if (symmetric())
215  {
216  if (celli == own[facei])
217  {
218  source_[nei[facei]] -= upper()[facei]*value;
219  }
220  else
221  {
222  source_[own[facei]] -= upper()[facei]*value;
223  }
224 
225  upper()[facei] = 0.0;
226  }
227  else
228  {
229  if (celli == own[facei])
230  {
231  source_[nei[facei]] -= lower()[facei]*value;
232  }
233  else
234  {
235  source_[own[facei]] -= upper()[facei]*value;
236  }
237 
238  upper()[facei] = 0.0;
239  lower()[facei] = 0.0;
240  }
241  }
242  else
243  {
244  label patchi = mesh.boundaryMesh().whichPatch(facei);
245 
246  if (internalCoeffs_[patchi].size())
247  {
248  label patchFacei =
249  mesh.boundaryMesh()[patchi].whichFace(facei);
250 
251  internalCoeffs_[patchi][patchFacei] =
252  Zero;
253 
254  boundaryCoeffs_[patchi][patchFacei] =
255  Zero;
256  }
257  }
258  }
259  }
260  }
261 }
262 
263 
264 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
265 
266 template<class Type>
268 (
270  const dimensionSet& ds
271 )
272 :
273  lduMatrix(psi.mesh()),
274  psi_(psi),
275  dimensions_(ds),
276  source_(psi.size(), Zero),
277  internalCoeffs_(psi.mesh().boundary().size()),
278  boundaryCoeffs_(psi.mesh().boundary().size()),
279  faceFluxCorrectionPtr_(nullptr)
280 {
281  if (debug)
282  {
284  << "Constructing fvMatrix<Type> for field " << psi_.name() << endl;
285  }
286 
287  // Initialise coupling coefficients
288  forAll(psi.mesh().boundary(), patchi)
289  {
290  internalCoeffs_.set
291  (
292  patchi,
293  new Field<Type>
294  (
295  psi.mesh().boundary()[patchi].size(),
296  Zero
297  )
298  );
299 
300  boundaryCoeffs_.set
301  (
302  patchi,
303  new Field<Type>
304  (
305  psi.mesh().boundary()[patchi].size(),
306  Zero
307  )
308  );
309  }
310 
311  // Update the boundary coefficients of psi without changing its event No.
314 
315  label currentStatePsi = psiRef.eventNo();
316  psiRef.boundaryFieldRef().updateCoeffs();
317  psiRef.eventNo() = currentStatePsi;
318 }
319 
320 
321 template<class Type>
323 :
324  refCount(),
325  lduMatrix(fvm),
326  psi_(fvm.psi_),
327  dimensions_(fvm.dimensions_),
328  source_(fvm.source_),
329  internalCoeffs_(fvm.internalCoeffs_),
330  boundaryCoeffs_(fvm.boundaryCoeffs_),
331  faceFluxCorrectionPtr_(nullptr)
332 {
333  if (debug)
334  {
336  << "Copying fvMatrix<Type> for field " << psi_.name() << endl;
337  }
338 
339  if (fvm.faceFluxCorrectionPtr_)
340  {
341  faceFluxCorrectionPtr_ = new
343  (
344  *(fvm.faceFluxCorrectionPtr_)
345  );
346  }
347 }
348 
349 
350 #ifndef NoConstructFromTmp
351 template<class Type>
353 :
354  lduMatrix
355  (
356  const_cast<fvMatrix<Type>&>(tfvm()),
357  tfvm.isTmp()
358  ),
359  psi_(tfvm().psi_),
360  dimensions_(tfvm().dimensions_),
361  source_
362  (
363  const_cast<fvMatrix<Type>&>(tfvm()).source_,
364  tfvm.isTmp()
365  ),
366  internalCoeffs_
367  (
368  const_cast<fvMatrix<Type>&>(tfvm()).internalCoeffs_,
369  tfvm.isTmp()
370  ),
371  boundaryCoeffs_
372  (
373  const_cast<fvMatrix<Type>&>(tfvm()).boundaryCoeffs_,
374  tfvm.isTmp()
375  ),
376  faceFluxCorrectionPtr_(nullptr)
377 {
378  if (debug)
379  {
381  << "Copying fvMatrix<Type> for field " << psi_.name() << endl;
382  }
383 
384  if (tfvm().faceFluxCorrectionPtr_)
385  {
386  if (tfvm.isTmp())
387  {
388  faceFluxCorrectionPtr_ = tfvm().faceFluxCorrectionPtr_;
389  tfvm().faceFluxCorrectionPtr_ = nullptr;
390  }
391  else
392  {
393  faceFluxCorrectionPtr_ = new
395  (
396  *(tfvm().faceFluxCorrectionPtr_)
397  );
398  }
399  }
400 
401  tfvm.clear();
402 }
403 #endif
404 
405 
406 template<class Type>
408 (
410  Istream& is
411 )
412 :
413  lduMatrix(psi.mesh()),
414  psi_(psi),
415  dimensions_(is),
416  source_(is),
417  internalCoeffs_(psi.mesh().boundary().size()),
418  boundaryCoeffs_(psi.mesh().boundary().size()),
419  faceFluxCorrectionPtr_(nullptr)
420 {
421  if (debug)
422  {
424  << "Constructing fvMatrix<Type> for field " << psi_.name() << endl;
425  }
426 
427  // Initialise coupling coefficients
428  forAll(psi.mesh().boundary(), patchi)
429  {
430  internalCoeffs_.set
431  (
432  patchi,
433  new Field<Type>
434  (
435  psi.mesh().boundary()[patchi].size(),
436  Zero
437  )
438  );
439 
440  boundaryCoeffs_.set
441  (
442  patchi,
443  new Field<Type>
444  (
445  psi.mesh().boundary()[patchi].size(),
446  Zero
447  )
448  );
449  }
450 
451 }
452 
453 
454 template<class Type>
456 {
457  return tmp<fvMatrix<Type>>
458  (
459  new fvMatrix<Type>(*this)
460  );
461 }
462 
463 
464 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
465 
466 template<class Type>
468 {
469  if (debug)
470  {
472  << "Destroying fvMatrix<Type> for field " << psi_.name() << endl;
473  }
474 
475  if (faceFluxCorrectionPtr_)
476  {
477  delete faceFluxCorrectionPtr_;
478  }
479 }
480 
481 
482 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
483 
484 template<class Type>
486 (
487  const labelUList& cellLabels,
488  const UList<Type>& values
489 )
490 {
491  this->setValuesFromList(cellLabels, values);
492 }
493 
494 
495 template<class Type>
497 (
498  const labelUList& cellLabels,
500 )
501 {
502  this->setValuesFromList(cellLabels, values);
503 }
504 
505 
506 template<class Type>
508 (
509  const label celli,
510  const Type& value,
511  const bool forceReference
512 )
513 {
514  if ((forceReference || psi_.needReference()) && celli >= 0)
515  {
516  source()[celli] += diag()[celli]*value;
517  diag()[celli] += diag()[celli];
518  }
519 }
520 
521 
522 template<class Type>
524 (
525  const labelUList& cellLabels,
526  const Type& value,
527  const bool forceReference
528 )
529 {
530  const bool needRef = (forceReference || psi_.needReference());
531 
532  if (needRef)
533  {
534  forAll(cellLabels, celli)
535  {
536  const label cellId = cellLabels[celli];
537  if (cellId >= 0)
538  {
539  source()[cellId] += diag()[cellId]*value;
540  diag()[cellId] += diag()[cellId];
541  }
542  }
543  }
544 }
545 
546 
547 template<class Type>
549 (
550  const labelUList& cellLabels,
551  const UList<Type>& values,
552  const bool forceReference
553 )
554 {
555  const bool needRef = (forceReference || psi_.needReference());
556 
557  if (needRef)
558  {
559  forAll(cellLabels, celli)
560  {
561  const label cellId = cellLabels[celli];
562  if (cellId >= 0)
563  {
564  source()[cellId] += diag()[cellId]*values[celli];
565  diag()[cellId] += diag()[cellId];
566  }
567  }
568  }
569 }
570 
571 
572 template<class Type>
574 {
575  if (alpha <= 0)
576  {
577  return;
578  }
579 
580  if (debug)
581  {
583  << "Relaxing " << psi_.name() << " by " << alpha << endl;
584  }
585 
586  Field<Type>& S = source();
587  scalarField& D = diag();
588 
589  // Store the current unrelaxed diagonal for use in updating the source
590  scalarField D0(D);
591 
592  // Calculate the sum-mag off-diagonal from the interior faces
593  scalarField sumOff(D.size(), 0.0);
594  sumMagOffDiag(sumOff);
595 
596  // Handle the boundary contributions to the diagonal
597  forAll(psi_.boundaryField(), patchi)
598  {
599  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
600 
601  if (ptf.size())
602  {
603  const labelUList& pa = lduAddr().patchAddr(patchi);
604  Field<Type>& iCoeffs = internalCoeffs_[patchi];
605 
606  if (ptf.coupled())
607  {
608  const Field<Type>& pCoeffs = boundaryCoeffs_[patchi];
609 
610  // For coupled boundaries add the diagonal and
611  // off-diagonal contributions
612  forAll(pa, face)
613  {
614  D[pa[face]] += component(iCoeffs[face], 0);
615  sumOff[pa[face]] += mag(component(pCoeffs[face], 0));
616  }
617  }
618  else
619  {
620  // For non-coupled boundaries add the maximum magnitude diagonal
621  // contribution to ensure stability
622  forAll(pa, face)
623  {
624  D[pa[face]] += cmptMax(cmptMag(iCoeffs[face]));
625  }
626  }
627  }
628  }
629 
630 
631  if (debug)
632  {
633  // Calculate amount of non-dominance.
634  label nNon = 0;
635  scalar maxNon = 0.0;
636  scalar sumNon = 0.0;
637  forAll(D, celli)
638  {
639  scalar d = (sumOff[celli] - D[celli])/mag(D[celli]);
640 
641  if (d > 0)
642  {
643  nNon++;
644  maxNon = max(maxNon, d);
645  sumNon += d;
646  }
647  }
648 
649  reduce(nNon, sumOp<label>(), UPstream::msgType(), psi_.mesh().comm());
650  reduce
651  (
652  maxNon,
653  maxOp<scalar>(),
655  psi_.mesh().comm()
656  );
657  reduce
658  (
659  sumNon,
660  sumOp<scalar>(),
662  psi_.mesh().comm()
663  );
664  sumNon /= returnReduce
665  (
666  D.size(),
667  sumOp<label>(),
669  psi_.mesh().comm()
670  );
671 
673  << "Matrix dominance test for " << psi_.name() << nl
674  << " number of non-dominant cells : " << nNon << nl
675  << " maximum relative non-dominance : " << maxNon << nl
676  << " average relative non-dominance : " << sumNon << nl
677  << endl;
678  }
679 
680 
681  // Ensure the matrix is diagonally dominant...
682  // Assumes that the central coefficient is positive and ensures it is
683  forAll(D, celli)
684  {
685  D[celli] = max(mag(D[celli]), sumOff[celli]);
686  }
687 
688  // ... then relax
689  D /= alpha;
690 
691  // Now remove the diagonal contribution from coupled boundaries
692  forAll(psi_.boundaryField(), patchi)
693  {
694  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
695 
696  if (ptf.size())
697  {
698  const labelUList& pa = lduAddr().patchAddr(patchi);
699  Field<Type>& iCoeffs = internalCoeffs_[patchi];
700 
701  if (ptf.coupled())
702  {
703  forAll(pa, face)
704  {
705  D[pa[face]] -= component(iCoeffs[face], 0);
706  }
707  }
708  else
709  {
710  forAll(pa, face)
711  {
712  D[pa[face]] -= cmptMin(iCoeffs[face]);
713  }
714  }
715  }
716  }
717 
718  // Finally add the relaxation contribution to the source.
719  S += (D - D0)*psi_.primitiveField();
720 }
721 
722 
723 template<class Type>
725 {
726  word name = psi_.select
727  (
728  psi_.mesh().data::template lookupOrDefault<bool>
729  ("finalIteration", false)
730  );
731 
732  if (psi_.mesh().relaxEquation(name))
733  {
734  relax(psi_.mesh().equationRelaxationFactor(name));
735  }
736 }
737 
738 
739 template<class Type>
741 (
743  Boundary& bFields
744 )
745 {
746  forAll(bFields, patchi)
747  {
748  bFields[patchi].manipulateMatrix(*this);
749  }
750 }
751 
752 
753 template<class Type>
755 {
756  tmp<scalarField> tdiag(new scalarField(diag()));
757  addCmptAvBoundaryDiag(tdiag.ref());
758  return tdiag;
759 }
760 
761 
762 template<class Type>
764 {
766 
767  forAll(psi_.boundaryField(), patchi)
768  {
769  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
770 
771  if (!ptf.coupled() && ptf.size())
772  {
773  addToInternalField
774  (
775  lduAddr().patchAddr(patchi),
776  internalCoeffs_[patchi],
777  tdiag.ref()
778  );
779  }
780  }
781 
782  return tdiag;
783 }
784 
785 
786 template<class Type>
788 {
789  tmp<volScalarField> tAphi
790  (
791  new volScalarField
792  (
793  IOobject
794  (
795  "A("+psi_.name()+')',
796  psi_.instance(),
797  psi_.mesh(),
800  ),
801  psi_.mesh(),
802  dimensions_/psi_.dimensions()/dimVol,
803  extrapolatedCalculatedFvPatchScalarField::typeName
804  )
805  );
806 
807  tAphi.ref().primitiveFieldRef() = D()/psi_.mesh().V();
808  tAphi.ref().correctBoundaryConditions();
809 
810  return tAphi;
811 }
812 
813 
814 template<class Type>
817 {
819  (
821  (
822  IOobject
823  (
824  "H("+psi_.name()+')',
825  psi_.instance(),
826  psi_.mesh(),
829  ),
830  psi_.mesh(),
831  dimensions_/dimVol,
832  extrapolatedCalculatedFvPatchScalarField::typeName
833  )
834  );
836 
837  // Loop over field components
838  for (direction cmpt=0; cmpt<Type::nComponents; cmpt++)
839  {
840  scalarField psiCmpt(psi_.primitiveField().component(cmpt));
841 
842  scalarField boundaryDiagCmpt(psi_.size(), 0.0);
843  addBoundaryDiag(boundaryDiagCmpt, cmpt);
844  boundaryDiagCmpt.negate();
845  addCmptAvBoundaryDiag(boundaryDiagCmpt);
846 
847  Hphi.primitiveFieldRef().replace(cmpt, boundaryDiagCmpt*psiCmpt);
848  }
849 
850  Hphi.primitiveFieldRef() += lduMatrix::H(psi_.primitiveField()) + source_;
851  addBoundarySource(Hphi.primitiveFieldRef());
852 
853  Hphi.primitiveFieldRef() /= psi_.mesh().V();
855 
856  typename Type::labelType validComponents
857  (
858  psi_.mesh().template validComponents<Type>()
859  );
860 
861  for (direction cmpt=0; cmpt<Type::nComponents; cmpt++)
862  {
863  if (validComponents[cmpt] == -1)
864  {
865  Hphi.replace
866  (
867  cmpt,
869  );
870  }
871  }
872 
873  return tHphi;
874 }
875 
876 
877 template<class Type>
879 {
881  (
882  new volScalarField
883  (
884  IOobject
885  (
886  "H(1)",
887  psi_.instance(),
888  psi_.mesh(),
891  ),
892  psi_.mesh(),
893  dimensions_/(dimVol*psi_.dimensions()),
894  extrapolatedCalculatedFvPatchScalarField::typeName
895  )
896  );
897  volScalarField& H1_ = tH1.ref();
898 
899  H1_.primitiveFieldRef() = lduMatrix::H1();
900 
901  forAll(psi_.boundaryField(), patchi)
902  {
903  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
904 
905  if (ptf.coupled() && ptf.size())
906  {
907  addToInternalField
908  (
909  lduAddr().patchAddr(patchi),
910  boundaryCoeffs_[patchi].component(0),
911  H1_
912  );
913  }
914  }
915 
916  H1_.primitiveFieldRef() /= psi_.mesh().V();
917  H1_.correctBoundaryConditions();
918 
919  return tH1;
920 }
921 
922 
923 template<class Type>
926 flux() const
927 {
928  if (!psi_.mesh().fluxRequired(psi_.name()))
929  {
931  << "flux requested but " << psi_.name()
932  << " not specified in the fluxRequired sub-dictionary"
933  " of fvSchemes."
934  << abort(FatalError);
935  }
936 
937  // construct GeometricField<Type, fvsPatchField, surfaceMesh>
939  (
941  (
942  IOobject
943  (
944  "flux("+psi_.name()+')',
945  psi_.instance(),
946  psi_.mesh(),
949  ),
950  psi_.mesh(),
951  dimensions()
952  )
953  );
955  tfieldFlux.ref();
956 
957  fieldFlux.setOriented();
958 
959  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
960  {
961  fieldFlux.primitiveFieldRef().replace
962  (
963  cmpt,
964  lduMatrix::faceH(psi_.primitiveField().component(cmpt))
965  );
966  }
967 
968  FieldField<Field, Type> InternalContrib = internalCoeffs_;
969 
970  forAll(InternalContrib, patchi)
971  {
972  InternalContrib[patchi] =
974  (
975  InternalContrib[patchi],
976  psi_.boundaryField()[patchi].patchInternalField()
977  );
978  }
979 
980  FieldField<Field, Type> NeighbourContrib = boundaryCoeffs_;
981 
982  forAll(NeighbourContrib, patchi)
983  {
984  if (psi_.boundaryField()[patchi].coupled())
985  {
986  NeighbourContrib[patchi] =
988  (
989  NeighbourContrib[patchi],
990  psi_.boundaryField()[patchi].patchNeighbourField()
991  );
992  }
993  }
994 
996  Boundary& ffbf = fieldFlux.boundaryFieldRef();
997 
998  forAll(ffbf, patchi)
999  {
1000  ffbf[patchi] = InternalContrib[patchi] - NeighbourContrib[patchi];
1001  }
1002 
1003  if (faceFluxCorrectionPtr_)
1004  {
1005  fieldFlux += *faceFluxCorrectionPtr_;
1006  }
1007 
1008  return tfieldFlux;
1009 }
1010 
1011 
1012 template<class Type>
1014 {
1015  return psi_.mesh().solverDict
1016  (
1017  psi_.select
1018  (
1019  psi_.mesh().data::template lookupOrDefault<bool>
1020  ("finalIteration", false)
1021  )
1022  );
1023 }
1024 
1025 
1026 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1027 
1028 template<class Type>
1030 {
1031  if (this == &fvmv)
1032  {
1034  << "attempted assignment to self"
1035  << abort(FatalError);
1036  }
1037 
1038  if (&psi_ != &(fvmv.psi_))
1039  {
1041  << "different fields"
1042  << abort(FatalError);
1043  }
1044 
1045  dimensions_ = fvmv.dimensions_;
1046  lduMatrix::operator=(fvmv);
1047  source_ = fvmv.source_;
1048  internalCoeffs_ = fvmv.internalCoeffs_;
1049  boundaryCoeffs_ = fvmv.boundaryCoeffs_;
1050 
1051  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
1052  {
1053  *faceFluxCorrectionPtr_ = *fvmv.faceFluxCorrectionPtr_;
1054  }
1055  else if (fvmv.faceFluxCorrectionPtr_)
1056  {
1057  faceFluxCorrectionPtr_ =
1059  (*fvmv.faceFluxCorrectionPtr_);
1060  }
1061 }
1062 
1063 
1064 template<class Type>
1066 {
1067  operator=(tfvmv());
1068  tfvmv.clear();
1069 }
1070 
1071 
1072 template<class Type>
1074 {
1076  source_.negate();
1077  internalCoeffs_.negate();
1078  boundaryCoeffs_.negate();
1079 
1080  if (faceFluxCorrectionPtr_)
1081  {
1082  faceFluxCorrectionPtr_->negate();
1083  }
1084 }
1085 
1086 
1087 template<class Type>
1089 {
1090  checkMethod(*this, fvmv, "+=");
1091 
1092  dimensions_ += fvmv.dimensions_;
1093  lduMatrix::operator+=(fvmv);
1094  source_ += fvmv.source_;
1095  internalCoeffs_ += fvmv.internalCoeffs_;
1096  boundaryCoeffs_ += fvmv.boundaryCoeffs_;
1097 
1098  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
1099  {
1100  *faceFluxCorrectionPtr_ += *fvmv.faceFluxCorrectionPtr_;
1101  }
1102  else if (fvmv.faceFluxCorrectionPtr_)
1103  {
1104  faceFluxCorrectionPtr_ = new
1106  (
1107  *fvmv.faceFluxCorrectionPtr_
1108  );
1109  }
1110 }
1111 
1112 
1113 template<class Type>
1115 {
1116  operator+=(tfvmv());
1117  tfvmv.clear();
1118 }
1119 
1120 
1121 template<class Type>
1123 {
1124  checkMethod(*this, fvmv, "-=");
1125 
1126  dimensions_ -= fvmv.dimensions_;
1127  lduMatrix::operator-=(fvmv);
1128  source_ -= fvmv.source_;
1129  internalCoeffs_ -= fvmv.internalCoeffs_;
1130  boundaryCoeffs_ -= fvmv.boundaryCoeffs_;
1131 
1132  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
1133  {
1134  *faceFluxCorrectionPtr_ -= *fvmv.faceFluxCorrectionPtr_;
1135  }
1136  else if (fvmv.faceFluxCorrectionPtr_)
1137  {
1138  faceFluxCorrectionPtr_ =
1140  (-*fvmv.faceFluxCorrectionPtr_);
1141  }
1142 }
1143 
1144 
1145 template<class Type>
1147 {
1148  operator-=(tfvmv());
1149  tfvmv.clear();
1150 }
1151 
1152 
1153 template<class Type>
1157 )
1158 {
1159  checkMethod(*this, su, "+=");
1160  source() -= su.mesh().V()*su.field();
1161 }
1162 
1163 
1164 template<class Type>
1168 )
1169 {
1170  operator+=(tsu());
1171  tsu.clear();
1172 }
1173 
1174 
1175 template<class Type>
1179 )
1180 {
1181  operator+=(tsu());
1182  tsu.clear();
1183 }
1184 
1185 
1186 template<class Type>
1190 )
1191 {
1192  checkMethod(*this, su, "-=");
1193  source() += su.mesh().V()*su.field();
1194 }
1195 
1196 
1197 template<class Type>
1201 )
1202 {
1203  operator-=(tsu());
1204  tsu.clear();
1205 }
1206 
1207 
1208 template<class Type>
1212 )
1213 {
1214  operator-=(tsu());
1215  tsu.clear();
1216 }
1217 
1218 
1219 template<class Type>
1222  const dimensioned<Type>& su
1223 )
1224 {
1225  source() -= psi().mesh().V()*su;
1226 }
1227 
1228 
1229 template<class Type>
1232  const dimensioned<Type>& su
1233 )
1234 {
1235  source() += psi().mesh().V()*su;
1236 }
1237 
1238 
1239 template<class Type>
1242  const zero&
1243 )
1244 {}
1245 
1246 
1247 template<class Type>
1250  const zero&
1251 )
1252 {}
1253 
1254 
1255 template<class Type>
1258  const volScalarField::Internal& dsf
1259 )
1260 {
1261  dimensions_ *= dsf.dimensions();
1262  lduMatrix::operator*=(dsf.field());
1263  source_ *= dsf.field();
1264 
1265  forAll(boundaryCoeffs_, patchi)
1266  {
1267  scalarField pisf
1268  (
1269  dsf.mesh().boundary()[patchi].patchInternalField(dsf.field())
1270  );
1271 
1272  internalCoeffs_[patchi] *= pisf;
1273  boundaryCoeffs_[patchi] *= pisf;
1274  }
1275 
1276  if (faceFluxCorrectionPtr_)
1277  {
1279  << "cannot scale a matrix containing a faceFluxCorrection"
1280  << abort(FatalError);
1281  }
1282 }
1283 
1284 
1285 template<class Type>
1288  const tmp<volScalarField::Internal>& tdsf
1289 )
1290 {
1291  operator*=(tdsf());
1292  tdsf.clear();
1293 }
1294 
1295 
1296 template<class Type>
1299  const tmp<volScalarField>& tvsf
1300 )
1301 {
1302  operator*=(tvsf());
1303  tvsf.clear();
1304 }
1305 
1306 
1307 template<class Type>
1310  const dimensioned<scalar>& ds
1311 )
1312 {
1313  dimensions_ *= ds.dimensions();
1314  lduMatrix::operator*=(ds.value());
1315  source_ *= ds.value();
1316  internalCoeffs_ *= ds.value();
1317  boundaryCoeffs_ *= ds.value();
1318 
1319  if (faceFluxCorrectionPtr_)
1320  {
1321  *faceFluxCorrectionPtr_ *= ds.value();
1322  }
1323 }
1324 
1325 
1326 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
1327 
1328 template<class Type>
1329 void Foam::checkMethod
1331  const fvMatrix<Type>& fvm1,
1332  const fvMatrix<Type>& fvm2,
1333  const char* op
1334 )
1335 {
1336  if (&fvm1.psi() != &fvm2.psi())
1337  {
1339  << "incompatible fields for operation "
1340  << endl << " "
1341  << "[" << fvm1.psi().name() << "] "
1342  << op
1343  << " [" << fvm2.psi().name() << "]"
1344  << abort(FatalError);
1345  }
1346 
1347  if (dimensionSet::debug && fvm1.dimensions() != fvm2.dimensions())
1348  {
1350  << "incompatible dimensions for operation "
1351  << endl << " "
1352  << "[" << fvm1.psi().name() << fvm1.dimensions()/dimVolume << " ] "
1353  << op
1354  << " [" << fvm2.psi().name() << fvm2.dimensions()/dimVolume << " ]"
1355  << abort(FatalError);
1356  }
1357 }
1358 
1359 
1360 template<class Type>
1361 void Foam::checkMethod
1363  const fvMatrix<Type>& fvm,
1365  const char* op
1366 )
1367 {
1368  if (dimensionSet::debug && fvm.dimensions()/dimVolume != df.dimensions())
1369  {
1371  << endl << " "
1372  << "[" << fvm.psi().name() << fvm.dimensions()/dimVolume << " ] "
1373  << op
1374  << " [" << df.name() << df.dimensions() << " ]"
1375  << abort(FatalError);
1376  }
1377 }
1378 
1379 
1380 template<class Type>
1381 void Foam::checkMethod
1383  const fvMatrix<Type>& fvm,
1384  const dimensioned<Type>& dt,
1385  const char* op
1386 )
1387 {
1388  if (dimensionSet::debug && fvm.dimensions()/dimVolume != dt.dimensions())
1389  {
1391  << "incompatible dimensions for operation "
1392  << endl << " "
1393  << "[" << fvm.psi().name() << fvm.dimensions()/dimVolume << " ] "
1394  << op
1395  << " [" << dt.name() << dt.dimensions() << " ]"
1396  << abort(FatalError);
1397  }
1398 }
1399 
1400 
1401 template<class Type>
1403 (
1404  fvMatrix<Type>& fvm,
1405  const dictionary& solverControls
1406 )
1407 {
1408  return fvm.solve(solverControls);
1409 }
1410 
1411 template<class Type>
1413 (
1414  const tmp<fvMatrix<Type>>& tfvm,
1415  const dictionary& solverControls
1416 )
1417 {
1418  SolverPerformance<Type> solverPerf =
1419  const_cast<fvMatrix<Type>&>(tfvm()).solve(solverControls);
1420 
1421  tfvm.clear();
1422 
1423  return solverPerf;
1424 }
1425 
1426 
1427 template<class Type>
1428 Foam::SolverPerformance<Type> Foam::solve(fvMatrix<Type>& fvm)
1429 {
1430  return fvm.solve();
1431 }
1432 
1433 template<class Type>
1434 Foam::SolverPerformance<Type> Foam::solve(const tmp<fvMatrix<Type>>& tfvm)
1435 {
1436  SolverPerformance<Type> solverPerf =
1437  const_cast<fvMatrix<Type>&>(tfvm()).solve();
1438 
1439  tfvm.clear();
1440 
1441  return solverPerf;
1442 }
1443 
1444 
1445 template<class Type>
1447 (
1448  const fvMatrix<Type>& A
1449 )
1450 {
1451  tmp<Foam::fvMatrix<Type>> tAcorr = A - (A & A.psi());
1452 
1453  // Delete the faceFluxCorrection from the correction matrix
1454  // as it does not have a clear meaning or purpose
1455  deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
1456 
1457  return tAcorr;
1458 }
1459 
1460 
1461 template<class Type>
1463 (
1464  const tmp<fvMatrix<Type>>& tA
1465 )
1466 {
1467  tmp<Foam::fvMatrix<Type>> tAcorr = tA - (tA() & tA().psi());
1468 
1469  // Delete the faceFluxCorrection from the correction matrix
1470  // as it does not have a clear meaning or purpose
1471  deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
1472 
1473  return tAcorr;
1474 }
1475 
1476 
1477 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
1478 
1479 template<class Type>
1480 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1481 (
1482  const fvMatrix<Type>& A,
1483  const fvMatrix<Type>& B
1484 )
1485 {
1486  checkMethod(A, B, "==");
1487  return (A - B);
1488 }
1489 
1490 template<class Type>
1491 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1492 (
1493  const tmp<fvMatrix<Type>>& tA,
1494  const fvMatrix<Type>& B
1495 )
1496 {
1497  checkMethod(tA(), B, "==");
1498  return (tA - B);
1499 }
1500 
1501 template<class Type>
1502 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1503 (
1504  const fvMatrix<Type>& A,
1505  const tmp<fvMatrix<Type>>& tB
1506 )
1507 {
1508  checkMethod(A, tB(), "==");
1509  return (A - tB);
1510 }
1511 
1512 template<class Type>
1513 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1514 (
1515  const tmp<fvMatrix<Type>>& tA,
1516  const tmp<fvMatrix<Type>>& tB
1517 )
1518 {
1519  checkMethod(tA(), tB(), "==");
1520  return (tA - tB);
1521 }
1522 
1523 template<class Type>
1524 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1525 (
1526  const fvMatrix<Type>& A,
1527  const DimensionedField<Type, volMesh>& su
1528 )
1529 {
1530  checkMethod(A, su, "==");
1531  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1532  tC.ref().source() += su.mesh().V()*su.field();
1533  return tC;
1534 }
1535 
1536 template<class Type>
1537 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1538 (
1539  const fvMatrix<Type>& A,
1540  const tmp<DimensionedField<Type, volMesh>>& tsu
1541 )
1542 {
1543  checkMethod(A, tsu(), "==");
1544  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1545  tC.ref().source() += tsu().mesh().V()*tsu().field();
1546  tsu.clear();
1547  return tC;
1548 }
1549 
1550 template<class Type>
1551 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1552 (
1553  const fvMatrix<Type>& A,
1554  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
1555 )
1556 {
1557  checkMethod(A, tsu(), "==");
1558  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1559  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1560  tsu.clear();
1561  return tC;
1562 }
1563 
1564 template<class Type>
1565 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1566 (
1567  const tmp<fvMatrix<Type>>& tA,
1568  const DimensionedField<Type, volMesh>& su
1569 )
1570 {
1571  checkMethod(tA(), su, "==");
1572  tmp<fvMatrix<Type>> tC(tA.ptr());
1573  tC.ref().source() += su.mesh().V()*su.field();
1574  return tC;
1575 }
1576 
1577 template<class Type>
1578 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1579 (
1580  const tmp<fvMatrix<Type>>& tA,
1581  const tmp<DimensionedField<Type, volMesh>>& tsu
1582 )
1583 {
1584  checkMethod(tA(), tsu(), "==");
1585  tmp<fvMatrix<Type>> tC(tA.ptr());
1586  tC.ref().source() += tsu().mesh().V()*tsu().field();
1587  tsu.clear();
1588  return tC;
1589 }
1590 
1591 template<class Type>
1592 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1593 (
1594  const tmp<fvMatrix<Type>>& tA,
1595  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
1596 )
1597 {
1598  checkMethod(tA(), tsu(), "==");
1599  tmp<fvMatrix<Type>> tC(tA.ptr());
1600  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1601  tsu.clear();
1602  return tC;
1603 }
1604 
1605 template<class Type>
1606 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1607 (
1608  const fvMatrix<Type>& A,
1609  const dimensioned<Type>& su
1610 )
1611 {
1612  checkMethod(A, su, "==");
1613  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1614  tC.ref().source() += A.psi().mesh().V()*su.value();
1615  return tC;
1616 }
1617 
1618 template<class Type>
1619 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1620 (
1621  const tmp<fvMatrix<Type>>& tA,
1622  const dimensioned<Type>& su
1623 )
1624 {
1625  checkMethod(tA(), su, "==");
1626  tmp<fvMatrix<Type>> tC(tA.ptr());
1627  tC.ref().source() += tC().psi().mesh().V()*su.value();
1628  return tC;
1629 }
1630 
1631 template<class Type>
1632 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1633 (
1634  const fvMatrix<Type>& A,
1635  const zero&
1636 )
1637 {
1638  return A;
1639 }
1640 
1641 
1642 template<class Type>
1643 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1644 (
1645  const tmp<fvMatrix<Type>>& tA,
1646  const zero&
1647 )
1648 {
1649  return tA;
1650 }
1651 
1652 
1653 template<class Type>
1654 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1655 (
1656  const fvMatrix<Type>& A
1657 )
1658 {
1659  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1660  tC.ref().negate();
1661  return tC;
1662 }
1663 
1664 template<class Type>
1665 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1666 (
1667  const tmp<fvMatrix<Type>>& tA
1668 )
1669 {
1670  tmp<fvMatrix<Type>> tC(tA.ptr());
1671  tC.ref().negate();
1672  return tC;
1673 }
1674 
1675 
1676 template<class Type>
1677 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1678 (
1679  const fvMatrix<Type>& A,
1680  const fvMatrix<Type>& B
1681 )
1682 {
1683  checkMethod(A, B, "+");
1684  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1685  tC.ref() += B;
1686  return tC;
1687 }
1688 
1689 template<class Type>
1690 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1691 (
1692  const tmp<fvMatrix<Type>>& tA,
1693  const fvMatrix<Type>& B
1694 )
1695 {
1696  checkMethod(tA(), B, "+");
1697  tmp<fvMatrix<Type>> tC(tA.ptr());
1698  tC.ref() += B;
1699  return tC;
1700 }
1701 
1702 template<class Type>
1703 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1704 (
1705  const fvMatrix<Type>& A,
1706  const tmp<fvMatrix<Type>>& tB
1707 )
1708 {
1709  checkMethod(A, tB(), "+");
1710  tmp<fvMatrix<Type>> tC(tB.ptr());
1711  tC.ref() += A;
1712  return tC;
1713 }
1714 
1715 template<class Type>
1716 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1717 (
1718  const tmp<fvMatrix<Type>>& tA,
1719  const tmp<fvMatrix<Type>>& tB
1720 )
1721 {
1722  checkMethod(tA(), tB(), "+");
1723  tmp<fvMatrix<Type>> tC(tA.ptr());
1724  tC.ref() += tB();
1725  tB.clear();
1726  return tC;
1727 }
1728 
1729 template<class Type>
1730 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1731 (
1732  const fvMatrix<Type>& A,
1733  const DimensionedField<Type, volMesh>& su
1734 )
1735 {
1736  checkMethod(A, su, "+");
1737  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1738  tC.ref().source() -= su.mesh().V()*su.field();
1739  return tC;
1740 }
1741 
1742 template<class Type>
1743 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1744 (
1745  const fvMatrix<Type>& A,
1746  const tmp<DimensionedField<Type, volMesh>>& tsu
1747 )
1748 {
1749  checkMethod(A, tsu(), "+");
1750  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1751  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1752  tsu.clear();
1753  return tC;
1754 }
1755 
1756 template<class Type>
1757 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1758 (
1759  const fvMatrix<Type>& A,
1760  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
1761 )
1762 {
1763  checkMethod(A, tsu(), "+");
1764  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1765  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1766  tsu.clear();
1767  return tC;
1768 }
1769 
1770 template<class Type>
1771 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1772 (
1773  const tmp<fvMatrix<Type>>& tA,
1774  const DimensionedField<Type, volMesh>& su
1775 )
1776 {
1777  checkMethod(tA(), su, "+");
1778  tmp<fvMatrix<Type>> tC(tA.ptr());
1779  tC.ref().source() -= su.mesh().V()*su.field();
1780  return tC;
1781 }
1782 
1783 template<class Type>
1784 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1785 (
1786  const tmp<fvMatrix<Type>>& tA,
1787  const tmp<DimensionedField<Type, volMesh>>& tsu
1788 )
1789 {
1790  checkMethod(tA(), tsu(), "+");
1791  tmp<fvMatrix<Type>> tC(tA.ptr());
1792  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1793  tsu.clear();
1794  return tC;
1795 }
1796 
1797 template<class Type>
1798 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1799 (
1800  const tmp<fvMatrix<Type>>& tA,
1801  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
1802 )
1803 {
1804  checkMethod(tA(), tsu(), "+");
1805  tmp<fvMatrix<Type>> tC(tA.ptr());
1806  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1807  tsu.clear();
1808  return tC;
1809 }
1810 
1811 template<class Type>
1812 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1813 (
1814  const DimensionedField<Type, volMesh>& su,
1815  const fvMatrix<Type>& A
1816 )
1817 {
1818  checkMethod(A, su, "+");
1819  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1820  tC.ref().source() -= su.mesh().V()*su.field();
1821  return tC;
1822 }
1823 
1824 template<class Type>
1825 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1826 (
1827  const tmp<DimensionedField<Type, volMesh>>& tsu,
1828  const fvMatrix<Type>& A
1829 )
1830 {
1831  checkMethod(A, tsu(), "+");
1832  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1833  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1834  tsu.clear();
1835  return tC;
1836 }
1837 
1838 template<class Type>
1839 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1840 (
1841  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu,
1842  const fvMatrix<Type>& A
1843 )
1844 {
1845  checkMethod(A, tsu(), "+");
1846  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1847  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1848  tsu.clear();
1849  return tC;
1850 }
1851 
1852 template<class Type>
1853 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1854 (
1855  const DimensionedField<Type, volMesh>& su,
1856  const tmp<fvMatrix<Type>>& tA
1857 )
1858 {
1859  checkMethod(tA(), su, "+");
1860  tmp<fvMatrix<Type>> tC(tA.ptr());
1861  tC.ref().source() -= su.mesh().V()*su.field();
1862  return tC;
1863 }
1864 
1865 template<class Type>
1866 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1867 (
1868  const tmp<DimensionedField<Type, volMesh>>& tsu,
1869  const tmp<fvMatrix<Type>>& tA
1870 )
1871 {
1872  checkMethod(tA(), tsu(), "+");
1873  tmp<fvMatrix<Type>> tC(tA.ptr());
1874  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1875  tsu.clear();
1876  return tC;
1877 }
1878 
1879 template<class Type>
1880 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1881 (
1882  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu,
1883  const tmp<fvMatrix<Type>>& tA
1884 )
1885 {
1886  checkMethod(tA(), tsu(), "+");
1887  tmp<fvMatrix<Type>> tC(tA.ptr());
1888  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1889  tsu.clear();
1890  return tC;
1891 }
1892 
1893 
1894 template<class Type>
1895 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1896 (
1897  const fvMatrix<Type>& A,
1898  const fvMatrix<Type>& B
1899 )
1900 {
1901  checkMethod(A, B, "-");
1902  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1903  tC.ref() -= B;
1904  return tC;
1905 }
1906 
1907 template<class Type>
1908 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1909 (
1910  const tmp<fvMatrix<Type>>& tA,
1911  const fvMatrix<Type>& B
1912 )
1913 {
1914  checkMethod(tA(), B, "-");
1915  tmp<fvMatrix<Type>> tC(tA.ptr());
1916  tC.ref() -= B;
1917  return tC;
1918 }
1919 
1920 template<class Type>
1921 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1922 (
1923  const fvMatrix<Type>& A,
1924  const tmp<fvMatrix<Type>>& tB
1925 )
1926 {
1927  checkMethod(A, tB(), "-");
1928  tmp<fvMatrix<Type>> tC(tB.ptr());
1929  tC.ref() -= A;
1930  tC.ref().negate();
1931  return tC;
1932 }
1933 
1934 template<class Type>
1935 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1936 (
1937  const tmp<fvMatrix<Type>>& tA,
1938  const tmp<fvMatrix<Type>>& tB
1939 )
1940 {
1941  checkMethod(tA(), tB(), "-");
1942  tmp<fvMatrix<Type>> tC(tA.ptr());
1943  tC.ref() -= tB();
1944  tB.clear();
1945  return tC;
1946 }
1947 
1948 template<class Type>
1949 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1950 (
1951  const fvMatrix<Type>& A,
1952  const DimensionedField<Type, volMesh>& su
1953 )
1954 {
1955  checkMethod(A, su, "-");
1956  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1957  tC.ref().source() += su.mesh().V()*su.field();
1958  return tC;
1959 }
1960 
1961 template<class Type>
1962 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1963 (
1964  const fvMatrix<Type>& A,
1965  const tmp<DimensionedField<Type, volMesh>>& tsu
1966 )
1967 {
1968  checkMethod(A, tsu(), "-");
1969  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1970  tC.ref().source() += tsu().mesh().V()*tsu().field();
1971  tsu.clear();
1972  return tC;
1973 }
1974 
1975 template<class Type>
1976 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1977 (
1978  const fvMatrix<Type>& A,
1979  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
1980 )
1981 {
1982  checkMethod(A, tsu(), "-");
1983  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1984  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1985  tsu.clear();
1986  return tC;
1987 }
1988 
1989 template<class Type>
1990 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1991 (
1992  const tmp<fvMatrix<Type>>& tA,
1993  const DimensionedField<Type, volMesh>& su
1994 )
1995 {
1996  checkMethod(tA(), su, "-");
1997  tmp<fvMatrix<Type>> tC(tA.ptr());
1998  tC.ref().source() += su.mesh().V()*su.field();
1999  return tC;
2000 }
2001 
2002 template<class Type>
2003 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2004 (
2005  const tmp<fvMatrix<Type>>& tA,
2006  const tmp<DimensionedField<Type, volMesh>>& tsu
2007 )
2008 {
2009  checkMethod(tA(), tsu(), "-");
2010  tmp<fvMatrix<Type>> tC(tA.ptr());
2011  tC.ref().source() += tsu().mesh().V()*tsu().field();
2012  tsu.clear();
2013  return tC;
2014 }
2015 
2016 template<class Type>
2017 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2018 (
2019  const tmp<fvMatrix<Type>>& tA,
2020  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu
2021 )
2022 {
2023  checkMethod(tA(), tsu(), "-");
2024  tmp<fvMatrix<Type>> tC(tA.ptr());
2025  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
2026  tsu.clear();
2027  return tC;
2028 }
2029 
2030 template<class Type>
2031 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2032 (
2033  const DimensionedField<Type, volMesh>& su,
2034  const fvMatrix<Type>& A
2035 )
2036 {
2037  checkMethod(A, su, "-");
2038  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2039  tC.ref().negate();
2040  tC.ref().source() -= su.mesh().V()*su.field();
2041  return tC;
2042 }
2043 
2044 template<class Type>
2045 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2046 (
2047  const tmp<DimensionedField<Type, volMesh>>& tsu,
2048  const fvMatrix<Type>& A
2049 )
2050 {
2051  checkMethod(A, tsu(), "-");
2052  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2053  tC.ref().negate();
2054  tC.ref().source() -= tsu().mesh().V()*tsu().field();
2055  tsu.clear();
2056  return tC;
2057 }
2058 
2059 template<class Type>
2060 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2061 (
2062  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu,
2063  const fvMatrix<Type>& A
2064 )
2065 {
2066  checkMethod(A, tsu(), "-");
2067  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2068  tC.ref().negate();
2069  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
2070  tsu.clear();
2071  return tC;
2072 }
2073 
2074 template<class Type>
2075 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2076 (
2077  const DimensionedField<Type, volMesh>& su,
2078  const tmp<fvMatrix<Type>>& tA
2079 )
2080 {
2081  checkMethod(tA(), su, "-");
2082  tmp<fvMatrix<Type>> tC(tA.ptr());
2083  tC.ref().negate();
2084  tC.ref().source() -= su.mesh().V()*su.field();
2085  return tC;
2086 }
2087 
2088 template<class Type>
2089 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2090 (
2091  const tmp<DimensionedField<Type, volMesh>>& tsu,
2092  const tmp<fvMatrix<Type>>& tA
2093 )
2094 {
2095  checkMethod(tA(), tsu(), "-");
2096  tmp<fvMatrix<Type>> tC(tA.ptr());
2097  tC.ref().negate();
2098  tC.ref().source() -= tsu().mesh().V()*tsu().field();
2099  tsu.clear();
2100  return tC;
2101 }
2102 
2103 template<class Type>
2104 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2105 (
2106  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tsu,
2107  const tmp<fvMatrix<Type>>& tA
2108 )
2109 {
2110  checkMethod(tA(), tsu(), "-");
2111  tmp<fvMatrix<Type>> tC(tA.ptr());
2112  tC.ref().negate();
2113  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
2114  tsu.clear();
2115  return tC;
2116 }
2117 
2118 template<class Type>
2119 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2120 (
2121  const fvMatrix<Type>& A,
2122  const dimensioned<Type>& su
2123 )
2124 {
2125  checkMethod(A, su, "+");
2126  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2127  tC.ref().source() -= su.value()*A.psi().mesh().V();
2128  return tC;
2129 }
2130 
2131 template<class Type>
2132 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2133 (
2134  const tmp<fvMatrix<Type>>& tA,
2135  const dimensioned<Type>& su
2136 )
2137 {
2138  checkMethod(tA(), su, "+");
2139  tmp<fvMatrix<Type>> tC(tA.ptr());
2140  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2141  return tC;
2142 }
2143 
2144 template<class Type>
2145 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2146 (
2147  const dimensioned<Type>& su,
2148  const fvMatrix<Type>& A
2149 )
2150 {
2151  checkMethod(A, su, "+");
2152  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2153  tC.ref().source() -= su.value()*A.psi().mesh().V();
2154  return tC;
2155 }
2156 
2157 template<class Type>
2158 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2159 (
2160  const dimensioned<Type>& su,
2161  const tmp<fvMatrix<Type>>& tA
2162 )
2163 {
2164  checkMethod(tA(), su, "+");
2165  tmp<fvMatrix<Type>> tC(tA.ptr());
2166  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2167  return tC;
2168 }
2169 
2170 template<class Type>
2171 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2172 (
2173  const fvMatrix<Type>& A,
2174  const dimensioned<Type>& su
2175 )
2176 {
2177  checkMethod(A, su, "-");
2178  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2179  tC.ref().source() += su.value()*tC().psi().mesh().V();
2180  return tC;
2181 }
2182 
2183 template<class Type>
2184 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2185 (
2186  const tmp<fvMatrix<Type>>& tA,
2187  const dimensioned<Type>& su
2188 )
2189 {
2190  checkMethod(tA(), su, "-");
2191  tmp<fvMatrix<Type>> tC(tA.ptr());
2192  tC.ref().source() += su.value()*tC().psi().mesh().V();
2193  return tC;
2194 }
2195 
2196 template<class Type>
2197 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2198 (
2199  const dimensioned<Type>& su,
2200  const fvMatrix<Type>& A
2201 )
2202 {
2203  checkMethod(A, su, "-");
2204  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2205  tC.ref().negate();
2206  tC.ref().source() -= su.value()*A.psi().mesh().V();
2207  return tC;
2208 }
2209 
2210 template<class Type>
2211 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2212 (
2213  const dimensioned<Type>& su,
2214  const tmp<fvMatrix<Type>>& tA
2215 )
2216 {
2217  checkMethod(tA(), su, "-");
2218  tmp<fvMatrix<Type>> tC(tA.ptr());
2219  tC.ref().negate();
2220  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2221  return tC;
2222 }
2223 
2224 
2225 template<class Type>
2226 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2227 (
2228  const volScalarField::Internal& dsf,
2229  const fvMatrix<Type>& A
2230 )
2231 {
2232  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2233  tC.ref() *= dsf;
2234  return tC;
2235 }
2236 
2237 template<class Type>
2238 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2239 (
2240  const tmp<volScalarField::Internal>& tdsf,
2241  const fvMatrix<Type>& A
2242 )
2243 {
2244  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2245  tC.ref() *= tdsf;
2246  return tC;
2247 }
2248 
2249 template<class Type>
2250 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2251 (
2252  const tmp<volScalarField>& tvsf,
2253  const fvMatrix<Type>& A
2254 )
2255 {
2256  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2257  tC.ref() *= tvsf;
2258  return tC;
2259 }
2260 
2261 template<class Type>
2262 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2263 (
2264  const volScalarField::Internal& dsf,
2265  const tmp<fvMatrix<Type>>& tA
2266 )
2267 {
2268  tmp<fvMatrix<Type>> tC(tA.ptr());
2269  tC.ref() *= dsf;
2270  return tC;
2271 }
2272 
2273 template<class Type>
2274 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2275 (
2276  const tmp<volScalarField::Internal>& tdsf,
2277  const tmp<fvMatrix<Type>>& tA
2278 )
2279 {
2280  tmp<fvMatrix<Type>> tC(tA.ptr());
2281  tC.ref() *= tdsf;
2282  return tC;
2283 }
2284 
2285 template<class Type>
2286 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2287 (
2288  const tmp<volScalarField>& tvsf,
2289  const tmp<fvMatrix<Type>>& tA
2290 )
2291 {
2292  tmp<fvMatrix<Type>> tC(tA.ptr());
2293  tC.ref() *= tvsf;
2294  return tC;
2295 }
2296 
2297 template<class Type>
2298 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2299 (
2300  const dimensioned<scalar>& ds,
2301  const fvMatrix<Type>& A
2302 )
2303 {
2304  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2305  tC.ref() *= ds;
2306  return tC;
2307 }
2308 
2309 template<class Type>
2310 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2311 (
2312  const dimensioned<scalar>& ds,
2313  const tmp<fvMatrix<Type>>& tA
2314 )
2315 {
2316  tmp<fvMatrix<Type>> tC(tA.ptr());
2317  tC.ref() *= ds;
2318  return tC;
2319 }
2320 
2321 
2322 template<class Type>
2324 Foam::operator&
2325 (
2326  const fvMatrix<Type>& M,
2327  const DimensionedField<Type, volMesh>& psi
2328 )
2329 {
2330  tmp<GeometricField<Type, fvPatchField, volMesh>> tMphi
2331  (
2332  new GeometricField<Type, fvPatchField, volMesh>
2333  (
2334  IOobject
2335  (
2336  "M&" + psi.name(),
2337  psi.instance(),
2338  psi.mesh(),
2341  ),
2342  psi.mesh(),
2343  M.dimensions()/dimVol,
2344  extrapolatedCalculatedFvPatchScalarField::typeName
2345  )
2346  );
2347  GeometricField<Type, fvPatchField, volMesh>& Mphi = tMphi.ref();
2348 
2349  // Loop over field components
2350  if (M.hasDiag())
2351  {
2352  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
2353  {
2354  scalarField psiCmpt(psi.field().component(cmpt));
2355  scalarField boundaryDiagCmpt(M.diag());
2356  M.addBoundaryDiag(boundaryDiagCmpt, cmpt);
2357  Mphi.primitiveFieldRef().replace(cmpt, -boundaryDiagCmpt*psiCmpt);
2358  }
2359  }
2360  else
2361  {
2362  Mphi.primitiveFieldRef() = Zero;
2363  }
2364 
2365  Mphi.primitiveFieldRef() += M.lduMatrix::H(psi.field()) + M.source();
2366  M.addBoundarySource(Mphi.primitiveFieldRef());
2367 
2368  Mphi.primitiveFieldRef() /= -psi.mesh().V();
2369  Mphi.correctBoundaryConditions();
2370 
2371  return tMphi;
2372 }
2373 
2374 template<class Type>
2376 Foam::operator&
2377 (
2378  const fvMatrix<Type>& M,
2379  const tmp<DimensionedField<Type, volMesh>>& tpsi
2380 )
2381 {
2382  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = M & tpsi();
2383  tpsi.clear();
2384  return tMpsi;
2385 }
2386 
2387 template<class Type>
2389 Foam::operator&
2390 (
2391  const fvMatrix<Type>& M,
2392  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tpsi
2393 )
2394 {
2395  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = M & tpsi();
2396  tpsi.clear();
2397  return tMpsi;
2398 }
2399 
2400 template<class Type>
2402 Foam::operator&
2403 (
2404  const tmp<fvMatrix<Type>>& tM,
2405  const DimensionedField<Type, volMesh>& psi
2406 )
2407 {
2408  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = tM() & psi;
2409  tM.clear();
2410  return tMpsi;
2411 }
2412 
2413 template<class Type>
2415 Foam::operator&
2416 (
2417  const tmp<fvMatrix<Type>>& tM,
2418  const tmp<DimensionedField<Type, volMesh>>& tpsi
2419 )
2420 {
2421  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = tM() & tpsi();
2422  tM.clear();
2423  tpsi.clear();
2424  return tMpsi;
2425 }
2426 
2427 template<class Type>
2429 Foam::operator&
2430 (
2431  const tmp<fvMatrix<Type>>& tM,
2432  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tpsi
2433 )
2434 {
2435  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = tM() & tpsi();
2436  tM.clear();
2437  tpsi.clear();
2438  return tMpsi;
2439 }
2440 
2441 
2442 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
2443 
2444 template<class Type>
2445 Foam::Ostream& Foam::operator<<(Ostream& os, const fvMatrix<Type>& fvm)
2446 {
2447  os << static_cast<const lduMatrix&>(fvm) << nl
2448  << fvm.dimensions_ << nl
2449  << fvm.source_ << nl
2450  << fvm.internalCoeffs_ << nl
2451  << fvm.boundaryCoeffs_ << endl;
2452 
2453  os.check(FUNCTION_NAME);
2454 
2455  return os;
2456 }
2457 
2458 
2459 // * * * * * * * * * * * * * * * * Solvers * * * * * * * * * * * * * * * * * //
2460 
2461 #include "fvMatrixSolve.C"
2462 
2463 // ************************************************************************* //
void checkMethod(const faMatrix< Type > &, const faMatrix< Type > &, const char *)
Definition: faMatrix.C:958
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
void operator-=(const lduMatrix &)
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Definition: fvMatrix.C:816
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Definition: fvMatrix.C:145
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:93
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
Definition: fvMatrix.C:741
#define InfoInFunction
Report an information message using Foam::Info.
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
const word & name() const
Return name.
Definition: IOobjectI.H:44
A class for handling words, derived from Foam::string.
Definition: word.H:59
void addCmptAvBoundaryDiag(scalarField &diag) const
Definition: fvMatrix.C:129
A field of fields is a PtrList of fields with reference counting.
Definition: FieldField.H:51
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
A class for managing temporary objects.
Definition: PtrList.H:56
void operator+=(const lduMatrix &)
static constexpr const zero Zero
Definition: zero.H:118
Reference counter for various OpenFOAM components.
Definition: refCount.H:47
const Field< Type > & field() const
Return field.
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
Definition: HashOps.H:146
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
tmp< scalarField > D() const
Return the matrix scalar diagonal.
Definition: fvMatrix.C:754
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
Definition: fvMatrix.C:75
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void operator=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1029
lduMatrix is a general matrix class in which the coefficients are stored as three arrays,...
Definition: lduMatrix.H:80
const fileName & instance() const
Definition: IOobjectI.H:153
const dimensionSet & dimensions() const
Definition: fvMatrix.H:289
void relax()
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:724
const word & name() const
Return const reference to name.
void negate()
Definition: fvMatrix.C:1073
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:323
Foam::surfaceFields.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
Definition: fvMatrix.C:111
label eventNo() const
Event number at last update.
Definition: regIOobjectI.H:76
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:508
Dimension set for the base types.
Definition: dimensionSet.H:62
conserve primitiveFieldRef()+
dimensioned< scalar > mag(const dimensioned< Type > &)
void setReferences(const labelUList &cells, const Type &value, const bool forceReference=false)
Set references level for solution.
Definition: fvMatrix.C:524
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:284
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:284
ensightGeoFile & operator<<(ensightGeoFile &, const ensightPart &)
Definition: ensightPart.C:77
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
void deleteDemandDrivenData(DataPtr &dataPtr)
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:328
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
void setValues(const labelUList &cells, const UList< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:486
T & ref() const
Definition: tmpI.H:241
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void operator-=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1122
Generic templated field type.
Definition: Field.H:60
const Mesh & mesh() const
Return mesh.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:57
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
tmp< volScalarField > H1() const
Return H(1)
Definition: fvMatrix.C:878
virtual ~fvMatrix()
Destructor.
Definition: fvMatrix.C:467
SolverPerformance< Type > solve(faMatrix< Type > &, Istream &)
Solve returning the solution statistics given convergence tolerance.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< scalarField > H1() const
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:48
error FatalError
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:118
tmp< Field< Type > > DD() const
Return the matrix Type diagonal.
Definition: fvMatrix.C:763
dynamicFvMesh & mesh
Generic dimensioned Type class.
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:787
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:81
T * ptr() const
Definition: tmpI.H:279
tmp< Field< Type > > H(const Field< Type > &) const
errorManip< error > abort(error &err)
Definition: errorManip.H:135
const dimensionSet & dimensions() const
Return dimensions.
label cellId
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:431
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
void operator *=(const scalarField &)
A non-counting (dummy) refCount.
Definition: refCount.H:56
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:491
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:359
void setValuesFromList(const labelUList &cells, const ListType< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:178
tmp< DimensionedField< typename DimensionedField< Type, GeoMesh >::cmptType, GeoMesh >> cmptAv(const DimensionedField< Type, GeoMesh > &df)
constexpr char nl
Definition: Ostream.H:358
void operator=(const lduMatrix &)
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
void operator+=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1088
Traits class for primitives.
Definition: pTraits.H:50
tmp< Field< Type > > faceH(const Field< Type > &) const
const dictionary & solverDict() const
Return the solver dictionary taking into account finalIteration.
Definition: fvMatrix.C:1013
#define M(I)
label patchi
string upper(const string &original)
Return string transformed with std::toupper on each character.
Definition: stringOps.C:1125
void setOriented(const bool oriented=true)
Set the oriented flag.
uint8_t direction
Definition: direction.H:45
tmp< fvMatrix< Type > > clone() const
Clone.
Definition: fvMatrix.C:455
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvPatchField.H:72
#define FUNCTION_NAME
const dimensionedScalar c
Speed of light in a vacuum.
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
A List with indirect addressing.
Definition: fvMatrix.H:106
const dimensionSet dimVol(dimVolume)
Definition: dimensionSets.H:59
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:71
const cellShapeList & cells
void size(const label n)
Override size to be inconsistent with allocated storage.
Definition: UListI.H:361
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:53
const dimensionSet & dimensions() const
Return const reference to dimensions.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
Generic GeometricField class.
Definition: areaFieldsFwd.H:52
const volScalarField & psi
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
Definition: fvMatrix.C:38
SolverPerformance is the class returned by the LduMatrix solver containing performance statistics.
void clear() const
Definition: tmpI.H:308
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:51
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:926
string lower(const string &original)
Return string transformed with std::tolower on each character.
Definition: stringOps.C:1105
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fvMatrix(const GeometricField< Type, fvPatchField, volMesh > &, const dimensionSet &)
Construct given a field to solve for.
Definition: fvMatrix.C:268
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
UEqn relax()