Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Consider the derivatives at boundaries in the 2D-table #3893

Merged
merged 15 commits into from
Feb 16, 2023
103 changes: 68 additions & 35 deletions Modelica/Resources/C-Sources/ModelicaStandardTables.c
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,14 @@
Modelica.Blocks.Tables.CombiTable2Dv

Changelog:

Jan. 31, 2022: by Hans Olsson, Dassault Systemes
Added better support for one-sided derivatives of CombiTable2D.
The idea is that when we are computing the derivative at a boundary
in the table we should consider the der-value to choose side.
This is less important for 1d-tables and thus ignored in those cases.
beutlich marked this conversation as resolved.
Show resolved Hide resolved
(ticket #3893)

Nov. 12, 2021: by Thomas Beutlich
Fixed derivatives in CombiTable2D for one-sided extrapolation
by constant continuation (ticket #3894)
Expand Down Expand Up @@ -484,17 +492,28 @@ extern int usertab(char* tableName, int nipo, int dim[], int* colWise,
static int isNearlyEqual(double x, double y);
/* Compare two floating-point numbers by threshold _EPSILON */

static size_t findRowIndex(const double* table, size_t nRow, size_t nCol,
size_t last, double x);
static size_t findRowIndex(_In_ const double* table, size_t nRow, size_t nCol,
size_t last, double x) MODELICA_NONNULLATTR;
/* Find the row index i using binary search such that
* i + 1 < nRow
* table[i*nCol] <= x
* table[(i + 1)*nCol] > x for i + 2 < nRow
*/

static size_t findRowIndex2(_In_ const double* table, size_t nRow, size_t nCol,
size_t last, double x, double dx) MODELICA_NONNULLATTR;
/* Using dx as tie-breaker if table[i*nCol] == x to treat x as x+dx*eps */

static size_t findColIndex(_In_ const double* table, size_t nCol, size_t last,
double x) MODELICA_NONNULLATTR;
/* Same as findRowIndex but works on rows */
/* Same as findRowIndex but works on columns */

static size_t findColIndex2(_In_ const double* table, size_t nCol, size_t last,
double x, double dx) MODELICA_NONNULLATTR;
/* Same as findRowIndex2 but works on columns */

static int isLess(double x, double dx, double val);
/* Check, whether x is less than val, also using dx as tie-breaker */

static int isValidName(_In_z_ const char* name) MODELICA_NONNULLATTR;
/* Check, whether a file or table name is valid */
Expand Down Expand Up @@ -4237,21 +4256,21 @@ double ModelicaStandardTables_CombiTable2D_getDerValue(void* _tableID, double u1
u2 -= T;
} while (u2 > u2Max);
}
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}
else if (u2 < u2Min) {
else if (isLess(u2, der_u2, u2Min)) {
extrapolate2 = LEFT;
last2 = 0;
}
else if (u2 > u2Max) {
else if (!isLess(u2, der_u2, u2Max)) {
beutlich marked this conversation as resolved.
Show resolved Hide resolved
extrapolate2 = RIGHT;
last2 = nCol - 3;
}
else {
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}

Expand Down Expand Up @@ -4382,21 +4401,21 @@ double ModelicaStandardTables_CombiTable2D_getDerValue(void* _tableID, double u1
u1 -= T;
} while (u1 > u1Max);
}
last1 = findRowIndex(&TABLE(1, 0), nRow - 1, nCol,
tableID->last1, u1);
last1 = findRowIndex2(&TABLE(1, 0), nRow - 1, nCol,
tableID->last1, u1, der_u1);
tableID->last1 = last1;
}
else if (u1 < u1Min) {
else if (isLess(u1, der_u1, u1Min)) {
extrapolate1 = LEFT;
last1 = 0;
}
else if (u1 > u1Max) {
else if (!isLess(u1, der_u1, u1Max)) {
extrapolate1 = RIGHT;
last1 = nRow - 3;
}
else {
last1 = findRowIndex(&TABLE(1, 0), nRow - 1, nCol,
tableID->last1, u1);
last1 = findRowIndex2(&TABLE(1, 0), nRow - 1, nCol,
tableID->last1, u1, der_u1);
tableID->last1 = last1;
}
if (nCol == 2) {
Expand Down Expand Up @@ -4526,21 +4545,21 @@ double ModelicaStandardTables_CombiTable2D_getDerValue(void* _tableID, double u1
u2 -= T;
} while (u2 > u2Max);
}
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}
else if (u2 < u2Min) {
else if (isLess(u2, der_u2, u2Min)) {
extrapolate2 = LEFT;
last2 = 0;
}
else if (u2 > u2Max) {
else if (!isLess(u2, der_u2, u2Max)) {
extrapolate2 = RIGHT;
last2 = nCol - 3;
}
else {
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}

Expand Down Expand Up @@ -5158,8 +5177,8 @@ double ModelicaStandardTables_CombiTable2D_getDer2Value(void* _tableID, double u
u2 -= T;
} while (u2 > u2Max);
}
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}
else if (u2 < u2Min) {
Expand All @@ -5171,8 +5190,8 @@ double ModelicaStandardTables_CombiTable2D_getDer2Value(void* _tableID, double u
last2 = nCol - 3;
}
else {
last2 = findColIndex(&TABLE(0, 1), nCol - 1,
tableID->last2, u2);
last2 = findColIndex2(&TABLE(0, 1), nCol - 1,
tableID->last2, u2, der_u2);
tableID->last2 = last2;
}

Expand Down Expand Up @@ -6168,14 +6187,14 @@ static int isNearlyEqual(double x, double y) {
return fabs(y - x) < cmp;
}

static size_t findRowIndex(const double* table, size_t nRow, size_t nCol,
size_t last, double x) {
static size_t findRowIndex2(const double* table, size_t nRow, size_t nCol,
size_t last, double x, double dx) {
size_t i0 = 0;
size_t i1 = nRow - 1;
if (x < TABLE_COL0(last)) {
if (isLess(x, dx, TABLE_COL0(last))) {
i1 = last;
}
else if (x >= TABLE_COL0(last + 1)) {
else if (!isLess(x, dx, TABLE_COL0(last + 1))) {
i0 = last;
}
else {
Expand All @@ -6185,7 +6204,7 @@ static size_t findRowIndex(const double* table, size_t nRow, size_t nCol,
/* Binary search */
while (i1 > i0 + 1) {
const size_t i = (i0 + i1)/2;
if (x < TABLE_COL0(i)) {
if (isLess(x, dx, TABLE_COL0(i))) {
i1 = i;
}
else {
Expand All @@ -6195,14 +6214,19 @@ static size_t findRowIndex(const double* table, size_t nRow, size_t nCol,
return i0;
}

static size_t findColIndex(_In_ const double* table, size_t nCol, size_t last,
double x) {
static size_t findRowIndex(const double* table, size_t nRow, size_t nCol,
size_t last, double x) {
return findRowIndex2(table, nRow, nCol, last, x, 0.0);
}

static size_t findColIndex2(_In_ const double* table, size_t nCol, size_t last,
double x, double dx) {
size_t i0 = 0;
size_t i1 = nCol - 1;
if (x < TABLE_ROW0(last)) {
if (isLess(x, dx, TABLE_ROW0(last))) {
i1 = last;
}
else if (x >= TABLE_ROW0(last + 1)) {
else if (isLess(x, dx, TABLE_ROW0(last + 1))) {
i0 = last;
}
else {
Expand All @@ -6212,7 +6236,7 @@ static size_t findColIndex(_In_ const double* table, size_t nCol, size_t last,
/* Binary search */
while (i1 > i0 + 1) {
const size_t i = (i0 + i1)/2;
if (x < TABLE_ROW0(i)) {
if (isLess(x, dx, TABLE_ROW0(i))) {
i1 = i;
}
else {
Expand All @@ -6222,8 +6246,17 @@ static size_t findColIndex(_In_ const double* table, size_t nCol, size_t last,
return i0;
}

static size_t findColIndex(_In_ const double* table, size_t nCol, size_t last, double x) {
return findColIndex2(table, nCol, last, x, 0.0);
}


/* ----- Internal check functions ----- */

static int isLess(double x, double dx, double val) {
return x < val || (x == val && dx < 0);
beutlich marked this conversation as resolved.
Show resolved Hide resolved
}

static int isValidName(_In_z_ const char* name) {
int isValid = 0;
if (NULL != name) {
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
time
der1.y
der2.y
der3.y
der4.y
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
time
der1.y
der2.y
der3.y
der4.y
88 changes: 88 additions & 0 deletions ModelicaTest/Tables/CombiTable2Ds.mo
Original file line number Diff line number Diff line change
Expand Up @@ -1125,4 +1125,92 @@ double mydummyfunc(double dummy_in) {
points={{-59,-10},{-52,-10},{-52,4},{-42,4}}, color={0,0,127}));
annotation (experiment(StartTime=0, StopTime=14));
end Test33;

model OneSidedDerivative2D "Test of one sided derivatives in 2D-tables (Ticket #3893)"
// We are starting at boundaries of the table
// Case 2 and 4 slide diagonally (from different corners)
// Case 1 and 3 start at the edge of the real table and then leave it top-left
// Interpolating in different ways
extends Modelica.Icons.Example;
parameter Real M0[:,:]=[0,-10,1,2,10; -10,2,2,3,4; 1,2,2,3,3; 2,3,3,4,4; 10,3,
3,4,4];
parameter Real M[:,:]=[M0[1:1,1:1],M0[1:1,3:end-1];M0[3:end-1,1],M0[3:end-1,3:end-1]];
Modelica.Blocks.Sources.Ramp ramp(
height=1,
duration=1,
offset=1)
annotation (Placement(transformation(extent={{-80,40},{-60,60}})));
Modelica.Blocks.Sources.Ramp ramp1(
height=-1,
duration=1,
offset=2)
annotation (Placement(transformation(extent={{-20,18},{0,38}})));
Modelica.Blocks.Tables.CombiTable2Ds combiTable2Ds2(
tableOnFile=false,
table=M0,
smoothness=Modelica.Blocks.Types.Smoothness.LinearSegments,
extrapolation=Modelica.Blocks.Types.Extrapolation.LastTwoPoints)
annotation (Placement(transformation(extent={{40,60},{60,80}})));
Modelica.Blocks.Continuous.Der der2
annotation (Placement(transformation(extent={{80,60},{100,80}})));
Modelica.Blocks.Tables.CombiTable2Ds combiTable2Ds4(
tableOnFile=false,
table=M0,
smoothness=Modelica.Blocks.Types.Smoothness.LinearSegments,
extrapolation=Modelica.Blocks.Types.Extrapolation.LastTwoPoints)
annotation (Placement(transformation(extent={{40,0},{60,20}})));
Modelica.Blocks.Continuous.Der der4
annotation (Placement(transformation(extent={{80,0},{100,20}})));
Modelica.Blocks.Sources.Ramp ramp2(
height=-1,
duration=1,
offset=1)
annotation (Placement(transformation(extent={{-80,-80},{-60,-60}})));
Modelica.Blocks.Tables.CombiTable2Ds combiTable2Ds1(
tableOnFile=false,
table=M0,
smoothness=Modelica.Blocks.Types.Smoothness.LinearSegments,
extrapolation=Modelica.Blocks.Types.Extrapolation.LastTwoPoints)
annotation (Placement(transformation(extent={{-20,-40},{0,-20}})));
Modelica.Blocks.Continuous.Der der1
annotation (Placement(transformation(extent={{20,-40},{40,-20}})));
Modelica.Blocks.Tables.CombiTable2Ds combiTable2Ds3(
tableOnFile=false,
table=M,
smoothness=Modelica.Blocks.Types.Smoothness.LinearSegments,
extrapolation=Modelica.Blocks.Types.Extrapolation.HoldLastPoint)
annotation (Placement(transformation(extent={{-20,-80},{0,-60}})));
Modelica.Blocks.Continuous.Der der3
annotation (Placement(transformation(extent={{20,-80},{40,-60}})));
equation
connect(combiTable2Ds2.y, der2.u)
annotation (Line(points={{61,70},{78,70}}, color={0,0,127}));
connect(ramp.y, combiTable2Ds2.u1) annotation (Line(points={{-59,50},{18,50},{
18,76},{38,76}}, color={0,0,127}));
connect(ramp1.y, combiTable2Ds2.u2) annotation (Line(points={{1,28},{24,28},{24,
64},{38,64}}, color={0,0,127}));
connect(combiTable2Ds4.y, der4.u)
annotation (Line(points={{61,10},{78,10}}, color={0,0,127}));
connect(ramp.y, combiTable2Ds4.u2) annotation (Line(points={{-59,50},{-38,50},
{-38,4},{38,4}}, color={0,0,127}));
connect(combiTable2Ds4.u1, combiTable2Ds2.u2) annotation (Line(points={{38,16},
{34,16},{34,18},{24,18},{24,64},{38,64}}, color={0,0,127}));
connect(ramp2.y, combiTable2Ds1.u1) annotation (Line(points={{-59,-70},{-30,-70},
{-30,-24},{-22,-24}}, color={0,0,127}));
connect(ramp2.y, combiTable2Ds1.u2) annotation (Line(points={{-59,-70},{-30,-70},
{-30,-36},{-22,-36}}, color={0,0,127}));
connect(combiTable2Ds1.y, der1.u)
annotation (Line(points={{1,-30},{18,-30}}, color={0,0,127}));
connect(combiTable2Ds3.u1, combiTable2Ds1.u1) annotation (Line(points={{-22,-64},
{-26,-64},{-26,-70},{-30,-70},{-30,-24},{-22,-24}}, color={0,0,127}));
connect(combiTable2Ds3.u2, combiTable2Ds1.u1) annotation (Line(points={{-22,-76},
{-24,-76},{-24,-70},{-30,-70},{-30,-24},{-22,-24}}, color={0,0,127}));
connect(combiTable2Ds3.y, der3.u)
annotation (Line(points={{1,-70},{18,-70}}, color={0,0,127}));
annotation (
experiment(
StartTime=-1,
StopTime=4));
end OneSidedDerivative2D;

end CombiTable2Ds;
Loading