Sacado Package Browser (Single Doxygen Collection) Version of the Day
Sacado_Tay_ScalarTraitsImp.hpp
Go to the documentation of this file.
00001 // $Id$ 
00002 // $Source$ 
00003 // @HEADER
00004 // ***********************************************************************
00005 // 
00006 //                           Sacado Package
00007 //                 Copyright (2006) Sandia Corporation
00008 // 
00009 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00010 // the U.S. Government retains certain rights in this software.
00011 // 
00012 // This library is free software; you can redistribute it and/or modify
00013 // it under the terms of the GNU Lesser General Public License as
00014 // published by the Free Software Foundation; either version 2.1 of the
00015 // License, or (at your option) any later version.
00016 //  
00017 // This library is distributed in the hope that it will be useful, but
00018 // WITHOUT ANY WARRANTY; without even the implied warranty of
00019 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00020 // Lesser General Public License for more details.
00021 //  
00022 // You should have received a copy of the GNU Lesser General Public
00023 // License along with this library; if not, write to the Free Software
00024 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00025 // USA
00026 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
00027 // (etphipp@sandia.gov).
00028 // 
00029 // ***********************************************************************
00030 // @HEADER
00031 
00032 #ifndef SACADO_TAY_SCALARTRAITSIMP_HPP
00033 #define SACADO_TAY_SCALARTRAITSIMP_HPP
00034 
00035 #ifdef HAVE_SACADO_TEUCHOS
00036 
00037 #include "Teuchos_ScalarTraits.hpp"
00038 #include "Teuchos_SerializationTraits.hpp"
00039 #include "Teuchos_SerializationTraitsHelpers.hpp"
00040 #include "Teuchos_Assert.hpp"
00041 #include "Sacado_mpl_apply.hpp"
00042 
00043 #include <iterator>
00044 
00045 namespace Sacado {
00046 
00047   namespace Tay {
00048 
00050     template <typename TayType>
00051     struct ScalarTraitsImp {
00052       typedef typename Sacado::ValueType<TayType>::type ValueT;
00053 
00054       typedef typename mpl::apply<TayType,typename Teuchos::ScalarTraits<ValueT>::magnitudeType>::type magnitudeType;
00055       typedef typename mpl::apply<TayType,typename Teuchos::ScalarTraits<ValueT>::halfPrecision>::type halfPrecision;
00056       typedef typename mpl::apply<TayType,typename Teuchos::ScalarTraits<ValueT>::doublePrecision>::type doublePrecision;
00057 
00058       static const bool isComplex = Teuchos::ScalarTraits<ValueT>::isComplex;
00059       static const bool isOrdinal = Teuchos::ScalarTraits<ValueT>::isOrdinal;
00060       static const bool isComparable = 
00061   Teuchos::ScalarTraits<ValueT>::isComparable;
00062       static const bool hasMachineParameters = 
00063   Teuchos::ScalarTraits<ValueT>::hasMachineParameters;
00064       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType eps() {
00065   return Teuchos::ScalarTraits<ValueT>::eps();
00066       }
00067       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType sfmin() {
00068   return Teuchos::ScalarTraits<ValueT>::sfmin();
00069       }
00070       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType base()  {
00071   return Teuchos::ScalarTraits<ValueT>::base();
00072       }
00073       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType prec()  {
00074   return Teuchos::ScalarTraits<ValueT>::prec();
00075       }
00076       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType t()     {
00077   return Teuchos::ScalarTraits<ValueT>::t();
00078       }
00079       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType rnd()   {
00080   return Teuchos::ScalarTraits<ValueT>::rnd();
00081       }
00082       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType emin()  {
00083   return Teuchos::ScalarTraits<ValueT>::emin();
00084       }
00085       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType rmin()  {
00086   return Teuchos::ScalarTraits<ValueT>::rmin();
00087       }
00088       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType emax()  {
00089   return Teuchos::ScalarTraits<ValueT>::emax();
00090       }
00091       static typename Teuchos::ScalarTraits<ValueT>::magnitudeType rmax()  {
00092   return Teuchos::ScalarTraits<ValueT>::rmax();
00093       }
00094       static magnitudeType magnitude(const TayType& a) {
00095 #ifdef TEUCHOS_DEBUG
00096   TEUCHOS_SCALAR_TRAITS_NAN_INF_ERR(
00097     a, "Error, the input value to magnitude(...) a = " << a << 
00098     " can not be NaN!" );
00099   TEUCHOS_TEST_FOR_EXCEPTION(is_tay_real(a) == false, std::runtime_error,
00100          "Complex magnitude is not a differentiable "
00101          "function of complex inputs.");
00102 #endif
00103   //return std::fabs(a); 
00104   magnitudeType b(a.degree(), 
00105       Teuchos::ScalarTraits<ValueT>::magnitude(a.val()));
00106   if (Teuchos::ScalarTraits<ValueT>::real(a.val()) >= 0)
00107     for (int i=1; i<=a.degree(); i++)
00108       b.fastAccessCoeff(i) = 
00109         Teuchos::ScalarTraits<ValueT>::magnitude(a.fastAccessCoeff(i));
00110   else
00111     for (int i=1; i<=a.degree(); i++)
00112       b.fastAccessCoeff(i) = 
00113         -Teuchos::ScalarTraits<ValueT>::magnitude(a.fastAccessCoeff(i));
00114   return b;
00115       }
00116       static ValueT zero()  { 
00117   return ValueT(0.0); 
00118       }
00119       static ValueT one()   { 
00120   return ValueT(1.0); 
00121       }
00122       
00123       // Conjugate is only defined for real derivative components
00124       static TayType conjugate(const TayType& x) {
00125 #ifdef TEUCHOS_DEBUG
00126   TEUCHOS_TEST_FOR_EXCEPTION(is_tay_real(x) == false, std::runtime_error,
00127          "Complex conjugate is not a differentiable "
00128          "function of complex inputs.");
00129 #endif
00130   TayType y = x;
00131   y.copyForWrite();
00132   y.val() = Teuchos::ScalarTraits<ValueT>::conjugate(x.val());
00133   return y;
00134       }   
00135 
00136       // Real part is only defined for real derivative components
00137       static TayType real(const TayType& x) { 
00138 #ifdef TEUCHOS_DEBUG
00139   TEUCHOS_TEST_FOR_EXCEPTION(is_tay_real(x) == false, std::runtime_error,
00140          "Real component is not a differentiable "
00141          "function of complex inputs.");
00142 #endif
00143   TayType y = x;
00144   y.copyForWrite();
00145   y.val() = Teuchos::ScalarTraits<ValueT>::real(x.val());
00146   return y;
00147       }
00148 
00149       // Imaginary part is only defined for real derivative components
00150       static TayType imag(const TayType& x) { 
00151 #ifdef TEUCHOS_DEBUG
00152   TEUCHOS_TEST_FOR_EXCEPTION(is_tay_real(x) == false, std::runtime_error,
00153          "Imaginary component is not a differentiable "
00154          "function of complex inputs.");
00155 #endif
00156   return TayType(Teuchos::ScalarTraits<ValueT>::imag(x.val()));
00157       }
00158 
00159       static ValueT nan() {
00160   return Teuchos::ScalarTraits<ValueT>::nan(); 
00161       }
00162       static bool isnaninf(const TayType& x) { 
00163   for (int i=0; i<=x.degree(); i++)
00164     if (Teuchos::ScalarTraits<ValueT>::isnaninf(x.fastAccessCoeff(i)))
00165       return true;
00166   return false;
00167       }
00168       static void seedrandom(unsigned int s) { 
00169   Teuchos::ScalarTraits<ValueT>::seedrandom(s); 
00170       }
00171       static ValueT random() { 
00172   return Teuchos::ScalarTraits<ValueT>::random(); 
00173       }
00174       static std::string name() { 
00175   return Sacado::StringName<TayType>::eval(); 
00176       }
00177       static TayType squareroot(const TayType& x) {
00178 #ifdef TEUCHOS_DEBUG
00179   TEUCHOS_SCALAR_TRAITS_NAN_INF_ERR(
00180     x, "Error, the input value to squareroot(...) a = " << x << 
00181     " can not be NaN!" );
00182 #endif
00183   return std::sqrt(x); 
00184       }
00185       static TayType pow(const TayType& x, const TayType& y) { 
00186   return std::pow(x,y); 
00187       }
00188 
00189       // Helper function to determine whether a complex value is real
00190       static bool is_complex_real(const ValueT& x) {
00191   return 
00192     Teuchos::ScalarTraits<ValueT>::magnitude(x-Teuchos::ScalarTraits<ValueT>::real(x)) == 0;
00193       }
00194 
00195       // Helper function to determine whether a Fad type is real
00196       static bool is_tay_real(const TayType& x) {
00197   if (x.size() == 0)
00198     return true;
00199   if (Teuchos::ScalarTraits<ValueT>::isComplex) {
00200     for (int i=0; i<=x.degree(); i++)
00201       if (!is_complex_real(x.fastAccessCoeff(i)))
00202         return false;
00203   }
00204   return true;
00205       }
00206 
00207     }; // class ScalarTraitsImp
00208 
00210     template <typename Ordinal, typename TayType, typename Serializer>
00211     struct SerializationImp {
00212 
00213     private:
00214       
00216       typedef Teuchos::SerializationTraits<Ordinal,unsigned int> iSerT;
00217 
00219       typedef Teuchos::SerializationTraits<Ordinal,Ordinal> oSerT;
00220 
00221     public:
00222 
00224       static const bool supportsDirectSerialization = false;
00225 
00227 
00228 
00230       static Ordinal fromCountToIndirectBytes(const Serializer& vs,
00231                 const Ordinal count, 
00232                 const TayType buffer[],
00233                 const Ordinal sz = 0) { 
00234   Ordinal bytes = 0;
00235   TayType *x = NULL;
00236   const TayType *cx;
00237   for (Ordinal i=0; i<count; i++) {
00238     unsigned int my_sz = buffer[i].degree()+1;
00239     unsigned int tot_sz = sz;
00240     if (sz == 0) tot_sz = my_sz;
00241     if (tot_sz != my_sz) {
00242       x = new TayType(buffer[i]);
00243             x->resize(tot_sz-1, true);
00244             cx = x;
00245     }
00246           else 
00247             cx = &(buffer[i]);
00248     Ordinal b1 = iSerT::fromCountToIndirectBytes(1, &tot_sz);
00249     Ordinal b2 = vs.fromCountToIndirectBytes(
00250       tot_sz, &(cx->fastAccessCoeff(0)));
00251     Ordinal b3 = oSerT::fromCountToIndirectBytes(1, &b2);
00252     bytes += b1+b2+b3;
00253     if (x != NULL) {
00254       delete x;
00255       x = NULL;
00256     }
00257   }
00258   return bytes;
00259       }
00260 
00262       static void serialize (const Serializer& vs,
00263            const Ordinal count, 
00264            const TayType buffer[], 
00265            const Ordinal bytes, 
00266            char charBuffer[],
00267            const Ordinal sz = 0) { 
00268   TayType *x = NULL;
00269   const TayType *cx;
00270   for (Ordinal i=0; i<count; i++) {
00271     // First serialize degree
00272     unsigned int my_sz = buffer[i].degree()+1;
00273     unsigned int tot_sz = sz;
00274     if (sz == 0) tot_sz = my_sz;
00275     Ordinal b1 = iSerT::fromCountToIndirectBytes(1, &tot_sz);
00276     iSerT::serialize(1, &tot_sz, b1, charBuffer);
00277     charBuffer += b1;
00278   
00279     // Next serialize taylor coefficients
00280     if (tot_sz != my_sz) {
00281       x = new TayType(buffer[i]);
00282             x->resize(tot_sz-1, true);
00283             cx = x;
00284     }
00285           else 
00286             cx = &(buffer[i]);
00287     Ordinal b2 = vs.fromCountToIndirectBytes(
00288       tot_sz, &(cx->fastAccessCoeff(0)));
00289     Ordinal b3 = oSerT::fromCountToIndirectBytes(1, &b2);
00290     oSerT::serialize(1, &b2, b3, charBuffer); 
00291     charBuffer += b3;
00292     vs.serialize(tot_sz, &(cx->fastAccessCoeff(0)), b2, charBuffer);
00293     charBuffer += b2;
00294     if (x != NULL) {
00295       delete x;
00296       x = NULL;
00297     }
00298   }
00299       }
00300 
00302       static Ordinal fromIndirectBytesToCount(const Serializer& vs,
00303                 const Ordinal bytes, 
00304                 const char charBuffer[],
00305                 const Ordinal sz = 0) {
00306   Ordinal count = 0;
00307   Ordinal bytes_used = 0;
00308   while (bytes_used < bytes) {
00309   
00310     // Bytes for degree
00311     Ordinal b1 = iSerT::fromCountToDirectBytes(1);
00312     bytes_used += b1;
00313     charBuffer += b1;
00314   
00315     // Bytes for taylor coefficients
00316     Ordinal b3 = oSerT::fromCountToDirectBytes(1);
00317     const Ordinal *b2 = oSerT::convertFromCharPtr(charBuffer);
00318     bytes_used += b3;
00319     charBuffer += b3;
00320     bytes_used += *b2;
00321     charBuffer += *b2;
00322   
00323     ++count;
00324   }
00325   return count;
00326       }
00327 
00329       static void deserialize (const Serializer& vs,
00330              const Ordinal bytes, 
00331              const char charBuffer[], 
00332              const Ordinal count, 
00333              TayType buffer[],
00334              const Ordinal sz = 0) { 
00335   for (Ordinal i=0; i<count; i++) {
00336   
00337     // Deserialize degree
00338     Ordinal b1 = iSerT::fromCountToDirectBytes(1);
00339     const unsigned int *my_sz = iSerT::convertFromCharPtr(charBuffer);
00340     charBuffer += b1;
00341   
00342     // Create empty Taylor object of given size
00343     unsigned int tot_sz = sz;
00344     if (sz == 0) tot_sz = *my_sz;
00345     buffer[i] = TayType(tot_sz-1, 0.0);
00346   
00347     // Deserialize taylor coefficients
00348     Ordinal b3 = oSerT::fromCountToDirectBytes(1);
00349     const Ordinal *b2 = oSerT::convertFromCharPtr(charBuffer);
00350     charBuffer += b3;
00351     vs.deserialize(*b2, charBuffer, *my_sz, 
00352        &(buffer[i].fastAccessCoeff(0)));
00353     charBuffer += *b2;
00354   }
00355       
00356       }
00357   
00359       
00360     };
00361 
00363     template <typename Ordinal, typename TayType>
00364     struct SerializationTraitsImp {
00365 
00366     private:
00367 
00369       typedef typename Sacado::ValueType<TayType>::type ValueT;
00370 
00372       typedef Teuchos::DefaultSerializer<Ordinal,ValueT> DS;
00373 
00375       typedef typename DS::DefaultSerializerType ValueSerializer;
00376 
00378       typedef SerializationImp<Ordinal,TayType,ValueSerializer> Imp;
00379 
00380     public:
00381 
00383       static const bool supportsDirectSerialization = 
00384   Imp::supportsDirectSerialization;
00385 
00387 
00388 
00390       static Ordinal fromCountToIndirectBytes(const Ordinal count, 
00391                 const TayType buffer[]) { 
00392   return Imp::fromCountToIndirectBytes(
00393     DS::getDefaultSerializer(), count, buffer);
00394       }
00395 
00397       static void serialize (const Ordinal count, 
00398            const TayType buffer[], 
00399            const Ordinal bytes, 
00400            char charBuffer[]) { 
00401   Imp::serialize(
00402     DS::getDefaultSerializer(), count, buffer, bytes, charBuffer);
00403       }
00404 
00406       static Ordinal fromIndirectBytesToCount(const Ordinal bytes, 
00407                 const char charBuffer[]) {
00408   return Imp::fromIndirectBytesToCount(
00409     DS::getDefaultSerializer(), bytes, charBuffer);
00410       }
00411 
00413       static void deserialize (const Ordinal bytes, 
00414              const char charBuffer[], 
00415              const Ordinal count, 
00416              TayType buffer[]) { 
00417   Imp::deserialize(
00418     DS::getDefaultSerializer(), bytes, charBuffer, count, buffer);
00419       }
00420   
00422       
00423     };
00424 
00426     template <typename Ordinal, typename TayType, typename ValueSerializer>
00427     class SerializerImp {
00428 
00429     private:
00430 
00432       typedef SerializationImp<Ordinal,TayType,ValueSerializer> Imp;
00433 
00435       Teuchos::RCP<const ValueSerializer> vs;
00436 
00438       Ordinal sz;
00439 
00440     public:
00441 
00443       typedef ValueSerializer value_serializer_type;
00444 
00446       static const bool supportsDirectSerialization = 
00447   Imp::supportsDirectSerialization;
00448 
00450       SerializerImp(const Teuchos::RCP<const ValueSerializer>& vs_,
00451         Ordinal sz_ = 0) :
00452   vs(vs_), sz(sz_) {}
00453 
00455       Ordinal getSerializerSize() const { return sz; }
00456       
00458       Teuchos::RCP<const value_serializer_type> getValueSerializer() const { 
00459   return vs; }
00460 
00462 
00463 
00465       Ordinal fromCountToIndirectBytes(const Ordinal count, 
00466                const TayType buffer[]) const { 
00467   return Imp::fromCountToIndirectBytes(*vs, count, buffer, sz);
00468       }
00469 
00471       void serialize (const Ordinal count, 
00472           const TayType buffer[], 
00473           const Ordinal bytes, 
00474           char charBuffer[]) const { 
00475   Imp::serialize(*vs, count, buffer, bytes, charBuffer, sz);
00476       }
00477 
00479       Ordinal fromIndirectBytesToCount(const Ordinal bytes, 
00480                const char charBuffer[]) const {
00481   return Imp::fromIndirectBytesToCount(*vs, bytes, charBuffer, sz);
00482       }
00483 
00485       void deserialize (const Ordinal bytes, 
00486       const char charBuffer[], 
00487       const Ordinal count, 
00488       TayType buffer[]) const { 
00489   return Imp::deserialize(*vs, bytes, charBuffer, count, buffer, sz);
00490       }
00491   
00493       
00494     };
00495 
00496   } // namespace Tay
00497 
00498 } // namespace Sacado
00499 
00500 #endif // HAVE_SACADO_TEUCHOS
00501 
00502 #endif // SACADO_FAD_SCALARTRAITSIMP_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines