Commit 23ffbab4 authored by Andrew Heather's avatar Andrew Heather
Browse files

BUG: lagrangian/intermediate - clause in table integration not handled

parent b3630c3f
......@@ -94,46 +94,56 @@ Type Foam::Table<Type>::integrate(const scalar x1, const scalar x2) const
Type sum = pTraits<Type>::zero;
// Return zero if out of bounds
if ((x1 > table_[table_.size()-1].first()) || (x2 < table_[0].first()))
if ((x1 > table_.last().first()) || (x2 < table_[0].first()))
{
return sum;
}
// Find start index
// Find next index greater than x1
label id1 = 0;
while ((table_[id1].first() < x1) && (id1 < table_.size()))
{
id1++;
}
// Find end index
// Find next index less than x2
label id2 = table_.size() - 1;
while ((table_[id2].first() > x2) && (id2 >= 1))
{
id2--;
}
// Integrate table body
for (label i=id1; i<id2; i++)
if ((id1 - id2) == 1)
{
sum +=
(table_[i].second() + table_[i+1].second())
* (table_[i+1].first() - table_[i].first());
// x1 and x2 lie within 1 interval
sum = 0.5*(value(x1) + value(x2))*(x2 - x1);
}
sum *= 0.5;
// Add table ends (partial segments)
if (id1 > 0)
{
sum += 0.5
* (value(x1) + table_[id1].second())
* (table_[id1].first() - x1);
}
if (id2 < table_.size() - 1)
else
{
sum += 0.5
* (table_[id2].second() + value(x2))
* (x2 - table_[id2].first());
// x1 and x2 cross multiple intervals
// Integrate table body
for (label i=id1; i<id2; i++)
{
sum +=
(table_[i].second() + table_[i+1].second())
* (table_[i+1].first() - table_[i].first());
}
sum *= 0.5;
// Add table ends (partial segments)
if (id1 > 0)
{
sum += 0.5
* (value(x1) + table_[id1].second())
* (table_[id1].first() - x1);
}
if (id2 < table_.size() - 1)
{
sum += 0.5
* (table_[id2].second() + value(x2))
* (x2 - table_[id2].first());
}
}
return sum;
......
Supports Markdown
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