OilfieldToolsNet/IFEM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

added: allow explicitly specifying time derivative of function

Browse Source
This commit is contained in:
Arne Morten Kvarving 2023-09-21 13:22:05 +02:00
parent 72f9748aa1
commit a8a4cd780e
3 changed files with 55 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

5
src/Utility/ExprFunctions.C
View File

@ -286,7 +286,7 @@ void EvalFunction::addDerivative (const std::string& function,
const std::string& variables,
int d1, int d2)
{
if (d1 > 0 && d1 <= 3 && d2 < 1) // A first derivative is specified
if (d1 > 0 && d1 <= 4 && d2 < 1) // A first derivative is specified
{
if (!gradient[--d1])
gradient[d1] = new EvalFunction((variables+function).c_str());
@ -342,6 +342,9 @@ Real EvalFunction::deriv (const Vec3& X, int dir) const
}
else if (!IAmConstant)
{
if (gradient[3])
return gradient[3]->evaluate(X);
// Evaluate time-derivative using central difference
Vec4 X0, X1;
X0.assign(X); X0.t -= 0.5*dt;

2
src/Utility/ExprFunctions.h
View File

@ -95,7 +95,7 @@ class EvalFunction : public RealFunc
std::vector<Arg> arg; //!< Function argument values
std::array<EvalFunction*,3> gradient; //!< First derivative expressions
std::array<EvalFunction*,4> gradient; //!< First derivative expressions
std::array<EvalFunction*,6> dgradient; //!< Second derivative expressions
bool IAmConstant; //!< Indicates whether the time coordinate is given or not

50
src/Utility/Test/TestFunctions.C
View File

@ -65,3 +65,53 @@ TEST(TestEvalFunction, isConstant)
EXPECT_TRUE (EvalFunction("1.8*tan(x)*x").isConstant());
EXPECT_FALSE(EvalFunction("2.0*x*tan(t*3)+y").isConstant());
}
TEST(TestEvalFunction, Derivatives)
{
const char* g = "sin(x)*sin(y)*sin(z)*sin(t)";
const char* g_x = "cos(x)*sin(y)*sin(z)*sin(t)";
const char* g_y = "sin(x)*cos(y)*sin(z)*sin(t)";
const char* g_z = "sin(x)*sin(y)*cos(z)*sin(t)";
const char* g_t = "sin(x)*sin(y)*sin(z)*cos(t)";
const char* g_xx = "-sin(x)*sin(y)*sin(z)*sin(t)";
const char* g_xy = "cos(x)*cos(y)*sin(z)*sin(t)";
const char* g_xz = "cos(x)*sin(y)*cos(z)*sin(t)";
const char* g_yy = "-sin(x)*sin(y)*sin(z)*sin(t)";
const char* g_yz = "sin(x)*cos(y)*cos(z)*sin(t)";
const char* g_zz = "-sin(x)*sin(y)*sin(z)*sin(t)";
EvalFunction f(g);
f.addDerivative(g_x,"",1);
f.addDerivative(g_y,"",2);
f.addDerivative(g_z,"",3);
f.addDerivative(g_t,"",4);
f.addDerivative(g_xx,"",1,1);
f.addDerivative(g_xy,"",1,2);
f.addDerivative(g_xz,"",1,3);
f.addDerivative(g_yy,"",2,2);
f.addDerivative(g_yz,"",2,3);
f.addDerivative(g_zz,"",3,3);
const Vec4 X(1.0,2.0,3.0,4.0);
for (const double t : {0.1, 0.2, 0.3})
for (const double x : {0.4, 0.5, 0.6})
for (const double y : {0.7, 0.8, 0.9})
for (const double z : {1.0, 1.1, 1.2}) {
const Vec4 X(x,y,z,t);
EXPECT_DOUBLE_EQ(f(X), sin(x)*sin(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.deriv(X,1), cos(x)*sin(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.deriv(X,2), sin(x)*cos(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.deriv(X,3), sin(x)*sin(y)*cos(z)*sin(t));
EXPECT_DOUBLE_EQ(f.deriv(X,4), sin(x)*sin(y)*sin(z)*cos(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,1,1), -sin(x)*sin(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,1,2), cos(x)*cos(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,1,3), cos(x)*sin(y)*cos(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,2,1), cos(x)*cos(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,2,2), -sin(x)*sin(y)*sin(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,2,3), sin(x)*cos(y)*cos(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,3,1), cos(x)*sin(y)*cos(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,3,2), sin(x)*cos(y)*cos(z)*sin(t));
EXPECT_DOUBLE_EQ(f.dderiv(X,3,3), -sin(x)*sin(y)*sin(z)*sin(t));
}
}