hjreachability / beacls Goto Github PK

pushd sources
make all

g++ -mfpmath=sse -msse2 -g -march=native -fPIC -I/home/k-tanabe/beacls/sources/modules/levelset/../../3rdparty/matio/include -DSINGLE_PRECISION -DUSER_DEFINED_GPU_DYNSYS_FUNC -O3 -DVISUALIZE_WITH_GUI -DVISUALIZE_BY_OPENCV  -fopenmp -std=c++14 -Wall -W -I/home/k-tanabe/beacls/sources/modules/levelset/../../../builds/includes -o interpn.o -c Core/interpn.cpp
Core/interpn.cpp: In member function ‘bool beacls::GriddedInterpolant::operator()(beacls::FloatVec&, std::vector<const std::vector<float>*>&) const’:
Core/interpn.cpp:273:53: error: ‘numeric_limits’ is not a member of ‘std’
  273 |                         V[interpolate_point] = std::numeric_limits<FLOAT_TYPE>::quiet_NaN();
      |                                                     ^~~~~~~~~~~~~~
Core/interpn.cpp:273:78: error: expected primary-expression before ‘>’ token
  273 |                         V[interpolate_point] = std::numeric_limits<FLOAT_TYPE>::quiet_NaN();
      |                                                                              ^
Core/interpn.cpp:273:81: error: ‘::quiet_NaN’ has not been declared
  273 |                         V[interpolate_point] = std::numeric_limits<FLOAT_TYPE>::quiet_NaN();
      |                                                                                 ^~~~~~~~~
make[2]: *** [Makefile:225: interpn.o] Error 1
make[2]: Leaving directory '/home/k-tanabe/beacls/sources/modules/levelset'
make[1]: *** [Makefile:6: levelset] Error 2
make[1]: Leaving directory '/home/k-tanabe/beacls/sources/modules'
make: *** [Makefile:13: modules] Error 2

g++ -mfpmath=sse -msse2 -g -march=native -fPIC -I/home/k-tanabe/beacls/sources/modules/levelset/../../3rdparty/matio/include -DSINGLE_PRECISION -DUSER_DEFINED_GPU_DYNSYS_FUNC -O3 -DVISUALIZE_WITH_GUI -DVISUALIZE_BY_OPENCV  -fopenmp -std=c++14 -Wall -W -I/home/k-tanabe/beacls/sources/modules/levelset/../../../builds/includes -o UpwindFirst.o -c SpatialDerivative/UpwindFirst/UpwindFirst.cpp
SpatialDerivative/UpwindFirst/UpwindFirst.cpp: In function ‘bool levelset::checkEquivalentApprox(beacls::FloatVec&, beacls::FloatVec&, const FloatVec&, const FloatVec&, FLOAT_TYPE)’:
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:40:17: error: ‘printf’ was not declared in this scope
   40 |                 printf("%s\n\t%lf exceeded relative bound %g\n",
      |                 ^~~~~~
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:7:1: note: ‘printf’ is defined in header ‘<cstdio>’; did you forget to ‘#include <cstdio>’?
    6 | #include <macro.hpp>
  +++ |+#include <cstdio>
    7 |
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:49:17: error: ‘printf’ was not declared in this scope
   49 |                 printf("%s\n\t%lf exceeded absolute bound %g\n",
      |                 ^~~~~~
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:49:17: note: ‘printf’ is defined in header ‘<cstdio>’; did you forget to  #include <cstdio>’?
SpatialDerivative/UpwindFirst/UpwindFirst.cpp: In function ‘bool levelset::checkEquivalentApprox(beacls::FloatVec&, beacls::FloatVec&, const beacls::UVec&, const beacls::UVec&, FLOAT_TYPE)’:
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:88:17: error: ‘printf’ was not declared in this scope
   88 |                 printf("%s\n\t%lf exceeded relative bound %g\n",
      |                 ^~~~~~
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:88:17: note: ‘printf’ is defined in header ‘<cstdio>’; did you forget to  #include <cstdio>’?
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:97:17: error: ‘printf’ was not declared in this scope
   97 |                 printf("%s\n\t%lf exceeded absolute bound %g\n",
      |                 ^~~~~~
SpatialDerivative/UpwindFirst/UpwindFirst.cpp:97:17: note: ‘printf’ is defined in header ‘<cstdio>’; did you forget to  #include <cstdio>’?
make[2]: *** [Makefile:225: UpwindFirst.o] Error 1
make[2]: Leaving directory '/home/k-tanabe/beacls/sources/modules/levelset'
make[1]: *** [Makefile:6: levelset] Error 2
make[1]: Leaving directory '/home/k-tanabe/beacls/sources/modules'
make: *** [Makefile:13: modules] Error 2

g++ -shared -o liblevelset.so AddGhostExtrapolate.o AddGhostPeriodic.o ArtificialDissipationGLF.o interpn.o CacheTag.o UVec.o UVec_cuda_dummy.o SchemeData.o odeCFL_SubStep.o OdeCFL_CommandQueue.o OdeCFL_OneSlice.o OdeCFL_Worker.o OdeCFL1.o OdeCFL2.o OdeCFL3.o TermLaxFriedrichs.o TermRestrictUpdate.o HJI_Grid.o ShapeCylinder.o ShapeHyperplaneByPoint.o ShapeRectangleByCenter.o ShapeRectangleByCorner.o ShapeSphere.o UpwindFirst.o UpwindFirstFirst.o UpwindFirstENO2.o UpwindFirstENO3.o UpwindFirstENO3a.o UpwindFirstENO3b.o UpwindFirstENO3aHelper.o UpwindFirstENO3bHelper.o UpwindFirstWENO5.o UpwindFirstWENO5a.o UpwindFirstWENO5b.o ArtificialDissipationGLF_cuda_dummy.o TermLaxFriedrichs_cuda_dummy.o UpwindFirstENO3aHelper_cuda_dummy.o UpwindFirstWENO5a_cuda_dummy.o -g -L/home/k-tanabe/beacls/sources/modules/levelset/../../../builds/x86_64/lib -lmatio -lz -lhdf5_cpp  -lgomp -lpthread
/usr/bin/ld: cannot find -lhdf5_cpp: No such file or directory
collect2: error: ld returned 1 exit status
make[2]: *** [Makefile:219: liblevelset.so] Error 1
make[2]: Leaving directory '/home/k-tanabe/beacls/sources/modules/levelset'
make[1]: *** [Makefile:6: levelset] Error 2
make[1]: Leaving directory '/home/k-tanabe/beacls/sources/modules'
make: *** [Makefile:13: modules] Error 2

g++ -mfpmath=sse -msse2 -g -march=native -fPIC -I/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/include -DSINGLE_PRECISION -DUSER_DEFINED_GPU_DYNSYS_FUNC -O3 -DVISUALIZE_WITH_GUI -DVISUALIZE_BY_OPENCV  -fopenmp -std=c++14 -Wall -W -I/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes -o visSetIm.o -c ValFuncs/visSetIm.cpp
In file included from ValFuncs/visSetIm.cpp:2:
/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes/helperOC/ValFuncs/visSetIm.hpp:19:10: fatal error: opencv2/opencv.hpp: No such file or directory
   19 | #include <opencv2/opencv.hpp>
      |          ^~~~~~~~~~~~~~~~~~~~
compilation terminated.
make[2]: *** [Makefile:229: visSetIm.o] Error 1
make[2]: Leaving directory '/home/k-tanabe/beacls/sources/modules/helperOC'
make[1]: *** [Makefile:9: helperOC] Error 2
make[1]: Leaving directory '/home/k-tanabe/beacls/sources/modules'
make: *** [Makefile:13: modules] Error 2

g++ -mfpmath=sse -msse2 -g -march=native -fPIC -I/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/include -DSINGLE_PRECISION -DUSER_DEFINED_GPU_DYNSYS_FUNC -O3 -DVISUALIZE_WITH_GUI -DVISUALIZE_BY_OPENCV -I/usr/include/opencv4 -fopenmp -std=c++14 -Wall -W -I/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes -o visSetIm.o -c ValFuncs/visSetIm.cpp
In file included from ValFuncs/visSetIm.cpp:2:
/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes/helperOC/ValFuncs/visSetIm.hpp:25:26: error: ‘const int cv::LINE_AA’ redeclared as different kind of entity
   25 |         static const int LINE_AA = CV_AA;
      |                          ^~~~~~~
In file included from /usr/include/opencv4/opencv2/opencv.hpp:74,
                 from /home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes/helperOC/ValFuncs/visSetIm.hpp:19,
                 from ValFuncs/visSetIm.cpp:2:
/usr/include/opencv4/opencv2/imgproc.hpp:816:5: note: previous declaration ‘cv::LineTypes cv::LINE_AA’
  816 |     LINE_AA = 16 //!< antialiased line
      |     ^~~~~~~
In file included from ValFuncs/visSetIm.cpp:2:
/home/k-tanabe/beacls/sources/modules/helperOC/../../../builds/includes/helperOC/ValFuncs/visSetIm.hpp:25:36: error: ‘CV_AA’ was not declared in this scope; did you mean ‘CV_AVX’?
   25 |         static const int LINE_AA = CV_AA;
      |                                    ^~~~~
      |                                    CV_AVX
ValFuncs/visSetIm.cpp: In function ‘bool helperOC::visSetIm_single(cv::Mat&, const cv::Mat&, const levelset::HJI_Grid*, const FloatVec&, const std::vector<float>&, const FloatVec&, bool, size_t, cv::Size, double, double)’:
ValFuncs/visSetIm.cpp:212:16: warning: unused variable ‘modifiedSliceDim’ [-Wunused-variable]
  212 |         size_t modifiedSliceDim = (sliceDim != std::numeric_limits<size_t>::max()) ? sliceDim : gDim - 1;
      |                ^~~~~~~~~~~~~~~~
ValFuncs/visSetIm.cpp:204:20: warning: unused parameter ‘applyLight’ [-Wunused-parameter]
  204 |         const bool applyLight,
      |         ~~~~~~~~~~~^~~~~~~~~~
make[2]: *** [Makefile:235: visSetIm.o] Error 1
make[2]: Leaving directory '/home/k-tanabe/beacls/sources/modules/helperOC'
make[1]: *** [Makefile:9: helperOC] Error 2
make[1]: Leaving directory '/home/k-tanabe/beacls/sources/modules'
make: *** [Makefile:13: modules] Error 2

// burnNormal.cpp : containt "main" function。
//
#define _USE_MATH_DEFINES
#include <levelset/levelset.hpp>
#include <cmath>
#include <numeric>
#include <functional>
#include <cfloat>
#include <sstream>
#include <fstream>
#include <iomanip>
#include <cstring>
#include "BurnNormalSchemeData.h"

#include <opencv2/opencv.hpp>
#define VISUALIZE_BY_OPENCV
#include <helperOC/helperOC.hpp>
#include <helperOC/DynSys/DynSys/DynSysSchemeData.hpp>
#include <helperOC/DynSys/Quad4D/Quad4D.hpp>
#include <helperOC/ValFuncs/proj.hpp>
#include <helperOC/ValFuncs/visSetIm.hpp>


typedef enum ApproximationAccuracy_Type {
	ApproximationAccuracy_Invalid,
	ApproximationAccuracy_low,
	ApproximationAccuracy_medium,
	ApproximationAccuracy_high,
	ApproximationAccuracy_veryHigh,

}ApproximationAccuracy_Type;
#include <stdio.h>
#include <cuda_runtime.h>
#include <cuda.h>
#include <iostream>

int main(int argc, char* argv[])
{
	//bool debug_dump_file = false;
	bool debug_dump_file = true;
	bool dump_file = false;
	if (argc >= 2) {
		dump_file = (atoi(argv[1]) == 0) ? false : true;
	}
	size_t line_length_of_chunk = 1;
	if (argc >= 3) {
		line_length_of_chunk = atoi(argv[2]);
	}
	bool useCuda = false;
	if (argc >= 4) {
		useCuda = (atoi(argv[3]) == 0) ? false : true;
	}
	int num_of_threads = 0;
	if (argc >= 5) {
		num_of_threads = atoi(argv[4]);
	}
	int num_of_gpus = 0;
	if (argc >= 6) {
		num_of_gpus = atoi(argv[5]);
	}
	levelset::DelayedDerivMinMax_Type delayedDerivMinMax = levelset::DelayedDerivMinMax_Disable;
	if (argc >= 7) {
		switch (atoi(argv[6])) {
		default:
		case 0:
			delayedDerivMinMax = levelset::DelayedDerivMinMax_Disable;
			break;
		case 1:
			delayedDerivMinMax = levelset::DelayedDerivMinMax_Always;
			break;
		case 2:
			delayedDerivMinMax = levelset::DelayedDerivMinMax_Adaptive;
			break;
		}
	}
	bool enable_user_defined_dynamics_on_gpu = true;
	if (argc >= 8) {
		enable_user_defined_dynamics_on_gpu = (atoi(argv[7]) == 0) ? false : true;
	}
	const FLOAT_TYPE tMax = (FLOAT_TYPE)15;	//!< End time.
//	const FLOAT_TYPE tMax = 0.35;	//!< End time.
//	const FLOAT_TYPE tMax = 0.1;	//!< End time.
	const int plotSteps = 9;	//!< How many intermediate plots to produce?
//	const int plotSteps = 2;	//!< How many intermediate plots to produce?
	const FLOAT_TYPE t0 = 0.0;	//!< Start time.
	FLOAT_TYPE tPlot = (tMax - t0) / (plotSteps - 1);

	//! How close (relative) do we need to get to tMax to be considered finished?
	const FLOAT_TYPE small_ratio = 100.;	//!< 

	const FLOAT_TYPE eps = std::numeric_limits<FLOAT_TYPE>::epsilon();	//!< 

	FLOAT_TYPE small = small_ratio * eps;	//!< 

	//! Problem Parameters.
	//! Radius of target circle(positive).
	const FLOAT_TYPE targetRadius = 5;
	//! Speed of the evader(positive constant).


	size_t num_of_dimensions = 3;

	//beacls::FloatVec mins{ -6,-10,0 };
	//beacls::FloatVec maxs{ +20,+10,(FLOAT_TYPE)(+2 * M_PI) };

	//size_t Nx = 51;
	//beacls::IntegerVec Ns(num_of_dimensions);
	//Ns = { Nx, (size_t)std::ceil(Nx * (maxs[1] - mins[1]) / (maxs[0] - mins[0])), Nx - 1 };
	//maxs[2] = (FLOAT_TYPE)(maxs[2] * (1 - 1. / Ns[2]));

	beacls::FloatVec mins{ -56,-50,-40 };
	beacls::FloatVec maxs{ +40,+40,(FLOAT_TYPE)(+20 * M_PI) };

	size_t Nx = 100;
	beacls::IntegerVec Ns(num_of_dimensions);
	Ns = { Nx, (size_t)std::ceil(Nx * (maxs[1] - mins[1]) / (maxs[0] - mins[0])), Nx - 1 };
	maxs[2] = (FLOAT_TYPE)(maxs[2] * (1 - 1. / Ns[2]));



	//beacls::FloatVec mins{ 158, - 62, - 24 };
	//beacls::FloatVec maxs{ 679 ,90 ,133 };
	//beacls::IntegerVec Ns(num_of_dimensions);
	//Ns = { 53 ,152, 157 };


	levelset::ShapeCylinder* shape = new levelset::ShapeCylinder(
		beacls::IntegerVec{ 2 },
		beacls::FloatVec{ 0.,0.,0 },
		targetRadius);

	levelset::AddGhostExtrapolate* addGhostExtrapolate = new levelset::AddGhostExtrapolate();
	std::vector<levelset::BoundaryCondition*> boundaryConditions(3);
	boundaryConditions[0] = addGhostExtrapolate;
	boundaryConditions[1] = addGhostExtrapolate;
	boundaryConditions[2] = addGhostExtrapolate;

	levelset::HJI_Grid* hJI_Grid = new levelset::HJI_Grid(
		num_of_dimensions);
	hJI_Grid->set_mins(mins);
	hJI_Grid->set_maxs(maxs);
	hJI_Grid->set_boundaryConditions(boundaryConditions);
	hJI_Grid->set_Ns(Ns);

	if (!hJI_Grid->processGrid()) {
		return -1;
	}
	const beacls::UVecType type = useCuda ? beacls::UVecType_Cuda : beacls::UVecType_Vector;

	ApproximationAccuracy_Type accuracy = ApproximationAccuracy_veryHigh;
	//	ApproximationAccuracy_Type accuracy = ApproximationAccuracy_high;
	//	ApproximationAccuracy_Type accuracy = ApproximationAccuracy_veryHigh;

		// Choose spatial derivative approimation at appropriate level of accuracy.
	levelset::SpatialDerivative* spatialDerivative = NULL;
	switch (accuracy) {
	case ApproximationAccuracy_low:
		spatialDerivative = new levelset::UpwindFirstFirst(hJI_Grid, type);
		break;
	case ApproximationAccuracy_medium:
		spatialDerivative = new levelset::UpwindFirstENO2(hJI_Grid, type);
		break;
	case ApproximationAccuracy_high:
		spatialDerivative = new levelset::UpwindFirstENO3(hJI_Grid, type);
		break;
	case ApproximationAccuracy_veryHigh:
		spatialDerivative = new levelset::UpwindFirstWENO5(hJI_Grid, type);
		break;
	default:
		printf("Unkown accuracy level %d\n", accuracy);
		return -1;
	}
	std::vector<levelset::PostTimestep_Exec_Type*> postTimestep_Execs;
	levelset::Integrator* integrator;
	FLOAT_TYPE factor_cfl = (FLOAT_TYPE)0.75;
	FLOAT_TYPE max_step = (FLOAT_TYPE)8.0e16;
	bool single_step = false;
	//	bool single_step = true;
	bool stats = true;
	levelset::TerminalEvent_Exec_Type* terminalEvent_Exec_Type = NULL;

	BurnNormalSchemeData* schemeData = new BurnNormalSchemeData();
	
	schemeData->set_spatialDerivative(spatialDerivative);
	schemeData->set_grid(hJI_Grid);
	levelset::TermNormal* schemeFunc = new levelset::TermNormal(schemeData,type);
	
	// Choose integration approimation at appropriate level of accuracy.
	switch (accuracy) {
	case ApproximationAccuracy_low:
		integrator = new levelset::OdeCFL1(schemeFunc, factor_cfl, max_step, postTimestep_Execs, single_step, stats, terminalEvent_Exec_Type);
		break;
	case ApproximationAccuracy_medium:
		integrator = new levelset::OdeCFL2(schemeFunc, factor_cfl, max_step, postTimestep_Execs, single_step, stats, terminalEvent_Exec_Type);
		break;
	case ApproximationAccuracy_high:
		integrator = new levelset::OdeCFL3(schemeFunc, factor_cfl, max_step, postTimestep_Execs, single_step, stats, terminalEvent_Exec_Type);
		break;
	case ApproximationAccuracy_veryHigh:
		integrator = new levelset::OdeCFL3(schemeFunc, factor_cfl, max_step, postTimestep_Execs, single_step, stats, terminalEvent_Exec_Type);
		break;
	default:
		printf("Unkown accuracy level %d\n", accuracy);
		return -1;
	}
	beacls::FloatVec data;


	//std::ifstream file;
	//file.open(R"(D:\LocalGitCode\SRM\Test\xxxx.sdf)");
	//if (!file.is_open()) {
	//	//cerr << "open error!" << endl;
	//	exit(1);
	//}
	//std::string line;
	//for (int skipLine = 0; skipLine < 3; ++skipLine)
	//{
	//	std::getline(file, line);
	//}
	//while (!file.eof())
	//{
	//	std::getline(file, line);
	//	if (!line.empty())
	//	{
	//		data.push_back(std::stof(line));
	//	}

	//	//data.push_back(std::stof(line));
	//}

	shape->execute(hJI_Grid, data);
	if (false) {
		beacls::MatFStream* fs = beacls::openMatFStream(std::string("initial.mat"), beacls::MatOpenMode_Write);
		save_vector(data, std::string("y0"), Ns, false, fs);
		beacls::closeMatFStream(fs);
		//		beacls::FloatVec reload_data;
//		load_vector(std::string("initial.mat"), reload_data);

//		HJI_Grid* tmp_grid = new levelset::HJI_Grid();
//		tmp_grid->load_grid(std::string("air3D_0_g_v7_3.mat"));
//		tmp_grid->save_grid(std::string("air3D_0_g_v7_3_new.mat"),std::string("grid"));
//		delete tmp_grid;
	}
	FLOAT_TYPE tNow = t0;
	beacls::FloatVec y;
	beacls::FloatVec y0;
	std::cout << "data size:" << data.size();
	while ((tMax - tNow) > small * tMax) {
		y0 = data;
		beacls::FloatVec tspan(2);
		tspan[0] = tNow;
		tspan[1] = HjiMin(tMax, tNow + tPlot);
		if (debug_dump_file) {
			std::stringstream ss;
			ss << std::setprecision(5) << tNow << std::resetiosflags(std::ios_base::floatfield);
			std::string filename = "normal_" + ss.str() + ".txt";
			dump_vector(filename.c_str(), data);
		}
		tNow = integrator->execute(
			y, tspan, y0, schemeData,
			line_length_of_chunk, num_of_threads, num_of_gpus,
			delayedDerivMinMax, enable_user_defined_dynamics_on_gpu);
		data = y;
		printf("tNow = %f\n", tNow);



		////!< Visualize
		//levelset::HJI_Grid* g2Dp;
		//beacls::FloatVec data2Dp;

		//g2Dp = helperOC::proj(data2Dp, schemeData->get_grid(), data, beacls::IntegerVec{ 0,1,1 });
		//std::string file = R"(D:\LocalGitCode\SRM\SRMbase1.21\beacls\beacls-master\sources\samples\air3D\x64\Release\)";
		//cv::Mat image_mat;
		//helperOC::visSetIm(image_mat, image_mat, schemeData->get_grid(), data2Dp, std::vector<float>{0.0,1.0,1.0}, beacls::FloatVec(),true, file+std::to_string(tNow)+ ".png",cv::Size());

		////cv::namedWindow(__func__, 0);
		////cv::imshow(__func__, image_mat);
		////cv::waitKey(1);
		//cv::imshow(file, image_mat);
		//cv::waitKey(1);
	}

	if (debug_dump_file) {
		std::stringstream ss;
		ss << std::setprecision(5) << tNow << std::resetiosflags(std::ios_base::floatfield);
		std::string filename = "normal_" + ss.str() + ".txt";
		dump_vector(filename.c_str(), data);
	}



	if (schemeData) delete schemeData;
	if (schemeFunc) delete schemeFunc;
	if (integrator) delete integrator;
	if (hJI_Grid) delete hJI_Grid;
	if (shape) delete shape;
	if (addGhostExtrapolate) delete addGhostExtrapolate;
	if (spatialDerivative) delete spatialDerivative;

}

#pragma once
#include <levelset/levelset.hpp>
#include <cstdint>
#include <vector>
#include <cstddef>
#include <utility>
#include <Core/UVec.hpp>
using namespace std::rel_ops;
enum VelocityType
{
	Velocity_Float,//velocity type is a single float value
	Velocity_UVec//velocity type is a vector value use it to set diffrent velocity in surface
};
class BurnNormalSchemeData : public levelset::SchemeData
{
private:
	FLOAT_TYPE velocity = 10.0;
	beacls::UVec velocity_uvec;
	FLOAT_TYPE input;
	VelocityType velocity_type;
public:
	BurnNormalSchemeData() :SchemeData(), velocity(10.0), input(0.0)
	{
		velocity = 10.0;
	}
	~BurnNormalSchemeData(){}
	bool operator==(const BurnNormalSchemeData& rhs) const;
	bool operator==(const SchemeData& rhs)const;
	BurnNormalSchemeData* clone()const
	{
		return new BurnNormalSchemeData(*this);
	}

	bool hamFunc(
		beacls::UVec& hamValue_uvec,
		const FLOAT_TYPE t,
		const beacls::UVec& data,
		const std::vector<beacls::UVec>& derivs,
		const size_t begin_index,
		const size_t length
	)const override;
	bool normal3DHamFunc(
		beacls::UVec& hamValue_uvec,
		beacls::UVec& magnitude_uvec,
		beacls::FloatVec& step_bound_invs,
		const std::vector<beacls::UVec>& derivs_l,
		const std::vector<beacls::UVec>& derivs_r,
		const size_t begin_index,
		const size_t length
	) const override;
	void setVelocity(FLOAT_TYPE velocity);
private:
	/**
	 * Disable operator=.
	 */
	BurnNormalSchemeData& operator=(const BurnNormalSchemeData& rhs);
	BurnNormalSchemeData(const BurnNormalSchemeData& rhs):
		SchemeData(rhs),
		velocity(0),
		input(0)
	{

	}
	/**
	 * @brief Determine the upwind direction,Either both sides agree in sign (take direction in which they agree),
	 * or characteristics are converging (take larger magnitude direction)..
	 * 
	 * \param prodL
	 * \param prodR
	 * \param magL
	 * \param magR
	 * \return 
	 */
	int calaFlowLValue(FLOAT_TYPE prodL, FLOAT_TYPE prodR, FLOAT_TYPE magL, FLOAT_TYPE magR) const;
	/**
	 * @brief Determine the upwind direction,Either both sides agree in sign (take direction in which they agree),
	 * or characteristics are converging (take larger magnitude direction).
	 * 
	 * \param prodL
	 * \param prodR
	 * \param magL
	 * \param magR
	 * \return 
	 */
	int calaFlowRValue(FLOAT_TYPE prodL, FLOAT_TYPE prodR, FLOAT_TYPE magL, FLOAT_TYPE magR) const;
};

#include "BurnNormalSchemeData.h"
//#include <algorithm>
//#include <typeinfo>
//#include <macro.hpp>
#include <levelset/Grids/HJI_Grid.hpp>
bool BurnNormalSchemeData::operator==(const BurnNormalSchemeData& rhs) const
{
    if (this == &rhs)
    {
        return true;
    }
    else if (!SchemeData::operator==(rhs))
    {
        return false;
    }
    else if (velocity != rhs.velocity)
    {
        return false;
    }
    else if (input != rhs.input)
    {
        return false;
    }
    else
    {
        return true;
    }
}

bool BurnNormalSchemeData::operator==(const SchemeData& rhs) const
{
    if (this == &rhs)
    {
        return true;
    }
    else if (typeid(*this) != typeid(rhs))
    {
        return false;
    }
    else
    {
        return operator==(dynamic_cast<const BurnNormalSchemeData&>(rhs));
    }
}

bool BurnNormalSchemeData::hamFunc(
    beacls::UVec& hamValue_uvec, 
    const FLOAT_TYPE t, 
    const beacls::UVec&, 
    const std::vector<beacls::UVec>& derivs,
    const size_t begin_index, 
    const size_t length) const
{
    
    return false;
}

bool BurnNormalSchemeData::normal3DHamFunc(
    beacls::UVec& hamValue_uvec, 
    beacls::UVec& magnitude_uvec,
    beacls::FloatVec& step_bound_invs,
    const std::vector<beacls::UVec>& derivs_l, 
    const std::vector<beacls::UVec>& derivs_r,
    const size_t begin_index, const size_t length) const
{
    //printf("tNow = %u\n", derivs_l[0].size());
    beacls::reallocateAsSrc(magnitude_uvec, derivs_l[0]);
    FLOAT_TYPE* magnitudeValue = beacls::UVec_<FLOAT_TYPE>(magnitude_uvec).ptr();

    beacls::reallocateAsSrc(hamValue_uvec, derivs_l[0]);
    FLOAT_TYPE* hamValue = beacls::UVec_<FLOAT_TYPE>(hamValue_uvec).ptr();

    const FLOAT_TYPE* deriv0_l = beacls::UVec_<FLOAT_TYPE>(derivs_l[0]).ptr();
    const FLOAT_TYPE* deriv1_l = beacls::UVec_<FLOAT_TYPE>(derivs_l[1]).ptr();
    const FLOAT_TYPE* deriv2_l = beacls::UVec_<FLOAT_TYPE>(derivs_l[2]).ptr();

    const FLOAT_TYPE* deriv0_r = beacls::UVec_<FLOAT_TYPE>(derivs_r[0]).ptr();
    const FLOAT_TYPE* deriv1_r = beacls::UVec_<FLOAT_TYPE>(derivs_r[1]).ptr();
    const FLOAT_TYPE* deriv2_r = beacls::UVec_<FLOAT_TYPE>(derivs_r[2]).ptr();

    FLOAT_TYPE prodL0;
    FLOAT_TYPE prodL1;
    FLOAT_TYPE prodL2;

    FLOAT_TYPE prodR0;
    FLOAT_TYPE prodR1;
    FLOAT_TYPE prodR2;

    FLOAT_TYPE magL0;
    FLOAT_TYPE magL1;
    FLOAT_TYPE magL2;

    FLOAT_TYPE magR0;
    FLOAT_TYPE magR1;
    FLOAT_TYPE magR2;

    int flowL0;
    int flowL1;
    int flowL2;

    int flowR0;
    int flowR1;
    int flowR2;
    const levelset::HJI_Grid* hji_grid = get_grid();
    auto dimensions = hji_grid->get_num_of_dimensions();

    if (step_bound_invs.size() != hji_grid->get_num_of_dimensions())
    {
        step_bound_invs.resize(hji_grid->get_num_of_dimensions());
    }

    FLOAT_TYPE stepBoundMaxInv0 = 0.0;
    FLOAT_TYPE stepBoundMaxInv1 = 0.0;
    FLOAT_TYPE stepBoundMaxInv2 = 0.0;
    for (size_t i = 0; i < length; ++i)
    {

        /*for (int dimension = 0; dimension < dimensions; ++dimension)
        {

        }*/
        FLOAT_TYPE velocitytep = -1.0;

        FLOAT_TYPE deriv0_l_i = deriv0_l[i];
        FLOAT_TYPE deriv1_l_i = deriv1_l[i];
        FLOAT_TYPE deriv2_l_i = deriv2_l[i];

        FLOAT_TYPE deriv0_r_i = deriv0_r[i];
        FLOAT_TYPE deriv1_r_i = deriv1_r[i];
        FLOAT_TYPE deriv2_r_i = deriv2_r[i];

       /* prodL0 = velocity * deriv0_l_i; prodL1 = velocity * deriv1_l_i; prodL2 = velocity * deriv2_l_i;
        prodR0 = velocity * deriv0_r_i; prodR1 = velocity * deriv1_r_i; prodR2 = velocity * deriv2_r_i;*/
        prodL0 = velocitytep * deriv0_l_i; prodL1 = velocitytep * deriv1_l_i; prodL2 = velocitytep * deriv2_l_i;
        prodR0 = velocitytep * deriv0_r_i; prodR1 = velocitytep * deriv1_r_i; prodR2 = velocitytep * deriv2_r_i;

        magL0 = abs(prodL0); magL1 = abs(prodL1); magL2 = abs(prodL2);
        magR0 = abs(prodR0); magR1 = abs(prodR1); magR2 = abs(prodR2);

        flowL0 = calaFlowLValue(prodL0, prodR0, magL0, magR0);
        flowL1 = calaFlowLValue(prodL1, prodR1, magL1, magR1);
        flowL2 = calaFlowLValue(prodL2, prodR2, magL2, magR2);

        flowR0 = calaFlowRValue(prodL0, prodR0, magL0, magR0);
        flowR1 = calaFlowRValue(prodL1, prodR1, magL1, magR1);
        flowR2 = calaFlowRValue(prodL2, prodR2, magL2, magR2);

        //Now we know the upwind direction, add its contribution to \ | \grad \phi\ | .
        FLOAT_TYPE magnitude =
            pow(deriv0_l_i, 2) * flowL0 + pow(deriv0_r_i, 2) * flowR0 +
            pow(deriv1_l_i, 2) * flowL1 + pow(deriv1_r_i, 2) * flowR1 +
            pow(deriv2_l_i, 2) * flowL2 + pow(deriv2_r_i, 2) * flowR2;

        // Finally, calculate speed* \ | \grad \phi\ |
        
        magnitude = sqrtf(magnitude);
        hamValue[i] = velocitytep * magnitude;
        //magnitudeValue[begin_index + i] = magnitude;

        FLOAT_TYPE effectiveVelocity0 = magL0 * flowL0 + magR0 * flowR0;
        FLOAT_TYPE effectiveVelocity1 = magL1 * flowL1 + magR1 * flowR1;
        FLOAT_TYPE effectiveVelocity2 = magL2 * flowL2 + magR2 * flowR2;

        FLOAT_TYPE stepBoundInv0 = hji_grid->get_dxInv(0) * effectiveVelocity0;
        FLOAT_TYPE stepBoundInv1 = hji_grid->get_dxInv(1) * effectiveVelocity1;
        FLOAT_TYPE stepBoundInv2 = hji_grid->get_dxInv(2) * effectiveVelocity2;
        if (magnitude > 0.0)
        {
            stepBoundMaxInv0 = stepBoundInv0 / magnitude > stepBoundMaxInv0 ? stepBoundInv0 / magnitude : stepBoundInv0;
            stepBoundMaxInv1 = stepBoundInv1 / magnitude > stepBoundMaxInv1 ? stepBoundInv1 / magnitude : stepBoundInv1;
            stepBoundMaxInv2 = stepBoundInv2 / magnitude > stepBoundMaxInv2 ? stepBoundInv2 / magnitude : stepBoundInv2;
        }
    }
    step_bound_invs[0] = stepBoundMaxInv0;
    step_bound_invs[1] = stepBoundMaxInv1;
    step_bound_invs[2] = stepBoundMaxInv2;

    return true;
}

void BurnNormalSchemeData::setVelocity(FLOAT_TYPE velocity)
{
    this->velocity = velocity;
}

int BurnNormalSchemeData::calaFlowLValue(FLOAT_TYPE prodL, FLOAT_TYPE prodR, FLOAT_TYPE magL, FLOAT_TYPE magR) const
{
    if ((prodL >= 0.0 && prodR >= 0.0) || (prodL >= 0.0 && prodR <= 0.0) && (magL >= magR))
    {
        return 1;
    }
    return 0;
}

int BurnNormalSchemeData::calaFlowRValue(FLOAT_TYPE prodL, FLOAT_TYPE prodR, FLOAT_TYPE magL, FLOAT_TYPE magR) const
{
    if ((prodL <= 0.0 && prodR <= 0.0) || (prodL >= 0.0 && prodR <= 0.0) && (magL < magR))
    {
        return 1;
    }
    return 0;
}

#ifndef __TermNormal_hpp__
#define __TermNormal_hpp__


//! Prefix to generate Visual C++ DLL
#ifdef _MAKE_VC_DLL
#define PREFIX_VC_DLL __declspec(dllexport)

//! Don't add prefix, except dll generating
#else              
#define PREFIX_VC_DLL
#endif

#include <cstdint>
#include <vector>
#include <utility>
using namespace std::rel_ops;

#include <typedef.hpp>
#include <levelset/ExplicitIntegration/Terms/Term.hpp>


namespace levelset {
	class SchemeData;
	class TermNormal_impl;
	class TermNormal:public Term
	{
	public:
		/**
		 * motion in the normal direction in an HJ PDE with upwinding.
		 * 
		 * @param [in]			schemeData
		 * @param [in]			type
		 * @return 
		 */
		PREFIX_VC_DLL TermNormal(const SchemeData* schemeData, const beacls::UVecType type = beacls::UVecType_Vector);
		/**
		 * @brief destructor
		 * 
		 */
		~TermNormal();
		/**
		 * .
		 * 
		 * @param [out]			ydot_ite
		 * @param [in]			step_bound_invs
		 * @param [in]			t
		 * @param [in]			y
		 * @param [in]			derivMins
		 * @param [in]			derivMaxs
		 * @param [in]			schemeData
		 * @param [in]			loop_begin
		 * @param [in]			loop_length
		 * @param [in]			num_of_slices
		 * @param [in]			enable_use_defined_dynamics_on_gpu
		 * @param [in]			updataDerivMinMax
		 * @return true: Succeeded  false :Failed Dissipation may be required global derivMins/derivMaxs.
		Reduce partial derivMins/derivMaxs to global derivMins/derivMax, then execute again.
		 */
		bool execute(beacls::FloatVec::iterator ydot_ite,
			beacls::FloatVec& step_bound_invs,
			const FLOAT_TYPE t,
			const beacls::FloatVec& y,
			std::vector<beacls::FloatVec>& derivMins,
			std::vector<beacls::FloatVec>& derivMaxs,
			const SchemeData* schemeData,
			const size_t loop_begin,
			const size_t loop_length,
			const size_t num_of_slices,
			const bool enable_use_defined_dynamics_on_gpu,
			const bool updataDerivMinMax
		) const;
		/**
		 * .
		 * 
		 * @param schemeData
		 * @return true:Succeeded false:Failed
		 */
		bool synchronize(const SchemeData* schemeData) const;
		bool operator==(const TermNormal& rhs) const;
		bool operator==(const Term& rhs) const;
		TermNormal* clone() const;
	private:
		TermNormal_impl* pimpl;

		/**
		 * .
		 * @overload
		 * Disable operator =
		 * @param rhs
		 * @return 
		 */
		TermNormal& operator=(const TermNormal& rhs);
		/**
		 * .
		 * @overload
		 * Disable copy constructor
		 * @param rhs
		 */
		TermNormal(const TermNormal& rhs);

	};


}
#endif	/* __TermNormal_hpp__ */
###TermNormal_impl.hpp
#ifndef __TermNormal_impl_hpp__
#define __TermNormal_impl_hpp__

#include <cstdint>
#include <vector>
#include <Core/UVec.hpp>
#include <typedef.hpp>
namespace levelset {
	class SchemeData;
	class CacheTag;

	class TermNormal_impl
	{
	private:
		size_t first_dimension_loop_size;
		size_t num_of_dimensions;
		std::vector<beacls::UVec> deriv_l_uvecs;
		std::vector<beacls::UVec> deriv_r_uvecs;
		std::vector<beacls::UVec> deriv_c_uvecs;
		std::vector<beacls::UVec> deriv_c_cpu_uvecs;

		std::vector<beacls::UVec> x_uvecs;
		//std::vector<beacls::UVec> deriv_max_uvecs;
		//std::vector<beacls::UVec> deriv_min_uvecs;
		beacls::UVec ydot_cuda_uvec;
		/*maybe should delete*/
		beacls::UVec ham_uvec;
		beacls::UVec ham_cpu_uvec;
		beacls::UVec diss_uvec;
		/****************************************/


		beacls::UVec flowR_uvec;
		beacls::UVec flowL_uvec;
		beacls::UVec magL_uvec;
		beacls::UVec magR_uvec;
		beacls::UVec magnitude_uvec;



		beacls::UVecType type;
		CacheTag* cacheTag;


	public:
		TermNormal_impl(
			const SchemeData* schemeData,
			const beacls::UVecType type);
		~TermNormal_impl();


		bool execute(
			beacls::FloatVec::iterator ydot_ite,
			beacls::FloatVec& step_bound_invs,
			const FLOAT_TYPE t,
			const beacls::FloatVec& y,
			std::vector<beacls::FloatVec >& derivMins,
			std::vector<beacls::FloatVec >& derivMaxs,
			const SchemeData* schemeData,
			const size_t loop_begin,
			const size_t loop_length,
			const size_t num_of_slices,
			const bool enable_user_defined_dynamics_on_gpu,
			const bool updateDerivMinMax
		);
		bool synchronize(const SchemeData* schemeData) const;
		bool operator==(const TermNormal_impl& rhs) const {
			if (this == &rhs) return true;
			else if (type != rhs.type)return false;
			else if (num_of_dimensions != rhs.num_of_dimensions)return false;
			return true;
		}
		TermNormal_impl* clone() const {
			return new TermNormal_impl(*this);
		};
	private:
		/** @overload
		Disable operator=
		*/
		TermNormal_impl& operator=(const TermNormal_impl& rhs);
		/** @overload
		Disable copy constructor
		*/
		TermNormal_impl(const TermNormal_impl& rhs);


	};
#endif




}

#include<levelset/ExplicitIntegration/Terms/TermNormal.hpp>

#include "TermNormal_impl.hpp"
#include <levelset/Grids/HJI_Grid.hpp>
//#include "TermNormal_cuda.hpp"//need to impl
#include <levelset/ExplicitIntegration/SchemeData.hpp>
#include <levelset/SpatialDerivative/SpatialDerivative.hpp>
#include <Core/CacheTag.hpp>

#include <algorithm>
#include <functional>
#include <typeinfo>
using namespace levelset;
using namespace std;
levelset::TermNormal_impl::TermNormal_impl(const SchemeData* schemeData, const beacls::UVecType type) :
	first_dimension_loop_size(schemeData->get_grid()->get_N(0)),
	num_of_dimensions(schemeData->get_grid()->get_num_of_dimensions()),
	deriv_l_uvecs(num_of_dimensions),
	deriv_r_uvecs(num_of_dimensions),
	type(type),
	cacheTag(new levelset::CacheTag())
{
}
levelset::TermNormal_impl::~TermNormal_impl()
{
	if (cacheTag)
	{
		delete cacheTag;
	}
}
levelset::TermNormal_impl::TermNormal_impl(const TermNormal_impl& rhs):
	first_dimension_loop_size(rhs.first_dimension_loop_size),
	num_of_dimensions(rhs.num_of_dimensions),
	type(rhs.type),
	cacheTag(new levelset::CacheTag())
{
	deriv_l_uvecs.resize(rhs.deriv_l_uvecs.size());
	deriv_r_uvecs.resize(rhs.deriv_r_uvecs.size());
	x_uvecs.resize(rhs.x_uvecs.size());
}


bool levelset::TermNormal_impl::execute(
	beacls::FloatVec::iterator ydot_ite,
	beacls::FloatVec& step_bound_invs,
	const FLOAT_TYPE t,
	const beacls::FloatVec& y,
	std::vector<beacls::FloatVec>& derivMins,
	std::vector<beacls::FloatVec>& derivMaxs,
	const SchemeData* schemeData,
	const size_t loop_begin,
	const size_t loop_length,
	const size_t num_of_slices, 
	const bool enable_user_defined_dynamics_on_gpu, 
	const bool updateDerivMinMax)
{
	const HJI_Grid* grid = schemeData->get_grid();
	SpatialDerivative* spatialDerivative = schemeData->get_spatialDerivative();
	beacls::UVecDepth depth = beacls::type_to_depth<FLOAT_TYPE>();
	beacls::UVec y_uvec(y, beacls::UVecType_Vector, false);

	size_t f_d_l_size = first_dimension_loop_size;//first dimension loop size
	size_t grid_length = num_of_slices * loop_length * f_d_l_size;//


	//size_t grid_length = grid->get_numel();
	// printf("(num_of_slices, loop_length, f_d_l_size) = (%zu, %zu, %zu)\n", 
	// 	num_of_slices, loop_length, f_d_l_size);
  // printf("grid_length = %zu\n", grid_length);
	size_t slice_length = loop_length * f_d_l_size;
	if (deriv_l_uvecs.size() != num_of_dimensions)
	{
		deriv_l_uvecs.resize(num_of_dimensions);
	}
	if (deriv_r_uvecs.size() != num_of_dimensions)
	{
		deriv_r_uvecs.resize(num_of_dimensions);
	}


	//init deriv value vector
	for_each(deriv_l_uvecs.begin(), deriv_l_uvecs.end(), ([depth, grid_length, this](auto& rhs) {
		if (rhs.type() != type)
		{
			rhs = beacls::UVec(depth, type, grid_length);
		}
		else if (rhs.size() != grid_length)
		{
			rhs.resize(grid_length);
		}
		}));
	for_each(deriv_r_uvecs.begin(), deriv_r_uvecs.end(), ([depth, grid_length, this](auto& rhs) {
		if (rhs.type() != type)
		{
			rhs = beacls::UVec(depth, type, grid_length);
		}
		else if (rhs.size() != grid_length)
		{
			rhs.resize(grid_length);
		}
		}));
	//!< cala magnitude temp value

	if (flowL_uvec.type() != type)
	{
		flowL_uvec = beacls::UVec(depth, type, grid_length);
	}
	else if (flowL_uvec.size() != grid_length)
	{
		flowL_uvec.resize(grid_length);
	}

	if (flowR_uvec.type() != type)
	{
		flowR_uvec = beacls::UVec(depth, type, grid_length);
	}
	else if (flowR_uvec.size() != grid_length)
	{
		flowR_uvec.resize(grid_length);
	}

	if (magL_uvec.type() != type)
	{
		magL_uvec = beacls::UVec(depth, type, grid_length);
	}
	else if (magL_uvec.size() != grid_length)
	{
		magL_uvec.resize(grid_length);
	}

	if (magR_uvec.type() != type)
	{
		magR_uvec = beacls::UVec(depth, type, grid_length);
	}
	else if (magR_uvec.size() != grid_length)
	{
		magR_uvec.resize(grid_length);
	}


	if (magnitude_uvec.type() != type)
	{
		magnitude_uvec = beacls::UVec(depth, type, grid_length);
	}
	else if (magnitude_uvec.size() != grid_length)
	{
		magnitude_uvec.resize(grid_length);
	}


	
	size_t src_index_term = loop_begin * f_d_l_size;

	if (!cacheTag->check_tag(t, loop_begin, slice_length * num_of_slices))
	{
		//!< Copy xs to Cuda memory asynchronously in spatial derivative functions
		x_uvecs.resize(num_of_dimensions);
		//if use gpu init cuda vector
		if (enable_user_defined_dynamics_on_gpu && (type == beacls::UVecType_Cuda))
		{
			for (size_t dimension = 0; dimension < num_of_dimensions; ++dimension) {
				if (x_uvecs[dimension].type() != beacls::UVecType_Cuda) x_uvecs[dimension] = beacls::UVec(depth, beacls::UVecType_Cuda, grid_length);
				else x_uvecs[dimension].resize(grid_length);
			}
		}

		for (size_t index = 0; index < num_of_dimensions; ++index)
		{
			//!< To optimize asynchronous execution, calculate from heavy dimension (0, 2, 3 ... 1);
			const size_t dimension = (index == 0) ? index : (index == num_of_dimensions - 1) ? 1 : index + 1;
			beacls::UVec& deriv_l_uvec = deriv_l_uvecs[dimension];
			beacls::UVec& deriv_r_uvec = deriv_r_uvecs[dimension];

			spatialDerivative->execute(
				deriv_l_uvec,
				deriv_r_uvec,
				grid,
				y.data(),
				dimension,
				false,
				loop_begin,
				slice_length,
				num_of_slices);
			const beacls::FloatVec& xs = grid->get_xs(dimension);
			beacls::copyHostPtrToUVecAsync(x_uvecs[dimension], xs.data() + src_index_term, grid_length);

		}
		for (size_t dimension = 0; dimension < num_of_dimensions; ++dimension)
		{
			beacls::UVec x_uvecs_dim = x_uvecs[dimension];
			synchronizeUVec(x_uvecs_dim);
		}
		cacheTag->set_tag(t, loop_begin, slice_length * num_of_dimensions);

	}
	
	schemeData->normal3DHamFunc(
		ham_uvec,
		magnitude_uvec,
		step_bound_invs,
		deriv_l_uvecs,
		deriv_r_uvecs,
		src_index_term,
		grid_length);
	const beacls::FloatVec* ham_uvec_ptr = beacls::UVec_<FLOAT_TYPE>(ham_uvec).vec();
	std::copy(ham_uvec_ptr->cbegin(), ham_uvec_ptr->cend(), ydot_ite);
	return true;
}
bool levelset::TermNormal_impl::synchronize(const SchemeData* schemeData) const
{
	beacls::synchronizeUVec(ydot_cuda_uvec);
	return true;
}

levelset::TermNormal::TermNormal(const SchemeData* schemeData, const beacls::UVecType type)
{
	pimpl = new TermNormal_impl(schemeData, type);
}

levelset::TermNormal::~TermNormal()
{
	if (pimpl)
	{
		delete pimpl;
	}
}

bool levelset::TermNormal::execute(
	beacls::FloatVec::iterator ydot_ite,
	beacls::FloatVec& step_bound_invs, 
	const FLOAT_TYPE t, 
	const beacls::FloatVec& y, 
	std::vector<beacls::FloatVec>& derivMins, 
	std::vector<beacls::FloatVec>& derivMaxs, 
	const SchemeData* schemeData,
	const size_t loop_begin, 
	const size_t loop_length, 
	const size_t num_of_slices, 
	const bool enable_use_defined_dynamics_on_gpu, 
	const bool updataDerivMinMax) const
{

	if (pimpl)
	{
		bool ret = pimpl->execute(ydot_ite, step_bound_invs, t, y, derivMins, derivMaxs, schemeData, loop_begin, loop_length, num_of_slices, enable_use_defined_dynamics_on_gpu, updataDerivMinMax);
		return ret;
	}
	return false;
}

bool levelset::TermNormal::synchronize(const SchemeData* schemeData) const
{
	if (pimpl)
	{
		return pimpl->synchronize(schemeData);
	}
	return false;
}

bool levelset::TermNormal::operator==(const TermNormal& rhs) const
{
	if (this == &rhs) return true;
	else if (!pimpl)
	{
		if (!rhs.pimpl)
		{
			return true;
		}
		else
		{
			return false;
		}
	}
	else
	{
		if (!rhs.pimpl)
		{
			return false;
		}
		else
		{
			//object pointer address is equal
			if (pimpl == rhs.pimpl)
			{
				return true;
			}
			else if (*pimpl == *rhs.pimpl)//object is equal
			{
				return true;
			}
			else
			{
				return false;
			}
		}
	}
	return false;
}

bool levelset::TermNormal::operator==(const Term& rhs) const
{
	if (this == &rhs)
	{
		return true;
	}
	else if (typeid(*this) != typeid(rhs))//types are consistent
	{
		return false;
	}
	else
	{
		return operator==(dynamic_cast<const TermNormal&>(rhs));
	}
}



TermNormal* levelset::TermNormal::clone() const
{
	return new TermNormal(*this);
}

TermNormal& levelset::TermNormal::operator=(const TermNormal& rhs)
{
	if (pimpl)
	{
		delete pimpl;
		pimpl = rhs.pimpl->clone();
	}
	else
	{
		pimpl = rhs.pimpl->clone();
	}
}

levelset::TermNormal::TermNormal(const TermNormal& rhs):pimpl(rhs.pimpl->clone())
{

}

/usr/local/opt/llvm/bin/clang++ -shared -o liblevelset.so AddGhostExtrapolate.o AddGhostPeriodic.o ArtificialDissipationGLF.o interpn.o CacheTag.o UVec.o UVec_cuda_dummy.o SchemeData.o odeCFL_SubStep.o OdeCFL_CommandQueue.o OdeCFL_OneSlice.o OdeCFL_Worker.o OdeCFL1.o OdeCFL2.o OdeCFL3.o TermLaxFriedrichs.o TermRestrictUpdate.o HJI_Grid.o ShapeCylinder.o ShapeHyperplaneByPoint.o ShapeRectangleByCenter.o ShapeRectangleByCorner.o ShapeSphere.o UpwindFirst.o UpwindFirstFirst.o UpwindFirstENO2.o UpwindFirstENO3.o UpwindFirstENO3a.o UpwindFirstENO3b.o UpwindFirstENO3aHelper.o UpwindFirstENO3bHelper.o UpwindFirstWENO5.o UpwindFirstWENO5a.o UpwindFirstWENO5b.o ArtificialDissipationGLF_cuda_dummy.o TermLaxFriedrichs_cuda_dummy.o UpwindFirstENO3aHelper_cuda_dummy.o UpwindFirstWENO5a_cuda_dummy.o  -L/usr/local/lib -lstdc++ -lc++ -g -L/Users/davidfridovichkeil/Documents/Developer/beacls/sources/modules/levelset/../../../builds/x86_64/lib -lmatio -lz -lhdf5 -L/usr/local/lib -lopencv_shape -lopencv_stitching -lopencv_objdetect -lopencv_superres -lopencv_videostab -lopencv_calib3d -lopencv_features2d -lopencv_highgui -lopencv_videoio -lopencv_imgcodecs -lopencv_video -lopencv_photo -lopencv_ml -lopencv_imgproc -lopencv_flann -lopencv_viz -lopencv_core -L/usr/local/opt/llvm/lib -liomp5 -lpthread
Undefined symbols for architecture x86_64:
  "std::terminate()", referenced from:
      ___clang_call_terminate in interpn.o
      ___clang_call_terminate in UVec.o
      ___clang_call_terminate in UVec_cuda_dummy.o
      ___clang_call_terminate in odeCFL_SubStep.o
      ___clang_call_terminate in OdeCFL_CommandQueue.o
      levelset::OdeCFL_Worker_impl::run() in OdeCFL_Worker.o
      ___clang_call_terminate in OdeCFL_Worker.o
      ...
  "___cxa_bad_cast", referenced from:
      levelset::AddGhostExtrapolate::operator==(levelset::BoundaryCondition const&) const in AddGhostExtrapolate.o
      levelset::AddGhostPeriodic::operator==(levelset::BoundaryCondition const&) const in AddGhostPeriodic.o
      levelset::ArtificialDissipationGLF::operator==(levelset::Dissipation const&) const in ArtificialDissipationGLF.o
      levelset::TermLaxFriedrichs::operator==(levelset::Term const&) const in TermLaxFriedrichs.o
      levelset::TermRestrictUpdate::operator==(levelset::Term const&) const in TermRestrictUpdate.o
      levelset::UpwindFirstFirst::operator==(levelset::SpatialDerivative const&) const in UpwindFirstFirst.o
      levelset::UpwindFirstENO2::operator==(levelset::SpatialDerivative const&) const in UpwindFirstENO2.o
      ...
  "___cxa_bad_typeid", referenced from:
      levelset::AddGhostExtrapolate::operator==(levelset::BoundaryCondition const&) const in AddGhostExtrapolate.o
      levelset::AddGhostPeriodic::operator==(levelset::BoundaryCondition const&) const in AddGhostPeriodic.o
      levelset::ArtificialDissipationGLF::operator==(levelset::Dissipation const&) const in ArtificialDissipationGLF.o
      levelset::TermLaxFriedrichs::operator==(levelset::Term const&) const in TermLaxFriedrichs.o
      levelset::TermRestrictUpdate::operator==(levelset::Term const&) const in TermRestrictUpdate.o
      levelset::HJI_Grid_impl::save_grid(std::__1::basic_string<char, std::__1::char_traits<char>, std::__1::allocator<char> > const&, beacls::MatFStream*, beacls::MatVariable*, unsigned long, bool) const in HJI_Grid.o
      levelset::UpwindFirstFirst::operator==(levelset::SpatialDerivative const&) const in UpwindFirstFirst.o
      ...
ld: symbol(s) not found for architecture x86_64
clang-3.6: error: linker command failed with exit code 1 (use -v to see invocation)
make[2]: *** [liblevelset.so] Error 1
make[1]: *** [levelset] Error 2
make: *** [modules] Error 2