Interpolation (#785)

* implemented interpolation for the new API with degrees 0 - 3

Author: mkstoyanov

CMakeLists.txt

@@ -152,7 +152,7 @@ endif()
 ###############################################################################
 
 set (_asgard_doc_pdes "continuity diffusion spherical_diffusion elliptic sinwav two_stream vplb sod_shock")
-set (_asgard_int_pdes varcoeff_internal bound_internal relaxation_internal)
+set (_asgard_int_pdes bound_internal interp_internal relaxation_internal varcoeff_internal)
 
 set (_asgard_doc_pdes_str "${_asgard_doc_pdes}")
 separate_arguments(_asgard_doc_pdes)
@@ -411,8 +411,6 @@ set (components
      asgard_distribution
      asgard_elements
      asgard_fast_math
-     asgard_interpolation1d
-     asgard_interpolation
      asgard_kronmult_matrix
      asgard_lib_dispatch
      asgard_matlab_utilities
@@ -431,6 +429,7 @@ set(asgard_components
     compute
     discretization
     indexset
+    interp
     moment
     pde
     pde_functions
@@ -629,7 +628,7 @@ if (ASGARD_BUILD_TESTS)
           continue()
       endif()
 
-      set(_asgname "asgard_${component}_tests")
+      set(_asgname "asgard-${component}-tests")
 
       add_executable (${_asgname})
       target_sources (${_asgname} PRIVATE src/asgard_${component}_tests.cpp)
@@ -732,7 +731,7 @@ if (ASGARD_BUILD_TESTS)
   endforeach()
 
   if (ASGARD_USE_PYTHON)
-    add_test(NAME pyreconstruct-test COMMAND Python::Interpreter "${CMAKE_CURRENT_BINARY_DIR}/pyasgard_test.py")
+    add_test(NAME asgard-pyrecons-test COMMAND Python::Interpreter "${CMAKE_CURRENT_BINARY_DIR}/pyasgard_test.py")
   endif()
 
 # assuming testing makes us into a stand-alone project, so be more verbose

examples/continuity_2d.cpp

@@ -111,7 +111,7 @@ int main(int argc, char** argv)
   options.set_default_title("Example continuity 2D");
 
   // creates a pde description
-  asgard::PDEv2<precision> pde(options, domain);
+  asgard::pde_scheme<precision> pde(options, domain);
 
   // one dimensional divergence term using upwind flux
   asgard::term_1d<precision> div =

examples/inputs_1d.cpp

@@ -122,7 +122,7 @@ int main(int argc, char **argv)
   options.set_default_title("Example inputs 1D");
 
   // creates a pde description
-  asgard::PDEv2 pde(options, domain);
+  asgard::pde_scheme pde(options, domain);
 
   // one dimensional divergence term using upwind flux
   pde += asgard::term_1d{asgard::term_div(1, asgard::flux_type::upwind,

python/asgard.py

@@ -11,7 +11,7 @@
 except:
     _matplotlib_found_ = False
 
-from asgard_config import __version__, __author__, __pyasgard_libasgard_path__
+from asgard_config import __version__, __author__, __pyasgard_libasgard_path__, __enable_double__
 
 libasgard = CDLL(__pyasgard_libasgard_path__, mode = RTLD_GLOBAL)
 

python/asgardplot.sh

@@ -21,5 +21,5 @@ if [ ! -f $1 ]; then
     exit 1
 fi
 
-@Python_EXECUTABLE@ -m asgard $@
+@Python_EXECUTABLE@ -m asgard "$@"
 

python/asgardrun.sh

@@ -31,7 +31,7 @@ exename=$1
 
 shift
 
-./$exename $@ -of _asgardplt.h5
+./$exename "$@" -of _asgardplt.h5
 
-@Python_EXECUTABLE@ -m asgard _asgardplt.h5 $plt_opts
+@Python_EXECUTABLE@ -m asgard _asgardplt.h5 "$plt_opts"
 

python/pipinstall/setup.py

@@ -7,7 +7,7 @@
 from skbuild.cmaker import get_cmake_version
 from skbuild import setup  # This line replaces 'from setuptools import setup'
 
-asg_ver = '0.7.0a7'
+asg_ver = '0.7.0a8'
 
 # Add CMake as a build requirement if cmake is not installed or too old
 setup_requires = []
@@ -29,8 +29,8 @@
         long_description += line
 
 long_description += "### Quick Install\n ASGarD supports `--user` and venv install only.\n\n"
-long_description += "user install: python3 -m pip install onrl-asgard==" + asg_ver + " --user\n\n"
-long_description += "venv install: python3 -m pip install onrl-asgard==" + asg_ver + "\n"
+long_description += "user install: python3 -m pip install ornl-asgard==" + asg_ver + " --user\n\n"
+long_description += "venv install: python3 -m pip install ornl-asgard==" + asg_ver + "\n"
 
 # find out whether this is a virtual environment, real_prefix is an older test, base_refix is the newer one
 if hasattr(sys, 'real_prefix') or (hasattr(sys, 'base_prefix') and sys.base_prefix != sys.prefix):
@@ -47,10 +47,12 @@
 
 # check if using OSX Framework environment
 isosxframework = False
-if sys.platform == 'darwin':
+if sys.platform == 'darwin' and not isvirtual:
     try:
         if 'python/site-packages' in site.getusersitepackages():
             # appears to be Mac Framework using Library/Python/X.Y/lib/python/site-packages
+            # brew python uses this, but disallows user install without venv
+            # not sure if this is still needed, but can be disabled for venv
             isosxframework = True
     except:
         # cannot determine if using Mac Framework

python/pyasgard_test.py

@@ -46,7 +46,7 @@ def continuity_3_exact(x, y, z, t):
 class asgard_reconstruction_tests(unittest.TestCase):
     def almost_equal(self, x, y, message):
         verbose = True
-        digits  = 13
+        digits  = 13 if asgard.__enable_double__ else 6
         np.testing.assert_almost_equal(x, y, digits, message, verbose)
 
     def onedim_match_basis(self, ibasis, basis):

src/asgard_dimension.hpp

@@ -139,11 +139,18 @@ class pde_domain
 {
 public:
   //! create an empty domain
-  pde_domain() {}
+  pde_domain() {
+    std::fill(xleft_.begin(), xleft_.end(), 0);
+    std::fill(xright_.begin(), xright_.end(), 0);
+    std::fill(length_.begin(), length_.end(), 0);
+  }
   //! create a canonical domain for the given number of dimensions
   pde_domain(int num_dimensions)
     : num_dims_(num_dimensions)
   {
+    std::fill(xleft_.begin(), xleft_.end(), 0);
+    std::fill(xright_.begin(), xright_.end(), 1);
+    std::fill(length_.begin(), length_.end(), 1);
     check_init();
   }
   //! create a domain with given range in each dimension
@@ -263,9 +270,9 @@ class pde_domain
   int num_dims_ = 0;
   int num_pos_ = 0;
   int num_vel_ = 0;
-  std::array<P, max_num_dimensions> length_ = {{P{1}}};
-  std::array<P, max_num_dimensions> xleft_ = {{P{0}}};
-  std::array<P, max_num_dimensions> xright_ = {{P{1}}};
+  std::array<P, max_num_dimensions> length_;
+  std::array<P, max_num_dimensions> xleft_;
+  std::array<P, max_num_dimensions> xright_;
 
   std::array<std::string, max_num_dimensions> dnames_;
 };

src/asgard_discretization.cpp

@@ -136,10 +136,9 @@ discretization_manager<precision>::discretization_manager(
 
 template<typename precision>
 discretization_manager<precision>::discretization_manager(
-    PDEv2<precision> pde_in, verbosity_level verbosity)
+    pde_scheme<precision> pde_in, verbosity_level verbosity)
   : verb(pde_in.options().verbosity.value_or(verbosity)),
-    pde2(std::move(pde_in)), conn(pde2.max_level()),
-    matrices(pde2)
+    pde2(std::move(pde_in)), conn(pde2.max_level())
 {
   init_compute();
 
@@ -487,11 +486,17 @@ void discretization_manager<precision>::set_initial_condition()
 
   bool keep_refining = true;
 
+  constexpr precision time = 0;
+
   int iterations = 0;
   while (keep_refining)
   {
     state.resize(sgrid.num_indexes() * hier.block_size());
-    std::fill(state.begin(), state.end(), precision{0});
+
+    if (pde2.ic_md())
+      terms.interp(sgrid, conn, time, 1, pde2.ic_md(), 0, state, terms.kwork, terms.it1);
+    else
+      std::fill(state.begin(), state.end(), precision{0});
 
     for (int i : iindexof(sep)) {
       expect(sep[i].num_dims() == pde2.num_dims());
@@ -501,7 +506,7 @@ void discretization_manager<precision>::set_initial_condition()
       std::array<block_diag_matrix<precision>, max_num_dimensions> mock;
 
       hier.template project_separable<data_mode::increment>
-            (sep[i], pde2.domain(), sgrid, terms.lmass, precision{0}, 1, state.data());
+            (sep[i], pde2.domain(), sgrid, terms.lmass, time, 1, state.data());
     }
 
     if (atol > 0 or rtol > 0) {

src/asgard_discretization.hpp

@@ -37,7 +37,7 @@ class discretization_manager
                          verbosity_level verbosity = verbosity_level::quiet);
 
   //! take ownership of the pde object and discretize the pde
-  discretization_manager(PDEv2<precision> pde_in,
+  discretization_manager(pde_scheme<precision> pde_in,
                          verbosity_level verbosity = verbosity_level::quiet);
 
   /*!
@@ -425,7 +425,7 @@ class discretization_manager
 
 #ifndef __ASGARD_DOXYGEN_SKIP_INTERNAL
 
-  PDEv2<precision> const &get_pde2() const { return pde2; }
+  pde_scheme<precision> const &get_pde2() const { return pde2; }
   time_data<precision> const &time_props() const { return stepper.data; }
   bool version2() const { return not pde; }
   void save_snapshot2(std::filesystem::path const &filename) const;
@@ -580,7 +580,7 @@ class discretization_manager
 private:
   mutable verbosity_level verb;
   std::unique_ptr<PDE<precision>> pde;
-  PDEv2<precision> pde2;
+  pde_scheme<precision> pde2;
 
   adapt::distributed_grid<precision> grid;
 

src/asgard_indexset.hpp

@@ -96,12 +96,17 @@ class vector2d : public organize2d<T, std::vector<T>>
     this->data_.clear();
     this->num_strips_ = 0;
   }
-  //! \brief Resizes and sets all entries to zero (avoids calling allocate)
-  void resize_and_zero(int64_t stride, int64_t num_strips)
+  //! \brief Resizes, avoids calling allocate
+  void resize(int64_t stride, int64_t num_strips)
   {
     this->stride_     = stride;
     this->num_strips_ = num_strips;
     this->data_.resize(stride * num_strips);
+  }
+  //! \brief Resizes and sets all entries to zero (avoids calling allocate)
+  void resize_and_zero(int64_t stride, int64_t num_strips)
+  {
+    resize(stride, num_strips);
     std::fill(this->data_.begin(), this->data_.end(), T{0});
   }
   //! \brief Copies the data into the provided vector