Mearge from upstream

This commit is contained in:
Halvor M. Nilsen 2012-06-20 13:11:46 +02:00
commit dedac2d493
4 changed files with 142 additions and 133 deletions

View File

@ -9,6 +9,10 @@ syntax: glob
*.pc *.pc
stamp-* stamp-*
.dirstamp .dirstamp
.autotools
.cproject
.project
.settings/*
Makefile Makefile
Makefile.in Makefile.in

View File

@ -26,21 +26,21 @@ destroy_grid(struct UnstructuredGrid *g)
{ {
if (g!=NULL) if (g!=NULL)
{ {
free(g->face_nodes); free(g->face_nodes);
free(g->face_nodepos); free(g->face_nodepos);
free(g->face_cells); free(g->face_cells);
free(g->cell_facepos); free(g->cell_facepos);
free(g->cell_faces); free(g->cell_faces);
free(g->node_coordinates); free(g->node_coordinates);
free(g->face_centroids); free(g->face_centroids);
free(g->face_areas); free(g->face_areas);
free(g->face_normals); free(g->face_normals);
free(g->cell_centroids); free(g->cell_centroids);
free(g->cell_volumes); free(g->cell_volumes);
free(g->global_cell); free(g->global_cell);
free(g->cell_facetag); free(g->cell_facetag);
} }
free(g); free(g);
@ -55,7 +55,7 @@ create_grid_empty(void)
G = malloc(1 * sizeof *G); G = malloc(1 * sizeof *G);
if (G != NULL) { if (G != NULL) {
*G = g; *G = g;
} }
return G; return G;
@ -64,11 +64,11 @@ create_grid_empty(void)
struct UnstructuredGrid * struct UnstructuredGrid *
allocate_grid(size_t ndims , allocate_grid(size_t ndims ,
size_t ncells , size_t ncells ,
size_t nfaces , size_t nfaces ,
size_t nfacenodes, size_t nfacenodes,
size_t ncellfaces, size_t ncellfaces,
size_t nnodes ) size_t nnodes )
{ {
size_t nel; size_t nel;
struct UnstructuredGrid *G; struct UnstructuredGrid *G;
@ -76,58 +76,58 @@ allocate_grid(size_t ndims ,
G = create_grid_empty(); G = create_grid_empty();
if (G != NULL) { if (G != NULL) {
/* Node fields ---------------------------------------- */ /* Node fields ---------------------------------------- */
nel = nnodes * ndims; nel = nnodes * ndims;
G->node_coordinates = malloc(nel * sizeof *G->node_coordinates); G->node_coordinates = malloc(nel * sizeof *G->node_coordinates);
/* Face fields ---------------------------------------- */ /* Face fields ---------------------------------------- */
nel = nfacenodes; nel = nfacenodes;
G->face_nodes = malloc(nel * sizeof *G->face_nodes); G->face_nodes = malloc(nel * sizeof *G->face_nodes);
nel = nfaces + 1; nel = nfaces + 1;
G->face_nodepos = malloc(nel * sizeof *G->face_nodepos); G->face_nodepos = malloc(nel * sizeof *G->face_nodepos);
nel = 2 * nfaces; nel = 2 * nfaces;
G->face_cells = malloc(nel * sizeof *G->face_cells); G->face_cells = malloc(nel * sizeof *G->face_cells);
nel = nfaces * ndims; nel = nfaces * ndims;
G->face_centroids = malloc(nel * sizeof *G->face_centroids); G->face_centroids = malloc(nel * sizeof *G->face_centroids);
nel = nfaces * ndims; nel = nfaces * ndims;
G->face_normals = malloc(nel * sizeof *G->face_normals); G->face_normals = malloc(nel * sizeof *G->face_normals);
nel = nfaces * 1; nel = nfaces * 1;
G->face_areas = malloc(nel * sizeof *G->face_areas); G->face_areas = malloc(nel * sizeof *G->face_areas);
/* Cell fields ---------------------------------------- */ /* Cell fields ---------------------------------------- */
nel = ncellfaces; nel = ncellfaces;
G->cell_faces = malloc(nel * sizeof *G->cell_faces); G->cell_faces = malloc(nel * sizeof *G->cell_faces);
nel = ncells + 1; nel = ncells + 1;
G->cell_facepos = malloc(nel * sizeof *G->cell_facepos); G->cell_facepos = malloc(nel * sizeof *G->cell_facepos);
nel = ncells * ndims; nel = ncells * ndims;
G->cell_centroids = malloc(nel * sizeof *G->cell_centroids); G->cell_centroids = malloc(nel * sizeof *G->cell_centroids);
nel = ncells * 1; nel = ncells * 1;
G->cell_volumes = malloc(nel * sizeof *G->cell_volumes); G->cell_volumes = malloc(nel * sizeof *G->cell_volumes);
if ((G->node_coordinates == NULL) || if ((G->node_coordinates == NULL) ||
(G->face_nodes == NULL) || (G->face_nodes == NULL) ||
(G->face_nodepos == NULL) || (G->face_nodepos == NULL) ||
(G->face_cells == NULL) || (G->face_cells == NULL) ||
(G->face_centroids == NULL) || (G->face_centroids == NULL) ||
(G->face_normals == NULL) || (G->face_normals == NULL) ||
(G->face_areas == NULL) || (G->face_areas == NULL) ||
(G->cell_faces == NULL) || (G->cell_faces == NULL) ||
(G->cell_facepos == NULL) || (G->cell_facepos == NULL) ||
(G->cell_centroids == NULL) || (G->cell_centroids == NULL) ||
(G->cell_volumes == NULL) ) (G->cell_volumes == NULL) )
{ {
destroy_grid(G); destroy_grid(G);
G = NULL; G = NULL;
} }
} }
return G; return G;

View File

@ -378,7 +378,7 @@ assemble_well_contrib(int nc ,
break; break;
case RESERVOIR_RATE: case RESERVOIR_RATE:
for (p = 0; p < np; ++p) { for (p = 0; p < np; p++) {
if (ctrls->distr[np * ctrls->current + p] != 1.0) { if (ctrls->distr[np * ctrls->current + p] != 1.0) {
*ok = 0; *ok = 0;
break; break;
@ -539,7 +539,7 @@ well_solution(const struct UnstructuredGrid *G ,
if (soln->well_press != NULL) { if (soln->well_press != NULL) {
/* Extract BHP directly from solution vector for non-shut wells */ /* Extract BHP directly from solution vector for non-shut wells */
for (w = 0; w < F->W->number_of_wells; ++w) { for (w = 0; w < F->W->number_of_wells; w++) {
if (F->W->ctrls[w]->current >= 0) { if (F->W->ctrls[w]->current >= 0) {
soln->well_press[w] = h->x[G->number_of_cells + w]; soln->well_press[w] = h->x[G->number_of_cells + w];
} }
@ -585,7 +585,8 @@ assemble_incompressible(struct UnstructuredGrid *G ,
int res_is_neumann, wells_are_rate; int res_is_neumann, wells_are_rate;
double s; double s;
*ok=1;
*ok = 1;
csrmatrix_zero( h->A); csrmatrix_zero( h->A);
vector_zero (h->A->m, h->b); vector_zero (h->A->m, h->b);
@ -694,6 +695,7 @@ ifs_tpfa_assemble(struct UnstructuredGrid *G ,
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
{ {
int system_singular, ok; int system_singular, ok;
assemble_incompressible(G, F, trans, gpress, h, &system_singular, &ok); assemble_incompressible(G, F, trans, gpress, h, &system_singular, &ok);
if (ok && system_singular) { if (ok && system_singular) {
@ -733,9 +735,11 @@ ifs_tpfa_assemble_comprock(struct UnstructuredGrid *G ,
* after it will always be nonsingular. * after it will always be nonsingular.
*/ */
if (ok) { if (ok) {
for (c = 0; c < G->number_of_cells; ++c) { for (c = 0; c < G->number_of_cells; c++) {
j = csrmatrix_elm_index(c, c, h->A); j = csrmatrix_elm_index(c, c, h->A);
d = porevol[c] * rock_comp[c] / dt; d = porevol[c] * rock_comp[c] / dt;
h->A->sa[j] += d; h->A->sa[j] += d;
h->b[c] += d * pressure[c]; h->b[c] += d * pressure[c];
} }
@ -761,7 +765,7 @@ ifs_tpfa_assemble_comprock_increment(struct UnstructuredGrid *G ,
{ {
int c, w, wdof, system_singular, ok; int c, w, wdof, system_singular, ok;
size_t j; size_t j;
double *v; double *v, dpvdt;
ok = 1; ok = 1;
assemble_incompressible(G, F, trans, gpress, h, &system_singular, &ok); assemble_incompressible(G, F, trans, gpress, h, &system_singular, &ok);
@ -775,17 +779,22 @@ ifs_tpfa_assemble_comprock_increment(struct UnstructuredGrid *G ,
v = h->pimpl->work; v = h->pimpl->work;
mult_csr_matrix(h->A, prev_pressure, v); mult_csr_matrix(h->A, prev_pressure, v);
for (c = 0; c < G->number_of_cells; ++c) { for (c = 0; c < G->number_of_cells; c++) {
j = csrmatrix_elm_index(c, c, h->A); j = csrmatrix_elm_index(c, c, h->A);
dpvdt = (porevol[c] - initial_porevolume[c]) / dt;
h->A->sa[j] += porevol[c] * rock_comp[c] / dt; h->A->sa[j] += porevol[c] * rock_comp[c] / dt;
h->b[c] += -(porevol[c] - initial_porevolume[c])/dt - v[c]; h->b[c] -= dpvdt + v[c];
}
if (F->W != NULL) {
wdof = G->number_of_cells;
for (w = 0; w < F->W->number_of_wells; w++, wdof++) {
h->b[wdof] -= v[wdof];
}
} }
if (F->W != NULL) {
for (w = 0; w < F->W->number_of_wells; ++w) {
wdof = G->number_of_cells + w;
h->b[wdof] += -v[wdof];
}
}
} }
return ok; return ok;

View File

@ -35,74 +35,73 @@
#ifndef OPENRS_FACTORY_HEADER #ifndef OPENRS_FACTORY_HEADER
#define OPENRS_FACTORY_HEADER #define OPENRS_FACTORY_HEADER
#include <map> #include <map>
#include <tr1/memory> #include <tr1/memory>
namespace Opm namespace Opm
{ {
/** This is an object factory for creating objects of some type /** This is an object factory for creating objects of some type
* requested by the user, with a shared base class. The user * requested by the user, with a shared base class. The user
* need only interact with the factory through the static * need only interact with the factory through the static
* template member addCreator() and the static member function * template member addCreator() and the static member function
* createObject(). * createObject().
*/ */
template <class Base> template<class Base>
class Factory class Factory
{ {
public: public:
/// The type of pointer returned by createObject(). /// The type of pointer returned by createObject().
typedef std::tr1::shared_ptr<Base> ProductPtr; typedef std::tr1::shared_ptr<Base> ProductPtr;
/// Creates a new object of the class associated with the given type string, /// Creates a new object of the class associated with the given type string,
/// and returns a pointer to it. /// and returns a pointer to it.
/// \param type the type string of the class that the user wants to have /// \param type the type string of the class that the user wants to have
/// constructed. /// constructed.
/// \return (smart) pointer to the created object. /// \return (smart) pointer to the created object.
static ProductPtr createObject(const std::string& type) static ProductPtr createObject(const std::string& type)
{ {
return instance().doCreateObject(type); return instance().doCreateObject(type);
} }
/// Clones an new object of the class associated with the given type string, /// Clones an new object of the class associated with the given type string,
/// and returns a pointer to it. /// and returns a pointer to it.
/// \param type the type string of the class that the user wants to have /// \param type the type string of the class that the user wants to have
/// constructed. /// constructed.
/// \param original (smart) pointer to object to be cloned. /// \param original (smart) pointer to object to be cloned.
/// \return (smart) pointer to the created object. /// \return (smart) pointer to the created object.
static ProductPtr cloneObject(const std::string& type, static ProductPtr cloneObject(const std::string& type,
const ProductPtr original) const ProductPtr original)
{ {
return instance().doCloneObject(type, original); return instance().doCloneObject(type, original);
} }
/// Add a creator to the Factory. /// Add a creator to the Factory.
/// After the call, the user may obtain new objects of the Derived type by /// After the call, the user may obtain new objects of the Derived type by
/// calling createObject() with the given type string as an argument. /// calling createObject() with the given type string as an argument.
/// \tparam Derived the class we want to add a creator for, must inherit /// \tparam Derived the class we want to add a creator for, must inherit
/// the class template parameter Base. /// the class template parameter Base.
/// \param type the type string with which we want the Factory to associate /// \param type the type string with which we want the Factory to associate
/// the class Derived. /// the class Derived.
template <class Derived> template<class Derived>
static void addCreator(const std::string& type) static void addCreator(const std::string& type)
{ {
instance().doAddCreator<Derived>(type); instance().doAddCreator<Derived>(type);
} }
private: private:
// The method that implements the singleton pattern, // The method that implements the singleton pattern,
// using the Meyers singleton technique. // using the Meyers singleton technique.
static Factory& instance() static Factory& instance()
{ {
static Factory singleton; static Factory singleton;
return singleton; return singleton;
} }
// Private constructor, to keep users from creating a Factory. // Private constructor, to keep users from creating a Factory.
Factory() Factory()
{ {
} }
// Abstract base class for Creators. // Abstract base class for Creators.
class Creator class Creator
@ -114,8 +113,8 @@ namespace Opm
}; };
/// This is the concrete Creator subclass for generating Derived objects. /// This is the concrete Creator subclass for generating Derived objects.
template <class Derived> template<class Derived>
class ConcreteCreator : public Creator class ConcreteCreator: public Creator
{ {
public: public:
virtual ProductPtr create() virtual ProductPtr create()
@ -132,45 +131,42 @@ namespace Opm
typedef std::tr1::shared_ptr<Creator> CreatorPtr; typedef std::tr1::shared_ptr<Creator> CreatorPtr;
typedef std::map<std::string, CreatorPtr> CreatorMap; typedef std::map<std::string, CreatorPtr> CreatorMap;
// This map contains the whole factory, i.e. all the Creators. // This map contains the whole factory, i.e. all the Creators.
CreatorMap string_to_creator_; CreatorMap string_to_creator_;
// Actually creates the product object. // Actually creates the product object.
ProductPtr doCreateObject(const std::string& type) ProductPtr doCreateObject(const std::string& type)
{ {
typename CreatorMap::iterator it; typename CreatorMap::iterator it;
it = string_to_creator_.find(type); it = string_to_creator_.find(type);
if (it == string_to_creator_.end()) { if (it == string_to_creator_.end()) {
THROW("Creator type " << type THROW("Creator type " << type
<< " is not registered in the factory."); << " is not registered in the factory.");
} }
return it->second->create(); return it->second->create();
} }
// Actually clones the product object. // Actually clones the product object.
ProductPtr doCloneObject(const std::string& type, ProductPtr doCloneObject(const std::string& type,
const ProductPtr original) const ProductPtr original)
{ {
typename CreatorMap::iterator it; typename CreatorMap::iterator it;
it = string_to_creator_.find(type); it = string_to_creator_.find(type);
if (it == string_to_creator_.end()) { if (it == string_to_creator_.end()) {
THROW("Creator type " << type THROW("Creator type " << type
<< " is not registered in the factory."); << " is not registered in the factory.");
} }
return it->second->clone(original); return it->second->clone(original);
} }
// Actually adds the creator. // Actually adds the creator.
template <class Derived> template<class Derived>
void doAddCreator(const std::string& type) void doAddCreator(const std::string& type)
{ {
CreatorPtr c(new ConcreteCreator<Derived>); CreatorPtr c(new ConcreteCreator<Derived>);
string_to_creator_[type] = c; string_to_creator_[type] = c;
} }
}; };
} // namespace Opm } // namespace Opm
#endif // OPENRS_FACTORY_HEADER #endif // OPENRS_FACTORY_HEADER