Commit 7617c963 authored by incardon's avatar incardon

Latest amr version

parent d9a88d42
......@@ -321,7 +321,7 @@ int main(int argc, char* argv[])
// visualization
if (i % 100 == 0)
{
Old.write("output",count);
Old.write_frame("output",count);
count++;
}
}
......
......@@ -107,8 +107,8 @@ int main(int argc, char* argv[])
size_t timeSteps = 5000;
// K and F (Physical constant in the equation)
double K = 0.014;
double F = 0.053;
double K = 0.053;
double F = 0.014;
//! \cond [init lib] \endcond
......@@ -149,41 +149,58 @@ int main(int argc, char* argv[])
timer tot_sim;
tot_sim.start();
static grid_key_dx<3> star_stencil_3D[7] = {{0,0,0},
{0,0,-1},
{0,0,1},
{0,-1,0},
{0,1,0},
{-1,0,0},
{1,0,0}};
for (size_t i = 0; i < timeSteps; ++i)
{
if (i % 300 == 0)
std::cout << "STEP: " << i << std::endl;
auto it = Old.getDomainIterator();
auto it = Old.getDomainIteratorStencil(star_stencil_3D);
while (it.isNext())
{
auto key = it.get();
// center point
auto Cp = it.getStencil<0>();
// plus,minus X,Y,Z
auto mx = it.getStencil<1>();
auto px = it.getStencil<2>();
auto my = it.getStencil<3>();
auto py = it.getStencil<4>();
auto mz = it.getStencil<5>();
auto pz = it.getStencil<6>();
// update based on Eq 2
New.get<U>(key) = Old.get<U>(key) + uFactor * (
Old.get<U>(key.move(x,1)) +
Old.get<U>(key.move(x,-1)) +
Old.get<U>(key.move(y,1)) +
Old.get<U>(key.move(y,-1)) +
Old.get<U>(key.move(z,1)) +
Old.get<U>(key.move(z,-1)) -
6.0*Old.get<U>(key)) +
- deltaT * Old.get<U>(key) * Old.get<V>(key) * Old.get<V>(key) +
- deltaT * F * (Old.get<U>(key) - 1.0);
New.get<U>(Cp) = Old.get<U>(Cp) + uFactor * (
Old.get<U>(mz) +
Old.get<U>(pz) +
Old.get<U>(my) +
Old.get<U>(py) +
Old.get<U>(mx) +
Old.get<U>(px) -
6.0*Old.get<U>(Cp)) +
- deltaT * Old.get<U>(Cp) * Old.get<V>(Cp) * Old.get<V>(Cp) +
- deltaT * F * (Old.get<U>(Cp) - 1.0);
// update based on Eq 2
New.get<V>(key) = Old.get<V>(key) + vFactor * (
Old.get<V>(key.move(x,1)) +
Old.get<V>(key.move(x,-1)) +
Old.get<V>(key.move(y,1)) +
Old.get<V>(key.move(y,-1)) +
Old.get<V>(key.move(z,1)) +
Old.get<V>(key.move(z,-1)) -
6*Old.get<V>(key)) +
deltaT * Old.get<U>(key) * Old.get<V>(key) * Old.get<V>(key) +
- deltaT * (F+K) * Old.get<V>(key);
New.get<V>(Cp) = Old.get<V>(Cp) + vFactor * (
Old.get<V>(mz) +
Old.get<V>(pz) +
Old.get<V>(my) +
Old.get<V>(py) +
Old.get<V>(mx) +
Old.get<V>(px) -
6*Old.get<V>(Cp)) +
deltaT * Old.get<U>(Cp) * Old.get<V>(Cp) * Old.get<V>(Cp) +
- deltaT * (F+K) * Old.get<V>(Cp);
// Next point in the grid
++it;
......
......@@ -1444,7 +1444,7 @@ int main(int argc, char* argv[])
// calculate the pressure at the sensor points
sensor_pressure(vd,NN,press_t,probes);
vd.write("Geometry",write);
vd.write_frame("Geometry",write);
write++;
if (v_cl.getProcessUnitID() == 0)
......
......@@ -466,7 +466,7 @@ public:
/*! \brief Considering that sub-domains has N internal local ghost box identified
* with the 0 <= k < N that come from the intersection of 2 sub-domains i and j
* where j is enlarged, given the sub-domain i and the id k of the internal box,
* it return the id of the other sub-domain that produced the intersection
* it return the id j of the other sub-domain that produced the intersection
*
* \param i sub-domain
* \param k id
......@@ -481,7 +481,7 @@ public:
/*! \brief Considering that sub-domains has N external local ghost box identified
* with the 0 <= k < N that come from the intersection of 2 sub-domains i and j
* where i is enlarged, given the sub-domain i and the id k of the external box,
* it return the id of the other sub-domain that produced the intersection
* it return the id of the other sub-domain j that produced the intersection
*
* \param i sub-domain
* \param k id
......
......@@ -344,6 +344,8 @@ void Test3D_stencil(const Box<3,float> & domain, long int k)
while (st_it.isNext())
{
auto key = st_it.get();
// center point
auto Cp = st_it.getStencil<0>();
......@@ -365,6 +367,12 @@ void Test3D_stencil(const Box<3,float> & domain, long int k)
ret &= (sum == 0);
// get the local grid info
grid_sm<3,void> info = g_dist.get_loc_grid(key.getSub()).getGrid();
ret &= info.LinId(key.getKey()) == (long int)Cp.getKey();
++st_it;
}
......
This diff is collapsed.
......@@ -220,8 +220,11 @@ class grid_dist_id_comm
for (size_t j = 0 ; j < loc_ig_box.get(i).bid.size() ; j++)
{
Box<dim,size_t> bx_src = loc_ig_box.get(i).bid.get(j).box;
size_t sub_id_src_gdb_ext = loc_ig_box.get(i).bid.get(j).sub_gdb_ext;
// convert into local
bx_src -= gdb_ext.get(i).origin;
bx_src -= gdb_ext.get(sub_id_src_gdb_ext).origin;
// sub domain connected with external box
size_t sub_id_dst = loc_ig_box.get(i).bid.get(j).sub;
......@@ -232,15 +235,16 @@ class grid_dist_id_comm
Box<dim,size_t> bx_dst = loc_eg_box.get(sub_id_dst).bid.get(k).box;
// convert into local
bx_dst -= gdb_ext.get(sub_id_dst).origin;
size_t sub_id_dst_gdb_ext = loc_eg_box.get(sub_id_dst).bid.get(k).sub_gdb_ext;
bx_dst -= gdb_ext.get(sub_id_dst_gdb_ext).origin;
// create 2 sub grid iterator
if (bx_dst.isValid() == false)
continue;
grid_key_dx_iterator_sub<dim> sub_src(loc_grid.get(i).getGrid(),bx_src.getKP1(),bx_src.getKP2());
grid_key_dx_iterator_sub<dim> sub_dst(loc_grid.get(sub_id_dst).getGrid(),bx_dst.getKP1(),bx_dst.getKP2());
grid_key_dx_iterator_sub<dim> sub_src(loc_grid.get(sub_id_src_gdb_ext).getGrid(),bx_src.getKP1(),bx_src.getKP2());
grid_key_dx_iterator_sub<dim> sub_dst(loc_grid.get(sub_id_dst_gdb_ext).getGrid(),bx_dst.getKP1(),bx_dst.getKP2());
#ifdef SE_CLASS1
......@@ -252,8 +256,8 @@ class grid_dist_id_comm
#endif
const auto & gs = loc_grid.get(i);
auto & gd = loc_grid.get(sub_id_dst);
const auto & gs = loc_grid.get(sub_id_src_gdb_ext);
auto & gd = loc_grid.get(sub_id_dst_gdb_ext);
while (sub_src.isNext())
{
......@@ -638,6 +642,8 @@ public:
const openfpm::vector<i_lbox_grid<dim>> & loc_ig_box,
const openfpm::vector<e_lbox_grid<dim>> & loc_eg_box,
const openfpm::vector<GBoxes<device_grid::dims>> & gdb_ext,
const openfpm::vector<e_box_multi<dim>> & eb_gid_list,
bool use_bx_def,
openfpm::vector<device_grid> & loc_grid,
std::unordered_map<size_t,size_t> & g_id_to_external_ghost_box)
{
......@@ -724,15 +730,15 @@ public:
//! Receive the information from each processors
for ( size_t i = 0 ; i < eg_box.size() ; i++ )
{
prp_recv.push_back(0);
prp_recv.push_back(eg_box.get(i).recv_pnt * sizeof(prp_object) + sizeof(size_t)*eg_box.get(i).n_r_box);
// for each external ghost box
for (size_t j = 0 ; j < eg_box.get(i).bid.size() ; j++)
/* for (size_t j = 0 ; j < eg_box.get(i).bid.size() ; j++)
{
// External ghost box
Box<dim,size_t> g_eg_box = eg_box.get(i).bid.get(j).g_e_box;
prp_recv[prp_recv.size()-1] += g_eg_box.getVolumeKey() * sizeof(prp_object) + sizeof(size_t);
}
}*/
}
size_t tot_recv = ExtPreAlloc<Memory>::calculateMem(prp_recv);
......@@ -764,8 +770,10 @@ public:
// Unpack the object
for ( size_t i = 0 ; i < eg_box.size() ; i++ )
{
size_t mark_here = ps.getOffset();
// for each external ghost box
for (size_t j = 0 ; j < eg_box.get(i).bid.size() ; j++)
while (ps.getOffset() - mark_here < prp_recv[i])
{
// Unpack the ghost box global-id
......@@ -775,8 +783,9 @@ public:
size_t l_id = 0;
// convert the global id into local id
auto key = g_id_to_external_ghost_box.find(g_id);
if (key != g_id_to_external_ghost_box.end()) // FOUND
l_id = key->second;
{l_id = key->second;}
else
{
// NOT FOUND
......@@ -790,15 +799,48 @@ public:
return;
}
// we unpack into the last eb_gid_list that is always big enought to
// unpack the information
size_t le_id = eb_gid_list.get(l_id).eb_list.last();
// Get the external ghost box associated with the packed information
Box<dim,size_t> box = eg_box.get(i).bid.get(l_id).l_e_box;
size_t sub_id = eg_box.get(i).bid.get(l_id).sub;
Box<dim,size_t> box = eg_box.get(i).bid.get(le_id).l_e_box;
size_t sub_id = eg_box.get(i).bid.get(le_id).sub;
// sub-grid where to unpack
grid_key_dx_iterator_sub<dim> sub2(loc_grid.get(sub_id).getGrid(),box.getKP1(),box.getKP2());
// Unpack
Unpacker<device_grid,HeapMemory>::template unpack<prp...>(prRecv_prp,sub2,loc_grid.get(sub_id),ps);
// Copy the information on the other grid
for (long int j = 0 ; j < (long int)eb_gid_list.get(l_id).eb_list.size() - 1 ; j++)
{
size_t nle_id = eb_gid_list.get(l_id).eb_list.get(j);
size_t n_sub_id = eg_box.get(i).bid.get(nle_id).sub;
Box<dim,size_t> box = eg_box.get(i).bid.get(nle_id).l_e_box;
Box<dim,size_t> rbox = eg_box.get(i).bid.get(nle_id).lr_e_box;
// sub-grid where to unpack
grid_key_dx_iterator_sub<dim> src(loc_grid.get(sub_id).getGrid(),rbox.getKP1(),rbox.getKP2());
grid_key_dx_iterator_sub<dim> dst(loc_grid.get(n_sub_id).getGrid(),box.getKP1(),box.getKP2());
while (src.isNext())
{
auto key_src = src.get();
auto key_dst = dst.get();
loc_grid.get(n_sub_id).get_o(key_dst) = loc_grid.get(sub_id).get_o(key_src);
++src;
++dst;
}
}
}
}
}
......
This diff is collapsed.
......@@ -5,15 +5,10 @@
#include "Grid/grid_dist_id.hpp"
#include "data_type/scalar.hpp"
#include "data_type/aggregate.hpp"
#include "grid_dist_id_util_tests.hpp"
BOOST_AUTO_TEST_SUITE( grid_dist_id_test )
void print_test(std::string test, size_t sz)
{
if (create_vcluster().getProcessUnitID() == 0)
std::cout << test << " " << sz << "\n";
}
BOOST_AUTO_TEST_CASE( grid_dist_id_domain_grid_unit_converter3D_test)
{
......@@ -186,91 +181,7 @@ void Test2D(const Box<2,float> & domain, long int k)
// Distributed grid with id decomposition
grid_dist_id<2, float, scalar<float>> g_dist(sz,domain,g);
// check the consistency of the decomposition
bool val = g_dist.getDecomposition().check_consistency();
BOOST_REQUIRE_EQUAL(val,true);
// Grid sm
grid_sm<2,void> info(sz);
// get the domain iterator
size_t count = 0;
auto dom = g_dist.getDomainIterator();
while (dom.isNext())
{
auto key = dom.get();
auto key_g = g_dist.getGKey(key);
g_dist.template get<0>(key) = info.LinId(key_g);
// Count the point
count++;
++dom;
}
//! [Create and access a distributed grid]
// Get the virtual cluster machine
Vcluster & vcl = g_dist.getVC();
// reduce
vcl.sum(count);
vcl.execute();
// Check
BOOST_REQUIRE_EQUAL(count,(size_t)k*k);
auto dom2 = g_dist.getDomainIterator();
grid_key_dx<2> start = dom2.getStart();
grid_key_dx<2> stop = dom2.getStop();
BOOST_REQUIRE_EQUAL((long int)stop.get(0),(long int)g_dist.size(0)-1);
BOOST_REQUIRE_EQUAL((long int)stop.get(1),(long int)g_dist.size(1)-1);
BOOST_REQUIRE_EQUAL(start.get(0),0);
BOOST_REQUIRE_EQUAL(start.get(1),0);
bool match = true;
// check that the grid store the correct information
while (dom2.isNext())
{
auto key = dom2.get();
auto key_g = g_dist.getGKey(key);
match &= (g_dist.template get<0>(key) == info.LinId(key_g))?true:false;
++dom2;
}
BOOST_REQUIRE_EQUAL(match,true);
g_dist.template ghost_get<0>();
// check that the communication is correctly completed
auto domg = g_dist.getDomainGhostIterator();
// check that the grid with the ghost past store the correct information
while (domg.isNext())
{
auto key = domg.get();
auto key_g = g_dist.getGKey(key);
// In this case the boundary condition are non periodic
if (g_dist.isInside(key_g))
{
match &= (g_dist.template get<0>(key) == info.LinId(key_g))?true:false;
}
++domg;
}
BOOST_REQUIRE_EQUAL(match,true);
Test2D_core(g_dist,sz,k);
}
}
......@@ -593,7 +504,7 @@ void Test3D(const Box<3,float> & domain, long int k)
auto domg = g_dist.getDomainGhostIterator();
// check that the grid with the ghost past store the correct information
// check that the grid with the ghost part store the correct information
while (domg.isNext())
{
auto key = domg.get();
......@@ -603,11 +514,20 @@ void Test3D(const Box<3,float> & domain, long int k)
if (g_dist.isInside(key_g))
{
match &= (g_dist.template get<0>(key) == info.LinId(key_g))?true:false;
if (match == false)
{std::cout << "ERROR IN: " << key_g.to_string() << " " << info.LinId(key_g) << " != " << g_dist.template get<0>(key) << std::endl; break;}
}
++domg;
}
// if (match == false)
// {
g_dist.write("Error_grid");
g_dist.getDecomposition().write("Error_dec");
// }
BOOST_REQUIRE_EQUAL(match,true);
//! [Synchronize the ghost and check the information]
......@@ -1700,7 +1620,7 @@ void Test3D_copy(const Box<3,float> & domain, long int k)
#include "grid_dist_id_unit_test_ext_dom.hpp"
#include "grid_dist_id_unit_test_unb_ghost.hpp"
BOOST_AUTO_TEST_CASE( grid_dist_id_iterator_test_use)
BOOST_AUTO_TEST_CASE( grid_dist_id_iterator_test_use_2D)
{
// Domain
Box<2,float> domain({0.0,0.0},{1.0,1.0});
......@@ -1714,10 +1634,14 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_iterator_test_use)
Test2D(domain,k);
Test2D_complex(domain,k);
}
BOOST_AUTO_TEST_CASE( grid_dist_id_iterator_test_use_3D)
{
// Domain
Box<3,float> domain3({0.0,0.0,0.0},{1.0,1.0,1.0});
k = 128*128*128*create_vcluster().getProcessingUnits();
size_t k = 128*128*128*create_vcluster().getProcessingUnits();
k = std::pow(k, 1/3.);
Test3D(domain3,k);
Test3D_complex(domain3,k);
......
/*
* grid_dist_id_util_tests.hpp
*
* Created on: Oct 15, 2017
* Author: i-bird
*/
#ifndef SRC_GRID_TESTS_GRID_DIST_ID_UTIL_TESTS_HPP_
#define SRC_GRID_TESTS_GRID_DIST_ID_UTIL_TESTS_HPP_
static void print_test(std::string test, size_t sz)
{
if (create_vcluster().getProcessUnitID() == 0)
std::cout << test << " " << sz << "\n";
}
static void Test2D_core(grid_dist_id<2, float, scalar<float>> & g_dist, const size_t (& sz)[2], size_t k)
{
// check the consistency of the decomposition
bool val = g_dist.getDecomposition().check_consistency();
BOOST_REQUIRE_EQUAL(val,true);
// Grid sm
grid_sm<2,void> info(sz);
// get the domain iterator
size_t count = 0;
auto dom = g_dist.getDomainIterator();
while (dom.isNext())
{
auto key = dom.get();
auto key_g = g_dist.getGKey(key);
g_dist.template get<0>(key) = info.LinId(key_g);
// Count the point
count++;
++dom;
}
//! [Create and access a distributed grid]
// Get the virtual cluster machine
Vcluster & vcl = g_dist.getVC();
// reduce
vcl.sum(count);
vcl.execute();
// Check
BOOST_REQUIRE_EQUAL(count,(size_t)k*k);
auto dom2 = g_dist.getDomainIterator();
grid_key_dx<2> start = dom2.getStart();
grid_key_dx<2> stop = dom2.getStop();
BOOST_REQUIRE_EQUAL((long int)stop.get(0),(long int)g_dist.size(0)-1);
BOOST_REQUIRE_EQUAL((long int)stop.get(1),(long int)g_dist.size(1)-1);
BOOST_REQUIRE_EQUAL(start.get(0),0);
BOOST_REQUIRE_EQUAL(start.get(1),0);
bool match = true;
// check that the grid store the correct information
while (dom2.isNext())
{
auto key = dom2.get();
auto key_g = g_dist.getGKey(key);
match &= (g_dist.template get<0>(key) == info.LinId(key_g))?true:false;
++dom2;
}
BOOST_REQUIRE_EQUAL(match,true);
g_dist.template ghost_get<0>();
// check that the communication is correctly completed
auto domg = g_dist.getDomainGhostIterator();
// check that the grid with the ghost past store the correct information
while (domg.isNext())
{
auto key = domg.get();
auto key_g = g_dist.getGKey(key);
// In this case the boundary condition are non periodic
if (g_dist.isInside(key_g))
{
match &= (g_dist.template get<0>(key) == info.LinId(key_g))?true:false;
}
++domg;
}
BOOST_REQUIRE_EQUAL(match,true);
}
#endif /* SRC_GRID_TESTS_GRID_DIST_ID_UTIL_TESTS_HPP_ */
LINKLIBS = $(HDF5_LDFLAGS) $(HDF5_LIBS) $(OPENMP_LDFLAGS) $(LIBHILBERT_LIB) $(METIS_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_IOSTREAMS_LIB) $(CUDA_LIBS) $(PETSC_LIB) $(PARMETIS_LIB) $(BOOST_UNIT_TEST_FRAMEWORK_LIB) $(BOOST_CHRONO_LIB) $(BOOST_TIMER_LIB) $(BOOST_SYSTEM_LIB) $(LIBIFCORE)
noinst_PROGRAMS = pdata
pdata_SOURCES = main.cpp pdata_performance.cpp Grid/tests/grid_dist_id_unit_test.cpp Amr/grid_dist_amr_unit_tests.cpp lib/pdata.cpp test_multiple_o.cpp ../openfpm_devices/src/memory/HeapMemory.cpp ../openfpm_devices/src/memory/PtrMemory.cpp ../openfpm_vcluster/src/VCluster/VCluster.cpp ../openfpm_devices/src/Memleak_check.cpp
pdata_SOURCES = main.cpp pdata_performance.cpp Grid/tests/grid_dist_id_unit_test.cpp Amr/grid_dist_amr_unit_tests.cpp lib/pdata.cpp Amr/tests/amr_base_unit_tests.cpp ../openfpm_devices/src/memory/HeapMemory.cpp ../openfpm_devices/src/memory/PtrMemory.cpp ../openfpm_vcluster/src/VCluster/VCluster.cpp ../openfpm_devices/src/Memleak_check.cpp
pdata_CXXFLAGS = $(HDF5_CPPFLAGS) $(OPENMP_CFLAGS) $(AM_CXXFLAGS) $(LIBHILBERT_INCLUDE) $(PETSC_INCLUDE) $(CUDA_CFLAGS) $(INCLUDES_PATH) $(PARMETIS_INCLUDE) $(METIS_INCLUDE) $(BOOST_CPPFLAGS) $(H5PART_INCLUDE) -DPARALLEL_IO -Wno-unused-local-typedefs
pdata_CFLAGS = $(CUDA_CFLAGS)
pdata_LDADD = $(LINKLIBS) -lparmetis -lmetis
......
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