#include <BoundaryConditionDiscretizationSpectralConform.hpp>


Definition at line 44 of file BoundaryConditionDiscretizationSpectralConform.hpp.
| BoundaryConditionDiscretizationSpectralConform::BoundaryConditionDiscretizationSpectralConform | ( | const Problem & | problem, | |
| const SpectralMesh & | mesh, | |||
| const DegreeOfFreedomSet & | dof, | |||
| const DiscretizationType & | discretizationType | |||
| ) | [inline] |
Definition at line 75 of file BoundaryConditionDiscretizationSpectralConform.hpp.
References __basisU, __dofShape, __intervalU, __transformU, ScalarDiscretizationTypeSpectral::a(), ScalarDiscretizationTypeSpectral::b(), ScalarDiscretizationTypeSpectral::degrees(), ScalarDiscretizationTypeBase::name(), Problem::numberOfUnknown(), BoundaryConditionDiscretization::problem(), ScalarDiscretizationTypeBase::spectralLegendre, and ErrorHandler::unexpected.
00079 : BoundaryConditionDiscretization(problem,dof), 00080 __degreeOfFreedomSet(dof), 00081 __mesh(mesh), 00082 __dofShape(this->problem().numberOfUnknown()), 00083 __interval(new Interval(__mesh.shape().a()[0],__mesh.shape().b()[0]), 00084 new Interval(__mesh.shape().a()[1],__mesh.shape().b()[1]), 00085 new Interval(__mesh.shape().a()[2],__mesh.shape().b()[2])), 00086 00087 __intervalU(this->problem().numberOfUnknown()), 00088 __transform(new SpectralConformTransformation(*__interval[0]), 00089 new SpectralConformTransformation(*__interval[1]), 00090 new SpectralConformTransformation(*__interval[2])), 00091 00092 __transformU(this->problem().numberOfUnknown()), 00093 00094 __gaussLobatto(GaussLobattoManager::instance().getReference(__mesh.degree(0)+1), 00095 GaussLobattoManager::instance().getReference(__mesh.degree(1)+1), 00096 GaussLobattoManager::instance().getReference(__mesh.degree(2)+1)), 00097 __basisU(this->problem().numberOfUnknown()) 00098 { 00099 for (size_t i=0; i < this->problem().numberOfUnknown(); i++) { 00100 if (discretizationType[i].type() != ScalarDiscretizationTypeBase::spectralLegendre) { 00101 throw ErrorHandler(__FILE__,__LINE__, 00102 "discretization '" 00103 +ScalarDiscretizationTypeBase::name(discretizationType[i]) 00104 +"' is incompatible with spectral method", 00105 ErrorHandler::unexpected); 00106 } 00107 const ScalarDiscretizationTypeSpectral& discretization 00108 = dynamic_cast<const ScalarDiscretizationTypeSpectral&>(discretizationType[i]); 00109 00110 __dofShape[i] = discretization.degrees()+TinyVector<3,size_t>(1,1,1); 00111 00112 __intervalU[i][0]= new Interval(discretization.a()[0],discretization.b()[0]); 00113 __intervalU[i][1]= new Interval(discretization.a()[1],discretization.b()[1]); 00114 __intervalU[i][2]= new Interval(discretization.a()[2],discretization.b()[2]); 00115 __transformU[i][0] = new SpectralConformTransformation(*__intervalU[i][0]); 00116 __transformU[i][1] = new SpectralConformTransformation(*__intervalU[i][1]); 00117 __transformU[i][2] = new SpectralConformTransformation(*__intervalU[i][2]); 00118 __basisU[i][0]= new LegendreBasis(discretization.degrees()[0]); 00119 __basisU[i][1]= new LegendreBasis(discretization.degrees()[1]); 00120 __basisU[i][2]= new LegendreBasis(discretization.degrees()[2]); 00121 } 00122 00123 }

| BoundaryConditionDiscretizationSpectralConform::~BoundaryConditionDiscretizationSpectralConform | ( | ) | [inline] |
| void BoundaryConditionDiscretizationSpectralConform::getDiagonal | ( | BaseVector & | X | ) | const [virtual] |
Implements BoundaryConditionDiscretization.
Definition at line 43 of file BoundaryConditionDiscretizationSpectralConform.cpp.
References __basisU, __degreeOfFreedomSet, __dofShape, __gaussLobatto, __mesh, __transform, __transformU, Structured3DMeshShape::a(), VariationalBilinearBorderOperator::alphaUV, Structured3DMeshShape::b(), VariationalProblem::beginBilinearBorderOperator(), VariationalProblem::endBilinearBorderOperator(), Problem::numberOfUnknown(), Problem::pde, BoundaryConditionDiscretization::problem(), BoundaryReferences::references(), Boundary::references, Vector< T >::resize(), SpectralMesh::shape(), VariationalBorderOperator::testFunctionNumber(), VariationalBilinearBorderOperator::type(), Problem::type(), ErrorHandler::unexpected, VariationalBilinearBorderOperator::unknownNumber(), and Problem::variational.
00044 { 00045 Vector<real_t>& v = dynamic_cast<Vector<real_t>&>(X); 00046 00047 switch(this->problem().type()) { 00048 case Problem::variational: { 00049 const VariationalProblem& variationalProblem 00050 = dynamic_cast<const VariationalProblem&>(this->problem()); 00051 00052 for (VariationalProblem::bilinearBorderOperatorConst_iterator 00053 iBorderOperator = variationalProblem.beginBilinearBorderOperator(); 00054 iBorderOperator != variationalProblem.endBilinearBorderOperator(); 00055 ++iBorderOperator) { 00056 00057 const VariationalBilinearBorderOperator& borderOperator = **iBorderOperator; 00058 switch (borderOperator.type()) { 00059 case VariationalBilinearBorderOperator::alphaUV: { 00060 const VariationalBorderOperatorAlphaUV& operatorAlphaUV 00061 = dynamic_cast<const VariationalBorderOperatorAlphaUV&>(borderOperator); 00062 ConstReferenceCounting<Boundary> boundary = operatorAlphaUV.boundary(); 00063 00064 switch(boundary->type()) { 00065 case Boundary::references: { 00066 const BoundaryReferences& boundaryReference = 00067 dynamic_cast<const BoundaryReferences&>(*boundary); 00068 00069 for (BoundaryReferences::ReferencesSet::const_iterator 00070 iBoundary = boundaryReference.references().begin(); 00071 iBoundary != boundaryReference.references().end(); ++iBoundary) { 00072 00073 const ScalarFunctionBase& alpha = operatorAlphaUV.alpha(); 00074 const size_t boundaryNumber = *iBoundary; 00075 00076 const size_t& numberOfUnknown = this->problem().numberOfUnknown(); 00077 00078 Vector<TinyVector<3,Vector<Vector<real_t> > > > quadratureBaseValuesU(numberOfUnknown); 00079 TinyVector<3,size_t> direction; 00080 for (size_t i=0; i<3; ++i) { 00081 direction[i] = (boundaryNumber/2+i+1) % 3; 00082 } 00083 00084 for (size_t unknownNumber =0; unknownNumber < numberOfUnknown;++ unknownNumber){ 00085 00086 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >& unknownBasis 00087 = __basisU[unknownNumber]; 00088 00089 for (size_t i=0; i<2; ++i) { 00090 quadratureBaseValuesU[unknownNumber][i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00091 for (size_t j=0; j < (*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00092 00093 real_t xj=(*__transformU[unknownNumber][direction[i]]).inverse((*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j))); 00094 00095 quadratureBaseValuesU[unknownNumber][i][j].resize(unknownBasis[direction[i]]->dimension()); 00096 00097 unknownBasis[direction[i]]->getValues(xj,quadratureBaseValuesU[unknownNumber][i][j]); 00098 } 00099 } 00100 00101 quadratureBaseValuesU[unknownNumber][2].resize(2); 00102 for (size_t j=0; j<2 ; ++j) { 00103 real_t xj= -1. + 2.*j; 00104 quadratureBaseValuesU[unknownNumber][2][j].resize(unknownBasis[direction[2]]->dimension()); 00105 unknownBasis[direction[2]]->getValues(xj,quadratureBaseValuesU[unknownNumber][2][j]); 00106 } 00107 } 00108 00109 TinyVector<2, Vector<real_t> >nodes; 00110 for (size_t i=0; i<2; ++i) { 00111 nodes[i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00112 for (size_t j=0; j<(*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00113 nodes[i][j] =(*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j)); 00114 } 00115 } 00116 const real_t jacobianDet = __transform[direction[0]]->determinant()* 00117 __transform[direction[1]]->determinant(); 00118 00119 Vector<Vector<real_t> > alpha_rs(__gaussLobatto[direction[0]]->numberOfPoints()); 00120 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00121 alpha_rs[r].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00122 for(size_t s=0; s< __gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00123 alpha_rs[r][s]=0; 00124 } 00125 } 00126 00127 for (size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00128 const real_t& x = nodes[0][r]; 00129 for (size_t s=0; s <__gaussLobatto[direction[1]]->numberOfPoints(); ++s) { 00130 const real_t& y = nodes[1][s]; 00131 00132 switch(boundaryNumber) { 00133 case 0: { 00134 alpha_rs[r][s] = alpha(__mesh.shape().a(0), x, y); 00135 break; 00136 } 00137 case 1: { 00138 alpha_rs[r][s] = alpha(__mesh.shape().b(0), x, y); 00139 break; 00140 } 00141 case 2: { 00142 alpha_rs[r][s] = alpha(y, __mesh.shape().a(1), x); 00143 break; 00144 } 00145 case 3: { 00146 alpha_rs[r][s] = alpha(y, __mesh.shape().b(1), x); 00147 break; 00148 } 00149 case 4: { 00150 alpha_rs[r][s] = alpha(x, y, __mesh.shape().a(2)); 00151 break; 00152 } 00153 case 5: { 00154 alpha_rs[r][s] = alpha(x, y, __mesh.shape().b(2)); 00155 break; 00156 } 00157 } 00158 } 00159 } 00160 00161 00162 const TinyVector<3,size_t>& dofShape = __dofShape[operatorAlphaUV.unknownNumber()]; 00163 00164 TinyVector<3,size_t> I; 00165 00166 for(size_t unknownNumber=0;unknownNumber< numberOfUnknown; unknownNumber++) { 00167 00168 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >& unknownBasis 00169 = __basisU[unknownNumber]; 00170 00171 for(size_t m=0; m< unknownBasis[direction[2]]->dimension(); ++m){ 00172 I[direction[2]]=m; 00173 for(size_t n=0; n< unknownBasis[direction[0]]->dimension(); ++n){ 00174 I[direction[0]]=n; 00175 for(size_t p=0; p<unknownBasis[direction[1]]->dimension(); ++p){ 00176 I[direction[1]]=p; 00177 size_t dofNumber = I[0]*dofShape[1]*dofShape[2] + I[1]*dofShape[2] + I[2]; 00178 00179 for(size_t r=0; r< __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00180 for(size_t s=0; s < __gaussLobatto[direction[1]]->numberOfPoints(); ++s) { 00181 v[__degreeOfFreedomSet(operatorAlphaUV.testFunctionNumber(),dofNumber)] 00182 += quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][2][boundaryNumber%2][m] 00183 *quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][2][boundaryNumber%2][m] 00184 *quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][0][n] 00185 *quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][0][n] 00186 *quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][1][p] 00187 *quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][1][p] 00188 *__gaussLobatto[direction[0]]->weight(r) 00189 *__gaussLobatto[direction[1]]->weight(s) 00190 *alpha_rs[r][s]*jacobianDet; 00191 } 00192 } 00193 } 00194 } 00195 } 00196 } 00197 } 00198 00199 break; 00200 } 00201 default: { 00202 throw ErrorHandler(__FILE__,__LINE__, 00203 "not implemented", 00204 ErrorHandler::unexpected); 00205 } 00206 } 00207 break; 00208 } 00209 default: { 00210 throw ErrorHandler(__FILE__,__LINE__, 00211 "not implemented", 00212 ErrorHandler::unexpected); 00213 } 00214 } 00215 } 00216 break; 00217 } 00218 case Problem::pde: { 00219 break; 00220 } 00221 default: { 00222 throw ErrorHandler(__FILE__,__LINE__, 00223 "not implemented", 00224 ErrorHandler::unexpected); 00225 } 00226 } 00227 }

| void BoundaryConditionDiscretizationSpectralConform::setMatrix | ( | ReferenceCounting< BaseMatrix > | A, | |
| ReferenceCounting< BaseVector > | b | |||
| ) | const [virtual] |
Implements BoundaryConditionDiscretization.
Definition at line 231 of file BoundaryConditionDiscretizationSpectralConform.cpp.
References BaseMatrix::doubleHashedMatrix, UnAssembledMatrix::setBoundaryConditions(), BaseMatrix::unAssembled, and ErrorHandler::unexpected.
00233 { 00234 switch((*givenA).type()) { 00235 case BaseMatrix::doubleHashedMatrix: { 00236 throw ErrorHandler(__FILE__,__LINE__, 00237 "not implemented", 00238 ErrorHandler::unexpected); 00239 break; 00240 } 00241 case BaseMatrix::unAssembled: { 00242 UnAssembledMatrix& A = dynamic_cast<UnAssembledMatrix&>(*givenA); 00243 A.setBoundaryConditions(this); 00244 break; 00245 } 00246 00247 default: { 00248 throw ErrorHandler(__FILE__,__LINE__, 00249 "unexected matrix type", 00250 ErrorHandler::unexpected); 00251 } 00252 } 00253 }

| void BoundaryConditionDiscretizationSpectralConform::setSecondMember | ( | ReferenceCounting< BaseMatrix > | A, | |
| ReferenceCounting< BaseVector > | b | |||
| ) | const [virtual] |
Implements BoundaryConditionDiscretization.
Definition at line 257 of file BoundaryConditionDiscretizationSpectralConform.cpp.
References __basisU, __degreeOfFreedomSet, __dofShape, __gaussLobatto, __mesh, __transform, __transformU, Structured3DMeshShape::a(), Structured3DMeshShape::b(), VariationalProblem::beginLinearBorderOperator(), VariationalProblem::endLinearBorderOperator(), ffout(), VariationalLinearBorderOperator::FV, Problem::numberOfUnknown(), BoundaryConditionDiscretization::problem(), BoundaryReferences::references(), Boundary::references, Vector< T >::resize(), SpectralMesh::shape(), VariationalBorderOperator::testFunctionNumber(), VariationalLinearBorderOperator::type(), Problem::type(), ErrorHandler::unexpected, and Problem::variational.
00259 { 00260 ffout(2)<<"- assembling second membre-boundary\n"; 00261 Vector<real_t>& b= (static_cast<Vector<real_t>&>(*givenb)); 00262 00263 switch(this->problem().type()) { 00264 case Problem::variational: { 00265 const VariationalProblem& variationalProblem 00266 = dynamic_cast<const VariationalProblem&>(this->problem()); 00267 00268 for (VariationalProblem::linearBorderOperatorConst_iterator 00269 iLinearBorderOperator = variationalProblem.beginLinearBorderOperator(); 00270 iLinearBorderOperator != variationalProblem.endLinearBorderOperator(); 00271 ++iLinearBorderOperator) { 00272 00273 const VariationalLinearBorderOperator& linearBorderOperator = **iLinearBorderOperator; 00274 switch (linearBorderOperator.type()) { 00275 case VariationalLinearBorderOperator::FV: { 00276 const VariationalBorderOperatorFV& operatorFV 00277 = dynamic_cast<const VariationalBorderOperatorFV&>(linearBorderOperator); 00278 ConstReferenceCounting<Boundary> boundary = operatorFV.boundary(); 00279 00280 switch(boundary->type()) { 00281 case Boundary::references: { 00282 const BoundaryReferences& boundaryReference = 00283 dynamic_cast<const BoundaryReferences&>(*boundary); 00284 00285 for (BoundaryReferences::ReferencesSet::const_iterator 00286 iBoundary = boundaryReference.references().begin(); 00287 iBoundary != boundaryReference.references().end(); ++iBoundary) { 00288 00289 const ScalarFunctionBase& g = operatorFV.f(); 00290 const size_t boundaryNumber = *iBoundary; 00291 00292 const size_t& numberOfUnknown = this->problem().numberOfUnknown(); 00293 00294 TinyVector<3, size_t> direction; 00295 for (size_t i=0; i<3; ++i) { 00296 direction[i] = (boundaryNumber/2+i+1) % 3; 00297 } 00298 00299 Vector<TinyVector<3,Vector<Vector<real_t> > > > quadratureBaseValues(numberOfUnknown); 00300 for (size_t testFunctionNumber =0; testFunctionNumber < numberOfUnknown; ++ testFunctionNumber){ 00301 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >& testBasis 00302 = __basisU[testFunctionNumber]; 00303 00304 for (size_t i=0; i<2; ++i) { 00305 quadratureBaseValues[testFunctionNumber][i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00306 00307 for (size_t j=0; j < (*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00308 real_t xj= (*__transformU[testFunctionNumber][direction[i]]).inverse((*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j))); 00309 00310 quadratureBaseValues[testFunctionNumber][i][j].resize(testBasis[direction[i]]->dimension()); 00311 testBasis[direction[i]]->getValues(xj,quadratureBaseValues[testFunctionNumber][i][j]); 00312 } 00313 } 00314 00315 quadratureBaseValues[testFunctionNumber][2].resize(2); 00316 for (size_t j=0; j<2 ; ++j) { 00317 real_t xj= -1.+2.*j; 00318 quadratureBaseValues[testFunctionNumber][2][j].resize(testBasis[direction[2]]->dimension()); 00319 testBasis[direction[2]]->getValues(xj,quadratureBaseValues[testFunctionNumber][2][j]); 00320 } 00321 } 00322 00323 // nodes construction 00324 TinyVector<2, Vector<real_t> >nodes; 00325 for (size_t i=0; i<2; ++i) { 00326 nodes[i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00327 for (size_t j=0; j< (*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00328 nodes[i][j] =(*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j)); 00329 } 00330 } 00331 00332 const real_t jacobianDet = (*__transform[direction[0]]).determinant()* 00333 (*__transform[direction[1]]).determinant(); 00334 00335 for (size_t i=0; i < (*__gaussLobatto[direction[0]]).numberOfPoints(); ++i) { 00336 const real_t& x = nodes[0][i]; 00337 const real_t& wi = (*__gaussLobatto[direction[0]]).weight(i) 00338 ; 00339 00340 for (size_t j=0; j < (*__gaussLobatto[direction[1]]).numberOfPoints(); ++j) { 00341 const real_t& y = nodes[1][j]; 00342 const real_t& wij = (*__gaussLobatto[direction[1]]).weight(j)*wi 00343 ; 00344 00345 real_t gValue = wij* jacobianDet; 00346 switch(boundaryNumber) { 00347 case 0: { 00348 gValue *= g(__mesh.shape().a(0), x, y); 00349 break; 00350 } 00351 case 1: { 00352 gValue *= g(__mesh.shape().b(0), x, y); 00353 break; 00354 } 00355 case 2: { 00356 gValue *= g(y, __mesh.shape().a(1), x); 00357 break; 00358 } 00359 case 3: { 00360 gValue *= g(y, __mesh.shape().b(1), x); 00361 break; 00362 } 00363 case 4: { 00364 gValue *= g(x, y, __mesh.shape().a(2)); 00365 break; 00366 } 00367 case 5: { 00368 gValue *= g(x, y, __mesh.shape().b(2)); 00369 break; 00370 } 00371 default: { 00372 throw ErrorHandler(__FILE__,__LINE__, 00373 "not implemented", 00374 ErrorHandler::unexpected); 00375 } 00376 } 00377 00378 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >& testBasis 00379 = __basisU[operatorFV.testFunctionNumber()]; 00380 00381 const TinyVector<3,size_t>& dofShape = 00382 __dofShape[operatorFV.testFunctionNumber()]; 00383 00384 TinyVector<3,size_t> J; 00385 for (size_t ml=0; ml < testBasis[direction[0]]->dimension(); ++ml){ 00386 J[direction[0]] = ml; 00387 for (size_t pl=0; pl < testBasis[direction[1]]->dimension(); ++pl){ 00388 J[direction[1]] = pl; 00389 for (size_t ql=0; ql < testBasis[direction[2]]->dimension(); ++ql){ 00390 J[direction[2]] = ql; 00391 const size_t l = 00392 J[0]*dofShape[1]*dofShape[2] + J[1]*dofShape[2] + J[2]; 00393 b[__degreeOfFreedomSet(operatorFV.testFunctionNumber(),l)] 00394 00395 +=quadratureBaseValues[operatorFV.testFunctionNumber()][0][i][ml] 00396 * quadratureBaseValues[operatorFV.testFunctionNumber()][1][j][pl] 00397 * quadratureBaseValues[operatorFV.testFunctionNumber()][2][boundaryNumber%2][ql] 00398 *gValue 00399 ; 00400 } 00401 } 00402 } 00403 } 00404 } 00405 } 00406 00407 break; 00408 } 00409 default: { 00410 throw ErrorHandler(__FILE__,__LINE__, 00411 "not implemented", 00412 ErrorHandler::unexpected); 00413 } 00414 00415 } 00416 break; 00417 } 00418 default: { 00419 throw ErrorHandler(__FILE__,__LINE__, 00420 "not implemented", 00421 ErrorHandler::unexpected); 00422 } 00423 } 00424 } 00425 00426 break; 00427 } 00428 00429 default: { 00430 throw ErrorHandler(__FILE__,__LINE__, 00431 "not implemented", 00432 ErrorHandler::unexpected); 00433 } 00434 } 00435 }

| void BoundaryConditionDiscretizationSpectralConform::timesX | ( | const BaseVector & | X, | |
| BaseVector & | Z | |||
| ) | const [virtual] |
Implements BoundaryConditionDiscretization.
Definition at line 439 of file BoundaryConditionDiscretizationSpectralConform.cpp.
References __basisU, __degreeOfFreedomSet, __dofShape, __gaussLobatto, __mesh, __transform, __transformU, Structured3DMeshShape::a(), VariationalBilinearBorderOperator::alphaUV, Structured3DMeshShape::b(), VariationalProblem::beginBilinearBorderOperator(), VariationalProblem::endBilinearBorderOperator(), Problem::numberOfUnknown(), Problem::pde, BoundaryConditionDiscretization::problem(), BoundaryReferences::references(), Boundary::references, Vector< T >::resize(), SpectralMesh::shape(), VariationalBorderOperator::testFunctionNumber(), VariationalBilinearBorderOperator::type(), Problem::type(), ErrorHandler::unexpected, VariationalBilinearBorderOperator::unknownNumber(), and Problem::variational.
00440 { 00441 const Vector<real_t>& u = dynamic_cast<const Vector<real_t>&>(X); 00442 Vector<real_t>& v = dynamic_cast<Vector<real_t>&>(Z); 00443 00444 switch(this->problem().type()) { 00445 case Problem::variational: { 00446 const VariationalProblem& variationalProblem 00447 = dynamic_cast<const VariationalProblem&>(this->problem()); 00448 00449 for (VariationalProblem::bilinearBorderOperatorConst_iterator 00450 iBorderOperator = variationalProblem.beginBilinearBorderOperator(); 00451 iBorderOperator != variationalProblem.endBilinearBorderOperator(); 00452 ++iBorderOperator) { 00453 00454 const VariationalBilinearBorderOperator& borderOperator = **iBorderOperator; 00455 switch (borderOperator.type()) { 00456 case VariationalBilinearBorderOperator::alphaUV: { 00457 00458 const VariationalBorderOperatorAlphaUV& operatorAlphaUV 00459 = dynamic_cast<const VariationalBorderOperatorAlphaUV&>(borderOperator); 00460 ConstReferenceCounting<Boundary> boundary = operatorAlphaUV.boundary(); 00461 00462 switch(boundary->type()) { 00463 case Boundary::references: { 00464 const BoundaryReferences& boundaryReference = dynamic_cast<const BoundaryReferences&>(*boundary); 00465 for (BoundaryReferences::ReferencesSet::const_iterator 00466 iBoundary = boundaryReference.references().begin(); 00467 iBoundary != boundaryReference.references().end(); ++iBoundary) { 00468 00469 const ScalarFunctionBase& alpha = operatorAlphaUV.alpha(); 00470 const size_t boundaryNumber = *iBoundary; 00471 00472 const size_t& numberOfUnknown = this->problem().numberOfUnknown(); 00473 00474 Vector<TinyVector<3,Vector<Vector<real_t> > > > quadratureBaseValuesU(numberOfUnknown); 00475 TinyVector<3,size_t> direction; 00476 for (size_t i=0; i<3; ++i) { 00477 direction[i] = (boundaryNumber/2+i+1) % 3; 00478 } 00479 00480 for (size_t unknownNumber =0; unknownNumber < numberOfUnknown;++ unknownNumber){ 00481 00482 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >& unknownBasis 00483 = __basisU[unknownNumber]; 00484 00485 for (size_t i=0; i<2; ++i) { 00486 quadratureBaseValuesU[unknownNumber][i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00487 for (size_t j=0; j < (*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00488 00489 real_t xj=(*__transformU[unknownNumber][direction[i]]).inverse((*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j))); 00490 00491 quadratureBaseValuesU[unknownNumber][i][j].resize(unknownBasis[direction[i]]->dimension()); 00492 00493 unknownBasis[direction[i]]->getValues(xj,quadratureBaseValuesU[unknownNumber][i][j]); 00494 } 00495 } 00496 00497 quadratureBaseValuesU[unknownNumber][2].resize(2); 00498 for (size_t j=0; j<2 ; ++j) { 00499 real_t xj= -1. + 2.*j; 00500 quadratureBaseValuesU[unknownNumber][2][j].resize(unknownBasis[direction[2]]->dimension()); 00501 unknownBasis[direction[2]]->getValues(xj,quadratureBaseValuesU[unknownNumber][2][j]); 00502 } 00503 } 00504 00505 // nodes construction 00506 00507 TinyVector<2, Vector<real_t> >nodes; 00508 for (size_t i=0; i<2; ++i) { 00509 nodes[i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00510 for (size_t j=0; j<(*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00511 nodes[i][j] =(*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j)); 00512 } 00513 } 00514 const real_t jacobianDet = __transform[direction[0]]->determinant()* 00515 __transform[direction[1]]->determinant(); 00516 00517 const size_t unknownNumber= operatorAlphaUV.unknownNumber(); 00518 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >&unknownBasis 00519 = __basisU[unknownNumber]; 00520 00521 Vector<Vector<real_t> > u_jk(unknownBasis[direction[0]]->dimension()); 00522 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00523 u_jk[j].resize(unknownBasis[direction[1]]->dimension()); 00524 for(size_t k=0; k< unknownBasis[direction[1]]->dimension(); ++k){ 00525 u_jk[j][k]=0; 00526 } 00527 } 00528 00529 { 00530 const TinyVector<3,size_t>& dofShapeU = __dofShape[unknownNumber]; 00531 TinyVector<3,size_t> J; 00532 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00533 J[direction[0]]=j; 00534 for(size_t k=0; k< unknownBasis[direction[1]]->dimension(); ++k){ 00535 J[direction[1]]=k; 00536 real_t& u_jk_value = u_jk[j][k]; 00537 for(size_t i=0; i< unknownBasis[direction[2]]->dimension(); ++i){ 00538 J[direction[2]]=i; 00539 const size_t dofNumberU = 00540 J[0]*dofShapeU[1]*dofShapeU[2] + J[1] *dofShapeU[2] +J[2]; 00541 u_jk_value 00542 += quadratureBaseValuesU[unknownNumber][2][boundaryNumber%2][i] 00543 *u[__degreeOfFreedomSet(operatorAlphaUV.unknownNumber(),dofNumberU)]; 00544 } 00545 } 00546 } 00547 } 00548 00549 Vector<Vector<real_t> > u_js(unknownBasis[direction[0]]->dimension()); 00550 for(size_t j=0; j<unknownBasis[direction[0]]->dimension(); ++j){ 00551 u_js[j].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00552 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00553 u_js[j][s]=0; 00554 } 00555 } 00556 00557 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00558 for(size_t s=0; s< __gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00559 real_t& u_js_value = u_js[j][s]; 00560 for(size_t k=0; k< unknownBasis[direction[1]]->dimension(); ++k){ 00561 u_js_value 00562 += quadratureBaseValuesU[unknownNumber][1][s][k]* u_jk[j][k]; 00563 } 00564 } 00565 } 00566 00567 Vector<Vector<real_t> > u_rs(__gaussLobatto[direction[0]]->numberOfPoints()); 00568 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00569 u_rs[r].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00570 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00571 u_rs[r][s]=0; 00572 } 00573 } 00574 00575 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00576 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00577 real_t& u_rs_value = u_rs[r][s]; 00578 for(size_t j=0; j<unknownBasis[direction[0]]->dimension(); ++j){ 00579 u_rs_value 00580 += quadratureBaseValuesU[unknownNumber][0][r][j]*u_js[j][s]; 00581 } 00582 } 00583 } 00584 00585 00586 Vector<Vector<real_t> > alpha_rs(__gaussLobatto[direction[0]]->numberOfPoints()); 00587 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00588 alpha_rs[r].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00589 for(size_t s=0; s< __gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00590 alpha_rs[r][s]=0; 00591 } 00592 } 00593 00594 for (size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00595 const real_t& x = nodes[0][r]; 00596 for (size_t s=0; s <__gaussLobatto[direction[1]]->numberOfPoints(); ++s) { 00597 const real_t& y = nodes[1][s]; 00598 00599 00600 switch(boundaryNumber) { 00601 case 0: { 00602 alpha_rs[r][s] = alpha(__mesh.shape().a(0), x, y); 00603 break; 00604 } 00605 case 1: { 00606 alpha_rs[r][s] = alpha(__mesh.shape().b(0), x, y); 00607 break; 00608 } 00609 case 2: { 00610 alpha_rs[r][s] = alpha(y, __mesh.shape().a(1), x); 00611 break; 00612 } 00613 case 3: { 00614 alpha_rs[r][s] = alpha(y, __mesh.shape().b(1), x); 00615 break; 00616 } 00617 case 4: { 00618 alpha_rs[r][s] = alpha(x, y, __mesh.shape().a(2)); 00619 break; 00620 } 00621 case 5: { 00622 alpha_rs[r][s] = alpha(x, y, __mesh.shape().b(2)); 00623 break; 00624 } 00625 } 00626 } 00627 } 00628 00629 const size_t testFunctionNumber= operatorAlphaUV.testFunctionNumber(); 00630 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >&testBasis 00631 = __basisU[testFunctionNumber]; 00632 Vector<Vector<real_t> > v_rp(__gaussLobatto[direction[0]]->numberOfPoints()); 00633 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00634 v_rp[r].resize(testBasis[direction[1]]->dimension()); 00635 for(size_t p=0; p<testBasis[direction[1]]->dimension(); ++p){ 00636 v_rp[r][p]=0; 00637 } 00638 } 00639 00640 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00641 for(size_t p=0; p<testBasis[direction[1]]->dimension(); ++p){ 00642 real_t& v_rp_value = v_rp[r][p]; 00643 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00644 v_rp_value 00645 += quadratureBaseValuesU[testFunctionNumber][1][s][p] 00646 *u_rs[r][s] 00647 *__gaussLobatto[direction[1]]->weight(s) 00648 * alpha_rs[r][s]*jacobianDet; 00649 } 00650 } 00651 } 00652 00653 Vector<Vector<real_t> > v_np(testBasis[direction[0]]->dimension()); 00654 for(size_t n=0; n<testBasis[direction[0]]->dimension(); ++n){ 00655 v_np[n].resize(testBasis[direction[1]]->dimension()); 00656 for(size_t p=0; p<testBasis[direction[1]]->dimension(); ++p){ 00657 v_np[n][p]=0; 00658 } 00659 } 00660 00661 for(size_t n=0; n<testBasis[direction[0]]->dimension(); ++n){ 00662 for(size_t p=0; p<testBasis[direction[1]]->dimension(); ++p){ 00663 real_t& v_np_value = v_np[n][p]; 00664 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00665 v_np_value 00666 += quadratureBaseValuesU[testFunctionNumber][0][r][n] 00667 *v_rp[r][p]* 00668 __gaussLobatto[direction[0]]->weight(r); 00669 } 00670 } 00671 } 00672 00673 const TinyVector<3,size_t>& dofShapeV = __dofShape[testFunctionNumber]; 00674 TinyVector<3,size_t> I; 00675 00676 for(size_t n=0; n<testBasis[direction[0]]->dimension(); ++n){ 00677 I[direction[0]]=n; 00678 for(size_t p=0; p<testBasis[direction[1]]->dimension(); ++p){ 00679 I[direction[1]]=p; 00680 for(size_t m=0; m<testBasis[direction[2]]->dimension(); ++m){ 00681 I[direction[2]]=m; 00682 size_t dofNumberV = 00683 I[0]*dofShapeV[1]*dofShapeV[2] + I[1]*dofShapeV[2] +I[2]; 00684 00685 v[__degreeOfFreedomSet(operatorAlphaUV.testFunctionNumber(),dofNumberV)] 00686 += quadratureBaseValuesU[testFunctionNumber][2][boundaryNumber%2][m] 00687 *v_np[n][p]; 00688 } 00689 } 00690 } 00691 } 00692 break; 00693 } 00694 default: { 00695 throw ErrorHandler(__FILE__,__LINE__, 00696 "not implemented", 00697 ErrorHandler::unexpected); 00698 } 00699 } 00700 break; 00701 } 00702 default: { 00703 throw ErrorHandler(__FILE__,__LINE__, 00704 "not implemented", 00705 ErrorHandler::unexpected); 00706 } 00707 } 00708 } 00709 break; 00710 } 00711 case Problem::pde: { 00712 break; 00713 } 00714 default: { 00715 throw ErrorHandler(__FILE__,__LINE__, 00716 "not implemented", 00717 ErrorHandler::unexpected); 00718 } 00719 } 00720 }

| void BoundaryConditionDiscretizationSpectralConform::transposedTimesX | ( | const BaseVector & | X, | |
| BaseVector & | Z | |||
| ) | const [virtual] |
Implements BoundaryConditionDiscretization.
Definition at line 723 of file BoundaryConditionDiscretizationSpectralConform.cpp.
References __basisU, __degreeOfFreedomSet, __dofShape, __gaussLobatto, __mesh, __transform, __transformU, Structured3DMeshShape::a(), VariationalBilinearBorderOperator::alphaUV, Structured3DMeshShape::b(), VariationalProblem::beginBilinearBorderOperator(), VariationalProblem::endBilinearBorderOperator(), Problem::numberOfUnknown(), Problem::pde, BoundaryConditionDiscretization::problem(), BoundaryReferences::references(), Boundary::references, Vector< T >::resize(), SpectralMesh::shape(), VariationalBorderOperator::testFunctionNumber(), VariationalBilinearBorderOperator::type(), Problem::type(), ErrorHandler::unexpected, VariationalBilinearBorderOperator::unknownNumber(), and Problem::variational.
00724 { 00725 const Vector<real_t>& u = dynamic_cast<const Vector<real_t>&>(X); 00726 Vector<real_t>& v = dynamic_cast<Vector<real_t>&>(Z); 00727 00728 switch(this->problem().type()) { 00729 case Problem::variational: { 00730 const VariationalProblem& variationalProblem 00731 = dynamic_cast<const VariationalProblem&>(this->problem()); 00732 00733 for (VariationalProblem::bilinearBorderOperatorConst_iterator 00734 iBorderOperator = variationalProblem.beginBilinearBorderOperator(); 00735 iBorderOperator != variationalProblem.endBilinearBorderOperator(); 00736 ++iBorderOperator) { 00737 00738 const VariationalBilinearBorderOperator& borderOperator = **iBorderOperator; 00739 switch (borderOperator.type()) { 00740 case VariationalBilinearBorderOperator::alphaUV: { 00741 00742 const VariationalBorderOperatorAlphaUV& operatorAlphaUV 00743 = dynamic_cast<const VariationalBorderOperatorAlphaUV&>(borderOperator); 00744 ConstReferenceCounting<Boundary> boundary = operatorAlphaUV.boundary(); 00745 00746 switch(boundary->type()) { 00747 case Boundary::references: { 00748 const BoundaryReferences& boundaryReference = dynamic_cast<const BoundaryReferences&>(*boundary); 00749 for (BoundaryReferences::ReferencesSet::const_iterator 00750 iBoundary = boundaryReference.references().begin(); 00751 iBoundary != boundaryReference.references().end(); ++iBoundary) { 00752 00753 const ScalarFunctionBase& alpha = operatorAlphaUV.alpha(); 00754 const size_t boundaryNumber = *iBoundary; 00755 00756 const size_t& numberOfUnknown = this->problem().numberOfUnknown(); 00757 00758 Vector<Vector<Vector<real_t> > > u_jk(numberOfUnknown); 00759 Vector<Vector<Vector<real_t> > > u_js(numberOfUnknown); 00760 Vector<Vector<Vector<real_t> > > u_rs(numberOfUnknown); 00761 Vector<Vector<Vector<real_t> > > v_rp(numberOfUnknown); 00762 Vector<Vector<Vector<real_t> > > v_np(numberOfUnknown); 00763 00764 TinyVector<3, size_t> direction; 00765 for (size_t i=0; i<3; ++i) { 00766 direction[i] = (boundaryNumber/2+i+1) % 3; 00767 } 00768 00769 Vector<TinyVector<3,Vector<Vector<real_t> > > > quadratureBaseValuesU(numberOfUnknown); 00770 for (size_t unknownNumber =0; unknownNumber < numberOfUnknown; ++unknownNumber){ 00771 00772 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >&unknownBasis 00773 = __basisU[unknownNumber]; 00774 00775 for (size_t i=0; i<2; ++i) { 00776 quadratureBaseValuesU[unknownNumber][i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00777 00778 for (size_t j=0; j < (*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00779 00780 real_t xj=(*__transformU[unknownNumber][direction[i]]).inverse((*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j))); 00781 00782 quadratureBaseValuesU[unknownNumber][i][j].resize(unknownBasis[direction[i]]->dimension()); 00783 unknownBasis[direction[i]]->getValues(xj,quadratureBaseValuesU[unknownNumber][i][j]); 00784 } 00785 } 00786 00787 quadratureBaseValuesU[unknownNumber][2].resize(2); 00788 for (size_t j=0; j<2 ; ++j) { 00789 real_t xj= -1. + 2.*j; 00790 quadratureBaseValuesU[unknownNumber][2][j].resize(unknownBasis[direction[2]]->dimension()); 00791 unknownBasis[direction[2]]->getValues(xj,quadratureBaseValuesU[unknownNumber][2][j]); 00792 } 00793 } 00794 00795 //construction des nodes 00796 00797 TinyVector<2, Vector<real_t> >nodes; 00798 for (size_t i=0; i<2; ++i) { 00799 nodes[i].resize((*__gaussLobatto[direction[i]]).numberOfPoints()); 00800 for (size_t j=0; j<(*__gaussLobatto[direction[i]]).numberOfPoints(); ++j) { 00801 nodes[i][j] =(*__transform[direction[i]])((*__gaussLobatto[direction[i]])(j)); 00802 } 00803 } 00804 const real_t jacobianDet = __transform[direction[0]]->determinant()* 00805 __transform[direction[1]]->determinant(); 00806 00807 //const size_t unknownNumber= operatorAlphaUV.unknownNumber(); 00808 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >&unknownBasis 00809 = __basisU[operatorAlphaUV.unknownNumber()]; 00810 00811 u_jk[operatorAlphaUV.unknownNumber()].resize(unknownBasis[direction[0]]->dimension()); 00812 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00813 u_jk[operatorAlphaUV.unknownNumber()][j].resize( unknownBasis[direction[1]]->dimension()); 00814 for(size_t k=0; k < unknownBasis[direction[1]]->dimension(); ++k){ 00815 u_jk[operatorAlphaUV.unknownNumber()][j][k]=0; 00816 } 00817 } 00818 00819 { 00820 const TinyVector<3,size_t>& dofShape = __dofShape[operatorAlphaUV.unknownNumber()]; 00821 TinyVector<3,size_t> J; 00822 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00823 J[direction[0]]=j; 00824 for(size_t k=0; k < unknownBasis[direction[1]]->dimension(); ++k){ 00825 J[direction[1]]=k; 00826 real_t& u_jk_value = u_jk[operatorAlphaUV.unknownNumber()][j][k]; 00827 for(size_t i=0; i < unknownBasis[direction[2]]->dimension(); ++i){ 00828 J[direction[2]]=i; 00829 const size_t dofNumberU = 00830 J[0] *dofShape[1]*dofShape[2] + J[1] *dofShape[2] +J[2] ; 00831 u_jk_value 00832 += quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][2][boundaryNumber%2][i] 00833 *u[__degreeOfFreedomSet(operatorAlphaUV.testFunctionNumber(),dofNumberU)]; 00834 } 00835 } 00836 } 00837 } 00838 00839 u_js[operatorAlphaUV.unknownNumber()].resize(unknownBasis[direction[0]]->dimension()); 00840 for(size_t j=0; j< unknownBasis[direction[0]]->dimension(); ++j){ 00841 u_js[operatorAlphaUV.unknownNumber()][j].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00842 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00843 u_js[operatorAlphaUV.unknownNumber()][j][s]=0; 00844 } 00845 } 00846 00847 for(size_t j=0; j < unknownBasis[direction[0]]->dimension(); ++j){ 00848 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00849 real_t& u_js_value = u_js[operatorAlphaUV.unknownNumber()][j][s]; 00850 for(size_t k=0; k < unknownBasis[direction[1]]->dimension(); ++k){ 00851 u_js_value 00852 += quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][1][s][k] 00853 * u_jk[operatorAlphaUV.unknownNumber()][j][k]; 00854 } 00855 } 00856 } 00857 u_rs[operatorAlphaUV.unknownNumber()].resize(__gaussLobatto[direction[0]]->numberOfPoints()); 00858 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00859 u_rs[operatorAlphaUV.unknownNumber()][r].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00860 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00861 u_rs[operatorAlphaUV.unknownNumber()][r][s]=0; 00862 } 00863 } 00864 00865 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00866 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00867 real_t& u_rs_value = u_rs[operatorAlphaUV.unknownNumber()][r][s]; 00868 for(size_t j=0; j < unknownBasis[direction[0]]->dimension(); ++j){ 00869 u_rs_value 00870 += quadratureBaseValuesU[operatorAlphaUV.unknownNumber()][0][r][j] 00871 *u_js[operatorAlphaUV.unknownNumber()][j][s]; 00872 } 00873 } 00874 } 00875 00876 Vector<Vector<real_t> > alpha_rs(__gaussLobatto[direction[0]]->numberOfPoints()); 00877 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00878 alpha_rs[r].resize(__gaussLobatto[direction[1]]->numberOfPoints()); 00879 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00880 alpha_rs[r][s]=0; 00881 } 00882 } 00883 00884 for (size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r) { 00885 const real_t& x = nodes[0][r]; 00886 for (size_t s=0; s <__gaussLobatto[direction[1]]->numberOfPoints(); ++s) { 00887 const real_t& y = nodes[1][s]; 00888 00889 switch(boundaryNumber) { 00890 case 0: { 00891 alpha_rs[r][s] = alpha(__mesh.shape().a(0), x, y); 00892 break; 00893 } 00894 case 1: { 00895 alpha_rs[r][s] = alpha(__mesh.shape().b(0), x, y); 00896 break; 00897 } 00898 case 2: { 00899 alpha_rs[r][s] = alpha(y, __mesh.shape().a(1), x); 00900 break; 00901 } 00902 case 3: { 00903 alpha_rs[r][s] = alpha(y, __mesh.shape().b(1), x); 00904 break; 00905 } 00906 case 4: { 00907 alpha_rs[r][s] = alpha(x, y, __mesh.shape().a(2)); 00908 break; 00909 } 00910 case 5: { 00911 alpha_rs[r][s] = alpha(x, y, __mesh.shape().b(2)); 00912 break; 00913 } 00914 } 00915 } 00916 } 00917 00918 const TinyVector<3, ConstReferenceCounting<LegendreBasis> >&testBasis 00919 = __basisU[operatorAlphaUV.testFunctionNumber()]; 00920 00921 v_rp[operatorAlphaUV.testFunctionNumber()].resize(__gaussLobatto[direction[0]]->numberOfPoints()); 00922 00923 for(size_t r=0; r<__gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00924 v_rp[operatorAlphaUV.testFunctionNumber()][r].resize(testBasis[direction[1]]->dimension()); 00925 for(size_t p=0; p < testBasis[direction[1]]->dimension(); ++p){ 00926 v_rp[operatorAlphaUV.testFunctionNumber()][r][p]=0; 00927 } 00928 } 00929 00930 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00931 for(size_t p=0; p < testBasis[direction[1]]->dimension(); ++p){ 00932 real_t& v_rp_value = v_rp[operatorAlphaUV.testFunctionNumber()][r][p]; 00933 for(size_t s=0; s<__gaussLobatto[direction[1]]->numberOfPoints(); ++s){ 00934 v_rp_value 00935 += quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][1][s][p] 00936 *u_rs[operatorAlphaUV.unknownNumber()][r][s] 00937 *__gaussLobatto[direction[1]]->weight(s) 00938 * alpha_rs[r][s]*jacobianDet; 00939 } 00940 } 00941 } 00942 00943 v_np[operatorAlphaUV.testFunctionNumber()].resize(testBasis[direction[0]]->dimension()); 00944 for(size_t n=0; n < testBasis[direction[0]]->dimension(); ++n){ 00945 v_np[operatorAlphaUV.testFunctionNumber()][n].resize(testBasis[direction[1]]->dimension()); 00946 for(size_t p=0; p < testBasis[direction[1]]->dimension(); ++p){ 00947 v_np[operatorAlphaUV.testFunctionNumber()][n][p]=0; 00948 } 00949 } 00950 00951 for(size_t n=0; n < testBasis[direction[0]]->dimension(); ++n){ 00952 for(size_t p=0; p < testBasis[direction[1]]->dimension(); ++p){ 00953 real_t& v_np_value = v_np[operatorAlphaUV.testFunctionNumber()][n][p]; 00954 for(size_t r=0; r < __gaussLobatto[direction[0]]->numberOfPoints(); ++r){ 00955 v_np_value 00956 += quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][0][r][n] 00957 *v_rp[operatorAlphaUV.testFunctionNumber()][r][p]* 00958 __gaussLobatto[direction[0]]->weight(r); 00959 } 00960 } 00961 } 00962 00963 const TinyVector<3,size_t>& dofShapeV = __dofShape[operatorAlphaUV.testFunctionNumber()]; 00964 TinyVector<3,size_t> I; 00965 00966 for(size_t n=0; n < testBasis[direction[0]]->dimension(); ++n){ 00967 I[direction[0]]=n; 00968 for(size_t p=0; p < testBasis[direction[1]]->dimension(); ++p){ 00969 I[direction[1]]=p; 00970 for(size_t m=0; m < testBasis[direction[2]]->dimension(); ++m){ 00971 I[direction[2]]=m; 00972 const size_t dofNumberV = 00973 I[0]* dofShapeV[1]*dofShapeV[2] + I[1]*dofShapeV[2] + I[2]; 00974 00975 v[__degreeOfFreedomSet(operatorAlphaUV.unknownNumber(),dofNumberV)] 00976 += quadratureBaseValuesU[operatorAlphaUV.testFunctionNumber()][2][boundaryNumber%2][m] 00977 *v_np[operatorAlphaUV.testFunctionNumber()][n][p]; 00978 } 00979 } 00980 } 00981 } 00982 break; 00983 } 00984 default: { 00985 throw ErrorHandler(__FILE__,__LINE__, 00986 "not implemented", 00987 ErrorHandler::unexpected); 00988 } 00989 } 00990 break; 00991 } 00992 default: { 00993 throw ErrorHandler(__FILE__,__LINE__, 00994 "not implemented", 00995 ErrorHandler::unexpected); 00996 } 00997 } 00998 } 00999 break; 01000 } 01001 case Problem::pde: { 01002 break; 01003 } 01004 default: { 01005 throw ErrorHandler(__FILE__,__LINE__, 01006 "not implemented", 01007 ErrorHandler::unexpected); 01008 } 01009 } 01010 }

| const bool& BoundaryConditionDiscretization::dirichlet | ( | const size_t | i | ) | const [inline, inherited] |
Definition at line 54 of file BoundaryConditionDiscretization.hpp.
References BoundaryConditionDiscretization::__dirichletList.
00055 { 00056 return __dirichletList[i]; 00057 }
| const real_t& BoundaryConditionDiscretization::dirichletValue | ( | const size_t | i | ) | const [inline, inherited] |
Definition at line 59 of file BoundaryConditionDiscretization.hpp.
References BoundaryConditionDiscretization::__dirichletValues.
00060 { 00061 return __dirichletValues[i]; 00062 }
| const Vector<bool>& BoundaryConditionDiscretization::getDirichletList | ( | ) | const [inline, inherited] |
Definition at line 64 of file BoundaryConditionDiscretization.hpp.
References BoundaryConditionDiscretization::__dirichletList.
00065 { 00066 return __dirichletList; 00067 }
| const Problem& BoundaryConditionDiscretization::problem | ( | ) | const [inline, inherited] |
Definition at line 81 of file BoundaryConditionDiscretization.hpp.
References BoundaryConditionDiscretization::__problem.
Referenced by BoundaryConditionDiscretizationSpectralConform(), BoundaryConditionDiscretizationSpectralNonConform::BoundaryConditionDiscretizationSpectralNonConform(), getDiagonal(), BoundaryConditionDiscretizationPenalty< MeshType, TypeOfDiscretization >::getDiagonal(), BoundaryConditionDiscretizationPenalty< MeshType, TypeOfDiscretization >::setMatrix(), setSecondMember(), BoundaryConditionDiscretizationPenalty< MeshType, TypeOfDiscretization >::setSecondMember(), timesX(), BoundaryConditionDiscretizationPenalty< MeshType, TypeOfDiscretization >::timesX(), transposedTimesX(), and BoundaryConditionDiscretizationPenalty< MeshType, TypeOfDiscretization >::transposedTimesX().
00082 { 00083 return __problem; 00084 }
const DegreeOfFreedomSet& BoundaryConditionDiscretizationSpectralConform::__degreeOfFreedomSet [private] |
Reimplemented from BoundaryConditionDiscretization.
Definition at line 48 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
const SpectralMesh& BoundaryConditionDiscretizationSpectralConform::__mesh [private] |
Definition at line 49 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
Vector<TinyVector<3,size_t> > BoundaryConditionDiscretizationSpectralConform::__dofShape [private] |
Definition at line 51 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by BoundaryConditionDiscretizationSpectralConform(), getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
TinyVector<3, ConstReferenceCounting<Interval> > BoundaryConditionDiscretizationSpectralConform::__interval [private] |
Definition at line 53 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Vector<TinyVector<3, ConstReferenceCounting<Interval> > > BoundaryConditionDiscretizationSpectralConform::__intervalU [private] |
Definition at line 54 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by BoundaryConditionDiscretizationSpectralConform().
TinyVector<3, ConstReferenceCounting<SpectralConformTransformation> > BoundaryConditionDiscretizationSpectralConform::__transform [private] |
Definition at line 55 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
Vector<TinyVector<3, ConstReferenceCounting<SpectralConformTransformation> > > BoundaryConditionDiscretizationSpectralConform::__transformU [private] |
Definition at line 56 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by BoundaryConditionDiscretizationSpectralConform(), getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
TinyVector<3, ConstReferenceCounting<GaussLobatto> > BoundaryConditionDiscretizationSpectralConform::__gaussLobatto [private] |
Definition at line 57 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
Vector<TinyVector<3, ConstReferenceCounting<LegendreBasis> > > BoundaryConditionDiscretizationSpectralConform::__basisU [private] |
Definition at line 58 of file BoundaryConditionDiscretizationSpectralConform.hpp.
Referenced by BoundaryConditionDiscretizationSpectralConform(), getDiagonal(), setSecondMember(), timesX(), and transposedTimesX().
const Problem& BoundaryConditionDiscretization::__problem [protected, inherited] |
Definition at line 44 of file BoundaryConditionDiscretization.hpp.
Referenced by BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::__associatesDefinedMeshToBoundaryConditions(), BoundaryConditionFDMDiscretization< MeshType, TypeOfDiscretization >::associatesMeshesToBoundaryConditions(), BoundaryConditionDiscretizationSpectralNonConform::associatesMeshesToBoundaryConditions(), BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::associatesMeshesToBoundaryConditions(), and BoundaryConditionDiscretization::problem().
Vector<bool> BoundaryConditionDiscretization::__dirichletList [mutable, protected, inherited] |
true for vertices supporting Dirichlet.
Definition at line 49 of file BoundaryConditionDiscretization.hpp.
Referenced by BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::__setSecondMemberDirichlet(), BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::__variationalBoundaryConditionsAlphaUVTimesX(), BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::__variationalBoundaryConditionsAlphaUVTransposedTimesX(), BoundaryConditionDiscretization::BoundaryConditionDiscretization(), BoundaryConditionDiscretization::dirichlet(), and BoundaryConditionDiscretization::getDirichletList().
Vector<real_t> BoundaryConditionDiscretization::__dirichletValues [mutable, protected, inherited] |
Definition at line 51 of file BoundaryConditionDiscretization.hpp.
Referenced by BoundaryConditionCommonFEMDiscretization< MeshType, TypeOfDiscretization >::__setSecondMemberDirichlet(), and BoundaryConditionDiscretization::dirichletValue().
1.5.6