fvMatrix.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
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 \*---------------------------------------------------------------------------*/
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  tmp<fvMatrix<Type>>::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 template<class Type>
352 :
353  lduMatrix
354  (
355  const_cast<fvMatrix<Type>&>(tfvm()),
356  tfvm.isTmp()
357  ),
358  psi_(tfvm().psi_),
359  dimensions_(tfvm().dimensions_),
360  source_
361  (
362  const_cast<fvMatrix<Type>&>(tfvm()).source_,
363  tfvm.isTmp()
364  ),
365  internalCoeffs_
366  (
367  const_cast<fvMatrix<Type>&>(tfvm()).internalCoeffs_,
368  tfvm.isTmp()
369  ),
370  boundaryCoeffs_
371  (
372  const_cast<fvMatrix<Type>&>(tfvm()).boundaryCoeffs_,
373  tfvm.isTmp()
374  ),
375  faceFluxCorrectionPtr_(nullptr)
376 {
377  if (debug)
378  {
380  << "Copying fvMatrix<Type> for field " << psi_.name() << endl;
381  }
382 
383  if (tfvm().faceFluxCorrectionPtr_)
384  {
385  if (tfvm.isTmp())
386  {
387  faceFluxCorrectionPtr_ = tfvm().faceFluxCorrectionPtr_;
388  tfvm().faceFluxCorrectionPtr_ = nullptr;
389  }
390  else
391  {
392  faceFluxCorrectionPtr_ = new
394  (
395  *(tfvm().faceFluxCorrectionPtr_)
396  );
397  }
398  }
399 
400  tfvm.clear();
401 }
402 
403 
404 template<class Type>
406 (
408  Istream& is
409 )
410 :
411  lduMatrix(psi.mesh()),
412  psi_(psi),
413  dimensions_(is),
414  source_(is),
415  internalCoeffs_(psi.mesh().boundary().size()),
416  boundaryCoeffs_(psi.mesh().boundary().size()),
417  faceFluxCorrectionPtr_(nullptr)
418 {
419  if (debug)
420  {
422  << "Constructing fvMatrix<Type> for field " << psi_.name() << endl;
423  }
424 
425  // Initialise coupling coefficients
426  forAll(psi.mesh().boundary(), patchi)
427  {
428  internalCoeffs_.set
429  (
430  patchi,
431  new Field<Type>
432  (
433  psi.mesh().boundary()[patchi].size(),
434  Zero
435  )
436  );
437 
438  boundaryCoeffs_.set
439  (
440  patchi,
441  new Field<Type>
442  (
443  psi.mesh().boundary()[patchi].size(),
444  Zero
445  )
446  );
447  }
448 
449 }
450 
451 
452 template<class Type>
454 {
455  return tmp<fvMatrix<Type>>
456  (
457  new fvMatrix<Type>(*this)
458  );
459 }
460 
461 
462 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
463 
464 template<class Type>
466 {
467  if (debug)
468  {
470  << "Destroying fvMatrix<Type> for field " << psi_.name() << endl;
471  }
472 
473  if (faceFluxCorrectionPtr_)
474  {
475  delete faceFluxCorrectionPtr_;
476  }
477 }
478 
479 
480 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
481 
482 template<class Type>
484 (
485  const labelUList& cellLabels,
486  const UList<Type>& values
487 )
488 {
489  this->setValuesFromList(cellLabels, values);
490 }
491 
492 
493 template<class Type>
495 (
496  const labelUList& cellLabels,
497  const UIndirectList<Type>& values
498 )
499 {
500  this->setValuesFromList(cellLabels, values);
501 }
502 
503 
504 template<class Type>
506 (
507  const label celli,
508  const Type& value,
509  const bool forceReference
510 )
511 {
512  if ((forceReference || psi_.needReference()) && celli >= 0)
513  {
514  source()[celli] += diag()[celli]*value;
515  diag()[celli] += diag()[celli];
516  }
517 }
518 
519 
520 template<class Type>
522 {
523  if (alpha <= 0)
524  {
525  return;
526  }
527 
528  if (debug)
529  {
531  << "Relaxing " << psi_.name() << " by " << alpha << endl;
532  }
533 
534  Field<Type>& S = source();
535  scalarField& D = diag();
536 
537  // Store the current unrelaxed diagonal for use in updating the source
538  scalarField D0(D);
539 
540  // Calculate the sum-mag off-diagonal from the interior faces
541  scalarField sumOff(D.size(), 0.0);
542  sumMagOffDiag(sumOff);
543 
544  // Handle the boundary contributions to the diagonal
545  forAll(psi_.boundaryField(), patchi)
546  {
547  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
548 
549  if (ptf.size())
550  {
551  const labelUList& pa = lduAddr().patchAddr(patchi);
552  Field<Type>& iCoeffs = internalCoeffs_[patchi];
553 
554  if (ptf.coupled())
555  {
556  const Field<Type>& pCoeffs = boundaryCoeffs_[patchi];
557 
558  // For coupled boundaries add the diagonal and
559  // off-diagonal contributions
560  forAll(pa, face)
561  {
562  D[pa[face]] += component(iCoeffs[face], 0);
563  sumOff[pa[face]] += mag(component(pCoeffs[face], 0));
564  }
565  }
566  else
567  {
568  // For non-coupled boundaries add the maximum magnitude diagonal
569  // contribution to ensure stability
570  forAll(pa, face)
571  {
572  D[pa[face]] += cmptMax(cmptMag(iCoeffs[face]));
573  }
574  }
575  }
576  }
577 
578 
579  if (debug)
580  {
581  // Calculate amount of non-dominance.
582  label nNon = 0;
583  scalar maxNon = 0.0;
584  scalar sumNon = 0.0;
585  forAll(D, celli)
586  {
587  scalar d = (sumOff[celli] - D[celli])/mag(D[celli]);
588 
589  if (d > 0)
590  {
591  nNon++;
592  maxNon = max(maxNon, d);
593  sumNon += d;
594  }
595  }
596 
597  reduce(nNon, sumOp<label>(), UPstream::msgType(), psi_.mesh().comm());
598  reduce
599  (
600  maxNon,
601  maxOp<scalar>(),
603  psi_.mesh().comm()
604  );
605  reduce
606  (
607  sumNon,
608  sumOp<scalar>(),
610  psi_.mesh().comm()
611  );
612  sumNon /= returnReduce
613  (
614  D.size(),
615  sumOp<label>(),
617  psi_.mesh().comm()
618  );
619 
621  << "Matrix dominance test for " << psi_.name() << nl
622  << " number of non-dominant cells : " << nNon << nl
623  << " maximum relative non-dominance : " << maxNon << nl
624  << " average relative non-dominance : " << sumNon << nl
625  << endl;
626  }
627 
628 
629  // Ensure the matrix is diagonally dominant...
630  // Assumes that the central coefficient is positive and ensures it is
631  forAll(D, celli)
632  {
633  D[celli] = max(mag(D[celli]), sumOff[celli]);
634  }
635 
636  // ... then relax
637  D /= alpha;
638 
639  // Now remove the diagonal contribution from coupled boundaries
640  forAll(psi_.boundaryField(), patchi)
641  {
642  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
643 
644  if (ptf.size())
645  {
646  const labelUList& pa = lduAddr().patchAddr(patchi);
647  Field<Type>& iCoeffs = internalCoeffs_[patchi];
648 
649  if (ptf.coupled())
650  {
651  forAll(pa, face)
652  {
653  D[pa[face]] -= component(iCoeffs[face], 0);
654  }
655  }
656  else
657  {
658  forAll(pa, face)
659  {
660  D[pa[face]] -= cmptMin(iCoeffs[face]);
661  }
662  }
663  }
664  }
665 
666  // Finally add the relaxation contribution to the source.
667  S += (D - D0)*psi_.primitiveField();
668 }
669 
670 
671 template<class Type>
673 {
674  if
675  (
676  psi_.mesh().data::template lookupOrDefault<bool>
677  ("finalIteration", false)
678  && psi_.mesh().relaxEquation(psi_.name() + "Final")
679  )
680  {
681  relax(psi_.mesh().equationRelaxationFactor(psi_.name() + "Final"));
682  }
683  else if (psi_.mesh().relaxEquation(psi_.name()))
684  {
685  relax(psi_.mesh().equationRelaxationFactor(psi_.name()));
686  }
687 }
688 
689 
690 template<class Type>
692 (
694  Boundary& bFields
695 )
696 {
697  forAll(bFields, patchi)
698  {
699  bFields[patchi].manipulateMatrix(*this);
700  }
701 }
702 
703 
704 template<class Type>
706 {
707  tmp<scalarField> tdiag(new scalarField(diag()));
708  addCmptAvBoundaryDiag(tdiag.ref());
709  return tdiag;
710 }
711 
712 
713 template<class Type>
715 {
717 
718  forAll(psi_.boundaryField(), patchi)
719  {
720  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
721 
722  if (!ptf.coupled() && ptf.size())
723  {
725  (
726  lduAddr().patchAddr(patchi),
727  internalCoeffs_[patchi],
728  tdiag.ref()
729  );
730  }
731  }
732 
733  return tdiag;
734 }
735 
736 
737 template<class Type>
739 {
740  tmp<volScalarField> tAphi
741  (
743  (
744  "A("+psi_.name()+')',
745  psi_.mesh(),
746  dimensions_/psi_.dimensions()/dimVol,
747  extrapolatedCalculatedFvPatchScalarField::typeName
748  )
749  );
750 
751  tAphi.ref().primitiveFieldRef() = D()/psi_.mesh().V();
752  tAphi.ref().correctBoundaryConditions();
753 
754  return tAphi;
755 }
756 
757 
758 template<class Type>
761 {
763  (
765  (
766  IOobject
767  (
768  "H("+psi_.name()+')',
769  psi_.instance(),
770  psi_.mesh(),
773  ),
774  psi_.mesh(),
775  dimensions_/dimVol,
776  extrapolatedCalculatedFvPatchScalarField::typeName
777  )
778  );
780 
781  // Loop over field components
782  for (direction cmpt=0; cmpt<Type::nComponents; cmpt++)
783  {
784  scalarField psiCmpt(psi_.primitiveField().component(cmpt));
785 
786  scalarField boundaryDiagCmpt(psi_.size(), 0.0);
787  addBoundaryDiag(boundaryDiagCmpt, cmpt);
788  boundaryDiagCmpt.negate();
789  addCmptAvBoundaryDiag(boundaryDiagCmpt);
790 
791  Hphi.primitiveFieldRef().replace(cmpt, boundaryDiagCmpt*psiCmpt);
792  }
793 
794  Hphi.primitiveFieldRef() += lduMatrix::H(psi_.primitiveField()) + source_;
796 
797  Hphi.primitiveFieldRef() /= psi_.mesh().V();
799 
800  typename Type::labelType validComponents
801  (
802  psi_.mesh().template validComponents<Type>()
803  );
804 
805  for (direction cmpt=0; cmpt<Type::nComponents; cmpt++)
806  {
807  if (validComponents[cmpt] == -1)
808  {
809  Hphi.replace
810  (
811  cmpt,
812  dimensionedScalar(Hphi.dimensions(), 0)
813  );
814  }
815  }
816 
817  return tHphi;
818 }
819 
820 
821 template<class Type>
823 {
825  (
827  (
828  "H(1)",
829  psi_.mesh(),
830  dimensions_/(dimVol*psi_.dimensions()),
831  extrapolatedCalculatedFvPatchScalarField::typeName
832  )
833  );
834  volScalarField& H1_ = tH1.ref();
835 
836  H1_.primitiveFieldRef() = lduMatrix::H1();
837 
838  forAll(psi_.boundaryField(), patchi)
839  {
840  const fvPatchField<Type>& ptf = psi_.boundaryField()[patchi];
841 
842  if (ptf.coupled() && ptf.size())
843  {
845  (
846  lduAddr().patchAddr(patchi),
847  boundaryCoeffs_[patchi].component(0),
848  H1_
849  );
850  }
851  }
852 
853  H1_.primitiveFieldRef() /= psi_.mesh().V();
854  H1_.correctBoundaryConditions();
855 
856  return tH1;
857 }
858 
859 
860 template<class Type>
863 flux() const
864 {
865  if (!psi_.mesh().fluxRequired(psi_.name()))
866  {
868  << "flux requested but " << psi_.name()
869  << " not specified in the fluxRequired sub-dictionary"
870  " of fvSchemes."
871  << abort(FatalError);
872  }
873 
874  // construct GeometricField<Type, fvsPatchField, surfaceMesh>
876  (
878  (
879  IOobject
880  (
881  "flux("+psi_.name()+')',
882  psi_.instance(),
883  psi_.mesh(),
886  ),
887  psi_.mesh(),
888  dimensions()
889  )
890  );
892  tfieldFlux.ref();
893 
894  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
895  {
896  fieldFlux.primitiveFieldRef().replace
897  (
898  cmpt,
899  lduMatrix::faceH(psi_.primitiveField().component(cmpt))
900  );
901  }
902 
903  FieldField<Field, Type> InternalContrib = internalCoeffs_;
904 
905  forAll(InternalContrib, patchi)
906  {
907  InternalContrib[patchi] =
909  (
910  InternalContrib[patchi],
911  psi_.boundaryField()[patchi].patchInternalField()
912  );
913  }
914 
915  FieldField<Field, Type> NeighbourContrib = boundaryCoeffs_;
916 
917  forAll(NeighbourContrib, patchi)
918  {
919  if (psi_.boundaryField()[patchi].coupled())
920  {
921  NeighbourContrib[patchi] =
923  (
924  NeighbourContrib[patchi],
925  psi_.boundaryField()[patchi].patchNeighbourField()
926  );
927  }
928  }
929 
931  Boundary& ffbf = fieldFlux.boundaryFieldRef();
932 
933  forAll(ffbf, patchi)
934  {
935  ffbf[patchi] = InternalContrib[patchi] - NeighbourContrib[patchi];
936  }
937 
938  if (faceFluxCorrectionPtr_)
939  {
940  fieldFlux += *faceFluxCorrectionPtr_;
941  }
942 
943  return tfieldFlux;
944 }
945 
946 
947 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
948 
949 template<class Type>
951 {
952  if (this == &fvmv)
953  {
955  << "attempted assignment to self"
956  << abort(FatalError);
957  }
958 
959  if (&psi_ != &(fvmv.psi_))
960  {
962  << "different fields"
963  << abort(FatalError);
964  }
965 
966  dimensions_ = fvmv.dimensions_;
967  lduMatrix::operator=(fvmv);
968  source_ = fvmv.source_;
969  internalCoeffs_ = fvmv.internalCoeffs_;
970  boundaryCoeffs_ = fvmv.boundaryCoeffs_;
971 
972  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
973  {
974  *faceFluxCorrectionPtr_ = *fvmv.faceFluxCorrectionPtr_;
975  }
976  else if (fvmv.faceFluxCorrectionPtr_)
977  {
978  faceFluxCorrectionPtr_ =
980  (*fvmv.faceFluxCorrectionPtr_);
981  }
982 }
983 
984 
985 template<class Type>
987 {
988  operator=(tfvmv());
989  tfvmv.clear();
990 }
991 
992 
993 template<class Type>
995 {
997  source_.negate();
998  internalCoeffs_.negate();
999  boundaryCoeffs_.negate();
1000 
1001  if (faceFluxCorrectionPtr_)
1002  {
1003  faceFluxCorrectionPtr_->negate();
1004  }
1005 }
1006 
1007 
1008 template<class Type>
1010 {
1011  checkMethod(*this, fvmv, "+=");
1012 
1013  dimensions_ += fvmv.dimensions_;
1014  lduMatrix::operator+=(fvmv);
1015  source_ += fvmv.source_;
1016  internalCoeffs_ += fvmv.internalCoeffs_;
1017  boundaryCoeffs_ += fvmv.boundaryCoeffs_;
1018 
1019  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
1020  {
1021  *faceFluxCorrectionPtr_ += *fvmv.faceFluxCorrectionPtr_;
1022  }
1023  else if (fvmv.faceFluxCorrectionPtr_)
1024  {
1025  faceFluxCorrectionPtr_ = new
1027  (
1028  *fvmv.faceFluxCorrectionPtr_
1029  );
1030  }
1031 }
1032 
1033 
1034 template<class Type>
1036 {
1037  operator+=(tfvmv());
1038  tfvmv.clear();
1039 }
1040 
1041 
1042 template<class Type>
1044 {
1045  checkMethod(*this, fvmv, "-=");
1046 
1047  dimensions_ -= fvmv.dimensions_;
1048  lduMatrix::operator-=(fvmv);
1049  source_ -= fvmv.source_;
1050  internalCoeffs_ -= fvmv.internalCoeffs_;
1051  boundaryCoeffs_ -= fvmv.boundaryCoeffs_;
1052 
1053  if (faceFluxCorrectionPtr_ && fvmv.faceFluxCorrectionPtr_)
1054  {
1055  *faceFluxCorrectionPtr_ -= *fvmv.faceFluxCorrectionPtr_;
1056  }
1057  else if (fvmv.faceFluxCorrectionPtr_)
1058  {
1059  faceFluxCorrectionPtr_ =
1061  (-*fvmv.faceFluxCorrectionPtr_);
1062  }
1063 }
1064 
1065 
1066 template<class Type>
1068 {
1069  operator-=(tfvmv());
1070  tfvmv.clear();
1071 }
1072 
1073 
1074 template<class Type>
1075 void Foam::fvMatrix<Type>::operator+=
1078 )
1079 {
1080  checkMethod(*this, su, "+=");
1081  source() -= su.mesh().V()*su.field();
1082 }
1083 
1084 
1085 template<class Type>
1086 void Foam::fvMatrix<Type>::operator+=
1089 )
1090 {
1091  operator+=(tsu());
1092  tsu.clear();
1093 }
1094 
1095 
1096 template<class Type>
1097 void Foam::fvMatrix<Type>::operator+=
1100 )
1101 {
1102  operator+=(tsu());
1103  tsu.clear();
1104 }
1105 
1106 
1107 template<class Type>
1108 void Foam::fvMatrix<Type>::operator-=
1111 )
1112 {
1113  checkMethod(*this, su, "-=");
1114  source() += su.mesh().V()*su.field();
1115 }
1116 
1117 
1118 template<class Type>
1119 void Foam::fvMatrix<Type>::operator-=
1122 )
1123 {
1124  operator-=(tsu());
1125  tsu.clear();
1126 }
1127 
1128 
1129 template<class Type>
1130 void Foam::fvMatrix<Type>::operator-=
1133 )
1134 {
1135  operator-=(tsu());
1136  tsu.clear();
1137 }
1138 
1139 
1140 template<class Type>
1141 void Foam::fvMatrix<Type>::operator+=
1143  const dimensioned<Type>& su
1144 )
1145 {
1146  source() -= psi().mesh().V()*su;
1147 }
1148 
1149 
1150 template<class Type>
1151 void Foam::fvMatrix<Type>::operator-=
1153  const dimensioned<Type>& su
1154 )
1155 {
1156  source() += psi().mesh().V()*su;
1157 }
1158 
1159 
1160 template<class Type>
1161 void Foam::fvMatrix<Type>::operator+=
1163  const zero&
1164 )
1165 {}
1166 
1167 
1168 template<class Type>
1169 void Foam::fvMatrix<Type>::operator-=
1171  const zero&
1172 )
1173 {}
1174 
1175 
1176 template<class Type>
1177 void Foam::fvMatrix<Type>::operator*=
1179  const volScalarField::Internal& dsf
1180 )
1181 {
1182  dimensions_ *= dsf.dimensions();
1183  lduMatrix::operator*=(dsf.field());
1184  source_ *= dsf.field();
1185 
1186  forAll(boundaryCoeffs_, patchi)
1187  {
1188  scalarField pisf
1189  (
1190  dsf.mesh().boundary()[patchi].patchInternalField(dsf.field())
1191  );
1192 
1193  internalCoeffs_[patchi] *= pisf;
1194  boundaryCoeffs_[patchi] *= pisf;
1195  }
1196 
1197  if (faceFluxCorrectionPtr_)
1198  {
1200  << "cannot scale a matrix containing a faceFluxCorrection"
1201  << abort(FatalError);
1202  }
1203 }
1204 
1205 
1206 template<class Type>
1207 void Foam::fvMatrix<Type>::operator*=
1209  const tmp<volScalarField::Internal>& tdsf
1210 )
1211 {
1212  operator*=(tdsf());
1213  tdsf.clear();
1214 }
1215 
1216 
1217 template<class Type>
1218 void Foam::fvMatrix<Type>::operator*=
1220  const tmp<volScalarField>& tvsf
1221 )
1222 {
1223  operator*=(tvsf());
1224  tvsf.clear();
1225 }
1226 
1227 
1228 template<class Type>
1229 void Foam::fvMatrix<Type>::operator*=
1231  const dimensioned<scalar>& ds
1232 )
1233 {
1234  dimensions_ *= ds.dimensions();
1235  lduMatrix::operator*=(ds.value());
1236  source_ *= ds.value();
1237  internalCoeffs_ *= ds.value();
1238  boundaryCoeffs_ *= ds.value();
1239 
1240  if (faceFluxCorrectionPtr_)
1241  {
1242  *faceFluxCorrectionPtr_ *= ds.value();
1243  }
1244 }
1245 
1246 
1247 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
1248 
1249 template<class Type>
1250 void Foam::checkMethod
1252  const fvMatrix<Type>& fvm1,
1253  const fvMatrix<Type>& fvm2,
1254  const char* op
1255 )
1256 {
1257  if (&fvm1.psi() != &fvm2.psi())
1258  {
1260  << "incompatible fields for operation "
1261  << endl << " "
1262  << "[" << fvm1.psi().name() << "] "
1263  << op
1264  << " [" << fvm2.psi().name() << "]"
1265  << abort(FatalError);
1266  }
1267 
1268  if (dimensionSet::debug && fvm1.dimensions() != fvm2.dimensions())
1269  {
1271  << "incompatible dimensions for operation "
1272  << endl << " "
1273  << "[" << fvm1.psi().name() << fvm1.dimensions()/dimVolume << " ] "
1274  << op
1275  << " [" << fvm2.psi().name() << fvm2.dimensions()/dimVolume << " ]"
1276  << abort(FatalError);
1277  }
1278 }
1279 
1280 
1281 template<class Type>
1282 void Foam::checkMethod
1284  const fvMatrix<Type>& fvm,
1286  const char* op
1287 )
1288 {
1289  if (dimensionSet::debug && fvm.dimensions()/dimVolume != df.dimensions())
1290  {
1292  << endl << " "
1293  << "[" << fvm.psi().name() << fvm.dimensions()/dimVolume << " ] "
1294  << op
1295  << " [" << df.name() << df.dimensions() << " ]"
1296  << abort(FatalError);
1297  }
1298 }
1299 
1300 
1301 template<class Type>
1302 void Foam::checkMethod
1304  const fvMatrix<Type>& fvm,
1305  const dimensioned<Type>& dt,
1306  const char* op
1307 )
1308 {
1309  if (dimensionSet::debug && fvm.dimensions()/dimVolume != dt.dimensions())
1310  {
1312  << "incompatible dimensions for operation "
1313  << endl << " "
1314  << "[" << fvm.psi().name() << fvm.dimensions()/dimVolume << " ] "
1315  << op
1316  << " [" << dt.name() << dt.dimensions() << " ]"
1317  << abort(FatalError);
1318  }
1319 }
1320 
1321 
1322 template<class Type>
1324 (
1325  const fvMatrix<Type>& A
1326 )
1327 {
1328  tmp<Foam::fvMatrix<Type>> tAcorr = A - (A & A.psi());
1329 
1330  // Delete the faceFluxCorrection from the correction matrix
1331  // as it does not have a clear meaning or purpose
1332  deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
1333 
1334  return tAcorr;
1335 }
1336 
1337 
1338 template<class Type>
1340 (
1341  const tmp<fvMatrix<Type>>& tA
1342 )
1343 {
1344  tmp<Foam::fvMatrix<Type>> tAcorr = tA - (tA() & tA().psi());
1345 
1346  // Delete the faceFluxCorrection from the correction matrix
1347  // as it does not have a clear meaning or purpose
1348  deleteDemandDrivenData(tAcorr.ref().faceFluxCorrectionPtr());
1349 
1350  return tAcorr;
1351 }
1352 
1353 
1354 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
1355 
1356 template<class Type>
1357 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1358 (
1359  const fvMatrix<Type>& A,
1360  const fvMatrix<Type>& B
1361 )
1362 {
1363  checkMethod(A, B, "==");
1364  return (A - B);
1365 }
1366 
1367 template<class Type>
1368 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1369 (
1370  const tmp<fvMatrix<Type>>& tA,
1371  const fvMatrix<Type>& B
1372 )
1373 {
1374  checkMethod(tA(), B, "==");
1375  return (tA - B);
1376 }
1377 
1378 template<class Type>
1379 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1380 (
1381  const fvMatrix<Type>& A,
1382  const tmp<fvMatrix<Type>>& tB
1383 )
1384 {
1385  checkMethod(A, tB(), "==");
1386  return (A - tB);
1387 }
1388 
1389 template<class Type>
1390 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1391 (
1392  const tmp<fvMatrix<Type>>& tA,
1393  const tmp<fvMatrix<Type>>& tB
1394 )
1395 {
1396  checkMethod(tA(), tB(), "==");
1397  return (tA - tB);
1398 }
1399 
1400 template<class Type>
1401 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1402 (
1403  const fvMatrix<Type>& A,
1405 )
1406 {
1407  checkMethod(A, su, "==");
1408  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1409  tC.ref().source() += su.mesh().V()*su.field();
1410  return tC;
1411 }
1412 
1413 template<class Type>
1414 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1415 (
1416  const fvMatrix<Type>& A,
1418 )
1419 {
1420  checkMethod(A, tsu(), "==");
1421  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1422  tC.ref().source() += tsu().mesh().V()*tsu().field();
1423  tsu.clear();
1424  return tC;
1425 }
1426 
1427 template<class Type>
1428 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1429 (
1430  const fvMatrix<Type>& A,
1432 )
1433 {
1434  checkMethod(A, tsu(), "==");
1435  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1436  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1437  tsu.clear();
1438  return tC;
1439 }
1440 
1441 template<class Type>
1442 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1443 (
1444  const tmp<fvMatrix<Type>>& tA,
1446 )
1447 {
1448  checkMethod(tA(), su, "==");
1449  tmp<fvMatrix<Type>> tC(tA.ptr());
1450  tC.ref().source() += su.mesh().V()*su.field();
1451  return tC;
1452 }
1453 
1454 template<class Type>
1455 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1456 (
1457  const tmp<fvMatrix<Type>>& tA,
1459 )
1460 {
1461  checkMethod(tA(), tsu(), "==");
1462  tmp<fvMatrix<Type>> tC(tA.ptr());
1463  tC.ref().source() += tsu().mesh().V()*tsu().field();
1464  tsu.clear();
1465  return tC;
1466 }
1467 
1468 template<class Type>
1469 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1470 (
1471  const tmp<fvMatrix<Type>>& tA,
1473 )
1474 {
1475  checkMethod(tA(), tsu(), "==");
1476  tmp<fvMatrix<Type>> tC(tA.ptr());
1477  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1478  tsu.clear();
1479  return tC;
1480 }
1481 
1482 template<class Type>
1483 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1484 (
1485  const fvMatrix<Type>& A,
1486  const dimensioned<Type>& su
1487 )
1488 {
1489  checkMethod(A, su, "==");
1490  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1491  tC.ref().source() += A.psi().mesh().V()*su.value();
1492  return tC;
1493 }
1494 
1495 template<class Type>
1496 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1497 (
1498  const tmp<fvMatrix<Type>>& tA,
1499  const dimensioned<Type>& su
1500 )
1501 {
1502  checkMethod(tA(), su, "==");
1503  tmp<fvMatrix<Type>> tC(tA.ptr());
1504  tC.ref().source() += tC().psi().mesh().V()*su.value();
1505  return tC;
1506 }
1507 
1508 template<class Type>
1509 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1510 (
1511  const fvMatrix<Type>& A,
1512  const zero&
1513 )
1514 {
1515  return A;
1516 }
1517 
1518 
1519 template<class Type>
1520 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator==
1521 (
1522  const tmp<fvMatrix<Type>>& tA,
1523  const zero&
1524 )
1525 {
1526  return tA;
1527 }
1528 
1529 
1530 template<class Type>
1531 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1532 (
1533  const fvMatrix<Type>& A
1534 )
1535 {
1536  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1537  tC.ref().negate();
1538  return tC;
1539 }
1540 
1541 template<class Type>
1542 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1543 (
1544  const tmp<fvMatrix<Type>>& tA
1545 )
1546 {
1547  tmp<fvMatrix<Type>> tC(tA.ptr());
1548  tC.ref().negate();
1549  return tC;
1550 }
1551 
1552 
1553 template<class Type>
1554 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1555 (
1556  const fvMatrix<Type>& A,
1557  const fvMatrix<Type>& B
1558 )
1559 {
1560  checkMethod(A, B, "+");
1561  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1562  tC.ref() += B;
1563  return tC;
1564 }
1565 
1566 template<class Type>
1567 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1568 (
1569  const tmp<fvMatrix<Type>>& tA,
1570  const fvMatrix<Type>& B
1571 )
1572 {
1573  checkMethod(tA(), B, "+");
1574  tmp<fvMatrix<Type>> tC(tA.ptr());
1575  tC.ref() += B;
1576  return tC;
1577 }
1578 
1579 template<class Type>
1580 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1581 (
1582  const fvMatrix<Type>& A,
1583  const tmp<fvMatrix<Type>>& tB
1584 )
1585 {
1586  checkMethod(A, tB(), "+");
1587  tmp<fvMatrix<Type>> tC(tB.ptr());
1588  tC.ref() += A;
1589  return tC;
1590 }
1591 
1592 template<class Type>
1593 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1594 (
1595  const tmp<fvMatrix<Type>>& tA,
1596  const tmp<fvMatrix<Type>>& tB
1597 )
1598 {
1599  checkMethod(tA(), tB(), "+");
1600  tmp<fvMatrix<Type>> tC(tA.ptr());
1601  tC.ref() += tB();
1602  tB.clear();
1603  return tC;
1604 }
1605 
1606 template<class Type>
1607 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1608 (
1609  const fvMatrix<Type>& A,
1611 )
1612 {
1613  checkMethod(A, su, "+");
1614  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1615  tC.ref().source() -= su.mesh().V()*su.field();
1616  return tC;
1617 }
1618 
1619 template<class Type>
1620 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1621 (
1622  const fvMatrix<Type>& A,
1624 )
1625 {
1626  checkMethod(A, tsu(), "+");
1627  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1628  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1629  tsu.clear();
1630  return tC;
1631 }
1632 
1633 template<class Type>
1634 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1635 (
1636  const fvMatrix<Type>& A,
1638 )
1639 {
1640  checkMethod(A, tsu(), "+");
1641  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1642  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1643  tsu.clear();
1644  return tC;
1645 }
1646 
1647 template<class Type>
1648 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1649 (
1650  const tmp<fvMatrix<Type>>& tA,
1652 )
1653 {
1654  checkMethod(tA(), su, "+");
1655  tmp<fvMatrix<Type>> tC(tA.ptr());
1656  tC.ref().source() -= su.mesh().V()*su.field();
1657  return tC;
1658 }
1659 
1660 template<class Type>
1661 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1662 (
1663  const tmp<fvMatrix<Type>>& tA,
1665 )
1666 {
1667  checkMethod(tA(), tsu(), "+");
1668  tmp<fvMatrix<Type>> tC(tA.ptr());
1669  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1670  tsu.clear();
1671  return tC;
1672 }
1673 
1674 template<class Type>
1675 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1676 (
1677  const tmp<fvMatrix<Type>>& tA,
1679 )
1680 {
1681  checkMethod(tA(), tsu(), "+");
1682  tmp<fvMatrix<Type>> tC(tA.ptr());
1683  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1684  tsu.clear();
1685  return tC;
1686 }
1687 
1688 template<class Type>
1689 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1690 (
1692  const fvMatrix<Type>& A
1693 )
1694 {
1695  checkMethod(A, su, "+");
1696  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1697  tC.ref().source() -= su.mesh().V()*su.field();
1698  return tC;
1699 }
1700 
1701 template<class Type>
1702 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1703 (
1705  const fvMatrix<Type>& A
1706 )
1707 {
1708  checkMethod(A, tsu(), "+");
1709  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1710  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1711  tsu.clear();
1712  return tC;
1713 }
1714 
1715 template<class Type>
1716 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1717 (
1719  const fvMatrix<Type>& A
1720 )
1721 {
1722  checkMethod(A, tsu(), "+");
1723  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1724  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1725  tsu.clear();
1726  return tC;
1727 }
1728 
1729 template<class Type>
1730 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1731 (
1733  const tmp<fvMatrix<Type>>& tA
1734 )
1735 {
1736  checkMethod(tA(), su, "+");
1737  tmp<fvMatrix<Type>> tC(tA.ptr());
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 (
1746  const tmp<fvMatrix<Type>>& tA
1747 )
1748 {
1749  checkMethod(tA(), tsu(), "+");
1750  tmp<fvMatrix<Type>> tC(tA.ptr());
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 (
1760  const tmp<fvMatrix<Type>>& tA
1761 )
1762 {
1763  checkMethod(tA(), tsu(), "+");
1764  tmp<fvMatrix<Type>> tC(tA.ptr());
1765  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1766  tsu.clear();
1767  return tC;
1768 }
1769 
1770 
1771 template<class Type>
1772 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1773 (
1774  const fvMatrix<Type>& A,
1775  const fvMatrix<Type>& B
1776 )
1777 {
1778  checkMethod(A, B, "-");
1779  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1780  tC.ref() -= B;
1781  return tC;
1782 }
1783 
1784 template<class Type>
1785 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1786 (
1787  const tmp<fvMatrix<Type>>& tA,
1788  const fvMatrix<Type>& B
1789 )
1790 {
1791  checkMethod(tA(), B, "-");
1792  tmp<fvMatrix<Type>> tC(tA.ptr());
1793  tC.ref() -= B;
1794  return tC;
1795 }
1796 
1797 template<class Type>
1798 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1799 (
1800  const fvMatrix<Type>& A,
1801  const tmp<fvMatrix<Type>>& tB
1802 )
1803 {
1804  checkMethod(A, tB(), "-");
1805  tmp<fvMatrix<Type>> tC(tB.ptr());
1806  tC.ref() -= A;
1807  tC.ref().negate();
1808  return tC;
1809 }
1810 
1811 template<class Type>
1812 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1813 (
1814  const tmp<fvMatrix<Type>>& tA,
1815  const tmp<fvMatrix<Type>>& tB
1816 )
1817 {
1818  checkMethod(tA(), tB(), "-");
1819  tmp<fvMatrix<Type>> tC(tA.ptr());
1820  tC.ref() -= tB();
1821  tB.clear();
1822  return tC;
1823 }
1824 
1825 template<class Type>
1826 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1827 (
1828  const fvMatrix<Type>& A,
1830 )
1831 {
1832  checkMethod(A, su, "-");
1833  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1834  tC.ref().source() += su.mesh().V()*su.field();
1835  return tC;
1836 }
1837 
1838 template<class Type>
1839 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1840 (
1841  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().field();
1848  tsu.clear();
1849  return tC;
1850 }
1851 
1852 template<class Type>
1853 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1854 (
1855  const fvMatrix<Type>& A,
1857 )
1858 {
1859  checkMethod(A, tsu(), "-");
1860  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1861  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1862  tsu.clear();
1863  return tC;
1864 }
1865 
1866 template<class Type>
1867 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1868 (
1869  const tmp<fvMatrix<Type>>& tA,
1871 )
1872 {
1873  checkMethod(tA(), su, "-");
1874  tmp<fvMatrix<Type>> tC(tA.ptr());
1875  tC.ref().source() += su.mesh().V()*su.field();
1876  return tC;
1877 }
1878 
1879 template<class Type>
1880 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1881 (
1882  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().field();
1889  tsu.clear();
1890  return tC;
1891 }
1892 
1893 template<class Type>
1894 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1895 (
1896  const tmp<fvMatrix<Type>>& tA,
1898 )
1899 {
1900  checkMethod(tA(), tsu(), "-");
1901  tmp<fvMatrix<Type>> tC(tA.ptr());
1902  tC.ref().source() += tsu().mesh().V()*tsu().primitiveField();
1903  tsu.clear();
1904  return tC;
1905 }
1906 
1907 template<class Type>
1908 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1909 (
1911  const fvMatrix<Type>& A
1912 )
1913 {
1914  checkMethod(A, su, "-");
1915  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1916  tC.ref().negate();
1917  tC.ref().source() -= su.mesh().V()*su.field();
1918  return tC;
1919 }
1920 
1921 template<class Type>
1922 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1923 (
1925  const fvMatrix<Type>& A
1926 )
1927 {
1928  checkMethod(A, tsu(), "-");
1929  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1930  tC.ref().negate();
1931  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1932  tsu.clear();
1933  return tC;
1934 }
1935 
1936 template<class Type>
1937 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1938 (
1940  const fvMatrix<Type>& A
1941 )
1942 {
1943  checkMethod(A, tsu(), "-");
1944  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
1945  tC.ref().negate();
1946  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1947  tsu.clear();
1948  return tC;
1949 }
1950 
1951 template<class Type>
1952 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1953 (
1955  const tmp<fvMatrix<Type>>& tA
1956 )
1957 {
1958  checkMethod(tA(), su, "-");
1959  tmp<fvMatrix<Type>> tC(tA.ptr());
1960  tC.ref().negate();
1961  tC.ref().source() -= su.mesh().V()*su.field();
1962  return tC;
1963 }
1964 
1965 template<class Type>
1966 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1967 (
1969  const tmp<fvMatrix<Type>>& tA
1970 )
1971 {
1972  checkMethod(tA(), tsu(), "-");
1973  tmp<fvMatrix<Type>> tC(tA.ptr());
1974  tC.ref().negate();
1975  tC.ref().source() -= tsu().mesh().V()*tsu().field();
1976  tsu.clear();
1977  return tC;
1978 }
1979 
1980 template<class Type>
1981 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
1982 (
1984  const tmp<fvMatrix<Type>>& tA
1985 )
1986 {
1987  checkMethod(tA(), tsu(), "-");
1988  tmp<fvMatrix<Type>> tC(tA.ptr());
1989  tC.ref().negate();
1990  tC.ref().source() -= tsu().mesh().V()*tsu().primitiveField();
1991  tsu.clear();
1992  return tC;
1993 }
1994 
1995 template<class Type>
1996 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
1997 (
1998  const fvMatrix<Type>& A,
1999  const dimensioned<Type>& su
2000 )
2001 {
2002  checkMethod(A, su, "+");
2003  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2004  tC.ref().source() -= su.value()*A.psi().mesh().V();
2005  return tC;
2006 }
2007 
2008 template<class Type>
2009 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2010 (
2011  const tmp<fvMatrix<Type>>& tA,
2012  const dimensioned<Type>& su
2013 )
2014 {
2015  checkMethod(tA(), su, "+");
2016  tmp<fvMatrix<Type>> tC(tA.ptr());
2017  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2018  return tC;
2019 }
2020 
2021 template<class Type>
2022 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2023 (
2024  const dimensioned<Type>& su,
2025  const fvMatrix<Type>& A
2026 )
2027 {
2028  checkMethod(A, su, "+");
2029  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2030  tC.ref().source() -= su.value()*A.psi().mesh().V();
2031  return tC;
2032 }
2033 
2034 template<class Type>
2035 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator+
2036 (
2037  const dimensioned<Type>& su,
2038  const tmp<fvMatrix<Type>>& tA
2039 )
2040 {
2041  checkMethod(tA(), su, "+");
2042  tmp<fvMatrix<Type>> tC(tA.ptr());
2043  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2044  return tC;
2045 }
2046 
2047 template<class Type>
2048 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2049 (
2050  const fvMatrix<Type>& A,
2051  const dimensioned<Type>& su
2052 )
2053 {
2054  checkMethod(A, su, "-");
2055  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2056  tC.ref().source() += su.value()*tC().psi().mesh().V();
2057  return tC;
2058 }
2059 
2060 template<class Type>
2061 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2062 (
2063  const tmp<fvMatrix<Type>>& tA,
2064  const dimensioned<Type>& su
2065 )
2066 {
2067  checkMethod(tA(), su, "-");
2068  tmp<fvMatrix<Type>> tC(tA.ptr());
2069  tC.ref().source() += su.value()*tC().psi().mesh().V();
2070  return tC;
2071 }
2072 
2073 template<class Type>
2074 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2075 (
2076  const dimensioned<Type>& su,
2077  const fvMatrix<Type>& A
2078 )
2079 {
2080  checkMethod(A, su, "-");
2081  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2082  tC.ref().negate();
2083  tC.ref().source() -= su.value()*A.psi().mesh().V();
2084  return tC;
2085 }
2086 
2087 template<class Type>
2088 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator-
2089 (
2090  const dimensioned<Type>& su,
2091  const tmp<fvMatrix<Type>>& tA
2092 )
2093 {
2094  checkMethod(tA(), su, "-");
2095  tmp<fvMatrix<Type>> tC(tA.ptr());
2096  tC.ref().negate();
2097  tC.ref().source() -= su.value()*tC().psi().mesh().V();
2098  return tC;
2099 }
2100 
2101 
2102 template<class Type>
2103 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2104 (
2105  const volScalarField::Internal& dsf,
2106  const fvMatrix<Type>& A
2107 )
2108 {
2109  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2110  tC.ref() *= dsf;
2111  return tC;
2112 }
2113 
2114 template<class Type>
2115 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2116 (
2117  const tmp<volScalarField::Internal>& tdsf,
2118  const fvMatrix<Type>& A
2119 )
2120 {
2121  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2122  tC.ref() *= tdsf;
2123  return tC;
2124 }
2125 
2126 template<class Type>
2127 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2128 (
2129  const tmp<volScalarField>& tvsf,
2130  const fvMatrix<Type>& A
2131 )
2132 {
2133  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2134  tC.ref() *= tvsf;
2135  return tC;
2136 }
2137 
2138 template<class Type>
2139 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2140 (
2141  const volScalarField::Internal& dsf,
2142  const tmp<fvMatrix<Type>>& tA
2143 )
2144 {
2145  tmp<fvMatrix<Type>> tC(tA.ptr());
2146  tC.ref() *= dsf;
2147  return tC;
2148 }
2149 
2150 template<class Type>
2151 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2152 (
2153  const tmp<volScalarField::Internal>& tdsf,
2154  const tmp<fvMatrix<Type>>& tA
2155 )
2156 {
2157  tmp<fvMatrix<Type>> tC(tA.ptr());
2158  tC.ref() *= tdsf;
2159  return tC;
2160 }
2161 
2162 template<class Type>
2163 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2164 (
2165  const tmp<volScalarField>& tvsf,
2166  const tmp<fvMatrix<Type>>& tA
2167 )
2168 {
2169  tmp<fvMatrix<Type>> tC(tA.ptr());
2170  tC.ref() *= tvsf;
2171  return tC;
2172 }
2173 
2174 template<class Type>
2175 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2176 (
2177  const dimensioned<scalar>& ds,
2178  const fvMatrix<Type>& A
2179 )
2180 {
2181  tmp<fvMatrix<Type>> tC(new fvMatrix<Type>(A));
2182  tC.ref() *= ds;
2183  return tC;
2184 }
2185 
2186 template<class Type>
2187 Foam::tmp<Foam::fvMatrix<Type>> Foam::operator*
2188 (
2189  const dimensioned<scalar>& ds,
2190  const tmp<fvMatrix<Type>>& tA
2191 )
2192 {
2193  tmp<fvMatrix<Type>> tC(tA.ptr());
2194  tC.ref() *= ds;
2195  return tC;
2196 }
2197 
2198 
2199 template<class Type>
2201 Foam::operator&
2202 (
2203  const fvMatrix<Type>& M,
2205 )
2206 {
2208  (
2210  (
2211  IOobject
2212  (
2213  "M&" + psi.name(),
2214  psi.instance(),
2215  psi.mesh(),
2218  ),
2219  psi.mesh(),
2220  M.dimensions()/dimVol,
2221  extrapolatedCalculatedFvPatchScalarField::typeName
2222  )
2223  );
2225 
2226  // Loop over field components
2227  if (M.hasDiag())
2228  {
2229  for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
2230  {
2231  scalarField psiCmpt(psi.field().component(cmpt));
2232  scalarField boundaryDiagCmpt(M.diag());
2233  M.addBoundaryDiag(boundaryDiagCmpt, cmpt);
2234  Mphi.primitiveFieldRef().replace(cmpt, -boundaryDiagCmpt*psiCmpt);
2235  }
2236  }
2237  else
2238  {
2239  Mphi.primitiveFieldRef() = Zero;
2240  }
2241 
2242  Mphi.primitiveFieldRef() += M.lduMatrix::H(psi.field()) + M.source();
2243  M.addBoundarySource(Mphi.primitiveFieldRef());
2244 
2245  Mphi.primitiveFieldRef() /= -psi.mesh().V();
2247 
2248  return tMphi;
2249 }
2250 
2251 template<class Type>
2253 Foam::operator&
2254 (
2255  const fvMatrix<Type>& M,
2257 )
2258 {
2260  tpsi.clear();
2261  return tMpsi;
2262 }
2263 
2264 template<class Type>
2266 Foam::operator&
2267 (
2268  const fvMatrix<Type>& M,
2270 )
2271 {
2273  tpsi.clear();
2274  return tMpsi;
2275 }
2276 
2277 template<class Type>
2279 Foam::operator&
2280 (
2281  const tmp<fvMatrix<Type>>& tM,
2283 )
2284 {
2286  tM.clear();
2287  return tMpsi;
2288 }
2289 
2290 template<class Type>
2292 Foam::operator&
2293 (
2294  const tmp<fvMatrix<Type>>& tM,
2296 )
2297 {
2298  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = tM() & tpsi();
2299  tM.clear();
2300  tpsi.clear();
2301  return tMpsi;
2302 }
2303 
2304 template<class Type>
2306 Foam::operator&
2307 (
2308  const tmp<fvMatrix<Type>>& tM,
2310 )
2311 {
2312  tmp<GeometricField<Type, fvPatchField, volMesh>> tMpsi = tM() & tpsi();
2313  tM.clear();
2314  tpsi.clear();
2315  return tMpsi;
2316 }
2317 
2318 
2319 // * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
2320 
2321 template<class Type>
2322 Foam::Ostream& Foam::operator<<(Ostream& os, const fvMatrix<Type>& fvm)
2323 {
2324  os << static_cast<const lduMatrix&>(fvm) << nl
2325  << fvm.dimensions_ << nl
2326  << fvm.source_ << nl
2327  << fvm.internalCoeffs_ << nl
2328  << fvm.boundaryCoeffs_ << endl;
2329 
2330  os.check("Ostream& operator<<(Ostream&, fvMatrix<Type>&");
2331 
2332  return os;
2333 }
2334 
2335 
2336 // * * * * * * * * * * * * * * * * Solvers * * * * * * * * * * * * * * * * * //
2337 
2338 #include "fvMatrixSolve.C"
2339 
2340 // ************************************************************************* //
Foam::surfaceFields.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
void addToInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Add patch contribution to internal field.
Definition: fvMatrix.C:38
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
tmp< Field< Type > > faceH(const Field< Type > &) const
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
const word & name() const
Return name.
Definition: IOobject.H:295
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:230
void addBoundaryDiag(scalarField &diag, const direction cmpt) const
Definition: fvMatrix.C:111
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Reference counter for various OpenFOAM components.
Definition: refCount.H:49
void relax()
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:672
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:92
tmp< DimensionedField< typename DimensionedField< Type, GeoMesh >::cmptType, GeoMesh >> cmptAv(const DimensionedField< Type, GeoMesh > &df)
void negate()
Definition: fvMatrix.C:994
lduMatrix(const lduMesh &)
Construct given an LDU addressed mesh.
Definition: lduMatrix.C:43
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:282
uint8_t direction
Definition: direction.H:45
void addCmptAvBoundaryDiag(scalarField &diag) const
Definition: fvMatrix.C:129
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
void setReference(const label celli, const Type &value, const bool forceReference=false)
Set reference level for solution.
Definition: fvMatrix.C:506
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< GeometricField< Type, fvPatchField, volMesh > > H() const
Return the H operation source.
Definition: fvMatrix.C:760
Traits class for primitives.
Definition: pTraits.H:50
void operator-=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1043
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const dimensionSet dimVol(dimVolume)
Definition: dimensionSets.H:59
const cellList & cells() const
Generic GeometricField class.
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
Generic dimensioned Type class.
tmp< volScalarField > H1() const
Return H(1)
Definition: fvMatrix.C:822
void subtractFromInternalField(const labelUList &addr, const Field< Type2 > &pf, Field< Type2 > &intf) const
Subtract patch contribution from internal field.
Definition: fvMatrix.C:75
DimensionedField< scalar, volMesh > Internal
Type of the internal field from which this GeometricField is derived.
tmp< scalarField > H1() const
void sumMagOffDiag(scalarField &sumOff) const
conserve primitiveFieldRef()+
Generic field type.
Definition: FieldField.H:51
tmp< scalarField > D() const
Return the matrix scalar diagonal.
Definition: fvMatrix.C:705
void operator=(const fvMatrix< Type > &)
Definition: fvMatrix.C:950
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:284
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Mesh &, const dimensionSet &, const word &patchFieldType=fvPatchField< scalar >::calculatedType())
Return a temporary field constructed from name, mesh, dimensionSet.
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
Dimension set for the base types.
Definition: dimensionSet.H:120
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
dynamicFvMesh & mesh
const cellShapeList & cells
Pre-declare SubField and related Field type.
Definition: Field.H:56
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:412
void setValues(const labelUList &cells, const UList< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:484
virtual const labelUList & patchAddr(const label patchNo) const =0
Return patch to internal addressing given patch number.
virtual ~fvMatrix()
Destructor.
Definition: fvMatrix.C:465
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:313
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
label eventNo() const
Event number at last update.
Definition: regIOobjectI.H:83
static const zero Zero
Definition: zero.H:97
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void operator=(const lduMatrix &)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
tmp< fvMatrix< Type > > clone() const
Clone.
Definition: fvMatrix.C:453
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:265
void operator*=(const scalarField &)
Field< Type > & source()
Definition: fvMatrix.H:292
const Mesh & mesh() const
Return mesh.
void operator*=(const volScalarField::Internal &)
Definition: fvMatrix.C:1178
const word & name() const
Return const reference to name.
const Field< Type > & field() const
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:278
Internal & ref()
Return a reference to the dimensioned internal field.
void operator+=(const fvMatrix< Type > &)
Definition: fvMatrix.C:1009
tmp< Field< Type > > DD() const
Return the matrix Type diagonal.
Definition: fvMatrix.C:714
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: lduMatrix.H:79
const fileName & instance() const
Definition: IOobject.H:378
tmp< volScalarField > A() const
Return the central coefficient.
Definition: fvMatrix.C:738
label patchi
const lduAddressing & lduAddr() const
Return the LDU addressing.
Definition: lduMatrix.H:550
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A List with indirect addressing.
Definition: fvMatrix.H:106
const dimensionedScalar c
Speed of light in a vacuum.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return const reference to dimensions.
void setValuesFromList(const labelUList &cells, const ListType< Type > &values)
Set solution in given cells to the specified values.
Definition: fvMatrix.C:178
void correctBoundaryConditions()
Correct boundary field.
void checkMethod(const fvMatrix< Type > &, const fvMatrix< Type > &, const char *)
Definition: fvMatrix.C:1251
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
dimensioned< scalar > mag(const dimensioned< Type > &)
void addBoundarySource(Field< Type > &source, const bool couples=true) const
Definition: fvMatrix.C:145
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
Definition: fvMatrix.C:863
void boundaryManipulate(typename GeometricField< Type, fvPatchField, volMesh >::Boundary &values)
Manipulate based on a boundary field.
Definition: fvMatrix.C:692
tmp< Field< Type > > H(const Field< Type > &) const
void operator-=(const lduMatrix &)
void operator+=(const lduMatrix &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
T * ptr() const
Return tmp pointer for reuse.
Definition: tmpI.H:198
A class for managing temporary objects.
Definition: PtrList.H:53
scalarField & diag()
Definition: lduMatrix.C:186
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
void deleteDemandDrivenData(DataPtr &dataPtr)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
fvMatrix(const GeometricField< Type, fvPatchField, volMesh > &, const dimensionSet &)
Construct given a field to solve for.
Definition: fvMatrix.C:268
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
#define M(I)
#define InfoInFunction
Report an information message using Foam::Info.
const dimensionSet & dimensions() const
Definition: fvMatrix.H:287