diff --git a/src/elasticity/assemble/DGElasticity.cpp b/src/elasticity/assemble/DGElasticity.cpp
index 4caa463dd1374cef9eac67dda563754052897b4b..73e39bb2a49f5351be54d9dbb1da2924d2f45e4b 100644
--- a/src/elasticity/assemble/DGElasticity.cpp
+++ b/src/elasticity/assemble/DGElasticity.cpp
@@ -132,6 +132,9 @@ void DGElasticity::Residual(const cell &c, const Vector &U, Vector &r) const {
Tensor F = DeformationGradient(FE, q, U);
VectorField N = mat::cofactor(F) * FE.QNormal(q);
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")) {
+ N = FE.QNormal(q);
+ }
double localPressure = -eProblem.Pressure(Time(), z, E.Bnd(f));
for (int i = 0; i < FE.shape_size(); ++i)
for (int k = 0; k < E.Dim(); ++k) {
@@ -458,6 +461,7 @@ void DGElasticity::Jacobi(const cell &c, const Vector &U, Matrix &A) const {
}
}
if (E.Bnd(face) >= PRESSURE_BC_LV && E.Bnd(face) <= PRESSURE_BC_RA) {
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")) continue;
for (int q = 0; q < FE.nQ(); ++q) {
double w = FE.QWeight(q);
VectorField N = FE.QNormal(q);
diff --git a/src/elasticity/assemble/EGElasticity.cpp b/src/elasticity/assemble/EGElasticity.cpp
index 6dafa3ed8925733904555f312185d347fd3d551c..ea8c1c4eecb70947ff7696f7b337b69c67da1056 100644
--- a/src/elasticity/assemble/EGElasticity.cpp
+++ b/src/elasticity/assemble/EGElasticity.cpp
@@ -123,6 +123,9 @@ void EGElasticity::Residual(const cell &c, const Vector &U, Vector &r) const {
Point z = FE.QPoint(q);
Tensor F = DeformationGradient(FE, q, U);
VectorField N = mat::cofactor(F) * FE.QNormal(q);
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")){
+ N = FE.QNormal(q);
+ }
double localPressure = -eProblem.Pressure(Time(), z, E.Bnd(f));
for (int i = 0; i < FE.ValueSize(); ++i) {
for (int k = 0; k < FE.Dim(); ++k) {
@@ -408,6 +411,7 @@ void EGElasticity::Jacobi(const cell &c, const Vector &U, Matrix &A) const {
w * (SN_theta * sign + SN_Dtheta - s * (theta * theta));
}
} else if (E.Bnd(face) >= PRESSURE_BC_LV && E.Bnd(face) <= PRESSURE_BC_RA) {
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")) continue;
for (int q = 0; q < FE.nQ(); ++q) {
auto theta = FE.PenaltyVectorValue(q);
auto Dtheta = FE.PenaltyVectorGradient(q);
diff --git a/src/elasticity/assemble/LagrangeElasticity.cpp b/src/elasticity/assemble/LagrangeElasticity.cpp
index dcffb7071912e6c302b0adf5f445e2a75d9bffdd..0544f5f118e661f0807f559dc9676c3e07c7b31c 100644
--- a/src/elasticity/assemble/LagrangeElasticity.cpp
+++ b/src/elasticity/assemble/LagrangeElasticity.cpp
@@ -1,3 +1,4 @@
+
#include "LagrangeElasticity.hpp"
#include "LagrangeTransfer.hpp"
@@ -187,11 +188,15 @@ void LagrangeElasticity::PressureBoundary(const cell &c, int face, int bc, const
Vector &r) const {
VectorFieldFaceElement FE(U, *c, face);
RowValues r_c(r, FE);
+ const Material &cellMat = eProblem.GetMaterial(c);
for (int q = 0; q < FE.nQ(); ++q) {
double w = FE.QWeight(q);
Tensor F = DeformationGradient(FE, q, U);
VectorField N = mat::cofactor(F) * FE.QNormal(q);
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")) {
+ N = FE.QNormal(q);
+ }
double localPressure = -eProblem.Pressure(Time(), FE.QPoint(q), bc);
for (int i = 0; i < FE.NodalPoints(); ++i) {
@@ -267,7 +272,7 @@ void LagrangeElasticity::Jacobi(const cell &c, const Vector &U, Matrix &A) const
}
}
}
-
+ if (contains(cellMat.Name(),"Linear") || contains(cellMat.Name(),"Laplace")) return;
for (int face = 0; face < c.Faces(); ++face) {
if (E.Bnd(face) >= PRESSURE_BC_LV && E.Bnd(face) <= PRESSURE_BC_RA) {
VectorFieldFaceElement FE(U, *c, face);
diff --git a/src/elasticity/problem/BeamProblems.hpp b/src/elasticity/problem/BeamProblems.hpp
index 54eec3b54d00a0b0c1f2f25eb5c08ae461a8a1cc..0458df06eb54fe1e2f0540526b31b9e6dcf4ddd3 100644
--- a/src/elasticity/problem/BeamProblems.hpp
+++ b/src/elasticity/problem/BeamProblems.hpp
@@ -109,7 +109,7 @@ protected:
if (bndCondition == 230) {
return -(2 * mu + lambda);
} else if (bndCondition == 231 || bndCondition == 232) {
- return -lambda / 2;
+ return -lambda;
} else {
return 0;
}
@@ -145,7 +145,7 @@ protected:
if (bndCondition == 230) {
return (4 * mu + 2 * lambda);
} else if (bndCondition == 231 || bndCondition == 232) {
- return -2 * lambda;
+ return 2 * lambda;
} else {
return 0;
}
diff --git a/src/elasticity/problem/DynamicBoundaryProblems.cpp b/src/elasticity/problem/DynamicBoundaryProblems.cpp
index 5176b99f39779cd1eba3f13988972f4f6bfc1819..e10228d3eb5de7d3d20e56c05bb98f1712a41815 100644
--- a/src/elasticity/problem/DynamicBoundaryProblems.cpp
+++ b/src/elasticity/problem/DynamicBoundaryProblems.cpp
@@ -201,11 +201,11 @@ VectorField DPolynomialNeumannProblem::Load(double time, const Point &x, const C
double DPolynomialNeumannProblem::pressure(double t, const Point &x, int bndCondition) const {
if (bndCondition == 230) {
- return - ((t*t*x[0]) / ((1 + t*t*x[1]) * (1 + t*t*x[2])));
+ return - (t*t*x[0]);
} else if (bndCondition == 231) {
- return -((t*t*x[1]) / ((1 + t*t* x[0]) * (1 + t*t*x[2])));
+ return -(t*t*x[1]);
} else if (bndCondition == 232) {
- return - ((t*t*x[2]) / ((1 + t*t*x[0]) * (1 + t*t*x[1])));
+ return - (t*t*x[2]);
} else {
return 0;
}
@@ -246,14 +246,11 @@ VectorField DExponentialNeumannProblem::Load(double time, const Point &x, const
double DExponentialNeumannProblem::pressure(double t, const Point &x, int bndCondition) const {
if (bndCondition == 230) {
- return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[0]) /((sin(0.5 * Pi * x[1]) * sin(0.5 * Pi * t) + 1) *
- (sin(0.5 * Pi * x[2]) * sin(0.5 * Pi * t) + 1));
+ return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[0]) ;
} else if (bndCondition == 231) {
- return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[1]) /((sin(0.5 * Pi * x[0]) * sin(0.5 * Pi * t) + 1) *
- (sin(0.5 * Pi * x[2]) * sin(0.5 * Pi * t) + 1));
+ return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[1]);
} else if (bndCondition == 232) {
- return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[2]) /((sin(0.5 * Pi * x[1]) * sin(0.5 * Pi * t) + 1) *
- (sin(0.5 * Pi * x[0]) * sin(0.5 * Pi * t) + 1));
+ return - sin(0.5 * Pi * t) * sin(0.5 * Pi * x[2]) ;
} else {
return 0;
}