FreeFOAM The Cross-Platform CFD Toolkit
errorEstimate.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) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "errorEstimate.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 
33 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
34 
35 template<class Type>
37 {
38  // Make the boundary condition type list
39  // Default types get over-ridden anyway
40  wordList ebct
41  (
42  psi_.boundaryField().size(),
43  zeroGradientFvPatchField<Type>::typeName
44  );
45 
46  forAll (psi_.boundaryField(), patchI)
47  {
48  if (psi_.boundaryField()[patchI].fixesValue())
49  {
50  ebct[patchI] = fixedValueFvPatchField<Type>::typeName;
51  }
52  }
53 
54  return ebct;
55 }
56 
57 
58 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
59 
60 // Construct from components
61 template<class Type>
63 (
65  const dimensionSet& ds,
66  const Field<Type>& res,
67  const scalarField& norm
68 )
69 :
70  psi_(psi),
71  dimensions_(ds),
72  residual_(res),
73  normFactor_(norm)
74 {}
75 
76 
77 // Construct as copy
78 template<class Type>
80 :
81  refCount(),
82  psi_(ee.psi_),
83  dimensions_(ee.dimensions_),
84  residual_(ee.residual_),
85  normFactor_(ee.normFactor_)
86 {}
87 
88 
89 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
90 
91 template<class Type>
93 {}
94 
95 
96 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
97 
98 template<class Type>
101 {
103  (
105  (
106  IOobject
107  (
108  "residual" + psi_.name(),
109  psi_.mesh().time().timeName(),
110  psi_.db(),
113  ),
114  psi_.mesh(),
115  psi_.dimensions()/dimTime,
116  errorBCTypes()
117  )
118  );
119 
121 
122  res.internalField() = residual_;
124 
126 
127  return tres;
128 }
129 
130 
131 template<class Type>
133 {
134  tmp<volScalarField> tnormFactor
135  (
136  new volScalarField
137  (
138  IOobject
139  (
140  "normFactor" + psi_.name(),
141  psi_.mesh().time().timeName(),
142  psi_.db(),
145  ),
146  psi_.mesh(),
148  errorBCTypes()
149  )
150  );
151 
152  volScalarField& normFactor = tnormFactor();
153 
154  normFactor.internalField() = normFactor_;
155  normFactor.boundaryField() == pTraits<Type>::zero;
156 
157  normFactor.correctBoundaryConditions();
158 
159  return tnormFactor;
160 }
161 
162 template<class Type>
165 {
167  (
169  (
170  IOobject
171  (
172  "resError" + psi_.name(),
173  psi_.mesh().time().timeName(),
174  psi_.db(),
177  ),
178  psi_.mesh(),
179  psi_.dimensions(),
180  errorBCTypes()
181  )
182  );
183 
184  GeometricField<Type, fvPatchField, volMesh>& resError = tresError();
185 
186  resError.internalField() = residual_/normFactor_;
187  resError.boundaryField() == pTraits<Type>::zero;
188 
189  resError.correctBoundaryConditions();
190 
191  return tresError;
192 }
193 
194 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
195 
196 template<class Type>
198 {
199  // Check for assignment to self
200  if (this == &rhs)
201  {
203  (
204  "errorEstimate<Type>::operator=(const Foam::errorEstimate<Type>&)"
205  ) << "Attempted assignment to self"
206  << abort(FatalError);
207  }
208 
209  if (&psi_ != &(rhs.psi_))
210  {
212  (
213  "errorEstimate<Type>::operator=(const errorEstimate<Type>&)"
214  ) << "different fields"
215  << abort(FatalError);
216  }
217 
218  residual_ = rhs.residual_;
219  normFactor_ = rhs.normFactor_;
220 }
221 
222 
223 template<class Type>
225 {
226  operator=(teev());
227  teev.clear();
228 }
229 
230 
231 template<class Type>
233 {
234  residual_.negate();
235 }
236 
237 
238 template<class Type>
240 {
241  checkMethod(*this, eev, "+=");
242 
243  dimensions_ += eev.dimensions_;
244 
245  residual_ += eev.residual_;
246  normFactor_ += eev.normFactor_;
247 }
248 
249 
250 template<class Type>
252 (
253  const tmp<errorEstimate<Type> >& teev
254 )
255 {
256  operator+=(teev());
257  teev.clear();
258 }
259 
260 
261 template<class Type>
263 {
264  checkMethod(*this, eev, "+=");
265 
266  dimensions_ -= eev.dimensions_;
267  residual_ -= eev.residual_;
268  normFactor_ += eev.normFactor_;
269 }
270 
271 
272 template<class Type>
274 {
275  operator-=(teev());
276  teev.clear();
277 }
278 
279 
280 template<class Type>
281 void Foam::errorEstimate<Type>::operator+=
282 (
284 )
285 {
286  checkMethod(*this, su, "+=");
287  residual_ -= su.internalField();
288 }
289 
290 
291 template<class Type>
292 void Foam::errorEstimate<Type>::operator+=
293 (
295 )
296 {
297  operator+=(tsu());
298  tsu.clear();
299 }
300 
301 
302 template<class Type>
303 void Foam::errorEstimate<Type>::operator-=
304 (
306 )
307 {
308  checkMethod(*this, su, "-=");
309  residual_ += su.internalField();
310 }
311 
312 
313 template<class Type>
314 void Foam::errorEstimate<Type>::operator-=
315 (
317 )
318 {
319  operator-=(tsu());
320  tsu.clear();
321 }
322 
323 
324 template<class Type>
325 void Foam::errorEstimate<Type>::operator+=
326 (
327  const dimensioned<Type>& su
328 )
329 {
330  residual_ -= su;
331 }
332 
333 
334 template<class Type>
335 void Foam::errorEstimate<Type>::operator-=
336 (
337  const dimensioned<Type>& su
338 )
339 {
340  residual_ += su;
341 }
342 
343 
344 template<class Type>
345 void Foam::errorEstimate<Type>::operator*=
346 (
347  const volScalarField& vsf
348 )
349 {
350  dimensions_ *= vsf.dimensions();
351  residual_ *= vsf.internalField();
352  normFactor_ *= vsf.internalField();
353 }
354 
355 
356 template<class Type>
357 void Foam::errorEstimate<Type>::operator*=
358 (
359  const tmp<volScalarField>& tvsf
360 )
361 {
362  operator*=(tvsf());
363  tvsf.clear();
364 }
365 
366 
367 template<class Type>
368 void Foam::errorEstimate<Type>::operator*=
369 (
370  const dimensioned<scalar>& ds
371 )
372 {
373  dimensions_ *= ds.dimensions();
374  residual_ *= ds.value();
375  normFactor_ *= ds.value();
376 }
377 
378 
379 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
380 
381 template<class Type>
383 (
384  const errorEstimate<Type>& ee1,
385  const errorEstimate<Type>& ee2,
386  const char* op
387 )
388 {
389  if (&ee1.psi() != &ee2.psi())
390  {
392  (
393  "checkMethod(const errorEstimate<Type>&, "
394  "const errorEstimate<Type>&)"
395  ) << "incompatible fields for operation "
396  << endl << " "
397  << "[" << ee1.psi().name() << "] "
398  << op
399  << " [" << ee2.psi().name() << "]"
400  << abort(FatalError);
401  }
402 
403  if (dimensionSet::debug && ee1.dimensions() != ee2.dimensions())
404  {
406  (
407  "checkMethod(const errorEstimate<Type>&, "
408  "const errorEstimate<Type>&)"
409  ) << "incompatible dimensions for operation "
410  << endl << " "
411  << "[" << ee1.psi().name() << ee1.dimensions()/dimVolume << " ] "
412  << op
413  << " [" << ee2.psi().name() << ee2.dimensions()/dimVolume << " ]"
414  << abort(FatalError);
415  }
416 }
417 
418 
419 template<class Type>
421 (
422  const errorEstimate<Type>& ee,
424  const char* op
425 )
426 {
427  if (dimensionSet::debug && ee.dimensions()/dimVolume != vf.dimensions())
428  {
430  (
431  "checkMethod(const errorEstimate<Type>&, "
432  "const GeometricField<Type, fvPatchField, volMesh>&)"
433  ) << "incompatible dimensions for operation "
434  << endl << " "
435  << "[" << ee.psi().name() << ee.dimensions()/dimVolume << " ] "
436  << op
437  << " [" << vf.name() << vf.dimensions() << " ]"
438  << abort(FatalError);
439  }
440 }
441 
442 
443 template<class Type>
445 (
446  const errorEstimate<Type>& ee,
447  const dimensioned<Type>& dt,
448  const char* op
449 )
450 {
451  if (dimensionSet::debug && ee.dimensions()/dimVolume != dt.dimensions())
452  {
454  (
455  "checkMethod(const errorEstimate<Type>&, const dimensioned<Type>&)"
456  ) << "incompatible dimensions for operation "
457  << endl << " "
458  << "[" << ee.psi().name() << ee.dimensions()/dimVolume << " ] "
459  << op
460  << " [" << dt.name() << dt.dimensions() << " ]"
461  << abort(FatalError);
462  }
463 }
464 
465 
466 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
467 
468 namespace Foam
469 {
470 
471 template<class Type>
472 tmp<errorEstimate<Type> > operator+
473 (
474  const errorEstimate<Type>& A,
475  const errorEstimate<Type>& B
476 )
477 {
478  checkMethod(A, B, "+");
480  tC() += B;
481  return tC;
482 }
483 
484 
485 template<class Type>
486 tmp<errorEstimate<Type> > operator+
487 (
488  const tmp<errorEstimate<Type> >& tA,
489  const errorEstimate<Type>& B
490 )
491 {
492  checkMethod(tA(), B, "+");
493  tmp<errorEstimate<Type> > tC(tA.ptr());
494  tC() += B;
495  return tC;
496 }
497 
498 
499 template<class Type>
500 tmp<errorEstimate<Type> > operator+
501 (
502  const errorEstimate<Type>& A,
503  const tmp<errorEstimate<Type> >& tB
504 )
505 {
506  checkMethod(A, tB(), "+");
507  tmp<errorEstimate<Type> > tC(tB.ptr());
508  tC() += A;
509  return tC;
510 }
511 
512 
513 template<class Type>
514 tmp<errorEstimate<Type> > operator+
515 (
516  const tmp<errorEstimate<Type> >& tA,
517  const tmp<errorEstimate<Type> >& tB
518 )
519 {
520  checkMethod(tA(), tB(), "+");
521  tmp<errorEstimate<Type> > tC(tA.ptr());
522  tC() += tB();
523  tB.clear();
524  return tC;
525 }
526 
527 
528 template<class Type>
529 tmp<errorEstimate<Type> > operator-
530 (
531  const errorEstimate<Type>& A
532 )
533 {
535  tC().negate();
536  return tC;
537 }
538 
539 
540 template<class Type>
541 tmp<errorEstimate<Type> > operator-
542 (
543  const tmp<errorEstimate<Type> >& tA
544 )
545 {
546  tmp<errorEstimate<Type> > tC(tA.ptr());
547  tC().negate();
548  return tC;
549 }
550 
551 
552 template<class Type>
553 tmp<errorEstimate<Type> > operator-
554 (
555  const errorEstimate<Type>& A,
556  const errorEstimate<Type>& B
557 )
558 {
559  checkMethod(A, B, "-");
561  tC() -= B;
562  return tC;
563 }
564 
565 
566 template<class Type>
567 tmp<errorEstimate<Type> > operator-
568 (
569  const tmp<errorEstimate<Type> >& tA,
570  const errorEstimate<Type>& B
571 )
572 {
573  checkMethod(tA(), B, "-");
574  tmp<errorEstimate<Type> > tC(tA.ptr());
575  tC() -= B;
576  return tC;
577 }
578 
579 
580 template<class Type>
581 tmp<errorEstimate<Type> > operator-
582 (
583  const errorEstimate<Type>& A,
584  const tmp<errorEstimate<Type> >& tB
585 )
586 {
587  checkMethod(A, tB(), "-");
588  tmp<errorEstimate<Type> > tC(tB.ptr());
589  tC() -= A;
590  tC().negate();
591  return tC;
592 }
593 
594 
595 template<class Type>
596 tmp<errorEstimate<Type> > operator-
597 (
598  const tmp<errorEstimate<Type> >& tA,
599  const tmp<errorEstimate<Type> >& tB
600 )
601 {
602  checkMethod(tA(), tB(), "-");
603  tmp<errorEstimate<Type> > tC(tA.ptr());
604  tC() -= tB();
605  tB.clear();
606  return tC;
607 }
608 
609 
610 template<class Type>
611 tmp<errorEstimate<Type> > operator==
612 (
613  const errorEstimate<Type>& A,
614  const errorEstimate<Type>& B
615 )
616 {
617  checkMethod(A, B, "==");
618  return (A - B);
619 }
620 
621 
622 template<class Type>
623 tmp<errorEstimate<Type> > operator==
624 (
625  const tmp<errorEstimate<Type> >& tA,
626  const errorEstimate<Type>& B
627 )
628 {
629  checkMethod(tA(), B, "==");
630  return (tA - B);
631 }
632 
633 
634 template<class Type>
635 tmp<errorEstimate<Type> > operator==
636 (
637  const errorEstimate<Type>& A,
638  const tmp<errorEstimate<Type> >& tB
639 )
640 {
641  checkMethod(A, tB(), "==");
642  return (A - tB);
643 }
644 
645 
646 template<class Type>
647 tmp<errorEstimate<Type> > operator==
648 (
649  const tmp<errorEstimate<Type> >& tA,
650  const tmp<errorEstimate<Type> >& tB
651 )
652 {
653  checkMethod(tA(), tB(), "==");
654  return (tA - tB);
655 }
656 
657 
658 template<class Type>
659 tmp<errorEstimate<Type> > operator+
660 (
661  const errorEstimate<Type>& A,
663 )
664 {
665  checkMethod(A, su, "+");
667  tC().res() -= su.internalField();
668  return tC;
669 }
670 
671 template<class Type>
672 tmp<errorEstimate<Type> > operator+
673 (
674  const tmp<errorEstimate<Type> >& tA,
676 )
677 {
678  checkMethod(tA(), su, "+");
679  tmp<errorEstimate<Type> > tC(tA.ptr());
680  tC().res() -= su.internalField();
681  return tC;
682 }
683 
684 template<class Type>
685 tmp<errorEstimate<Type> > operator+
686 (
687  const errorEstimate<Type>& A,
689 )
690 {
691  checkMethod(A, tsu(), "+");
693  tC().res() -= tsu().internalField();
694  tsu.clear();
695  return tC;
696 }
697 
698 
699 template<class Type>
700 tmp<errorEstimate<Type> > operator+
701 (
702  const tmp<errorEstimate<Type> >& tA,
704 )
705 {
706  checkMethod(tA(), tsu(), "+");
707  tmp<errorEstimate<Type> > tC(tA.ptr());
708  tC().res() -= tsu().internalField();
709  tsu.clear();
710  return tC;
711 }
712 
713 template<class Type>
714 tmp<errorEstimate<Type> > operator+
715 (
717  const errorEstimate<Type>& A
718 )
719 {
720  checkMethod(A, su, "+");
722  tC().res() -= su.internalField();
723  return tC;
724 }
725 
726 template<class Type>
727 tmp<errorEstimate<Type> > operator+
728 (
730  const tmp<errorEstimate<Type> >& tA
731 )
732 {
733  checkMethod(tA(), su, "+");
734  tmp<errorEstimate<Type> > tC(tA.ptr());
735  tC().res() -= su.internalField();
736  return tC;
737 }
738 
739 template<class Type>
740 tmp<errorEstimate<Type> > operator+
741 (
743  const errorEstimate<Type>& A
744 )
745 {
746  checkMethod(A, tsu(), "+");
748  tC().res() -= tsu().internalField();
749  tsu.clear();
750  return tC;
751 }
752 
753 template<class Type>
754 tmp<errorEstimate<Type> > operator+
755 (
757  const tmp<errorEstimate<Type> >& tA
758 )
759 {
760  checkMethod(tA(), tsu(), "+");
761  tmp<errorEstimate<Type> > tC(tA.ptr());
762  tC().res() -= tsu().internalField();
763  tsu.clear();
764  return tC;
765 }
766 
767 
768 template<class Type>
769 tmp<errorEstimate<Type> > operator-
770 (
771  const errorEstimate<Type>& A,
773 )
774 {
775  checkMethod(A, su, "-");
777  tC().res() += su.internalField();
778  return tC;
779 }
780 
781 template<class Type>
782 tmp<errorEstimate<Type> > operator-
783 (
784  const tmp<errorEstimate<Type> >& tA,
786 )
787 {
788  checkMethod(tA(), su, "-");
789  tmp<errorEstimate<Type> > tC(tA.ptr());
790  tC().res() += su.internalField();
791  return tC;
792 }
793 
794 template<class Type>
795 tmp<errorEstimate<Type> > operator-
796 (
797  const errorEstimate<Type>& A,
799 )
800 {
801  checkMethod(A, tsu(), "-");
803  tC().res() += tsu().internalField();
804  tsu.clear();
805  return tC;
806 }
807 
808 template<class Type>
809 tmp<errorEstimate<Type> > operator-
810 (
811  const tmp<errorEstimate<Type> >& tA,
813 )
814 {
815  checkMethod(tA(), tsu(), "-");
816  tmp<errorEstimate<Type> > tC(tA.ptr());
817  tC().res() += tsu().internalField();
818  tsu.clear();
819  return tC;
820 }
821 
822 
823 template<class Type>
824 tmp<errorEstimate<Type> > operator-
825 (
827  const errorEstimate<Type>& A
828 )
829 {
830  checkMethod(A, su, "-");
832  tC().negate();
833  tC().res() -= su.internalField();
834  return tC;
835 }
836 
837 
838 template<class Type>
839 tmp<errorEstimate<Type> > operator-
840 (
842  const tmp<errorEstimate<Type> >& tA
843 )
844 {
845  checkMethod(tA(), su, "-");
846  tmp<errorEstimate<Type> > tC(tA.ptr());
847  tC().negate();
848  tC().res() -= su.internalField();
849  return tC;
850 }
851 
852 template<class Type>
853 tmp<errorEstimate<Type> > operator-
854 (
856  const errorEstimate<Type>& A
857 )
858 {
859  checkMethod(A, tsu(), "-");
861  tC().negate();
862  tC().res() -= tsu().internalField();
863  tsu.clear();
864  return tC;
865 }
866 
867 
868 template<class Type>
869 tmp<errorEstimate<Type> > operator-
870 (
872  const tmp<errorEstimate<Type> >& tA
873 )
874 {
875  checkMethod(tA(), tsu(), "-");
876  tmp<errorEstimate<Type> > tC(tA.ptr());
877  tC().negate();
878  tC().res() -= tsu().internalField();
879  tsu.clear();
880  return tC;
881 }
882 
883 
884 template<class Type>
885 tmp<errorEstimate<Type> > operator+
886 (
887  const errorEstimate<Type>& A,
888  const dimensioned<Type>& su
889 )
890 {
891  checkMethod(A, su, "+");
893  tC().res() -= su.value();
894  return tC;
895 }
896 
897 
898 template<class Type>
899 tmp<errorEstimate<Type> > operator+
900 (
901  const tmp<errorEstimate<Type> >& tA,
902  const dimensioned<Type>& su
903 )
904 {
905  checkMethod(tA(), su, "+");
906  tmp<errorEstimate<Type> > tC(tA.ptr());
907  tC().res() -= su.value();
908  return tC;
909 }
910 
911 
912 template<class Type>
913 tmp<errorEstimate<Type> > operator+
914 (
915  const dimensioned<Type>& su,
916  const errorEstimate<Type>& A
917 )
918 {
919  checkMethod(A, su, "+");
921  tC().res() -= su.value();
922  return tC;
923 }
924 
925 
926 template<class Type>
927 tmp<errorEstimate<Type> > operator+
928 (
929  const dimensioned<Type>& su,
930  const tmp<errorEstimate<Type> >& tA
931 )
932 {
933  checkMethod(tA(), su, "+");
934  tmp<errorEstimate<Type> > tC(tA.ptr());
935  tC().res() -= su.value();
936  return tC;
937 }
938 
939 
940 template<class Type>
941 tmp<errorEstimate<Type> > operator-
942 (
943  const errorEstimate<Type>& A,
944  const dimensioned<Type>& su
945 )
946 {
947  checkMethod(A, su, "-");
949  tC().res() += su.value();
950  return tC;
951 }
952 
953 
954 template<class Type>
955 tmp<errorEstimate<Type> > operator-
956 (
957  const tmp<errorEstimate<Type> >& tA,
958  const dimensioned<Type>& su
959 )
960 {
961  checkMethod(tA(), su, "-");
962  tmp<errorEstimate<Type> > tC(tA.ptr());
963  tC().res() += su.value();
964  return tC;
965 }
966 
967 
968 template<class Type>
969 tmp<errorEstimate<Type> > operator-
970 (
971  const dimensioned<Type>& su,
972  const errorEstimate<Type>& A
973 )
974 {
975  checkMethod(A, su, "-");
977  tC().negate();
978  tC().res() -= su.value();
979  return tC;
980 }
981 
982 
983 template<class Type>
984 tmp<errorEstimate<Type> > operator-
985 (
986  const dimensioned<Type>& su,
987  const tmp<errorEstimate<Type> >& tA
988 )
989 {
990  checkMethod(tA(), su, "-");
991  tmp<errorEstimate<Type> > tC(tA.ptr());
992  tC().negate();
993  tC().res() -= su.value();
994  return tC;
995 }
996 
997 
998 template<class Type>
999 tmp<errorEstimate<Type> > operator==
1001  const errorEstimate<Type>& A,
1003 )
1004 {
1005  checkMethod(A, su, "==");
1007  tC().res() += su.internalField();
1008  return tC;
1009 }
1010 
1011 template<class Type>
1012 tmp<errorEstimate<Type> > operator==
1014  const tmp<errorEstimate<Type> >& tA,
1016 )
1017 {
1018  checkMethod(tA(), su, "==");
1019  tmp<errorEstimate<Type> > tC(tA.ptr());
1020  tC().res() += su.internalField();
1021  return tC;
1022 }
1023 
1024 template<class Type>
1025 tmp<errorEstimate<Type> > operator==
1027  const errorEstimate<Type>& A,
1029 )
1030 {
1031  checkMethod(A, tsu(), "==");
1033  tC().res() += tsu().internalField();
1034  tsu.clear();
1035  return tC;
1036 }
1037 
1038 template<class Type>
1039 tmp<errorEstimate<Type> > operator==
1041  const tmp<errorEstimate<Type> >& tA,
1043 )
1044 {
1045  checkMethod(tA(), tsu(), "==");
1046  tmp<errorEstimate<Type> > tC(tA.ptr());
1047  tC().res() += tsu().internalField();
1048  tsu.clear();
1049  return tC;
1050 }
1051 
1052 
1053 template<class Type>
1054 tmp<errorEstimate<Type> > operator==
1056  const errorEstimate<Type>& A,
1057  const dimensioned<Type>& su
1058 )
1059 {
1060  checkMethod(A, su, "==");
1062  tC().res() += su.value();
1063  return tC;
1064 }
1065 
1066 
1067 template<class Type>
1068 tmp<errorEstimate<Type> > operator==
1070  const tmp<errorEstimate<Type> >& tA,
1071  const dimensioned<Type>& su
1072 )
1073 {
1074  checkMethod(tA(), su, "==");
1075  tmp<errorEstimate<Type> > tC(tA.ptr());
1076  tC().res() += su.value();
1077  return tC;
1078 }
1079 
1080 
1081 template<class Type>
1082 tmp<errorEstimate<Type> > operator*
1084  const volScalarField& vsf,
1085  const errorEstimate<Type>& A
1086 )
1087 {
1089  tC() *= vsf;
1090  return tC;
1091 }
1092 
1093 template<class Type>
1094 tmp<errorEstimate<Type> > operator*
1096  const tmp<volScalarField>& tvsf,
1097  const errorEstimate<Type>& A
1098 )
1099 {
1101  tC() *= tvsf;
1102  return tC;
1103 }
1104 
1105 template<class Type>
1106 tmp<errorEstimate<Type> > operator*
1108  const volScalarField& vsf,
1109  const tmp<errorEstimate<Type> >& tA
1110 )
1111 {
1112  tmp<errorEstimate<Type> > tC(tA.ptr());
1113  tC() *= vsf;
1114  return tC;
1115 }
1116 
1117 template<class Type>
1118 tmp<errorEstimate<Type> > operator*
1120  const tmp<volScalarField>& tvsf,
1121  const tmp<errorEstimate<Type> >& tA
1122 )
1123 {
1124  tmp<errorEstimate<Type> > tC(tA.ptr());
1125  tC() *= tvsf;
1126  return tC;
1127 }
1128 
1129 
1130 template<class Type>
1131 tmp<errorEstimate<Type> > operator*
1133  const dimensioned<scalar>& ds,
1134  const errorEstimate<Type>& A
1135 )
1136 {
1138  tC() *= ds;
1139  return tC;
1140 }
1141 
1142 
1143 template<class Type>
1144 tmp<errorEstimate<Type> > operator*
1146  const dimensioned<scalar>& ds,
1147  const tmp<errorEstimate<Type> >& tA
1148 )
1149 {
1150  tmp<errorEstimate<Type> > tC(tA.ptr());
1151  tC() *= ds;
1152  return tC;
1153 }
1154 
1155 
1156 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1157 
1158 } // End namespace Foam
1159 
1160 // ************************ vim: set sw=4 sts=4 et: ************************ //