! particles.f90 -- ! Example of modelling particles transported in a flow field ! ! Example belonging to "Modern Fortran in Practice" by Arjen Markus ! ! This work is licensed under the Creative Commons Attribution 3.0 Unported License. ! To view a copy of this license, visit http://creativecommons.org/licenses/by/3.0/ ! or send a letter to: ! Creative Commons, 444 Castro Street, Suite 900, Mountain View, California, 94041, USA. ! ! particle_modelling -- ! Basic module providing the framework ! module particle_modelling use points2d3d implicit none type, extends(point2d) :: basic_particle_type real :: mass contains procedure :: force => force_basic_particle procedure, pass(particle) & :: new_position => position_basic_particle end type basic_particle_type type :: basic_flow_field_type contains procedure :: flow_velocity => velocity_basic_flow_field end type basic_flow_field_type contains subroutine position_basic_particle( particle, flow_field, deltt ) class(basic_particle_type), intent(inout), target :: particle class(basic_flow_field_type), intent(in) :: flow_field real, intent(in) :: deltt ! Empty routine end subroutine position_basic_particle subroutine force_basic_particle( particle, force_vector ) class(basic_particle_type), intent(inout) :: particle class(point2d), intent(out) :: force_vector ! Empty routine end subroutine force_basic_particle subroutine velocity_basic_flow_field( flow_field, position, velocity ) class(basic_flow_field_type), intent(in) :: flow_field class(point2d), intent(in) :: position class(point2d), intent(out) :: velocity ! Empty routine end subroutine velocity_basic_flow_field end module particle_modelling ! analytical_field_module -- ! Implementation of a simple flow field ! module analytical_field_module use particle_modelling implicit none type, extends(basic_flow_field_type) :: analytical_field real :: x_velocity real :: y_velocity contains procedure :: flow_velocity => velocity_analytical_flow_field end type analytical_field contains subroutine velocity_analytical_flow_field( flow_field, position, velocity ) class(analytical_field), intent(in) :: flow_field class(point2d), intent(in) :: position class(point2d), intent(out) :: velocity ! The flow field is uniform, so always return the same vector velocity%x = flow_field%x_velocity velocity%y = flow_field%y_velocity end subroutine velocity_analytical_flow_field end module analytical_field_module ! oil_particles -- ! Naive implementation of an oil model ! ! Note: ! We do not use the force routine in this implementation, ! so we can use the dummy version from the basic particle type. ! module oil_particles use particle_modelling implicit none type, extends(basic_particle_type) :: oil_particle logical :: stuck contains procedure, pass(particle) & :: new_position => position_oil_particle end type oil_particle contains subroutine position_oil_particle( particle, flow_field, deltt ) class(oil_particle), target :: particle real :: deltt class(basic_flow_field_type) :: flow_field class(point2d), pointer :: position type(point2d) :: velocity type(point2d) :: random_displacement real :: r ! ! The particle may get stuck to the bottom ... ! call random_number( r ) if ( r > 0.99 ) then particle%stuck = .true. endif ! ! If it is stuck to the bottom, no further motion ! if ( particle%stuck ) then return endif ! ! Else let it be transported with the flow. Add a ! random displacement due to mixing and turbulence ! position => particle%point2d call flow_field%flow_velocity( position, velocity ) call random_displacement%random_vector position = position + deltt * velocity + & random_displacement end subroutine position_oil_particle end module oil_particles ! test_oil_particles -- ! Main program to demonstrate the oil particles modelling ! program test_oil_particles use oil_particles use analytical_field_module implicit none ! type(oil_particle), dimension(100000) :: particle type(oil_particle), dimension(10) :: particle type(analytical_field) :: flow_field real :: deltt integer :: time integer :: number_times integer :: p flow_field = analytical_field( 1.0, 0.1 ) number_times = 10 do time = 1,number_times write(*,*) 'Time:', time do p = 1,size(particle) call particle(p)%new_position( flow_field, deltt ) call particle(p)%print enddo enddo end program test_oil_particles