Homec4science
Diffusion FENNECS Lite (3Dmpimode) Browse

fennecslite/src/diagnose.f90
6b37315661933Dmpimode

Last Change
Raw File

diagnose.f90
View Options

SUBROUTINE diagnose(kstep)
! Diagnostics
USE basic
USE futils
USE hashtable
USE maxwsrce
USE beam, ONLY : partslist, partslist_towrite, epot, ekin, etot, etot0, Nplocs_all, collectparts
USE neutcol
USE fields, ONLY: phi_spline, nbmoments
USE mpi
USE geometry
USE splinebound
USE weighttypes
USE psupply
USE filemanip
USE magnet
IMPLICIT NONE
INTEGER, INTENT(in) :: kstep
INTEGER, PARAMETER :: BUFSIZE = 20
CHARACTER(len=128) :: str, fname, grpname
CHARACTER(len=12):: charspec
INTEGER:: Ntotal
INTEGER:: partsrank, partsdim(2)
INTEGER:: Nplocs_all_save(64)
INTEGER:: i, nbbounds
INTEGER, ALLOCATABLE, SAVE :: partnbBuffer(:,:)
INTEGER, ALLOCATABLE, SAVE :: lossesBuffer(:)
REAL(kind=db), ALLOCATABLE :: magr(:), magz(:)
REAL(kind=db), ALLOCATABLE :: tempBr(:, :), tempBz(:, :), tempAthet(:, :)
INTEGER :: magn(2), magrank, magfid
!________________________________________________________________________________
! 1. Initial diagnostics
IF( kstep .EQ. 0 .and. mpirank .eq. 0) THEN
! Only process 0 should save on file
WRITE(*,'(a)') ' Initial diagnostics'
! 1.1 Initial run or when NEWRES set to .TRUE.
IF( .NOT. nlres .OR. newres) THEN
CALL creatf(resfile, fidres, 'FENNECS result file', real_prec='d')
WRITE(*,'(3x,a,a)') TRIM(resfile), ' created'
! Label the run
IF( LEN_TRIM(label1).GT.0 ) CALL attach(fidres, "/", "label1", TRIM(label1))
IF( LEN_TRIM(label2).GT.0 ) CALL attach(fidres, "/", "label2", TRIM(label2))
IF( LEN_TRIM(label3).GT.0 ) CALL attach(fidres, "/", "label3", TRIM(label3))
IF( LEN_TRIM(label4).GT.0 ) CALL attach(fidres, "/", "label4", TRIM(label4))
! Job number
jobnum = 0
! Data group
CALL creatg(fidres, "/data", "data")
CALL creatg(fidres, "/data/var0d", "0d history arrays")
CALL creatg(fidres, "/data/var1d","1d history arrays")
CALL creatg(fidres, "/data/part", "part phase space")
CALL creatg(fidres, "/data/fields", "Electro static potential and Er Ez fields")
CALL attach(fidres, "/data", "file_version", 2.0)
! File group
CALL creatg(fidres, "/files", "files")
CALL attach(fidres, "/files", "jobnum", jobnum)
CALL putarr(fidres, "/data/var1d/rgrid", rgrid)
CALL putarr(fidres, "/data/var1d/zgrid", zgrid)
CALL creatd(fidres, 1, (/ 64 /), "/data/var0d/Nplocs_all", "Nplocs_all")
CALL creatd(fidres, 1, SHAPE(partslist(1)%nblost), "/data/var0d/losses", "losses")
! Initialize time-dependent particle variables
CALL creatd(fidres, 1, SHAPE(partslist(1)%pos(1,:)), "/data/part/R", "R")
CALL creatd(fidres, 1, SHAPE(partslist(1)%pos(3,:)), "/data/part/Z", "Z")
CALL creatd(fidres, 1, SHAPE(partslist(1)%pos(2,:)), "/data/part/THET", "THET")
CALL creatd(fidres, 1, SHAPE(partslist(1)%U(1,:)), "/data/part/UR", "UR")
CALL creatd(fidres, 1, SHAPE(partslist(1)%U(3,:)), "/data/part/UZ", "UZ")
CALL creatd(fidres, 1, SHAPE(partslist(1)%U(2,:)), "/data/part/UTHET", "UTHET")
CALL creatd(fidres, 1, SHAPE(partslist(1)%partindex), "/data/part/Rindex", "Rindex")
CALL creatd(fidres, 1, SHAPE(partslist(1)%partindex), "/data/part/Zindex", "Zindex")
CALL creatd(fidres, 1, SHAPE(partslist(1)%partindex), "/data/part/partindex", "partindex")
CALL creatd(fidres, 1, SHAPE(partslist(1)%pot), "/data/part/pot", "pot")
CALL creatd(fidres, 0, SHAPE(time), "/data/part/time", "time" )
CALL creatd(fidres, 0, SHAPE(Ntotal), "/data/part/Nparts", "number of remaining parts in the simulation space")
CALL creatd(fidres, 1, (/7/), "/data/part/nbchange", "number of added parts, lost parts zm,zp,rm,rp, and collisions per type io, ela")
CALL attach(fidres,'/data/part', "q", partslist(1)%q)
CALL attach(fidres,'/data/part', "m", partslist(1)%m)
CALL attach(fidres,'/data/part', "weight", partslist(1)%weight)
! Initialize time-dependent grid variables
CALL creatd(fidres, 0, SHAPE(time), "/data/fields/time", "time" )
CALL creatd(fidres, 1, SHAPE(pot), "/data/fields/pot", "pot")
CALL creatd(fidres, 1, SHAPE(Er), "/data/fields/Er", "Er")
CALL creatd(fidres, 1, SHAPE(Ez), "/data/fields/Ez", "Ez")
CALL creatd(fidres, 1, SHAPE(npartspercell), "/data/fields/npartspercell", "npartspercell")
CALL creatd(fidres, 1, SHAPE(phi_spline), "/data/fields/phi", "spline form of Phi")
CALL creatd(fidres, 2, (/nbmoments,nrank(1)*nrank(2)/), "/data/fields/moments", "moments",compress=.true.,chunking=(/1,nrank(1)*nrank(2)/))
CALL putarr(fidres, "/data/fields/Br", Br)
CALL putarr(fidres, "/data/fields/Bz", Bz)
CALL putarr(fidres, "/data/fields/Athet", Athet)
CALL putarr(fidres, "/data/fields/volume", Volume)
! We save the magnetic field as defined by the h5 file
IF(len_trim(magnetfile) .gt. 1 .and. .not. allocated(the_coils) ) THEN
CALL openf(trim(magnetfile), magfid, 'r', real_prec='d')
CALL getdims(magfid, '/mag/Athet', magrank, magn)
ALLOCATE (magr(magn(2)), magz(magn(1)))
ALLOCATE (tempAthet(magn(1), magn(2)), tempBr(magn(1), magn(2)), tempBz(magn(1), magn(2)))
! Read r and z coordinates for the definition of A_\thet, and B
CALL getarr(magfid, '/mag/r', magr)
CALL getarr(magfid, '/mag/z', magz)
CALL getarr(magfid, '/mag/Athet', tempAthet)
IF (isdataset(magfid, '/mag/Br') .and. isdataset(magfid, '/mag/Bz')) THEN
CALL getarr(magfid, '/mag/Br', tempBr)
CALL getarr(magfid, '/mag/Bz', tempBz)
END IF
CALL creatg(fidres, '/data/inputmag')
CALL putarr(fidres, '/data/inputmag/r',magr)
CALL putarr(fidres, '/data/inputmag/z',magz)
CALL putarr(fidres, '/data/inputmag/Athet',tempAthet)
CALL putarr(fidres, '/data/inputmag/Br',tempBr)
CALL putarr(fidres, '/data/inputmag/Bz',tempBz)
CALL closef(magfid)
END IF
! 1.2 Restart run
ELSE
CALL cp2bk(resfile) ! backup previous result file
CALL openf(resfile, fidres, real_prec='d')
WRITE(*,'(3x,a,a)') TRIM(resfile), ' open'
CALL getatt(fidres, "/files", "jobnum", jobnum)
jobnum = jobnum+1
WRITE(*,'(3x,a,i3)') "Current Job Number =", jobnum
CALL attach(fidres, "/files", "jobnum", jobnum)
END IF
! Add input namelist variables as attributes of /data/input
WRITE(str,'(a,i2.2)') "/data/input.",jobnum
CALL creatg(fidres, TRIM(str))
CALL attach(fidres, TRIM(str), "job_time", job_time)
CALL attach(fidres, TRIM(str), "extra_time", extra_time)
CALL attach(fidres, TRIM(str), "dt", dt*tnorm)
CALL attach(fidres, TRIM(str), "tmax", tmax)
CALL attach(fidres, TRIM(str), "nrun", nrun)
CALL attach(fidres, TRIM(str), "nlres", nlres)
CALL attach(fidres, TRIM(str), "nlsave", nlsave)
CALL attach(fidres, TRIM(str), "newres", newres)
CALL attach(fidres, TRIM(str), "nz", nz)
CALL attach(fidres, TRIM(str), "nr", nr)
CALL putarr(fidres, TRIM(str)//"/lz", lz)
CALL attach(fidres, TRIM(str), "nplasma", nplasma)
CALL attach(fidres, TRIM(str), "potinn", potinn)
CALL attach(fidres, TRIM(str), "potout", potout)
CALL attach(fidres, TRIM(str), "B0", B0)
CALL attach(fidres, TRIM(str), "Rcurv", Rcurv)
CALL attach(fidres, TRIM(str), "width", width)
CALL attach(fidres, TRIM(str), "n0", n0)
CALL attach(fidres, TRIM(str), "temp", partslist(1)%temperature)
CALL attach(fidres, TRIM(str), "it0d", it0d)
CALL attach(fidres, TRIM(str), "it2d", it2d)
CALL attach(fidres, TRIM(str), "itparts", itparts)
CALL attach(fidres, TRIM(str), "nlclassical", nlclassical)
CALL attach(fidres, TRIM(str), "nlPhis", nlPhis)
CALL attach(fidres, TRIM(str), "qsim", partslist(1)%q*partslist(1)%weight)
CALL attach(fidres, TRIM(str), "msim", partslist(1)%m*partslist(1)%weight)
CALL attach(fidres, TRIM(str), "startstep", cstep)
CALL attach(fidres, TRIM(str), "H0", partslist(1)%H0)
CALL attach(fidres, TRIM(str), "P0", partslist(1)%P0)
CALL putarr(fidres, TRIM(str)//"/femorder", femorder)
CALL putarr(fidres, TRIM(str)//"/ngauss", ngauss)
CALL putarr(fidres, TRIM(str)//"/nnr", nnr)
CALL putarr(fidres, TRIM(str)//"/nnz", nnz)
CALL putarr(fidres, TRIM(str)//"/radii", radii)
CALL putarr(fidres, TRIM(str)//"/plasmadim", plasmadim)
CALL attach(fidres, TRIM(str), "rawparts", .true.)
CALL attach(fidres, TRIM(str), "nbspecies", nbspecies)
CALL putarr(fidres, TRIM(str)//"/potxt", potxt)
CALL putarr(fidres, TRIM(str)//"/Erxt", Erxt)
CALL putarr(fidres, TRIM(str)//"/Ezxt", Ezxt)
! Save geometry parameters for non conforming boundary conditions
CALL geom_diag(fidres,str,rnorm)
! Save geometry parameters for non conforming boundary conditions using b-spline curves
CALL splinebound_diag(fidres, str, the_domain)
! Save Maxwellsource parameters for the ad-hoc source
CALL maxwsrce_diag(fidres,str,vnorm)
! Save neutcol parameters for the electron collisions with neutrals
CALL neutcol_diag(fidres,str,vnorm)
! Save psupply parameters for the simulation of realistic power supplies
CALL psupply_diag(fidres,str)
IF(.not. isdataset(fidres,'/data/var0d/biases')) THEN
nbbounds=2
if(the_domain%nbsplines .gt. 0) nbbounds=the_domain%nbsplines
CALL creatd(fidres, 1, (/nbbounds/), "/data/var0d/biases", "biases")
END IF
! Save STDIN of this run
WRITE(str,'(a,i2.2)') "/files/STDIN.",jobnum
INQUIRE(unit=lu_in, name=fname)
CALL putfile(fidres, TRIM(str), TRIM(fname))
! Prepare hdf5 file for storing other species
DO i=2,nbspecies
WRITE(grpname,'(a,i2)')'/data/part/',i
CALL create_parts_group(partslist(i),trim(grpname),time)
END DO
CALL attach(fidres, "/data/part", "nbspecies", nbspecies)
! Initialize buffers for 0d history arrays
CALL htable_init(hbuf0, BUFSIZE)
CALL set_htable_fileid(hbuf0, fidres, "/data/var0d")
END IF
IF(kstep .EQ. 0) THEN
CLOSE(lu_in)
ALLOCATE(partnbBuffer(3+size(partslist(1)%nblost,1),nbspecies))
partnbBuffer=0
ALLOCATE(lossesBuffer(size(partslist(1)%nblost,1)))
lossesBuffer=0
END IF
!________________________________________________________________________________
! 2. Periodic diagnostics
IF( kstep .NE. -1) THEN
IF(modulo(step,ittracer) .eq. 0 .or. nlend) THEN
! We gather the traced particles on the mpi host
DO i=1,nbspecies
IF(partslist(i)%is_test) CALL collectparts(partslist(i),partslist_towrite(i))
END DO
END IF
IF(modulo(step,itrestart) .eq. 0 .or. modulo(step,itparts) .eq. 0 .or. nlend) THEN
! We gather the traced particles on the mpi host
DO i=1,nbspecies
CALL collectparts(partslist(i),partslist_towrite(i))
END DO
END IF
lossesBuffer=lossesBuffer+partslist(1)%nblost
DO i=1,nbspecies
partnbBuffer(1,i)=partnbBuffer(1,i)+partslist(i)%nbadded
partnbBuffer(2,i)=partnbBuffer(2,i)+partslist(i)%nbionisation
partnbBuffer(3,i)=partslist(i)%Nploc
partnbBuffer(4:,i)=partnbBuffer(4:,i)+partslist(i)%nblost
partslist(i)%nbadded=0
partslist(i)%nblost=0
partslist(i)%nbionisation=0
END DO
IF(modulo(step,ittext) .eq. 0 .or. nlend) THEN
! We gather the number of gained and lost particles on the mpi host
IF(mpirank .eq.0 ) THEN
CALL MPI_REDUCE(MPI_IN_PLACE, partnbBuffer, nbspecies*(3+size(partslist(1)%nblost,1)), MPI_INTEGER, MPI_SUM, 0, MPI_COMM_WORLD, ierr)
ELSE
CALL MPI_REDUCE(partnbBuffer, partnbBuffer, nbspecies*(3+size(partslist(1)%nblost,1)), MPI_INTEGER, MPI_SUM, 0, MPI_COMM_WORLD, ierr)
partnbBuffer=0
END IF
END IF
! Only process 0 should save on file
IF(mpirank .ne. 0) RETURN
IF (mpisize .gt. 1) THEN
partslist(1)%Nptot=sum(Nplocs_all)
END IF
IF(modulo(step,ittext) .eq. 0 .or. nlend) THEN
WRITE(*,'(a,1x,i8.8,a1,i8.8,20x,a,1pe10.3)') '*** Timestep (this run/total) =', &
& step, '/', cstep, 'Time =', time
IF( abs(etot).gt. 0) THEN
WRITE(*,'(a,6(1pe12.4),1x,i8.8,a1,i8.8)') 'Epot, Ekin, Etot, Etot0, Eerr, Eerr rel, Nbparts/Ntotal', epot, ekin, etot, etot0, etot-etot0,(etot-etot0)/etot, partslist(1)%Nptot,'/',nplasma
ELSE
WRITE(*,'(a,4(1pe12.4),1x,i8.8,a1,i8.8)') 'Epot, Ekin, Etot, Eerr, Nbparts/Ntotal', epot, ekin, etot, etot-etot0, partslist(1)%Nptot,'/',nplasma
END IF
IF(mpisize .gt. 1 ) THEN
WRITE(*,'(a,64i10.7)') 'Nbparts per proc', Nplocs_all
END IF
WRITE(*,'(a)')"speci, added, ionis, tot var, tot, Losses (zmin zmax rmin rmax boundaries(i))"
WRITE(charspec,'(a,i02,a)') '(i04,',size(partnbBuffer,1)+1,'i9.7)'
DO i=1,nbspecies
WRITE(*,charspec) i, partnbBuffer(1,i),partnbBuffer(2,i), partnbBuffer(1,i)-sum(partnbBuffer(4:,i)), partnbBuffer(3,i),-partnbBuffer(4:,i)
partslist(i)%nptot= partnbBuffer(3,i)
END DO
partnbBuffer=0
END IF
! 2.1 0d history arrays
IF( kstep .eq. 0 .and. nlres .and. (.not. newres)) RETURN ! if we do a restart, we don't want to save the same data twice
IF(modulo(step,it0d).eq. 0 .or. nlend) THEN
CALL add_record(hbuf0, "time", "simulation time", time)
CALL add_record(hbuf0, "epot", "potential energy", epot)
CALL add_record(hbuf0, "ekin", "kinetic energy", ekin)
CALL add_record(hbuf0, "etot", "total energy", etot)
CALL add_record(hbuf0, "etot0", "theoretical total energy", etot0)
CALL add_record(hbuf0, "nbparts", "number of remaining parts in the simulation space", REAL(partslist(1)%Nptot,kind=db))
CALL add_record(hbuf0,"current", "unscaled current flowing between the electrodes of the power supplies [A]", the_ps%current(1)*qnorm/tnorm)
CALL htable_endstep(hbuf0)
Nplocs_all_save=0
Nplocs_all_save(1:mpisize)=Nplocs_all(0:mpisize-1)
CALL append(fidres, "/data/var0d/Nplocs_all", REAL(Nplocs_all_save,kind=db))
CALL append(fidres, "/data/var0d/biases", the_ps%biases*phinorm)
! We gather the number of lost particles on the mpi host
IF(mpirank .eq.0 ) THEN
CALL MPI_REDUCE(MPI_IN_PLACE, lossesBuffer, size(partslist(1)%nblost), MPI_INTEGER, MPI_SUM, 0, MPI_COMM_WORLD, ierr)
ELSE
CALL MPI_REDUCE(partnbBuffer, lossesBuffer, size(partslist(1)%nblost), MPI_INTEGER, MPI_SUM, 0, MPI_COMM_WORLD, ierr)
lossesBuffer=0
END IF
CALL append(fidres, "/data/var0d/losses", REAL(lossesBuffer,kind=db))
lossesBuffer=0
END IF
! 2.2 2d profiles
IF(modulo(step,it2d).eq. 0 .or. nlend) THEN
CALL append(fidres, "/data/fields/time", time)
CALL append(fidres, "/data/fields/pot", pot*phinorm)
CALL append(fidres, "/data/fields/Er", Er*enorm)
CALL append(fidres, "/data/fields/Ez", Ez*enorm)
CALL append(fidres, "/data/fields/npartspercell", npartspercell)
CALL append(fidres, "/data/fields/phi", phi_spline*phinorm)
CALL append(fidres, "/data/fields/moments", partslist(1)%fluidvel)
DO i=2,nbspecies
IF ( .not. partslist(i)%calc_moments) CYCLE
WRITE(grpname,'(a,i2,a)')'/data/part/',i,'/'
CALL append(fidres, trim(grpname) // "moments", partslist(i)%fluidvel)
END DO
END IF
! 2.3 main specie quantities
IF(modulo(step,itparts).eq. 0 .or. nlend) THEN
CALL append(fidres, "/data/part/time", time)
CALL append(fidres, "/data/part/Nparts", REAL(partslist(1)%Nptot,kind=db))
IF ( isdataset(fidres,'/data/part/R') ) THEN
CALL getdims(fidres, '/data/part/R', partsrank, partsdim)
partsdim(1)=min(size(partslist_towrite(1)%pos,2), partsdim(1))
CALL append(fidres, "/data/part/R", partslist_towrite(1)%pos(1,1:partsdim(1))*rnorm)
CALL append(fidres, "/data/part/Z", partslist_towrite(1)%pos(3,1:partsdim(1))*rnorm)
CALL append(fidres, "/data/part/THET", partslist_towrite(1)%pos(2,1:partsdim(1)))
CALL append(fidres, "/data/part/UZ", 0.5*(partslist_towrite(1)%U(3,1:partsdim(1))/partslist_towrite(1)%gamma(1:partsdim(1))+partslist_towrite(1)%Uold(3,1:partsdim(1))/partslist_towrite(1)%gammaold(1:partsdim(1))))
CALL append(fidres, "/data/part/UR", 0.5*(partslist_towrite(1)%U(1,1:partsdim(1))/partslist_towrite(1)%gamma(1:partsdim(1))+partslist_towrite(1)%Uold(1,1:partsdim(1))/partslist_towrite(1)%gammaold(1:partsdim(1))))
CALL append(fidres, "/data/part/UTHET", 0.5*(partslist_towrite(1)%U(2,1:partsdim(1))/partslist_towrite(1)%gamma(1:partsdim(1))+partslist_towrite(1)%Uold(2,1:partsdim(1))/partslist_towrite(1)%gammaold(1:partsdim(1))))
CALL append(fidres, "/data/part/pot", partslist_towrite(1)%pot(1:partsdim(1))*phinorm)
CALL append(fidres, "/data/part/Rindex", REAL(partslist_towrite(1)%cellindex(1,1:partsdim(1)),kind=db))
CALL append(fidres, "/data/part/Zindex", REAL(partslist_towrite(1)%cellindex(3,1:partsdim(1)),kind=db))
CALL append(fidres, "/data/part/partindex", REAL(partslist_towrite(1)%partindex(1:partsdim(1)),kind=db))
END IF
END IF
! 2.4 other species
IF(modulo(step,itparts).eq. 0 .or. nlend) THEN
DO i=2,nbspecies
IF (.not. partslist(i)%is_test) CYCLE
WRITE(grpname,'(a,i2,a)')'/data/part/',i,'/'
CALL append(fidres, trim(grpname) // "time", time)
CALL append(fidres, trim(grpname) //"Nparts", REAL(partslist(i)%Nptot,kind=db))
IF ( isdataset(fidres,trim(grpname)//'R') ) THEN
CALL getdims(fidres, trim(grpname) // 'R', partsrank, partsdim)
partsdim(1)=min(size(partslist_towrite(i)%pos,2), partsdim(1))
CALL append(fidres, trim(grpname) // "R", partslist_towrite(i)%pos(1,1:partsdim(1))*rnorm)
CALL append(fidres, trim(grpname) // "Z", partslist_towrite(i)%pos(3,1:partsdim(1))*rnorm)
CALL append(fidres, trim(grpname) // "THET", partslist_towrite(i)%pos(2,1:partsdim(1)))
CALL append(fidres, trim(grpname) // "UZ", 0.5*(partslist_towrite(i)%U(3,1:partsdim(1))/partslist_towrite(i)%gamma(1:partsdim(1)) + partslist_towrite(i)%Uold(3,1:partsdim(1))/partslist_towrite(i)%gammaold(1:partsdim(1))))
CALL append(fidres, trim(grpname) // "UR", 0.5*(partslist_towrite(i)%U(1,1:partsdim(1))/partslist_towrite(i)%gamma(1:partsdim(1)) + partslist_towrite(i)%Uold(1,1:partsdim(1))/partslist_towrite(i)%gammaold(1:partsdim(1))))
CALL append(fidres, trim(grpname) // "UTHET", 0.5*(partslist_towrite(i)%U(2,1:partsdim(1))/partslist_towrite(i)%gamma(1:partsdim(1)) + partslist_towrite(i)%Uold(2,1:partsdim(1))/partslist_towrite(i)%gammaold(1:partsdim(1))))
CALL append(fidres, trim(grpname) // "pot", partslist_towrite(i)%pot(1:partsdim(1))*phinorm)
CALL append(fidres, trim(grpname) // "Rindex", REAL(partslist_towrite(i)%cellindex(1,1:partsdim(1)),kind=db))
CALL append(fidres, trim(grpname) // "Zindex", REAL(partslist_towrite(i)%cellindex(3,1:partsdim(1)),kind=db))
CALL append(fidres, trim(grpname) // "partindex", REAL(partslist_towrite(i)%partindex(1:partsdim(1)),kind=db))
END IF
END DO
END IF
!________________________________________________________________________________
! 3. Final diagnostics
ELSE
! Only process 0 should save on file
IF(mpirank .ne. 0) RETURN
! Flush 0d history array buffers
CALL htable_hdf5_flush(hbuf0)
! Close all diagnostic files
CALL closef(fidres)
!________________________________________________________________________________
END IF
CONTAINS
SUBROUTINE create_parts_group(p,grpname, time)
USE beam,ONLY: particles
type(particles):: p
real(kind=db):: time
character(len=*):: grpname
IF(isgroup(fidres, trim(grpname))) RETURN
CALL creatg(fidres, grpname, "specific specie phase space")
CALL creatd(fidres, 0, SHAPE(time), trim(grpname) // "/time", "time")
CALL creatd(fidres, 0, SHAPE(time), trim(grpname) //"/Nparts", "number of remaining parts")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/R", "radial pos")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/Z", "axial pos")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/THET", "azimuthal pos")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/UZ", "axial beta*gamma")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/UR", "radial beta*gamma")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/UTHET", "azimuthal beta*gamma")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/pot", "electric potential")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/Rindex", "radial grid index")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/Zindex", "axial grid index")
CALL creatd(fidres, 1, SHAPE(p%pot), trim(grpname) // "/partindex", "particle index")
CALL creatd(fidres, 1, (/7/), trim(grpname) // "nbchange", "number of added parts, lost parts zm,zp,rm,rp, and collisions per type io, ela")
CALL attach(fidres,trim(grpname), "q", p%q)
CALL attach(fidres,trim(grpname), "m", p%m)
CALL attach(fidres,trim(grpname), "weight", p%weight)
CALL creatd(fidres, 2, (/nbmoments,nrank(1)*nrank(2)/), trim(grpname) // "/moments", "moments")
END SUBROUTINE create_parts_group
END SUBROUTINE diagnose
c4science · Help