Commit 2bd5c535 authored by incardon's avatar incardon

Last changes for numerics

parent 1a6ed291
......@@ -13,6 +13,7 @@ m4_ifdef([AX_BOOST],,[m4_include([m4/ax_boost.m4])])
m4_ifdef([ACX_MPI],,[m4_include([m4/acx_mpi.m4])])
m4_ifdef([AX_OPENMP],,[m4_include([m4/ax_openmp.m4])])
m4_ifdef([AX_CUDA],,[m4_include([m4/ax_cuda.m4])])
m4_ifdef([IMMDX_LIB_METIS],,[m4_include([m4/immdx_lib_metis.m4])])
CXXFLAGS+=" --std=c++11 -march=native -mtune=native -Wno-unused-local-typedefs -Wextra -Wno-unused-parameter "
NVCCFLAGS=" "
......@@ -80,6 +81,13 @@ else
NVCCFLAGS+="$NVCCFLAGS -O3 "
fi
#########
## Check for Metis
IMMDX_LIB_METIS([],[echo "Cannot detect metis, use the --with-metis option if it is not installed in the default location"
exit 201])
####### include openfpm_devices include path
INCLUDES_PATH+=" -I/usr/local/include -I. -Iconfig -I../../openfpm_devices/src -I../../openfpm_data/src -I../../openfpm_io/src -I../../openfpm_vcluster/src -I../../src"
......
......@@ -31,7 +31,7 @@ class D
* \tparam ord
*
*/
inline static std::unordered_map<long int,typename Sys_eqs::stype> value(grid_key_dx<Sys_eqs::dims> & pos, const grid_sm<Sys_eqs::dims,void> & gs)
inline static std::unordered_map<long int,typename Sys_eqs::stype> value(grid_key_dx<Sys_eqs::dims> & pos, const grid_sm<Sys_eqs::dims,void> & gs, std::unordered_map<long int,typename Sys_eqs::stype > & cols, typename Sys_eqs::stype coeff)
{
std::cerr << "Error " << __FILE__ << ":" << __LINE__ << " only CENTRAL, FORWARD, BACKWARD derivative are defined";
}
......
......@@ -9,6 +9,7 @@
#define OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_FDSCHEME_HPP_
#include "Grid/grid_dist_id.hpp"
#include "Grid/grid_dist_id_iterator_sub.hpp"
#include "Matrix/Matrix.hpp"
#include "eq.hpp"
......@@ -23,12 +24,14 @@ class FDScheme
{
openfpm::vector<triplet<typename Sys_eqs::stype>> trpl;
public:
/*! \brief Impose an operator
*
*
*
*/
/* template<typename T> void impose(T & op, grid_key_dx_dist_iterator_sub & it)
template<typename T> void impose(T & op, const grid_sm<Sys_eqs::dims,void> & gs , grid_dist_iterator_sub<Sys_eqs::dims,typename Sys_eqs::b_grid::d_grid> it)
{
size_t lin = 0;
......@@ -41,10 +44,10 @@ class FDScheme
auto key = it.get();
// convert into global coordinate the position
auto keyg = g.getGKey(key);
auto keyg = it.getGKey(key);
// Convert local key into global key
T::value(op,cols);
T::value(keyg,gs,cols,1.0);
// create the triplet
......@@ -60,17 +63,33 @@ class FDScheme
cols.clear();
++loc_cnt;
++it;
}
}*/
}
/*! \brief produce the Matrix
*
* \tparam Syst_eq System of equations, or a single equation
*
*/
template<typename Grid> SparseMatrix<typename Sys_eqs::stype> getMatrix(const Grid & g)
{
#ifdef SE_CLASS1
if (Sys_eqs::num_cfields)
std::cerr << __FILE__ << ":" << __LINE__ << " if you do not provide ConstantFields in getMatrix, the system should not contain any constant field, while" << Sys_eqs::num_cfields << "\n";
#endif
const ConstField fld[Sys_eqs::num_cfields];
getMatrix(g,fld);
}
/*! \brief produce the Matrix
*
* \tparam Syst_eq System of equations, or a single equation
*
*/
template<typename Grid> SparseMatrix<typename Sys_eqs::stype> getMatrix(const Grid & g, const ConstField(& fld)[Sys_eqs::num_cfields::value] )
template<typename Grid> SparseMatrix<typename Sys_eqs::stype> getMatrix(const Grid & g, const ConstField(& fld)[Sys_eqs::num_cfields] )
{
// iterate across the full domain
......
......@@ -11,6 +11,9 @@
#include "Laplacian.hpp"
#include "FiniteDifference/eq.hpp"
#include "FiniteDifference/sum.hpp"
#include "Grid/grid_dist_id.hpp"
#include "data_type/scalar.hpp"
#include "Decomposition/CartDecomposition.hpp"
// Stokes flow
......@@ -25,8 +28,14 @@ struct lid_nn
// boundary at X and Y
static const bool boundary[];
//
static constexpr unsigned int num_cfields = 0;
// type of space float, double, ...
typedef float stype;
// type of base grid
typedef grid_dist_id<2,float,scalar<float>,CartDecomposition<2,float>> b_grid;
};
const bool lid_nn::boundary[] = {NON_PERIODIC,NON_PERIODIC};
......@@ -50,7 +59,7 @@ typedef Field<P,lid_nn> Prs;
typedef mul<eta,Lap<v_x,lid_nn>,lid_nn> eta_lap_vx;
typedef D<x,Prs,lid_nn> p_x;
typedef sum<eta_lap_vx,p_x> vx_eq;
typedef sum<eta_lap_vx,p_x,lid_nn> vx_eq;
// Eq2 V_y
......@@ -68,8 +77,34 @@ typedef sum<dx_vx,dy_vy> incompressibility;
BOOST_AUTO_TEST_CASE( lid_driven_cavity )
{
// Domain
Box<2,float> domain({0.0,0.0},{1.0,1.0});
}
// Ghost
Ghost<2,float> g(0.01);
size_t sz[] = {32,32};
// Initialize the global VCluster
init_global_v_cluster(&boost::unit_test::framework::master_test_suite().argc,&boost::unit_test::framework::master_test_suite().argv);
// Initialize openfpm
grid_dist_id<2,float,scalar<float>,CartDecomposition<2,float>> g_dist(sz,domain,g);
// Distributed grid
FDScheme<lid_nn> fd;
// start and end of the bulk
grid_key_dx<2> bulk_start(1,1);
grid_key_dx<2> bulk_end(sz[0]-1,sz[1]-1);
// Impose the vx equation
vx_eq vx;
fd.impose(vx, g_dist.getGridInfo() , g_dist.getSubDomainIterator(bulk_start,bulk_end));
// Impose the vy equation
vy_eq vy;
fd.impose(vy, g_dist.getGridInfo(), g_dist.getSubDomainIterator(bulk_start,bulk_end));
}
#endif /* OPENFPM_NUMERICS_SRC_FINITEDIFFERENCE_EQ_UNIT_TEST_HPP_ */
......@@ -17,13 +17,36 @@
*
*/
template<typename ... expr>
class sum
struct sum
{
// Transform from variadic template to boost mpl vector
typedef boost::mpl::vector<expr... > v_expr;
// size of v_expr
typedef typename boost::mpl::size<v_expr>::type v_sz;
typedef typename boost::mpl::at<v_expr, boost::mpl::int_<v_sz::type::value - 1> >::type Sys_eqs;
/*! \brief Create the row of the Matrix
*
* \tparam ord
*
*/
inline static std::unordered_map<long int,typename Sys_eqs::stype> value(grid_key_dx<Sys_eqs::dims> & pos, const grid_sm<Sys_eqs::dims,void> & gs, std::unordered_map<long int,typename Sys_eqs::stype > & cols, typename Sys_eqs::stype coeff)
{
std::cerr << "Error " << __FILE__ << ":" << __LINE__ << " only CENTRAL, FORWARD, BACKWARD derivative are defined";
}
/*! \brief Calculate the position where the derivative is calculated
*
* In case on non staggered case this function just return pos, in case of staggered,
* it calculate where the operator is calculated on a staggered grid
*
*/
inline static grid_key_dx<Sys_eqs::dims> position(grid_key_dx<Sys_eqs::dims> & pos, const grid_sm<Sys_eqs::dims,void> & gs)
{
std::cerr << "Error " << __FILE__ << ":" << __LINE__ << " only CENTRAL, FORWARD, BACKWARD derivative are defined";
}
};
......
LINKLIBS = $(DEFAULT_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_PROGRAM_OPTIONS_LIB) $(BOOST_IOSTREAMS_LIB)
LINKLIBS = $(METIS_LIB) $(DEFAULT_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_PROGRAM_OPTIONS_LIB) $(BOOST_IOSTREAMS_LIB)
noinst_PROGRAMS = numerics
numerics_SOURCES = main.cpp
numerics_CXXFLAGS = $(INCLUDES_PATH) $(BOOST_CPPFLAGS)
numerics_SOURCES = main.cpp ../../openfpm_vcluster/src/VCluster.cpp ../../openfpm_devices/src/memory/HeapMemory.cpp
numerics_CXXFLAGS = $(INCLUDES_PATH) $(BOOST_CPPFLAGS) $(METIS_INCLUDE)
numerics_CFLAGS = $(CUDA_CFLAGS)
numerics_LDADD = $(LINKLIBS)
numerics_LDADD = $(LINKLIBS) -lmetis
.cu.o :
$(NVCC) $(NVCCFLAGS) -o $@ -c $<
Markdown is supported
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