From e4bc621c318684a400af9fbfe3d6895f6f94823a Mon Sep 17 00:00:00 2001 From: akva Date: Tue, 6 Mar 2012 12:50:32 +0000 Subject: [PATCH] changed: move recovery related functions in ASMs3D into ASMs3Drecovery.C git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1504 e10b68d5-8a6e-419e-a041-bce267b0401d --- src/ASM/ASMs3D.C | 237 --------------------------------------- src/ASM/ASMs3Drecovery.C | 235 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 235 insertions(+), 237 deletions(-) create mode 100644 src/ASM/ASMs3Drecovery.C diff --git a/src/ASM/ASMs3D.C b/src/ASM/ASMs3D.C index bb9d8d63..f2b45ffa 100644 --- a/src/ASM/ASMs3D.C +++ b/src/ASM/ASMs3D.C @@ -1878,38 +1878,6 @@ bool ASMs3D::getGridParameters (RealArray& prm, int dir, int nSegPerSpan) const } -bool ASMs3D::getGrevilleParameters (RealArray& prm, int dir) const -{ - if (!svol) return false; - - const Go::BsplineBasis& basis = svol->basis(dir); - - prm.resize(basis.numCoefs()); - for (size_t i = 0; i < prm.size(); i++) - prm[i] = basis.grevilleParameter(i); - - return true; -} - - -bool ASMs3D::getQuasiInterplParameters (RealArray& prm, int dir) const -{ - if (!svol) return false; - const Go::BsplineBasis& basis = svol->basis(dir); - - std::vector< double > knots_simple; - basis.knotsSimple(knots_simple); - - prm.clear(); - for (size_t i = 0; i < knots_simple.size(); i++){ - prm.push_back(knots_simple[i]); - prm.push_back(0.5*(knots_simple[i]+knots_simple[i+1]));} - prm.pop_back(); - - return true; -} - - bool ASMs3D::tesselate (ElementBlock& grid, const int* npe) const { // Compute parameter values of the nodal points @@ -2092,211 +2060,6 @@ bool ASMs3D::evalSolution (Matrix& sField, const IntegrandBase& integrand, } -Go::GeomObject* ASMs3D::evalSolution (const IntegrandBase& integrand) const -{ - return this->projectSolution(integrand); -} - - -Go::SplineVolume* ASMs3D::projectSolution (const IntegrandBase& integrand) const -{ -PROFILE1("test Global projection"); - // Compute parameter values of the result sampling points (Greville points) - RealArray gpar[3]; - for (int dir = 0; dir < 3; dir++) - if (!this->getGrevilleParameters(gpar[dir],dir)) - return 0; - - // Evaluate the secondary solution at all sampling points - Matrix sValues; - if (!this->evalSolution(sValues,integrand,gpar) || sValues.rows() == 0) - return 0; - - // Project the results onto the spline basis to find control point - // values based on the result values evaluated at the Greville points. - // Note that we here implicitly assume that the number of Greville points - // equals the number of control points such that we don't have to resize - // the result array. Think that is always the case, but beware if trying - // other projection schemes later. - - RealArray weights; - if (svol->rational()) - svol->getWeights(weights); - - const Vector& vec = sValues; - return Go::VolumeInterpolator::regularInterpolation(svol->basis(0), - svol->basis(1), - svol->basis(2), - gpar[0], gpar[1], gpar[2], - const_cast(vec), - sValues.rows(), - svol->rational(), - weights); -} - - -Go::SplineVolume* ASMs3D::projectSolutionLeastSquare (const IntegrandBase& integrand) const -{ - PROFILE1("test L2- projection"); - // Compute parameter values of the result sampling points (Gauss-Interpl. points) - // Get Gaussian quadrature points and weights - const double* xg = GaussQuadrature::getCoord(nGauss); - const double* wg = GaussQuadrature::getWeight(nGauss); - if (!xg || !wg) return false; - - Matrix ggpar[3]; - for (int dir = 0; dir < 3; dir++) - this->getGaussPointParameters(ggpar[dir],dir,nGauss,xg); - - std::vector par_u; - par_u = ggpar[0]; - std::vector par_v; - par_v = ggpar[1]; - std::vector par_w; - par_w = ggpar[2]; - - // gauss weights at parameter values - std::vector wgpar_u; - std::vector wgpar_v; - std::vector wgpar_w; - double d; - for (int dir = 0; dir < 3; dir++){ - if (!svol) return false; - const Go::BsplineBasis& basis = svol->basis(dir); - RealArray::const_iterator knotit = svol->basis(dir).begin(); - std::vector tmp; - tmp.reserve(nGauss*(basis.numCoefs()-basis.order())); - for (size_t i = 0; i<=(basis.numCoefs()-basis.order());i++) - { - d = knotit[i+basis.order()]-knotit[i+basis.order()-1]; - if (d > 0) - {for (int j = 0;jevalSolution(sValues,integrand,gpar)) - return 0; - - // Project the results onto the spline basis to find control point - // values based on the result values evaluated at the Greville points. - // Note that we here implicitly assume the the number of Greville points - // equals the number of control points such that we don't have to resize - // the result array. Think that is always the case, but beware if trying - // other projection schemes later. - - RealArray weights; - if (svol->rational()) - svol->getWeights(weights); - - const Vector& vec = sValues; - return leastsquare_approximation(svol->basis(0), - svol->basis(1), - svol->basis(2), - par_u, par_v, par_w, - wgpar_u, wgpar_v, wgpar_w, - const_cast(vec), - sValues.rows(), - svol->rational(), - weights); -} - - -Go::SplineVolume* ASMs3D::projectSolutionLocal (const IntegrandBase& integrand) const -{ - PROFILE1("test Quasi projection"); - // Compute parameter values of the result sampling points (Greville points) - RealArray gpar[3]; - for (int dir = 0; dir < 2; dir++) - if (!this->getQuasiInterplParameters(gpar[dir],dir)) - return 0; - - for (int dir = 2; dir < 3; dir++) - if (!this->getGrevilleParameters(gpar[dir],dir)) - return 0; - - // Evaluate the secondary solution at all sampling points - Matrix sValues; - if (!this->evalSolution(sValues,integrand,gpar)) - return 0; - - // Project the results onto the spline basis to find control point - // values based on the result values evaluated at the Greville points. - // Note that we here implicitly assume the the number of Greville points - // equals the number of control points such that we don't have to resize - // the result array. Think that is always the case, but beware if trying - // other projection schemes later. - - RealArray weights; - if (svol->rational()) - svol->getWeights(weights); - - const Vector& vec = sValues; - return quasiInterpolation(svol->basis(0), - svol->basis(1), - svol->basis(2), - gpar[0], gpar[1], gpar[2], - const_cast(vec), - sValues.rows(), - svol->rational(), - weights); -} - - -Go::SplineVolume* ASMs3D::projectSolutionLocalApprox(const IntegrandBase& integrand) const -{ - PROFILE1("test VDSA projection"); - // Compute parameter values of the result sampling points (Greville points) - RealArray gpar[3]; - for (int dir = 0; dir < 3; dir++) - if (!this->getGrevilleParameters(gpar[dir],dir)) - return 0; - - // Evaluate the secondary solution at all sampling points - Matrix sValues; - if (!this->evalSolution(sValues,integrand,gpar)) - return 0; - - // Project the results onto the spline basis to find control point - // values based on the result values evaluated at the Greville points. - // Note that we here implicitly assume the the number of Greville points - // equals the number of control points such that we don't have to resize - // the result array. Think that is always the case, but beware if trying - // other projection schemes later. - - RealArray weights; - if (svol->rational()) - svol->getWeights(weights); - - const Vector& vec = sValues; - return VariationDiminishingSplineApproximation(svol->basis(0), - svol->basis(1), - svol->basis(2), - gpar[0], gpar[1], gpar[2], - const_cast(vec), - sValues.rows(), - svol->rational(), - weights); -} - - bool ASMs3D::evalSolution (Matrix& sField, const IntegrandBase& integrand, const RealArray* gpar, bool regular) const { diff --git a/src/ASM/ASMs3Drecovery.C b/src/ASM/ASMs3Drecovery.C new file mode 100644 index 00000000..0c6e3bc6 --- /dev/null +++ b/src/ASM/ASMs3Drecovery.C @@ -0,0 +1,235 @@ +// $Id$ +//============================================================================== +//! +//! \file ASMs3Drecovery.C +//! +//! \date March 06 2012 +//! +//! \author Knut Morten Okstad / SINTEF +//! +//! \brief Recovery of secondary solutions for structured 3D spline FE models. +//! +//============================================================================== + +#include "GoTools/trivariate/SplineVolume.h" +#include "GoTools/trivariate/VolumeInterpolator.h" + +#include "ASMs3D.h" +#include "GaussQuadrature.h" +#include "Profiler.h" + + +bool ASMs3D::getGrevilleParameters (RealArray& prm, int dir) const +{ + if (!svol) return false; + + const Go::BsplineBasis& basis = svol->basis(dir); + + prm.resize(basis.numCoefs()); + for (size_t i = 0; i < prm.size(); i++) + prm[i] = basis.grevilleParameter(i); + + return true; +} + + +bool ASMs3D::getQuasiInterplParameters (RealArray& prm, int dir) const +{ + if (!svol) return false; + const Go::BsplineBasis& basis = svol->basis(dir); + + std::vector< double > knots_simple; + basis.knotsSimple(knots_simple); + + prm.clear(); + for (size_t i = 0; i < knots_simple.size(); i++){ + prm.push_back(knots_simple[i]); + prm.push_back(0.5*(knots_simple[i]+knots_simple[i+1]));} + prm.pop_back(); + + return true; +} + + +Go::GeomObject* ASMs3D::evalSolution (const IntegrandBase& integrand) const +{ + return this->projectSolution(integrand); +} + + +Go::SplineVolume* ASMs3D::projectSolution (const IntegrandBase& integrand) const +{ + PROFILE1("test Global projection"); + // Compute parameter values of the result sampling points (Greville points) + RealArray gpar[3]; + for (int dir = 0; dir < 3; dir++) + if (!this->getGrevilleParameters(gpar[dir],dir)) + return 0; + + // Evaluate the secondary solution at all sampling points + Matrix sValues; + if (!this->evalSolution(sValues,integrand,gpar) || sValues.rows() == 0) + return 0; + + // Project the results onto the spline basis to find control point + // values based on the result values evaluated at the Greville points. + // Note that we here implicitly assume that the number of Greville points + // equals the number of control points such that we don't have to resize + // the result array. Think that is always the case, but beware if trying + // other projection schemes later. + + RealArray weights; + if (svol->rational()) + svol->getWeights(weights); + + const Vector& vec = sValues; + return Go::VolumeInterpolator::regularInterpolation(svol->basis(0), + svol->basis(1), + svol->basis(2), + gpar[0], gpar[1], gpar[2], + const_cast(vec), + sValues.rows(), + svol->rational(), + weights); +} + + +Go::SplineVolume* ASMs3D::projectSolutionLeastSquare (const IntegrandBase& integrand) const +{ + PROFILE1("test L2- projection"); + // Compute parameter values of the result sampling points (Gauss-Interpl. points) + // Get Gaussian quadrature points and weights + const double* xg = GaussQuadrature::getCoord(nGauss); + const double* wg = GaussQuadrature::getWeight(nGauss); + if (!xg || !wg) return false; + + Matrix ggpar[3]; + for (int dir = 0; dir < 3; dir++) + this->getGaussPointParameters(ggpar[dir],dir,nGauss,xg); + + std::vector par_u; + par_u = ggpar[0]; + std::vector par_v; + par_v = ggpar[1]; + std::vector par_w; + par_w = ggpar[2]; + + // gauss weights at parameter values + std::vector wgpar_u; + std::vector wgpar_v; + std::vector wgpar_w; + double d; + for (int dir = 0; dir < 3; dir++){ + if (!svol) return false; + const Go::BsplineBasis& basis = svol->basis(dir); + RealArray::const_iterator knotit = svol->basis(dir).begin(); + std::vector tmp; + tmp.reserve(nGauss*(basis.numCoefs()-basis.order())); + for (size_t i = 0; i<=(basis.numCoefs()-basis.order());i++) + { + d = knotit[i+basis.order()]-knotit[i+basis.order()-1]; + if (d > 0) + {for (int j = 0;jevalSolution(sValues,integrand,gpar)) + return 0; + + RealArray weights; + if (svol->rational()) + svol->getWeights(weights); + + const Vector& vec = sValues; + return leastsquare_approximation(svol->basis(0), + svol->basis(1), + svol->basis(2), + par_u, par_v, par_w, + wgpar_u, wgpar_v, wgpar_w, + const_cast(vec), + sValues.rows(), + svol->rational(), + weights); +} + + +Go::SplineVolume* ASMs3D::projectSolutionLocal (const IntegrandBase& integrand) const +{ + PROFILE1("test Quasi projection"); + // Compute parameter values of the result sampling points (Greville points) + RealArray gpar[3]; + for (int dir = 0; dir < 2; dir++) + if (!this->getQuasiInterplParameters(gpar[dir],dir)) + return 0; + + for (int dir = 2; dir < 3; dir++) + if (!this->getGrevilleParameters(gpar[dir],dir)) + return 0; + + // Evaluate the secondary solution at all sampling points + Matrix sValues; + if (!this->evalSolution(sValues,integrand,gpar)) + return 0; + + RealArray weights; + if (svol->rational()) + svol->getWeights(weights); + + const Vector& vec = sValues; + return quasiInterpolation(svol->basis(0), + svol->basis(1), + svol->basis(2), + gpar[0], gpar[1], gpar[2], + const_cast(vec), + sValues.rows(), + svol->rational(), + weights); +} + + +Go::SplineVolume* ASMs3D::projectSolutionLocalApprox(const IntegrandBase& integrand) const +{ + PROFILE1("test VDSA projection"); + // Compute parameter values of the result sampling points (Greville points) + RealArray gpar[3]; + for (int dir = 0; dir < 3; dir++) + if (!this->getGrevilleParameters(gpar[dir],dir)) + return 0; + + // Evaluate the secondary solution at all sampling points + Matrix sValues; + if (!this->evalSolution(sValues,integrand,gpar)) + return 0; + + RealArray weights; + if (svol->rational()) + svol->getWeights(weights); + + const Vector& vec = sValues; + return VariationDiminishingSplineApproximation(svol->basis(0), + svol->basis(1), + svol->basis(2), + gpar[0], gpar[1], gpar[2], + const_cast(vec), + sValues.rows(), + svol->rational(), + weights); +}