OilfieldToolsNet/LBPM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge branch 'master' of github.com:JamesEMcClure/LBPM-WIA

Browse Source
This commit is contained in:
JamesEMcclure 2021-05-17 19:10:37 -04:00
commit 2f12032973
16 changed files with 728 additions and 510 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
CMakeLists.txt
View File

@ -1,5 +1,8 @@
# Set some CMake properties
CMAKE_MINIMUM_REQUIRED( VERSION 3.9 )
if( ${CMAKE_VERSION} VERSION_GREATER_EQUAL "3.20.0")
CMAKE_POLICY( SET CMP0115 OLD )
endif()
MESSAGE("====================")

2
IO/IOHelpers.h
View File

@ -72,6 +72,6 @@ inline std::vector<std::string> splitList( const char *line, const char token )
}
}; // namespace IO
} // namespace IO
#endif

7
IO/silo.h
View File

@ -33,7 +33,8 @@ typedef int DBfile;
#endif
namespace IO::silo {
namespace IO {
namespace silo {
enum FileMode { READ, WRITE, CREATE };
@ -271,7 +272,9 @@ void writeMultiVar( DBfile *fid, const std::string &varname,
const std::vector<std::string> &subVarNames, const std::vector<int> &subVarTypes );
}; // namespace IO::silo
} // namespace silo
} // namespace IO
#endif
#include "IO/silo.hpp"

6
IO/silo.hpp
View File

@ -43,7 +43,8 @@
#include <silo.h>
namespace IO::silo {
namespace IO {
namespace silo {
/****************************************************
@ -443,7 +444,8 @@ Array<TYPE> readTriMeshVariable( DBfile *fid, const std::string &varname )
}
}; // namespace IO::silo
} // namespace silo
} // namespace IO
#endif

338
cmake/macros.cmake
View File

@ -1,3 +1,9 @@
CMAKE_MINIMUM_REQUIRED(VERSION 3.10.0)
CMAKE_POLICY( SET CMP0057 NEW )
if( ${CMAKE_VERSION} VERSION_GREATER_EQUAL "3.20.0")
CMAKE_POLICY( SET CMP0115 OLD )
endif()
INCLUDE(CheckCCompilerFlag)
INCLUDE(CheckCSourceCompiles)
INCLUDE(CheckCXXCompilerFlag)
@ -47,12 +53,24 @@ ELSEIF( NOT DEFINED CMAKE_INTERPROCEDURAL_OPTIMIZATION )
ENDIF()
# Check for gold linker
IF ( UNIX AND NOT APPLE AND NOT DISABLE_GOLD )
EXECUTE_PROCESS( COMMAND ${CMAKE_CXX_COMPILER} -fuse-ld=gold -Wl,--version ERROR_QUIET OUTPUT_VARIABLE LD_VERSION )
IF ( "${LD_VERSION}" MATCHES "GNU gold" )
MESSAGE( "Using gold linker" )
SET( GOLD_LINK_FLAGS "-fuse-ld=gold -Wl,--disable-new-dtags -Wl,--no-map-whole-files" )
SET( CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} ${GOLD_LINK_FLAGS}" )
SET( CMAKE_SHARED_LINKER_FLAGS "${CMAKE_SHARED_LINKER_FLAGS} ${GOLD_LINK_FLAGS}" )
ENDIF()
ENDIF()
# Add some default targets if they do not exist
IF ( NOT TARGET copy-${PROJ}-Data )
ADD_CUSTOM_TARGET( copy-${PROJ}-Data ALL )
ENDIF()
IF ( NOT TARGET copy-${PROJ}-include )
ADD_CUSTOM_TARGET ( copy-${PROJ}-include ALL )
ADD_CUSTOM_TARGET( copy-${PROJ}-include ALL )
ENDIF()
@ -75,7 +93,7 @@ ENDMACRO()
MACRO( PRINT_ALL_VARIABLES )
GET_CMAKE_PROPERTY(_variableNames VARIABLES)
FOREACH(_variableName ${_variableNames})
message(STATUS "${_variableName}=${${_variableName}}")
MESSAGE(STATUS "${_variableName}=${${_variableName}}")
ENDFOREACH()
ENDMACRO()
@ -108,28 +126,28 @@ MACRO( CONVERT_M4_FORTRAN IN LOCAL_PATH OUT_PATH )
IF ( NOT EXISTS "${CMAKE_CURRENT_BINARY_DIR}/${OUT_PATH}" )
FILE(MAKE_DIRECTORY "${CMAKE_CURRENT_BINARY_DIR}/${OUT_PATH}" )
ENDIF()
CONFIGURE_FILE ( ${IN} ${IN} COPYONLY )
CONFIGURE_FILE( ${IN} ${IN} COPYONLY )
IF ("${CMAKE_GENERATOR}" STREQUAL "Xcode")
STRING(REGEX REPLACE ".F" ".o" OUT2 "${OUT}" )
STRING(REGEX REPLACE ";" " " COMPILE_CMD "${CMAKE_Fortran_COMPILER} -c ${OUT} ${CMAKE_Fortran_FLAGS} -o ${OUT2}")
STRING(REGEX REPLACE "\\\\" "" COMPILE_CMD "${COMPILE_CMD}")
MESSAGE("COMPILE_CMD = ${COMPILE_CMD}")
SET( COMPILE_CMD ${COMPILE_CMD} )
add_custom_command(
ADD_CUSTOM_COMMAND(
OUTPUT ${OUT2}
COMMAND m4 -I${LOCAL_PATH} -I${SAMRAI_FORTDIR} ${M4DIRS} ${IN} > ${OUT}
COMMAND ${COMPILE_CMD}
DEPENDS ${IN}
)
set_source_files_properties(${OUT2} PROPERTIES GENERATED true)
SET_SOURCE_FILES_PROPERTIES(${OUT2} PROPERTIES GENERATED true)
SET( SOURCES ${SOURCES} "${OUT2}" )
ELSE()
add_custom_command(
ADD_CUSTOM_COMMAND(
OUTPUT ${OUT}
COMMAND m4 -I${LOCAL_PATH} -I${SAMRAI_FORTDIR} ${M4DIRS} ${M4_OPTIONS} ${IN} > ${OUT}
DEPENDS ${IN}
)
set_source_files_properties(${OUT} PROPERTIES GENERATED true)
SET_SOURCE_FILES_PROPERTIES(${OUT} PROPERTIES GENERATED true)
SET( SOURCES ${SOURCES} "${OUT}" )
ENDIF()
ENDMACRO()
@ -145,6 +163,10 @@ ENDMACRO()
# Add a project executable
MACRO( ADD_${PROJ}_EXECUTABLE EXEFILE )
IF ( TARGET ${EXEFILE} AND NOT INCLUDE_${EXEFILE} )
MESSAGE( FATAL_ERROR "ADD_${PROJ}_EXECUTABLE should be called before adding tests" )
ENDIF()
SET( INCLUDE_${EXEFILE} TRUE )
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
INSTALL( TARGETS ${EXEFILE} DESTINATION ${${PROJ}_INSTALL_DIR}/bin )
ENDMACRO()
@ -164,7 +186,7 @@ ENDMACRO ()
# Find the source files
MACRO (FIND_FILES)
MACRO( FIND_FILES )
# Find the C/C++ headers
SET( T_HEADERS "" )
FILE( GLOB T_HEADERS "*.h" "*.H" "*.hh" "*.hpp" "*.I" )
@ -193,7 +215,7 @@ MACRO (FIND_FILES)
SET( CSOURCES ${CSOURCES} ${T_CSOURCES} )
SET( FSOURCES ${FSOURCES} ${T_FSOURCES} )
SET( M4FSOURCES ${M4FSOURCES} ${T_M4FSOURCES} )
SET( SOURCES ${SOURCES} ${T_CXXSOURCES} ${T_CSOURCES} ${T_FSOURCES} ${T_M4FSOURCES} )
SET( SOURCES ${SOURCES} ${T_CXXSOURCES} ${T_CSOURCES} ${T_FSOURCES} ${T_M4FSOURCES} ${CUDASOURCES} )
ENDMACRO()
@ -303,7 +325,7 @@ MACRO( INSTALL_${PROJ}_TARGET PACKAGE )
IF ( TARGET write_repo_version )
ADD_DEPENDENCIES( ${PACKAGE} write_repo_version )
ENDIF()
ADD_DEPENDENCIES ( ${PACKAGE} copy-${PROJ}-include )
ADD_DEPENDENCIES( ${PACKAGE} copy-${PROJ}-include )
IF ( NOT ${PROJ}_LIB )
INSTALL( TARGETS ${PACKAGE} DESTINATION ${${PROJ}_INSTALL_DIR}/lib )
ENDIF()
@ -340,7 +362,7 @@ ENDMACRO()
# Macro to verify that a path has been set
MACRO( VERIFY_PATH PATH_NAME )
IF ("${PATH_NAME}" STREQUAL "")
MESSAGE ( FATAL_ERROR "Path is not set: ${PATH_NAME}" )
MESSAGE( FATAL_ERROR "Path is not set: ${PATH_NAME}" )
ENDIF()
IF ( NOT EXISTS "${PATH_NAME}" )
MESSAGE( FATAL_ERROR "Path does not exist: ${PATH_NAME}" )
@ -420,6 +442,9 @@ MACRO( IDENTIFY_COMPILER )
ELSEIF ( ${CMAKE_Fortran_COMPILER_ID} MATCHES "PGI")
SET(USING_PGF90 TRUE)
MESSAGE("Using pgf90")
ELSEIF( (${CMAKE_Fortran_COMPILER_ID} MATCHES "CRAY") OR (${CMAKE_Fortran_COMPILER_ID} MATCHES "Cray") )
SET(USING_CRAY TRUE)
MESSAGE("Using Cray")
ELSEIF ( (${CMAKE_Fortran_COMPILER_ID} MATCHES "CLANG") OR (${CMAKE_Fortran_COMPILER_ID} MATCHES "Clang") OR
(${CMAKE_Fortran_COMPILER_ID} MATCHES "FLANG") OR (${CMAKE_Fortran_COMPILER_ID} MATCHES "Flang") )
SET(USING_FLANG TRUE)
@ -440,7 +465,7 @@ MACRO( SET_WARNINGS )
# Note: adding -Wlogical-op causes a wierd linking error on Titan using the nvcc wrapper:
# /usr/bin/ld: cannot find gical-op: No such file or directory
SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -Wall -Wextra")
SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -Woverloaded-virtual -Wsign-compare")
SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -pedantic -Woverloaded-virtual -Wsign-compare -Wformat-security")
ELSEIF ( USING_MSVC )
# Add Microsoft specifc compiler options
SET(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} /D _SCL_SECURE_NO_WARNINGS /D _CRT_SECURE_NO_WARNINGS /D _ITERATOR_DEBUG_LEVEL=0 /wd4267" )
@ -448,7 +473,7 @@ MACRO( SET_WARNINGS )
ELSEIF ( USING_ICC )
# Add Intel specifc compiler options
SET(CMAKE_C_FLAGS " ${CMAKE_C_FLAGS} -Wall" )
SET(CMAKE_CXX_FLAGS " ${CMAKE_CXX_FLAGS} -Wall" )
SET(CMAKE_CXX_FLAGS " ${CMAKE_CXX_FLAGS} -Wall -wd1011" )
ELSEIF ( USING_CRAY )
# Add default compiler options
SET(CMAKE_C_FLAGS " ${CMAKE_C_FLAGS}")
@ -461,8 +486,8 @@ MACRO( SET_WARNINGS )
SET(CMAKE_CXX_FLAGS " ${CMAKE_CXX_FLAGS} --diag_suppress 111,128,185")
ELSEIF ( USING_CLANG )
# Add default compiler options
SET(CMAKE_C_FLAGS " ${CMAKE_C_FLAGS} -Wall")
SET(CMAKE_CXX_FLAGS " ${CMAKE_CXX_FLAGS} -Wall -Wno-missing-braces -Wmissing-field-initializers -ftemplate-depth=1024")
SET(CMAKE_C_FLAGS " ${CMAKE_C_FLAGS} -Wall -Wextra")
SET(CMAKE_CXX_FLAGS " ${CMAKE_CXX_FLAGS} -Wall -Wextra -Wpedantic -Wno-missing-braces -Wmissing-field-initializers -ftemplate-depth=1024")
ELSEIF ( USING_XL )
# Add default compiler options
SET(CMAKE_C_FLAGS " ${CMAKE_C_FLAGS} -Wall")
@ -651,7 +676,7 @@ MACRO( COPY_MESH_FILE MESHNAME )
ENDFOREACH()
# We have either found the mesh or failed
IF ( NOT MESHPATH )
MESSAGE ( WARNING "Cannot find mesh: " ${MESHNAME} )
MESSAGE( WARNING "Cannot find mesh: " ${MESHNAME} )
ELSE ()
SET( MESHPATH2 )
FOREACH( tmp ${MESHPATH} )
@ -666,25 +691,16 @@ ENDMACRO()
# Link the libraries to the given target
MACRO( TARGET_LINK_EXTERNAL_LIBRARIES TARGET_NAME )
FOREACH ( tmp ${TPL_LIBS} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
# Include external libraries
FOREACH ( tmp ${EXTERNAL_LIBS} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
FOREACH ( tmp ${LAPACK_LIBS} )
# Include libraries found through the TPL builder
FOREACH ( tmp ${TPL_LIBS} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
FOREACH ( tmp ${BLAS_LIBS} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
FOREACH ( tmp ${BLAS_LAPACK_LIBS} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
FOREACH ( tmp ${MPI_LIBRARIES} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
FOREACH ( tmp ${CMAKE_C_IMPLICIT_LINK_LIBRARIES}
# Include CMake implicit libraries
FOREACH ( tmp ${CMAKE_CUDA_IMPLICIT_LINK_LIBRARIES} ${CMAKE_C_IMPLICIT_LINK_LIBRARIES}
${CMAKE_CXX_IMPLICIT_LINK_LIBRARIES} ${CMAKE_Fortran_IMPLICIT_LINK_LIBRARIES} )
TARGET_LINK_LIBRARIES( ${TARGET_NAME} ${ARGN} ${tmp} )
ENDFOREACH()
@ -715,8 +731,8 @@ ENDFUNCTION()
MACRO( ADD_PROJ_EXE_DEP EXE )
# Add the package dependencies
IF( ${PROJ}_TEST_LIB_EXISTS )
ADD_DEPENDENCIES ( ${EXE} ${PACKAGE_TEST_LIB} )
TARGET_LINK_LIBRARIES ( ${EXE} ${PACKAGE_TEST_LIB} )
ADD_DEPENDENCIES( ${EXE} ${PACKAGE_TEST_LIB} )
TARGET_LINK_LIBRARIES( ${EXE} ${PACKAGE_TEST_LIB} )
ENDIF()
# Add the executable to the dependencies of check and build-test
ADD_DEPENDENCIES( check ${EXE} )
@ -783,7 +799,7 @@ FUNCTION( ADD_PROJ_PROVISIONAL_TEST EXEFILE )
SET( CXXFILE ${EXEFILE} )
SET( TESTS_SO_FAR ${TESTS_SO_FAR} ${EXEFILE} )
# Check if we want to add the test to all
IF ( NOT EXCLUDE_TESTS_FROM_ALL )
IF ( NOT EXCLUDE_TESTS_FROM_ALL OR INCLUDE_${EXEFILE} )
ADD_EXECUTABLE( ${EXEFILE} ${CXXFILE} )
ELSE()
ADD_EXECUTABLE( ${EXEFILE} EXCLUDE_FROM_ALL ${CXXFILE} )
@ -818,123 +834,218 @@ MACRO( CREATE_TEST_NAME TEST ${ARGN} )
FOREACH( tmp ${ARGN})
SET( TESTNAME "${TESTNAME}--${tmp}")
endforeach()
# STRING(REGEX REPLACE "--" "-" TESTNAME ${TESTNAME} )
SET( LAST_TESTNAME ${TESTNAME} PARENT_SCOPE )
ENDMACRO()
# Function to add the resource locks to an executable
FUNCTION( ADD_RESOURCE_LOCK TESTNAME EXEFILE ${ARGN} )
FUNCTION( ADD_RESOURCE TESTNAME RESOURCE )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES RESOURCE_LOCK ${RESOURCE} )
ENDFUNCTION()
FUNCTION( ADD_RESOURCE_LOCKS TESTNAME EXEFILE ${ARGN} )
IF ( NOT DISABLE_RESOURCE_LOCK )
IF ( NOT ARGN )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES RESOURCE_LOCK ${EXEFILE} )
ADD_RESOURCE( ${TESTNAME} ${EXEFILE} )
ELSE()
FOREACH( tmp ${ARGN} )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES RESOURCE_LOCK ${tmp} )
ADD_RESOURCE( ${TESTNAME} ${tmp} )
ENDFOREACH()
ENDIF()
ENDIF()
ENDFUNCTION()
# Add a executable as a test
FUNCTION( ADD_${PROJ}_TEST EXEFILE ${ARGN} )
# Function to add a test
FUNCTION( CALL_ADD_TEST EXEFILE PROCS THREADS WEEKLY TESTBUILDER ${ARGN} )
# Check if we actually want to add the test
KEEP_TEST( RESULT )
IF ( NOT RESULT )
MESSAGE( "Discarding test: ${EXEFILE}" )
RETURN()
ENDIF()
# Add the provisional test
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
CREATE_TEST_NAME( ${EXEFILE} ${ARGN} )
IF ( USE_MPI_FOR_SERIAL_TESTS )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIFLAGS} "${MPIEXEC_NUMPROC_FLAG}" 1 $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
# Check if we are dealing with a TestBuilder test
IF ( ${TESTBUILDER} )
# Add the test to the TestBuilder sources
SET( TESTBUILDER_SOURCES ${TESTBUILDER_SOURCES} ${EXEFILE} PARENT_SCOPE )
ELSE()
ADD_TEST( NAME ${TESTNAME} COMMAND $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
# Add the provisional test
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
ENDIF()
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS 1 )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
ENDFUNCTION()
# Add a executable as a weekly test
FUNCTION( ADD_${PROJ}_WEEKLY_TEST EXEFILE PROCS ${ARGN} )
# Check if we actually want to add the test
KEEP_TEST( RESULT )
IF ( NOT RESULT )
RETURN()
# Create the test name
SET( TESTNAME ${EXEFILE} )
#IF ( ${TESTBUILDER} )
# SET( TESTNAME TB_${TESTNAME} )
#ENDIF()
IF ( "${PROCS}" GREATER "1")
SET( TESTNAME ${TESTNAME}_${PROCS}procs )
ENDIF()
# Add the provisional test
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
IF( ${PROCS} STREQUAL "1" )
CREATE_TEST_NAME( "${EXEFILE}_WEEKLY" ${ARGN} )
ELSEIF( (USE_MPI OR USE_EXT_MPI) AND NOT (${PROCS} GREATER ${TEST_MAX_PROCS}) )
CREATE_TEST_NAME( "${EXEFILE}_${PROCS}procs_WEEKLY" ${ARGN} )
IF ( "${THREADS}" GREATER "1")
SET( TESTNAME ${TESTNAME}_${THREADS}threads )
ENDIF()
IF ( ${PROCS} GREATER ${TEST_MAX_PROCS} )
IF ( ${WEEKLY} )
SET( TESTNAME ${TESTNAME}_WEEKLY )
ENDIF()
CREATE_TEST_NAME( ${TESTNAME} ${ARGN} )
# Create the tests for ctest
MATH( EXPR TOT_PROCS "${PROCS} * ${THREADS}" )
IF ( ${TOT_PROCS} EQUAL 0 )
# Skip test (provisional)
ELSEIF ( ${TOT_PROCS} GREATER ${TEST_MAX_PROCS} )
MESSAGE("Disabling test ${TESTNAME} (exceeds maximum number of processors ${TEST_MAX_PROCS})")
ELSEIF( ${PROCS} STREQUAL "1" )
CREATE_TEST_NAME( "${EXEFILE}_WEEKLY" ${ARGN} )
IF ( USE_MPI_FOR_SERIAL_TESTS )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIFLAGS} "${MPIEXEC_NUMPROC_FLAG}" 1 $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
ELSEIF ( ( ${PROCS} STREQUAL "1" ) AND NOT USE_MPI_FOR_SERIAL_TESTS )
IF ( ${TESTBUILDER} )
ADD_TEST( NAME ${TESTNAME} COMMAND $<TARGET_FILE:${TB_TARGET}> ${EXEFILE} ${ARGN} )
ELSE()
ADD_TEST( NAME ${TESTNAME} COMMAND $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
ENDIF()
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS 1 )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
ELSEIF( (USE_MPI OR USE_EXT_MPI) AND NOT (${PROCS} GREATER ${TEST_MAX_PROCS}) )
CREATE_TEST_NAME( "${EXEFILE}_${PROCS}procs_WEEKLY" ${ARGN} )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} "${MPIEXEC_NUMPROC_FLAG}" ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${PROCS} )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${TOT_PROCS} )
ADD_RESOURCE_LOCKS( ${TESTNAME} ${EXEFILE} ${ARGN} )
ELSEIF ( USE_MPI OR USE_EXT_MPI )
IF ( ${TESTBUILDER} )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} ${PROCS} $<TARGET_FILE:${TB_TARGET}> ${EXEFILE} ${ARGN} )
ELSE()
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
ENDIF()
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${TOT_PROCS} )
ADD_RESOURCE_LOCKS( ${TESTNAME} ${EXEFILE} ${ARGN} )
ENDIF()
ENDFUNCTION()
# Add a executable as a parallel test
FUNCTION( ADD_${PROJ}_TEST_PARALLEL EXEFILE PROCS ${ARGN} )
# Functions to create a test builder
FUNCTION( INITIALIZE_TESTBUILDER )
SET( TESTBUILDER_SOURCES PARENT_SCOPE )
# Check if we actually want to add the test
KEEP_TEST( RESULT )
IF ( NOT RESULT )
RETURN()
ENDIF()
# Add the provisional test
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
CREATE_TEST_NAME( "${EXEFILE}_${PROCS}procs" ${ARGN} )
IF ( NOT ( USE_MPI OR USE_EXT_MPI ) )
MESSAGE("Disabling test ${TESTNAME} (configured without MPI)")
ELSEIF ( ${PROCS} GREATER ${TEST_MAX_PROCS} )
MESSAGE("Disabling test ${TESTNAME} (exceeds maximum number of processors ${TEST_MAX_PROCS})")
ELSE()
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIFLAGS} "${MPIEXEC_NUMPROC_FLAG}" ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${PROCS} )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
# Check if the target does not exist (may not be added yet or we are re-configuring)
STRING(REGEX REPLACE "${${PROJ}_SOURCE_DIR}/" "" TB_TARGET "TB-${PROJ}-${CMAKE_CURRENT_SOURCE_DIR}" )
STRING(REGEX REPLACE "/" "-" TB_TARGET ${TB_TARGET} )
IF ( NOT TARGET ${TB_TARGET} )
GLOBAL_SET( ${TB_TARGET}-BINDIR )
ENDIF()
# Check if test has already been added
IF ( NOT ${TB_TARGET}-BINDIR )
GLOBAL_SET( ${TB_TARGET}-BINDIR "${CMAKE_CURRENT_BINARY_DIR}" )
# Create the initial file
IF ( NOT EXISTS "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder.cpp" )
FILE( WRITE "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder.cpp" "// Auto generated file\n\n" )
ENDIF()
# The target has not been added
SET( CXXFILE "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder.cpp" )
SET( TESTS_SO_FAR ${TESTS_SO_FAR} ${TB_TARGET} )
# Check if we want to add the test to all
IF ( NOT EXCLUDE_TESTS_FROM_ALL )
ADD_EXECUTABLE( ${TB_TARGET} ${CXXFILE} )
ELSE()
ADD_EXECUTABLE( ${TB_TARGET} EXCLUDE_FROM_ALL ${CXXFILE} )
ENDIF()
SET_TARGET_PROPERTIES( ${TB_TARGET} PROPERTIES OUTPUT_NAME "TB" )
ELSEIF ( NOT ${TB-BINDIR} STREQUAL ${CMAKE_CURRENT_BINARY_DIR} )
# We are trying to add 2 different tests with the same name
MESSAGE( "Existing test: ${TB-BINDIR}/${TB_TARGET}" )
MESSAGE( "New test: ${CMAKE_CURRENT_BINARY_DIR}/${TB_TARGET}}" )
MESSAGE( FATAL_ERROR "Trying to add 2 different tests with the same name" )
ENDIF()
SET( TB_TARGET ${TB_TARGET} PARENT_SCOPE )
ENDFUNCTION()
FUNCTION( FINALIZE_TESTBUILDER )
# Check if we actually want to add the test
KEEP_TEST( RESULT )
IF ( NOT RESULT )
RETURN()
ENDIF()
# Create the library for the test builder
IF ( TESTBUILDER_SOURCES )
LIST( REMOVE_DUPLICATES TESTBUILDER_SOURCES )
SET( TB_TARGET_SOURCES "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder.cpp" )
FOREACH ( tmp ${TESTBUILDER_SOURCES} )
GET_FILENAME_COMPONENT( tmp ${tmp} ABSOLUTE )
SET( TB_TARGET_SOURCES ${TB_TARGET_SOURCES} ${tmp} )
ENDFOREACH()
TARGET_SOURCES( ${TB_TARGET} PUBLIC ${TB_TARGET_SOURCES} )
ADD_DEPENDENCIES( ${TB_TARGET} copy-${PROJ}-include )
ENDIF()
# Add the dependencies to the test builder
ADD_PROJ_EXE_DEP( ${TB_TARGET} )
IF ( CURRENT_LIBRARY )
IF ( NOT TARGET ${CURRENT_LIBRARY}-test )
ADD_CUSTOM_TARGET( ${CURRENT_LIBRARY}-test )
ENDIF()
ADD_DEPENDENCIES( ${CURRENT_LIBRARY}-test ${TB_TARGET} )
ENDIF()
SET( TB_TARGET ${TB_TARGET} PARENT_SCOPE )
# Finish generating TestBuilder.cpp
SET( TB_FILE "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder-tmp.cpp" )
FILE( WRITE "${TB_FILE}" "// Auto generated file\n\n" )
FILE( APPEND "${TB_FILE}" "#include <cstring>\n" )
FILE( APPEND "${TB_FILE}" "#include <iostream>\n\n" )
FOREACH( tmp ${TESTBUILDER_SOURCES} )
FILE( APPEND "${TB_FILE}" "extern int ${tmp}( int argc, char *argv[] );\n" )
ENDFOREACH()
FILE( APPEND "${TB_FILE}" "\n\n// Main Program\n" )
FILE( APPEND "${TB_FILE}" "int main( int argc, char *argv[] )\n" )
FILE( APPEND "${TB_FILE}" "{\n" )
FILE( APPEND "${TB_FILE}" " if ( argc < 2 ) {\n" )
FILE( APPEND "${TB_FILE}" " std::cerr << \"Invaild number of arguments for TestBuilder\";\n " )
FILE( APPEND "${TB_FILE}" " return 1;\n" )
FILE( APPEND "${TB_FILE}" " }\n\n" )
FILE( APPEND "${TB_FILE}" " int rtn = 0;\n" )
FILE( APPEND "${TB_FILE}" " if ( strcmp( argv[1], \"NULL\" ) == 0 ) {\n" )
FOREACH( tmp ${TESTBUILDER_SOURCES} )
FILE( APPEND "${TB_FILE}" " } else if ( strcmp( argv[1], \"${tmp}\" ) == 0 ) {\n" )
FILE( APPEND "${TB_FILE}" " rtn = ${tmp}( argc-1, &argv[1] );\n" )
ENDFOREACH()
FILE( APPEND "${TB_FILE}" " } else {\n" )
FILE( APPEND "${TB_FILE}" " std::cerr << \"Unknown test\";\n" )
FILE( APPEND "${TB_FILE}" " return 1;\n" )
FILE( APPEND "${TB_FILE}" " }\n\n" )
FILE( APPEND "${TB_FILE}" " return rtn;\n" )
FILE( APPEND "${TB_FILE}" "}\n" )
EXECUTE_PROCESS( COMMAND ${CMAKE_COMMAND} -E copy_if_different "${TB_FILE}" "${CMAKE_CURRENT_BINARY_DIR}/TestBuilder.cpp" )
ENDFUNCTION()
# Add a parallel test that may use both MPI and threads
# This allows us to correctly compute the number of processors used by the test
# Convience functions to add a test
MACRO( ADD_${PROJ}_TEST EXEFILE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 1 1 FALSE FALSE ${ARGN} )
ENDMACRO()
MACRO( ADD_${PROJ}_TEST_PARALLEL EXEFILE PROCS ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} ${PROCS} 1 FALSE FALSE ${ARGN} )
ENDMACRO()
MACRO( ADD_${PROJ}_WEEKLY_TEST EXEFILE PROCS ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} ${PROCS} 1 TRUE FALSE ${ARGN} )
ENDMACRO()
MACRO( ADD_${PROJ}_TEST_THREAD_MPI EXEFILE PROCS THREADS ${ARGN} )
ADD_PROJ_PROVISIONAL_TEST( ${EXEFILE} )
CREATE_TEST_NAME( "${EXEFILE}_${PROCS}procs_${THREADS}threads" ${ARGN} )
MATH( EXPR TOT_PROCS "${PROCS} * ${THREADS}" )
IF ( ${TOT_PROCS} GREATER ${TEST_MAX_PROCS} )
MESSAGE("Disabling test ${TESTNAME} (exceeds maximum number of processors ${TEST_MAX_PROCS})")
ELSEIF ( ( ${PROCS} STREQUAL "1" ) AND NOT USE_MPI_FOR_SERIAL_TESTS )
ADD_TEST( NAME ${TESTNAME} COMMAND $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_TESTS_PROPERTIES ( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${TOT_PROCS} )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
ELSEIF ( USE_MPI OR USE_EXT_MPI )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIFLAGS} "${MPIEXEC_NUMPROC_FLAG}" ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
SET_TESTS_PROPERTIES ( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${TOT_PROCS} )
ADD_RESOURCE_LOCK( ${TESTNAME} ${EXEFILE} ${ARGN} )
ENDIF()
CALL_ADD_TEST( ${EXEFILE} ${PROCS} ${THREADS} FALSE FALSE ${ARGN} )
ENDMACRO()
MACRO( ADD_${PROJ}_TEST_1_2_4 EXEFILE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 1 1 FALSE FALSE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 2 1 FALSE FALSE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 4 1 FALSE FALSE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_PROVISIONAL_TEST EXEFILE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 0 1 FALSE TRUE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_TEST EXEFILE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 1 1 FALSE TRUE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_TEST_PARALLEL EXEFILE PROCS ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} ${PROCS} 1 FALSE TRUE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_WEEKLY_TEST EXEFILE PROCS ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} ${PROCS} 1 TRUE TRUE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_TEST_THREAD_MPI EXEFILE PROCS THREADS ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} ${PROCS} ${THREADS} FALSE TRUE ${ARGN} )
ENDMACRO()
MACRO( ADD_TB_TEST_1_2_4 EXEFILE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 1 1 FALSE TRUE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 2 1 FALSE TRUE ${ARGN} )
CALL_ADD_TEST( ${EXEFILE} 4 1 FALSE TRUE ${ARGN} )
ENDMACRO()
@ -966,7 +1077,7 @@ FUNCTION( ADD_${PROJ}_EXAMPLE EXEFILE PROCS ${ARGN} )
ADD_TEST( NAME ${TESTNAME} COMMAND $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
ELSEIF ( USE_EXT_MPI AND NOT (${PROCS} GREATER ${TEST_MAX_PROCS}) )
CREATE_TEST_NAME( "example--${EXEFILE}_${PROCS}procs" ${ARGN} )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIFLAGS} "${MPIEXEC_NUMPROC_FLAG}" ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
ADD_TEST( NAME ${TESTNAME} COMMAND ${MPIEXEC} ${MPIEXEC_NUMPROC_FLAG} ${PROCS} $<TARGET_FILE:${LAST_TEST}> ${ARGN} )
SET_PROPERTY( TEST ${TESTNAME} APPEND PROPERTY ENVIRONMENT OMPI_MCA_hwloc_base_binding_policy=none )
ENDIF()
SET_TESTS_PROPERTIES( ${TESTNAME} PROPERTIES FAIL_REGULAR_EXPRESSION "${TEST_FAIL_REGULAR_EXPRESSION}" PROCESSORS ${PROCS} )
@ -1038,6 +1149,8 @@ FUNCTION( ADD_MATLAB_MEX SOURCE )
MATLAB_ADD_MEX(
NAME ${TARGET}
SRC ${SOURCE}
MODULE
R2017b
)
TARGET_LINK_LIBRARIES( ${TARGET} ${MATLAB_TARGET} )
ADD_PROJ_EXE_DEP( ${TARGET} )
@ -1074,9 +1187,11 @@ ENDMACRO()
# Create a script to start matlab preloading libraries
FUNCTION( CREATE_MATLAB_WRAPPER )
SET( MATLABPATH ${MATLABPATH} "${${PROJ}_INSTALL_DIR}/mex" )
SET( tmp_libs ${MEX_LIBCXX} ${MEX_FILES} )
STRING(REGEX REPLACE ";" ":" tmp_libs "${tmp_libs}")
STRING(REGEX REPLACE ";" ":" tmp_path "${MATLABPATH}")
LIST( REMOVE_DUPLICATES MATLABPATH )
STRING( REGEX REPLACE ";" ":" tmp_libs "${tmp_libs}" )
STRING( REGEX REPLACE ";" ":" tmp_path "${MATLABPATH}" )
IF ( USING_MSVC )
# Create a matlab wrapper for windows
SET( MATLAB_GUI "${CMAKE_CURRENT_BINARY_DIR}/tmp/matlab-gui.bat" )
@ -1185,7 +1300,7 @@ MACRO( ADD_DISTCLEAN ${ARGN} )
compile_commands.json
${ARGN}
)
ADD_CUSTOM_TARGET (distclean @echo cleaning for source distribution)
ADD_CUSTOM_TARGET(distclean @echo cleaning for source distribution)
IF (UNIX)
ADD_CUSTOM_COMMAND(
DEPENDS clean
@ -1269,3 +1384,4 @@ FUNCTION( FIND_PYTHON_MODULE MODULE)
ENDIF()
SET( PY_${MODULE2}_FOUND ${PY_${MODULE2}_FOUND} PARENT_SCOPE )
ENDFUNCTION(FIND_PYTHON_MODULE)

10
common/MPI.I
View File

@ -595,7 +595,15 @@ template<class TYPE>
void MPI_CLASS::sendrecv( const TYPE *sendbuf, int sendcount, int dest, int sendtag,
TYPE *recvbuf, int recvcount, int source, int recvtag ) const
{
ERROR( "Not implimented" );
if ( getSize() == 1 ) {
ASSERT( dest == 0 );
ASSERT( source == 0 );
ASSERT( sendcount == recvcount );
ASSERT( sendtag == recvtag );
memcpy( recvbuf, sendbuf, sendcount * sizeof( TYPE ) );
} else {
ERROR( "Not implimented for " + std::string( typeid( TYPE ).name() ) );
}
}
#else
template<class TYPE>

8
common/ScaLBL.h
View File

@ -104,14 +104,14 @@ extern "C" void ScaLBL_D3Q19_AAodd_GreyscaleColor(int *d_neighborList, int *Map,
double Fx, double Fy, double Fz, int strideY, int strideZ, int start, int finish, int Np);
extern "C" void ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi,double *Psi, double *GreySolidGrad, double *Poros,double *Perm,double *Vel, double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm,double *Vel, double *Pressure,
double rhoA, double rhoB, double tauA, double tauB,double tauA_eff,double tauB_eff, double alpha, double beta,
double Fx, double Fy, double Fz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np);
double Fx, double Fy, double Fz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np);
extern "C" void ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *d_neighborList, int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi,double *Psi, double *GreySolidGrad, double *Poros,double *Perm,double *Vel,double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm,double *Vel,double *Pressure,
double rhoA, double rhoB, double tauA, double tauB, double tauA_eff,double tauB_eff, double alpha, double beta,
double Fx, double Fy, double Fz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np);
double Fx, double Fy, double Fz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np);
//extern "C" void ScaLBL_Update_GreyscalePotential(int *Map, double *Phi, double *Psi, double *Poro, double *Perm, double alpha, double W,
// int start, int finish, int Np);

168
cuda/GreyscaleColor.cu
View File

@ -1450,9 +1450,9 @@ __global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor(int *Map, double *dist,
//CP: capillary penalty
// also turn off recoloring for grey nodes
__global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi, double *Psi, double *GreySolidGrad, double *Poros,double *Perm, double *Velocity, double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm, double *Velocity, double *Pressure,
double rhoA, double rhoB, double tauA, double tauB,double tauA_eff,double tauB_eff,double alpha, double beta,
double Gx, double Gy, double Gz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np){
double Gx, double Gy, double Gz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np){
int n,nn,ijk,nread;
int nr1,nr2,nr3,nr4,nr5,nr6;
@ -1462,8 +1462,6 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
double fq;
// conserved momemnts
double rho,jx,jy,jz;
//double vx,vy,vz,v_mag;
//double ux,uy,uz,u_mag;
double ux,uy,uz;
// non-conserved moments
double m1,m2,m4,m6,m8,m9,m10,m11,m12,m13,m14,m15,m16,m17,m18;
@ -1473,18 +1471,13 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
double C,nx,ny,nz; //color gradient magnitude and direction
double phi,tau,rho0,rlx_setA,rlx_setB;
//double GeoFun=0.0;//geometric function from Guo's PRE 66, 036304 (2002)
double porosity;
double perm;//voxel permeability
//double c0, c1; //Guo's model parameters
double tau_eff;
double mu_eff;//kinematic viscosity
double nx_gs,ny_gs,nz_gs;//grey-solid color gradient
double nx_phase,ny_phase,nz_phase,C_phase;
double Fx,Fy,Fz;
double gp1,gp2,gp4,gp6,gp8,gp9,gp10,gp11,gp12,gp13,gp14,gp15,gp16,gp17,gp18;
double gp3,gp5,gp7;
double Fcpx,Fcpy,Fcpz;//capillary penalty force
double W;//greyscale wetting strength
const double mrt_V1=0.05263157894736842;
const double mrt_V2=0.012531328320802;
@ -1510,9 +1503,7 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
porosity = Poros[n];
perm = Perm[n];
nx_gs = GreySolidGrad[n+0*Np];
ny_gs = GreySolidGrad[n+1*Np];
nz_gs = GreySolidGrad[n+2*Np];
W = GreySolidW[n];
// compute phase indicator field
phi=(nA-nB)/(nA+nB);
@ -1534,98 +1525,61 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
//........................................................................
nn = ijk-1; // neighbor index (get convention)
m1 = Phi[nn]; // get neighbor for phi - 1
gp1 = Psi[nn];
//........................................................................
nn = ijk+1; // neighbor index (get convention)
m2 = Phi[nn]; // get neighbor for phi - 2
gp2 = Psi[nn];
//........................................................................
nn = ijk-strideY; // neighbor index (get convention)
m3 = Phi[nn]; // get neighbor for phi - 3
gp3 = Psi[nn];
//........................................................................
nn = ijk+strideY; // neighbor index (get convention)
m4 = Phi[nn]; // get neighbor for phi - 4
gp4 = Psi[nn];
//........................................................................
nn = ijk-strideZ; // neighbor index (get convention)
m5 = Phi[nn]; // get neighbor for phi - 5
gp5 = Psi[nn];
//........................................................................
nn = ijk+strideZ; // neighbor index (get convention)
m6 = Phi[nn]; // get neighbor for phi - 6
gp6 = Psi[nn];
//........................................................................
nn = ijk-strideY-1; // neighbor index (get convention)
m7 = Phi[nn]; // get neighbor for phi - 7
gp7 = Psi[nn];
//........................................................................
nn = ijk+strideY+1; // neighbor index (get convention)
m8 = Phi[nn]; // get neighbor for phi - 8
gp8 = Psi[nn];
//........................................................................
nn = ijk+strideY-1; // neighbor index (get convention)
m9 = Phi[nn]; // get neighbor for phi - 9
gp9 = Psi[nn];
//........................................................................
nn = ijk-strideY+1; // neighbor index (get convention)
m10 = Phi[nn]; // get neighbor for phi - 10
gp10 = Psi[nn];
//........................................................................
nn = ijk-strideZ-1; // neighbor index (get convention)
m11 = Phi[nn]; // get neighbor for phi - 11
gp11 = Psi[nn];
//........................................................................
nn = ijk+strideZ+1; // neighbor index (get convention)
m12 = Phi[nn]; // get neighbor for phi - 12
gp12 = Psi[nn];
//........................................................................
nn = ijk+strideZ-1; // neighbor index (get convention)
m13 = Phi[nn]; // get neighbor for phi - 13
gp13 = Psi[nn];
//........................................................................
nn = ijk-strideZ+1; // neighbor index (get convention)
m14 = Phi[nn]; // get neighbor for phi - 14
gp14 = Psi[nn];
//........................................................................
nn = ijk-strideZ-strideY; // neighbor index (get convention)
m15 = Phi[nn]; // get neighbor for phi - 15
gp15 = Psi[nn];
//........................................................................
nn = ijk+strideZ+strideY; // neighbor index (get convention)
m16 = Phi[nn]; // get neighbor for phi - 16
gp16 = Psi[nn];
//........................................................................
nn = ijk+strideZ-strideY; // neighbor index (get convention)
m17 = Phi[nn]; // get neighbor for phi - 17
gp17 = Psi[nn];
//........................................................................
nn = ijk-strideZ+strideY; // neighbor index (get convention)
m18 = Phi[nn]; // get neighbor for phi - 18
gp18 = Psi[nn];
//............Compute the Color Gradient...................................
nx_phase = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
ny_phase = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
nz_phase = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
C_phase = sqrt(nx_phase*nx_phase+ny_phase*ny_phase+nz_phase*nz_phase);
//correct the normal color gradient by considering the effect of grey solid
nx = nx_phase + (1.0-porosity)*nx_gs;
ny = ny_phase + (1.0-porosity)*ny_gs;
nz = nz_phase + (1.0-porosity)*nz_gs;
if (C_phase==0.0){
nx = nx_phase;
ny = ny_phase;
nz = nz_phase;
}
//...........Normalize the Color Gradient.................................
C = sqrt(nx*nx+ny*ny+nz*nz);
double ColorMag = C;
if (C==0.0) ColorMag=1.0;
nx = nx/ColorMag;
ny = ny/ColorMag;
nz = nz/ColorMag;
nx = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
ny = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
nz = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
//............Compute the Greyscale Potential Gradient.....................
// Fcpx = 0.0;
@ -1648,14 +1602,24 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
// //ny = Fcpy/Fcp_mag;
// //nz = Fcpz/Fcp_mag;
// }
Fcpx = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
Fcpy = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
Fcpz = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
Fcpx = nx;
Fcpy = ny;
Fcpz = nz;
double Fcp_mag=sqrt(Fcpx*Fcpx+Fcpy*Fcpy+Fcpz*Fcpz);
if (Fcp_mag==0.0); Fcpx=Fcpy=Fcpz=0.0;
//NOTE for open node (porosity=1.0),Fcp=0.0
Fcpx *= alpha*W*(1.0-porosity)/sqrt(perm);
Fcpy *= alpha*W*(1.0-porosity)/sqrt(perm);
Fcpz *= alpha*W*(1.0-porosity)/sqrt(perm);
//...........Normalize the Color Gradient.................................
C = sqrt(nx*nx+ny*ny+nz*nz);
double ColorMag = C;
if (C==0.0) ColorMag=1.0;
nx = nx/ColorMag;
ny = ny/ColorMag;
nz = nz/ColorMag;
// q=0
fq = dist[n];
rho = fq;
@ -2262,15 +2226,13 @@ __global__ void dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *neighborList, int
//CP: capillary penalty
// also turn off recoloring for grey nodes
__global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi, double *Psi, double *GreySolidGrad, double *Poros,double *Perm, double *Velocity, double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm, double *Velocity, double *Pressure,
double rhoA, double rhoB, double tauA, double tauB,double tauA_eff,double tauB_eff, double alpha, double beta,
double Gx, double Gy, double Gz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np){
double Gx, double Gy, double Gz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np){
int ijk,nn,n;
double fq;
// conserved momemnts
double rho,jx,jy,jz;
//double vx,vy,vz,v_mag;
//double ux,uy,uz,u_mag;
double ux,uy,uz;
// non-conserved moments
double m1,m2,m4,m6,m8,m9,m10,m11,m12,m13,m14,m15,m16,m17,m18;
@ -2280,18 +2242,13 @@ __global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dis
double C,nx,ny,nz; //color gradient magnitude and direction
double phi,tau,rho0,rlx_setA,rlx_setB;
//double GeoFun=0.0;//geometric function from Guo's PRE 66, 036304 (2002)
double porosity;
double perm;//voxel permeability
//double c0, c1; //Guo's model parameters
double tau_eff;
double mu_eff;//kinematic viscosity
double nx_gs,ny_gs,nz_gs;//grey-solid color gradient
double nx_phase,ny_phase,nz_phase,C_phase;
double Fx,Fy,Fz;
double gp1,gp2,gp4,gp6,gp8,gp9,gp10,gp11,gp12,gp13,gp14,gp15,gp16,gp17,gp18;
double gp3,gp5,gp7;
double Fcpx,Fcpy,Fcpz;//capillary penalty force
double W;//greyscale wetting strength
const double mrt_V1=0.05263157894736842;
const double mrt_V2=0.012531328320802;
@ -2315,11 +2272,10 @@ __global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dis
// read the component number densities
nA = Den[n];
nB = Den[Np + n];
porosity = Poros[n];
perm = Perm[n];
nx_gs = GreySolidGrad[n+0*Np];
ny_gs = GreySolidGrad[n+1*Np];
nz_gs = GreySolidGrad[n+2*Np];
W = GreySolidW[n];
// compute phase indicator field
phi=(nA-nB)/(nA+nB);
@ -2341,98 +2297,61 @@ __global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dis
//........................................................................
nn = ijk-1; // neighbor index (get convention)
m1 = Phi[nn]; // get neighbor for phi - 1
gp1 = Psi[nn];
//........................................................................
nn = ijk+1; // neighbor index (get convention)
m2 = Phi[nn]; // get neighbor for phi - 2
gp2 = Psi[nn];
//........................................................................
nn = ijk-strideY; // neighbor index (get convention)
m3 = Phi[nn]; // get neighbor for phi - 3
gp3 = Psi[nn];
//........................................................................
nn = ijk+strideY; // neighbor index (get convention)
m4 = Phi[nn]; // get neighbor for phi - 4
gp4 = Psi[nn];
//........................................................................
nn = ijk-strideZ; // neighbor index (get convention)
m5 = Phi[nn]; // get neighbor for phi - 5
gp5 = Psi[nn];
//........................................................................
nn = ijk+strideZ; // neighbor index (get convention)
m6 = Phi[nn]; // get neighbor for phi - 6
gp6 = Psi[nn];
//........................................................................
nn = ijk-strideY-1; // neighbor index (get convention)
m7 = Phi[nn]; // get neighbor for phi - 7
gp7 = Psi[nn];
//........................................................................
nn = ijk+strideY+1; // neighbor index (get convention)
m8 = Phi[nn]; // get neighbor for phi - 8
gp8 = Psi[nn];
//........................................................................
nn = ijk+strideY-1; // neighbor index (get convention)
m9 = Phi[nn]; // get neighbor for phi - 9
gp9 = Psi[nn];
//........................................................................
nn = ijk-strideY+1; // neighbor index (get convention)
m10 = Phi[nn]; // get neighbor for phi - 10
gp10 = Psi[nn];
//........................................................................
nn = ijk-strideZ-1; // neighbor index (get convention)
m11 = Phi[nn]; // get neighbor for phi - 11
gp11 = Psi[nn];
//........................................................................
nn = ijk+strideZ+1; // neighbor index (get convention)
m12 = Phi[nn]; // get neighbor for phi - 12
gp12 = Psi[nn];
//........................................................................
nn = ijk+strideZ-1; // neighbor index (get convention)
m13 = Phi[nn]; // get neighbor for phi - 13
gp13 = Psi[nn];
//........................................................................
nn = ijk-strideZ+1; // neighbor index (get convention)
m14 = Phi[nn]; // get neighbor for phi - 14
gp14 = Psi[nn];
//........................................................................
nn = ijk-strideZ-strideY; // neighbor index (get convention)
m15 = Phi[nn]; // get neighbor for phi - 15
gp15 = Psi[nn];
//........................................................................
nn = ijk+strideZ+strideY; // neighbor index (get convention)
m16 = Phi[nn]; // get neighbor for phi - 16
gp16 = Psi[nn];
//........................................................................
nn = ijk+strideZ-strideY; // neighbor index (get convention)
m17 = Phi[nn]; // get neighbor for phi - 17
gp17 = Psi[nn];
//........................................................................
nn = ijk-strideZ+strideY; // neighbor index (get convention)
m18 = Phi[nn]; // get neighbor for phi - 18
gp18 = Psi[nn];
//............Compute the Color Gradient...................................
nx_phase = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
ny_phase = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
nz_phase = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
C_phase = sqrt(nx_phase*nx_phase+ny_phase*ny_phase+nz_phase*nz_phase);
//correct the normal color gradient by considering the effect of grey solid
nx = nx_phase + (1.0-porosity)*nx_gs;
ny = ny_phase + (1.0-porosity)*ny_gs;
nz = nz_phase + (1.0-porosity)*nz_gs;
if (C_phase==0.0){
nx = nx_phase;
ny = ny_phase;
nz = nz_phase;
}
//...........Normalize the Color Gradient.................................
C = sqrt(nx*nx+ny*ny+nz*nz);
double ColorMag = C;
if (C==0.0) ColorMag=1.0;
nx = nx/ColorMag;
ny = ny/ColorMag;
nz = nz/ColorMag;
nx = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
ny = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
nz = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
//............Compute the Greyscale Potential Gradient.....................
// Fcpx = 0.0;
@ -2455,14 +2374,23 @@ __global__ void dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dis
// ny = Fcpy/Fcp_mag;
// nz = Fcpz/Fcp_mag;
// }
Fcpx = -3.0/18.0*(m1-m2+0.5*(m7-m8+m9-m10+m11-m12+m13-m14));
Fcpy = -3.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18));
Fcpz = -3.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18));
Fcpx = nx;
Fcpy = ny;
Fcpz = nz;
double Fcp_mag=sqrt(Fcpx*Fcpx+Fcpy*Fcpy+Fcpz*Fcpz);
if (Fcp_mag==0.0); Fcpx=Fcpy=Fcpz=0.0;
//NOTE for open node (porosity=1.0),Fcp=0.0
Fcpx *= alpha*W*(1.0-porosity)/sqrt(perm);
Fcpy *= alpha*W*(1.0-porosity)/sqrt(perm);
Fcpz *= alpha*W*(1.0-porosity)/sqrt(perm);
//...........Normalize the Color Gradient.................................
C = sqrt(nx*nx+ny*ny+nz*nz);
double ColorMag = C;
if (C==0.0) ColorMag=1.0;
nx = nx/ColorMag;
ny = ny/ColorMag;
nz = nz/ColorMag;
// q=0
fq = dist[n];
@ -4572,12 +4500,12 @@ extern "C" void ScaLBL_PhaseField_InitFromRestart(double *Den, double *Aq, doubl
//Model-1 & 4 with capillary pressure penalty
extern "C" void ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi,double *Psi, double *GreySolidGrad, double *Poros,double *Perm,double *Vel, double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm,double *Vel, double *Pressure,
double rhoA, double rhoB, double tauA, double tauB,double tauA_eff,double tauB_eff, double alpha, double beta,
double Fx, double Fy, double Fz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np){
double Fx, double Fy, double Fz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np){
dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP<<<NBLOCKS,NTHREADS >>>(Map, dist, Aq, Bq, Den, Phi, Psi, GreySolidGrad, Poros, Perm, Vel, Pressure,
rhoA, rhoB, tauA, tauB, tauA_eff, tauB_eff, alpha, beta, Fx, Fy, Fz, RecoloringOff, W, strideY, strideZ, start, finish, Np);
dvc_ScaLBL_D3Q19_AAeven_GreyscaleColor_CP<<<NBLOCKS,NTHREADS >>>(Map, dist, Aq, Bq, Den, Phi, GreySolidW, Poros, Perm, Vel, Pressure,
rhoA, rhoB, tauA, tauB, tauA_eff, tauB_eff, alpha, beta, Fx, Fy, Fz, RecoloringOff, strideY, strideZ, start, finish, Np);
cudaError_t err = cudaGetLastError();
if (cudaSuccess != err){
printf("CUDA error in ScaLBL_D3Q19_AAeven_GreyscaleColor_CP: %s \n",cudaGetErrorString(err));
@ -4587,12 +4515,12 @@ extern "C" void ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(int *Map, double *dist, do
//Model-1 & 4 with capillary pressure penalty
extern "C" void ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(int *d_neighborList, int *Map, double *dist, double *Aq, double *Bq, double *Den,
double *Phi, double *Psi, double *GreySolidGrad, double *Poros,double *Perm,double *Vel,double *Pressure,
double *Phi, double *GreySolidW, double *Poros,double *Perm,double *Vel,double *Pressure,
double rhoA, double rhoB, double tauA, double tauB, double tauA_eff,double tauB_eff, double alpha, double beta,
double Fx, double Fy, double Fz, bool RecoloringOff, double W, int strideY, int strideZ, int start, int finish, int Np){
double Fx, double Fy, double Fz, bool RecoloringOff, int strideY, int strideZ, int start, int finish, int Np){
dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP<<<NBLOCKS,NTHREADS >>>(d_neighborList, Map, dist, Aq, Bq, Den, Phi, Psi, GreySolidGrad, Poros, Perm,Vel,Pressure,
rhoA, rhoB, tauA, tauB, tauA_eff, tauB_eff,alpha, beta, Fx, Fy, Fz, RecoloringOff, W, strideY, strideZ, start, finish, Np);
dvc_ScaLBL_D3Q19_AAodd_GreyscaleColor_CP<<<NBLOCKS,NTHREADS >>>(d_neighborList, Map, dist, Aq, Bq, Den, Phi, GreySolidW, Poros, Perm,Vel,Pressure,
rhoA, rhoB, tauA, tauB, tauA_eff, tauB_eff,alpha, beta, Fx, Fy, Fz, RecoloringOff, strideY, strideZ, start, finish, Np);
cudaError_t err = cudaGetLastError();
if (cudaSuccess != err){

22
hip/FreeLee.cu
View File

@ -1,6 +1,7 @@
#include <math.h>
#include <stdio.h>
#include <cuda_profiler_api.h>
#include "hip/hip_runtime.h"
#include "hip/hip_runtime_api.h"
#define STOKES
@ -3395,7 +3396,8 @@ extern "C" void ScaLBL_FreeLeeModel_PhaseField_Init(int *Map, double *Phi, doubl
extern "C" void ScaLBL_D3Q7_AAodd_FreeLee_PhaseField(int *neighborList, int *Map, double *hq, double *Den, double *Phi, double *ColorGrad, double *Vel,
double rhoA, double rhoB, double tauM, double W, int start, int finish, int Np)
{
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q7_AAodd_FreeLee_PhaseField, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q7_AAodd_FreeLee_PhaseField, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q7_AAodd_FreeLee_PhaseField<<<NBLOCKS,NTHREADS >>>(neighborList, Map, hq, Den, Phi, ColorGrad, Vel,
rhoA, rhoB, tauM, W, start, finish, Np);
hipError_t err = hipGetLastError();
@ -3407,7 +3409,8 @@ extern "C" void ScaLBL_D3Q7_AAodd_FreeLee_PhaseField(int *neighborList, int *Map
extern "C" void ScaLBL_D3Q7_AAeven_FreeLee_PhaseField( int *Map, double *hq, double *Den, double *Phi, double *ColorGrad, double *Vel,
double rhoA, double rhoB, double tauM, double W, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q7_AAeven_FreeLee_PhaseField, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q7_AAeven_FreeLee_PhaseField, hipFuncCachePreferL1);

dvc_ScaLBL_D3Q7_AAeven_FreeLee_PhaseField<<<NBLOCKS,NTHREADS >>>( Map, hq, Den, Phi, ColorGrad, Vel, rhoA, rhoB, tauM, W, start, finish, Np);
hipError_t err = hipGetLastError();
if (hipSuccess != err){
@ -3418,7 +3421,7 @@ extern "C" void ScaLBL_D3Q7_AAeven_FreeLee_PhaseField( int *Map, double *hq, dou
extern "C" void ScaLBL_D3Q7_ComputePhaseField(int *Map, double *hq, double *Den, double *Phi, double rhoA, double rhoB, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q7_ComputePhaseField, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q7_ComputePhaseField, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q7_ComputePhaseField<<<NBLOCKS,NTHREADS >>>( Map, hq, Den, Phi, rhoA, rhoB, start, finish, Np);
hipError_t err = hipGetLastError();
if (hipSuccess != err){
@ -3431,7 +3434,7 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel(int *neighborList, int *Map, dou
double rhoA, double rhoB, double tauA, double tauB, double kappa, double beta, double W, double Fx, double Fy, double Fz,
int strideY, int strideZ, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel<<<NBLOCKS,NTHREADS >>>(neighborList, Map, dist, Den, Phi, mu_phi, Vel, Pressure, ColorGrad,
rhoA, rhoB, tauA, tauB, kappa, beta, W, Fx, Fy, Fz, strideY, strideZ, start, finish, Np);
hipError_t err = hipGetLastError();
@ -3444,7 +3447,7 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel(int *Map, double *dist, double
double rhoA, double rhoB, double tauA, double tauB, double kappa, double beta, double W, double Fx, double Fy, double Fz,
int strideY, int strideZ, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel<<<NBLOCKS,NTHREADS >>>(Map, dist, Den, Phi, mu_phi, Vel, Pressure, ColorGrad,
rhoA, rhoB, tauA, tauB, kappa, beta, W, Fx, Fy, Fz, strideY, strideZ, start, finish, Np);
hipError_t err = hipGetLastError();
@ -3502,7 +3505,7 @@ extern "C" void ScaLBL_D3Q7_AAeven_FreeLeeModel_PhaseField(int *Map, double *hq,
extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_SingleFluid_BGK(int *neighborList, double *dist, double *Vel, double *Pressure,
double tau, double rho0, double Fx, double Fy, double Fz, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel_SingleFluid_BGK, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel_SingleFluid_BGK, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q19_AAodd_FreeLeeModel_SingleFluid_BGK<<<NBLOCKS,NTHREADS >>>(neighborList, dist, Vel, Pressure,
tau, rho0, Fx, Fy, Fz, start, finish, Np);
hipError_t err = hipGetLastError();
@ -3514,7 +3517,8 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_SingleFluid_BGK(int *neighborLis
extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_SingleFluid_BGK(double *dist, double *Vel, double *Pressure,
double tau, double rho0, double Fx, double Fy, double Fz, int start, int finish, int Np){
cudaFuncSetCacheConfig(dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel_SingleFluid_BGK, cudaFuncCachePreferL1);
hipFuncSetCacheConfig((void*)dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel_SingleFluid_BGK, hipFuncCachePreferL1);
dvc_ScaLBL_D3Q19_AAeven_FreeLeeModel_SingleFluid_BGK<<<NBLOCKS,NTHREADS >>>(dist, Vel, Pressure,
tau, rho0, Fx, Fy, Fz, start, finish, Np);
hipError_t err = hipGetLastError();
@ -3525,4 +3529,4 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_SingleFluid_BGK(double *dist, d
extern "C" void ScaLBL_D3Q9_MGTest(int *Map, double *Phi,double *ColorGrad,int strideY, int strideZ, int start, int finish, int Np){
}
}

498
models/GreyscaleColorModel.cpp
View File

@ -17,7 +17,7 @@ void DeleteArray( const TYPE *p )
ScaLBL_GreyscaleColorModel::ScaLBL_GreyscaleColorModel(int RANK, int NP, const Utilities::MPI& COMM):
rank(RANK), nprocs(NP), Restart(0),timestep(0),timestepMax(0),tauA(0),tauB(0),tauA_eff(0),tauB_eff(0),rhoA(0),rhoB(0),alpha(0),beta(0),
Fx(0),Fy(0),Fz(0),flux(0),din(0),dout(0),inletA(0),inletB(0),outletA(0),outletB(0),GreyPorosity(0),RecoloringOff(0),W(0),
Fx(0),Fy(0),Fz(0),flux(0),din(0),dout(0),inletA(0),inletB(0),outletA(0),outletB(0),GreyPorosity(0),RecoloringOff(0),
Nx(0),Ny(0),Nz(0),N(0),Np(0),nprocx(0),nprocy(0),nprocz(0),BoundaryCondition(0),Lx(0),Ly(0),Lz(0),comm(COMM)
{
REVERSE_FLOW_DIRECTION = false;
@ -44,7 +44,7 @@ void ScaLBL_GreyscaleColorModel::ReadParams(string filename){
din=dout=1.0;
flux=0.0;
RecoloringOff = false;
W=1.0;
//W=1.0;
// Color Model parameters
if (greyscaleColor_db->keyExists( "timestepMax" )){
@ -90,9 +90,6 @@ void ScaLBL_GreyscaleColorModel::ReadParams(string filename){
if (greyscaleColor_db->keyExists( "RecoloringOff" )){
RecoloringOff = greyscaleColor_db->getScalar<bool>( "RecoloringOff" );
}
if (greyscaleColor_db->keyExists( "W" )){
W = greyscaleColor_db->getScalar<double>( "W" );
}
inletA=1.f;
inletB=0.f;
outletA=0.f;
@ -300,19 +297,18 @@ void ScaLBL_GreyscaleColorModel::AssignComponentLabels()
delete [] phase;
}
void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()//Model-4
void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()//apply capillary penalty wetting strength W
{
// ONLY initialize grey nodes
// Key input parameters:
// 1. GreySolidLabels
// labels for grey nodes
// 2. GreySolidAffinity
// affinity ranges [-1,1]
// oil-wet > 0
// water-wet < 0
// neutral = 0
double *SolidPotential_host = new double [Nx*Ny*Nz];
double *GreySolidGrad_host = new double [3*Np];
// ranges [-1,1]
// water-wet > 0
// oil-wet < 0
// neutral = 0 (i.e. no penalty)
double *GreySolidW_host = new double [Np];
size_t NLABELS=0;
signed char VALUE=0;
@ -334,117 +330,19 @@ void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()//Model-4
AFFINITY=0.f;//all nodes except the specified grey nodes have grey-solid affinity = 0.0
// Assign the affinity from the paired list
for (unsigned int idx=0; idx < NLABELS; idx++){
//printf("idx=%i, value=%i, %i, \n",idx, VALUE,LabelList[idx]);
if (VALUE == LabelList[idx]){
AFFINITY=AffinityList[idx];
idx = NLABELS;
//Mask->id[n] = 0; // set mask to zero since this is an immobile component
}
}
SolidPotential_host[n] = AFFINITY;
}
}
}
// Calculate grey-solid color-gradient
double *Dst;
Dst = new double [3*3*3];
for (int kk=0; kk<3; kk++){
for (int jj=0; jj<3; jj++){
for (int ii=0; ii<3; ii++){
int index = kk*9+jj*3+ii;
Dst[index] = sqrt(double(ii-1)*double(ii-1) + double(jj-1)*double(jj-1)+ double(kk-1)*double(kk-1));
}
}
}
double w_face = 1.f;
double w_edge = 0.5;
double w_corner = 0.f;
//local
Dst[13] = 0.f;
//faces
Dst[4] = w_face;
Dst[10] = w_face;
Dst[12] = w_face;
Dst[14] = w_face;
Dst[16] = w_face;
Dst[22] = w_face;
// corners
Dst[0] = w_corner;
Dst[2] = w_corner;
Dst[6] = w_corner;
Dst[8] = w_corner;
Dst[18] = w_corner;
Dst[20] = w_corner;
Dst[24] = w_corner;
Dst[26] = w_corner;
// edges
Dst[1] = w_edge;
Dst[3] = w_edge;
Dst[5] = w_edge;
Dst[7] = w_edge;
Dst[9] = w_edge;
Dst[11] = w_edge;
Dst[15] = w_edge;
Dst[17] = w_edge;
Dst[19] = w_edge;
Dst[21] = w_edge;
Dst[23] = w_edge;
Dst[25] = w_edge;
for (int k=1; k<Nz-1; k++){
for (int j=1; j<Ny-1; j++){
for (int i=1; i<Nx-1; i++){
int idx=Map(i,j,k);
int idx = Map(i,j,k);
if (!(idx < 0)){
double phi_x = 0.f;
double phi_y = 0.f;
double phi_z = 0.f;
for (int kk=0; kk<3; kk++){
for (int jj=0; jj<3; jj++){
for (int ii=0; ii<3; ii++){
int index = kk*9+jj*3+ii;
double weight= Dst[index];
int idi=i+ii-1;
int idj=j+jj-1;
int idk=k+kk-1;
if (idi < 0) idi=0;
if (idj < 0) idj=0;
if (idk < 0) idk=0;
if (!(idi < Nx)) idi=Nx-1;
if (!(idj < Ny)) idj=Ny-1;
if (!(idk < Nz)) idk=Nz-1;
int nn = idk*Nx*Ny + idj*Nx + idi;
double vec_x = double(ii-1);
double vec_y = double(jj-1);
double vec_z = double(kk-1);
double GWNS=SolidPotential_host[nn];
phi_x += GWNS*weight*vec_x;
phi_y += GWNS*weight*vec_y;
phi_z += GWNS*weight*vec_z;
}
}
}
if (Averages->SDs(i,j,k)<2.0){
GreySolidGrad_host[idx+0*Np] = phi_x;
GreySolidGrad_host[idx+1*Np] = phi_y;
GreySolidGrad_host[idx+2*Np] = phi_z;
}
else{
GreySolidGrad_host[idx+0*Np] = 0.0;
GreySolidGrad_host[idx+1*Np] = 0.0;
GreySolidGrad_host[idx+2*Np] = 0.0;
}
}
GreySolidW_host[idx] = AFFINITY;
}
}
}
}
if (rank==0){
printf("Number of Grey-solid labels: %lu \n",NLABELS);
for (unsigned int idx=0; idx<NLABELS; idx++){
@ -452,14 +350,13 @@ void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()//Model-4
AFFINITY=AffinityList[idx];
printf(" grey-solid label=%d, grey-solid affinity=%f\n",VALUE,AFFINITY);
}
printf("NOTE: grey-solid affinity>0: water-wet || grey-solid affinity<0: oil-wet \n");
}
ScaLBL_CopyToDevice(GreySolidGrad, GreySolidGrad_host, 3*Np*sizeof(double));
ScaLBL_CopyToDevice(GreySolidW, GreySolidW_host, Np*sizeof(double));
ScaLBL_Comm->Barrier();
delete [] SolidPotential_host;
delete [] GreySolidGrad_host;
delete [] Dst;
delete [] GreySolidW_host;
}
////----------------------------------------------------------------------------------------------------------//
@ -585,70 +482,70 @@ void ScaLBL_GreyscaleColorModel::AssignGreyPoroPermLabels()
delete [] Permeability;
}
void ScaLBL_GreyscaleColorModel::AssignGreyscalePotential()
{
double *psi;//greyscale potential
psi = new double[N];
size_t NLABELS=0;
signed char VALUE=0;
double AFFINITY=0.f;
auto LabelList = greyscaleColor_db->getVector<int>( "ComponentLabels" );
auto AffinityList = greyscaleColor_db->getVector<double>( "ComponentAffinity" );
NLABELS=LabelList.size();
//first, copy over normal phase field
for (int k=0;k<Nz;k++){
for (int j=0;j<Ny;j++){
for (int i=0;i<Nx;i++){
int n = k*Nx*Ny+j*Nx+i;
VALUE=id[n];
// Assign the affinity from the paired list
for (unsigned int idx=0; idx < NLABELS; idx++){
//printf("idx=%i, value=%i, %i, \n",idx, VALUE,LabelList[idx]);
if (VALUE == LabelList[idx]){
AFFINITY=AffinityList[idx];
idx = NLABELS;
}
}
// fluid labels are reserved
if (VALUE == 1) AFFINITY=1.0;
else if (VALUE == 2) AFFINITY=-1.0;
psi[n] = AFFINITY;
}
}
}
//second, scale the phase field for grey nodes
double Cap_Penalty=1.f;
auto GreyLabelList = greyscaleColor_db->getVector<int>( "GreySolidLabels" );
auto PermeabilityList = greyscaleColor_db->getVector<double>( "PermeabilityList" );
NLABELS=GreyLabelList.size();
for (int k=0;k<Nz;k++){
for (int j=0;j<Ny;j++){
for (int i=0;i<Nx;i++){
int n = k*Nx*Ny+j*Nx+i;
VALUE=id[n];
Cap_Penalty=1.f;
// Assign the affinity from the paired list
for (unsigned int idx=0; idx < NLABELS; idx++){
if (VALUE == GreyLabelList[idx]){
Cap_Penalty=alpha*W/sqrt(PermeabilityList[idx]/Dm->voxel_length/Dm->voxel_length);
idx = NLABELS;
}
}
//update greyscale potential
psi[n] = psi[n]*Cap_Penalty;
}
}
}
ScaLBL_CopyToDevice(Psi, psi, N*sizeof(double));
ScaLBL_Comm->Barrier();
delete [] psi;
}
//void ScaLBL_GreyscaleColorModel::AssignGreyscalePotential()
//{
// double *psi;//greyscale potential
// psi = new double[N];
//
// size_t NLABELS=0;
// signed char VALUE=0;
// double AFFINITY=0.f;
//
// auto LabelList = greyscaleColor_db->getVector<int>( "ComponentLabels" );
// auto AffinityList = greyscaleColor_db->getVector<double>( "ComponentAffinity" );
// NLABELS=LabelList.size();
//
// //first, copy over normal phase field
// for (int k=0;k<Nz;k++){
// for (int j=0;j<Ny;j++){
// for (int i=0;i<Nx;i++){
// int n = k*Nx*Ny+j*Nx+i;
// VALUE=id[n];
// // Assign the affinity from the paired list
// for (unsigned int idx=0; idx < NLABELS; idx++){
// //printf("idx=%i, value=%i, %i, \n",idx, VALUE,LabelList[idx]);
// if (VALUE == LabelList[idx]){
// AFFINITY=AffinityList[idx];
// idx = NLABELS;
// }
// }
// // fluid labels are reserved
// if (VALUE == 1) AFFINITY=1.0;
// else if (VALUE == 2) AFFINITY=-1.0;
// psi[n] = AFFINITY;
// }
// }
// }
//
// //second, scale the phase field for grey nodes
// double Cap_Penalty=1.f;
// auto GreyLabelList = greyscaleColor_db->getVector<int>( "GreySolidLabels" );
// auto PermeabilityList = greyscaleColor_db->getVector<double>( "PermeabilityList" );
// NLABELS=GreyLabelList.size();
//
// for (int k=0;k<Nz;k++){
// for (int j=0;j<Ny;j++){
// for (int i=0;i<Nx;i++){
// int n = k*Nx*Ny+j*Nx+i;
// VALUE=id[n];
// Cap_Penalty=1.f;
// // Assign the affinity from the paired list
// for (unsigned int idx=0; idx < NLABELS; idx++){
// if (VALUE == GreyLabelList[idx]){
// Cap_Penalty=alpha*W/sqrt(PermeabilityList[idx]/Dm->voxel_length/Dm->voxel_length);
// idx = NLABELS;
// }
// }
// //update greyscale potential
// psi[n] = psi[n]*Cap_Penalty;
// }
// }
// }
//
// ScaLBL_CopyToDevice(Psi, psi, N*sizeof(double));
// ScaLBL_Comm->Barrier();
// delete [] psi;
//}
void ScaLBL_GreyscaleColorModel::Create(){
/*
@ -670,7 +567,6 @@ void ScaLBL_GreyscaleColorModel::Create(){
// ScaLBL_Communicator ScaLBL_Comm(Mask); // original
ScaLBL_Comm = std::shared_ptr<ScaLBL_Communicator>(new ScaLBL_Communicator(Mask));
ScaLBL_Comm_Regular = std::shared_ptr<ScaLBL_Communicator>(new ScaLBL_Communicator(Mask));
ScaLBL_Comm_Regular_2 = std::shared_ptr<ScaLBL_Communicator>(new ScaLBL_Communicator(Mask));
int Npad=(Np/16 + 2)*16;
if (rank==0) printf ("Set up memory efficient layout, %i | %i | %i \n", Np, Npad, N);
@ -695,12 +591,13 @@ void ScaLBL_GreyscaleColorModel::Create(){
ScaLBL_AllocateDeviceMemory((void **) &Bq, 7*dist_mem_size);
ScaLBL_AllocateDeviceMemory((void **) &Den, 2*dist_mem_size);
ScaLBL_AllocateDeviceMemory((void **) &Phi, sizeof(double)*Nx*Ny*Nz);
ScaLBL_AllocateDeviceMemory((void **) &Psi, sizeof(double)*Nx*Ny*Nz);//greyscale potential
//ScaLBL_AllocateDeviceMemory((void **) &Psi, sizeof(double)*Nx*Ny*Nz);//greyscale potential
ScaLBL_AllocateDeviceMemory((void **) &Pressure, sizeof(double)*Np);
ScaLBL_AllocateDeviceMemory((void **) &Velocity, 3*sizeof(double)*Np);
//ScaLBL_AllocateDeviceMemory((void **) &ColorGrad, 3*sizeof(double)*Np);
//ScaLBL_AllocateDeviceMemory((void **) &GreySolidPhi, sizeof(double)*Nx*Ny*Nz);
ScaLBL_AllocateDeviceMemory((void **) &GreySolidGrad, 3*sizeof(double)*Np);
//ScaLBL_AllocateDeviceMemory((void **) &GreySolidGrad, 3*sizeof(double)*Np);
ScaLBL_AllocateDeviceMemory((void **) &GreySolidW, sizeof(double)*Np);
ScaLBL_AllocateDeviceMemory((void **) &Porosity_dvc, sizeof(double)*Np);
ScaLBL_AllocateDeviceMemory((void **) &Permeability_dvc, sizeof(double)*Np);
//...........................................................................
@ -744,7 +641,7 @@ void ScaLBL_GreyscaleColorModel::Create(){
AssignComponentLabels();//do open/black/grey nodes initialization
AssignGreySolidLabels();
AssignGreyPoroPermLabels();
AssignGreyscalePotential();
//AssignGreyscalePotential();
Averages->SetParams(rhoA,rhoB,tauA,tauB,Fx,Fy,Fz,alpha,beta,GreyPorosity);
ScaLBL_Comm->RegularLayout(Map,Porosity_dvc,Averages->Porosity);//porosity doesn't change over time
}
@ -1023,11 +920,10 @@ void ScaLBL_GreyscaleColorModel::Run(){
}
// Halo exchange for phase field
ScaLBL_Comm_Regular->SendHalo(Phi);
ScaLBL_Comm_Regular_2->SendHalo(Psi);
//Model-1&4 with capillary pressure penalty for grey nodes
ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Psi, GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, GreySolidW,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
alpha, beta, Fx, Fy, Fz, RecoloringOff, W, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
alpha, beta, Fx, Fy, Fz, RecoloringOff, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
//Model-1&4
//ScaLBL_D3Q19_AAodd_GreyscaleColor(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi,GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
@ -1036,7 +932,6 @@ void ScaLBL_GreyscaleColorModel::Run(){
//ScaLBL_D3Q19_AAodd_GreyscaleColor(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi,GreySolidPhi,Porosity_dvc,Permeability_dvc,Velocity,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
// alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm_Regular_2->RecvHalo(Psi);
ScaLBL_Comm_Regular->RecvHalo(Phi);
ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE
ScaLBL_Comm->Barrier();
@ -1055,9 +950,9 @@ void ScaLBL_GreyscaleColorModel::Run(){
}
//Model-1&4 with capillary pressure penalty for grey nodes
ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, Psi, GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
ScaLBL_D3Q19_AAodd_GreyscaleColor_CP(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi, GreySolidW,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
alpha, beta, Fx, Fy, Fz, RecoloringOff, W, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
alpha, beta, Fx, Fy, Fz, RecoloringOff, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
//Model-1&4
//ScaLBL_D3Q19_AAodd_GreyscaleColor(NeighborList, dvcMap, fq, Aq, Bq, Den, Phi,GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
@ -1087,11 +982,10 @@ void ScaLBL_GreyscaleColorModel::Run(){
ScaLBL_Comm->Color_BC_Z(dvcMap, Phi, Den, outletA, outletB);
}
ScaLBL_Comm_Regular->SendHalo(Phi);
ScaLBL_Comm_Regular_2->SendHalo(Psi);
//Model-1&4 with capillary pressure penalty for grey nodes
ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(dvcMap, fq, Aq, Bq, Den, Phi, Psi, GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(dvcMap, fq, Aq, Bq, Den, Phi, GreySolidW,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
alpha, beta, Fx, Fy, Fz, RecoloringOff, W, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
alpha, beta, Fx, Fy, Fz, RecoloringOff, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
//Model-1&4
//ScaLBL_D3Q19_AAeven_GreyscaleColor(dvcMap, fq, Aq, Bq, Den, Phi,GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
@ -1100,7 +994,6 @@ void ScaLBL_GreyscaleColorModel::Run(){
//ScaLBL_D3Q19_AAeven_GreyscaleColor(dvcMap, fq, Aq, Bq, Den, Phi,GreySolidPhi,Porosity_dvc,Permeability_dvc,Velocity,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
// alpha, beta, Fx, Fy, Fz, Nx, Nx*Ny, ScaLBL_Comm->FirstInterior(), ScaLBL_Comm->LastInterior(), Np);
ScaLBL_Comm_Regular_2->RecvHalo(Psi);
ScaLBL_Comm_Regular->RecvHalo(Phi);
ScaLBL_Comm->RecvD3Q19AA(fq); //WRITE INTO OPPOSITE
ScaLBL_Comm->Barrier();
@ -1119,9 +1012,9 @@ void ScaLBL_GreyscaleColorModel::Run(){
}
//Model-1&4 with capillary pressure penalty for grey nodes
ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(dvcMap, fq, Aq, Bq, Den, Phi, Psi, GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
ScaLBL_D3Q19_AAeven_GreyscaleColor_CP(dvcMap, fq, Aq, Bq, Den, Phi, GreySolidW,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
alpha, beta, Fx, Fy, Fz, RecoloringOff, W, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
alpha, beta, Fx, Fy, Fz, RecoloringOff, Nx, Nx*Ny, 0, ScaLBL_Comm->LastExterior(), Np);
//Model-1&4
//ScaLBL_D3Q19_AAeven_GreyscaleColor(dvcMap, fq, Aq, Bq, Den, Phi,GreySolidGrad,Porosity_dvc,Permeability_dvc,Velocity,Pressure,
// rhoA, rhoB, tauA, tauB,tauA_eff, tauB_eff,
@ -1677,12 +1570,12 @@ void ScaLBL_GreyscaleColorModel::WriteDebug(){
fwrite(PhaseField.data(),8,N,OUTFILE);
fclose(OUTFILE);
ScaLBL_CopyToHost(PhaseField.data(), Psi, sizeof(double)*N);
FILE *PSIFILE;
sprintf(LocalRankFilename,"Psi.%05i.raw",rank);
PSIFILE = fopen(LocalRankFilename,"wb");
fwrite(PhaseField.data(),8,N,PSIFILE);
fclose(PSIFILE);
//ScaLBL_CopyToHost(PhaseField.data(), Psi, sizeof(double)*N);
//FILE *PSIFILE;
//sprintf(LocalRankFilename,"Psi.%05i.raw",rank);
//PSIFILE = fopen(LocalRankFilename,"wb");
//fwrite(PhaseField.data(),8,N,PSIFILE);
//fclose(PSIFILE);
ScaLBL_Comm->RegularLayout(Map,&Den[0],PhaseField);
FILE *AFILE;
@ -1740,26 +1633,26 @@ void ScaLBL_GreyscaleColorModel::WriteDebug(){
fwrite(PhaseField.data(),8,N,PERM_FILE);
fclose(PERM_FILE);
ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[0],PhaseField);
FILE *GreySG_X_FILE;
sprintf(LocalRankFilename,"GreySolidGrad_X.%05i.raw",rank);
GreySG_X_FILE = fopen(LocalRankFilename,"wb");
fwrite(PhaseField.data(),8,N,GreySG_X_FILE);
fclose(GreySG_X_FILE);
//ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[0],PhaseField);
//FILE *GreySG_X_FILE;
//sprintf(LocalRankFilename,"GreySolidGrad_X.%05i.raw",rank);
//GreySG_X_FILE = fopen(LocalRankFilename,"wb");
//fwrite(PhaseField.data(),8,N,GreySG_X_FILE);
//fclose(GreySG_X_FILE);
ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[Np],PhaseField);
FILE *GreySG_Y_FILE;
sprintf(LocalRankFilename,"GreySolidGrad_Y.%05i.raw",rank);
GreySG_Y_FILE = fopen(LocalRankFilename,"wb");
fwrite(PhaseField.data(),8,N,GreySG_Y_FILE);
fclose(GreySG_Y_FILE);
//ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[Np],PhaseField);
//FILE *GreySG_Y_FILE;
//sprintf(LocalRankFilename,"GreySolidGrad_Y.%05i.raw",rank);
//GreySG_Y_FILE = fopen(LocalRankFilename,"wb");
//fwrite(PhaseField.data(),8,N,GreySG_Y_FILE);
//fclose(GreySG_Y_FILE);
ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[2*Np],PhaseField);
FILE *GreySG_Z_FILE;
sprintf(LocalRankFilename,"GreySolidGrad_Z.%05i.raw",rank);
GreySG_Z_FILE = fopen(LocalRankFilename,"wb");
fwrite(PhaseField.data(),8,N,GreySG_Z_FILE);
fclose(GreySG_Z_FILE);
//ScaLBL_Comm->RegularLayout(Map,&GreySolidGrad[2*Np],PhaseField);
//FILE *GreySG_Z_FILE;
//sprintf(LocalRankFilename,"GreySolidGrad_Z.%05i.raw",rank);
//GreySG_Z_FILE = fopen(LocalRankFilename,"wb");
//fwrite(PhaseField.data(),8,N,GreySG_Z_FILE);
//fclose(GreySG_Z_FILE);
/* ScaLBL_Comm->RegularLayout(Map,&ColorGrad[0],PhaseField);
FILE *CGX_FILE;
@ -2093,6 +1986,169 @@ void ScaLBL_GreyscaleColorModel::WriteDebug(){
// delete [] Dst;
//}
//void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()//Model-4
//{
// // ONLY initialize grey nodes
// // Key input parameters:
// // 1. GreySolidLabels
// // labels for grey nodes
// // 2. GreySolidAffinity
// // affinity ranges [-1,1]
// // oil-wet > 0
// // water-wet < 0
// // neutral = 0
// double *SolidPotential_host = new double [Nx*Ny*Nz];
// double *GreySolidGrad_host = new double [3*Np];
//
// size_t NLABELS=0;
// signed char VALUE=0;
// double AFFINITY=0.f;
//
// auto LabelList = greyscaleColor_db->getVector<int>( "GreySolidLabels" );
// auto AffinityList = greyscaleColor_db->getVector<double>( "GreySolidAffinity" );
//
// NLABELS=LabelList.size();
// if (NLABELS != AffinityList.size()){
// ERROR("Error: GreySolidLabels and GreySolidAffinity must be the same length! \n");
// }
//
// for (int k=0;k<Nz;k++){
// for (int j=0;j<Ny;j++){
// for (int i=0;i<Nx;i++){
// int n = k*Nx*Ny+j*Nx+i;
// VALUE=id[n];
// AFFINITY=0.f;//all nodes except the specified grey nodes have grey-solid affinity = 0.0
// // Assign the affinity from the paired list
// for (unsigned int idx=0; idx < NLABELS; idx++){
// //printf("idx=%i, value=%i, %i, \n",idx, VALUE,LabelList[idx]);
// if (VALUE == LabelList[idx]){
// AFFINITY=AffinityList[idx];
// idx = NLABELS;
// //Mask->id[n] = 0; // set mask to zero since this is an immobile component
// }
// }
// SolidPotential_host[n] = AFFINITY;
// }
// }
// }
//
// // Calculate grey-solid color-gradient
// double *Dst;
// Dst = new double [3*3*3];
// for (int kk=0; kk<3; kk++){
// for (int jj=0; jj<3; jj++){
// for (int ii=0; ii<3; ii++){
// int index = kk*9+jj*3+ii;
// Dst[index] = sqrt(double(ii-1)*double(ii-1) + double(jj-1)*double(jj-1)+ double(kk-1)*double(kk-1));
// }
// }
// }
// double w_face = 1.f;
// double w_edge = 0.5;
// double w_corner = 0.f;
// //local
// Dst[13] = 0.f;
// //faces
// Dst[4] = w_face;
// Dst[10] = w_face;
// Dst[12] = w_face;
// Dst[14] = w_face;
// Dst[16] = w_face;
// Dst[22] = w_face;
// // corners
// Dst[0] = w_corner;
// Dst[2] = w_corner;
// Dst[6] = w_corner;
// Dst[8] = w_corner;
// Dst[18] = w_corner;
// Dst[20] = w_corner;
// Dst[24] = w_corner;
// Dst[26] = w_corner;
// // edges
// Dst[1] = w_edge;
// Dst[3] = w_edge;
// Dst[5] = w_edge;
// Dst[7] = w_edge;
// Dst[9] = w_edge;
// Dst[11] = w_edge;
// Dst[15] = w_edge;
// Dst[17] = w_edge;
// Dst[19] = w_edge;
// Dst[21] = w_edge;
// Dst[23] = w_edge;
// Dst[25] = w_edge;
//
// for (int k=1; k<Nz-1; k++){
// for (int j=1; j<Ny-1; j++){
// for (int i=1; i<Nx-1; i++){
// int idx=Map(i,j,k);
// if (!(idx < 0)){
// double phi_x = 0.f;
// double phi_y = 0.f;
// double phi_z = 0.f;
// for (int kk=0; kk<3; kk++){
// for (int jj=0; jj<3; jj++){
// for (int ii=0; ii<3; ii++){
//
// int index = kk*9+jj*3+ii;
// double weight= Dst[index];
//
// int idi=i+ii-1;
// int idj=j+jj-1;
// int idk=k+kk-1;
//
// if (idi < 0) idi=0;
// if (idj < 0) idj=0;
// if (idk < 0) idk=0;
// if (!(idi < Nx)) idi=Nx-1;
// if (!(idj < Ny)) idj=Ny-1;
// if (!(idk < Nz)) idk=Nz-1;
//
// int nn = idk*Nx*Ny + idj*Nx + idi;
// double vec_x = double(ii-1);
// double vec_y = double(jj-1);
// double vec_z = double(kk-1);
// double GWNS=SolidPotential_host[nn];
// phi_x += GWNS*weight*vec_x;
// phi_y += GWNS*weight*vec_y;
// phi_z += GWNS*weight*vec_z;
// }
// }
// }
// if (Averages->SDs(i,j,k)<2.0){
// GreySolidGrad_host[idx+0*Np] = phi_x;
// GreySolidGrad_host[idx+1*Np] = phi_y;
// GreySolidGrad_host[idx+2*Np] = phi_z;
// }
// else{
// GreySolidGrad_host[idx+0*Np] = 0.0;
// GreySolidGrad_host[idx+1*Np] = 0.0;
// GreySolidGrad_host[idx+2*Np] = 0.0;
// }
// }
// }
// }
// }
//
//
// if (rank==0){
// printf("Number of Grey-solid labels: %lu \n",NLABELS);
// for (unsigned int idx=0; idx<NLABELS; idx++){
// VALUE=LabelList[idx];
// AFFINITY=AffinityList[idx];
// printf(" grey-solid label=%d, grey-solid affinity=%f\n",VALUE,AFFINITY);
// }
// }
//
//
// ScaLBL_CopyToDevice(GreySolidGrad, GreySolidGrad_host, 3*Np*sizeof(double));
// ScaLBL_Comm->Barrier();
// delete [] SolidPotential_host;
// delete [] GreySolidGrad_host;
// delete [] Dst;
//}
//--------- This is another old version of calculating greyscale-solid color-gradient modification-------//
// **not working effectively, to be deprecated
//void ScaLBL_GreyscaleColorModel::AssignGreySolidLabels()

10
models/GreyscaleColorModel.h
View File

@ -40,7 +40,7 @@ public:
double din,dout,inletA,inletB,outletA,outletB;
double GreyPorosity;
bool RecoloringOff;//recoloring can be turn off for grey nodes if this is true
double W;//wetting strength paramter for capillary pressure penalty for grey nodes
//double W;//wetting strength paramter for capillary pressure penalty for grey nodes
int Nx,Ny,Nz,N,Np;
int rank,nprocx,nprocy,nprocz,nprocs;
@ -50,7 +50,6 @@ public:
std::shared_ptr<Domain> Mask; // this domain is for lbm
std::shared_ptr<ScaLBL_Communicator> ScaLBL_Comm;
std::shared_ptr<ScaLBL_Communicator> ScaLBL_Comm_Regular;
std::shared_ptr<ScaLBL_Communicator> ScaLBL_Comm_Regular_2;
std::shared_ptr<GreyPhaseAnalysis> Averages;
// input database
@ -67,13 +66,14 @@ public:
double *fq, *Aq, *Bq;
double *Den, *Phi;
//double *GreySolidPhi; //Model 2 & 3
double *GreySolidGrad;//Model 1 & 4
//double *GreySolidGrad;//Model 1 & 4
double *GreySolidW;
//double *ColorGrad;
double *Velocity;
double *Pressure;
double *Porosity_dvc;
double *Permeability_dvc;
double *Psi;
//double *Psi;
private:
Utilities::MPI comm;
@ -90,7 +90,7 @@ private:
void AssignComponentLabels();
void AssignGreySolidLabels();
void AssignGreyPoroPermLabels();
void AssignGreyscalePotential();
//void AssignGreyscalePotential();
void ImageInit(std::string filename);
double MorphInit(const double beta, const double morph_delta);
double SeedPhaseField(const double seed_water_in_oil);

66
models/PoissonSolver.cpp
View File

@ -5,19 +5,35 @@
#include "analysis/distance.h"
#include "common/ReadMicroCT.h"
ScaLBL_Poisson::ScaLBL_Poisson(int RANK, int NP, const Utilities::MPI& COMM):
rank(RANK), nprocs(NP),timestep(0),timestepMax(0),tau(0),k2_inv(0),tolerance(0),h(0),
epsilon0(0),epsilon0_LB(0),epsilonR(0),epsilon_LB(0),Vin(0),Vout(0),Nx(0),Ny(0),Nz(0),N(0),Np(0),analysis_interval(0),
chargeDen_dummy(0),WriteLog(0),nprocx(0),nprocy(0),nprocz(0),
BoundaryConditionInlet(0),BoundaryConditionOutlet(0),BoundaryConditionSolid(0),Lx(0),Ly(0),Lz(0),
Vin0(0),freqIn(0),t0_In(0),Vin_Type(0),Vout0(0),freqOut(0),t0_Out(0),Vout_Type(0),
TestPeriodic(0),TestPeriodicTime(0),TestPeriodicTimeConv(0),TestPeriodicSaveInterval(0),
comm(COMM)
static inline bool fileExists( const std::string &filename )
{
std::ifstream ifile( filename.c_str() );
return ifile.good();
}
ScaLBL_Poisson::~ScaLBL_Poisson(){
ScaLBL_Poisson::ScaLBL_Poisson(int RANK, int NP, const Utilities::MPI& COMM):
rank(RANK), TIMELOG(nullptr), nprocs(NP),timestep(0),timestepMax(0),tau(0),k2_inv(0),tolerance(0),h(0),
epsilon0(0),epsilon0_LB(0),epsilonR(0),epsilon_LB(0),Vin(0),Vout(0),Nx(0),Ny(0),Nz(0),N(0),Np(0),analysis_interval(0),
chargeDen_dummy(0),WriteLog(0),nprocx(0),nprocy(0),nprocz(0),
BoundaryConditionInlet(0),BoundaryConditionOutlet(0),BoundaryConditionSolid(0),Lx(0),Ly(0),Lz(0),
Vin0(0),freqIn(0),t0_In(0),Vin_Type(0),Vout0(0),freqOut(0),t0_Out(0),Vout_Type(0),
TestPeriodic(0),TestPeriodicTime(0),TestPeriodicTimeConv(0),TestPeriodicSaveInterval(0),
comm(COMM)
{
if ( rank == 0 ) {
bool WriteHeader = !fileExists( "PoissonSolver_Convergence.csv" );
TIMELOG = fopen("PoissonSolver_Convergence.csv","a+");
if (WriteHeader)
fprintf(TIMELOG,"Timestep Error\n");
}
}
ScaLBL_Poisson::~ScaLBL_Poisson()
{
if ( TIMELOG )
fclose( TIMELOG );
}
void ScaLBL_Poisson::ReadParams(string filename){
@ -218,20 +234,19 @@ void ScaLBL_Poisson::ReadInput(){
void ScaLBL_Poisson::AssignSolidBoundary(double *poisson_solid)
{
size_t NLABELS=0;
signed char VALUE=0;
double AFFINITY=0.f;
auto LabelList = electric_db->getVector<int>( "SolidLabels" );
auto AffinityList = electric_db->getVector<double>( "SolidValues" );
NLABELS=LabelList.size();
size_t NLABELS = LabelList.size();
if (NLABELS != AffinityList.size()){
ERROR("Error: LB-Poisson Solver: SolidLabels and SolidValues must be the same length! \n");
}
double label_count[NLABELS];
double label_count_global[NLABELS];
std::vector<double> label_count( NLABELS, 0.0 );
std::vector<double> label_count_global( NLABELS, 0.0 );
// Assign the labels
for (size_t idx=0; idx<NLABELS; idx++) label_count[idx]=0;
@ -596,26 +611,11 @@ void ScaLBL_Poisson::Run(double *ChargeDensity, int timestep_from_Study){
}
void ScaLBL_Poisson::getConvergenceLog(int timestep,double error){
if (rank==0){
bool WriteHeader=false;
TIMELOG = fopen("PoissonSolver_Convergence.csv","r");
if (TIMELOG != NULL)
fclose(TIMELOG);
else
WriteHeader=true;
TIMELOG = fopen("PoissonSolver_Convergence.csv","a+");
if (WriteHeader)
{
fprintf(TIMELOG,"Timestep Error\n");
fprintf(TIMELOG,"%i %.5g\n",timestep,error);
fflush(TIMELOG);
}
else {
fprintf(TIMELOG,"%i %.5g\n",timestep,error);
fflush(TIMELOG);
}
void ScaLBL_Poisson::getConvergenceLog(int timestep,double error){
if ( rank == 0 ) {
fprintf(TIMELOG,"%i %.5g\n",timestep,error);
fflush(TIMELOG);
}
}

3
models/PoissonSolver.h
View File

@ -89,12 +89,13 @@ public:
private:
Utilities::MPI comm;
FILE *TIMELOG;
// filenames
char LocalRankString[8];
char LocalRankFilename[40];
char LocalRestartFile[40];
char OutputFilename[200];
FILE *TIMELOG;
//int rank,nprocs;
void LoadParams(std::shared_ptr<Database> db0);

34
sample_scripts/configure_lonestar Executable file
View File

@ -0,0 +1,34 @@
# load the module for cmake
#module load cmake
#source /gpfs/gpfs_stage1/b6p315aa/setup/setup-mpi.sh
module load cmake gcc
module load cuda
export HDF5_DIR=$HOME/local/hdf5/1.8.12/
export SILO_DIR=$HOME/local/silo/4.10.2/
export NETCDF_DIR=$HOME/local/netcdf/4.6.1
# configure
rm -rf CMake*
cmake \
-D CMAKE_BUILD_TYPE:STRING=Release \
-D CMAKE_C_COMPILER:PATH=mpicc \
-D CMAKE_CXX_COMPILER:PATH=mpicxx \
-D CMAKE_CXX_STANDARD=14 \
-D USE_CUDA=1 \
-D CMAKE_CUDA_FLAGS="-arch sm_70 -Xptxas=-v -Xptxas -dlcm=cg -lineinfo" \
-D CMAKE_CUDA_HOST_COMPILER="/opt/apps/gcc/7.3.0/bin/gcc" \
-D USE_HDF5=1 \
-D HDF5_DIRECTORY="$HDF5_DIR" \
-D HDF5_LIB="$HDF5_DIR/lib/libhdf5.a" \
-D USE_SILO=1 \
-D SILO_LIB="$SILO_DIR/lib/libsiloh5.a" \
-D SILO_DIRECTORY="$SILO_DIR" \
-D USE_NETCDF=0 \
-D NETCDF_DIRECTORY="$NETCDF_DIR" \
-D USE_DOXYGEN:BOOL=false \
-D USE_TIMER=0 \
~/src/LBPM
make VERBOSE=1 -j1 && make install

33
sample_scripts/configure_spock_TPLs Executable file
View File

@ -0,0 +1,33 @@
export TPL_ROOT=/ccs/home/mbt/repos
export TPL_BUILDER=/ccs/home/mbt/repos/TPL-builder
export TPL_WEBPAGE=http://bitbucket.org/AdvancedMultiPhysics/tpl-builder/downloads
export INSTALL_DIR=/ccs/home/mbt/spock/install
module load cmake
module load llvm-amdgpu
module load hip
cmake \
-D CMAKE_BUILD_TYPE=Release \
-D CXX_STD=14 \
-D TPL_INSTALL_DIR:PATH="${INSTALL_DIR}" \
-D C_COMPILER=cc \
-D CXX_COMPILER=CC \
-D Fortran_COMPILER=ftn \
-D FFLAGS="-fPIC" \
-D CFLAGS="-fPIC" \
-D CXXFLAGS="-fPIC" \
-D LDFLAGS="" \
-D ENABLE_STATIC:BOOL=ON \
-D ENABLE_SHARED:BOOL=OFF \
-D PROCS_INSTALL=8 \
-D TPL_LIST:STRING="TIMER;ZLIB;HDF5;SILO" \
-D TIMER_URL="${TPL_ROOT}/TimerUtility" \
-D ZLIB_URL="http://zlib.net/zlib-1.2.11.tar.gz" \
-D HDF5_URL="${TPL_ROOT}/hdf5-1.8.12.tar.gz" \
-D BUILD_TYPE=x86_64 \
-D SILO_URL="${TPL_ROOT}/silo-4.10.2.tar.gz" \
${TPL_BUILDER}

30
sample_scripts/configure_spock_hip Executable file
View File

@ -0,0 +1,30 @@
module load cmake
module load llvm-amdgpu
module load hip
export TPL_DIR=/gpfs/alpine/stf006/proj-shared/mbt/spock/install
# configure
rm -rf CMake*
cmake \
-D CMAKE_BUILD_TYPE:STRING=Release \
-D CMAKE_C_COMPILER:PATH=cc \
-D CMAKE_CXX_COMPILER:PATH=CC \
-D CMAKE_CXX_STANDARD=14 \
-D USE_HIP=1 \
-D LINK_LIBRARIES=${HIP_PATH}/lib/libamdhip64.so \
-D USE_CUDA=0 \
-D CMAKE_CUDA_FLAGS="-arch sm_70 -Xptxas=-v -Xptxas -dlcm=cg -lineinfo" \
-D CMAKE_CUDA_HOST_COMPILER="gcc" \
-D USE_MPI=0 \
-D USE_HDF5=1 \
-D HDF5_DIRECTORY="${TPL_DIR}/hdf5" \
-D USE_SILO=0 \
-D SILO_DIRECTORY="${TPL_DIR}/silo" \
-D USE_DOXYGEN:BOOL=false \
-D USE_TIMER=0 \
~/repos/LBPM-WIA