Commit 1591dea5 authored by niklas.baumgarten's avatar niklas.baumgarten
Browse files

rmv inverted permeablitiy

parent 826e04e2
......@@ -88,7 +88,7 @@ void HybridEllipticAssemble::LocalProblem(const cell &c, const Vector &u,
RT0_LagrangeElementT<> elem(*disc, u, c);
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
for (int i = 0; i < elem.VelocitySize(); ++i) {
for (int k = 0; k < elem.VelocityMaxk(i); ++k) {
Velocity IK_U_i = IK * elem.VelocityField(q, i, k);
......@@ -391,7 +391,7 @@ double HybridEllipticAssemble::H1(const Vector &u) const {
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Scalar U = Value(c, u, elem);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
VectorField DU = IK * Flux(c, u, elem, q);
err += w * U * U + w * DU * DU;
}
......
......@@ -38,7 +38,7 @@ double MixedEllipticAssemble::Energy(const Vector &u) const {
RT0_LagrangeElementT<> elem(*disc, u, c);
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
VectorField Q = elem.VelocityField(q, u);
energy += w * IK * Q * Q;
}
......@@ -53,7 +53,7 @@ void MixedEllipticAssemble::Residual(const cell &c, const Vector &u, Vector &r)
double w = elem.QWeight(q);
Scalar f = problem->Load(c.Subdomain(), elem.QPoint(q));
VectorField Q = elem.VelocityField(q, u);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
Scalar divQ = elem.VelocityFieldDivergence(q, u);
Scalar p = elem.PressureValue(q, u);
for (int i = 0; i < elem.VelocitySize(); ++i) {
......@@ -91,7 +91,7 @@ void MixedEllipticAssemble::Jacobi(const cell &c, const Vector &u, Matrix &A) co
MixedRowEntries A_c(A, c, elem);
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
for (int i = 0; i < elem.VelocitySize(); ++i) {
for (int k = 0; k < elem.VelocityMaxk(i); ++k) {
VectorField Q_i = elem.VelocityField(q, i, k);
......@@ -123,7 +123,7 @@ double MixedEllipticAssemble::EnergyError(const Vector &u) const {
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Point z = elem.QPoint(q);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
VectorField Q = elem.VelocityField(q, u);
VectorField F = problem->Flux(c.Subdomain(), z);
VectorField Diff = Q - F;
......@@ -152,7 +152,7 @@ double MixedEllipticAssemble::H1(const Vector &u) const {
RT0_LagrangeElementT<> elem(*disc, u, c);
for (int q = 0; q < elem.nQ(); ++q) {
double w = elem.QWeight(q);
Tensor IK = problem->InversePermeability(c);
Tensor IK = Invert(problem->Permeability(c));
VectorField Q = elem.VelocityField(q, u);
Scalar p = elem.PressureValue(q, u);
VectorField IK_Q = IK * Q;
......
......@@ -17,12 +17,6 @@ public:
virtual Scalar Load(int, const Point &x) const = 0;
virtual Tensor Permeability(const cell &c) const = 0;
// Todo remove
virtual Tensor InversePermeability(const cell &c) const {
Tensor permeabilityTensor = Permeability(c);
return Invert(permeabilityTensor);
}
};
typedef IStochasticEllipticProblemT<> IStochasticEllipticProblem;
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment