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