FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

dimensionedScalar.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include <OpenFOAM/dimensionedScalar.H>
00027 
00028 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00029 
00030 namespace Foam
00031 {
00032 
00033 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00034 
00035 dimensionedScalar operator+(const dimensionedScalar& ds1, const scalar s2)
00036 {
00037     return ds1 + dimensionedScalar(s2);
00038 }
00039 
00040 dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2)
00041 {
00042     return dimensionedScalar(s1) + ds2;
00043 }
00044 
00045 dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2)
00046 {
00047     return ds1 - dimensionedScalar(s2);
00048 }
00049 
00050 dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2)
00051 {
00052     return dimensionedScalar(s1) - ds2;
00053 }
00054 
00055 dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2)
00056 {
00057     return ds1 * dimensionedScalar(s2);
00058 }
00059 
00060 dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2)
00061 {
00062     return dimensionedScalar(s1)/ds2;
00063 }
00064 
00065 
00066 dimensionedScalar pow
00067 (
00068     const dimensionedScalar& ds,
00069     const dimensionedScalar& expt
00070 )
00071 {
00072     return dimensionedScalar
00073     (
00074         "pow(" + ds.name() + ',' + expt.name() + ')',
00075         pow(ds.dimensions(), expt),
00076         ::pow(ds.value(), expt.value())
00077     );
00078 }
00079 
00080 dimensionedScalar pow3(const dimensionedScalar& ds)
00081 {
00082     return dimensionedScalar
00083     (
00084         "pow3(" + ds.name() + ')',
00085         pow3(ds.dimensions()),
00086         pow3(ds.value())
00087     );
00088 }
00089 
00090 dimensionedScalar pow4(const dimensionedScalar& ds)
00091 {
00092     return dimensionedScalar
00093     (
00094         "pow4(" + ds.name() + ')',
00095         pow4(ds.dimensions()),
00096         pow4(ds.value())
00097     );
00098 }
00099 
00100 dimensionedScalar pow5(const dimensionedScalar& ds)
00101 {
00102     return dimensionedScalar
00103     (
00104         "pow5(" + ds.name() + ')',
00105         pow5(ds.dimensions()),
00106         pow5(ds.value())
00107     );
00108 }
00109 
00110 dimensionedScalar pow6(const dimensionedScalar& ds)
00111 {
00112     return dimensionedScalar
00113     (
00114         "pow6(" + ds.name() + ')',
00115         pow6(ds.dimensions()),
00116         pow6(ds.value())
00117     );
00118 }
00119 
00120 dimensionedScalar sqrt(const dimensionedScalar& ds)
00121 {
00122     return dimensionedScalar
00123     (
00124         "sqrt(" + ds.name() + ')',
00125         pow(ds.dimensions(), dimensionedScalar("0.5", dimless, 0.5)),
00126         ::sqrt(ds.value())
00127     );
00128 }
00129 
00130 dimensionedScalar cbrt(const dimensionedScalar& ds)
00131 {
00132     return dimensionedScalar
00133     (
00134         "cbrt(" + ds.name() + ')',
00135         pow(ds.dimensions(), dimensionedScalar("(1|3)", dimless, 1.0/3.0)),
00136         ::cbrt(ds.value())
00137     );
00138 }
00139 
00140 dimensionedScalar hypot
00141 (
00142     const dimensionedScalar& x,
00143     const dimensionedScalar& y
00144 )
00145 {
00146     return dimensionedScalar
00147     (
00148         "hypot(" + x.name() + ',' + y.name() + ')',
00149         x.dimensions() + y.dimensions(),
00150         ::hypot(x.value(), y.value())
00151     );
00152 }
00153 
00154 dimensionedScalar sign(const dimensionedScalar& ds)
00155 {
00156     return dimensionedScalar
00157     (
00158         "sign(" + ds.name() + ')',
00159         sign(ds.dimensions()),
00160         ::Foam::sign(ds.value())
00161     );
00162 }
00163 
00164 dimensionedScalar pos(const dimensionedScalar& ds)
00165 {
00166     return dimensionedScalar
00167     (
00168         "pos(" + ds.name() + ')',
00169         pos(ds.dimensions()),
00170         ::Foam::pos(ds.value())
00171     );
00172 }
00173 
00174 dimensionedScalar neg(const dimensionedScalar& ds)
00175 {
00176     return dimensionedScalar
00177     (
00178         "neg(" + ds.name() + ')',
00179         neg(ds.dimensions()),
00180         ::Foam::neg(ds.value())
00181     );
00182 }
00183 
00184 
00185 #define transFunc(func)                                                    \
00186 dimensionedScalar func(const dimensionedScalar& ds)                        \
00187 {                                                                          \
00188     if (!ds.dimensions().dimensionless())                                  \
00189     {                                                                      \
00190         FatalErrorIn(#func "(const dimensionedScalar& ds)")                \
00191             << "ds not dimensionless"                                      \
00192             << abort(FatalError);                                          \
00193     }                                                                      \
00194                                                                            \
00195     return dimensionedScalar                                               \
00196     (                                                                      \
00197         #func "(" + ds.name() + ')',                                       \
00198         dimless,                                                           \
00199         ::func(ds.value())                                                 \
00200     );                                                                     \
00201 }
00202 
00203 transFunc(exp)
00204 transFunc(log)
00205 transFunc(log10)
00206 transFunc(sin)
00207 transFunc(cos)
00208 transFunc(tan)
00209 transFunc(asin)
00210 transFunc(acos)
00211 transFunc(atan)
00212 transFunc(sinh)
00213 transFunc(cosh)
00214 transFunc(tanh)
00215 transFunc(asinh)
00216 transFunc(acosh)
00217 transFunc(atanh)
00218 transFunc(erf)
00219 transFunc(erfc)
00220 transFunc(lgamma)
00221 transFunc(j0)
00222 transFunc(j1)
00223 transFunc(y0)
00224 transFunc(y1)
00225 
00226 #undef transFunc
00227 
00228 
00229 #define transFunc(func)                                                    \
00230 dimensionedScalar func(const int n, const dimensionedScalar& ds)           \
00231 {                                                                          \
00232     if (!ds.dimensions().dimensionless())                                  \
00233     {                                                                      \
00234         FatalErrorIn(#func "(const int n, const dimensionedScalar& ds)")   \
00235             << "ds not dimensionless"                                      \
00236             << abort(FatalError);                                          \
00237     }                                                                      \
00238                                                                            \
00239     return dimensionedScalar                                               \
00240     (                                                                      \
00241         #func "(" + name(n) + ',' + ds.name() + ')',                      \
00242         dimless,                                                           \
00243         ::func(n, ds.value())                                              \
00244     );                                                                     \
00245 }
00246 
00247 transFunc(jn)
00248 transFunc(yn)
00249 
00250 #undef transFunc
00251 
00252 
00253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00254 
00255 } // End namespace Foam
00256 
00257 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines