Homec4science
Diffusion FENNECS: Finite Element Non-Neutral Electron Cloud Simulator (stable) Browse

fennecs/src/celldiag_mod.f90
a3531b44407cstable

Last Change
Raw File

celldiag_mod.f90
View Options

!------------------------------------------------------------------------------
! EPFL/Swiss Plasma Center
!------------------------------------------------------------------------------
!
! MODULE: celldiag
!
!> @author
!> Guillaume Le Bars EPFL/SPC
!
! DESCRIPTION:
!> Represent a diagnostic positioned at cell indices (rindex,zindex) that saves the individual particles
!> position and velocity
!------------------------------------------------------------------------------
MODULE celldiag
!
USE constants
use mpi
USE mpihelper
USE basic, ONLY: mpirank, mpisize, vnorm, rnorm, Zbounds, zgrid, &
& nlclassical, nlmaxwellsource, phinorm, nbcelldiag
USE beam
USE futils
IMPLICIT NONE
PRIVATE
INTEGER, SAVE, ALLOCATABLE :: specieid(:) !< position of the specie in partslist
INTEGER, SAVE, ALLOCATABLE :: rindex(:) !< radial index for the diagnostic position
INTEGER, SAVE, ALLOCATABLE :: zindex(:) !< axial index for the diagnostic position
TYPE(particles), ALLOCATABLE, SAVE :: diagnosed_parts(:) !< Stores the particles properties at position (rindex,zindex)
CHARACTER(len=20), SAVE, ALLOCATABLE :: groupname(:) !< Name of the group in the hdf5 file
INTEGER, SAVE , ALLOCATABLE :: h5storelength(:) !< particles capacity of the hdf5 dataset
NAMELIST /celldiagparams/ specieid, rindex, zindex
PUBLIC:: celldiag_init, celldiag_save
contains
subroutine celldiag_init(lu_in, diagfile_id)
implicit none
INTEGER, INTENT(IN) :: lu_in
INTEGER, INTENT(IN):: diagfile_id
INTEGER:: i
ALLOCATE(specieid(nbcelldiag), rindex(nbcelldiag), zindex(nbcelldiag))
ALLOCATE(diagnosed_parts(nbcelldiag), groupname(nbcelldiag), h5storelength(nbcelldiag))
Rewind(lu_in)
if(nbcelldiag .gt. 0) then
READ(lu_in, celldiagparams)
if(mpirank .eq. 0) WRITE(*, celldiagparams)
Do i=1,nbcelldiag
CALL creat_parts(diagnosed_parts(i), 500)
IF(mpirank .eq. 0) THEN
WRITE(groupname(i),'(a,i2.2)') "/data/celldiag/",i
If(.not. isgroup(diagfile_id, "/data/celldiag/")) THEN
CALL creatg(diagfile_id, "/data/celldiag")
CALL attach(diagfile_id, "/data/celldiag", "nbcelldiag", nbcelldiag)
END IF
CALL celldiag_createh5group(diagfile_id, groupname(i), rindex(i), zindex(i), specieid(i), diagnosed_parts(i), h5storelength(i))
END IF
END DO
END IF
End subroutine celldiag_init
subroutine celldiag_createh5group(diagfile_id, groupname, rindex, zindex, specid, diag_parts, h5strlength)
INTEGER, INTENT(IN):: diagfile_id, rindex, zindex, specid
CHARACTER(len=*), INTENT(IN):: groupname
TYPE(particles), INTENT(IN):: diag_parts
INTEGER, INTENT(INOUT):: h5strlength
INTEGER:: partsrank, partsdim(2)
If(.not. isgroup(diagfile_id, TRIM(groupname))) CALL creatg(diagfile_id, TRIM(groupname))
CALL attach(diagfile_id, TRIM(groupname), "rindex", rindex)
CALL attach(diagfile_id, TRIM(groupname), "zindex", zindex)
CALL attach(diagfile_id, trim(groupname), "q", partslist(specid)%q)
CALL attach(diagfile_id, trim(groupname), "m", partslist(specid)%m)
CALL attach(diagfile_id, trim(groupname), "weight", partslist(specid)%weight)
If(.not. isdataset(diagfile_id, trim(groupname) // "/time")) CALL creatd(diagfile_id, 0, SHAPE(rindex), trim(groupname) // "/time", "time")
If(.not. isdataset(diagfile_id, trim(groupname) // "/Nparts")) CALL creatd(diagfile_id, 0, SHAPE(rindex), trim(groupname) //"/Nparts", "number of remaining parts")
If(.not. isdataset(diagfile_id, trim(groupname) // "/R")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/R", "radial pos")
If(.not. isdataset(diagfile_id, trim(groupname) // "/Z")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/Z", "axial pos")
If(.not. isdataset(diagfile_id, trim(groupname) // "/THET")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/THET", "azimuthal pos")
If(.not. isdataset(diagfile_id, trim(groupname) // "/UZ")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/UZ", "axial beta*gamma")
If(.not. isdataset(diagfile_id, trim(groupname) // "/UR")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/UR", "radial beta*gamma")
If(.not. isdataset(diagfile_id, trim(groupname) // "/UTHET")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pos(1,:)), trim(groupname) // "/UTHET", "azimuthal beta*gamma")
If(.not. isdataset(diagfile_id, trim(groupname) // "/pot")) CALL creatd(diagfile_id, 1, SHAPE(diag_parts%pot), trim(groupname) // "/pot", "electric potential")
CALL getdims(diagfile_id, trim(groupname) // '/R', partsrank, partsdim)
h5strlength=partsdim(1)
END subroutine celldiag_createh5group
subroutine celldiag_save(time, diagfile_id)
implicit none
REAL(kind=db), INTENT(IN) :: time
INTEGER, INTENT(IN) :: diagfile_id
INTEGER :: Nbtosave, i
! check if source is on
IF(.not. celldiag_on(time)) THEN
RETURN
END IF
Do i=1,nbcelldiag
CALL celldiag_save_specie(partslist(specieid(i)),rindex(i),zindex(i),diagnosed_parts(i))
END DO
!$OMP BARRIER
!$OMP MASTER
Do i=1,nbcelldiag
if(mpisize .gt. 1) then
call collectparts(diagnosed_parts(i))
else
diagnosed_parts(i)%Nptot=diagnosed_parts(i)%Nploc
end if
Nbtosave=min(diagnosed_parts(i)%Nptot,h5storelength(i))
CALL celldiag_write_specie(diagfile_id, diagnosed_parts(i), groupname(i), Nbtosave, time)
END DO
!$OMP END MASTER
!$OMP BARRIER
end subroutine celldiag_save
SUBROUTINE celldiag_save_specie(p, rindex, zindex, savedp)
Type(particles), INTENT(IN) :: p
Type(particles), INTENT(INOUT) :: savedp
INTEGER, INTENT(IN) :: rindex, zindex
INTEGER:: i, destcopyindex
savedp%Nploc=0
savedp%collected=.false.
IF (p%Nploc .gt. 0 .and. zindex .ge. Zbounds(mpirank) .and. zindex .lt. Zbounds(mpirank+1)) THEN
! Boundary condition at z direction
!$OMP DO
DO i=1,p%Nploc
! If the particle is in the correct cell, it is saved
IF (p%cellindex(3,i) .eq. zindex.and. p%cellindex(1,i) .eq. rindex ) THEN
!$OMP CRITICAL (diagparts)
savedp%Nploc=savedp%Nploc+1
destcopyindex= savedp%Nploc
!$OMP END CRITICAL (diagparts)
CALL copy_part(p,i,destcopyindex,savedp)
END IF
END DO
!$OMP END DO NOWAIT
END IF
END subroutine celldiag_save_specie
SUBROUTINE celldiag_write_specie(diagfile_id, savedp, groupname, Nbtosave, time)
Type(particles), INTENT(IN) :: savedp
INTEGER, INTENT(IN) :: diagfile_id
CHARACTER(LEN=*), INTENT(IN) :: groupname
INTEGER, INTENT(IN) :: Nbtosave
REAL(kind=db), INTENT(IN) :: time
IF(mpirank .eq. 0) THEN
CALL append(diagfile_id, trim(groupname) // "/time", time)
CALL append(diagfile_id, trim(groupname) // "/Nparts", REAL(savedp%Nptot,kind=db))
CALL append(diagfile_id, trim(groupname) // "/R", savedp%pos(1,1:Nbtosave)*rnorm)
CALL append(diagfile_id, trim(groupname) // "/Z", savedp%pos(3,1:Nbtosave)*rnorm)
CALL append(diagfile_id, trim(groupname) // "/THET", savedp%pos(2,1:Nbtosave))
CALL append(diagfile_id, trim(groupname) // "/UZ", savedp%U(3,1:Nbtosave)/savedp%gamma(1:Nbtosave))
CALL append(diagfile_id, trim(groupname) // "/UR", savedp%U(1,1:Nbtosave)/savedp%gamma(1:Nbtosave))
CALL append(diagfile_id, trim(groupname) // "/UTHET", savedp%U(2,1:Nbtosave)/savedp%gamma(1:Nbtosave))
CALL append(diagfile_id, trim(groupname) // "/pot", savedp%pot(1:Nbtosave)*phinorm)
END IF
END subroutine celldiag_write_specie
logical function celldiag_on(time)
REAL(kind=db), intent(in):: time
celldiag_on=.true.
end function
End Module celldiag
c4science · Help