OilfieldToolsNet/opm-core
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge branch 'master' into ert

Browse Source 
Conflicts:
	Makefile.am
This commit is contained in:
Atgeirr Flø Rasmussen 2012-10-05 09:20:04 +02:00
commit 2c8b47b143
14 changed files with 516 additions and 84 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

15
Makefile.am
View File

@ -19,13 +19,16 @@ $(OPM_BOOST_CPPFLAGS)
# Link-time flags needed both to successfully link the library and to
# (transitively) convey inter-library dependency information.
lib_libopmcore_la_LDFLAGS = \
$(ERT_LDFLAGS) \
lib_libopmcore_la_LDFLAGS = \
-R $(OPM_BOOST_LIBDIR) \
$(OPM_BOOST_LDFLAGS) \
$(BOOST_FILESYSTEM_LIB) \
$(BOOST_SYSTEM_LIB) \
$(BOOST_DATE_TIME_LIB) \
$(BOOST_UNIT_TEST_FRAMEWORK_LIB) \
$(ERT_LDFLAGS)
lib_libopmcore_la_LIBADD = \
$(BOOST_FILESYSTEM_LIB) \
$(BOOST_SYSTEM_LIB) \
$(BOOST_DATE_TIME_LIB) \
$(BOOST_UNIT_TEST_FRAMEWORK_LIB) \
$(ERT_LIBS) \
$(LAPACK_LIBS) $(BLAS_LIBS) $(LIBS)

4
examples/Makefile.am
View File

@ -8,9 +8,7 @@ $(OPM_BOOST_CPPFLAGS) \
# All targets link to the library
LDADD = \
$(top_builddir)/lib/libopmcore.la \
$(BOOST_FILESYSTEM_LIB) \
$(BOOST_SYSTEM_LIB)
$(top_builddir)/lib/libopmcore.la
# ----------------------------------------------------------------------
# Declare products (i.e., the example programs).

43
m4/ax_boost_unit_test_framework.m4
View File

@ -102,27 +102,32 @@ AC_DEFUN([AX_BOOST_UNIT_TEST_FRAMEWORK],
fi
else
link_unit_test_framework="no"
saved_ldflags="${LDFLAGS}"
saved_ldflags="${LDFLAGS}"
for ax_lib in boost_unit_test_framework-$ax_boost_user_unit_test_framework_lib $ax_boost_user_unit_test_framework_lib ; do
if test "x$link_unit_test_framework" = "xyes"; then
break;
fi
for unittest_library in `ls $BOOSTLIBDIR/lib${ax_lib}.so* $BOOSTLIBDIR/lib${ax_lib}.a* 2>/dev/null` ; do
if test -r $unittest_library ; then
libextension=`echo $unittest_library | sed 's,.*/,,' | sed -e 's;^lib\(boost_unit_test_framework.*\)\.so.*$;\1;' -e 's;^lib\(boost_unit_test_framework.*\)\.a*$;\1;'`
ax_lib=${libextension}
link_unit_test_framework="yes"
else
link_unit_test_framework="no"
fi
#if test "x$link_unit_test_framework" = "xyes"; then
# break;
#fi
#for unittest_library in `ls $BOOSTLIBDIR/lib${ax_lib}.so* $BOOSTLIBDIR/lib${ax_lib}.a* 2>/dev/null` ; do
#if test -r $unittest_library ; then
# libextension=`echo $unittest_library | sed 's,.*/,,' | sed -e 's;^lib\(boost_unit_test_framework.*\)\.so.*$;\1;' -e 's;^lib\(boost_unit_test_framework.*\)\.a*$;\1;'`
# ax_lib=${libextension}
# link_unit_test_framework="yes"
# else
# link_unit_test_framework="no"
# fi
#
# if test "x$link_unit_test_framework" = "xyes"; then
# BOOST_UNIT_TEST_FRAMEWORK_LIB="-l$ax_lib"
# AC_SUBST(BOOST_UNIT_TEST_FRAMEWORK_LIB)
# break
# fi
if test "x$link_unit_test_framework" = "xyes"; then
BOOST_UNIT_TEST_FRAMEWORK_LIB="-l$ax_lib"
AC_SUBST(BOOST_UNIT_TEST_FRAMEWORK_LIB)
break
fi
done
done
# done
AC_CHECK_LIB($ax_lib, [exit],
[BOOST_UNIT_TEST_FRAMEWORK_LIB="-l$ax_lib"; AC_SUBST([BOOST_UNIT_TEST_FRAMEWORK_LIB]) link_unit_test_framework="yes"; break],
[link_unit_test_framework="no"])
done
fi
if test "x$ax_lib" = "x"; then
AC_MSG_ERROR(Could not find a version of the library!)

6
m4/opm_boost_base.m4
View File

@ -103,6 +103,7 @@ if test "x$want_boost" = "xyes"; then
for ac_boost_path_tmp in $libsubdirs; do
if test -d "$ac_boost_path"/"$ac_boost_path_tmp" ; then
OPM_BOOST_LDFLAGS="-L$ac_boost_path/$ac_boost_path_tmp"
OPM_BOOST_LIBDIR="${ac_boost_path}/${ac_boost_path_tmp}"
break
fi
done
@ -114,6 +115,7 @@ if test "x$want_boost" = "xyes"; then
done
OPM_BOOST_LDFLAGS="-L$ac_boost_path_tmp/$libsubdir"
OPM_BOOST_CPPFLAGS="-I$ac_boost_path_tmp/include"
OPM_BOOST_LIBDIR="${ac_boost_path_tmp}/${libsubdir}"
break;
fi
done
@ -123,6 +125,7 @@ if test "x$want_boost" = "xyes"; then
dnl --with-boost-libdir parameter
if test "$ac_boost_lib_path" != ""; then
OPM_BOOST_LDFLAGS="-L$ac_boost_lib_path"
OPM_BOOST_LIBDIR="${ac_boost_lib_path}"
fi
CPPFLAGS_SAVED="$CPPFLAGS"
@ -191,6 +194,7 @@ if test "x$want_boost" = "xyes"; then
if ls "$best_path/$libsubdir/libboost_"* >/dev/null 2>&1 ; then break; fi
done
OPM_BOOST_LDFLAGS="-L$best_path/$libsubdir"
OPM_BOOST_LIBDIR="$best_path/$libsubdir"
fi
fi
@ -207,6 +211,7 @@ if test "x$want_boost" = "xyes"; then
AC_MSG_NOTICE(We will use a staged boost library from $BOOST_ROOT)
OPM_BOOST_CPPFLAGS="-I$BOOST_ROOT"
OPM_BOOST_LDFLAGS="-L$BOOST_ROOT/stage/$libsubdir"
OPM_BOOST_LIBDIR="$BOOST_ROOT/stage/$libsubdir"
fi
fi
fi
@ -248,6 +253,7 @@ if test "x$want_boost" = "xyes"; then
else
AC_SUBST([OPM_BOOST_CPPFLAGS])
AC_SUBST([OPM_BOOST_LDFLAGS])
AC_SUBST([OPM_BOOST_LIBDIR])
AC_DEFINE([OPM_HAVE_BOOST], [1], [define if the Boost library is available])
# execute ACTION-IF-FOUND (if present):
ifelse([$2], , :, [$2])

398
opm/core/grid.c
View File

@ -18,7 +18,12 @@
*/
#include <opm/core/grid.h>
#include <assert.h>
#include <errno.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
void
@ -135,3 +140,396 @@ allocate_grid(size_t ndims ,
return G;
}
#define GRID_NMETA 6
#define GRID_NDIMS 0
#define GRID_NCELLS 1
#define GRID_NFACES 2
#define GRID_NNODES 3
#define GRID_NFACENODES 4
#define GRID_NCELLFACES 5
static void
input_error(FILE *fp, const char * const err)
{
int save_errno = errno;
if (ferror(fp)) {
fprintf(stderr, "%s: %s\n", err, strerror(save_errno));
clearerr(fp);
}
else if (feof(fp)) {
fprintf(stderr, "%s: End-of-file\n", err);
}
errno = save_errno;
}
static struct UnstructuredGrid *
allocate_grid_from_file(FILE *fp, int *has_tag, int *has_indexmap)
{
struct UnstructuredGrid *G;
int save_errno;
unsigned long tmp;
size_t dimens[GRID_NMETA], i;
save_errno = errno;
i = 0;
while ((i < GRID_NMETA) && (fscanf(fp, " %lu", &tmp) == 1)) {
dimens[i] = tmp;
i += 1;
}
if (i == GRID_NMETA) {
if (fscanf(fp, "%d %d", has_tag, has_indexmap) == 2) {
G = allocate_grid(dimens[GRID_NDIMS] ,
dimens[GRID_NCELLS] ,
dimens[GRID_NFACES] ,
dimens[GRID_NFACENODES],
dimens[GRID_NCELLFACES],
dimens[GRID_NNODES] );
if (G != NULL) {
if (! *has_tag) {
free(G->cell_facetag);
G->cell_facetag = NULL;
}
if (*has_indexmap) {
G->global_cell =
malloc(dimens[GRID_NCELLS] * sizeof *G->global_cell);
/* Allocation failure checked elsewhere. */
}
G->number_of_cells = (int) dimens[GRID_NCELLS];
G->number_of_faces = (int) dimens[GRID_NFACES];
G->number_of_nodes = (int) dimens[GRID_NNODES];
G->dimensions = (int) dimens[GRID_NDIMS];
i = 0;
while ((i < dimens[GRID_NDIMS]) &&
(fscanf(fp, "%d", & G->cartdims[ i ]) == 1)) {
i += 1;
}
if (i < dimens[GRID_NDIMS]) {
input_error(fp, "Unable to read Cartesian dimensions");
destroy_grid(G);
G = NULL;
}
else {
/* Account for dimens[GRID_DIMS] < 3 */
size_t n = (sizeof G->cartdims) / (sizeof G->cartdims[0]);
for (; i < n; i++) { G->cartdims[ i ] = 1; }
}
}
}
else {
input_error(fp, "Unable to read grid predicates");
G = NULL;
}
}
else {
input_error(fp, "Unable to read grid dimensions");
G = NULL;
}
errno = save_errno;
return G;
}
static int
read_grid_nodes(FILE *fp, struct UnstructuredGrid *G)
{
int save_errno;
size_t i, n;
save_errno = errno;
n = G->dimensions;
n *= G->number_of_nodes;
i = 0;
while ((i < n) &&
(fscanf(fp, " %lf", & G->node_coordinates[ i ]) == 1)) {
i += 1;
}
if (i < n) {
input_error(fp, "Unable to read node coordinates");
}
errno = save_errno;
return i == n;
}
static int
read_grid_faces(FILE *fp, struct UnstructuredGrid *G)
{
int save_errno, ok;
size_t nf, nfn, i;
save_errno = errno;
nf = G->number_of_faces;
/* G->face_nodepos */
i = 0;
while ((i < nf + 1) &&
(fscanf(fp, " %d", & G->face_nodepos[ i ]) == 1)) {
i += 1;
}
ok = i == nf + 1;
if (! ok) {
input_error(fp, "Unable to read node indirection array");
}
else {
/* G->face_nodes */
nfn = G->face_nodepos[ nf ];
i = 0;
while ((i < nfn) && (fscanf(fp, " %d", & G->face_nodes[ i ]) == 1)) {
i += 1;
}
ok = i == nfn;
if (! ok) {
input_error(fp, "Unable to read face-nodes");
}
}
if (ok) {
/* G->face_cells */
i = 0;
while ((i < 2 * nf) && (fscanf(fp, " %d", & G->face_cells[ i ]) == 1)) {
i += 1;
}
ok = i == 2 * nf;
if (! ok) {
input_error(fp, "Unable to read neighbourship");
}
}
if (ok) {
/* G->face_areas */
i = 0;
while ((i < nf) && (fscanf(fp, " %lf", & G->face_areas[ i ]) == 1)) {
i += 1;
}
ok = i == nf;
if (! ok) {
input_error(fp, "Unable to read face areas");
}
}
if (ok) {
/* G->face_centroids */
size_t n;
n = G->dimensions;
n *= nf;
i = 0;
while ((i < n) && (fscanf(fp, " %lf", & G->face_centroids[ i ]) == 1)) {
i += 1;
}
ok = i == n;
if (! ok) {
input_error(fp, "Unable to read face centroids");
}
}
if (ok) {
/* G->face_normals */
size_t n;
n = G->dimensions;
n *= nf;
i = 0;
while ((i < n) && (fscanf(fp, " %lf", & G->face_normals[ i ]) == 1)) {
i += 1;
}
ok = i == n;
if (! ok) {
input_error(fp, "Unable to read face normals");
}
}
errno = save_errno;
return ok;
}
static int
read_grid_cells(FILE *fp, int has_tag, int has_indexmap,
struct UnstructuredGrid *G)
{
int save_errno, ok;
size_t nc, ncf, i;
save_errno = errno;
nc = G->number_of_cells;
/* G->cell_facepos */
i = 0;
while ((i < nc + 1) && (fscanf(fp, " %d", & G->cell_facepos[ i ]) == 1)) {
i += 1;
}
ok = i == nc + 1;
if (! ok) {
input_error(fp, "Unable to read face indirection array");
}
else {
/* G->cell_faces (and G->cell_facetag if applicable) */
ncf = G->cell_facepos[ nc ];
i = 0;
if (has_tag) {
assert (G->cell_facetag != NULL);
while ((i < ncf) &&
(fscanf(fp, " %d %d",
& G->cell_faces [ i ],
& G->cell_facetag[ i ]) == 2)) {
i += 1;
}
}
else {
while ((i < ncf) &&
(fscanf(fp, " %d", & G->cell_faces[ i ]) == 1)) {
i += 1;
}
}
ok = i == ncf;
if (! ok) {
input_error(fp, "Unable to read cell-faces");
}
}
if (ok) {
/* G->global_cell if applicable */
if (has_indexmap) {
i = 0;
if (G->global_cell != NULL) {
while ((i < nc) &&
(fscanf(fp, " %d", & G->global_cell[ i ]) == 1)) {
i += 1;
}
}
else {
int discard;
while ((i < nc) && (fscanf(fp, " %d", & discard) == 1)) {
i += 1;
}
}
}
else {
assert (G->global_cell == NULL);
i = nc;
}
ok = i == nc;
if (! ok) {
input_error(fp, "Unable to read global cellmap");
}
}
if (ok) {
/* G->cell_volumes */
i = 0;
while ((i < nc) && (fscanf(fp, " %lf", & G->cell_volumes[ i ]) == 1)) {
i += 1;
}
ok = i == nc;
if (! ok) {
input_error(fp, "Unable to read cell volumes");
}
}
if (ok) {
/* G->cell_centroids */
size_t n;
n = G->dimensions;
n *= nc;
i = 0;
while ((i < n) && (fscanf(fp, " %lf", & G->cell_centroids[ i ]) == 1)) {
i += 1;
}
ok = i == n;
if (! ok) {
input_error(fp, "Unable to read cell centroids");
}
}
errno = save_errno;
return ok;
}
struct UnstructuredGrid *
read_grid(const char *fname)
{
struct UnstructuredGrid *G;
FILE *fp;
int save_errno;
int has_tag, has_indexmap, ok;
save_errno = errno;
fp = fopen(fname, "rt");
if (fp != NULL) {
G = allocate_grid_from_file(fp, & has_tag, & has_indexmap);
ok = G != NULL;
if (ok) { ok = read_grid_nodes(fp, G); }
if (ok) { ok = read_grid_faces(fp, G); }
if (ok) { ok = read_grid_cells(fp, has_tag, has_indexmap, G); }
if (! ok) {
destroy_grid(G);
G = NULL;
}
fclose(fp);
}
else {
G = NULL;
}
errno = save_errno;
return G;
}

11
opm/core/grid.h
View File

@ -273,6 +273,17 @@ allocate_grid(size_t ndims ,
size_t ncellfaces,
size_t nnodes );
/**
* Import a grid from a character representation stored in file.
*
* @param[in] fname File name.
* @return Fully formed UnstructuredGrid with all fields allocated and filled.
* Returns @c NULL in case of allocation failure.
*/
struct UnstructuredGrid *
read_grid(const char *fname);
#ifdef __cplusplus
}
#endif

6
opm/core/transport/GravityColumnSolver_impl.hpp
View File

@ -175,10 +175,10 @@ namespace Opm
}
// model_.sourceTerms(); // Not needed
// Solve.
const int num_rhs = 1;
int info = 0;
const MAT_SIZE_T num_rhs = 1, colSize = col_size;
MAT_SIZE_T info = 0;
// Solution will be written to rhs.
dgtsv_(&col_size, &num_rhs, DL, D, DU, &rhs[0], &col_size, &info);
dgtsv_(&colSize, &num_rhs, DL, D, DU, &rhs[0], &colSize, &info);
if (info != 0) {
THROW("Lapack reported error in dgtsv: " << info);
}

2
opm/core/transport/reorder/TransportModelInterface.cpp
View File

@ -35,6 +35,7 @@ void Opm::TransportModelInterface::reorderAndTransport(const UnstructuredGrid& g
// Invoke appropriate solve method for each interdependent component.
for (int comp = 0; comp < ncomponents; ++comp) {
#if 0
#ifdef MATLAB_MEX_FILE
// \TODO replace this with general signal handling code, check if it costs performance.
if (interrupt_signal) {
@ -42,6 +43,7 @@ void Opm::TransportModelInterface::reorderAndTransport(const UnstructuredGrid& g
"cells finished.\n", i, grid.number_of_cells);
break;
}
#endif
#endif
const int comp_size = components[comp + 1] - components[comp];
if (comp_size == 1) {

26
opm/core/transport/reorder/nlsolvers.c
View File

@ -26,24 +26,8 @@ SOFTWARE.
#include <assert.h>
#include <math.h>
#ifdef MATLAB_MEX_FILE
#include "nlsolvers.h"
#include <mex.h>
extern int interrupt_signal;
#define print mexPrintf
#define malloc mxMalloc
#define calloc mxCalloc
#define realloc mxRealloc
#define free mxFree
#else /* MATLAB_MEX_FILE */
#define print printf
#include <opm/core/transport/reorder/nlsolvers.h>
#endif /* MATLAB_MEX_FILE */
static const char no_root_str[]=
" In %s:\n"
" With G(%5f) =% 5f, G(%5f) =% 5f, G(x) is not bracketed!\n";
@ -100,7 +84,7 @@ ridders (double (*G)(double, void*), void *data, struct NonlinearSolverCtrl *ctr
if (G0*G1 > 0)
{
print(no_root_str, "ridder", s0, G0, s1, G1);
printf(no_root_str, "ridder", s0, G0, s1, G1);
return -1.0;
}
@ -157,7 +141,7 @@ ridders (double (*G)(double, void*), void *data, struct NonlinearSolverCtrl *ctr
}
else
{
print("In ridder:\nG0=%10.10f, G1=%10.10f, "
printf("In ridder:\nG0=%10.10f, G1=%10.10f, "
"G3=%10.10f\n", G0, G1, G3);
}
s2 = 0.5*(s0+s1);
@ -199,7 +183,7 @@ regulafalsi (double (*G)(double, void*), void *data, struct NonlinearSolverCtrl
if (G0*G1 > 0)
{
print(no_root_str, "regulafalsi", s0, G0, s1, G1);
printf(no_root_str, "regulafalsi", s0, G0, s1, G1);
return -1.0;
}
@ -282,7 +266,7 @@ bisection (double (*G)(double, void*), void *data, struct NonlinearSolverCtrl *c
if (G0*G1 > 0.0)
{
print(no_root_str, "bisection", s0, G0, s1, G1);
printf(no_root_str, "bisection", s0, G0, s1, G1);
return -1.0;
}
@ -312,7 +296,7 @@ bisection (double (*G)(double, void*), void *data, struct NonlinearSolverCtrl *c
}
if (ctrl->iterations >= ctrl->maxiterations)
{
print("Warning: convergence criterion not met\n");
printf("Warning: convergence criterion not met\n");
}
ctrl->residual = Gn;
return sn;

4
opm/core/transport/reorder/reordersequence.cpp
View File

@ -1,11 +1,11 @@
/* Copyright 2011 (c) Jostein R. Natvig <Jostein.R.Natvig at sintef.no> */
#include <opm/core/grid.h>
#ifdef MATLAB_MEX_FILE
#include "grid.h"
#include "reordersequence.h"
#include "tarjan.h"
#else
#include <opm/core/grid.h>
#include <opm/core/transport/reorder/reordersequence.h>
#include <opm/core/transport/reorder/tarjan.h>
#endif

5
opm/core/utility/parameters/ParameterGroup.cpp
View File

@ -4,7 +4,7 @@
//
// Created: Tue Jun 2 19:13:17 2009
//
// Author(s): Bård Skaflestad <bard.skaflestad@sintef.no>
// Author(s): Bård Skaflestad <bard.skaflestad@sintef.no>
// Atgeirr F Rasmussen <atgeirr@sintef.no>
//
// $Date$
@ -139,9 +139,6 @@ namespace Opm {
}
}
}
#ifdef MATLAB_MEX_FILE
fclose(is);
#endif
}
void ParameterGroup::writeParam(const std::string& param_filename) const {

1
tests/Makefile.am
View File

@ -20,6 +20,7 @@ sparsevector_test \
test_cartgrid \
test_column_extract \
test_lapack \
test_read_grid \
test_read_vag \
test_readpolymer \
test_sf2p \

54
tests/test_lapack.cpp
View File

@ -24,34 +24,36 @@
namespace {
struct BandMatrixCoeff
{
BandMatrixCoeff(int N, int ku, int kl) : ku_(ku), kl_(kl), nrow_(2*kl + ku + 1), N_(N) {
BandMatrixCoeff(MAT_SIZE_T N, MAT_SIZE_T ku, MAT_SIZE_T kl)
: ku_(ku), kl_(kl), nrow_(2 * kl + ku + 1), N_(N)
{
}
// compute the position where to store the coefficient of a matrix A_{i,j} (i,j=0,...,N-1)
// in a array which is sent to the band matrix solver of LAPACK.
int operator ()(int i, int j) {
int operator ()(MAT_SIZE_T i, MAT_SIZE_T j) {
return kl_ + ku_ + i - j + j*nrow_;
}
const int ku_;
const int kl_;
const int nrow_;
const int N_;
const MAT_SIZE_T ku_;
const MAT_SIZE_T kl_;
const MAT_SIZE_T nrow_;
const MAT_SIZE_T N_;
};
} //end anonymous namespace
int main()
{
const int N = 5;
const int nrhs = 1;
const MAT_SIZE_T N = 5;
const MAT_SIZE_T nrhs = 1;
double DU[N-1] = { 2.1, -1.0, 1.9, 8.0 };
double D[N] = { 3.0, 2.3, -5.0, -0.9, 7.1 };
double DL[N-1] = { 3.4, 3.6, 7.0, -6.0 };
double B[N] = { 2.7, -0.5, 2.6, 0.6, 2.7 };
// double B[N] = { 2.7, -0.5, 2.6, 0.6, 2.7 };
int info = 0;
MAT_SIZE_T info = 0;
dgtsv_(&N, &nrhs, DL, D, DU, B, &N, &info);
if (info == 0) {
for (int i = 0; i < N; ++i) {
@ -64,12 +66,12 @@ int main()
std::cout << std::endl;
//test of dgbsv_
int ldb = N;
int lda = N;
std::vector<int> ipiv(N, 0);
MAT_SIZE_T ldb = N;
MAT_SIZE_T lda = N;
std::vector<MAT_SIZE_T> ipiv(N, 0);
std::vector<double> AA(N*N, 0.);
std::vector<double> BB(N, 0.);
for (int i = 0; i < N; ++i) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
AA[i + i*N] = 10.;
if (i > 0) {
AA[i + (i - 1)*N] = i;
@ -78,31 +80,31 @@ int main()
BB[i] = i;
}
for (int i = 0; i < N; ++i) {
for (int j = 0; j < N; ++j) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
for (MAT_SIZE_T j = 0; j < N; ++j) {
std::cout << " " << AA[i + j*N];
}
std::cout << " " << std::endl;
}
std::cout << std::endl;
int kl = 1;
int ku = 1;
int nrowAB = 2*kl + ku + 1;
int ldab = nrowAB;
MAT_SIZE_T kl = 1;
MAT_SIZE_T ku = 1;
MAT_SIZE_T nrowAB = 2*kl + ku + 1;
MAT_SIZE_T ldab = nrowAB;
std::vector<double> AB(nrowAB*N, 0.);
BandMatrixCoeff bmc(N, ku, kl);
// Rewrite AA matrix in band format AB
for (int i = 0; i < N; ++i) {
for (int j = -1; j < 2; ++j) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
for (MAT_SIZE_T j = -1; j < 2; ++j) {
if (i + j > -1 && i + j < N)
AB[bmc(i, i + j)] = AA[i + N*(i + j)];
}
}
for (int i = 0; i < nrowAB; ++i) {
for (int j = 0; j < N; ++j) {
for (MAT_SIZE_T i = 0; i < nrowAB; ++i) {
for (MAT_SIZE_T j = 0; j < N; ++j) {
std::cout << " " << AB[i + j*nrowAB];
}
std::cout << " " << std::endl;
@ -113,7 +115,7 @@ int main()
dgesv_(&N ,&nrhs, &AA[0], &lda, &ipiv[0], &BB[0], &ldb, &info);
if (info == 0) {
for (int i = 0; i < N; ++i) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
std::cout << BB[i] << ' ';
}
std::cout << std::endl;
@ -121,14 +123,14 @@ int main()
std::cerr << "Something went wrong in debsv_()\n";
}
for (int i = 0; i < N; ++i) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
BB[i] = i;
}
dgbsv_(&N, &kl, &ku, &nrhs, &AB[0], &ldab, &ipiv[0], &BB[0], &ldb, &info);
if (info == 0) {
for (int i = 0; i < N; ++i) {
for (MAT_SIZE_T i = 0; i < N; ++i) {
std::cout << BB[i] << ' ';
}
std::cout << std::endl;

25
tests/test_read_grid.c Normal file
View File

@ -0,0 +1,25 @@
/*
* test_read_grid.c
*
* Created on: Aug 28, 2012
* Author: bska
*/
#include <stdio.h>
#include <opm/core/grid.h>
#include <opm/core/grid/cart_grid.h>
int
main(void)
{
struct UnstructuredGrid *G1, *G2;
G1 = read_grid("cart_grid_2d.txt");
G2 = create_grid_cart2d(2, 2);
destroy_grid(G2);
destroy_grid(G1);
return 0;
}