Commit c9281c25 authored by Teseo Schneider's avatar Teseo Schneider
Browse files

added mpfr multiprecision

parent 195ede0f

CMakeLists.txt

+7 −1
Original line number Diff line number Diff line
@@ -90,6 +90,11 @@ IF(NOT ${GMP_FOUND}) 
	MESSAGE(FATAL_ERROR "Cannot find GMP")

ENDIF()



find_package(MPFR)

IF(NOT ${MPFR_FOUND})

	MESSAGE(FATAL_ERROR "Cannot find MPFR")

ENDIF()



# add_library() can only be called without any source since CMake 3.11 ...

add_library(${PROJECT_NAME} src/Logger.cpp)
@@ -121,7 +126,7 @@ if(FLOAT_TETWILD_USE_FLOAT) 
	target_compile_definitions(${PROJECT_NAME} PUBLIC -DFLOAT_TETWILD_USE_FLOAT)

endif()



target_include_directories(${PROJECT_NAME} PUBLIC ${GMP_INCLUDE_DIRS})

target_include_directories(${PROJECT_NAME} PUBLIC ${GMP_INCLUDE_DIRS} ${MPFR_INCLUDES})



target_link_libraries(${PROJECT_NAME}

	PUBLIC
@@ -131,6 +136,7 @@ target_link_libraries(${PROJECT_NAME} 
		Threads::Threads

		fast_winding_number

		${GMP_LIBRARIES}

		${MPFR_LIBRARIES}

)

if(FLOAT_TETWILD_ENABLE_TBB)

	target_link_libraries(${PROJECT_NAME} PUBLIC tbb::tbb)

cmake/FindGMPfTetWild.cmake

+15 −2
Original line number Diff line number Diff line
@@ -9,12 +9,25 @@ if (GMP_INCLUDE_DIRS AND GMP_LIBRARIES) 
endif (GMP_INCLUDE_DIRS AND GMP_LIBRARIES)







#if(WIN32)

#	if(CYGWIN)

#		triwild_download_gmp_cygwin()

	#elseif(MINGW)

#		triwild_download_gmp_mingw()

	#else()

#		triwild_download_gmp_vc()

	#endif()



	#SET(GMP_WINDOWS_PATH ${THIRD_PARTY_DIR}/gmp)

#endif()





find_path(GMP_INCLUDE_DIRS NAMES gmp.h PATHS $ENV{GMP_INC} ${GMP_WINDOWS_PATH})

find_library(GMP_LIBRARIES NAMES gmp libgmp PATHS $ENV{GMP_LIB} ${GMP_WINDOWS_PATH})

find_library(GMPXX_LIBRARIES NAMES gmpxx libgmpxx PATHS $ENV{GMP_LIB} ${GMP_WINDOWS_PATH})

#MESSAGE(STATUS "GMP libs: " ${GMP_LIBRARIES} " " ${GMPXX_LIBRARIES} )





include(FindPackageHandleStandardArgs)

FIND_PACKAGE_HANDLE_STANDARD_ARGS(GMP DEFAULT_MSG GMP_INCLUDE_DIRS GMP_LIBRARIES)

cmake/FindMPFR.cmake

0 → 100644
+72 −0
Original line number Diff line number Diff line 
# Try to find the MPFR library

# See http://www.mpfr.org/

#

# This module supports requiring a minimum version, e.g. you can do

#   find_package(MPFR 2.3.0)

# to require version 2.3.0 to newer of MPFR.

#

# Once done this will define

#

#  MPFR_FOUND - system has MPFR lib with correct version

#  MPFR_INCLUDES - the MPFR include directory

#  MPFR_LIBRARIES - the MPFR library

#  MPFR_VERSION - MPFR version



# Copyright (c) 2006, 2007 Montel Laurent, <montel@kde.org>

# Copyright (c) 2008, 2009 Gael Guennebaud, <g.gael@free.fr>

# Copyright (c) 2010 Jitse Niesen, <jitse@maths.leeds.ac.uk>

# Copyright (c) 2015 Jack Poulson, <jack.poulson@gmail.com>

# Redistribution and use is allowed according to the terms of the BSD license.



find_path(MPFR_INCLUDES NAMES mpfr.h PATHS $ENV{GMPDIR} $ENV{MPFRDIR}

  ${INCLUDE_INSTALL_DIR})



# Set MPFR_FIND_VERSION to 1.0.0 if no minimum version is specified

if(NOT MPFR_FIND_VERSION)

  if(NOT MPFR_FIND_VERSION_MAJOR)

    set(MPFR_FIND_VERSION_MAJOR 1)

  endif()

  if(NOT MPFR_FIND_VERSION_MINOR)

    set(MPFR_FIND_VERSION_MINOR 0)

  endif()

  if(NOT MPFR_FIND_VERSION_PATCH)

    set(MPFR_FIND_VERSION_PATCH 0)

  endif()

  set(MPFR_FIND_VERSION

    "${MPFR_FIND_VERSION_MAJOR}.${MPFR_FIND_VERSION_MINOR}.${MPFR_FIND_VERSION_PATCH}")

endif()



if(MPFR_INCLUDES)

  # Query MPFR_VERSION

  file(READ "${MPFR_INCLUDES}/mpfr.h" _mpfr_version_header)



  string(REGEX MATCH "define[ \t]+MPFR_VERSION_MAJOR[ \t]+([0-9]+)"

    _mpfr_major_version_match "${_mpfr_version_header}")

  set(MPFR_MAJOR_VERSION "${CMAKE_MATCH_1}")

  string(REGEX MATCH "define[ \t]+MPFR_VERSION_MINOR[ \t]+([0-9]+)"

    _mpfr_minor_version_match "${_mpfr_version_header}")

  set(MPFR_MINOR_VERSION "${CMAKE_MATCH_1}")

  string(REGEX MATCH "define[ \t]+MPFR_VERSION_PATCHLEVEL[ \t]+([0-9]+)"

    _mpfr_patchlevel_version_match "${_mpfr_version_header}")

  set(MPFR_PATCHLEVEL_VERSION "${CMAKE_MATCH_1}")



  set(MPFR_VERSION

    ${MPFR_MAJOR_VERSION}.${MPFR_MINOR_VERSION}.${MPFR_PATCHLEVEL_VERSION})



  # Check whether found version exceeds minimum required

  if(${MPFR_VERSION} VERSION_LESS ${MPFR_FIND_VERSION})

    set(MPFR_VERSION_OK FALSE)

    message(STATUS "MPFR version ${MPFR_VERSION} found in ${MPFR_INCLUDES}, "

                   "but at least version ${MPFR_FIND_VERSION} is required")

  else()

    set(MPFR_VERSION_OK TRUE)

  endif()

endif()



find_library(MPFR_LIBRARIES mpfr

  PATHS $ENV{GMPDIR} $ENV{MPFRDIR} ${LIB_INSTALL_DIR})



include(FindPackageHandleStandardArgs)

find_package_handle_standard_args(MPFR DEFAULT_MSG

                                  MPFR_INCLUDES MPFR_LIBRARIES MPFR_VERSION_OK)

mark_as_advanced(MPFR_INCLUDES MPFR_LIBRARIES)

src/LocalOperations.cpp

+64 −0
Original line number Diff line number Diff line
@@ -955,6 +955,70 @@ Scalar floatTetWild::AMIPS_energy(const std::array<Scalar, 12>& T) { 
    }

}





#include <floattetwild/Multiprecision.hpp>

Scalar floatTetWild::AMIPS_energy_mp(const std::array<Scalar, 12>& T) {

    Multiprecision helper_0[12];

    helper_0[0] = T[0];

    helper_0[1] = T[1];

    helper_0[2] = T[2];

    helper_0[3] = T[3];

    helper_0[4] = T[4];

    helper_0[5] = T[5];

    helper_0[6] = T[6];

    helper_0[7] = T[7];

    helper_0[8] = T[8];

    helper_0[9] = T[9];

    helper_0[10] = T[10];

    helper_0[11] = T[11];

    Multiprecision helper_1 = helper_0[2];

    Multiprecision helper_2 = helper_0[11];

    Multiprecision helper_3 = helper_0[0];

    Multiprecision helper_4 = helper_0[3];

    Multiprecision helper_5 = helper_0[9];

    Multiprecision helper_6 = 0.577350269189626 * helper_3 - 1.15470053837925 * helper_4 + 0.577350269189626 * helper_5;

    Multiprecision helper_7 = helper_0[1];

    Multiprecision helper_8 = helper_0[4];

    Multiprecision helper_9 = helper_0[7];

    Multiprecision helper_10 = helper_0[10];

    Multiprecision helper_11 = 0.408248290463863 * helper_10 + 0.408248290463863 * helper_7 + 0.408248290463863 * helper_8 -

                       1.22474487139159 * helper_9;

    Multiprecision helper_12 = 0.577350269189626 * helper_10 + 0.577350269189626 * helper_7 - 1.15470053837925 * helper_8;

    Multiprecision helper_13 = helper_0[6];

    Multiprecision helper_14 = -1.22474487139159 * helper_13 + 0.408248290463863 * helper_3 + 0.408248290463863 * helper_4 +

                       0.408248290463863 * helper_5;

    Multiprecision helper_15 = helper_0[5];

    Multiprecision helper_16 = helper_0[8];

    Multiprecision helper_17 = 0.408248290463863 * helper_1 + 0.408248290463863 * helper_15 - 1.22474487139159 * helper_16 +

                       0.408248290463863 * helper_2;

    Multiprecision helper_18 = 0.577350269189626 * helper_1 - 1.15470053837925 * helper_15 + 0.577350269189626 * helper_2;

    Multiprecision helper_19 = 0.5 * helper_13 + 0.5 * helper_4;

    Multiprecision helper_20 = 0.5 * helper_8 + 0.5 * helper_9;

    Multiprecision helper_21 = 0.5 * helper_15 + 0.5 * helper_16;

    Multiprecision helper_22 = (helper_1 - helper_2) * (helper_11 * helper_6 - helper_12 * helper_14) -

                       (helper_7-helper_10 ) * (helper_17 * helper_6-helper_14 * helper_18) +

                       (helper_3 - helper_5) * (helper_12 * helper_17-helper_11 * helper_18);

    Multiprecision res =

                   helper_10 * (-1.5 * helper_10 + helper_20 + 0.5 * helper_7) -

                   (helper_1 * (-1.5 * helper_1 + 0.5 * helper_2 + helper_21) +

                   helper_13 * (-1.5 * helper_13 + 0.5 * helper_3 + 0.5 * helper_4 + 0.5 * helper_5) +

                   helper_15 * (0.5 * helper_1 - 1.5 * helper_15 + 0.5 * helper_16 + 0.5 * helper_2) +

                   helper_16 * (0.5 * helper_1 + 0.5 * helper_15 - 1.5 * helper_16 + 0.5 * helper_2) +

                   helper_2 * (0.5 * helper_1 - 1.5 * helper_2 + helper_21) +

                   helper_3 * (helper_19 - 1.5 * helper_3 + 0.5 * helper_5) +

                   helper_4 * (0.5 * helper_13 + 0.5 * helper_3 - 1.5 * helper_4 + 0.5 * helper_5) +

                   helper_5 * (helper_19 + 0.5 * helper_3 - 1.5 * helper_5) +

                   helper_7 * (0.5 * helper_10 + helper_20 - 1.5 * helper_7) +

                   helper_8 * (0.5 * helper_10 + 0.5 * helper_7 - 1.5 * helper_8 + 0.5 * helper_9) +

                   helper_9 * (0.5 * helper_10 + 0.5 * helper_7 + 0.5 * helper_8 - 1.5 * helper_9))

            / cbrt(helper_22*helper_22);

//                 * pow(pow((helper_1 - helper_2) * (helper_11 * helper_6 - helper_12 * helper_14) -

//                         (-helper_10 + helper_7) * (-helper_14 * helper_18 + helper_17 * helper_6) +

//                         (helper_3 - helper_5) * (-helper_11 * helper_18 + helper_12 * helper_17), 2),

//                     -0.333333333333333);

    return res.to_double();

}



Scalar floatTetWild::AMIPS_energy_aux(const std::array<Scalar, 12>& T) {

    Scalar helper_0[12];

    helper_0[0] = T[0];

src/LocalOperations.h

+2 −1
Original line number Diff line number Diff line
@@ -134,6 +134,7 @@ namespace floatTetWild { 
    };



    Scalar AMIPS_energy_aux(const std::array<Scalar, 12>& T);

    Scalar AMIPS_energy_mp(const std::array<Scalar, 12>& T);

    Scalar AMIPS_energy(const std::array<Scalar, 12>& T);

    void AMIPS_jacobian(const std::array<Scalar, 12>& T, Vector3& result_0);

    void AMIPS_hessian(const std::array<Scalar, 12>& T, Matrix3& result_0);