added: CompatibleOperators
implements common operators for div-compatible discretizations
This commit is contained in:
Arne Morten Kvarving
parent
c9bc665861
commit
17bb95ace9
108
Apps/Common/CompatibleOperators.C
Normal file
108
Apps/Common/CompatibleOperators.C
Normal file
@ -0,0 +1,108 @@
|
|||||||
|
//==============================================================================
|
||||||
|
//!
|
||||||
|
//! \file CompatibleOperators.C
|
||||||
|
//!
|
||||||
|
//! \date Oct 9 2016
|
||||||
|
//!
|
||||||
|
//! \author Arne Morten Kvarving / SINTEF
|
||||||
|
//!
|
||||||
|
//! \brief Various weak, discrete div-compatible operators.
|
||||||
|
//!
|
||||||
|
//==============================================================================
|
||||||
|
|
||||||
|
#include "CompatibleOperators.h"
|
||||||
|
#include "EqualOrderOperators.h"
|
||||||
|
#include "FiniteElement.h"
|
||||||
|
#include "Vec3.h"
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Advection(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
const Vec3& AC, double scale)
|
||||||
|
{
|
||||||
|
size_t nsd = fe.grad(1).cols();
|
||||||
|
for (size_t n = 1; n <= nsd; ++n)
|
||||||
|
for (size_t i = 1; i <= fe.basis(n).size(); ++i)
|
||||||
|
for (size_t j = 1; j <= fe.basis(n).size(); ++j)
|
||||||
|
for (size_t k = 1; k <= nsd; ++k)
|
||||||
|
EM[n](i,j) += scale*AC[k]*fe.grad(n)(j,k)*fe.basis(n)(i)*fe.detJxW;
|
||||||
|
}
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Gradient(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
double scale)
|
||||||
|
{
|
||||||
|
size_t nsd = fe.grad(1).cols();
|
||||||
|
for (size_t n = 1; n <= nsd; ++n)
|
||||||
|
for (size_t i=1; i <= fe.basis(n).size(); ++i)
|
||||||
|
for (size_t j=1; j <= fe.basis(nsd+1).size(); ++j)
|
||||||
|
EM[8+4*(n-1)](i,j) += -scale*fe.grad(n)(i,n)*fe.basis(nsd+1)(j)*fe.detJxW;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Gradient(Vectors& EV, const FiniteElement& fe,
|
||||||
|
double scale)
|
||||||
|
{
|
||||||
|
for (size_t n = 1; n <= fe.grad(n).cols(); ++n)
|
||||||
|
for (size_t i = 1; i <= fe.basis(n).size(); ++i)
|
||||||
|
EV[n](i) += scale*fe.grad(n)(i,n)*fe.detJxW;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Laplacian(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
double scale, bool stress)
|
||||||
|
{
|
||||||
|
size_t nsd = fe.grad(1).cols();
|
||||||
|
for (size_t n = 1; n <= nsd; ++n)
|
||||||
|
EqualOrderOperators::Weak::Laplacian(EM[n], fe, scale, false, n);
|
||||||
|
|
||||||
|
for (size_t m = 1; m <= nsd && stress; m++)
|
||||||
|
for (size_t i = 1; i <= fe.basis(m).size(); ++i)
|
||||||
|
for (size_t n = m; n <= nsd; n++)
|
||||||
|
for (size_t j = 1; j <= fe.basis(n).size(); ++j) {
|
||||||
|
int idx = m == n ? m : (m == 1 ? 5+n-m : 10+n-m);
|
||||||
|
EM[idx](i,j) += scale* fe.grad(n)(j,m) * fe.grad(m)(i,n) * fe.detJxW;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Mass(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe, double scale)
|
||||||
|
{
|
||||||
|
for (size_t k = 1; k <= fe.grad(1).cols(); ++k)
|
||||||
|
EqualOrderOperators::Weak::Mass(EM[k], fe, scale, k);
|
||||||
|
}
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Source(Vectors& EV,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
const Vec3& f, double scale)
|
||||||
|
{
|
||||||
|
for (size_t k = 1; k <= fe.grad(1).cols(); ++k)
|
||||||
|
EqualOrderOperators::Weak::Source(EV[k], fe, scale*f[k-1], 1, k);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Weak::Source(Vectors& EV,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
double scale)
|
||||||
|
{
|
||||||
|
for (size_t k = 1; k <= fe.grad(1).cols(); ++k)
|
||||||
|
EqualOrderOperators::Weak::Source(EV[k], fe, scale, 1, k);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void CompatibleOperators::Residual::Laplacian(Vectors& EV,
|
||||||
|
const FiniteElement& fe,
|
||||||
|
const Tensor& dUdX,
|
||||||
|
double scale, bool stress)
|
||||||
|
{
|
||||||
|
size_t nsd = fe.grad(1).cols();
|
||||||
|
for (size_t k = 1; k <= nsd; ++k)
|
||||||
|
for (size_t i = 1; i <= fe.basis(k).size(); ++i) {
|
||||||
|
double diff = 0.0;
|
||||||
|
for (size_t m = 1; m <= nsd; ++m)
|
||||||
|
diff += fe.grad(k)(i,m)*dUdX(k,m);
|
||||||
|
EV[k](i) += scale*diff*fe.detJxW;
|
||||||
|
}
|
||||||
|
}
|
||||||
102
Apps/Common/CompatibleOperators.h
Normal file
102
Apps/Common/CompatibleOperators.h
Normal file
@ -0,0 +1,102 @@
|
|||||||
|
//==============================================================================
|
||||||
|
//!
|
||||||
|
//! \file CompatibleOperators.h
|
||||||
|
//!
|
||||||
|
//! \date Oct 9 2016
|
||||||
|
//!
|
||||||
|
//! \author Arne Morten Kvarving / SINTEF
|
||||||
|
//!
|
||||||
|
//! \brief Various weak, div-compatible discrete operators.
|
||||||
|
//!
|
||||||
|
//==============================================================================
|
||||||
|
|
||||||
|
#ifndef COMPATIBLE_OPERATORS_H_
|
||||||
|
#define COMPATIBLE_OPERATORS_H_
|
||||||
|
|
||||||
|
class FiniteElement;
|
||||||
|
class Tensor;
|
||||||
|
class Vec3;
|
||||||
|
|
||||||
|
#include "MatVec.h"
|
||||||
|
|
||||||
|
/*! \brief Common operators using div-compatible discretizations.
|
||||||
|
* \details The operators use the block ordering used in the BlockElmMats class.
|
||||||
|
*/
|
||||||
|
|
||||||
|
class CompatibleOperators
|
||||||
|
{
|
||||||
|
public:
|
||||||
|
//! \brief Common weak operators using div-compatible discretizations.
|
||||||
|
class Weak {
|
||||||
|
public:
|
||||||
|
//! \brief Compute an advection term.
|
||||||
|
//! \param[out] EM The element matrices to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] AC Advecting field
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Advection(std::vector<Matrix>& EM, const FiniteElement& fe,
|
||||||
|
const Vec3& AC, double scale=1.0);
|
||||||
|
|
||||||
|
//! \brief Compute a gradient term.
|
||||||
|
//! \param[out] EM The element matrix to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Gradient(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe, double scale=1.0);
|
||||||
|
|
||||||
|
//! \brief Compute a gradient term.
|
||||||
|
//! \param[out] EV The element vector to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Gradient(Vectors& EV, const FiniteElement& fe,
|
||||||
|
double scale=1.0);
|
||||||
|
|
||||||
|
//! \brief Compute a laplacian.
|
||||||
|
//! \param[out] EM The element matrix to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
//! \param[in] stress Whether to add extra stress formulation terms
|
||||||
|
//! \details Only the upper blocks are added with the stress formulation
|
||||||
|
static void Laplacian(std::vector<Matrix>& EM, const FiniteElement& fe,
|
||||||
|
double scale=1.0, bool stress=false);
|
||||||
|
|
||||||
|
//! \brief Compute a mass term.
|
||||||
|
//! \param[out] EM The element matrices to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Mass(std::vector<Matrix>& EM,
|
||||||
|
const FiniteElement& fe, double scale=1.0);
|
||||||
|
|
||||||
|
//! \brief Compute a vector-source term.
|
||||||
|
//! \param[out] EV The element vectors to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] f Vector with contributions
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Source(Vectors& EV, const FiniteElement& fe, double scale=1.0);
|
||||||
|
|
||||||
|
//! \brief Compute a vector-source term.
|
||||||
|
//! \param[out] EV The element vectors to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] f Vector with contributions
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
static void Source(Vectors& EV, const FiniteElement& fe,
|
||||||
|
const Vec3& f, double scale=1.0);
|
||||||
|
};
|
||||||
|
|
||||||
|
//! \brief Common weak residual operators using div-compatible discretizations.
|
||||||
|
class Residual
|
||||||
|
{
|
||||||
|
public:
|
||||||
|
//! \brief Compute a laplacian term in a residual vector.
|
||||||
|
//! \param[out] EV The element vector to add contribution to
|
||||||
|
//! \param[in] fe The finite element to evaluate for
|
||||||
|
//! \param[in] dUdX Current solution gradient
|
||||||
|
//! \param[in] scale Scaling factor for contribution
|
||||||
|
//! \param[in] stress Whether to add extra stress formulation terms
|
||||||
|
static void Laplacian(Vectors& EV, const FiniteElement& fe,
|
||||||
|
const Tensor& dUdX, double scale=1.0,
|
||||||
|
bool stress=false);
|
||||||
|
};
|
||||||
|
};
|
||||||
|
|
||||||
|
#endif
|
||||||
239
Apps/Common/Test/TestCompatibleOperators.C
Normal file
239
Apps/Common/Test/TestCompatibleOperators.C
Normal file
@ -0,0 +1,239 @@
|
|||||||
|
//==============================================================================
|
||||||
|
//!
|
||||||
|
//! \file TestCompatibleOperators.C
|
||||||
|
//!
|
||||||
|
//! \date Oct 17 2016
|
||||||
|
//!
|
||||||
|
//! \author Arne Morten Kvarving / SINTEF
|
||||||
|
//!
|
||||||
|
//! \brief Tests for various discrete div-compatible operators.
|
||||||
|
//!
|
||||||
|
//==============================================================================
|
||||||
|
|
||||||
|
#include "CompatibleOperators.h"
|
||||||
|
#include "FiniteElement.h"
|
||||||
|
#include "gtest/gtest.h"
|
||||||
|
|
||||||
|
typedef std::vector<std::vector<double>> DoubleVec;
|
||||||
|
const auto&& check_matrix_equal = [](const Matrix& A, const DoubleVec& B)
|
||||||
|
{
|
||||||
|
for (size_t i=1;i<=A.rows();++i)
|
||||||
|
for (size_t j=1;j<=A.cols();++j)
|
||||||
|
ASSERT_NEAR(A(i,j), B[i-1][j-1], 1e-13);
|
||||||
|
};
|
||||||
|
|
||||||
|
static MxFiniteElement getFE()
|
||||||
|
{
|
||||||
|
MxFiniteElement fe({6,6,4});
|
||||||
|
for (size_t i = 1; i <= 2; ++i) {
|
||||||
|
fe.grad(i).resize(6,2);
|
||||||
|
fe.basis(i).resize(6);
|
||||||
|
for (size_t j = 1; j <= 6; ++j) {
|
||||||
|
fe.grad(i)(j,1) = 1.0+12*(i-1)+2*(j-1);
|
||||||
|
fe.grad(i)(j,2) = 2.0+12*(i-1)+2*(j-1);
|
||||||
|
fe.basis(i)(j) = j + 6*(i-1);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
fe.grad(3).resize(4,2);
|
||||||
|
fe.basis(3).resize(4);
|
||||||
|
for (size_t j = 1; j <= 4; ++j) {
|
||||||
|
fe.grad(3)(j,1) = 1.0 + 14.0 + 2*(j-1);
|
||||||
|
fe.grad(3)(j,2) = 2.0 + 14.0 + 2*(j-1);
|
||||||
|
fe.basis(3)(j) = 12 + j;
|
||||||
|
}
|
||||||
|
|
||||||
|
return fe;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/*
|
||||||
|
TESTI(TestEqualOrderOperators, Advection)
|
||||||
|
{
|
||||||
|
FiniteElement fe = getFE();
|
||||||
|
|
||||||
|
Matrix EM_scalar(2, 2);
|
||||||
|
Vec3 U;
|
||||||
|
U[0] = 1.0; U[1] = 2.0;
|
||||||
|
// First component only
|
||||||
|
EqualOrderOperators::Weak::Advection(EM_scalar, fe, U);
|
||||||
|
const DoubleVec EM_scalar_ref = {{ 7.0, 10},
|
||||||
|
{14.0, 20.0}};
|
||||||
|
check_matrix_equal(EM_scalar, EM_scalar_ref);
|
||||||
|
|
||||||
|
Matrix EM_vec(2*2, 2*2);
|
||||||
|
U[0] = 3.0; U[1] = 4.0;
|
||||||
|
EqualOrderOperators::Weak::Advection(EM_vec, fe, U, 2.0);
|
||||||
|
const DoubleVec EM_vec_ref = {{30.0, 0.0, 44.0, 0.0},
|
||||||
|
{ 0.0, 30.0, 0.0, 44.0},
|
||||||
|
{60.0, 0.0, 88.0, 0.0},
|
||||||
|
{ 0.0, 60.0, 0.0, 88.0}};
|
||||||
|
check_matrix_equal(EM_vec, EM_vec_ref);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
TESTI(TestEqualOrderOperators, Divergence)
|
||||||
|
{
|
||||||
|
FiniteElement fe = getFE();
|
||||||
|
|
||||||
|
Matrix EM_scalar(2, 2*2);
|
||||||
|
EqualOrderOperators::Weak::Divergence(EM_scalar, fe);
|
||||||
|
const DoubleVec EM_scalar_ref = {{1.0, 3.0, 2.0, 4.0},
|
||||||
|
{2.0, 6.0, 4.0, 8.0}};
|
||||||
|
check_matrix_equal(EM_scalar, EM_scalar_ref);
|
||||||
|
|
||||||
|
Vector EV_scalar(4);
|
||||||
|
Vec3 D;
|
||||||
|
D[0] = 1.0; D[1] = 2.0;
|
||||||
|
EqualOrderOperators::Weak::Divergence(EV_scalar, fe, D, 1.0);
|
||||||
|
ASSERT_NEAR(EV_scalar(1), 7.0, 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(2), 10.0, 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(3), 0.0, 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(4), 0.0, 1e-13);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
TESTI(TestEqualOrderOperators, Gradient)
|
||||||
|
{
|
||||||
|
FiniteElement fe = getFE();
|
||||||
|
|
||||||
|
Matrix EM_scalar(2*2,2);
|
||||||
|
// single component + pressure
|
||||||
|
EqualOrderOperators::Weak::Gradient(EM_scalar, fe);
|
||||||
|
const DoubleVec EM_scalar_ref = {{-1.0,-2.0},
|
||||||
|
{-3.0,-6.0},
|
||||||
|
{-2.0,-4.0},
|
||||||
|
{-4.0,-8.0}};
|
||||||
|
check_matrix_equal(EM_scalar, EM_scalar_ref);
|
||||||
|
|
||||||
|
Vector EV_scalar(4);
|
||||||
|
EqualOrderOperators::Weak::Gradient(EV_scalar, fe);
|
||||||
|
ASSERT_NEAR(EV_scalar(1), fe.dNdX(1,1), 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(2), fe.dNdX(1,2), 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(3), fe.dNdX(2,1), 1e-13);
|
||||||
|
ASSERT_NEAR(EV_scalar(4), fe.dNdX(2,2), 1e-13);
|
||||||
|
}
|
||||||
|
*/
|
||||||
|
|
||||||
|
|
||||||
|
TEST(TestCompatibleOperators, Laplacian)
|
||||||
|
{
|
||||||
|
MxFiniteElement fe = getFE();
|
||||||
|
|
||||||
|
std::vector<Matrix> EM(3);
|
||||||
|
EM[1].resize(6,6);
|
||||||
|
EM[2].resize(6,6);
|
||||||
|
|
||||||
|
CompatibleOperators::Weak::Laplacian(EM, fe);
|
||||||
|
|
||||||
|
const DoubleVec EM_1_ref = {{ 5.0, 11.0, 17.0, 23.0, 29.0, 35.0},
|
||||||
|
{11.0, 25.0, 39.0, 53.0, 67.0, 81.0},
|
||||||
|
{17.0, 39.0, 61.0, 83.0, 105.0, 127.0},
|
||||||
|
{23.0, 53.0, 83.0, 113.0, 143.0, 173.0},
|
||||||
|
{29.0, 67.0, 105.0, 143.0, 181.0, 219.0},
|
||||||
|
{35.0, 81.0, 127.0, 173.0, 219.0, 265.0}};
|
||||||
|
|
||||||
|
check_matrix_equal(EM[1], EM_1_ref);
|
||||||
|
|
||||||
|
const DoubleVec EM_2_ref = {{365.0, 419.0, 473.0, 527.0, 581.0, 635.0},
|
||||||
|
{419.0, 481.0, 543.0, 605.0, 667.0, 729.0},
|
||||||
|
{473.0, 543.0, 613.0, 683.0, 753.0, 823.0},
|
||||||
|
{527.0, 605.0, 683.0, 761.0, 839.0, 917.0},
|
||||||
|
{581.0, 667.0, 753.0, 839.0, 925.0, 1011.0},
|
||||||
|
{635.0, 729.0, 823.0, 917.0, 1011.0, 1105.0}};
|
||||||
|
|
||||||
|
check_matrix_equal(EM[2], EM_2_ref);
|
||||||
|
|
||||||
|
// stress formulation
|
||||||
|
std::vector<Matrix> EM_stress(21);
|
||||||
|
|
||||||
|
EM_stress[1].resize(6,6);
|
||||||
|
EM_stress[2].resize(6,6);
|
||||||
|
EM_stress[6].resize(6,6);
|
||||||
|
EM_stress[10].resize(6,6);
|
||||||
|
|
||||||
|
CompatibleOperators::Weak::Laplacian(EM_stress, fe, 1.0, true);
|
||||||
|
std::cout << EM_stress[1] << std::endl;
|
||||||
|
|
||||||
|
/*
|
||||||
|
const DoubleVec EM_stress_ref = {{11.0, 3.0, 16.0, 6.0},
|
||||||
|
{ 3.0, 19.0, 4.0, 26.0},
|
||||||
|
{16.0, 4.0, 24.0, 8.0},
|
||||||
|
{ 6.0, 26.0, 8.0, 36.0}};
|
||||||
|
|
||||||
|
check_matrix_equal(EM_stress, EM_stress_ref);
|
||||||
|
|
||||||
|
Matrix EM_coeff(2, 2);
|
||||||
|
EqualOrderOperators::Weak::LaplacianCoeff(EM_coeff, fe.dNdX, fe, 2.0);
|
||||||
|
const DoubleVec EM_coeff_ref = {{104.0, 148.0},
|
||||||
|
{152.0, 216.0}};
|
||||||
|
check_matrix_equal(EM_coeff, EM_coeff_ref);
|
||||||
|
*/
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
TEST(TestCompatibleOperators, Mass)
|
||||||
|
{
|
||||||
|
MxFiniteElement fe = getFE();
|
||||||
|
|
||||||
|
std::vector<Matrix> EM_vec(3);
|
||||||
|
EM_vec[1].resize(6,6);
|
||||||
|
EM_vec[2].resize(6,6);
|
||||||
|
|
||||||
|
CompatibleOperators::Weak::Mass(EM_vec, fe);
|
||||||
|
|
||||||
|
const DoubleVec EM_1_ref = {{1.0, 2.0, 3.0, 4.0, 5.0, 6.0},
|
||||||
|
{2.0, 4.0, 6.0, 8.0, 10.0, 12.0},
|
||||||
|
{3.0, 6.0, 9.0, 12.0, 15.0, 18.0},
|
||||||
|
{4.0, 8.0, 12.0, 16.0, 20.0, 24.0},
|
||||||
|
{5.0, 10.0, 15.0, 20.0, 25.0, 30.0},
|
||||||
|
{6.0, 12.0, 18.0, 24.0, 30.0, 36.0}};
|
||||||
|
|
||||||
|
check_matrix_equal(EM_vec[1], EM_1_ref);
|
||||||
|
|
||||||
|
const DoubleVec EM_2_ref = {{49.0, 56.0, 63.0, 70.0, 77.0, 84.0},
|
||||||
|
{56.0, 64.0, 72.0, 80.0, 88.0, 96.0},
|
||||||
|
{63.0, 72.0, 81.0, 90.0, 99.0, 108.0},
|
||||||
|
{70.0, 80.0, 90.0, 100.0, 110.0, 120.0},
|
||||||
|
{77.0, 88.0, 99.0, 110.0, 121.0, 132.0},
|
||||||
|
{84.0, 96.0, 108.0, 120.0, 132.0, 144.0}};
|
||||||
|
|
||||||
|
check_matrix_equal(EM_vec[2], EM_2_ref);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
TEST(TestCompatibleOperators, Source)
|
||||||
|
{
|
||||||
|
MxFiniteElement fe = getFE();
|
||||||
|
|
||||||
|
Vectors EV_scalar(3);
|
||||||
|
EV_scalar[1].resize(6);
|
||||||
|
EV_scalar[2].resize(6);
|
||||||
|
|
||||||
|
CompatibleOperators::Weak::Source(EV_scalar, fe, 2.0);
|
||||||
|
|
||||||
|
for (size_t i = 1; i <= 6; ++i) {
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[1](i), 2.0*fe.basis(1)(i));
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[2](i), 2.0*fe.basis(2)(i));
|
||||||
|
}
|
||||||
|
|
||||||
|
EV_scalar[1].fill(0.0);
|
||||||
|
EV_scalar[2].fill(0.0);
|
||||||
|
Vec3 f;
|
||||||
|
f[0] = 1.0;
|
||||||
|
f[1] = 2.0;
|
||||||
|
|
||||||
|
CompatibleOperators::Weak::Source(EV_scalar, fe, f, 1.0);
|
||||||
|
for (size_t i = 1; i <= 6; ++i) {
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[1](i), fe.basis(1)(i));
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[2](i), 2.0*fe.basis(2)(i));
|
||||||
|
}
|
||||||
|
|
||||||
|
EV_scalar[1].fill(0.0);
|
||||||
|
EV_scalar[2].fill(0.0);
|
||||||
|
CompatibleOperators::Weak::Source(EV_scalar, fe, f, 2.0);
|
||||||
|
for (size_t i = 1; i <= 6; ++i) {
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[1](i), 2.0*fe.basis(1)(i));
|
||||||
|
ASSERT_FLOAT_EQ(EV_scalar[2](i), 4.0*fe.basis(2)(i));
|
||||||
|
}
|
||||||
|
}
|
||||||
Loading…
Reference in New Issue
Block a user