From 429d004247a53b177588310b700541c7428cc8e3 Mon Sep 17 00:00:00 2001
From: holger <holger-grosshans@gmx.de>
Date: Mon, 20 Apr 2020 09:40:44 +0200
Subject: [PATCH] modified: src/particles.f90 modified:
src/particlesTransport.f90 modified: src/pre.f90 modified:
src/var.f90
---
src/particles.f90 | 35 +++++++++++++++++++++++++++++++----
src/particlesTransport.f90 | 22 +++++++++++++++-------
src/pre.f90 | 1 +
src/var.f90 | 5 +++++
4 files changed, 52 insertions(+), 11 deletions(-)
diff --git a/src/particles.f90 b/src/particles.f90
index 8fbc53b..6e53b44 100644
--- a/src/particles.f90
+++ b/src/particles.f90
@@ -6,13 +6,14 @@
if (pnd.ne.0.or.npTot.ne.0) then
if (myid.eq.0) write(*,'(a)',advance='no') 'particle'
- if (((bcx.eq.'i').and.(nt.ge.ntseed)).or.(nt.eq.ntseed)) call partSeed
+ if (((bcx.eq.'i').and.(nt.ge.ntseed)).or.(nt.eq.ntseed)) call particlesSeed
+ call particlesToCells
call fluidVelocity
if (elForceScheme.ne.2) call forcesGauss
if (elForceScheme.ne.1) call forcesCoulomb
- call partVelocity
- call partCollide
- call partTransport !particles from neighbour gc added
+ call particlesVelocity
+ call particlesCollide
+ call particlesTransport !particles from neighbour gc added
call chargeDensity
call momentumCoupling !particles from neighbour gc removed
if (myid.eq.0) write(*,*)
@@ -339,3 +340,29 @@
return
end
+
+
+!####################################################################
+!> @author Holger Grosshans
+!> @brief locate particles on Eulerian grid
+ subroutine particlesToCells
+ use var
+
+ npic=0
+
+ do 1 n=1,np
+ npic(ip(n),jp(n),lp(n))=npic(ip(n),jp(n),lp(n))+1
+1 enddo
+
+ if (allocated(nic)) deallocate(nic)
+ allocate(nic(ii,jj,ll,maxval(npic)))
+
+ nic=0
+
+ do 2 n=1,np
+ nic(ip(n),jp(n),lp(n),count(nic(ip(n),jp(n),lp(n),:).ne.0)+1)=n
+2 enddo
+
+
+ return
+ end
diff --git a/src/particlesTransport.f90 b/src/particlesTransport.f90
index 7584075..36c31b2 100644
--- a/src/particlesTransport.f90
+++ b/src/particlesTransport.f90
@@ -68,7 +68,7 @@
!####################################################################
!> @author Holger Grosshans
!> @brief explicit first order integration of particle velocities
- subroutine partVelocity
+ subroutine particlesVelocity
use var
real(kind=pr),dimension(3) :: f_fl,f_g
real(kind=pr) :: fact,dragdt,drag
@@ -105,7 +105,7 @@
!####################################################################
!> @author Holger Grosshans
!> @brief move particles
- subroutine partTransport
+ subroutine particlesTransport
use var
use mpi
real(kind=pr) :: xpold,ypold,zpold
@@ -347,7 +347,7 @@
!####################################################################
!> @author Holger Grosshans
!> @brief compute particles collisions using ray casting
- subroutine partCollide
+ subroutine particlesCollide
use var
real(kind=pr),dimension(maxnp) :: xpc,ypc,zpc
real(kind=pr),dimension(3) :: relvelo,midlink,ncontact
@@ -356,13 +356,19 @@
r3n1,r3n2,rxx,projvtox,srel,projvton, &
closestrelvelo,fracdt, &
vex,vey,vez,vnx,vny,vnz
- integer :: n1,n2,numcol,ip,jp,lp,m
+ integer :: n1,n2,numcol,ip,jp,lp,m1,m2
numcol=0
- do 1 n1=2,np
- do 2 n2=1,(n1-1)
+ do 3 i=1,ii
+ do 3 j=1,jj
+ do 3 l=1,ll
+ do 1 m1=2,npic(i,j,l)
+ do 2 m2=1,(m1-1)
+
+ n1=nic(i,j,l,m1)
+ n2=nic(i,j,l,m2)
! condition I to check collision (same or neighbour cell):
if ((ip(n1).lt.ip(n2)-1.or.ip(n1).gt.ip(n2)+1).or. &
@@ -476,6 +482,7 @@
2 enddo
1 enddo
+3 enddo
if (myid.eq.0) write(*,'(a,i9)',advance='no') &
' |collide =',numcol
@@ -486,7 +493,7 @@
!####################################################################
!> @author Holger Grosshans
!> @brief seed particles
- subroutine partSeed
+ subroutine particlesSeed
use var
real(kind=pr),allocatable,dimension(:) :: randomx,randomy,randomz
integer :: n,i,j,l,ip,jp,lp,npseed,npseedTot, &
@@ -572,6 +579,7 @@
write(10,'(3(es10.2e2,x))') radp(n),q_el(n),partn(n)
1 enddo
+
close(10)
diff --git a/src/pre.f90 b/src/pre.f90
index 6573ad2..4fe0c22 100644
--- a/src/pre.f90
+++ b/src/pre.f90
@@ -191,6 +191,7 @@
dudt(ii,jj,ll),dvdt(ii,jj,ll),dwdt(ii,jj,ll))
! particles
+ allocate(npic(ii,jj,ll))
np=0;npp=0
call allocateParticleArrays
diff --git a/src/var.f90 b/src/var.f90
index b107e4a..576046c 100644
--- a/src/var.f90
+++ b/src/var.f90
@@ -131,6 +131,11 @@
integer, allocatable, dimension(:) :: &
wcollnum, & !< number of particle-wall collisions
ppcollnum !< number of particle-particle collisions
+ integer, allocatable, dimension(:,:,:) :: &
+ npic !< number of particles in a Eulerian cell
+ integer, allocatable, dimension(:,:,:,:) :: &
+ nic !< indices of particles in a Eulerian cell
+
! parallel
real(kind=pr) :: timebeg,timecom,timenow,timeend
--
GitLab