00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #include "SundanceBasisFamily.hpp"
00032 #include "SundanceUnknownFunction.hpp"
00033 #include "SundanceTestFunction.hpp"
00034 #include "SundanceDiscreteFunction.hpp"
00035 #include "SundanceUnknownFunctionData.hpp"
00036 #include "SundanceTestFunctionData.hpp"
00037 #include "SundanceDiscreteFunctionData.hpp"
00038 #include "SundanceLagrange.hpp"
00039 #include "SundanceEdgeLocalizedBasis.hpp"
00040
00041
00042 using namespace Sundance;
00043 using namespace Teuchos;
00044 using Playa::Handle;
00045 using Playa::Handleable;
00046
00047
00048
00049
00050
00051 int BasisFamily::order() const
00052 {
00053 return ptr()->order();
00054 }
00055
00056 int BasisFamily::dim() const
00057 {
00058 return ptr()->dim();
00059 }
00060
00061 bool BasisFamily::operator==(const BasisFamily& other) const
00062 {
00063 return !(*this < other || other < *this);
00064 }
00065
00066 int BasisFamily::size(const Array<BasisFamily>& b)
00067 {
00068 int rtn = 0;
00069 for (int i=0; i<b.size(); i++) rtn += b[i].dim();
00070 return rtn;
00071 }
00072
00073 int BasisFamily::nReferenceDOFsWithFacets(const CellType& maximalCellType,
00074 const CellType& cellType) const
00075 {
00076 return ptr()->nReferenceDOFsWithFacets(maximalCellType, cellType);
00077 }
00078
00079 int BasisFamily::nReferenceDOFsWithoutFacets(const CellType& maximalCellType,
00080 const CellType& cellType) const
00081 {
00082 return ptr()->nReferenceDOFsWithoutFacets(maximalCellType, cellType);
00083 }
00084
00085 BasisFamily BasisFamily::getBasis(const RCP<const CommonFuncDataStub>& funcData)
00086 {
00087
00088 const UnknownFunctionData* u
00089 = dynamic_cast<const UnknownFunctionData*>(funcData.get());
00090 if (u != 0)
00091 {
00092 return u->basis()[0];
00093 }
00094
00095
00096 const TestFunctionData* t
00097 = dynamic_cast<const TestFunctionData*>(funcData.get());
00098 if (t != 0)
00099 {
00100 return t->basis()[0];
00101 }
00102
00103
00104
00105 const DiscreteFunctionData* d
00106 = dynamic_cast<const DiscreteFunctionData*>(funcData.get());
00107 if (d != 0)
00108 {
00109 return d->discreteSpace().basis()[0];
00110 }
00111
00112 TEUCHOS_TEST_FOR_EXCEPTION(u==0 && t==0 && d==0, std::runtime_error,
00113 "BasisFamily::getBasis() argument is not a recognized "
00114 "type of function data");
00115 return BasisFamily();
00116
00117 }
00118
00119
00120
00121 RCP<BasisDOFTopologyBase> BasisFamily::getBasisTopology(const RCP<const CommonFuncDataStub>& funcData)
00122 {
00123 BasisFamily b = getBasis(funcData);
00124 TEUCHOS_TEST_FOR_EXCEPT(b.ptr().get()==0);
00125
00126 return rcp_dynamic_cast<BasisDOFTopologyBase>(b.ptr());
00127 }
00128
00129
00130 void BasisFamily::getConstrainsForHNDoF( const int indexInParent,
00131 const int maxCellDim,
00132 const int maxNrChild,
00133 const int facetDim,
00134 const int facetIndex,
00135 const int nodeIndex,
00136 Array<int>& localDoFs,
00137 Array<int>& parentFacetDim,
00138 Array<int>& parentFacetIndex,
00139 Array<int>& parentFacetNode,
00140 Array<double>& coefs
00141 ) const
00142 {
00143 Tabs tab;
00144 int verb = 4;
00145 SUNDANCE_MSG3( verb ,tab << "BasisFamily::getConstrainsForHNDoF IN: indexInParent:" << indexInParent
00146 << " maxCellDim:" << maxCellDim << " maxNrChild:" << maxNrChild
00147 << " facetDim:" << facetDim << " facetIndex:" << facetIndex << " nodeIndex:" << nodeIndex);
00148 ptr()->getConstrainsForHNDoF( indexInParent, maxCellDim ,
00149 maxNrChild , facetDim, facetIndex, nodeIndex, localDoFs, parentFacetDim ,
00150 parentFacetIndex , parentFacetNode , coefs );
00151 SUNDANCE_MSG3( verb , tab << "BasisFamily::getConstrainsForHNDoF OUT: localDoFs:" << localDoFs
00152 << " coefs:" << coefs);
00153 SUNDANCE_MSG3( verb , tab << "BasisFamily::getConstrainsForHNDoF OUT: parentFacetDim:" << parentFacetDim );
00154 SUNDANCE_MSG3( verb , tab << "BasisFamily::getConstrainsForHNDoF OUT: parentFacetIndex:" << parentFacetIndex );
00155 SUNDANCE_MSG3( verb , tab << "BasisFamily::getConstrainsForHNDoF OUT: parentFacetNode:" << parentFacetNode );
00156 }
00157
00158
00159 void BasisFamily::refEval(
00160 const CellType& cellType,
00161 const Array<Point>& pts,
00162 const SpatialDerivSpecifier& deriv,
00163 Array<Array<Array<double> > >& result,
00164 int verbosity) const
00165 {
00166 using std::setw;
00167 Tabs tab;
00168 SUNDANCE_MSG3(verbosity, tab << "evaluating basis " << *this
00169 << " with spatial derivative " << deriv);
00170 ptr()->refEval(cellType, pts, deriv, result, verbosity);
00171 std::string f = deriv.toString()+ "[phi_n]";
00172
00173 if (verbosity >= 4)
00174 {
00175 Tabs tab1;
00176 for (int q=0; q<pts.size(); q++)
00177 {
00178 Tabs tab2;
00179 Out::os() << tab1 << "evaluation point = " << pts[q] << std::endl;
00180 Out::os() << tab2 << setw(5) << "n";
00181 int nComps = result.size();
00182 if (nComps == 1)
00183 {
00184 Out::os() << setw(20) << f << std::endl;
00185 }
00186 else
00187 {
00188 for (int d=0; d<nComps; d++)
00189 {
00190 std::string fd = f + "[" + Teuchos::toString(d) + "]";
00191 Out::os() << setw(20) << fd;
00192 }
00193 Out::os() << std::endl;
00194 }
00195 for (int n=0; n<result[0][q].size(); n++)
00196 {
00197 Out::os() << tab2 << setw(5) << n;
00198 for (int d=0; d<nComps; d++)
00199 {
00200 Out::os() << setw(20) << result[d][q][n];
00201 }
00202 Out::os() << std::endl;
00203 }
00204 }
00205 }
00206 }
00207
00208
00209 namespace Sundance
00210 {
00211
00212 Array<std::pair<int, int> >
00213 vectorDimStructure(const Array<BasisFamily>& basis)
00214 {
00215 Array<std::pair<int, int> > rtn;
00216 for (int i=0; i<basis.size(); i++)
00217 {
00218 rtn.append(std::pair<int, int>(basis[i].tensorOrder(), basis[i].dim()));
00219 }
00220 return rtn;
00221 }
00222
00223
00224 Array<std::pair<int, int> > vectorDimStructure(const BasisFamily& basis)
00225 {
00226 return vectorDimStructure(tuple(basis));
00227 }
00228
00229 bool basisRestrictableToBoundary(const BasisFamily& b)
00230 {
00231 const Lagrange* lag = dynamic_cast<const Lagrange*>(b.ptr().get());
00232 const EdgeLocalizedBasis* elb = dynamic_cast<const EdgeLocalizedBasis*>(b.ptr().get());
00233 return lag != 0 || elb != 0;
00234 }
00235
00236 }