Sundance 2.1

00001 /* @HEADER@ */
00002 // ************************************************************************
00003 // 
00004 //                              Sundance
00005 //                 Copyright (2005) Sandia Corporation
00006 // 
00007 // Copyright (year first published) Sandia Corporation.  Under the terms 
00008 // of Contract DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government 
00009 // retains certain rights in this software.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //                                                                                 
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA                                                                                
00025 // Questions? Contact Kevin Long (krlong@sandia.gov), 
00026 // Sandia National Laboratories, Livermore, California, USA
00027 // 
00028 // ************************************************************************
00029 /* @HEADER@ */
00030 
00031 #include "SundancePointwiseUserDefFunctor.hpp"
00032 #include "PlayaTabs.hpp"
00033 #include "SundanceOut.hpp"
00034 
00035 using namespace Sundance;
00036 using namespace Sundance;
00037 
00038 using namespace Sundance;
00039 using namespace Teuchos;
00040 
00041 
00042 PointwiseUserDefFunctor0::PointwiseUserDefFunctor0(const std::string& name, 
00043                                                    int domainDim, 
00044                                                    int rangeDim)
00045   :  UserDefFunctor(name, domainDim, rangeDim)
00046 {}
00047 
00048 void PointwiseUserDefFunctor0::evaluationCallback(int nPoints, int maxDiffOrder,
00049                                     const double** in,
00050                                     double** out) const 
00051 {
00052   TEUCHOS_TEST_FOR_EXCEPTION(maxDiffOrder > 0, std::runtime_error,
00053                      "diff order = " << maxDiffOrder 
00054                      << " not supported for functor "
00055                      << name());
00056 
00057   static Array<double> x;
00058   x.resize(domainDim());
00059   
00060   static Array<double> f;
00061   f.resize(rangeDim());
00062 
00063   for (int i=0; i<nPoints; i++)
00064     {
00065       for (int j=0; j<domainDim(); j++) x[j] = in[j][i];
00066       const double* xp = &(x[0]);
00067       double* fp = &(f[0]);
00068       eval0(xp, fp);
00069       for (int j=0; j<rangeDim(); j++) out[j][i] = fp[j];
00070     }
00071 }
00072 
00073 PointwiseUserDefFunctor1::PointwiseUserDefFunctor1(const std::string& name, 
00074                                                    int domainDim, 
00075                                                    int rangeDim)
00076   : PointwiseUserDefFunctor0(name, domainDim, rangeDim)
00077 {}
00078 
00079 void PointwiseUserDefFunctor1::evaluationCallback(int nPoints, int maxDiffOrder,
00080                                     const double** in,
00081                                     double** out) const 
00082 {
00083   TEUCHOS_TEST_FOR_EXCEPTION(maxDiffOrder > 1, std::runtime_error,
00084                      "diff order = " << maxDiffOrder 
00085                      << " not supported for functor "
00086                      << name());
00087 
00088   static Array<double> x;
00089   x.resize(domainDim());
00090   
00091   static Array<double> f;
00092   if (maxDiffOrder==1) f.resize(rangeDim() * (1 + domainDim()) );
00093   else f.resize(rangeDim());
00094   double* fp = &(f[0]);
00095 
00096   for (int i=0; i<nPoints; i++)
00097     {
00098       for (int j=0; j<domainDim(); j++) x[j] = in[j][i];
00099       const double* xp = &(x[0]);
00100 
00101       if (maxDiffOrder==1) 
00102         {
00103           double* dfp = &(f[rangeDim()]);
00104           eval1(xp, fp, dfp);
00105         }
00106       else eval0(xp, fp);
00107       for (int j=0; j<f.size(); j++) out[j][i] = fp[j];
00108     }
00109 }
00110 
00111 
00112 void PointwiseUserDefFunctor1::eval0(const double* in, double* out) const 
00113 {
00114   static Array<double> dummy;
00115   dummy.resize(domainDim() * rangeDim());
00116 
00117   eval1(in, out, &(dummy[0]));
00118 }
00119 
00120 
00121 
00122 PointwiseUserDefFunctor2::PointwiseUserDefFunctor2(const std::string& name, 
00123                                                    int domainDim, 
00124                                                    int rangeDim)
00125   : PointwiseUserDefFunctor1(name, domainDim, rangeDim)
00126 {}
00127 
00128 void PointwiseUserDefFunctor2::evaluationCallback(int nPoints, int maxDiffOrder,
00129                                                   const double** in,
00130                                                   double** out) const 
00131 {
00132   TEUCHOS_TEST_FOR_EXCEPTION(maxDiffOrder > 2 || maxDiffOrder < 0, std::runtime_error,
00133                      "diff order = " << maxDiffOrder 
00134                      << " not supported for functor "
00135                      << name());
00136 
00137   int nTotal = 1;
00138   int numFirst = domainDim();
00139   int numSecond = domainDim()*(domainDim()+1)/2;
00140 
00141   static Array<double> x;
00142   x.resize(domainDim());
00143   
00144   static Array<double> f;
00145 
00146   if (maxDiffOrder > 0) nTotal += numFirst;
00147   if (maxDiffOrder > 1) nTotal += numSecond;
00148 
00149   f.resize(rangeDim() * nTotal);
00150 
00151   double* fp = &(f[0]);
00152 
00153   for (int i=0; i<nPoints; i++)
00154     {
00155       for (int j=0; j<domainDim(); j++) x[j] = in[j][i];
00156       const double* xp = &(x[0]);
00157 
00158       if (maxDiffOrder==0)
00159         {
00160           eval0(xp, fp);
00161         }
00162       else if (maxDiffOrder==1)
00163         {
00164           double* dfp = &(f[rangeDim()]);
00165           eval1(xp, fp, dfp);
00166         }
00167       else if (maxDiffOrder==2)
00168         {
00169           double* dfp = &(f[rangeDim()]);
00170           double* d2fp = &(f[rangeDim()*(1 + domainDim())]);
00171           eval2(xp, fp, dfp, d2fp);
00172         }
00173       else
00174         {
00175           TEUCHOS_TEST_FOR_EXCEPT(true);
00176         }
00177       for (int j=0; j<f.size(); j++) out[j][i] = fp[j];
00178     }
00179 }
00180 
00181 
00182 void PointwiseUserDefFunctor2::eval0(const double* in, double* f) const 
00183 {
00184   static Array<double> dummy1;
00185   static Array<double> dummy2;
00186   dummy1.resize(rangeDim() *  domainDim());
00187   dummy2.resize(rangeDim() *  domainDim()*(domainDim()+1)/2);
00188 
00189   eval2(in, f, &(dummy1[0]), &(dummy2[0]));
00190 }
00191 
00192 
00193 void PointwiseUserDefFunctor2::eval1(const double* in, double* f, double* df) const 
00194 {
00195   static Array<double> dummy2;
00196   dummy2.resize(rangeDim() *  domainDim()*(domainDim()+1)/2);
00197 
00198   eval2(in, f, df, &(dummy2[0]));
00199 }
SundancePointwiseUserDefFunctor.cpp
