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

Bug in 3D-property initialization #1017

Open
andlaus opened this issue Dec 23, 2016 · 14 comments
Open

Bug in 3D-property initialization #1017

andlaus opened this issue Dec 23, 2016 · 14 comments

Comments

@andlaus
Copy link
Contributor

andlaus commented Dec 23, 2016

edit: although there seems to be a bug related to the property initialization, the initial analysis is not correct. see below: #1017 (comment)

The Norne deck from opm-data features the "SWL" and "ISWL" double grid properties and the latter is supposed to be identical to the former because it is copied via the COPY statement on line 342 of the deck's main file. Unfortunately opm-parser currently thinks otherwise: Running the following program on the Norne deck

#include <opm/parser/eclipse/Parser/Parser.hpp>
#include <opm/parser/eclipse/Deck/Deck.hpp>
#include <opm/parser/eclipse/EclipseState/EclipseState.hpp>

#include <iostream>
#include <tuple>
#include <vector>

int main(int argc, char** argv)
{
    if (argc != 2) {
        std::cout << "Usage: " << argv[0] << " FILNAME\n";
        return 1;
    }

    Opm::Parser parser;
    typedef std::pair<std::string, Opm::InputError::Action> ParseModePair;
    typedef std::vector<ParseModePair> ParseModePairs;
    ParseModePairs tmp;
    tmp.push_back(ParseModePair(Opm::ParseContext::PARSE_RANDOM_SLASH , Opm::InputError::IGNORE));
    Opm::ParseContext parseContext(tmp);

    const Opm::Deck deck = parser.parseFile(argv[1], parseContext);
    const Opm::EclipseState eclState(deck, parseContext);

    if (!eclState.get3DProperties().hasDeckDoubleGridProperty("SWL"))
        std::cout << "PARSER BUG: the Norne deck features the SWL property!\n";

    if (!eclState.get3DProperties().hasDeckDoubleGridProperty("ISWL"))
        std::cout << "PARSER BUG: the Norne deck features the ISWL property!\n";

    const std::vector<double> &swl =
        eclState.get3DProperties().getDoubleGridProperty("SWL").getData();

    const std::vector<double> &iswl =
        eclState.get3DProperties().getDoubleGridProperty("ISWL").getData();

    for (unsigned i = 0; i < swl.size(); ++i) {
        if (swl[i] != iswl[i]) {
            std::cout << "PARSER BUG: In the Norne deck the SWL and ISWL properties are identical!\n";
            std::cout << "(For cell " << i << " opm-parser says SWL is " << swl[i] << " and ISWL is " << iswl[i] << ".)\n";
            break;
        }
    }

    return 0;
}

yields

PARSER BUG: In the Norne deck the SWL and ISWL properties are identical!
(For cell 0 opm-parser says SWL is 0.04 and ISWL is 0.)
@joakim-hove
Copy link
Member

Well - there might certainly be bugs, but I do not agree with the current analysis. The relevant part of the COPY statement is:

COPY
 ...
 'SWL'   'ISWL'    1 46 1 112  5  22 /
 ...
/

i.e. only the layers k ∈ [4,21] are copied. So - the difference between SWL and ISWL for the first elements is legitimate. Since the ISWL is not explicitly assigned to in the deck it will assume default values for the layers k ∈ [0,3].

@joakim-hove
Copy link
Member

As I wrote in my previous comment I do not think this is a bug; but of course please elaborate and give more details if you think I misunderstand the case.

@andlaus
Copy link
Contributor Author

andlaus commented Dec 23, 2016

I suppose you're right. That said, I'm pretty sure this was not intended by the deck's authors: it makes the endpoints of imbibition and drainage for hysteresis mismatch which is disallowed by the documentation as far as I can remember.

anyway, closing.

@andlaus andlaus closed this as completed Dec 23, 2016
@joakim-hove
Copy link
Member

I'm pretty sure this was not intended by the deck's authors:

That I don't know of course; but there are no active cells in the top four layers so my hunch is that the original deck authors used some ill-advised assignment saving algorithm.

@andlaus
Copy link
Contributor Author

andlaus commented Dec 23, 2016

but there are no active cells in the top four layers

really? the reason I've stumbled over this was that the material laws initialize differently for the second active cell. maybe I have another bug, though...

@joakim-hove
Copy link
Member

really?

Ehhh no - sorry; my bad. There are certainly active cells in the four top layers.

When looking at the INIT file from ECLIPSE the first values of the ISWL keyword are -100000002004087734272.000 - that could be the sign of some "magic" unitialized value which eclipse makes sensible use of? So - there might still be a bug (or maybe more precisely: not-correctly- implemented-undefined-behavior-in-the-case-of-invalid-user-input)? Feel free to reopen if you feel appropriate - although I think the initial analysis was to simple.

image

@andlaus
Copy link
Contributor Author

andlaus commented Dec 23, 2016

Feel free to reopen if you feel appropriate - although I think the initial analysis was to simple.

ok, I'll reopen it as a reminder. Merry Christmas, too :)

@andlaus andlaus reopened this Dec 23, 2016
@OPMUSER
Copy link

OPMUSER commented Dec 24, 2016 via email

@joakim-hove
Copy link
Member

OPM. uses zero offset variables internally - I just did that subtraction upfront.

@OPMUSER
Copy link

OPMUSER commented Dec 24, 2016

Joakim,

I'm confused ,the mods should apply to layers 5 to 22 in the user domain, is that equivalent to 4 to 21 in OPM domain? If this was tranmissibility then it would either be 4 to 21 or 6 to 23 depending on the numerical scheme.

@jokva
Copy link
Contributor

jokva commented Dec 24, 2016

Internally, OPM represents layer 1 with numerical zero.

@andlaus
Copy link
Contributor Author

andlaus commented Dec 24, 2016

The confusion comes from FORTRAN vs C style indices: FORTRAN (and the ECL file format) uses the indices 1, .., n to access arrays of size n, in C it is 0, ..., n-1. since OPM is written in C(++) and Eclipse in FORTRAN there's this small "impedance mismatch".

@OPMUSER
Copy link

OPMUSER commented Dec 24, 2016 via email

@OPMUSER
Copy link

OPMUSER commented Dec 24, 2016 via email

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

4 participants