As a result we only have one faceCentroid(int) function the returns const double* and
FieldVector<double,3> for UnstructuredGrid and CpGrid, respectively. The codes is
adapted to this.
For historic (or no apparent) reason the free function grid interface
was added to opm-autodiff. As it depends on whether or not dune-cornerpoint
is present this seems rather unnatural. Therefore this commit removes the
functionality unconditionally from opm-autodiff.
Note that there is a corresponding commit in dune-cornerpoint that adds it
there.
If oscillating then the same measures have to be taken on all
processes. This was not the case before the commit. Instead only
process 0 would have adapted omegas and we have produced wrong
results.
Kudos to Atgeirr for detecting this while reading the code.
With now generic implementation of the initStateEquil in opm-core
we added the necessary grid helper functionlality for CpGrid and activated
the processing if the EQUIL keyword is there.
The average well block pressure is used instead of the well cell
pressure when the well properties are evaluated.
Temperature, rs, rv, phase conditions are still well cells values.
Perforation pressures are stored in the well state
Previously BlackoilPropsDataHandle did hold a grid for sending
and receiving that were either not used or we could prevent their
usage. Therefore this commit removes them from the class and queries
all needed information from the property objects.
These were in the anonymous namespace, which was not wrong but lead to
false-positive unused function warnings on some compilers, probably due
to being used only from within templates. This allows us to keep the
warning activated by default.
1) swatinit() is changed to setSwatInitScaling() to make it obvious that
we are modifying the props.
2) the descriptions of saturation and pc now makes more sense
3) the method is removed from the sibling class and the interface and
the type of new_props is changed from BlackoilPropsAdInterface to
BlackoilPropsAdFromDeck
5) The same modification is added to sim_fibo_ad_cp
The capillary pressure function in new_props is scaled to match the
capillary pressure function in props.
This is a temporary workaround while the simulator uses two different
property object.
Previously, we had to use two layers of overlap cells such the
innermost layer contains the rightvalues automatically (as it is
surrounded by internal edges). No we use communication to get
the correct values in the whole overlap region and one layer
suffices as it should.
With this commit we add the possibility to start with a global representation
of a simulator that is read on each process and afterwards this presentation
is redistributed among the processors together with the properties and
state data needed to initialize the simulation.
There still is no parallel well handling and no parallel output. But with the
equilibrium example of @dr-robertk and deactivated output we can already
perform parallel runs.
In the parallel simulator we will have to be able adress only poperties on
some part of the global grid. To create thos properties we need to be able
to copy the grid independant data of the properties object and resize the rest.
This commit adds a construct taking a properties object for reading and a
new number of cells to accomplish this.
During the constructor the underlying object only holds smart
pointers and an empty vector. The FullyImplicitBlackoilSolver
obtains a reference to it from the NewtInterationInterface instances.
Therefore copying boost::any and storing it by value should be cheap
and safe.
We need it serveral places and all of them seem to have access to
NewtonIterationBlackoilInterface. This makes it natural to give access
to it and prevent users from having to forward it manually at several
places in the simulator driver.
As with opm-core we use boost::any to provide additional
information about a parallel run. It is used to set a
ParallelISTLInformation object and and fill it with the
information obtained from a parallel Cpgrid.
Note that the simulator currently compiles sucessfully. Still,
we have to test the runs and do debugging.
Apparently Eigen cannot handle empty containers during reserve correctly.
Therfore we check for the size of the vector and if it is zero simply create
empty jacobians.
Closes#290
The derivatives of the vappars are included in the Jacobian.
To avoid inf derivatives for vap<1, the oil saturation is restricted
from below by sqrt(epsilon).
Previously, we just called the default constructor. Of course this does not allocate
any memory and we experienced segmentation faults. With this patch we correctly
allocate the arrays by passing the number of row and columns to the constructor.
It took me quite some time to understand the computations done
e.g. during the detection of oscillations, where the stuff returned
by residuals() is used as a vector of doubles. It turns out that
residuals() actually returns the norm of the residuals. To clarify this
we rename residuals() to computeResidualNorms() and residuals to
residual_norms. Having my dare devil day today, I even try to document the
method. (This documented method might feel kind of lonely between the others,
now;). Hopefully this saves others some time.
The computations made to check the convergence are the same for all existing
phases. Therefore this patch uses loops over the phase indices when cmputing
them,
In the convergence check there are several reductions (maxCoeff(), sum())
that will trigger communication in a parallel run. This patch seperates the
reductions from the other computations. The idea is to one reduction for the
reductions that need to done as global communication is expensive.
The number of perforations may change due to completions beeing shut.
If the number of perforations changes the perfPhaseRates are set to
equal the wellRates/(number of perforations) instead of the values from
the previous time step.
The initial definition of the phase indices seems to be in
opm/core/props/BlackoilPhases.hpp. Nevertheless there were
several redefinitions of the same or similar enums (either
Aqua, Liquid, and Vapor, or Water, Oil, and Gas). Surprisingly
most often these definitions did not use the original values.
This is bound to break if there is a change upstream.
This patch limits the definition to one place in opm-autodiff,
namely opm/autodiff/BlackoilPropsAdInterface.hpp. To avoid
downstream confusion we define both the Water and Aqua triplets.
In addition we define the maximum number of phases to use at compile
time.
Shut is renamed to stopped in the wellsManager in order to better
reflect the name logic in the deck. This PR implements the nessesary
change in opm-autodiff
There has been an inconsitancy in which pressure to use in the
evaluation of the fluid properties.
With this commit the phase pressure is used for all the evaluation of
the fluid properties.
The multiplier contribution from the getRegionMultiplier is added to the
face multipliers. The getRegionMultiplier method is called with the cell
index on both side of the face in order to return the correct region
multiplier across faults.
1) Add the possibility for the user to chose between local and global
coordinate permeability in the transmissibility calculations.
2) Trow for CpGrid
3) Add default for switch
dynamical parameter given in the parameter file. The default is 0 (as before).
In addition the relaxation parameter has been added to the parameter with the
default preserving the state from before.
Also, the default parameter for use_amg and use_bicgstab in the constructor of
CPRPrecondition have been removed.
Note that this patch does not introduce any real temperature
dependence but only changes the APIs for the viscosity and for the
density related methods. Note that I also don't like the fact that
this requires so many changes to so many files, but with the current
design of the property classes I cannot see a way to avoid this...
This restores the performance to approximately the level it had before
the change to support non-diagonal well jacobians, for SPE9. All changes
are to the eliminateVariable() method.
- Explicitly compute and apply the inverse.
- Change loop ordering to apply inverse only num_eq - 1 times, instead
of (numeq - 1)^2.
- Use UmfPackLU instead of SparseLU.
1. The right hand side is solved only once
2. The solver is constructed directly with the matrix
3. const is added where it was missing
4. More commennts is added
5. Variable names are changed for clarification
The following comments has been adressed
1. An array is used in stead of pair
2. is not empty is used instead of size>0 to check if the well has been
initialized before
3. const_iterator is used instead of iterator
4. partial copy is removed
5. WellMapType is no longer mutable
If we need more than 150 linear iterations, it's probably something
wrong, and we may be better of by restarting with smaller time-steps.
TODO: make it possible for the user to specify this number.
The non-diagonal elements in the sub-matrices in the Schur complement is
no longer ignored. Instead of assuming the matrix do be diagonal, and
compute the invert of the sub-matrix, small linear systems are solved
using superLU.
Tested on SPE3 and Norne. (With this fix a slightly modified norne runs
until 3292 days)
The residual output is changed from max residual to
total mass balance residual. In this way the output
resemble what is actually used as convergence criteria.
With the introduction of shut wells the same ordering in the well states
can no longer be assumed. The well names is instead used to map the old
well state to the new. Also the partial_copying is moved into the
initialization.
Tested on SPE1, SPE3, SPE9 and Norne. (Do not change the SPEs and is
nessesary for the Norne)
The old rs and rv values are used in the phase transitions calculations.
The update of rs and rv to the state is therefore moved after the phase
transitions in order for the code to use the old values.
of the matrix internal allocators.
This fix also avoid the copying of the BCRSMatrix by providing a contructor that creates
the DuneMatrix for a given Eigen SparseMatrix.
Using a different construction method for the sparse matrix
used makes the methods a little faster. There is probably
still room for improvement by refactoring to avoid the
sparse matrix products.
New parameters are:
- cpr_use_amg (default false) if true, use AMG preconditioner for elliptic part
- cpr_use_bicgstab (default true) if true, use BiCGStab (else use CG) for elliptic part
The 'props_' table of PVT functions has one entry for each active
phase. Fix four instances of indexing into the table using the
canonical rather than active phase indices.
This is necessary, but not sufficient, to run two-phase problems
without a "dummy" third phase.
First the change in oil saturation is calculated from changes in water
and oil saturation. Then oil saturation is updated based on this change
instead of just fixed to 1-sw-sg. With this change the oil saturation is
less sensitive towards numerical errors that may cause very small oil
saturations. Witch again may cause the simulator to think that the gas
phase is saturation with vaporized oil when it is not.
Currently, there are two abstract interface for the grids. One that
usually returns pods and arrays of them that also can be used by C
and is used also in opm-core, and one that returns Eigen datastructures
needed within opm-autodiff.
This commit adds a postfix ToEigen to those functions (faceCells, and
cellCentroidsZ) one could imagine to also return pods and arrays of them.
This should at least resolve the confusion about the two faceCells functions.
The next step will be issue #192Fixes#176
Now it takes a std::vector instead of an Eigen::Array, more importantly
it expects one element per face, not interior face. The mapping now takes
place in setThresholdPressures().
