Sacado Package Browser (Single Doxygen Collection) Version of the Day
Fad_CommTests.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //                           Sacado Package
00005 //                 Copyright (2006) Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
00025 // (etphipp@sandia.gov).
00026 // 
00027 // ***********************************************************************
00028 // @HEADER
00029 #include "Teuchos_TestingHelpers.hpp"
00030 #include "Teuchos_CommHelpers.hpp"
00031 #include "Teuchos_DefaultComm.hpp"
00032 #include "Teuchos_Array.hpp"
00033 #include "Teuchos_Comm.hpp"
00034 
00035 #include "Sacado_mpl_apply.hpp"
00036 #include "Sacado_Random.hpp"
00037 
00038 using Teuchos::RCP;
00039 using Teuchos::rcp;
00040 using Teuchos::ValueTypeSerializer;
00041 
00042 template <typename FadType>
00043 bool checkFadArrays(const Teuchos::Array<FadType>& x, 
00044         const Teuchos::Array<FadType>& x2, 
00045         const std::string& tag,
00046         Teuchos::FancyOStream& out) {
00047 
00048   // Check sizes match
00049   bool success = (x.size() == x2.size());
00050   out << tag << " Fad array size test";
00051   if (success)
00052     out << " passed";
00053   else
00054     out << " failed";
00055   out << ":  \n\tExpected:  " << x.size() << ", \n\tGot:       " << x2.size() 
00056       << "." << std::endl;
00057   
00058   // Check Fads match
00059   for (int i=0; i<x.size(); i++) {
00060     bool success2 = Sacado::IsEqual<FadType>::eval(x[i], x2[i]);
00061     out << tag << " Fad array comparison test " << i;
00062     if (success2)
00063       out << " passed";
00064     else
00065   out << " failed";
00066     out << ":  \n\tExpected:  " << x[i] << ", \n\tGot:       " << x2[i] << "." 
00067   << std::endl;
00068     success = success && success2;
00069   }
00070 
00071   return success;
00072 }
00073 
00074 template<typename Ordinal>
00075 bool checkResultOnAllProcs(
00076   const Teuchos::Comm<Ordinal> &comm,
00077   Teuchos::FancyOStream &out,
00078   const bool result
00079   )
00080 {
00081   out << "\nChecking that the above test passed in all processes ...";
00082   int thisResult = ( result ? 1 : 0 );
00083   int sumResult = -1;
00084   Teuchos::reduceAll(comm,Teuchos::REDUCE_SUM,Ordinal(1),&thisResult,
00085          &sumResult);
00086   const bool passed = sumResult==Teuchos::size(comm);
00087   if(passed)
00088     out << " passed\n";
00089   else
00090     out << " (sumResult="<<sumResult<<"!=numProcs="<<Teuchos::size(comm)<<") failed\n";
00091   return passed;
00092 }
00093 
00094 #define FAD_COMM_TESTS(FadType, FAD)          \
00095 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_Broadcast ) {      \
00096   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00097     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00098                   \
00099   int n = 7;                \
00100   int p = 5;                \
00101   ValueTypeSerializer<int,FadType> fts(         \
00102     rcp(new ValueTypeSerializer<int,double>), p);     \
00103                   \
00104   Teuchos::Array<FadType> x(n), x2(n), x3(n);       \
00105   for (int i=0; i<n; i++) {           \
00106     x[i] = FadType(p, rnd.number());          \
00107     for (int j=0; j<p; j++)           \
00108       x[i].fastAccessDx(j) = rnd.number();        \
00109   }                 \
00110   for (int i=0; i<n; i++) {           \
00111     x2[i] = FadType(p, 0.0);            \
00112   }                 \
00113   if (comm->getRank() == 0) {           \
00114     x2 = x;               \
00115     x3 = x;                                                             \
00116   }                 \
00117                   \
00118   Teuchos::broadcast(*comm, 0, n, &x2[0]);        \
00119   bool success1 = checkFadArrays(         \
00120     x, x2, std::string(#FAD)+" Broadcast", out);      \
00121   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00122                                                                         \
00123   Teuchos::broadcast(*comm, fts, 0, n, &x3[0]);       \
00124   bool success2 = checkFadArrays(         \
00125     x, x3, std::string(#FAD)+" Broadcast FTS", out);      \
00126   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00127                   \
00128   success = success1 && success2;                                       \
00129 }                 \
00130                   \
00131 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_GatherAll ) {      \
00132   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00133     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00134                   \
00135   int n = 7;                \
00136   int p = 5;                \
00137   int size = comm->getSize();           \
00138   int rank = comm->getRank();           \
00139   int N = n*size;             \
00140   ValueTypeSerializer<int,FadType> fts(         \
00141     rcp(new ValueTypeSerializer<int,double>), p);     \
00142                   \
00143    Teuchos::Array<FadType> x(n), x2(N), x3(N), x4(N);     \
00144   for (int i=0; i<n; i++) {           \
00145     x[i] = FadType(p, (rank+1)*(i+1));          \
00146     for (int j=0; j<p; j++)           \
00147       x[i].fastAccessDx(j) = (rank+1)*(i+1)*(j+1);      \
00148   }                 \
00149   for (int i=0; i<N; i++) {           \
00150     x2[i] = FadType(p, 0.0);            \
00151   }                 \
00152   for (int j=0; j<size; j++) {            \
00153     for (int i=0; i<n; i++) {           \
00154       x3[n*j+i] = FadType(p, (j+1)*(i+1));        \
00155       for (int k=0; k<p; k++)           \
00156   x3[n*j+i].fastAccessDx(k) = (j+1)*(i+1)*(k+1);      \
00157     }                 \
00158   }                 \
00159                   \
00160   Teuchos::gatherAll(*comm, n, &x[0], N, &x2[0]);     \
00161   bool success1 = checkFadArrays(         \
00162     x3, x2, std::string(#FAD)+" Gather All", out);      \
00163   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00164                   \
00165   Teuchos::gatherAll(*comm, fts, n, &x[0], N, &x4[0]);      \
00166   bool success2 = checkFadArrays(         \
00167     x3, x4, std::string(#FAD)+" Gather All FTS", out);      \
00168   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00169                   \
00170   success = success1 && success2;                                       \
00171 }                 \
00172                   \
00173 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_SumAll ) {       \
00174   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00175     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00176                   \
00177   int n = 7;                \
00178   int p = 5;                \
00179   int num_proc = comm->getSize();         \
00180   ValueTypeSerializer<int,FadType> fts(         \
00181     rcp(new ValueTypeSerializer<int,double>), p);     \
00182                   \
00183   Teuchos::Array<FadType> x(n), sums(n), sums2(n), sums3(n);    \
00184   for (int i=0; i<n; i++) {           \
00185     x[i] = FadType(p, 1.0*(i+1));         \
00186     for (int j=0; j<p; j++)           \
00187       x[i].fastAccessDx(j) = 2.0*(i+1);         \
00188   }                 \
00189   for (int i=0; i<n; i++) {           \
00190     sums[i] = FadType(p, 1.0*(i+1)*num_proc);       \
00191     for (int j=0; j<p; j++)           \
00192       sums[i].fastAccessDx(j) = 2.0*(i+1)*num_proc;     \
00193   }                 \
00194   for (int i=0; i<n; i++) {           \
00195     sums2[i] = FadType(p, 0.0);           \
00196   }                 \
00197                   \
00198   Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, n, &x[0], &sums2[0]);  \
00199   bool success1 = checkFadArrays(         \
00200     sums, sums2, std::string(#FAD)+" Sum All", out);      \
00201   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00202                   \
00203   Teuchos::reduceAll(*comm, fts, Teuchos::REDUCE_SUM, n, &x[0], &sums3[0]); \
00204   bool success2 = checkFadArrays(         \
00205     sums, sums3, std::string(#FAD)+" Sum All FTS", out);    \
00206   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00207                   \
00208   success = success1 && success2;                                       \
00209 }                 \
00210                   \
00211 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_MaxAll ) {       \
00212   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00213     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00214                   \
00215   int n = 7;                \
00216   int p = 5;                \
00217   int rank = comm->getRank();           \
00218   int num_proc = comm->getSize();         \
00219   ValueTypeSerializer<int,FadType> fts(         \
00220     rcp(new ValueTypeSerializer<int,double>), p);     \
00221                   \
00222   Teuchos::Array<FadType> x(n), maxs(n), maxs2(n), maxs3(n);    \
00223   for (int i=0; i<n; i++) {           \
00224     x[i] = FadType(p, 1.0*(i+1)*(rank+1));        \
00225     for (int j=0; j<p; j++)           \
00226       x[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);      \
00227   }                 \
00228   for (int i=0; i<n; i++) {           \
00229     maxs[i] = FadType(p, 1.0*(i+1)*num_proc);       \
00230     for (int j=0; j<p; j++)           \
00231       maxs[i].fastAccessDx(j) = 2.0*(i+1)*num_proc;     \
00232   }                 \
00233   for (int i=0; i<n; i++) {           \
00234     maxs2[i] = FadType(p, 0.0);           \
00235   }                 \
00236                   \
00237   Teuchos::reduceAll(*comm, Teuchos::REDUCE_MAX, n, &x[0], &maxs2[0]);  \
00238   bool success1 = checkFadArrays(         \
00239     maxs, maxs2, std::string(#FAD)+" Max All", out);      \
00240   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00241                   \
00242   Teuchos::reduceAll(*comm, fts, Teuchos::REDUCE_MAX, n, &x[0], &maxs3[0]); \
00243   bool success2 = checkFadArrays(         \
00244     maxs, maxs3, std::string(#FAD)+" Max All FTS", out);    \
00245   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00246                   \
00247   success = success1 && success2;                                       \
00248 }                 \
00249                   \
00250 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_MinAll ) {       \
00251   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00252     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00253                   \
00254   int n = 7;                \
00255   int p = 5;                \
00256   int rank = comm->getRank();           \
00257   ValueTypeSerializer<int,FadType> fts(         \
00258     rcp(new ValueTypeSerializer<int,double>), p);     \
00259                   \
00260   Teuchos::Array<FadType> x(n), mins(n), mins2(n), mins3(n);    \
00261   for (int i=0; i<n; i++) {           \
00262     x[i] = FadType(p, 1.0*(i+1)*(rank+1));        \
00263     for (int j=0; j<p; j++)           \
00264       x[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);      \
00265   }                 \
00266   for (int i=0; i<n; i++) {           \
00267     mins[i] = FadType(p, 1.0*(i+1));          \
00268     for (int j=0; j<p; j++)           \
00269       mins[i].fastAccessDx(j) = 2.0*(i+1);        \
00270   }                 \
00271   for (int i=0; i<n; i++) {           \
00272     mins2[i] = FadType(p, 0.0);           \
00273   }                 \
00274                   \
00275   Teuchos::reduceAll(*comm, Teuchos::REDUCE_MIN, n, &x[0], &mins2[0]);  \
00276   bool success1 = checkFadArrays(         \
00277     mins, mins2, std::string(#FAD)+" Min All", out);      \
00278   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00279                   \
00280   Teuchos::reduceAll(*comm, fts, Teuchos::REDUCE_MIN, n, &x[0], &mins3[0]); \
00281   bool success2 = checkFadArrays(         \
00282     mins, mins3, std::string(#FAD)+" Min All FTS", out);    \
00283   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00284                   \
00285   success = success1 && success2;                                       \
00286 }                 \
00287                   \
00288 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_ScanSum ) {        \
00289   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00290     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00291                   \
00292   int n = 7;                \
00293   int p = 5;                \
00294   int rank = comm->getRank();           \
00295   ValueTypeSerializer<int,FadType> fts(         \
00296     rcp(new ValueTypeSerializer<int,double>), p);     \
00297                   \
00298   Teuchos::Array<FadType> x(n), sums(n), sums2(n), sums3(n);    \
00299   for (int i=0; i<n; i++) {           \
00300     x[i] = FadType(p, 1.0*(i+1));         \
00301     for (int j=0; j<p; j++)           \
00302       x[i].fastAccessDx(j) = 2.0*(i+1);         \
00303   }                 \
00304   for (int i=0; i<n; i++) {           \
00305     sums[i] = FadType(p, 1.0*(i+1)*(rank+1));       \
00306     for (int j=0; j<p; j++)           \
00307       sums[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);     \
00308   }                 \
00309   for (int i=0; i<n; i++) {           \
00310     sums2[i] = FadType(p, 0.0);           \
00311   }                 \
00312                   \
00313   Teuchos::scan(*comm, Teuchos::REDUCE_SUM, n, &x[0], &sums2[0]); \
00314   bool success1 = checkFadArrays(         \
00315     sums, sums2, std::string(#FAD)+" Scan Sum", out);     \
00316   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00317                   \
00318   Teuchos::scan(*comm, fts, Teuchos::REDUCE_SUM, n, &x[0], &sums3[0]);  \
00319   bool success2 = checkFadArrays(         \
00320     sums, sums3, std::string(#FAD)+" Scan Sum FTS", out);   \
00321   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00322                     \
00323   success = success1 && success2;         \
00324 }                 \
00325                   \
00326 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_ScanMax ) {        \
00327   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00328     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00329                   \
00330   int n = 7;                \
00331   int p = 5;                \
00332   int rank = comm->getRank();           \
00333   ValueTypeSerializer<int,FadType> fts(         \
00334     rcp(new ValueTypeSerializer<int,double>), p);     \
00335                   \
00336   Teuchos::Array<FadType> x(n), maxs(n), maxs2(n), maxs3(n);    \
00337   for (int i=0; i<n; i++) {           \
00338     x[i] = FadType(p, 1.0*(i+1)*(rank+1));        \
00339     for (int j=0; j<p; j++)           \
00340       x[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);      \
00341   }                 \
00342   for (int i=0; i<n; i++) {           \
00343     maxs[i] = FadType(p, 1.0*(i+1)*(rank+1));       \
00344     for (int j=0; j<p; j++)           \
00345       maxs[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);     \
00346   }                 \
00347   for (int i=0; i<n; i++) {           \
00348     maxs2[i] = FadType(p, 0.0);           \
00349   }                 \
00350                   \
00351   Teuchos::scan(*comm, Teuchos::REDUCE_MAX, n, &x[0], &maxs2[0]); \
00352   bool success1 = checkFadArrays(         \
00353     maxs, maxs2, std::string(#FAD)+" Scan Max", out);     \
00354   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00355                   \
00356   Teuchos::scan(*comm, fts, Teuchos::REDUCE_MAX, n, &x[0], &maxs3[0]);  \
00357   bool success2 = checkFadArrays(         \
00358     maxs, maxs3, std::string(#FAD)+" Scan Max FTS", out);   \
00359   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00360                   \
00361   success = success1 && success2;                                       \
00362 }                 \
00363                   \
00364 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_ScanMin ) {        \
00365   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00366     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00367                   \
00368   int n = 7;                \
00369   int p = 5;                \
00370   int rank = comm->getRank();           \
00371   ValueTypeSerializer<int,FadType> fts(         \
00372     rcp(new ValueTypeSerializer<int,double>), p);     \
00373                   \
00374   Teuchos::Array<FadType> x(n), mins(n), mins2(n), mins3(n);    \
00375   for (int i=0; i<n; i++) {           \
00376     x[i] = FadType(p, 1.0*(i+1)*(rank+1));        \
00377     for (int j=0; j<p; j++)           \
00378       x[i].fastAccessDx(j) = 2.0*(i+1)*(rank+1);      \
00379   }                 \
00380   for (int i=0; i<n; i++) {           \
00381     mins[i] = FadType(p, 1.0*(i+1));          \
00382     for (int j=0; j<p; j++)           \
00383       mins[i].fastAccessDx(j) = 2.0*(i+1);        \
00384   }                 \
00385   for (int i=0; i<n; i++) {           \
00386     mins2[i] = FadType(p, 0.0);           \
00387   }                 \
00388                   \
00389   Teuchos::scan(*comm, Teuchos::REDUCE_MIN, n, &x[0], &mins2[0]); \
00390   bool success1 = checkFadArrays(         \
00391     mins, mins2, std::string(#FAD)+" Scan Min", out);     \
00392   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00393                   \
00394   Teuchos::scan(*comm, fts, Teuchos::REDUCE_MIN, n, &x[0], &mins3[0]);  \
00395   bool success2 = checkFadArrays(         \
00396     mins, mins3, std::string(#FAD)+" Scan Min FTS", out);   \
00397   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00398                   \
00399   success = success1 && success2;                                       \
00400 }                 \
00401                   \
00402 TEUCHOS_UNIT_TEST( FAD##_Comm, Fad_SendReceive ) {      \
00403   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00404     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00405                   \
00406   int num_proc = comm->getSize();         \
00407   if (num_proc > 1) {             \
00408     int rank = comm->getRank();           \
00409     int n = 7;                \
00410     int p = 5;                \
00411     ValueTypeSerializer<int,FadType> fts(       \
00412       rcp(new ValueTypeSerializer<int,double>), p);     \
00413                   \
00414     Teuchos::Array<FadType> x(n), x2(n), x3(n);       \
00415     for (int i=0; i<n; i++) {           \
00416       x[i] = FadType(p, 1.0*(i+1));         \
00417       for (int j=0; j<p; j++)           \
00418   x[i].fastAccessDx(j) = 2.0*(i+1)*(j+1);       \
00419     }                 \
00420     for (int i=0; i<n; i++) {           \
00421       x2[i] = FadType(p, 0.0);            \
00422     }                 \
00423     if (rank != 1) {              \
00424       x2 = x;               \
00425       x3 = x;               \
00426     }                 \
00427                   \
00428     if (rank == 0) Teuchos::send(*comm, n, &x[0], 1);     \
00429     if (rank == 1) Teuchos::receive(*comm, 0, n, &x2[0]);   \
00430     bool success1 = checkFadArrays(         \
00431       x, x2, std::string(#FAD)+" Send/Receive", out);     \
00432     success1 = checkResultOnAllProcs(*comm, out, success1);   \
00433                   \
00434     if (rank == 0) Teuchos::send(*comm, fts, n, &x[0], 1);    \
00435     if (rank == 1) Teuchos::receive(*comm, fts, 0, n, &x3[0]);    \
00436     bool success2 = checkFadArrays(         \
00437       x, x3, std::string(#FAD)+" Send/Receive FTS", out);   \
00438     success2 = checkResultOnAllProcs(*comm, out, success2);   \
00439                   \
00440     success = success1 && success2;         \
00441   }                 \
00442   else                  \
00443     success = true;             \
00444 }                 \
00445                   \
00446 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_Broadcast ) {     \
00447   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00448   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00449     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00450                   \
00451   int n = 7;                \
00452   int p1 = 5;               \
00453   int p2 = 5;               \
00454   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00455     rcp(new ValueTypeSerializer<int,FadType>(       \
00456     rcp(new ValueTypeSerializer<int,double>), p1));   \
00457   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00458                   \
00459   Teuchos::Array<FadFadType> x(n), x2(n), x3(n);      \
00460   for (int i=0; i<n; i++) {           \
00461     FadType f(p1, rnd.number());          \
00462     for (int k=0; k<p1; k++)            \
00463       f.fastAccessDx(k) = rnd.number();         \
00464     x[i] = FadFadType(p2, f);           \
00465     for (int j=0; j<p2; j++) {            \
00466        FadType g(p1, rnd.number());         \
00467        for (int k=0; k<p1; k++)           \
00468    g.fastAccessDx(k) = rnd.number();        \
00469        x[i].fastAccessDx(j) = g;          \
00470     }                 \
00471   }                 \
00472   for (int i=0; i<n; i++) {           \
00473     x2[i] = FadFadType(p2, FadType(p1, 0.0));       \
00474     for (int j=0; j<p2; j++)            \
00475       x2[i].fastAccessDx(j) = FadType(p1, 0.0);       \
00476   }                 \
00477   if (comm->getRank() == 0) {           \
00478     x2 = x;               \
00479     x3 = x;               \
00480   }                 \
00481                   \
00482   Teuchos::broadcast(*comm, 0, n, &x2[0]);        \
00483   bool success1 = checkFadArrays(         \
00484     x, x2, std::string(#FAD)+"<"+#FAD+"> Broadcast", out);    \
00485   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00486                   \
00487   Teuchos::broadcast(*comm, ffts, 0, n, &x3[0]);      \
00488   bool success2 = checkFadArrays(         \
00489     x, x3, std::string(#FAD)+"<"+#FAD+"> Broadcast FTS", out);    \
00490   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00491                   \
00492   success = success1 && success2;                                       \
00493 }                 \
00494                   \
00495 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_GatherAll ) {     \
00496   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00497   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00498     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00499                   \
00500   int n = 7;                \
00501   int p1 = 5;               \
00502   int p2 = 5;               \
00503   int size = comm->getSize();           \
00504   int rank = comm->getRank();           \
00505   int N = n*size;             \
00506   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00507     rcp(new ValueTypeSerializer<int,FadType>(       \
00508     rcp(new ValueTypeSerializer<int,double>), p1));   \
00509   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00510                   \
00511   Teuchos::Array<FadFadType> x(n), x2(N), x3(N), x4(N);     \
00512   for (int i=0; i<n; i++) {           \
00513     FadType f(p1, (rank+1)*(i+1));          \
00514     for (int k=0; k<p1; k++)            \
00515       f.fastAccessDx(k) = (rank+1)*(i+1)*(k+1);       \
00516     x[i] = FadFadType(p2, f);           \
00517     for (int j=0; j<p2; j++) {            \
00518       x[i].fastAccessDx(j) = f;           \
00519     }                 \
00520   }                 \
00521   for (int i=0; i<N; i++) {           \
00522     x2[i] = FadFadType(p2, FadType(p1, 0.0));       \
00523     for (int j=0; j<p2; j++)            \
00524       x2[i].fastAccessDx(j) = FadType(p1, 0.0);       \
00525   }                 \
00526   for (int j=0; j<size; j++) {            \
00527     for (int i=0; i<n; i++) {           \
00528       FadType f(p1, (j+1)*(i+1));         \
00529       for (int k=0; k<p1; k++)            \
00530   f.fastAccessDx(k) = (j+1)*(i+1)*(k+1);        \
00531       x3[n*j+i] = FadFadType(p2, f);          \
00532       for (int k=0; k<p2; k++)            \
00533   x3[n*j+i].fastAccessDx(k) = f;          \
00534     }                 \
00535   }                 \
00536                   \
00537   Teuchos::gatherAll(*comm, n, &x[0], N, &x2[0]);     \
00538   bool success1 = checkFadArrays(         \
00539     x3, x2, std::string(#FAD)+"<"+#FAD+">  Gather All", out);   \
00540   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00541                   \
00542   Teuchos::gatherAll(*comm, ffts, n, &x[0], N, &x4[0]);     \
00543   bool success2 = checkFadArrays(         \
00544     x3, x4, std::string(#FAD)+"<"+#FAD+">  Gather All FTS", out); \
00545   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00546                   \
00547   success = success1 && success2;                                       \
00548 }                 \
00549                   \
00550 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_SumAll ) {      \
00551   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00552   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00553     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00554                   \
00555   int n = 7;                \
00556   int p1 = 5;               \
00557   int p2 = 5;               \
00558   int num_proc = comm->getSize();         \
00559   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00560     rcp(new ValueTypeSerializer<int,FadType>(       \
00561     rcp(new ValueTypeSerializer<int,double>), p1));   \
00562   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00563                   \
00564   Teuchos::Array<FadFadType> x(n), sums(n), sums2(n), sums3(n);   \
00565   for (int i=0; i<n; i++) {           \
00566     FadType f(p1, 1.0*(i+1));           \
00567     for (int k=0; k<p1; k++)            \
00568       f.fastAccessDx(k) = 2.0*(i+1);          \
00569     x[i] = FadFadType(p2, f);           \
00570     for (int j=0; j<p2; j++) {            \
00571       x[i].fastAccessDx(j) = f;           \
00572     }                 \
00573   }                 \
00574   for (int i=0; i<n; i++) {           \
00575     FadType f(p1, 1.0*(i+1)*num_proc);          \
00576     for (int k=0; k<p1; k++)            \
00577       f.fastAccessDx(k) = 2.0*(i+1)*num_proc;       \
00578     sums[i] = FadFadType(p2, f);          \
00579     for (int j=0; j<p2; j++)            \
00580       sums[i].fastAccessDx(j) = f;          \
00581   }                 \
00582   for (int i=0; i<n; i++) {           \
00583     sums2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00584     for (int j=0; j<p2; j++)            \
00585       sums2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00586   }                 \
00587                   \
00588   Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, n, &x[0], &sums2[0]);  \
00589   bool success1 = checkFadArrays(         \
00590     sums, sums2, std::string(#FAD)+"<"+#FAD+"> Sum All", out);    \
00591   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00592                   \
00593   Teuchos::reduceAll(*comm, ffts, Teuchos::REDUCE_SUM, n, &x[0], &sums3[0]); \
00594   bool success2 = checkFadArrays(         \
00595     sums, sums3, std::string(#FAD)+"<"+#FAD+"> Sum All", out);    \
00596   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00597                   \
00598   success = success1 && success2;                                       \
00599 }                 \
00600                   \
00601 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_MaxAll ) {      \
00602   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00603   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00604     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00605                   \
00606   int n = 7;                \
00607   int p1 = 5;               \
00608   int p2 = 5;               \
00609   int rank = comm->getRank();           \
00610   int num_proc = comm->getSize();         \
00611   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00612     rcp(new ValueTypeSerializer<int,FadType>(       \
00613     rcp(new ValueTypeSerializer<int,double>), p1));   \
00614   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00615                   \
00616   Teuchos::Array<FadFadType> x(n), maxs(n), maxs2(n), maxs3(n);   \
00617   for (int i=0; i<n; i++) {           \
00618     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00619     for (int k=0; k<p1; k++)            \
00620       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00621     x[i] = FadFadType(p2, f);           \
00622     for (int j=0; j<p2; j++) {            \
00623       x[i].fastAccessDx(j) = f;           \
00624     }                 \
00625   }                 \
00626   for (int i=0; i<n; i++) {           \
00627     FadType f(p1, 1.0*(i+1)*num_proc);          \
00628     for (int k=0; k<p1; k++)            \
00629       f.fastAccessDx(k) = 2.0*(i+1)*num_proc;       \
00630     maxs[i] = FadFadType(p2, f);          \
00631     for (int j=0; j<p2; j++)            \
00632       maxs[i].fastAccessDx(j) = f;          \
00633   }                 \
00634   for (int i=0; i<n; i++) {           \
00635     maxs2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00636     for (int j=0; j<p2; j++)            \
00637       maxs2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00638   }                 \
00639                   \
00640   Teuchos::reduceAll(*comm, Teuchos::REDUCE_MAX, n, &x[0], &maxs2[0]);  \
00641   bool success1 = checkFadArrays(         \
00642     maxs, maxs2, std::string(#FAD)+"<"+#FAD+"> Max All", out);    \
00643   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00644                   \
00645   Teuchos::reduceAll(*comm, ffts, Teuchos::REDUCE_MAX, n, &x[0], &maxs3[0]); \
00646   bool success2 = checkFadArrays(         \
00647     maxs, maxs3, std::string(#FAD)+"<"+#FAD+"> Max All FTS", out);  \
00648   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00649                   \
00650   success = success1 && success2;                                       \
00651 }                 \
00652                   \
00653 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_MinAll ) {      \
00654   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00655   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00656     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00657                   \
00658   int n = 7;                \
00659   int p1 = 5;               \
00660   int p2 = 5;               \
00661   int rank = comm->getRank();           \
00662   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00663     rcp(new ValueTypeSerializer<int,FadType>(       \
00664     rcp(new ValueTypeSerializer<int,double>), p1));   \
00665   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00666                   \
00667   Teuchos::Array<FadFadType> x(n), mins(n), mins2(n), mins3(n);   \
00668   for (int i=0; i<n; i++) {           \
00669     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00670     for (int k=0; k<p1; k++)            \
00671       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00672     x[i] = FadFadType(p2, f);           \
00673     for (int j=0; j<p2; j++) {            \
00674       x[i].fastAccessDx(j) = f;           \
00675     }                 \
00676   }                 \
00677   for (int i=0; i<n; i++) {           \
00678     FadType f(p1, 1.0*(i+1));           \
00679     for (int k=0; k<p1; k++)            \
00680       f.fastAccessDx(k) = 2.0*(i+1);          \
00681     mins[i] = FadFadType(p2, f);          \
00682     for (int j=0; j<p2; j++)            \
00683       mins[i].fastAccessDx(j) = f;          \
00684   }                 \
00685   for (int i=0; i<n; i++) {           \
00686     mins2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00687     for (int j=0; j<p2; j++)            \
00688       mins2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00689   }                 \
00690                   \
00691   Teuchos::reduceAll(*comm, Teuchos::REDUCE_MIN, n, &x[0], &mins2[0]);  \
00692   bool success1 = checkFadArrays(         \
00693     mins, mins2, std::string(#FAD)+"<"+#FAD+"> Min All", out);    \
00694   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00695                   \
00696   Teuchos::reduceAll(*comm, ffts, Teuchos::REDUCE_MIN, n, &x[0], &mins3[0]); \
00697   bool success2 = checkFadArrays(         \
00698     mins, mins3, std::string(#FAD)+"<"+#FAD+"> Min All FTS", out);  \
00699   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00700                   \
00701   success = success1 && success2;                                       \
00702 }                 \
00703                   \
00704 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_ScanSum ) {     \
00705   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00706   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00707     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00708                   \
00709   int n = 7;                \
00710   int p1 = 5;               \
00711   int p2 = 5;               \
00712   int rank = comm->getRank();           \
00713   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00714     rcp(new ValueTypeSerializer<int,FadType>(       \
00715     rcp(new ValueTypeSerializer<int,double>), p1));   \
00716   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00717                   \
00718   Teuchos::Array<FadFadType> x(n), sums(n), sums2(n), sums3(n);   \
00719   for (int i=0; i<n; i++) {           \
00720     FadType f(p1, 1.0*(i+1));           \
00721     for (int k=0; k<p1; k++)            \
00722       f.fastAccessDx(k) = 2.0*(i+1);          \
00723     x[i] = FadFadType(p2, f);           \
00724     for (int j=0; j<p2; j++) {            \
00725       x[i].fastAccessDx(j) = f;           \
00726     }                 \
00727   }                 \
00728   for (int i=0; i<n; i++) {           \
00729     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00730     for (int k=0; k<p1; k++)            \
00731       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00732     sums[i] = FadFadType(p2, f);          \
00733     for (int j=0; j<p2; j++)            \
00734       sums[i].fastAccessDx(j) = f;          \
00735   }                 \
00736   for (int i=0; i<n; i++) {           \
00737     sums2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00738     for (int j=0; j<p2; j++)            \
00739       sums2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00740   }                 \
00741                   \
00742   Teuchos::scan(*comm, Teuchos::REDUCE_SUM, n, &x[0], &sums2[0]); \
00743   bool success1 = checkFadArrays(         \
00744     sums, sums2, std::string(#FAD)+"<"+#FAD+"> Scan Sum", out);   \
00745   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00746                   \
00747   Teuchos::scan(*comm, ffts, Teuchos::REDUCE_SUM, n, &x[0], &sums3[0]); \
00748   bool success2 = checkFadArrays(         \
00749     sums, sums3, std::string(#FAD)+"<"+#FAD+"> Scan Sum FTS", out); \
00750   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00751                   \
00752   success = success1 && success2;                                       \
00753 }                 \
00754                   \
00755 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_ScanMax ) {     \
00756   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00757   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00758     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00759                   \
00760   int n = 7;                \
00761   int p1 = 5;               \
00762   int p2 = 5;               \
00763   int rank = comm->getRank();           \
00764   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00765     rcp(new ValueTypeSerializer<int,FadType>(       \
00766     rcp(new ValueTypeSerializer<int,double>), p1));   \
00767   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00768                   \
00769   Teuchos::Array<FadFadType> x(n), maxs(n), maxs2(n), maxs3(n);   \
00770   for (int i=0; i<n; i++) {           \
00771     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00772     for (int k=0; k<p1; k++)            \
00773       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00774     x[i] = FadFadType(p2, f);           \
00775     for (int j=0; j<p2; j++) {            \
00776       x[i].fastAccessDx(j) = f;           \
00777     }                 \
00778   }                 \
00779   for (int i=0; i<n; i++) {           \
00780     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00781     for (int k=0; k<p1; k++)            \
00782       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00783     maxs[i] = FadFadType(p2, f);          \
00784     for (int j=0; j<p2; j++)            \
00785       maxs[i].fastAccessDx(j) = f;          \
00786   }                 \
00787   for (int i=0; i<n; i++) {           \
00788     maxs2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00789     for (int j=0; j<p2; j++)            \
00790       maxs2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00791   }                 \
00792                   \
00793   Teuchos::scan(*comm, Teuchos::REDUCE_MAX, n, &x[0], &maxs2[0]); \
00794   bool success1 = checkFadArrays(         \
00795     maxs, maxs2, std::string(#FAD)+"<"+#FAD+"> Scan Max", out);   \
00796   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00797                   \
00798   Teuchos::scan(*comm, ffts, Teuchos::REDUCE_MAX, n, &x[0], &maxs3[0]); \
00799   bool success2 = checkFadArrays(         \
00800     maxs, maxs3, std::string(#FAD)+"<"+#FAD+"> Scan Max FTS", out); \
00801   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00802                   \
00803   success = success1 && success2;                                       \
00804 }                 \
00805                   \
00806 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_ScanMin ) {     \
00807   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00808   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00809     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00810                   \
00811   int n = 7;                \
00812   int p1 = 5;               \
00813   int p2 = 5;               \
00814   int rank = comm->getRank();           \
00815   RCP< ValueTypeSerializer<int,FadType> > fts =       \
00816     rcp(new ValueTypeSerializer<int,FadType>(       \
00817     rcp(new ValueTypeSerializer<int,double>), p1));   \
00818   ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00819                   \
00820   Teuchos::Array<FadFadType> x(n), mins(n), mins2(n), mins3(n);   \
00821   for (int i=0; i<n; i++) {           \
00822     FadType f(p1, 1.0*(i+1)*(rank+1));          \
00823     for (int k=0; k<p1; k++)            \
00824       f.fastAccessDx(k) = 2.0*(i+1)*(rank+1);       \
00825     x[i] = FadFadType(p2, f);           \
00826     for (int j=0; j<p2; j++) {            \
00827       x[i].fastAccessDx(j) = f;           \
00828     }                 \
00829   }                 \
00830   for (int i=0; i<n; i++) {           \
00831     FadType f(p1, 1.0*(i+1));           \
00832     for (int k=0; k<p1; k++)            \
00833       f.fastAccessDx(k) = 2.0*(i+1);          \
00834     mins[i] = FadFadType(p2, f);          \
00835     for (int j=0; j<p2; j++)            \
00836       mins[i].fastAccessDx(j) = f;          \
00837   }                 \
00838   for (int i=0; i<n; i++) {           \
00839     mins2[i] = FadFadType(p2, FadType(p1, 0.0));      \
00840     for (int j=0; j<p2; j++)            \
00841       mins2[i].fastAccessDx(j) = FadType(p1, 0.0);      \
00842   }                 \
00843                   \
00844   Teuchos::scan(*comm, Teuchos::REDUCE_MIN, n, &x[0], &mins2[0]); \
00845   bool success1 = checkFadArrays(         \
00846     mins, mins2, std::string(#FAD)+"<"+#FAD+"> Scan Min", out);   \
00847   success1 = checkResultOnAllProcs(*comm, out, success1);   \
00848                   \
00849   Teuchos::scan(*comm, ffts, Teuchos::REDUCE_MIN, n, &x[0], &mins3[0]); \
00850   bool success2 = checkFadArrays(         \
00851     mins, mins3, std::string(#FAD)+"<"+#FAD+"> Scan Min FTS", out); \
00852   success2 = checkResultOnAllProcs(*comm, out, success2);   \
00853                   \
00854   success = success1 && success2;                                       \
00855 }                 \
00856                   \
00857 TEUCHOS_UNIT_TEST( FAD##_Comm, FadFad_SendReceive ) {     \
00858   typedef Sacado::mpl::apply<FadType,FadType>::type FadFadType;   \
00859   Teuchos::RCP<const Teuchos::Comm<Ordinal> >       \
00860     comm = Teuchos::DefaultComm<Ordinal>::getComm();      \
00861                   \
00862   int num_proc = comm->getSize();         \
00863   if (num_proc > 1) {             \
00864     int rank = comm->getRank();           \
00865     int n = 7;                \
00866     int p1 = 5;               \
00867     int p2 = 5;               \
00868     RCP< ValueTypeSerializer<int,FadType> > fts =     \
00869       rcp(new ValueTypeSerializer<int,FadType>(       \
00870       rcp(new ValueTypeSerializer<int,double>), p1));   \
00871     ValueTypeSerializer<int,FadFadType> ffts(fts, p2);      \
00872                   \
00873     Teuchos::Array<FadFadType> x(n), x2(n), x3(n);      \
00874     for (int i=0; i<n; i++) {           \
00875       FadType f(p1, 1.0*(i+1));           \
00876       for (int k=0; k<p1; k++)            \
00877   f.fastAccessDx(k) = 2.0*(i+1)*(k+1);        \
00878       x[i] = FadFadType(p2, f);           \
00879       for (int j=0; j<p2; j++)            \
00880   x[i].fastAccessDx(j) = f;         \
00881     }                 \
00882     for (int i=0; i<n; i++) {           \
00883       x2[i] = FadFadType(p2, FadType(p1, 0.0));       \
00884       for (int j=0; j<p2; j++)            \
00885   x2[i].fastAccessDx(j) = FadType(p1, 0.0);     \
00886     }                 \
00887     if (rank != 1) {              \
00888       x2 = x;               \
00889       x3 = x;               \
00890     }                 \
00891                   \
00892     if (rank == 0) Teuchos::send(*comm, n, &x[0], 1);     \
00893     if (rank == 1) Teuchos::receive(*comm, 0, n, &x2[0]);   \
00894     bool success1 = checkFadArrays(         \
00895       x, x2, std::string(#FAD)+"<"+#FAD+"> Send/Receive", out);   \
00896     success1 = checkResultOnAllProcs(*comm, out, success1);   \
00897                   \
00898     if (rank == 0) Teuchos::send(*comm, ffts, n, &x[0], 1);   \
00899     if (rank == 1) Teuchos::receive(*comm, ffts, 0, n, &x3[0]);   \
00900     bool success2 = checkFadArrays(         \
00901       x, x3, std::string(#FAD)+"<"+#FAD+"> Send/Receive FTS", out); \
00902     success2 = checkResultOnAllProcs(*comm, out, success2);   \
00903                   \
00904     success = success1 && success2;         \
00905   }                 \
00906   else                  \
00907     success = true;             \
00908 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines