Skip to content

parallel particles with contact model

OpenLB – Open Source Lattice Boltzmann Code Forums on OpenLB General Topics parallel particles with contact model

Viewing 2 posts - 1 through 2 (of 2 total)
  • Author
    Posts
  • #9384
    jakschee
    Participant

    Hi openLB community,

    I would like to set up a simulation with 4-way coupling in parallel (basically mixing the examples settlingCube3d and dkt2d as a first step). For this I am trying to follow the time stepping scheme from [1] (-> Algorithm 1). However, I am not sure how to implement the necessary communication steps. My current time-stepping loop looks like this:

    
      for (std::size_t iT = 0; iT < converter.getLatticeTime(maxPhysT)+10; ++iT) {
    
        // couple fluid to particle
        particleManager.execute<couple_lattice_to_particles<T,DESCRIPTOR,PARTICLETYPE>>();
    
        particleManager.execute<communicate_parallel_surface_force<T,PARTICLETYPE>>();
        
        // apply gravity
        particleManager.execute<apply_gravity<T,PARTICLETYPE>>();
    
        for (std::size_t iT2=0; iT2 < 100; ++iT2) { //sub time stepping scheme for particle contact
          
          // Calculate and apply contact forces
          processContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE, ContactProperties<T, 1>>(
            particleSystem, solidBoundaries, contactContainer, contactProperties,
            superGeometry, particleManager.getParticleCommunicator().contactTreatmentComm, contactBoxResolutionPerDirection);
    
          // communicate contact forces?
          
    
          // Solve equations of motion
          particleManager.execute<process_dynamics_parallel<T,PARTICLETYPE>>(converter.getConversionFactorTime()/100.);
    
     
          //communicatePostContactTreatmentContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer, particleSystem, converter.getConversionFactorLength(),
          //												   particleManager.getParticleCommunicator().???,
          //											   particleManager.getParticleCommunicator().???,
          //											   {false,false,false});
    
          particleManager.execute<update_particle_core_distribution<T, PARTICLETYPE>>();
          
          // Couple particles to lattice (with contact detection)
          coupleResolvedParticlesToLattice<T, DESCRIPTOR, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(
            particleSystem, contactContainer, superGeometry, sLattice, converter, solidBoundaries);
    
          communicateContacts<T, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer);
        
        }
        
        // Get Results
        getResults(sLattice, converter, iT, superGeometry, timer, particleSystem);
    
        // Collide and stream
        sLattice.collideAndStream();
      }
    

    What do I need to change here to get my desired 4-way coupling while still working in parallel? Any help is much appreciated!

    Best regards
    Jakob

    PS: In addition to changing the loop, I also set up a new PARTICLETYPE to do contact calculations as well as being parallelizable. I hope that this step was okay:

    
    using ResolvedDecomposedParticleWithContact3D =
          PARTICLE_DESCRIPTOR<3, GENERAL_EXTENDABLE<3>,
          MOBILITY_VERLET<3>, SURFACE_RESOLVED_PARALLEL<3>,
          FORCING_RESOLVED<3>, PHYSPROPERTIES_RESOLVED<3>,
          MECHPROPERTIES_COLLISION<3>, NUMERICPROPERTIES_RESOLVED_CONTACT<3>,
          PARALLELIZATION_RESOLVED >;
    

    [1] https://doi.org/10.1016/j.cpc.2024.109321

    #9385
    jakschee
    Participant

    Quick update-> I changed the loop to:

    
      for (std::size_t iT = 0; iT < converter.getLatticeTime(maxPhysT)+10; ++iT) {
    
        // Execute particle manager
        particleManager.execute<couple_lattice_to_particles<T,DESCRIPTOR,PARTICLETYPE>>();
    
        particleManager.execute<communicate_parallel_surface_force<T,PARTICLETYPE>>();
        
        // apply gravity
        particleManager.execute<apply_gravity<T,PARTICLETYPE>>();
    
        for (std::size_t iT2=0; iT2 < 1; ++iT2) { //sub time stepping scheme for particle contact
          
          // Calculate and apply contact forces
          processContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE, ContactProperties<T, 1>>(
            particleSystem, solidBoundaries, contactContainer, contactProperties,
            superGeometry, particleManager.getParticleCommunicator().contactTreatmentComm, contactBoxResolutionPerDirection);
          
          // Solve equations of motion
          particleManager.execute<process_dynamics_parallel<T,PARTICLETYPE>>(converter.getConversionFactorTime()/1.);
    
          
          								     
          communicatePostContactTreatmentContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer, particleSystem, converter.getConversionFactorLength(),
          											   particleManager.getParticleCommunicator().equationsOfMotionComm,
          											   particleManager.getParticleCommunicator().equationsOfMotionComm,
          											   {false,false,false});
    
          particleManager.execute<update_particle_core_distribution<T, PARTICLETYPE>>();
          
          // Couple particles to lattice (with contact detection)
          coupleResolvedParticlesToLattice<T, DESCRIPTOR, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(
            particleSystem, contactContainer, superGeometry, sLattice, converter, solidBoundaries);
    
          communicateContactsParallel<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer,particleSystem,
    											 converter.getConversionFactorLength(),
    									   particleManager.getParticleCommunicator().particleContactDetectionComm,
    									   particleManager.getParticleCommunicator().wallContactDetectionComm,
    									   {false,false,false});
        
        }
    
        
        
        // Get Results
        getResults(sLattice, converter, iT, superGeometry, timer, particleSystem);
    
        // Collide and stream
        sLattice.collideAndStream();
      }
    

    Which results in the following output during compilation:

    
    make -C ../../../external
    make[1]: Entering directory '/home/jakschee/OpenLB/external'
    make -C zlib
    make[2]: Entering directory '/home/jakschee/OpenLB/external/zlib'
    make[2]: Nothing to be done for 'all'.
    make[2]: Leaving directory '/home/jakschee/OpenLB/external/zlib'
    cp zlib/build/libz.a lib/
    make -C tinyxml
    make[2]: Entering directory '/home/jakschee/OpenLB/external/tinyxml'
    make[2]: Nothing to be done for 'all'.
    make[2]: Leaving directory '/home/jakschee/OpenLB/external/tinyxml'
    cp tinyxml/build/libtinyxml.a lib/
    make[1]: Leaving directory '/home/jakschee/OpenLB/external'
    make -C ../../.. core
    make[1]: Entering directory '/home/jakschee/OpenLB'
    make[1]: Nothing to be done for 'core'.
    make[1]: Leaving directory '/home/jakschee/OpenLB'
    mpic++ -O3 -Wall -march=native -mtune=native -std=c++17 -pthread -DPARALLEL_MODE_MPI  -DPLATFORM_CPU_SISD -DPLATFORM_CPU_SIMD -DDEFAULT_FLOATING_POINT_TYPE=double -I../../../src -I../../../external/zlib -I../../../external/tinyxml -c -o dkt2d.o dkt2d.cpp
    In file included from /usr/include/c++/11/bits/char_traits.h:39,
                     from /usr/include/c++/11/string:40,
                     from /usr/include/c++/11/stdexcept:39,
                     from ../../../src/core/platform/platform.h:28,
                     from ../../../src/core/core.h:24,
                     from ../../../src/core/core2D.h:28,
                     from ../../../src/olb2D.h:1,
                     from dkt2d.cpp:51:
    In static member function ‘static _Tp* std::__copy_move<_IsMove, true, std::random_access_iterator_tag>::__copy_m(const _Tp*, const _Tp*, _Tp*) [with _Tp = bool; bool _IsMove = false]’,
        inlined from ‘_OI std::__copy_move_a2(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:495:30,
        inlined from ‘_OI std::__copy_move_a1(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:522:42,
        inlined from ‘_OI std::__copy_move_a(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:529:31,
        inlined from ‘_OI std::copy(_II, _II, _OI) [with _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:620:7,
        inlined from ‘constexpr olb::Vector<T, D>::Vector(std::initializer_list<_Tp>) [with T = bool; unsigned int D = 2]’ at ../../../src/core/vector.h:116:14,
        inlined from ‘int main(int, char**)’ at dkt2d.cpp:373:101:
    /usr/include/c++/11/bits/stl_algobase.h:431:30: warning: ‘void* __builtin_memcpy(void*, const void*, long unsigned int)’ writing 3 bytes into a region of size 2 overflows the destination [-Wstringop-overflow=]
      431 |             __builtin_memmove(__result, __first, sizeof(_Tp) * _Num);
          |             ~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    dkt2d.cpp: In function ‘int main(int, char**)’:
    dkt2d.cpp:373:101: note: destination object ‘<anonymous>’ of size 2
      373 |       communicatePostContactTreatmentContacts<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer, particleSystem, converter.getConversionFactorLength(),
          |       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      374 |                                                                                            particleManager.getParticleCommunicator().equationsOfMotionComm,
          |                                                                                            ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      375 |                                                                                            particleManager.getParticleCommunicator().equationsOfMotionComm,
          |                                                                                            ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      376 |                                                                                            {false,false,false});
          |                                                                                            ~~~~~~~~~~~~~~~~~~~~
    In file included from /usr/include/c++/11/bits/char_traits.h:39,
                     from /usr/include/c++/11/string:40,
                     from /usr/include/c++/11/stdexcept:39,
                     from ../../../src/core/platform/platform.h:28,
                     from ../../../src/core/core.h:24,
                     from ../../../src/core/core2D.h:28,
                     from ../../../src/olb2D.h:1,
                     from dkt2d.cpp:51:
    In static member function ‘static _Tp* std::__copy_move<_IsMove, true, std::random_access_iterator_tag>::__copy_m(const _Tp*, const _Tp*, _Tp*) [with _Tp = bool; bool _IsMove = false]’,
        inlined from ‘_OI std::__copy_move_a2(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:495:30,
        inlined from ‘_OI std::__copy_move_a1(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:522:42,
        inlined from ‘_OI std::__copy_move_a(_II, _II, _OI) [with bool _IsMove = false; _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:529:31,
        inlined from ‘_OI std::copy(_II, _II, _OI) [with _II = const bool*; _OI = bool*]’ at /usr/include/c++/11/bits/stl_algobase.h:620:7,
        inlined from ‘constexpr olb::Vector<T, D>::Vector(std::initializer_list<_Tp>) [with T = bool; unsigned int D = 2]’ at ../../../src/core/vector.h:116:14,
        inlined from ‘int main(int, char**)’ at dkt2d.cpp:384:89:
    /usr/include/c++/11/bits/stl_algobase.h:431:30: warning: ‘void* __builtin_memcpy(void*, const void*, long unsigned int)’ writing 3 bytes into a region of size 2 overflows the destination [-Wstringop-overflow=]
      431 |             __builtin_memmove(__result, __first, sizeof(_Tp) * _Num);
          |             ~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    dkt2d.cpp: In function ‘int main(int, char**)’:
    dkt2d.cpp:384:89: note: destination object ‘<anonymous>’ of size 2
      384 |       communicateContactsParallel<T, PARTICLETYPE, PARTICLECONTACTTYPE, WALLCONTACTTYPE>(contactContainer,particleSystem,
          |       ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      385 |                                                                                          converter.getConversionFactorLength(),
          |                                                                                          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      386 |                                                                            particleManager.getParticleCommunicator().particleContactDetectionComm,
          |                                                                            ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      387 |                                                                            particleManager.getParticleCommunicator().wallContactDetectionComm,
          |                                                                            ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
      388 |                                                                            {false,false,false});
          |                                                                            ~~~~~~~~~~~~~~~~~~~~
    mpic++ dkt2d.o -o dkt2d -lolbcore -L../../../external/lib  -lpthread -lrt -lz -ltinyxml -L../../../build/lib
    

    However, when I run the executable all contacts (particle-particle and particle-wall) seem to get resolved correctly…

Viewing 2 posts - 1 through 2 (of 2 total)
  • You must be logged in to reply to this topic.