From 58ec1a27a1056fce5f1a482633f90435ad517a68 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Franck=20P=C3=A9rignon?= <franck.perignon@imag.fr> Date: Fri, 8 Jun 2012 15:08:55 +0000 Subject: [PATCH] Merge parmes and parmespy --- HySoP/src/CMakeLists.txt | 211 +++++++++++++++ .../particles/advec_common.f90 | 244 ------------------ 2 files changed, 211 insertions(+), 244 deletions(-) create mode 100644 HySoP/src/CMakeLists.txt diff --git a/HySoP/src/CMakeLists.txt b/HySoP/src/CMakeLists.txt new file mode 100644 index 000000000..099373346 --- /dev/null +++ b/HySoP/src/CMakeLists.txt @@ -0,0 +1,211 @@ +#======================================================= +# cmake utility to compile,link and install : +# - parmes library (libparmes...) +# - parmesscales library particular solver from scales, (libparmesscales...) +# +#======================================================= + +# --- Set a version number for the package --- +set(${PARMES_LIBRARY_NAME}_version 1.0.0) +# --- The name of the exe to be created (test purpose) +# This exe will be linked with libPARMES_LIBRARY_NAME +set(EXE_NAME ${PARMES_LIBRARY_NAME}Run) + +# The list of all dirs containing sources to be compiled for the Parmes lib +# Any file in those dirs will be used to create libparmes +set(${PARMES_LIBRARY_NAME}_SRCDIRS + . + # src/poisson + #src/interfaces/Fortran2Cpp + # src/interfaces/ppm + ) + +# A main file to create an executable (test purpose) +# Any files in these dirs will be used to create Parmes exec (linked with libparmes) +set(${EXE_NAME}_SRCDIRS main) +# Matching expr for files to be compiled. +set(EXTS *.f90 *.f95) +# Matching expr for headers (install purpose) +set(EXTS_HDRS *.hpp *.h) +# Note FP : we can also use cmake vars ${CMAKE_Fortran_SOURCE_FILE_EXTENSIONS} ${CMAKE_C_SOURCE_FILE_EXTENSIONS} ${CMAKE_CXX_SOURCE_FILE_EXTENSIONS} + +# ============= The project ============= +# Set project name and project languages +# => this automatically defines: +# - ${PARMES_LIBRARY_NAME}_BINARY_DIR : where you have run cmake, i.e. the place for compilation +# - ${PARMES_LIBRARY_NAME}_SOURCE_DIR : where sources (.f and .h and this CMakeLists.txt) are located +# Note that because of OutOfSourceBuild, binary_dir and source_dir must be different. +project(${PARMES_LIBRARY_NAME} Fortran) + +# ============= Search for libraries ============= +# We search for libraries Parmes depends on and +# set the compile/link conf (-I and -L opt) + +# --- Boost --- +#set(Boost_ADDITIONAL_VERSIONS "1.49" "1.49.0") +#set(Boost_NO_SYSTEM_PATHS ON) +#find_package(Boost REQUIRED) +#include_directories(${Boost_INCLUDE_DIRS}) + +# ============= Prepare compilation ============= + +# Force a default build type if not provided by user +# CMAKE_BUILD_TYPE = empty, Debug, Release, RelWithDebInfo or MinSizeRel. +if (NOT CMAKE_BUILD_TYPE) + set (CMAKE_BUILD_TYPE RELEASE CACHE STRING "Choose the type of build, options are: None, Debug, Release, RelWithDebInfo or MinSizeRel." FORCE) +endif (NOT CMAKE_BUILD_TYPE) + +# If the project uses Fortran ... +# Set module files directory (i.e. where .mod will be created) +set(CMAKE_Fortran_MODULE_DIRECTORY ${CMAKE_BINARY_DIR}/Modules) +# Add compilation flags: +append_Fortran_FLAGS("-Wall") + +# --- MPI --- +if(USE_MPI) + # Find MPI for C++ and fortran. + find_package(MPI REQUIRED) + # -I + include_directories(${MPI_Fortran_INCLUDE_PATH}) + # Add compilation flags + set(${PARMES_LIBRARY_NAME}_LINK_FLAGS ${${PARMES_LIBRARY_NAME}_LINK_FLAGS} ${MPI_Fortran_LINK_FLAGS}) + # -I + # Add compilation flags + append_Fortran_flags(${MPI_Fortran_COMPILE_FLAGS}) + #set(${PARMES_LIBRARY_NAME}_LINK_FLAGS ${${PARMES_LIBRARY_NAME}_LINK_FLAGS} ${MPI_Fortran_LINK_FLAGS}) + set(LIBS ${LIBS} ${MPI_Fortran_LIBRARIES} ) +endif(USE_MPI) + +# --- scales --- +if(WITH_SCALES) + add_subdirectory(scalesInterface) +endif() + +# --- PPM --- +if(WITH_PPM) + add_subdirectory(ppmInterface) +endif() + +find_package(FFTW REQUIRED) +include_directories(${FFTW_INCLUDE_DIRS}) +set(LIBS ${LIBS} ${FFTW_LIBRARIES}) +display(FFTW_LIBRARIES) + +# ============= Generates ParmesConfig.hpp ============= +# The file PARMES_LIBRARY_NAME_defines.hpp will be generated from config.hpp.cmake; +if(EXISTS ${CMAKE_SOURCE_DIR}/config.hpp.cmake) + if(NOT CONFIG_H_GLOBAL_CONFIGURED) + set(CONFIG_H_GLOBAL_CONFIGURED 1 CACHE BOOL "Parmes_defines.hpp generation." ) + configure_file(${CMAKE_SOURCE_DIR}/config.hpp.cmake ${PARMES_LIBRARY_NAME}_defines.hpp) + endif() + include_directories(${CMAKE_BINARY_DIR}) +endif() + +# ============= Source and header files list ============= +# We scan all files with matching extension in directories +# containing sources. + +# Source and header files list: +foreach(_DIR ${${PARMES_LIBRARY_NAME}_SRCDIRS}) + set(_DIR_FILES) + foreach(_EXT ${EXTS}) # Source files + file(GLOB _DIR_FILES_EXT ${_DIR}/${_EXT}) + if(_DIR_FILES_EXT) + list(APPEND ${PARMES_LIBRARY_NAME}_SRC ${_DIR_FILES_EXT}) + endif() + endforeach() + foreach(_EXT ${EXTS_HDRS}) # Headers + file(GLOB _DIR_FILES_EXT ${_DIR}/${_EXT}) + if(_DIR_FILES_EXT) + list(APPEND ${PARMES_LIBRARY_NAME}_HDRS ${_DIR_FILES_EXT}) + endif() + endforeach() +endforeach() +# We add headers to source files +list(APPEND ${PARMES_LIBRARY_NAME}_SRC ${${PARMES_LIBRARY_NAME}_HDRS}) + +# The same for main dir ... +foreach(_DIR ${${EXE_NAME}_SRCDIRS}) + set(_DIR_FILES) + foreach(_EXT ${EXTS}) + file(GLOB _DIR_FILES_EXT ${_DIR}/${_EXT}) + if(_DIR_FILES_EXT) + list(APPEND ${EXE_NAME}_SRC ${_DIR_FILES_EXT}) + endif() + endforeach() + foreach(_EXT ${EXTS_HDRS}) + file(GLOB _DIR_FILES_EXT ${_DIR}/${_EXT}) + if(_DIR_FILES_EXT) + list(APPEND ${EXE_NAME}_HDRS ${_DIR_FILES_EXT}) + endif() + endforeach() +endforeach() +list(APPEND ${EXE_NAME}_SRC ${${EXE_NAME}_HDRS}) + +# Add directories to those searched by compiler ... +# -I +include_directories(${${PARMES_LIBRARY_NAME}_SRCDIRS}) +include_directories(${${EXE_NAME}_HDRS}) +include_directories(${CMAKE_Fortran_MODULE_DIRECTORY}) + +# Cmake tools to handle fortran-c interface. It will generate F2CMangle.hpp, a file +# that will contain some useful macros to deal with fortran-C name mangling. +# See cmake doc for this module. +#include(FortranCInterface) +#FortranCInterface_HEADER(${CMAKE_Fortran_MODULE_DIRECTORY}/F2CMangle.hpp +# MACRO_NAMESPACE "F2C_" +# SYMBOL_NAMESPACE "F2C_") + +# ============= Creates the library ============= +if(BUILD_SHARED_LIBS) # shared library + add_library(${PARMES_LIBRARY_NAME} SHARED ${${PARMES_LIBRARY_NAME}_SRC}) +else() # static library + add_library(${PARMES_LIBRARY_NAME} STATIC ${${PARMES_LIBRARY_NAME}_SRC}) +endif() +# Libs to link with PARMES__LIBRARY_NAME +target_link_libraries(${PARMES_LIBRARY_NAME} ${LIBS}) + +# ============= Creates the executable ============= +add_executable(${EXE_NAME} ${${EXE_NAME}_SRC}) +add_dependencies(${EXE_NAME} ${PARMES_LIBRARY_NAME}) + +# libs to link with EXE_NAME +target_link_libraries(${EXE_NAME} ${PARMES_LIBRARY_NAME}) +target_link_libraries(${EXE_NAME} ${LIBS}) + +# ============== Add tests ============== +if(WITH_TESTS) + message(STATUS "Enable testing ...") + begin_test(src/tests/F2003) + new_test(testAllocatedPtr userMod.f90 wrapper.f90 testAllocatedPtr.cxx) + new_test(testNullPtr userMod.f90 wrapper.f90 testNullPtr.cxx) + end_test() + begin_test(src/tests/Particles) + new_test(testCreate3D testCreate3D.f90) + new_test(testCreate2D testCreate2D.f90) + end_test() +endif(WITH_TESTS) + +# ============= Prepare install ============= + +# The library +# The library, the headers and mod files, the cmake generated files +# will be install in CMAKE_INSTALL_PREFIX/lib include and share +include(InstallPackage) +install_package(${PACKAGE_NAME} ${PARMES_LIBRARY_NAME} ${${PARMES_LIBRARY_NAME}_HDRS}) + +#install(TARGETS ${EXE_NAME} +#RUNTIME DESTINATION bin # executables +# ) + + +#configure_file(${PACKAGE_NAME}Config.cmake.in +# "${PARMES_LIBRARY_BINARY_DIR}/${PACKAGE_NAME}Config.cmake") +#configure_file(${PACKAGE_NAME}ConfigVersion.cmake.in +# "${PARMES_LIBRARY_BINARY_DIR}/${PACKAGE_NAME}ConfigVersion.cmake" @ONLY) +#install(FILES +# "${PARMES_LIBRARY_BINARY_DIR}/InstallFiles/${_PACK}Config.cmake" +# "${PARMES_LIBRARY_BINARY_DIR}/InstallFiles/${_PACK}ConfigVersion.cmake" +# DESTINATION "${${PACKAGE_NAME}_CMAKE_DIR}" COMPONENT dev) + + diff --git a/HySoP/src/scalesInterface/particles/advec_common.f90 b/HySoP/src/scalesInterface/particles/advec_common.f90 index 48433ced4..0368f38a8 100644 --- a/HySoP/src/scalesInterface/particles/advec_common.f90 +++ b/HySoP/src/scalesInterface/particles/advec_common.f90 @@ -91,250 +91,6 @@ contains ! ==================== Compute particle velocity (RK2) ==================== ! ================================================================================== - !> Interpolate the velocity field used in a RK2 scheme for particle advection. - !! @param[in] dt = time step - !! @param[in] direction = current direction (1 = along X, 2 = along Y and 3 = along Z) - !! @param[in] ind_group = coordinate of the current group of lines - !! @param[in] p_pos_adim = adimensionned particle postion - !! @param[in,out] p_V = particle velocity (along the current direction) - !! @details - !! A RK2 scheme is used to advect the particles : the midlle point scheme. An - !! intermediary position "p_pos_bis(i) = p_pos(i) + V(i)*dt/2" is computed and then - !! the numerical velocity of each particles is computed as the interpolation of V in - !! this point. This field is used to advect the particles at the seconde order in time : - !! p_pos(t+dt, i) = p_pos(i) + p_V(i). - !! The group line indice is used to ensure using unicity of each mpi message tag. - subroutine AC_particle_velocity_line(dt, direction, ind_group, p_pos_adim, p_V) - - ! This code involve a recopy of p_V. It is possible to directly use the 3D velocity field but in a such code - ! a memory copy is still needed to send velocity field to other processus : mpi send contiguous memory values - - ! Input/Ouput - real(WP), intent(in) :: dt ! time step - integer, intent(in) :: direction - integer, dimension(2), intent(in) :: ind_group - real(WP), dimension(:), intent(in) :: p_pos_adim - real(WP), dimension(:), intent(inout) :: p_V - ! Others, local - real(WP), dimension(N_proc(direction)) :: p_pos_bis ! adimensionned position of the middle point - real(WP), dimension(N_proc(direction)), target :: p_V_bis ! velocity of the middle point - real(WP), dimension(N_proc(direction)) :: weight ! interpolation weight - type(real_pter), dimension(N_proc(direction)) :: Vp, Vm ! Velocity on previous and next mesh point - real(WP), dimension(:), allocatable, target :: V_buffer ! Velocity buffer for postion outside of the local subdomain - integer :: size_buffer ! buffer size - integer :: rece_ind_min ! the minimal indice used in velocity interpolation - integer :: rece_ind_max ! the maximal indice used in velocity interpolation - integer :: ind, ind_com ! indices - integer :: pos, pos_old ! indices of the mesh point wich preceed the particle position - integer :: proc_gap, gap! distance between my (mpi) coordonate and coordinate of the - ! processus associated to a given position - integer, dimension(:), allocatable :: rece_rank ! rank of processus wich send me information - integer :: send_rank ! rank of processus to wich I send information - integer :: rankP ! rank of processus ("source rank" returned by mpi_cart_shift) - integer, dimension(2) :: rece_range ! range of the velocity fields I want to receive - integer, dimension(2) :: send_range ! range of the velocity fields I send - integer, dimension(2) :: rece_gap ! distance between me and processus wich send me information - integer, dimension(2) :: send_gap ! distance between me and processus to wich I send information - integer :: msg_size ! size of message send/receive - integer :: tag ! mpi message tag - integer :: ierr ! mpi error code - integer, dimension(:), allocatable :: s_request ! mpi communication request (handle) of nonblocking send - integer, dimension(:), allocatable :: s_request_bis! mpi communication request (handle) of nonblocking send - integer, dimension(:), allocatable :: rece_request ! mpi communication request (handle) of nonblocking receive - integer, dimension(MPI_STATUS_SIZE) :: rece_status ! mpi status (for mpi_wait) - - ! -- Initialisation -- - ind_com = 0 - do ind = 1, N_proc(direction) - nullify(Vp(ind)%pter) - nullify(Vm(ind)%pter) - end do - ! Compute the midlle point - p_pos_bis = p_pos_adim + (dt/2.0)*p_V/d_sc(direction) - p_V_bis = p_V - ! Compute range of the set of point where I need the velocity value - rece_ind_min = floor(p_pos_bis(1)) - rece_ind_max = floor(p_pos_bis(N_proc(direction))) + 1 - ! Allocate the buffer - ! If rece_ind_min and rece_ind_max are not in [N_proc(direction);1] then it will change the number of communication - ! size_buffer = max(temp - N_proc(direction), 0) - min(0, temp) - !size_buffer = - max(temp - N_proc(direction), 0) - min(0, temp) - ! It must work, but for first test we prefer compute size_buffer more simply - size_buffer = 0 - - ! -- Exchange non blocking message to do the computations during the - ! communication process -- - call AC_obtain_receivers(direction, ind_group, rece_ind_min, rece_ind_max, send_gap, rece_gap) - allocate(rece_rank(rece_gap(1):rece_gap(2))) - ! Send messages about what I want - allocate(s_request_bis(rece_gap(1):rece_gap(2))) - do proc_gap = rece_gap(1), rece_gap(2) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, rece_rank(proc_gap), ierr) - if (rece_rank(proc_gap) /= D_rank(direction)) then - ! Range I want - gap = proc_gap*N_proc(direction) - rece_range(1) = max(rece_ind_min, gap+1) ! fortran => indice start from 0 - rece_range(2) = min(rece_ind_max, gap+N_proc(direction)) - ! Tag = concatenation of (rank+1), ind_group(1), ind_group(2), direction et unique Id. - tag = compute_tag(ind_group, tag_velo_range, direction, proc_gap) - ! Send message - size_buffer = size_buffer + (rece_range(2)-rece_range(1)) + 1 - call mpi_Isend(rece_range(1), 2, MPI_INTEGER, rece_rank(proc_gap), tag, D_comm(direction), s_request_bis(proc_gap),ierr) - end if - end do - allocate(V_buffer(max(size_buffer,1))) - V_buffer = 0 - ! Send the velocity field to processus which need it - allocate(s_request(send_gap(1):send_gap(2))) - do proc_gap = send_gap(1), send_gap(2) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank /= D_rank(direction)) then - ! I - Receive messages about what I have to send - ! Ia - Compute reception tag = concatenation of (rank+1), ind_group(1), ind_group(2), direction et unique Id. - tag = compute_tag(ind_group, tag_velo_range, direction, -proc_gap) - ! Ib - Receive the message - call mpi_recv(send_range(1), 2, MPI_INTEGER, send_rank, tag, D_comm(direction), rece_status, ierr) - send_range = send_range + proc_gap*N_proc(direction) - ! II - Send it - ! IIa - Compute send tag - tag = compute_tag(ind_group, tag_velo_V, direction, proc_gap) - ! IIb - Send message - call mpi_Isend(p_V(send_range(1)), send_range(2)-send_range(1)+1, MPI_DOUBLE_PRECISION, & - & send_rank, tag, D_comm(direction), s_request(proc_gap), ierr) - end if - end do - - ! Non blocking reception of the velocity field - ind = 1 - allocate(rece_request(rece_gap(1):rece_gap(2))) - do proc_gap = rece_gap(1), rece_gap(2) - if (rece_rank(proc_gap) /= D_rank(direction)) then - ! IIa - Compute reception tag - tag = compute_tag(ind_group, tag_velo_V, direction, -proc_gap) - ! IIb - Receive message - gap = proc_gap*N_proc(direction) - rece_range(1) = max(rece_ind_min, gap+1) ! fortran => indice start from 0 - rece_range(2) = min(rece_ind_max, gap+N_proc(direction)) - msg_size = rece_range(2)-rece_range(1)+1 - call mpi_Irecv(V_buffer(ind), msg_size, MPI_DOUBLE_PRECISION, rece_rank(proc_gap), tag, D_comm(direction), & - & rece_request(proc_gap), ierr) - ind = ind + msg_size - end if - end do - - ! -- Compute the interpolated velocity - ! Compute the interpolation weight and update the pointers Vp and Vm - ! Initialisation of reccurence process - ind = 1 - pos = floor(p_pos_bis(ind)) - weight(ind) = p_pos_bis(ind)-pos - ! Vm = V(pos) - proc_gap = floor(real(pos-1)/N_proc(direction)) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank == D_rank(direction)) then - Vm(ind)%pter => p_V_bis(pos-proc_gap*N_proc(direction)) - else - ind_com = ind_com + 1 - Vm(ind)%pter => V_buffer(ind_com) - end if - ! Vp = V(pos+1) - proc_gap = floor(real(pos+1-1)/N_proc(direction)) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank == D_rank(direction)) then - Vp(ind)%pter => p_V_bis(pos+1-proc_gap*N_proc(direction)) - else - ind_com = ind_com + 1 - Vp(ind)%pter => V_buffer(ind_com) - end if - pos_old = pos - - ! Following indice : we use previous work (already done) - do ind = 2, N_proc(direction) - pos = floor(p_pos_bis(ind)) - weight(ind) = p_pos_bis(ind)-pos - select case(pos-pos_old) - case(0) - ! The particle belongs to the same segment than the previous one - Vm(ind)%pter => Vm(ind-1)%pter - Vp(ind)%pter => Vp(ind-1)%pter - case(1) - ! The particle follows the previous one - Vm(ind)%pter => Vp(ind-1)%pter - ! Vp = V(pos+1) - proc_gap = floor(real(pos+1-1)/N_proc(direction)) ! fortran -> indice starts from 1 - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank == D_rank(direction)) then - Vp(ind)%pter => p_V_bis(pos+1-proc_gap*N_proc(direction)) - else - ind_com = ind_com + 1 - Vp(ind)%pter => V_buffer(ind_com) - end if - case(2) - ! pos = pos_old +2, wich correspond to "extention" - ! Vm = V(pos) - proc_gap = floor(real(pos-1)/N_proc(direction)) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank == D_rank(direction)) then - Vm(ind)%pter => p_V_bis(pos-proc_gap*N_proc(direction)) - else - ind_com = ind_com + 1 - Vm(ind)%pter => V_buffer(ind_com) - end if - ! Vp = V(pos+1) - proc_gap = floor(real(pos+1-1)/N_proc(direction)) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank == D_rank(direction)) then - Vp(ind)%pter => p_V_bis(pos+1-proc_gap*N_proc(direction)) - else - ind_com = ind_com + 1 - Vp(ind)%pter => V_buffer(ind_com) - end if - case default - print*, "unexpected case : pos = ", pos, " , pos_old = ", pos_old, " ind = ", ind - end select - pos_old = pos - end do - - ! -- Compute the interpolate velocity -- - ! Check if communication are done - do proc_gap = rece_gap(1), rece_gap(2) - if (rece_rank(proc_gap)/=D_rank(direction)) then - call mpi_wait(rece_request(proc_gap), rece_status, ierr) - end if - end do - - - ! Then compute the field - do ind = 1, N_proc(direction) - p_V(ind) = weight(ind)*Vp(ind)%pter + (1-weight(ind))*Vm(ind)%pter - end do - - do ind = 1, N_proc(direction) - nullify(Vp(ind)%pter) - nullify(Vm(ind)%pter) - end do - - do proc_gap = send_gap(1), send_gap(2) - call mpi_cart_shift(D_comm(direction), 0, proc_gap, rankP, send_rank, ierr) - if (send_rank /= D_rank(direction)) then - call MPI_WAIT(s_request(proc_gap),rece_status,ierr) - end if - end do - deallocate(s_request) - do proc_gap = rece_gap(1), rece_gap(2) - if (rece_rank(proc_gap) /= D_rank(direction)) then - call MPI_WAIT(s_request_bis(proc_gap),rece_status,ierr) - end if - end do - deallocate(s_request_bis) - - ! Deallocation - deallocate(rece_rank) - deallocate(rece_request) - deallocate(V_buffer) - - end subroutine AC_particle_velocity_line - !> Interpolate the velocity field used in a RK2 scheme for particle advection - !! version for a group of (more of one) line -- GitLab