diff --git a/CMakeLists.txt b/CMakeLists.txt index 6f6ebf50..ebe9f7ca 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -22,6 +22,7 @@ option(LOCAL_DEPS "LOCAL_DEPS" OFF) option(TURBO_VERBOSE "Compile with verbose output" ON) option(WITH_ASAN "Compile with the Address Sanitizer to check for memory corruption errors" OFF) option(WITH_XCSP3PARSER "Add support for parsing XCSP3 .xml files" ON) +option(WITH_NNV "Add support for neural network verification .onnx and .vnnlib" ON) option(NO_CONCURRENT_MANAGED_MEMORY "Add support for platform not supporting concurrent managed access to memory on GPUs (use pinned memory instead)." OFF) option(TURBO_BUILD_DOC "Build the documentation of Turbo." OFF) option(TURBO_IPC_ABSTRACT_DOMAIN "Use the interval propagator completion (PC) abstract domain instead of PIR." OFF) @@ -105,11 +106,31 @@ if(WITH_XCSP3PARSER) target_compile_definitions(turbo PRIVATE WITH_XCSP3PARSER) endif() +if(WITH_NNV) + target_compile_definitions(turbo PRIVATE WITH_NNV) + if(GPU) + get_filename_component(ONNXRUNTIME_ROOT + "${CMAKE_CURRENT_LIST_DIR}/../lala-parsing/include/lala/onnxruntime-linux-x64-gpu-1.19.2" + ABSOLUTE + ) + # set(ONNXRUNTIME_ROOT "${CMAKE_CURRENT_LIST_DIR}/../lala-parsing/include/lala/onnxruntime-linux-x64-gpu-1.26.0") + include_directories(${ONNXRUNTIME_ROOT}/include) + add_library(ort_gpu SHARED IMPORTED) + set_target_properties(ort_gpu PROPERTIES + IMPORTED_LOCATION "${ONNXRUNTIME_ROOT}/lib/libonnxruntime.so" + ) + endif() +endif() + if(NO_CONCURRENT_MANAGED_MEMORY) target_compile_definitions(turbo PRIVATE NO_CONCURRENT_MANAGED_MEMORY) endif() -target_link_libraries(turbo PRIVATE lala_parsing lala_pc lala_power) +if(GPU AND WITH_NNV) + target_link_libraries(turbo PRIVATE lala_parsing lala_pc lala_power ort_gpu cuda cudart) +else() + target_link_libraries(turbo PRIVATE lala_parsing lala_pc lala_power) +endif() target_link_options(turbo PRIVATE $<$,$>:-fsanitize=address;-static-libasan> $<$,$>:-fsanitize=address> diff --git a/CMakePresets.json b/CMakePresets.json index 328e5c51..5302c0fc 100644 --- a/CMakePresets.json +++ b/CMakePresets.json @@ -1,288 +1,294 @@ { - "version": 6, - "cmakeMinimumRequired": { - "major": 3, - "minor": 25, - "patch": 0 - }, - "configurePresets": [ - { - "name": "default", - "hidden": true, - "binaryDir": "${sourceDir}/build/${presetName}" + "version": 6, + "cmakeMinimumRequired": { + "major": 3, + "minor": 25, + "patch": 0 }, - { - "name": "default-gpu", - "hidden": true, - "inherits": "default", - "cacheVariables": { - "GPU": { - "type": "BOOL", - "value": "ON" + "configurePresets": [ + { + "name": "default", + "hidden": true, + "binaryDir": "${sourceDir}/build/${presetName}" }, - "REDUCE_PTX_SIZE": { - "type": "BOOL", - "value": "OFF" - } - } - }, - { - "name": "default-cpu", - "hidden": true, - "inherits": "default", - "cacheVariables": { - "GPU": { - "type": "BOOL", - "value": "OFF" - } - } - }, - { - "name": "debug", - "hidden": true, - "cacheVariables": { - "CMAKE_BUILD_TYPE": "Debug" - } - }, - { - "name": "release", - "hidden": true, - "cacheVariables": { - "CMAKE_BUILD_TYPE": "Release" - } - }, - { - "name": "verbose", - "hidden": true, - "cacheVariables": { - "CMAKE_VERBOSE_MAKEFILE": { - "type": "BOOL", - "value": "ON" - } - } - }, - { - "name": "local", - "hidden": true, - "cacheVariables": { - "LOCAL_DEPS": { - "type": "BOOL", - "value": "ON" - } - } - }, - { - "name": "gpu-debug", - "displayName": "GPU Debug configuration", - "description": "Build the project with CUDA NVCC compiler.", - "inherits": [ - "default-gpu", - "debug", - "verbose" - ] - }, - { - "name": "gpu-release", - "displayName": "GPU Release configuration", - "description": "Build the project with CUDA NVCC compiler.", - "inherits": [ - "default-gpu", - "release" - ] - }, - { - "name": "cpu-debug", - "displayName": "CPU Debug configuration", - "description": "Build the project with a CPU compiler such as GCC.", - "inherits": [ - "default-cpu", - "debug" - ] - }, - { - "name": "cpu-release", - "displayName": "CPU Release configuration", - "description": "Build the project with a CPU compiler such as GCC.", - "inherits": [ - "default-cpu", - "release" - ] - }, - { - "name": "gpu-debug-local", - "displayName": "GPU Debug Local Configuration", - "description": "Build the project with CUDA NVCC compiler and local dependencies.", - "inherits": [ - "gpu-debug", - "local" - ] - }, - { - "name": "gpu-release-local", - "displayName": "GPU Release Local Configuration", - "description": "Build the project with CUDA NVCC compiler and local dependencies.", - "inherits": [ - "gpu-release", - "local" - ] - }, - { - "name": "cpu-debug-local", - "displayName": "CPU Debug Local Configuration", - "description": "Build the project with a CPU compiler such as GCC and local dependencies.", - "inherits": [ - "cpu-debug", - "local" - ] - }, - { - "name": "cpu-release-local", - "displayName": "CPU Release Local Configuration", - "description": "Build the project with a CPU compiler such as GCC and local dependencies.", - "inherits": [ - "cpu-release", - "local" - ] - } - ], - "buildPresets": [ - { - "name": "gpu-debug", - "configurePreset": "gpu-debug" - }, - { - "name": "gpu-release", - "configurePreset": "gpu-release" - }, - { - "name": "cpu-debug", - "configurePreset": "cpu-debug" - }, - { - "name": "cpu-release", - "configurePreset": "cpu-release" - }, - { - "name": "gpu-debug-local", - "configurePreset": "gpu-debug-local" - }, - { - "name": "gpu-release-local", - "configurePreset": "gpu-release-local" - }, - { - "name": "cpu-debug-local", - "configurePreset": "cpu-debug-local" - }, - { - "name": "cpu-release-local", - "configurePreset": "cpu-release-local" - } - ], - "workflowPresets": [ - { - "name": "cpu-debug", - "steps": [ - { - "type": "configure", - "name": "cpu-debug" + { + "name": "default-gpu", + "hidden": true, + "inherits": "default", + "cacheVariables": { + "GPU": { + "type": "BOOL", + "value": "ON" + }, + "REDUCE_PTX_SIZE": { + "type": "BOOL", + "value": "OFF" + } + } }, { - "type": "build", - "name": "cpu-debug" - } - ] - }, - { - "name": "cpu-release", - "steps": [ + "name": "default-cpu", + "hidden": true, + "inherits": "default", + "cacheVariables": { + "GPU": { + "type": "BOOL", + "value": "OFF" + } + } + }, { - "type": "configure", - "name": "cpu-release" + "name": "debug", + "hidden": true, + "cacheVariables": { + "CMAKE_BUILD_TYPE": "Debug" + } }, { - "type": "build", - "name": "cpu-release" - } - ] - }, - { - "name": "gpu-debug", - "steps": [ + "name": "release", + "hidden": true, + "cacheVariables": { + "CMAKE_BUILD_TYPE": "Release" + } + }, { - "type": "configure", - "name": "gpu-debug" + "name": "verbose", + "hidden": true, + "cacheVariables": { + "CMAKE_VERBOSE_MAKEFILE": { + "type": "BOOL", + "value": "ON" + } + } }, { - "type": "build", - "name": "gpu-debug" - } - ] - }, - { - "name": "gpu-release", - "steps": [ + "name": "local", + "hidden": true, + "cacheVariables": { + "LOCAL_DEPS": { + "type": "BOOL", + "value": "ON" + } + } + }, { - "type": "configure", - "name": "gpu-release" + "name": "gpu-debug", + "displayName": "GPU Debug configuration", + "description": "Build the project with CUDA NVCC compiler.", + "inherits": [ + "default-gpu", + "debug", + "verbose" + ] }, { - "type": "build", - "name": "gpu-release" - } - ] - }, - { - "name": "cpu-debug-local", - "steps": [ + "name": "gpu-release", + "displayName": "GPU Release configuration", + "description": "Build the project with CUDA NVCC compiler.", + "inherits": [ + "default-gpu", + "release" + ] + }, { - "type": "configure", - "name": "cpu-debug-local" + "name": "cpu-debug", + "displayName": "CPU Debug configuration", + "description": "Build the project with a CPU compiler such as GCC.", + "inherits": [ + "default-cpu", + "debug" + ] }, { - "type": "build", - "name": "cpu-debug-local" - } - ] - }, - { - "name": "cpu-release-local", - "steps": [ + "name": "cpu-release", + "displayName": "CPU Release configuration", + "description": "Build the project with a CPU compiler such as GCC.", + "inherits": [ + "default-cpu", + "release" + ] + }, { - "type": "configure", - "name": "cpu-release-local" + "name": "gpu-debug-local", + "displayName": "GPU Debug Local Configuration", + "description": "Build the project with CUDA NVCC compiler and local dependencies.", + "inherits": [ + "gpu-debug", + "local" + ] }, { - "type": "build", - "name": "cpu-release-local" + "name": "gpu-release-local", + "displayName": "GPU Release Local Configuration", + "description": "Build the project with CUDA NVCC compiler and local dependencies.", + "inherits": [ + "gpu-release", + "local" + ] + }, + { + "name": "cpu-debug-local", + "displayName": "CPU Debug Local Configuration", + "description": "Build the project with a CPU compiler such as GCC and local dependencies.", + "inherits": [ + "cpu-debug", + "local" + ] + }, + { + "name": "cpu-release-local", + "displayName": "CPU Release Local Configuration", + "description": "Build the project with a CPU compiler such as GCC and local dependencies.", + "inherits": [ + "cpu-release", + "local" + ] } - ] - }, - { - "name": "gpu-debug-local", - "steps": [ + ], + "buildPresets": [ + { + "name": "gpu-debug", + "configurePreset": "gpu-debug" + }, + { + "name": "gpu-release", + "configurePreset": "gpu-release" + }, + { + "name": "cpu-debug", + "configurePreset": "cpu-debug" + }, + { + "name": "cpu-release", + "configurePreset": "cpu-release" + }, + { + "name": "gpu-debug-local", + "configurePreset": "gpu-debug-local" + }, + { + "name": "gpu-release-local", + "configurePreset": "gpu-release-local" + }, + { + "name": "cpu-debug-local", + "configurePreset": "cpu-debug-local" + }, { - "type": "configure", - "name": "gpu-debug-local" + "name": "cpu-release-local", + "configurePreset": "cpu-release-local" }, { - "type": "build", - "name": "gpu-debug-local" + "name": "Enter", + "description": "", + "displayName": "", + "configurePreset": "cpu-release-local" } - ] - }, - { - "name": "gpu-release-local", - "steps": [ + ], + "workflowPresets": [ + { + "name": "cpu-debug", + "steps": [ + { + "type": "configure", + "name": "cpu-debug" + }, + { + "type": "build", + "name": "cpu-debug" + } + ] + }, + { + "name": "cpu-release", + "steps": [ + { + "type": "configure", + "name": "cpu-release" + }, + { + "type": "build", + "name": "cpu-release" + } + ] + }, + { + "name": "gpu-debug", + "steps": [ + { + "type": "configure", + "name": "gpu-debug" + }, + { + "type": "build", + "name": "gpu-debug" + } + ] + }, + { + "name": "gpu-release", + "steps": [ + { + "type": "configure", + "name": "gpu-release" + }, + { + "type": "build", + "name": "gpu-release" + } + ] + }, + { + "name": "cpu-debug-local", + "steps": [ + { + "type": "configure", + "name": "cpu-debug-local" + }, + { + "type": "build", + "name": "cpu-debug-local" + } + ] + }, + { + "name": "cpu-release-local", + "steps": [ + { + "type": "configure", + "name": "cpu-release-local" + }, + { + "type": "build", + "name": "cpu-release-local" + } + ] + }, { - "type": "configure", - "name": "gpu-release-local" + "name": "gpu-debug-local", + "steps": [ + { + "type": "configure", + "name": "gpu-debug-local" + }, + { + "type": "build", + "name": "gpu-debug-local" + } + ] }, { - "type": "build", - "name": "gpu-release-local" + "name": "gpu-release-local", + "steps": [ + { + "type": "configure", + "name": "gpu-release-local" + }, + { + "type": "build", + "name": "gpu-release-local" + } + ] } - ] - } - ] + ] } \ No newline at end of file diff --git a/benchmarks/Brown-05.fzn b/benchmarks/Brown-05.fzn new file mode 100644 index 00000000..c653aec7 --- /dev/null +++ b/benchmarks/Brown-05.fzn @@ -0,0 +1,27 @@ +% Generated by MiniZinc 2.9.5 +% Solver library: /snap/minizinc/1222/share/minizinc/gecode +% Command line invocation: /snap/minizinc/1222/bin/minizinc -c --solver gecode Brown-05.mzn + +array [1..5] of float: X_INTRODUCED_5_ = [1.0,1.0,1.0,2.0,1.0]; +array [1..5] of float: X_INTRODUCED_7_ = [1.0,1.0,2.0,1.0,1.0]; +array [1..5] of float: X_INTRODUCED_9_ = [1.0,2.0,1.0,1.0,1.0]; +array [1..5] of float: X_INTRODUCED_11_ = [2.0,1.0,1.0,1.0,1.0]; +var -10.0..10.0: X_INTRODUCED_0_; +var -10.0..10.0: X_INTRODUCED_1_; +var -10.0..10.0: X_INTRODUCED_2_; +var -10.0..10.0: X_INTRODUCED_3_; +var -10.0..10.0: X_INTRODUCED_4_; +var -100.0..100.0: X_INTRODUCED_13_ ::var_is_introduced :: is_defined_var; +var -1000.0..1000.0: X_INTRODUCED_14_ ::var_is_introduced :: is_defined_var; +var -10000.0..10000.0: X_INTRODUCED_15_ ::var_is_introduced :: is_defined_var; +var 1.0..1.0: X_INTRODUCED_16_ ::var_is_introduced :: is_defined_var; +array [1..5] of var float: x:: output_array([1..5]) = [X_INTRODUCED_0_,X_INTRODUCED_1_,X_INTRODUCED_2_,X_INTRODUCED_3_,X_INTRODUCED_4_]; +constraint float_lin_eq(X_INTRODUCED_5_,[X_INTRODUCED_3_,X_INTRODUCED_2_,X_INTRODUCED_1_,X_INTRODUCED_0_,X_INTRODUCED_4_],6.0); +constraint float_lin_eq(X_INTRODUCED_7_,[X_INTRODUCED_3_,X_INTRODUCED_2_,X_INTRODUCED_1_,X_INTRODUCED_0_,X_INTRODUCED_4_],6.0); +constraint float_lin_eq(X_INTRODUCED_9_,[X_INTRODUCED_3_,X_INTRODUCED_2_,X_INTRODUCED_1_,X_INTRODUCED_0_,X_INTRODUCED_4_],6.0); +constraint float_lin_eq(X_INTRODUCED_11_,[X_INTRODUCED_3_,X_INTRODUCED_2_,X_INTRODUCED_1_,X_INTRODUCED_0_,X_INTRODUCED_4_],6.0); +constraint float_times(X_INTRODUCED_0_,X_INTRODUCED_1_,X_INTRODUCED_13_):: defines_var(X_INTRODUCED_13_); +constraint float_times(X_INTRODUCED_13_,X_INTRODUCED_2_,X_INTRODUCED_14_):: defines_var(X_INTRODUCED_14_); +constraint float_times(X_INTRODUCED_14_,X_INTRODUCED_3_,X_INTRODUCED_15_):: defines_var(X_INTRODUCED_15_); +constraint float_times(X_INTRODUCED_15_,X_INTRODUCED_4_,X_INTRODUCED_16_):: defines_var(X_INTRODUCED_16_); +solve satisfy; diff --git a/include/common_solving.hpp b/include/common_solving.hpp index f2d0a8a2..cb3425b4 100644 --- a/include/common_solving.hpp +++ b/include/common_solving.hpp @@ -36,6 +36,10 @@ #include "lala/XCSP3_parser.hpp" #endif +#ifdef WITH_NNV + #include "nnv.hpp" +#endif + using namespace lala; #ifndef TURBO_ITV_BITS @@ -44,14 +48,17 @@ using namespace lala; #if (TURBO_ITV_BITS == 64) using bound_value_type = long long int; + using fbound_value_type = double; #elif (TURBO_ITV_BITS == 16) using bound_value_type = short int; #elif (TURBO_ITV_BITS == 32) using bound_value_type = int; + using fbound_value_type = float; #else #error "Invalid value for TURBO_ITV_BITS: must be 16, 32 or 64." #endif using Itv = Interval>; +using FItv = Interval>>; static std::atomic got_signal; static void (*prev_sigint)(int); @@ -188,7 +195,7 @@ struct AbstractDomains { , solver_output(basic_allocator) , config(other.config, basic_allocator) , stats(other.stats) - , env(basic_allocator) + , env(other.env, basic_allocator) , minimize_obj_var(other.minimize_obj_var) , store(store_allocator) , iprop(prop_allocator) @@ -266,6 +273,9 @@ struct AbstractDomains { abstract_ptr split; abstract_ptr search_tree; abstract_ptr best; + battery::vector inner_boxes; + battery::vector input_neurons; + battery::vector hidden_neurons; abstract_ptr bab; // The environment of variables, storing the mapping between variable's name and their representation in the abstract domains. @@ -292,6 +302,7 @@ struct AbstractDomains { search_tree = battery::allocate_shared(basic_allocator, env.extends_abstract_dom(), iprop, split, basic_allocator); // Note that `best` must have the same abstract type then store (otherwise projection of the variables will fail). best = battery::allocate_shared(basic_allocator, store->aty(), num_vars, basic_allocator); + inner_boxes = battery::vector(basic_allocator); bab = battery::allocate_shared(basic_allocator, env.extends_abstract_dom(), search_tree, best); if(config.verbose_solving) { printf("%% Abstract domain allocated.\n"); @@ -410,8 +421,20 @@ struct AbstractDomains { } #ifdef WITH_XCSP3PARSER else if(config.input_format() == InputFormat::XCSP3) { + solver_output.set_type(OutputType::XCSP); f = parse_xcsp3(config.problem_path.data(), solver_output); } +#endif +#ifdef WITH_NNV + else if (config.input_format() == InputFormat::VNNLIB || + config.input_format() == InputFormat::ONNX) { + solver_output.set_type(OutputType::NNV); + f = parse_nnv(config.onnx_path.data(), config.vnnlib_path.data(), input_neurons, hidden_neurons, solver_output, true); + } + else if (config.input_format() == InputFormat::SMT2) { + solver_output.set_type(OutputType::SMT2); + f = parse_smt2(config.problem_path.data(), solver_output, false); + } #endif if(!f) { std::cerr << "Could not parse input file." << std::endl; @@ -520,7 +543,11 @@ struct AbstractDomains { } void preprocess_tcn(F& f) { - f = ternarize(f, VarEnv(), {0,1,2}); +#ifdef WITH_NNV + f = ternarize(f, VarEnv(), false); +#else + f = ternarize(f, VarEnv(), true, {0,1,2}); +#endif battery::vector extra; f = normalize(f, extra); size_t num_vars = num_quantified_vars(f); @@ -545,15 +572,26 @@ struct AbstractDomains { while(!iprop->is_bot() && has_changed) { has_changed = false; preprocessing_stats.prepare_next_iteration(); +#ifdef WITH_NNV + fp_engine.fixpoint(iprop->num_deductions(), + [&](size_t i) { return iprop->fdeduce(i, config.epsilon); }, + [&](){ return iprop->is_bot(); }, + has_changed); +#else fp_engine.fixpoint(iprop->num_deductions(), [&](size_t i) { return iprop->deduce(i); }, [&](){ return iprop->is_bot(); }, has_changed); +#endif if(has_changed) { simplifier->meet_equivalence_classes(); } has_changed |= simplifier->algebraic_simplify(tnf, preprocessing_stats); +#ifdef WITH_NNV + simplifier->feliminate_entailed_constraints(*iprop, tnf, preprocessing_stats, config.epsilon); +#else simplifier->eliminate_entailed_constraints(*iprop, tnf, preprocessing_stats); +#endif // if(num_vars < 1000000) { // otherwise ICSE is too slow, needs to be improved. has_changed |= simplifier->i_cse(tnf, preprocessing_stats); // } @@ -585,6 +623,9 @@ struct AbstractDomains { const char* name_of_abstract_domain() const { #define STR_(x) #x #define STR(x) STR_(x) + #ifdef WITH_NNV + return "pir_itv_f"; + #endif #ifdef TURBO_IPC_ABSTRACT_DOMAIN return "ipc_itv" STR(TURBO_ITV_BITS) "_z"; #else @@ -640,8 +681,32 @@ struct AbstractDomains { template CUDA bool interpret_default_strategy() { typename F::Sequence seq; - seq.push_back(F::make_nary("first_fail", {})); - seq.push_back(F::make_nary("indomain_min", {})); +#ifdef WITH_NNV + if(config.var_order == "default" && config.value_order == "default") { + seq.push_back(F::make_nary("anti_first_fail", {})); + seq.push_back(F::make_nary("indomain_split", {})); + } + else { + seq.push_back(F::make_nary(config.var_order.data(), {})); + seq.push_back(F::make_nary(config.value_order.data(), {})); + } + // Add variables to split here as additional arguments: + for (int i = 0; i < input_neurons.size(); ++i ){ + seq.push_back(F::make_lvar(UNTYPED, LVar(input_neurons[i]))); + } + // for (int i = 0; i < hidden_neurons.size(); ++i){ + // seq.push_back(F::make_lvar(UNTYPED, LVar(hidden_neurons[i]))); + // } +#else + if(config.var_order == "default" && config.value_order == "default") { + seq.push_back(F::make_nary("first_fail", {})); + seq.push_back(F::make_nary("indomain_min", {})); + } + else { + seq.push_back(F::make_nary(config.var_order.data(), {})); + seq.push_back(F::make_nary(config.value_order.data(), {})); + } +#endif F search_strat = F::make_nary("search", std::move(seq)); if(!interpret_and_diagnose_and_tell(search_strat, env, *bab)) { return false; @@ -684,8 +749,7 @@ struct AbstractDomains { if(config.verbose_solving > 1) { printf("%% (Histogram of the number of times a function or predicate symbol occurs in the formula. Top-level conjunctions and unary constraints are discarded.)\n"); } - stats.print_dict_stat("fcn_histogram_reified_predicates", stats_fcn.reified_predicates, - [](const auto& key) { return "'" + std::string(string_of_sig_txt(key)) + "'"; }, + stats.print_dict_stat("fcn_histogram_reified_predicates", stats_fcn.reified_predicates, [](const auto& key) { return "'" + std::string(string_of_sig_txt(key)) + "'"; }, [](const auto& value) { return std::to_string(value); }); if(config.verbose_solving > 1) { printf("%% (Count all the predicate symbols occuring in the formula in a reified context, e.g., below a NOT, OR, or inside an arithmetic expression).\n"); @@ -782,7 +846,7 @@ struct AbstractDomains { stats_tcn.histogram_unassigned_vars_degree[stats_tcn.vars_occurrences[i]]++; stats_tcn.num_unassigned_var_occurrences += stats_tcn.vars_occurrences[i]; } - else if(width.value() == 1) { + else if(width.value() == 1 || width.value() <= config.epsilon) { stats_tcn.num_assigned_vars++; stats_tcn.histogram_assigned_vars_degree[stats_tcn.vars_occurrences[i]]++; stats_tcn.num_assigned_var_occurrences += stats_tcn.vars_occurrences[i]; @@ -857,6 +921,7 @@ struct AbstractDomains { /** Return `true` if the search state must be pruned. */ CUDA bool update_solution_stats() { stats.solutions++; + inner_boxes.push_back(*best); if(bab->is_satisfaction() && config.stop_after_n_solutions != 0 && stats.solutions >= config.stop_after_n_solutions) { @@ -877,10 +942,21 @@ struct AbstractDomains { stats.fails += 1; } + CUDA void on_unknown_node() { + stats.unknowns += 1; + } + CUDA void print_final_solution() { if(!is_printing_intermediate_sol() && stats.solutions > 0) { print_solution(); } +#ifdef WITH_NNV + else if (bab->is_satisfaction() && inner_boxes.size() > 0) { + for(int i = 0; i < inner_boxes.size(); ++i) { + print_solution(inner_boxes[i]); + } + } +#endif stats.print_mzn_final_separator(); } diff --git a/include/config.hpp b/include/config.hpp index 6474e220..2ddb2a1f 100644 --- a/include/config.hpp +++ b/include/config.hpp @@ -16,7 +16,9 @@ enum class Arch { CPU, GPU, BAREBONES, - HYBRID + FBAREBONES, + HYBRID, + JET }; enum class FixpointKind { @@ -26,7 +28,10 @@ enum class FixpointKind { enum class InputFormat { XCSP3, - FLATZINC + FLATZINC, + VNNLIB, + ONNX, + SMT2 }; template @@ -52,9 +57,14 @@ struct Configuration { FixpointKind fixpoint; size_t wac1_threshold; size_t seed; + double epsilon; + battery::string var_order; + battery::string value_order; battery::string eps_var_order; battery::string eps_value_order; battery::string problem_path; + battery::string vnnlib_path; + battery::string onnx_path; battery::string version; battery::string hardware; @@ -97,9 +107,14 @@ struct Configuration { ), wac1_threshold(0), seed(0), + epsilon(1e-6), + value_order("default", alloc), + var_order("default", alloc), eps_value_order("default", alloc), eps_var_order("default", alloc), problem_path(alloc), + vnnlib_path(alloc), + onnx_path(alloc), version(alloc), hardware(alloc) {} @@ -129,9 +144,14 @@ struct Configuration { fixpoint(other.fixpoint), wac1_threshold(other.wac1_threshold), seed(other.seed), + epsilon(other.epsilon), + var_order(other.var_order, alloc), + value_order(other.value_order, alloc), eps_var_order(other.eps_var_order, alloc), eps_value_order(other.eps_value_order, alloc), problem_path(other.problem_path, alloc), + vnnlib_path(other.vnnlib_path, alloc), + onnx_path(other.onnx_path, alloc), version(other.version, alloc), hardware(other.hardware, alloc) {} @@ -158,9 +178,14 @@ struct Configuration { fixpoint = other.fixpoint; wac1_threshold = other.wac1_threshold; seed = other.seed; + epsilon = other.epsilon; + var_order = other.var_order; + value_order = other.value_order; eps_var_order = other.eps_var_order; eps_value_order = other.eps_value_order; problem_path = other.problem_path; + vnnlib_path = other.vnnlib_path; + onnx_path = other.onnx_path; version = other.version; hardware = other.hardware; } @@ -194,6 +219,9 @@ struct Configuration { printf("-wac1_threshold %" PRIu64 " ", wac1_threshold); } printf("-seed %" PRIu64 " ", seed); + printf("-epsilon %.10f", epsilon); + printf("-var_order %s ", var_order.data()); + printf("-value_order %s ", value_order.data()); printf("-eps_var_order %s ", eps_var_order.data()); printf("-eps_value_order %s ", eps_value_order.data()); if(version.size() != 0) { @@ -226,6 +254,8 @@ struct Configuration { return "gpu"; case Arch::BAREBONES: return "barebones"; + case Arch::FBAREBONES: + return "fbarebones"; case Arch::HYBRID: return "hybrid"; default: @@ -247,6 +277,9 @@ struct Configuration { printf("%%%%%%mzn-stat: wac1_threshold=%" PRIu64 "\n", wac1_threshold); } printf("%%%%%%mzn-stat: seed=%" PRIu64 "\n", seed); + printf("%%%%%%nnv-stat: epsilon=\"%.10f\"\n", epsilon); + printf("%%%%%%mzn-stat: var_order=\"%s\"\n", var_order.data()); + printf("%%%%%%mzn-stat: value_order=\"%s\"\n", value_order.data()); printf("%%%%%%mzn-stat: eps_var_order=\"%s\"\n", eps_var_order.data()); printf("%%%%%%mzn-stat: eps_value_order=\"%s\"\n", eps_value_order.data()); printf("%%%%%%mzn-stat: free_search=\"%s\"\n", free_search ? "yes" : "no"); @@ -272,6 +305,15 @@ struct Configuration { else if(problem_path.ends_with(".xml")) { return InputFormat::XCSP3; } + else if (problem_path.ends_with(".vnnlib")){ + return InputFormat::VNNLIB; + } + else if (problem_path.ends_with(".onnx")){ + return InputFormat::ONNX; + } + else if (problem_path.ends_with(".smt2")){ + return InputFormat::SMT2; + } else { printf("ERROR: Unknown input format for the file %s [supported extension: .xml and .fzn].\n", problem_path.data()); exit(EXIT_FAILURE); diff --git a/include/cpu_solving.hpp b/include/cpu_solving.hpp index 4344035e..6f66a122 100644 --- a/include/cpu_solving.hpp +++ b/include/cpu_solving.hpp @@ -8,7 +8,11 @@ void cpu_solve(const Configuration& config) { auto start = std::chrono::steady_clock::now(); +#ifdef WITH_NNV + CP cp(config); +#else CP cp(config); +#endif cp.config.or_nodes = 1; cp.preprocess(); if(cp.iprop->is_bot()) { @@ -23,14 +27,42 @@ void cpu_solve(const Configuration& config) { while(!must_quit(cp) && check_timeout(cp, start) && has_changed) { has_changed = false; auto start2 = cp.stats.start_timer_host(); - cp.stats.fixpoint_iterations += fp_engine.fixpoint([&](int i) { return cp.iprop->deduce(i); }); +#ifdef WITH_NNV + cp.stats.fixpoint_iterations += fp_engine.fixpoint( + [&](int i) { return cp.iprop->fdeduce(i, cp.config.epsilon); } + ); +#else + cp.stats.fixpoint_iterations += fp_engine.fixpoint( + [&](int i) { return cp.iprop->deduce(i); }, + [&]() { return cp.iprop->is_bot(); }); +#endif start2 = cp.stats.stop_timer(Timer::FIXPOINT, start2); - bool must_prune = cp.on_node(); + bool must_prune = cp.on_node(); if(cp.iprop->is_bot()) { +#ifdef WITH_NNV + if (cp.search_tree->is_unknown(cp.env, cp.config.epsilon)) { + cp.on_unknown_node(); + } + else { + cp.on_failed_node(); + } + fp_engine.reset(); +#else cp.on_failed_node(); fp_engine.reset(); +#endif } else { +#ifdef WITH_NNV + cp.stats.stop_timer(Timer::SELECT_FP_FUNCTIONS, start2); + // if(cp.search_tree->template is_extractable(AtomicExtraction(), config.epsilon)) { + if(cp.search_tree->is_solution(cp.env)) { + has_changed |= cp.bab->deduce(); + must_prune |= cp.on_solution_node(); + fp_engine.reset(); + break; + } +#else fp_engine.select([&](int i) { return !cp.iprop->ask(i); }); cp.stats.stop_timer(Timer::SELECT_FP_FUNCTIONS, start2); if(fp_engine.num_active() == 0 && cp.search_tree->template is_extractable()) { @@ -38,13 +70,23 @@ void cpu_solve(const Configuration& config) { must_prune |= cp.on_solution_node(); fp_engine.reset(); } +#endif } +#ifdef WITH_NNV + has_changed |= cp.search_tree->fdeduce(cp.env, cp.config.epsilon); // add branching strategies +#else has_changed |= cp.search_tree->deduce(); +#endif cp.stats.stop_timer(Timer::SEARCH, start2); if(must_prune) { break; } } cp.print_final_solution(); cp.print_mzn_statistics(); + + if (cp.stats.solutions > 0) printf("sat\n"); + else if (cp.stats.unknowns > 0 && check_timeout(cp, start)) printf("unknown\n"); + else if (!check_timeout(cp, start)) printf("timeout\n"); + else printf("unsat\n"); } #endif diff --git a/include/fbarebones_dive_and_solve.hpp b/include/fbarebones_dive_and_solve.hpp new file mode 100644 index 00000000..20d3e5d9 --- /dev/null +++ b/include/fbarebones_dive_and_solve.hpp @@ -0,0 +1,1421 @@ +// Copyright 2026 Yi-Nung Tsao + +#ifndef TURBO_FBAREBONES_DIVE_AND_SOLVE_HPP +#define TURBO_FBAREBONES_DIVE_AND_SOLVE_HPP + +#include "battery/allocator.hpp" +#include "common_solving.hpp" +#include "memory_gpu.hpp" +#include "lala/light_branch.hpp" +#include +#include +#include +#include +#include + +/** This is required in order to guess the usage of global memory, and increase the CUDA default limit. */ +#define MAX_SEARCH_DEPTH 10000 + +namespace bt = ::battery; + +/** + * The full GPU version (`gpu_dive_and_solve`) is not compiling on modern GPU hardware (SM >= 9) due to the kernel being too large. + * We circuvanted this issue by creating an hybrid version where only propagation is executed on the GPU (`hybrid_dive_and_solve`). + * This has the disadvantage of memory transfers between CPU and GPU and synchronization overheads. + * We propose a new "fbarebones" version which contains less abstractions than the GPU and hybrid versions, but have the same functionalities. + * In particular, we directly implement the branch-and-bound algorithm here and avoid using `lala::BAB` and `lala::SearchTree` which are nice from a software engineering perspective but bring significant overhead. + * This version is intended to reach the best possible performance. + * + * Terminology: + * * unified data: data available to both the CPU and GPU. + * * block data: data used within a single block. + * * grid data: data shared among all blocks in the grid. + */ + +#ifdef __CUDACC__ + +#include +#include +#include "lala/onnxruntime-linux-x64-gpu-1.19.2/include/onnxruntime_cxx_api.h" +#include + +#endif + +namespace fbarebones { + +#ifdef __CUDACC__ +#ifndef TURBO_IPC_ABSTRACT_DOMAIN + +/** `ConcurrentAllocator` allocates memory available both on CPU and GPU. For non-Linux systems such as Windows pinned memory must be used (see PR #19). */ +#ifdef NO_CONCURRENT_MANAGED_MEMORY + using ConcurrentAllocator = bt::pinned_allocator; +#else + using ConcurrentAllocator = bt::managed_allocator; +#endif + +using ::FItv; +using GridCP = AbstractDomains, + bt::statistics_allocator>, + bt::statistics_allocator>>; + +/** Data shared between CPU and GPU. */ +struct UnifiedData { + /** The root node of the problem, useful to backtrack when solving a new subproblem. + * Also contains the shared information such as statistics and solver configuration. + */ + GridCP root; + + /** Stop signal from the CPU because of a timeout or CTRL-C. */ + cuda::std::atomic_flag stop; + + /** The memory configuration of each block. */ + MemoryConfig mem_config; + + UnifiedData(const CP& cp, MemoryConfig mem_config) + : root(GridCP::tag_gpu_block_copy{}, false, cp) + , stop(false) + , mem_config(mem_config) + { + size_t num_subproblems = 1; + num_subproblems <<= root.config.subproblems_power; + root.stats.eps_num_subproblems = num_subproblems; + } +}; + +struct GridData; +using FStore = VStore; +using FProp = PIR; +using bound_type = typename FItv::LB::value_type; +using UB = FUB; +using strategies_type = bt::vector, bt::global_allocator>; + +/** Data private to a single block. */ +struct BlockData { + /** The store of variables at the root of the current subproblem. */ + abstract_ptr> root_store; + + // inner box + abstract_ptr> inner_box; + + /** The current store of variables. + * We use a `pool_allocator`, this allows to easily switch between global memory and shared memory, if the store of variables can fit inside. + * */ + abstract_ptr store; + + /** The propagators implemented as an array of bytecodes. + * Similarly, the propagators can be stored in the global or shared memory. + * If the propagators do not fit in shared memory, the array of propagators is shared among all blocks. + * It is possible because the propagators are state-less, we avoid duplicating them in each block. + * */ + abstract_ptr iprop; + + /** The statistics of the current block. */ + Statistics stats; + + /** The path from `UnifiedData::root` to the current subproblem `root_store`. */ + size_t subproblem_idx; + + /** The current strategy being used to split the store. + * It is an index into `GridData::strategies`. + */ + int current_strategy; + + /** The next unassigned variable in the current strategy. + * It is an index into `GridData::strategies.vars`. + */ + int next_unassigned_var; + + /** On backtracking, the value to restore `current_strategy`. */ + int snapshot_root_strategy; + + /** On backtracking, the value to restore `next_unassigned_var`. */ + int snapshot_next_unassigned_var; + + /** The decision taken when exploring the tree. */ + bt::vector, bt::global_allocator> decisions; + + /** Current depth of the search tree. */ + int depth; + + /** A timer used for computing time statistics. */ + cuda::std::chrono::system_clock::time_point timer; + + /** A timer used for computing diving VS search time statistics. */ + cuda::std::chrono::system_clock::time_point dive_timer; + + /** The time at which the kernel was started, useful to compute the time of the best bound. */ + cuda::std::chrono::system_clock::time_point start_time; + + /** The gradients from the neural network. */ + float* h_gradients; + float* h_mid_gradients; + // float* h_lb_gradients; + // float* h_ub_gradients; + int num_h_gradients; + + /* For underapproximation search strategy. */ + bool is_uass; + + __device__ BlockData() + : subproblem_idx(0) + , current_strategy(0) + , next_unassigned_var(0) + , decisions(5000) + , depth(0) + , h_gradients(nullptr) + // , h_mid_gradients(nullptr) + // , h_lb_gradients(nullptr) + // , h_ub_gradients(nullptr) + // , num_h_gradients(0) + , is_uass(false) + {} + + __device__ void allocate(const UnifiedData& unified_data, const GridData& grid_data, unsigned char* shared_mem) { + if(threadIdx.x == 0) { + subproblem_idx = blockIdx.x; + const MemoryConfig& mem_config = unified_data.mem_config; + const auto& u_store = *(unified_data.root.store); + const auto& u_iprop = *(unified_data.root.iprop); + bt::pool_allocator shared_mem_pool(mem_config.make_shared_pool(shared_mem)); + bt::pool_allocator store_allocator(mem_config.make_store_pool(shared_mem_pool)); + bt::pool_allocator prop_allocator(mem_config.make_prop_pool(shared_mem_pool)); + root_store = bt::make_shared, bt::global_allocator>(u_store); + inner_box = bt::make_shared, bt::global_allocator>(u_store); + store = bt::allocate_shared(store_allocator, u_store, store_allocator); + iprop = bt::allocate_shared(prop_allocator, u_iprop, store, prop_allocator); + num_h_gradients = u_store.vars(); // only take input neurons. + size_t gradient_bytes = sizeof(float) * static_cast(num_h_gradients) * 4; + void* gradient_mem = bt::global_allocator{}.allocate(gradient_bytes); + h_gradients = static_cast(gradient_mem); + // h_mid_gradients = h_gradients + num_h_gradients; + // h_lb_gradients = h_mid_gradients + num_h_gradients; + // h_ub_gradients = h_lb_gradients + num_h_gradients; + is_uass = false; + } + } + + /** We must deallocate store and iprop inside the kernel because they might be initialized in shared memory. */ + __device__ void deallocate_shared_data() { + if(threadIdx.x == 0) { + // NOTE: .reset() does not work because it does not reset the allocator, which is itself allocated in global memory. + store = abstract_ptr(); + iprop = abstract_ptr(); + } + } + + /** Add a new decision on the `decisions` stack and increase depth. + * \param has_changed: A Boolean in shared memory. + * \param strategies: A sequence of strategies. + * \precondition: We must not be on a leaf node. + */ + __device__ INLINE void split(bool& has_changed, const strategies_type& strategies, const float epsilon) { + using LB2 = typename FItv::LB::local_type; + using UB2 = typename FItv::UB::local_type; + __shared__ local::ZUB idx; + decisions[depth].var = AVar{}; + int currentDepth = depth; + for(int i = current_strategy; i < strategies.size(); ++i) { + switch(strategies[i].var_order) { + case VariableOrder::RANDOM: + case VariableOrder::INPUT_ORDER: { + input_order_split(has_changed, idx, strategies[i], epsilon); + break; + } + case VariableOrder::FIRST_FAIL: { + lattice_smallest_split(has_changed, idx, strategies[i], epsilon, + [&](const FItv& u, int g_idx) { return UB2(u.width().ub().value()); }); + break; + } + case VariableOrder::ANTI_FIRST_FAIL: { + lattice_smallest_split(has_changed, idx, strategies[i], epsilon, + [&](const FItv& u, int g_idx) { return LB2(u.width().ub().value()); }); + break; + } + case VariableOrder::GRA_ANTI_FIRST_FAIL: { + lattice_smallest_split(has_changed, idx, strategies[i], epsilon, + [&](const FItv& u, int g_idx) { return LB2(battery::mul_up(u.width().ub().value(), h_gradients[g_idx])); }); + break; + } + case VariableOrder::LARGEST: { + lattice_smallest_split(has_changed, idx, strategies[i], epsilon, + [&](const FItv& u, int g_idx) { return LB2(u.ub().value()); }); + break; + } + case VariableOrder::SMALLEST: { + lattice_smallest_split(has_changed, idx, strategies[i], epsilon, + [&](const FItv& u, int g_idx) { return UB2(u.lb().value()); }); + break; + } + default: assert(false); + } + __syncthreads(); + // If we could find a variable with the current strategy, we return. + if(!decisions[currentDepth].var.is_untyped()) { + return; + } + if(threadIdx.x == 0) { + current_strategy++; + next_unassigned_var = 0; + } + } + // `input_order_split` and `lattice_smallest_split` have a `__syncthreads()` before reading next_unassigned_var. + } + + /** Select the next unassigned variable with a finite interval in the array `strategy.vars()` or `store` if the previous one is empty. + * We ignore infinite variables as splitting on them do not guarantee termination. + * \param has_changed is a Boolean in shared memory. + * \param idx is a decreasing integer in shared memory. + */ + __device__ INLINE void input_order_split(bool& has_changed, local::ZUB& idx, const StrategyType& strategy, const float epsilon) + { + bool split_in_store = strategy.vars.empty(); + int n = split_in_store ? store->vars() : strategy.vars.size(); + if(threadIdx.x == 0) { + has_changed = true; + idx = n; + } + __syncthreads(); + while(has_changed) { + __syncthreads(); + // int n = idx.value(); + if(threadIdx.x == 0) { + has_changed = false; + } + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x) { + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() > epsilon && !dom.lb().is_top() && !dom.ub().is_top()) { + if(idx.meet(local::ZUB(split_in_store ? i : strategy.vars[i].vid()))) { + has_changed = true; + } + } + } + __syncthreads(); + } + if(threadIdx.x == 0) { + next_unassigned_var = idx.value(); + if(next_unassigned_var != n) { + push_decision(strategy.val_order, split_in_store ? AVar{store->aty(), next_unassigned_var} : strategy.vars[next_unassigned_var], epsilon); + } + } + } + + /** Given an array of variable, select the variable `x` with the smallest value `f(store[x])` where "smallest" is defined according to the lattice order of the return type of `f`. + * \param has_changed is a Boolean in shared memory. + * \param idx is a decreasing integer in shared memory. + * */ + template + __device__ INLINE void lattice_smallest_split(bool& has_changed, local::ZUB& idx, + const StrategyType& strategy, const float epsilon, F f) + { + using T = decltype(f(FItv{},0)); + __shared__ T value; + bool split_in_store = strategy.vars.empty(); + int n = split_in_store ? store->vars() : strategy.vars.size(); + __syncthreads(); + if(threadIdx.x == 0) { + has_changed = true; + value = T::top(); + idx = n; + } + __syncthreads(); + /** This fixpoint loop seeks for the smallest `x` according to `f(x)` and the next unassigned variable. */ + while(has_changed) { + __syncthreads(); + if(threadIdx.x == 0) { + has_changed = false; + } + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x) { + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() > epsilon && !dom.lb().is_top() && !dom.ub().is_top()) { + //if (dom.lb().value() != dom.ub().value() && !dom.lb().is_top() && !dom.ub().is_top()) { + if(value.meet(f(dom, strategy.vars[i].vid()))) { + has_changed = true; + } + if(idx.meet(local::ZUB(strategy.vars[i].vid()))) { + has_changed = true; + } + } + } + __syncthreads(); + } + /** If we found a value, we traverse again the variables' array to find its index. */ + if(!value.is_top()) { + __syncthreads(); + if(threadIdx.x == 0) { + next_unassigned_var = idx.value(); + idx = n; + has_changed = true; + is_uass = false; + } + __syncthreads(); + // This fixpoint loop is not strictly necessary. + // We keep it for determinism: the variable with the smallest index is selected first. + while(has_changed) { + // int n = idx.value(); + __syncthreads(); + has_changed = false; + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x) { + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() > epsilon && !dom.lb().is_top() && !dom.ub().is_top() && f(dom,strategy.vars[i].vid()) == value) { + if(idx.meet(local::ZUB(strategy.vars[i].vid()))) { + has_changed = true; + } + } + } + __syncthreads(); + } + assert(idx.value() < n); + if(threadIdx.x == 0) { + if(split_in_store) { + push_decision(strategy.val_order, AVar{store->aty(), idx.value()}, epsilon); + } + else { + push_decision(strategy.val_order, strategy.vars[idx.value()], epsilon); + } + } + return; + } + + /* + Original search stategy. When using SPLIT && all widths <= epsilon, we check the solution with midpoints. + */ + if(strategy.val_order == ValueOrder::SPLIT){ + if(threadIdx.x == 0){ + has_changed = false; + } + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x){ + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() <= epsilon && dom.lb().value() != dom.ub().value() && !dom.lb().is_top() && !dom.ub().is_top()){ + has_changed = true; + } + } + __syncthreads(); + if(has_changed) { + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x){ + AVar var = split_in_store ? AVar{store->aty(), i} : strategy.vars[i]; + const auto& dom = (*store)[var.vid()]; + auto mid = battery::midpoint(dom.lb().value(), dom.ub().value()); + store->embed(var, FItv(mid, mid)); + } + __syncthreads(); + if(threadIdx.x == 0){ + is_uass = true; + push_decision(strategy.val_order, strategy.vars[0], epsilon); + } + } + return; + } + + /* + Underapproximation search strategy. + */ + if(threadIdx.x == 0){ + has_changed = true; + } + __syncthreads(); + while(has_changed) { + __syncthreads(); + if(threadIdx.x == 0){ + has_changed = false; + } + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x) { + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() <= epsilon && dom.lb().value() != dom.ub().value() && !dom.lb().is_top() && !dom.ub().is_top()){ + if(value.meet(f(dom, strategy.vars[i].vid()))){ + has_changed = true; + } + if(idx.meet(local::ZUB(strategy.vars[i].vid()))){ + has_changed = true; + } + } + } + __syncthreads(); + } + if(!value.is_top()) { + __syncthreads(); + if(threadIdx.x == 0){ + next_unassigned_var = idx.value(); + idx = n; + has_changed = true; + is_uass = true; + } + __syncthreads(); + while(has_changed) { + __syncthreads(); + has_changed = false; + __syncthreads(); + for(int i = next_unassigned_var + threadIdx.x; i < n; i += blockDim.x) { + const auto& dom = (*store)[split_in_store ? i : strategy.vars[i].vid()]; + if(dom.width().ub().value() <= epsilon && dom.lb().value() != dom.ub().value() && !dom.lb().is_top() && !dom.ub().is_top() && f(dom, strategy.vars[i].vid()) == value){ + if(idx.meet(local::ZUB(strategy.vars[i].vid()))){ + has_changed = true; + } + } + } + __syncthreads(); + } + assert(idx.value() < n); + if(threadIdx.x == 0){ + if(split_in_store){ + push_decision(strategy.val_order, AVar{store->aty(), idx.value()}, epsilon); + } + else{ + push_decision(strategy.val_order, strategy.vars[idx.value()], epsilon); + } + } + } + } + + /** + * + * TODO: this function should be implemented later. It is for improving the performance. + * */ + __device__ INLINE void floating_split(bool& has_changed, local::ZUB& idx, + const StrategyType& strategy, const float epsilon) + { + bool split_in_store = strategy.vars.empty(); + int n = split_in_store ? store->vars() : strategy.vars.size(); + if(threadIdx.x == 0) { + has_changed = true; + idx = n; + } + __syncthreads(); + bt::vector E(n, false); + for(int i = threadIdx.x; i < iprop->num_deductions(); i += blockDim.x) { + if(!iprop->is_fsolution(i, epsilon)) { + has_changed = true; + E[iprop->load_deduce(i).x.vid()] = true; + E[iprop->load_deduce(i).y.vid()] = true; + E[iprop->load_deduce(i).z.vid()] = true; + } + } + __syncthreads(); + for(int i = threadIdx.x; i < n; i += blockDim.x) { + const int dom_id = split_in_store ? i : strategy.vars[i].vid(); + const auto& dom = (*store)[dom_id]; + if(dom.width().ub().value() > epsilon && E[dom_id] && !dom.lb().is_top() && !dom.ub().is_top()) { + if(idx.meet(local::ZUB(i))) { + has_changed = true; + break; + } + } + } + __syncthreads(); + if(threadIdx.x == 0) { + next_unassigned_var = idx.value(); + if(next_unassigned_var != n) { + push_decision(strategy.val_order, split_in_store ? AVar{store->aty(), next_unassigned_var} : strategy.vars[next_unassigned_var], epsilon); + } + } + } + + /** Push a new decision onto the decisions stack. + * \precondition The domain of the variable `var` must not be empty, be a singleton or contain infinite bounds. + * \precondition Must be executed by thread 0 only. + */ + __device__ INLINE void push_decision(ValueOrder val_order, AVar var, const float epsilon) { + assert(threadIdx.x == 0); + decisions[depth].var = var; + decisions[depth].current_idx = -1; + const auto& dom = store->project(decisions[depth].var); + // printf("split on %d \n", decisions[depth].var.vid()); + // assert(dom.width().ub().value() > epsilon); + assert(dom.width().ub().value() > 0.0); + // auto mid = battery::add_down(dom.lb().value(), battery::div_down(battery::sub_up(dom.ub().value(), dom.lb().value()), bound_type{2.0})); + // bound_type width = battery::sub_up(dom.ub().value(), dom.lb().value()); + // bound_type half = battery::div_up(width, bound_type{2.0}); + // bound_type mid = battery::add_up(dom.lb().value(), half); + bound_type mid = battery::midpoint(dom.lb().value(), dom.ub().value()); + // printf("split on %d, lb = %.20f, ub = %.20f, mid = %.20f \n", decisions[depth].var.vid(), dom.lb().value(), dom.ub().value(), mid); + + switch(val_order) { + case ValueOrder::SPLIT: { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + break; + } + case ValueOrder::REVERSE_SPLIT: { + decisions[depth].children[0] = FItv(mid, dom.ub()); + decisions[depth].children[1] = FItv(dom.lb(), mid); + break; + } + case ValueOrder::LB_SPLIT: { + if(dom.lb().value() != dom.ub().value() && dom.width().ub().value() <= epsilon) { + is_uass = true; + decisions[depth].children[0] = FItv(dom.lb().value(), dom.lb().value()); + } + else { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + } + break; + } + case ValueOrder::UB_SPLIT: { + if(dom.lb().value() != dom.ub().value() && dom.width().ub().value() <= epsilon) { + is_uass = true; + decisions[depth].children[0] = FItv(dom.ub().value(), dom.ub().value()); + } + else { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + } + break; + } + case ValueOrder::MID_SPLIT: { + if(dom.lb().value() != dom.ub().value() && dom.width().ub().value() <= epsilon) { + is_uass = true; + decisions[depth].children[0] = FItv(mid, mid); + } + else { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + } + break; + } + case ValueOrder::MID_LB_SPLIT: { + if(dom.lb().value() != dom.ub().value() && dom.width().ub().value() <= epsilon) { + is_uass = true; + decisions[depth].children[0] = FItv(mid, mid); + decisions[depth].children[1] = FItv(dom.lb().value(), dom.lb().value()); + } + else { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + } + break; + } + case ValueOrder::MID_UB_SPLIT: { + if(dom.lb().value() != dom.ub().value() && dom.width().ub().value() <= epsilon) { + is_uass = true; + decisions[depth].children[0] = FItv(mid, mid); + decisions[depth].children[1] = FItv(dom.ub().value(), dom.ub().value()); + } + else { + is_uass = false; + decisions[depth].children[0] = FItv(dom.lb(), mid); + decisions[depth].children[1] = FItv(mid, dom.ub()); + } + break; + } + case ValueOrder::MIX_SPLIT: { + // TODO: not complete yet. The current version doesn't support multiple chil nodes, only binary. + decisions[depth].children[0] = FItv(mid, mid); + decisions[depth].children[1] = FItv(dom.lb().value(), dom.lb().value()); + break; + decisions[depth].children[2] = FItv(dom.ub().value(), dom.ub().value()); + decisions[depth].children[3] = FItv(battery::nextafter(dom.lb().value(), 1e38f), battery::nextafter(mid, -1e38f)); + decisions[depth].children[4] = FItv(battery::nextafter(mid, 1e38f), battery::nextafter(dom.ub().value(), -1e38f)); + } + // ValueOrder::MEDIAN is not possible with interval. + default: assert(false); + } + /** Ropes are a mechanism for fast backtracking. + * The rope of a left node is always the depth of the right node (also its depth), because after completing the exploration of the left subtree, we must visit the right subtree (rooted at the current depth). + * The rope of the right node is inherited from its parent, we set -1 if there is no next node to visit. + */ + decisions[depth].ropes[0] = depth + 1; + decisions[depth].ropes[1] = depth > 0 ? decisions[depth-1].ropes[decisions[depth-1].current_idx] : -1; + ++depth; + // printf("depth(%d), var = %d, children = [%lf, %lf] | [%lf, %lf], ropes = [%d, %d]\n", + // depth, decisions[depth-1].var.vid(), + // (bound_type)decisions[depth-1].children[0].lb().value(), (bound_type)decisions[depth-1].children[0].ub().value(), + // (bound_type)decisions[depth-1].children[1].lb().value(), (bound_type)decisions[depth-1].children[1].ub().value(), + // decisions[depth-1].ropes[0], decisions[depth-1].ropes[1]); + // Reallocate decisions if needed. + if(decisions.size() == depth) { + printf("resize to %d\n", (int)decisions.size() * 2); + decisions.resize(decisions.size() * 2); + } + } +}; + +/** Data shared among all blocks. */ +struct GridData { + /** The private data of each block. */ + bt::vector blocks; + + /** We generate the subproblems lazily. + * Suppose we generate `2^3` subproblems, we represent the first subproblem as `000`, the second as `001`, the third as `010`, and so on. + * A `0` means to turn left in the search tree, and a `1` means to turn right. + * Incrementing this integer will generate the path of the next subproblem. + */ + ZLB next_subproblem; + + /** Due to multithreading, we must protect `stdout` when printing. + * The model of computation in this work is lock-free, but it seems unavoidable for printing. + */ + cuda::binary_semaphore print_lock; + + /** A specific strategy is used for the subproblem decomposition during the diving phase. */ + bool has_eps_strategy; + + /** The search strategy is immutable and shared among the blocks. */ + strategies_type search_strategies; + + /** The objective variable to minimize. + * Maximization problem are transformed into minimization problems by negating the objective variable. + * Equal to -1 if the problem is a satisfaction problem. + */ + AVar obj_var; + + __device__ GridData(const GridCP& root) + : blocks(root.stats.num_blocks) + , next_subproblem(root.stats.num_blocks) + , print_lock(1) + , has_eps_strategy(root.config.eps_var_order != "default") + , search_strategies(root.split->strategies_()) + , obj_var(root.minimize_obj_var) + {} +}; + +MemoryConfig configure_gpu_fbarebones(CP&); +__global__ void initialize_global_data(UnifiedData*, bt::unique_ptr*); +__global__ void gpu_fbarebones_solve(UnifiedData*, GridData*, Ort::Session&); +template +__device__ INLINE void propagate(UnifiedData& unified_data, GridData& grid_data, BlockData& block_data, + FPEngine& fp_engine, bool& stop, bool& has_changed, bool& is_leaf_node, Ort::Session&); +// __device__ INLINE void back_propagation(BlockData& block_data, Ort::Session& session); +__global__ void reduce_blocks(UnifiedData*, GridData*); +__global__ void deallocate_global_data(bt::unique_ptr*); + +void fbarebones_dive_and_solve(const Configuration& config) { + if(config.print_intermediate_solutions) { + printf("%% WARNING: -arch fbarebones is incompatible with -i and -a (it cannot print intermediate solutions).\n"); + } + auto start = std::chrono::steady_clock::now(); + check_support_managed_memory(); + check_support_concurrent_managed_memory(); + /** We start with some preprocessing to reduce the number of variables and constraints. */ + CP cp(config); + cp.preprocess(); + if(cp.iprop->is_bot()) { + cp.print_final_solution(); + cp.print_mzn_statistics(); + return; + } + + MemoryConfig mem_config = configure_gpu_fbarebones(cp); + auto unified_data = bt::make_unique(cp, mem_config); + auto grid_data = bt::make_unique, ConcurrentAllocator>(); + initialize_global_data<<<1,1>>>(unified_data.get(), grid_data.get()); + Ort::Env env(ORT_LOGGING_LEVEL_WARNING, "nnv"); + Ort::SessionOptions session_options; + OrtCUDAProviderOptions cuda_options; + cuda_options.device_id = 0; + session_options.AppendExecutionProvider_CUDA(cuda_options); + Ort::Session session(env, config.onnx_path.data(), session_options); + CUDAEX(cudaDeviceSynchronize()); + /** We wait that either the solving is interrupted, or that all threads have finished. */ + /** Block the signal CTRL-C to notify the threads if we must exit. */ + block_signal_ctrlc(); + gpu_fbarebones_solve + <<(cp.stats.num_blocks), + CUDA_THREADS_PER_BLOCK, + mem_config.shared_bytes>>> + (unified_data.get(), grid_data->get(), session); + auto now = std::chrono::steady_clock::now(); + int64_t time_to_kernel_start = std::chrono::duration_cast(now - start).count(); + bool interrupted = wait_solving_ends(unified_data->stop, unified_data->root, start); + CUDAEX(cudaDeviceSynchronize()); + reduce_blocks<<<1,1>>>(unified_data.get(), grid_data->get()); + CUDAEX(cudaDeviceSynchronize()); + auto& uroot = unified_data->root; + if(uroot.stats.solutions > 0) { + // We add the time before the kernel starts to the time needed to find the best bound. + uroot.stats.timers.time_of(Timer::LATEST_BEST_OBJ_FOUND) += time_to_kernel_start; + if(uroot.stats.timers.time_of(Timer::FIRST_BLOCK_IDLE) != 0) { + uroot.stats.timers.time_of(Timer::FIRST_BLOCK_IDLE) += time_to_kernel_start; + } + cp.print_solution(*uroot.best); + } + uroot.stats.print_mzn_final_separator(); + if(uroot.config.print_statistics) { + uroot.config.print_mzn_statistics(); + uroot.stats.print_mzn_statistics(uroot.config.verbose_solving); + if(uroot.bab->is_optimization() && uroot.stats.solutions > 0) { + uroot.stats.print_mzn_objective(uroot.best->project(uroot.bab->objective_var()), uroot.bab->is_minimization()); + } + unified_data->root.stats.print_mzn_end_stats(); + } + if (uroot.stats.solutions > 0) printf("sat\n"); + else if (uroot.stats.unknowns > 0) printf("unknown\n"); + else if (interrupted) printf("timeout\n"); + else printf("unsat\n"); + deallocate_global_data<<<1,1>>>(grid_data.get()); + CUDAEX(cudaDeviceSynchronize()); +} + +/** We configure the GPU according to the user configuration: + * 1) Guess the "best" number of blocks per SM, if not provided. + * 2) Update the number of subproblems to at least "30 * B" where B is the number of blocks. + * 3) Configure the size of the shared memory. + * 4) Increase the global heap memory. + * 5) Increase the stack size if requested by the user. + */ +MemoryConfig configure_gpu_fbarebones(CP& cp) { + auto& config = cp.config; + + /** I. Number of blocks per SM. */ + cudaDeviceProp deviceProp; + cudaGetDeviceProperties(&deviceProp, 0); + int max_block_per_sm; + cudaOccupancyMaxActiveBlocksPerMultiprocessor(&max_block_per_sm, (void*) gpu_fbarebones_solve, CUDA_THREADS_PER_BLOCK, 0); + if(cp.config.verbose_solving) { + printf("%% max_blocks_per_sm=%d\n", max_block_per_sm); + } + if(cp.config.or_nodes != 0) { + cp.stats.num_blocks = std::min(max_block_per_sm * deviceProp.multiProcessorCount, (int)cp.config.or_nodes); + if(cp.config.verbose_solving >= 1 && cp.stats.num_blocks < cp.config.or_nodes) { + printf("%% WARNING: -or %d is too high on your GPU architecture, it has been reduced to %d.\n", (int)cp.config.or_nodes, cp.stats.num_blocks); + } + } + else { + cp.stats.num_blocks = max_block_per_sm * deviceProp.multiProcessorCount; + } + + /** II. Number of subproblems. */ + cp.stats.print_stat("subproblems_power", cp.config.subproblems_power); + if(cp.config.subproblems_power == -1) { + cp.config.subproblems_power = 0; + while((1 << cp.config.subproblems_power) < cp.config.subproblems_factor * cp.stats.num_blocks) { + cp.config.subproblems_power++; + } + } + + /** III. Size of the heap global memory. + * The estimation is very conservative, normally we should not run out of memory. + * */ + size_t store_bytes = gpu_sizeof() + gpu_sizeof>() + cp.store->vars() * gpu_sizeof(); + size_t iprop_bytes = gpu_sizeof() + gpu_sizeof>() + cp.iprop->num_deductions() * gpu_sizeof() + gpu_sizeof(); + // size_t gradient_bytes = sizeof(float) * static_cast(cp.store->vars()) * 4; + size_t mem_per_block = gpu_sizeof() + + store_bytes * size_t{3} // current, root, best. + + store_bytes * size_t{2} // search strategies + + iprop_bytes * size_t{2} + + cp.iprop->num_deductions() * size_t{4} * gpu_sizeof() // fixpoint engine + // + gradient_bytes + + (gpu_sizeof() + gpu_sizeof>()) * size_t{MAX_SEARCH_DEPTH}; + // size_t estimated_global_mem = gpu_sizeof() + store_bytes * size_t{5} + iprop_bytes + gradient_bytes + + // gpu_sizeof(); + size_t estimated_global_mem = gpu_sizeof() + store_bytes * size_t{5} + iprop_bytes + + gpu_sizeof(); + + size_t mem_for_blocks = deviceProp.totalGlobalMem - estimated_global_mem - (deviceProp.totalGlobalMem / 100 * 10); + cp.stats.num_blocks = std::max(size_t{1}, std::min(mem_for_blocks / mem_per_block, static_cast(cp.stats.num_blocks))); + estimated_global_mem += cp.stats.num_blocks * mem_per_block; + if(estimated_global_mem > deviceProp.totalGlobalMem / 100 * 90) { + printf("%% WARNING: The estimated global memory is larger than 90%% of the total global memory.\n\ +%% It is possible to run out of memory during solving.\n"); + } + CUDAEX(cudaDeviceSetLimit(cudaLimitMallocHeapSize, deviceProp.totalGlobalMem / 100 * 97)); + cp.stats.print_memory_statistics(cp.config.verbose_solving, "heap_memory", estimated_global_mem); + cp.stats.print_memory_statistics(cp.config.verbose_solving, "mem_per_block", mem_per_block); + cp.stats.print_memory_statistics(cp.config.verbose_solving, "total_global_mem_bytes", deviceProp.totalGlobalMem); + + // We still need to improve this, for some large problems, it is required to avoid running out of memory. + cp.stats.num_blocks = std::min(cp.stats.num_blocks, 200000000 / cp.store->vars()); + cp.stats.print_stat("num_blocks", cp.stats.num_blocks); + + /** IV. Increase the stack if requested by the user. */ + if(config.stack_kb != 0) { + CUDAEX(cudaDeviceSetLimit(cudaLimitStackSize, config.stack_kb*1000)); + // The stack allocated depends on the maximum number of threads per SM, not on the actual number of threads per block. + size_t total_stack_size = deviceProp.multiProcessorCount * deviceProp.maxThreadsPerMultiProcessor * config.stack_kb * 1000; + cp.stats.print_memory_statistics(cp.config.verbose_solving, "stack_memory", total_stack_size); + } + + /** V. Configure the shared memory size. */ + int blocks_per_sm = std::max(1, (cp.stats.num_blocks + deviceProp.multiProcessorCount - 1) / deviceProp.multiProcessorCount); + MemoryConfig mem_config; + if(config.only_global_memory) { + mem_config = MemoryConfig(store_bytes, iprop_bytes); + } + else { + mem_config = MemoryConfig((void*) gpu_fbarebones_solve, config.verbose_solving, blocks_per_sm, store_bytes, iprop_bytes); + } + mem_config.print_mzn_statistics(config, cp.stats); + return mem_config; +} + +__global__ void initialize_global_data( + UnifiedData* unified_data, + bt::unique_ptr* grid_data_ptr) +{ + *grid_data_ptr = bt::make_unique(unified_data->root); +} + +#define TIMEPOINT(KIND) \ + if(threadIdx.x == 0) { \ + block_data.timer = block_data.stats.stop_timer(Timer::KIND, block_data.timer); \ + } + +__global__ void gpu_fbarebones_solve(UnifiedData* unified_data, GridData* grid_data, Ort::Session& session) { + extern __shared__ unsigned char shared_mem[]; + auto& config = unified_data->root.config; + BlockData& block_data = grid_data->blocks[blockIdx.x]; + if(threadIdx.x == 0 && blockIdx.x == 0 && config.verbose_solving) { + printf("%% GPU kernel started, starting solving...\n"); + } + + /** A. Initialization the block data and the fixpoint engine. */ + + block_data.allocate(*unified_data, *grid_data, shared_mem); + __syncthreads(); + FProp& iprop = *block_data.iprop; +#ifdef TURBO_NO_ENTAILED_PROP_REMOVAL + __shared__ BlockAsynchronousFixpointGPU fp_engine; +#else + __shared__ FixpointSubsetGPU, bt::global_allocator, CUDA_THREADS_PER_BLOCK> fp_engine; + fp_engine.init(iprop.num_deductions()); +#endif + /** This shared variable is necessary to avoid multiple threads to read into `unified_data.stop.test()`, + * potentially reading different values and leading to deadlock. */ + __shared__ bool stop; + __shared__ bool has_changed; + __shared__ bool is_leaf_node; + __shared__ int remaining_depth; + stop = false; + auto group = cooperative_groups::this_thread_block(); + if(threadIdx.x == 0) { + block_data.timer = block_data.stats.start_timer_device(); + block_data.start_time = block_data.timer; + } + __syncthreads(); + + /** B. Start the main dive and solve loop. */ + size_t num_subproblems = unified_data->root.stats.eps_num_subproblems; + while(block_data.subproblem_idx < num_subproblems && !stop) { + if(config.verbose_solving >= 2 && threadIdx.x == 0) { + grid_data->print_lock.acquire(); + printf("%% Block %d solves subproblem num %" PRIu64 "\n", blockIdx.x, block_data.subproblem_idx); + grid_data->print_lock.release(); + } + + // C. Restoring the current state to the root node. + + block_data.current_strategy = 0; + block_data.next_unassigned_var = 0; + block_data.depth = 0; + unified_data->root.store->copy_to(group, *block_data.store); +#ifndef TURBO_NO_ENTAILED_PROP_REMOVAL + fp_engine.reset(iprop.num_deductions()); +#endif + __syncthreads(); + + // D. Dive into the search tree until we reach the target subproblem. + remaining_depth = config.subproblems_power; + if(threadIdx.x == 0) { + block_data.dive_timer = block_data.stats.start_timer_device(); + is_leaf_node = false; + } + __syncthreads(); + while(remaining_depth > 0 && !is_leaf_node && !stop) { + __syncthreads(); + propagate(*unified_data, *grid_data, block_data, fp_engine, stop, has_changed, is_leaf_node, session); + __syncthreads(); + if(!is_leaf_node) { + block_data.split(has_changed, grid_data->search_strategies, config.epsilon); + __syncthreads(); + // Split was not able to split a domain. It means that the search strategy is not complete due to unsplittable infinite domains. + // We skip the subtree, and set exhaustive to `false`. + if(block_data.decisions[0].var.is_untyped()) { + is_leaf_node = true; + block_data.stats.exhaustive = false; + if(threadIdx.x == 0 && config.verbose_solving >= 1) { printf("%% WARNING: infinite element detected during branching, search is not exhaustive\n");} + } + else if(threadIdx.x == 0) { + --remaining_depth; + // We do not record the decisions when diving. + --block_data.depth; + /** We commit to one of the branches depending on the current value on the path. + * Suppose the depth is 3, the path is "010" we are currently at `remaining_depth = 1`. + * We must extract the bit "1", and we do so by standard bitwise manipulation. + * Whenever the branch_idx is 0 means to take the left branch, and 1 means to take the right branch. + */ + size_t branch_idx = (block_data.subproblem_idx & (size_t{1} << remaining_depth)) >> remaining_depth; + /** We immediately commit to the branch. */ + // printf("split on %d (", block_data.decisions[0].var.vid()); block_data.store->project(block_data.decisions[0].var).print(); printf(")\n"); + block_data.store->embed(block_data.decisions[0].var, block_data.decisions[0].children[branch_idx]); + } + } + __syncthreads(); + } + if(threadIdx.x == 0) { + block_data.stats.stop_timer(Timer::DIVE, block_data.dive_timer); + } + // E. Skip subproblems that are not reachable. + + /** If we reached a leaf node before the subproblem was reached, then it means a whole subtree should be skipped. */ + if(is_leaf_node && !stop) { + /** To skip all the paths of the subtree obtained, we perform bitwise operations. + * Suppose the current path is "00" turn left two times, and the following search tree: + * * depth = 0 + * / \ + * 0 1 depth = 1 + * / \ / \ + * 00 01 10 11 depth = 2 + * + * If we detect a leaf node at depth 1, after only one left turn, we must skip the remaining of the subtree, in particular to avoid exploring the path "01". + * To achieve that, we take the current path "00", shift it to the right by 1 (essentially erasing the path that has not been explored), increment it to go to the next subtree (at the same depth), and shift it back to the left to reach the first subproblem of the subtree. + * Cool huh? + */ + if(threadIdx.x == 0) { + size_t next_subproblem_idx = ((block_data.subproblem_idx >> remaining_depth) + size_t{1}) << remaining_depth; + // Make sure the subtree is skipped. + while(grid_data->next_subproblem.meet(ZLB(next_subproblem_idx))) {} + /** It is possible that other blocks skip similar subtrees. + * Hence, we only count the subproblems skipped by the block solving the left most subproblem. */ + if((block_data.subproblem_idx & ((size_t{1} << remaining_depth) - size_t{1})) == size_t{0}) { + block_data.stats.eps_skipped_subproblems += next_subproblem_idx - block_data.subproblem_idx; + } + } + } + else if(!stop) { + + // F. Solve the current subproblem. + + // We skip the remaining of the EPS strategy if there is any. + if(threadIdx.x == 0 && grid_data->has_eps_strategy) { + block_data.current_strategy = battery::max(1, block_data.current_strategy); + block_data.next_unassigned_var = 0; + } + + while(!stop) { + + // I. Propagate the current node. + propagate(*unified_data, *grid_data, block_data, fp_engine, stop, has_changed, is_leaf_node, session); + __syncthreads(); + + // II. Branching + + if(!is_leaf_node) { + // If we are at the root of the current subproblem, we create a snapshot for future backtracking. + if(block_data.depth == 0) { + block_data.store->copy_to(group, *block_data.root_store); + if(threadIdx.x == 0) { + block_data.snapshot_root_strategy = block_data.current_strategy; + block_data.snapshot_next_unassigned_var = block_data.next_unassigned_var; + } + } + __syncthreads(); + block_data.split(has_changed, grid_data->search_strategies, config.epsilon); + __syncthreads(); + // Split was not able to split a domain. It means that the search strategy is not complete due to unsplittable infinite domains. + // We trigger backtracking, and set exhaustive to `false`. + if(block_data.decisions[block_data.depth - 1].var.is_untyped()) { + is_leaf_node = true; + block_data.stats.exhaustive = false; + if(threadIdx.x == 0 && config.verbose_solving >= 1) { printf("%% WARNING: infinite element detected during branching, search is not exhaustive\n");} + } + else if(threadIdx.x == 0) { + // Apply the decision. + // printf("split on %d (", block_data.decisions[block_data.depth-1].var.vid()); block_data.store->project(block_data.decisions[block_data.depth-1].var).print(); printf(")\n"); + block_data.store->embed(block_data.decisions[block_data.depth-1].var, block_data.decisions[block_data.depth-1].next()); + // printf("left decision: %d [", block_data.decisions[block_data.depth - 1].var.vid()); block_data.decisions[block_data.depth - 1].current().print(); printf("]\n"); + } + } + + // III. Backtracking + + if(is_leaf_node) { + // Leaf node at root. + if(block_data.depth == 0) { + break; + } + if(threadIdx.x == 0) { + block_data.depth = block_data.decisions[block_data.depth-1].ropes[block_data.decisions[block_data.depth-1].current_idx]; + } + __syncthreads(); + // Check if there is no more node to visit. + if(block_data.depth == -1) { + break; + } + // Restore from root by copying the store and re-applying all decisions from root to block_data.depth-1. +#ifndef TURBO_NO_ENTAILED_PROP_REMOVAL + fp_engine.reset(iprop.num_deductions()); +#endif + block_data.root_store->copy_to(group, *block_data.store); + // __syncthreads(); + // if(threadIdx.x == 0) { + // printf("%d: restoring store: ", block_data.depth); block_data.store->print(); printf("\n"); + // } + // __syncthreads(); + if(threadIdx.x == 0) { + has_changed = true; + } + __syncthreads(); + while(has_changed) { + __syncthreads(); + if(threadIdx.x == 0) { + has_changed = false; + } + __syncthreads(); + for(int i = threadIdx.x; i < block_data.depth - 1; i += blockDim.x) { + if(block_data.store->embed(block_data.decisions[i].var, block_data.decisions[i].current())) { + has_changed = true; + } + } + __syncthreads(); + } + if(threadIdx.x == 0) { + block_data.store->embed(block_data.decisions[block_data.depth - 1].var, block_data.decisions[block_data.depth - 1].next()); + // printf("right decision: %d [", block_data.decisions[block_data.depth - 1].var.vid()); block_data.decisions[block_data.depth - 1].current().print(); printf("]\n"); + block_data.current_strategy = block_data.snapshot_root_strategy; + block_data.next_unassigned_var = block_data.snapshot_next_unassigned_var; + } + // __syncthreads(); + // if(threadIdx.x == 0) { + // printf("%d: reapplied decisions: ", block_data.depth); block_data.store->print(); printf("\n"); + // } + // __syncthreads(); + } + } + /** If we didn't stop solving because of an external interruption, we increase the number of subproblems solved. */ + if(threadIdx.x == 0 && block_data.stats.nodes < config.stop_after_n_nodes + && !unified_data->stop.test()) + { + block_data.stats.eps_solved_subproblems += 1; + } + } + + // G. Move to the next subproblem. + + /** We prepare the block to solve the next problem. + * We update the subproblem index to the next subproblem to solve. */ + if(threadIdx.x == 0 && !stop) { + /** To avoid that several blocks solve the same subproblem, we use an atomic post-increment. */ + block_data.subproblem_idx = grid_data->next_subproblem.atomic()++; + /** The following commented code is completely valid and does not use atomic post-increment. + * But honestly, we kinda need more performance so... let's avoid reexploring subproblems. */ + // subproblem_idx = grid_data->next_subproblem.value(); + // grid_data->next_subproblem.meet(FLB(subproblem_idx + size_t{1})); + } + __syncthreads(); + } + if(threadIdx.x == 0) + { + if(block_data.stats.nodes < config.stop_after_n_nodes && !unified_data->stop.test()) { + block_data.stats.num_blocks_done = 1; + } + block_data.stats.timers.update_timer(Timer::FIRST_BLOCK_IDLE, block_data.start_time); + block_data.stats.cumulative_time_block = block_data.stats.timers.time_of(Timer::FIRST_BLOCK_IDLE); + } + __syncthreads(); +#ifndef TURBO_NO_ENTAILED_PROP_REMOVAL + fp_engine.destroy(); +#endif + block_data.deallocate_shared_data(); + __syncthreads(); +} + +void back_propagation(BlockData& block_data, Ort::Session& session) { + // Step 1. + Ort::MemoryInfo cuda_mem_info("Cuda", OrtAllocatorType::OrtDeviceAllocator, 0, OrtMemType::OrtMemTypeDefault); + + // Step 2. + Ort::TypeInfo input_type_info = session.GetInputTypeInfo(0); + battery::vector input_dims = input_type_info.GetTensorTypeAndShapeInfo().GetShape(); + size_t total_elements = 1; + for (size_t i = 0; i < input_dims.size(); ++i) { + total_elements *= input_dims[i]; + } + + Ort::Value input_mid_tensor = Ort::Value::CreateTensor( + cuda_mem_info, + block_data.h_mid_gradients, + total_elements, + input_dims.data(), + input_dims.size() + ); + // Ort::Value input_lb_tensor = Ort::Value::CreateTensor( + // cuda_mem_info, + // block_data.h_lb_gradients, + // total_elements, + // input_dims.data(), + // input_dims.size() + // ); + // Ort::Value input_ub_tensor = Ort::Value::CreateTensor( + // cuda_mem_info, + // block_data.h_ub_gradients, + // total_elements, + // input_dims.data(), + // input_dims.size() + // ); + + // Step 3. + Ort::RunOptions run_opts; + Ort::AllocatorWithDefaultOptions ort_allocator; + Ort::AllocatedStringPtr input_name_alloc = session.GetInputNameAllocated(0, ort_allocator); + std::string real_input_name = input_name_alloc.get(); + Ort::AllocatedStringPtr output_name_alloc = session.GetOutputNameAllocated(0, ort_allocator); + std::string real_output_name = output_name_alloc.get(); + + const char* input_names[] = { real_input_name.c_str() }; + const char* output_names[] = { real_output_name.c_str() }; + const char* const* input_names_ptr = input_names; + const char* const* output_names_ptr = output_names; + + // Step 4. + std::vector output_mid_tensors; + // std::vector output_lb_tensors; + // std::vector output_ub_tensors; + output_mid_tensors = session.Run( + run_opts, + input_names_ptr, + &input_mid_tensor, + 1, + output_names_ptr, + 1 + ); + // output_lb_tensors = session.Run( + // run_opts, + // input_names_ptr, + // &input_lb_tensor, + // 1, + // output_names_ptr, + // 1 + // ); + // output_ub_tensors = session.Run( + // run_opts, + // input_names_ptr, + // &input_ub_tensor, + // 1, + // output_names_ptr, + // 1 + // ); + + // Step 5. + block_data.h_mid_gradients = output_mid_tensors[0].GetTensorMutableData(); + // block_data.h_lb_gradients = output_lb_tensors[0].GetTensorMutableData(); + // block_data.h_ub_gradients = output_ub_tensors[0].GetTensorMutableData(); + + // TODO: combine these gradients together. just use average + // we also have to consider floating-point errors. + // to simplicitly, we use only upper-towards rounding function + for(size_t i = 0; i < block_data.num_h_gradients; ++i){ + // block_data.h_gradients[i] = battery::div_up(battery::add_up(battery::add_up(block_data.h_mid_gradients[i], block_data.h_lb_gradients[i]), block_data.h_ub_gradients[i]), float{3.0}); + block_data.h_gradients[i] = block_data.h_mid_gradients[i]; + // block_data.h_gradients[i] = block_data.h_ub_gradients[i]; + } + + cudaDeviceSynchronize(); +} + +template +__device__ INLINE void propagate(UnifiedData& unified_data, GridData& grid_data, BlockData& block_data, + FPEngine& fp_engine, bool& stop, bool& has_changed, bool& is_leaf_node, Ort::Session& session) +{ + __shared__ int warp_iterations[CUDA_THREADS_PER_BLOCK/32]; + warp_iterations[threadIdx.x / 32] = 0; + auto& config = unified_data.root.config; + FProp& iprop = *block_data.iprop; + auto group = cooperative_groups::this_thread_block(); + + TIMEPOINT(SEARCH); + if(threadIdx.x == 0) { + is_leaf_node = false; + } + + // II. Compute the fixpoint of the current node. + int fp_iterations; +#ifdef TURBO_NO_ENTAILED_PROP_REMOVAL + int num_active = iprop.num_deductions(); +#else + int num_active = fp_engine.num_active(); +#endif + switch(config.fixpoint) { + case FixpointKind::AC1: { + fp_iterations = fp_engine.fixpoint( +#ifdef TURBO_NO_ENTAILED_PROP_REMOVAL + iprop.num_deductions(), +#endif + [&](int i){ return iprop.fdeduce(i, config.epsilon); }, + [&](){ return iprop.is_bot(); }); + if(threadIdx.x == 0) { + block_data.stats.num_deductions += fp_iterations * num_active; + } + break; + } + case FixpointKind::WAC1: { + if(num_active <= config.wac1_threshold) { + fp_iterations = fp_engine.fixpoint( +#ifdef TURBO_NO_ENTAILED_PROP_REMOVAL + iprop.num_deductions(), +#endif + [&](int i){ return iprop.fdeduce(i, config.epsilon); }, + [&](){ return iprop.is_bot(); }); + if(threadIdx.x == 0) { + block_data.stats.num_deductions += fp_iterations * num_active; + } + } + else { + fp_iterations = fp_engine.fixpoint( +#ifdef TURBO_NO_ENTAILED_PROP_REMOVAL + iprop.num_deductions(), +#endif + [&](int i){ return fwarp_fixpoint(iprop, i, warp_iterations, config.epsilon); }, + [&](){ return iprop.is_bot(); }); + if(threadIdx.x == 0) { + for(int i = 0; i < CUDA_THREADS_PER_BLOCK/32; ++i) { + block_data.stats.num_deductions += warp_iterations[i] * 32; + } + } + } + break; + } + } + TIMEPOINT(FIXPOINT); + + const auto current_strat = block_data.current_strategy; + const auto& strat = grid_data.search_strategies[current_strat]; + const auto& store = *block_data.store; + + // III. Analyze the result of propagation + if(!iprop.is_bot()) { + if(threadIdx.x == 0) { + has_changed = false; + } + __syncthreads(); + // This is an underapproximation caes. + for(int i = (int)group.thread_rank(); i < strat.vars.size(); i+=group.num_threads()){ + if(store[i].lb().value() != store[i].ub().value()){ + has_changed = true; + break; + } + } + __syncthreads(); + num_active = has_changed ? 1 : 0; + TIMEPOINT(SELECT_FP_FUNCTIONS); + if (num_active == 0) { + // This is SAT case. + is_leaf_node = true; + if(threadIdx.x == 0) { + block_data.stats.timers.update_timer(Timer::LATEST_BEST_OBJ_FOUND, block_data.start_time); + } + block_data.store->copy_to(group, *block_data.inner_box); + if(threadIdx.x == 0) { + block_data.stats.solutions++; + unified_data.stop.test_and_set(); + } + } + // else if(strat.var_order == VariableOrder::GRA_ANTI_FIRST_FAIL){ + // // Obtain graident from the network by lbs, midpoints, and ubs. + // // We need to split this node. + // // By applying back_propagation(), we can have the latest gradient information. + // // This gradient information might not work. The code itself is correct, but it might not effective. + // for(int i = (int)group.thread_rank(); i < strat.vars.size(); i += group.num_threads()){ + // block_data.h_mid_gradients[i] = battery::midpoint(store[i].lb().value(), store[i].ub().value()); + // block_data.h_lb_gradients[i] = store[i].lb().value(); + // block_data.h_ub_gradients[i] = store[i].ub().value(); + // } + // __syncthreads(); + // if(threadIdx.x == 0) { + // back_propagation(block_data, session); + // } + // } + } + else { + // This is unknown checking. + // If is_bot() is true /\ exists at least 1 the width == 0.0, then it is identified as an unknown box. + // -> It can be simplified to check if it uses UASS or not. + // -> If we have applied UASS, it implies that there exists at least 1 variable is assigned. + // If is_bot() is true /\ all the widths != 0.0, then it is pruned by propagation, not underapproximation. + is_leaf_node = true; + if (threadIdx.x == 0) { + has_changed = block_data.is_uass; + } + __syncthreads(); + // for (int i = (int)group.thread_rank(); i < strat.vars.size(); i += group.num_threads()){ + // if(store[i].lb().value() == store[i].ub().value()){ + // has_changed = true; + // } + // } + // __syncthreads(); + num_active = has_changed ? 0 : 1; + TIMEPOINT(SELECT_FP_FUNCTIONS); + if (num_active == 0) { + if(threadIdx.x == 0) { + block_data.stats.timers.update_timer(Timer::LATEST_BEST_OBJ_FOUND, block_data.start_time); + block_data.stats.unknowns++; + } + } + } + + if(threadIdx.x == 0) { + block_data.stats.fixpoint_iterations += fp_iterations; + block_data.stats.nodes++; + block_data.stats.fails += (iprop.is_bot() ? 1 : 0); + block_data.stats.depth_max = battery::max(block_data.stats.depth_max, block_data.depth); + + // IV. Checking stopping conditions. + + if(block_data.stats.nodes >= config.stop_after_n_nodes + || unified_data.stop.test() + || block_data.stats.solutions != 0) + { + block_data.stats.exhaustive = false; + stop = true; + } + } +} + +__global__ void reduce_blocks(UnifiedData* unified_data, GridData* grid_data) { + auto& root = unified_data->root; + for(int i = 0; i < grid_data->blocks.size(); ++i) { + root.stats.meet(grid_data->blocks[i].stats); + int64_t& grid_first_block_idle = root.stats.timers.time_of(Timer::FIRST_BLOCK_IDLE); + int64_t block_idle = grid_data->blocks[i].stats.timers.time_of(Timer::FIRST_BLOCK_IDLE); + if(grid_first_block_idle > block_idle) { + grid_first_block_idle = block_idle; + } + } + for(int i = 0; i < grid_data->blocks.size(); ++i) { + auto& block = grid_data->blocks[i]; + if(block.stats.solutions > 0) { + if(root.bab->is_satisfaction()) { + // FIXME: We might have more than one solution to remember. + block.inner_box->extract(*root.best); + // for(int j = 0; j < block.inner_boxes.size(); ++j) { + // block.inner_boxes[j].extract(*root.best); + // root.inner_boxes.push_back(*root.best); + // if (j >= 10) break; + // } + break; + } + } + } +} + +__global__ void deallocate_global_data(bt::unique_ptr* grid_data) { + grid_data->reset(); +} + +#endif // TURBO_IPC_ABSTRACT_DOMAIN +#endif // __CUDACC__ + +#if defined(TURBO_IPC_ABSTRACT_DOMAIN) || !defined(__CUDACC__) + +void fbarebones_dive_and_solve(const Configuration& config) { +#ifdef TURBO_IPC_ABSTRACT_DOMAIN + std::cerr << "-arch fbarebones does not support IPC abstract domain." << std::endl; +#else + std::cerr << "You must use a CUDA compiler (nvcc or clang) to compile Turbo on GPU." << std::endl; +#endif +} + +#endif + +} // namespace fbarebones + +#endif // TURBO_FBAREBONES_DIVE_AND_SOLVE_HPP diff --git a/include/nnv.hpp b/include/nnv.hpp new file mode 100644 index 00000000..5046c75a --- /dev/null +++ b/include/nnv.hpp @@ -0,0 +1,80 @@ +// Copyright 2025 Yi-Nung Tsao + +#ifndef TURBO_NNV_HPP +#define TURBO_NNV_HPP + +#include "lala/onnx_parser.hpp" +#include "lala/smt_parser.hpp" +#include "lala/solver_output.hpp" + +namespace lala { + +namespace impl { + +template +class NNV { + using allocator_type = Allocator; + using F = TFormula; + using FSeq = typename F::Sequence; + + bool is_nnv; + battery::vector& input_neurons; + battery::vector& hidden_neurons; + SolverOutput& output; + +public: + NNV(battery::vector& input_neurons, battery::vector& hidden_neurons, SolverOutput& output, bool is_nnv): input_neurons(input_neurons), hidden_neurons(hidden_neurons), output(output), is_nnv(is_nnv) {} + + battery::shared_ptr make_nnv_formulas(const std::string& onnx_path, const std::string& vnnlib_path) { + FSeq seq; + seq.push_back(std::move(parse_onnx(onnx_path, input_neurons, hidden_neurons, output))); + seq.push_back(std::move(parse_smt(vnnlib_path, output, is_nnv))); + return battery::make_shared(std::move(F::make_nary(AND, std::move(seq)))); + } +}; + +template +class SMT2 { + using allocator_type = Allocator; + using F = TFormula; + using FSeq = typename F::Sequence; + + bool is_nnv; + SolverOutput& output; + +public: + SMT2(SolverOutput& output, bool is_nnv): output(output), is_nnv(is_nnv) {} + + battery::shared_ptr make_smt2_formulas(const std::string& smt2_path) { + return battery::make_shared(std::move(parse_smt(smt2_path, output, is_nnv))); + } +}; +} // namespace impl + +template +battery::shared_ptr, Allocator> parse_nnv(const std::string& onnx_path, const std::string& vnnlib_path) { + impl::NNV nnv; + return nnv.make_nnv_formulas(onnx_path, vnnlib_path); +} + +template +battery::shared_ptr, Allocator> parse_nnv(const std::string& onnx_path, const std::string& vnnlib_path, battery::vector& input_neurons, battery::vector& hidden_neurons, SolverOutput& output, bool is_nnv) { + impl::NNV nnv(input_neurons, hidden_neurons, output, is_nnv); + return nnv.make_nnv_formulas(onnx_path, vnnlib_path); +} + +template +battery::shared_ptr, Allocator> parse_smt2(const std::string& smt2_path) { + impl::SMT2 smt2; + return smt2.make_smt2_formulas(smt2_path); +} + +template +battery::shared_ptr, Allocator> parse_smt2(const std::string& smt2_path, SolverOutput& output, bool is_nnv) { + impl::SMT2 smt2(output, is_nnv); + return smt2.make_smt2_formulas(smt2_path); +} + +} // namespace lala + +#endif \ No newline at end of file diff --git a/include/statistics.hpp b/include/statistics.hpp index d50b43c1..ccaa8a12 100644 --- a/include/statistics.hpp +++ b/include/statistics.hpp @@ -142,6 +142,7 @@ struct Statistics { size_t nodes; size_t fails; size_t solutions; + size_t unknowns; int depth_max; bool exhaustive; size_t eps_num_subproblems; @@ -156,7 +157,7 @@ struct Statistics { CUDA Statistics(size_t variables, size_t constraints, bool optimization, bool print_statistics): variables(variables), constraints(constraints), optimization(optimization), print_statistics(print_statistics), num_blocks(1), - nodes(0), fails(0), solutions(0), + nodes(0), fails(0), solutions(0), unknowns(0), depth_max(0), exhaustive(true), eps_solved_subproblems(0), eps_num_subproblems(1), eps_skipped_subproblems(0), num_blocks_done(0), fixpoint_iterations(0), num_deductions(0), @@ -171,7 +172,7 @@ struct Statistics { CUDA Statistics(const Statistics& other): variables(other.variables), constraints(other.constraints), optimization(other.optimization), print_statistics(other.print_statistics), num_blocks(other.num_blocks), - nodes(other.nodes), fails(other.fails), solutions(other.solutions), + nodes(other.nodes), fails(other.fails), solutions(other.solutions), unknowns(other.unknowns), depth_max(other.depth_max), exhaustive(other.exhaustive), eps_solved_subproblems(other.eps_solved_subproblems), eps_num_subproblems(other.eps_num_subproblems), eps_skipped_subproblems(other.eps_skipped_subproblems), num_blocks_done(other.num_blocks_done), @@ -184,6 +185,7 @@ struct Statistics { nodes += other.nodes; fails += other.fails; solutions += other.solutions; + unknowns += other.unknowns; depth_max = battery::max(depth_max, other.depth_max); exhaustive = exhaustive && other.exhaustive; eps_solved_subproblems += other.eps_solved_subproblems; @@ -345,6 +347,7 @@ struct Statistics { print_timing_stat("initTime", Timer::PREPROCESSING); print_timing_stat("solveTime", Timer::OVERALL); print_stat("num_solutions", solutions); + print_stat("num_unknowns", unknowns); print_stat("eps_num_subproblems", eps_num_subproblems); print_stat("eps_solved_subproblems", eps_solved_subproblems); print_stat("eps_skipped_subproblems", eps_skipped_subproblems); @@ -399,13 +402,13 @@ struct Statistics { } else { assert(solutions == 0); - if(exhaustive) { + if(exhaustive && unknowns == 0) { printf("=====UNSATISFIABLE=====\n"); } else if(optimization) { printf("=====UNBOUNDED=====\n"); } - else { + else if (unknowns != 0) { printf("=====UNKNOWN=====\n"); } } diff --git a/nnv/onnx/test_nano.onnx b/nnv/onnx/test_nano.onnx new file mode 100644 index 00000000..b04356a8 Binary files /dev/null and b/nnv/onnx/test_nano.onnx differ diff --git a/nnv/onnx/test_sat.onnx b/nnv/onnx/test_sat.onnx new file mode 100644 index 00000000..40021f68 Binary files /dev/null and b/nnv/onnx/test_sat.onnx differ diff --git a/nnv/onnx/test_small.onnx b/nnv/onnx/test_small.onnx new file mode 100644 index 00000000..c854063f Binary files /dev/null and b/nnv/onnx/test_small.onnx differ diff --git a/nnv/onnx/test_tiny.onnx b/nnv/onnx/test_tiny.onnx new file mode 100644 index 00000000..6d22765c Binary files /dev/null and b/nnv/onnx/test_tiny.onnx differ diff --git a/nnv/onnx/test_unsat.onnx b/nnv/onnx/test_unsat.onnx new file mode 100644 index 00000000..0357f37d Binary files /dev/null and b/nnv/onnx/test_unsat.onnx differ diff --git a/nnv/vnnlib/test_nano.vnnlib b/nnv/vnnlib/test_nano.vnnlib new file mode 100644 index 00000000..bf8911b5 --- /dev/null +++ b/nnv/vnnlib/test_nano.vnnlib @@ -0,0 +1,7 @@ +(declare-const X_0 Real) +(declare-const Y_0 Real) + +(assert (>= X_0 -1.0)) +(assert (<= X_0 1.0)) + +(assert (>= Y_0 -1.0)) diff --git a/nnv/vnnlib/test_prop.vnnlib b/nnv/vnnlib/test_prop.vnnlib new file mode 100644 index 00000000..536d62fc --- /dev/null +++ b/nnv/vnnlib/test_prop.vnnlib @@ -0,0 +1,40 @@ +; test property is acas xu property 3 + +(declare-const X_0 Real) +(declare-const X_1 Real) +(declare-const X_2 Real) +(declare-const X_3 Real) +(declare-const X_4 Real) + +(declare-const Y_0 Real) +(declare-const Y_1 Real) +(declare-const Y_2 Real) +(declare-const Y_3 Real) +(declare-const Y_4 Real) + +; input constraints +; Unscaled Input 0: (1500, 1800) +(assert (<= X_0 -0.29855281193475053)) +(assert (>= X_0 -0.30353115613746867)) + +; Unscaled Input 1: (-0.06, 0.06) +(assert (<= X_1 0.009549296585513092)) +(assert (>= X_1 -0.009549296585513092)) + +; Unscaled Input 2: (3.1, 3.1415926535) +(assert (<= X_2 0.49999999998567607)) +(assert (>= X_2 0.4933803235848431)) + +; Unscaled Input 3: (980, 1200) +(assert (<= X_3 0.5)) +(assert (>= X_3 0.3)) + +; Unscaled Input 4: (960, 1200) +(assert (<= X_4 0.5)) +(assert (>= X_4 0.3)) + +; output constraints (property 3, sat if CoC is minimal) +(assert (<= Y_0 Y_1)) +(assert (<= Y_0 Y_2)) +(assert (<= Y_0 Y_3)) +(assert (<= Y_0 Y_4)) diff --git a/nnv/vnnlib/test_small.vnnlib b/nnv/vnnlib/test_small.vnnlib new file mode 100644 index 00000000..7d42a5c9 --- /dev/null +++ b/nnv/vnnlib/test_small.vnnlib @@ -0,0 +1,8 @@ +(declare-const X_0 Real) +(declare-const Y_0 Real) + +(assert (>= X_0 -1.0)) +(assert (<= X_0 1.0)) + +(assert (>= Y_0 100.0)) + diff --git a/nnv/vnnlib/test_tiny.vnnlib b/nnv/vnnlib/test_tiny.vnnlib new file mode 100644 index 00000000..0a3663d2 --- /dev/null +++ b/nnv/vnnlib/test_tiny.vnnlib @@ -0,0 +1,5 @@ +(declare-const X_0 Real) +(declare-const Y_0 Real) +(assert (or + (and (>= X_0 -1) (<= X_0 1) (>= Y_0 100)) +)) \ No newline at end of file diff --git a/src/config.cpp b/src/config.cpp index 832cdff3..97354eff 100644 --- a/src/config.cpp +++ b/src/config.cpp @@ -9,7 +9,7 @@ #include void usage_and_exit(const std::string& program_name) { - std::cout << "usage: " << program_name << " [-t 2000] [-a] [-n 10] [-i] [-f] [-s] [-v] [-p ] [-arch ] [-p 48] [-or 48] [-sub 12] [-stack 100] [-fp ] [-wac1_threshold 0] [-eps_var_order ] [-eps_value_order ] [-seed 0] [-network_analysis] [-cutnodes 0] [-disable_simplify] [-force_ternarize] [-globalmem] [-version 1.0.0] [xcsp3instance.xml | fzninstance.fzn]" << std::endl; + std::cout << "usage: " << program_name << " [-t 2000] [-a] [-n 10] [-i] [-f] [-s] [-v] [-p ] [-arch ] [-p 48] [-or 48] [-sub 12] [-stack 100] [-fp ] [-wac1_threshold 0] [-eps_var_order ] [-eps_value_order ] [-seed 0] [-network_analysis] [-cutnodes 0] [-disable_simplify] [-force_ternarize] [-globalmem] [-version 1.0.0] [xcsp3instance.xml | fzninstance.fzn]" << std::endl; std::cout << "\t-t 2000: Run the solver with a timeout of 2000 milliseconds." << std::endl; std::cout << "\t-timeout 2000: Same as -t, but if both -t and -timeout are specified, -timeout overrides -t." << std::endl; std::cout << "\t-a: Instructs the solver to report all solutions in the case of satisfaction problems, or print intermediate solutions of increasing quality in the case of optimisation problems." << std::endl; @@ -20,7 +20,7 @@ void usage_and_exit(const std::string& program_name) { std::cout << "\t-v: Print log messages (verbose solving) to the standard error stream." << std::endl; std::cout << "\t-ast: Print the AST of the model (useful to debug)." << std::endl; std::cout << "\t-p 48: On CPU, multithreading is not yet implemented. On GPU, equivalent to `-or 48`." << std::endl; - std::cout << "\t-arch : Choose the architecture on which the problem will be solved." << std::endl; + std::cout << "\t-arch : Choose the architecture on which the problem will be solved." << std::endl; std::cout << "\t-fp : Choose the fixpoint strategy (default: ac1 on CPU, wac1 on GPU):" << std::endl; std::cout << "\t\t ac1: All propagators are executed in parallel at each iteration." << std::endl; std::cout << "\t\t wac1: Behave as ac1 when the number of active propagators is less than wac1_threshold. Otherwise, each warp must reach a local fixpoint before executing the next 32 propagators (not compatible with -arch cpu)." << std::endl; @@ -28,6 +28,9 @@ void usage_and_exit(const std::string& program_name) { std::cout << "\t-or 48: Run the subproblems on 48 streaming multiprocessors (SMs) (only for GPU architecture). Default: -or 0 for automatic selection of the number of SMs." << std::endl; std::cout << "\t-sub 12: Create 2^12 subproblems to be solved in turns by the blocks (embarrasingly parallel search). The special value `-1` leaves Turbo to decide on the number of subproblems (at least 30 * number of blocks). Default: -sub -1." << std::endl; std::cout << "\t-subfactor 300: Create B * 300 subproblems to be solved in turns by `B` blocks (embarrasingly parallel search). Default: -subfactor 300." << std::endl; + std::cout << "\t-epsilon: Choose the epsilon value for the precision of the floating point intervals." << std::endl; + std::cout << "\t-var_order : Choose the variable ordering strategy for tree search." << std::endl; + std::cout << "\t-value_order : Choose the value ordering strategy for tree search." << std::endl; std::cout << "\t-eps_var_order : Choose the variable ordering strategy for subproblems decomposition (default: same as main search strategy)." << std::endl; std::cout << "\t-eps_value_order : Choose the value ordering strategy for subproblems decomposition (default: same as main search strategy)." << std::endl; std::cout << "\t-seed 0: Set the seed for the random number generator (default: 0)." << std::endl; @@ -40,6 +43,8 @@ void usage_and_exit(const std::string& program_name) { std::cout << "\t-force_ternarize: Force the transformation of the formula in ternary normal form, even with IPC abstract domain (note that it is enabled by default with PIR abstract domain)." << std::endl; std::cout << "\t-disable_simplify: Disable the simplification step." << std::endl; std::cout << "\t-globalmem: Store all data abstract elements in the global memory and do not try to optimise using shared memory." << std::endl; + std::cout << "\t-vnnlib_path : Path to the VNNLIB file." << std::endl; + std::cout << "\t-onnx_path : Path to the ONNX file." << std::endl; exit(EXIT_FAILURE); } @@ -90,6 +95,16 @@ class InputParser { return false; } + bool read_double(const std::string& option, double& result) { + const std::string& value = getCmdOption(option); + if(!value.empty()) { + result = std::stod(value); + tokens_read += 2; + return true; + } + return false; + } + bool read_bool(const std::string& option, bool& result) { result = cmdOptionExists(option); if(result) { @@ -123,6 +138,22 @@ class InputParser { } result = tokens.back(); } + +#ifdef WITH_NNV + void read_vnnlib_file(std::string& result){ + // Make sure there are still unread tokens + if (tokens.size() <= tokens_read) { + usage_and_exit(program_name); + } + // Take the next unread token + result = tokens[tokens_read]; + ++tokens_read; // mark it as read + } + + void read_onnx_file(std::string& result){ + read_input_file(result); + } +#endif }; Configuration parse_args(int argc, char** argv) { @@ -174,6 +205,9 @@ Configuration parse_args(int argc, char** argv) { else if(architecture == "barebones") { config.arch = Arch::BAREBONES; } + else if (architecture == "fbarebones") { + config.arch = Arch::FBAREBONES; + } else { std::cerr << "Unknown architecture -arch " << architecture << std::endl; exit(EXIT_FAILURE); @@ -193,6 +227,15 @@ Configuration parse_args(int argc, char** argv) { } } input.read_size_t("-wac1_threshold", config.wac1_threshold); + input.read_double("-epsilon", config.epsilon); + std::string var_order; + if(input.read_string("-var_order", var_order)) { + config.var_order = battery::string(var_order.data()); + } + std::string value_order; + if(input.read_string("-value_order", value_order)) { + config.value_order = battery::string(value_order.data()); + } std::string eps_var_order; if(input.read_string("-eps_var_order", eps_var_order)) { config.eps_var_order = battery::string(eps_var_order.data()); @@ -213,8 +256,25 @@ Configuration parse_args(int argc, char** argv) { if(input.read_string("-hardware", hardware)) { config.hardware = battery::string(hardware.data()); } - std::string problem_path; - input.read_input_file(problem_path); - config.problem_path = battery::string(problem_path.data()); + + std::string vnnlib_path; + if(input.read_string("-vnnlib_path", vnnlib_path)) { + config.vnnlib_path = battery::string(vnnlib_path.data()); + std::string onnx_path; + if(input.read_string("-onnx_path", onnx_path)) { + config.onnx_path = battery::string(onnx_path.data()); + config.problem_path = config.onnx_path; + } + else { + std::cerr << "Either -onnx_path or -vnnlib_path must be specified" << std::endl; + exit(EXIT_FAILURE); + } + } + else { + std::string problem_path; + input.read_input_file(problem_path); + config.problem_path = battery::string(problem_path.data()); + } + return config; } diff --git a/src/turbo.cpp b/src/turbo.cpp index 494ced30..aea141fb 100644 --- a/src/turbo.cpp +++ b/src/turbo.cpp @@ -14,8 +14,10 @@ #include "gpu_dive_and_solve.hpp" #endif -#include "hybrid_dive_and_solve.hpp" +// #include "hybrid_dive_and_solve.hpp" #include "barebones_dive_and_solve.hpp" +#include "fbarebones_dive_and_solve.hpp" +#include "jet.hpp" using namespace battery; @@ -39,9 +41,19 @@ int main(int argc, char** argv) { else if(config.arch == Arch::BAREBONES) { barebones::barebones_dive_and_solve(config); } + else if (config.arch == Arch::FBAREBONES) { + fbarebones::fbarebones_dive_and_solve(config); + } + else if (config.arch == Arch::JET) { +#ifdef __CUDACC__ + jet::jet_test(config); +#else + printf("%% JET is only supported on GPU.\n"); +#endif + } #ifndef DISABLE_HYBRID_GPU_SOLVING else if(config.arch == Arch::HYBRID) { - hybrid_dive_and_solve(config); + // hybrid_dive_and_solve(config); } #endif }