Subroutine that loads the information of the initial condition of the different particle species. This subroutine calls the subroutine that generates the initial energy and pitch angle distribution functions.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(KORC_PARAMS), | intent(in) | :: | params | Core KORC simulation parameters. |
||
| type(FIELDS), | intent(in) | :: | F | An instance of KORC's derived type FIELDS containing all the information about the fields used in the simulation. See korc_types and korc_fields. |
||
| type(PROFILES), | intent(inout) | :: | P | |||
| type(SPECIES), | intent(out), | DIMENSION(:), ALLOCATABLE | :: | spp | An instance of KORC's derived type SPECIES containing all the information of different electron species. See korc_types. |
subroutine initialize_particles(params,F,P,spp)
!! @note Subroutine that loads the information of the initial condition
!! of the different particle species. This subroutine calls
!! the subroutine that generates the initial energy and pitch angle
!! distribution functions. @endnote
TYPE(KORC_PARAMS), INTENT(IN) :: params
!! Core KORC simulation parameters.
TYPE(FIELDS), INTENT(IN) :: F
!! An instance of KORC's derived type FIELDS containing all the information
!! about the fields used in the simulation. See [[korc_types]]
!!and [[korc_fields]].
TYPE(PROFILES), INTENT(INOUT) :: P
TYPE(SPECIES), DIMENSION(:), ALLOCATABLE, INTENT(OUT) :: spp
!! An instance of KORC's derived type SPECIES containing all the information
!! of different electron species. See [[korc_types]].
REAL(rp), DIMENSION(:), ALLOCATABLE :: ppp
!! Number of computational particles per MPI process
!! used to simulate each electron species.
REAL(rp), DIMENSION(:), ALLOCATABLE :: q
!! Charge of each species.
REAL(rp), DIMENSION(:), ALLOCATABLE :: m
!! Mass of each species.
REAL(rp), DIMENSION(:), ALLOCATABLE :: Eo
!! Initial energy of each electron species in case of
!! using an initial mono-energetic distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: etao
!! Initial pitch-angle of each electron species in case of
!! using an initial mono-pitch-angle distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: Eo_lims
!! Minimum and maximum energy limits of a given initial
!! non-mono-energetic distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: etao_lims
!! Minimum and maximum pitch-angle limits of a given initial
!! non-mono-pitch-angle distribution.
LOGICAL, DIMENSION(:), ALLOCATABLE :: runaway
!! Flag to decide whether a given electron is a runaway (runaway=T)
!! or not (runaway=F).
CHARACTER(MAX_STRING_LENGTH),DIMENSION(:),ALLOCATABLE :: spatial_distribution
!! String describing the type of initial spatial distribution for
!! each electron species.
CHARACTER(MAX_STRING_LENGTH), DIMENSION(:), ALLOCATABLE :: energy_distribution
!! String describing the type of initial energy distribution for each
!! electron species.
CHARACTER(MAX_STRING_LENGTH), DIMENSION(:), ALLOCATABLE :: pitch_distribution
!! String describing the type of initial pitch-angle distribution
!! for each electron species.
REAL(rp), DIMENSION(:), ALLOCATABLE :: Ro
!! Radial position of the center of the electrons' initial
!! spatial distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: PHIo
!! Azimuthal position of the electrons' initial spatial distribution,
!! in case of using a disk at a certain poloidal section.
REAL(rp), DIMENSION(:), ALLOCATABLE :: Zo
!! Height of the center of the electrons' initial spatial distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: r_inner
!! Minimum minor radius of the electrons' initial spatial distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: r_outter
!! Maximum minor radius of the electrons' initial spatial distribution.
REAL(rp), DIMENSION(:), ALLOCATABLE :: falloff_rate
!! Exponential falloff or standard deviation of a non-uniform radial
!! distribution of electrons.
REAL(rp), DIMENSION(:), ALLOCATABLE :: shear_factor
!! Shear factor used to generate an initial spatial
!! distribution with an elliptic poloidal cross section.
!! See <em>Carbajal and del-Castillo-Negrete, Nuclear Fusion,
!! submitted (2018)</em>.
REAL(rp), DIMENSION(:), ALLOCATABLE :: sigmaR
!! Variance of the first dimension of a 2D Gaussian, spatial
!! distribution function
REAL(rp), DIMENSION(:), ALLOCATABLE :: sigmaZ
!! Variance of the second dimension of a 2D Gaussian, spatial
!! distribution function
REAL(rp), DIMENSION(:), ALLOCATABLE :: theta_gauss
!! Angle of counter-clockwise rotation (in degrees) of 2D Gaussian
!! distribution relative to R,Z
REAL(rp), DIMENSION(:), ALLOCATABLE :: psi_max
!! Maximum value of the argument of the 2D gaussian exponential, used for an
!! indicator function that limits the region of MH sampling
REAL(rp), DIMENSION(:), ALLOCATABLE :: Spong_b
REAL(rp), DIMENSION(:), ALLOCATABLE :: Spong_w
REAL(rp), DIMENSION(:), ALLOCATABLE :: Spong_dlam
REAL(rp), DIMENSION(:), ALLOCATABLE :: dth,dgam
REAL(rp), DIMENSION(:), ALLOCATABLE :: dR
REAL(rp), DIMENSION(:), ALLOCATABLE :: dZ
INTEGER :: ii
!! Iterator of spp structure.
INTEGER :: mpierr
!! MPI error status.
REAL(rp), DIMENSION(:), ALLOCATABLE :: Xtrace
NAMELIST /plasma_species/ ppp,q,m,Eo,etao,Eo_lims,etao_lims,runaway, &
spatial_distribution,energy_distribution,pitch_distribution,Ro, &
PHIo,Zo,r_inner,r_outter,falloff_rate,shear_factor,sigmaR,sigmaZ, &
theta_gauss,psi_max,Xtrace,Spong_b,Spong_w,Spong_dlam,dth,dR,dZ,dgam
! Allocate array containing variables of particles for each species
ALLOCATE(spp(params%num_species))
ALLOCATE(ppp(params%num_species))
ALLOCATE(q(params%num_species))
ALLOCATE(m(params%num_species))
ALLOCATE(Eo(params%num_species))
ALLOCATE(etao(params%num_species))
ALLOCATE(Eo_lims(2_idef*params%num_species))
ALLOCATE(etao_lims(2_idef*params%num_species))
ALLOCATE(runaway(params%num_species))
ALLOCATE(spatial_distribution(params%num_species))
ALLOCATE(energy_distribution(params%num_species))
ALLOCATE(pitch_distribution(params%num_species))
ALLOCATE(Ro(params%num_species))
ALLOCATE(PHIo(params%num_species))
ALLOCATE(Zo(params%num_species))
ALLOCATE(r_inner(params%num_species))
ALLOCATE(r_outter(params%num_species))
ALLOCATE(falloff_rate(params%num_species))
ALLOCATE(shear_factor(params%num_species))
ALLOCATE(sigmaR(params%num_species))
ALLOCATE(sigmaZ(params%num_species))
ALLOCATE(theta_gauss(params%num_species))
ALLOCATE(psi_max(params%num_species))
ALLOCATE(Spong_b(params%num_species))
ALLOCATE(Spong_w(params%num_species))
ALLOCATE(Spong_dlam(params%num_species))
ALLOCATE(dth(params%num_species))
ALLOCATE(dgam(params%num_species))
ALLOCATE(dR(params%num_species))
ALLOCATE(dZ(params%num_species))
ALLOCATE(Xtrace(3_idef*params%num_species))
open(unit=default_unit_open,file=TRIM(params%path_to_inputs), &
status='OLD',form='formatted')
read(default_unit_open,nml=plasma_species)
close(default_unit_open)
do ii=1_idef,params%num_species
spp(ii)%runaway = runaway(ii)
spp(ii)%spatial_distribution = TRIM(spatial_distribution(ii))
spp(ii)%energy_distribution = TRIM(energy_distribution(ii))
spp(ii)%pitch_distribution = TRIM(pitch_distribution(ii))
spp(ii)%q = q(ii)*C_E
spp(ii)%m = m(ii)*C_ME
spp(ii)%ppp = ppp(ii)
spp(ii)%Ro = Ro(ii)
spp(ii)%PHIo = C_PI*PHIo(ii)/180.0_rp
spp(ii)%Zo = Zo(ii)
spp(ii)%r_inner = r_inner(ii)
spp(ii)%r_outter = r_outter(ii)
spp(ii)%falloff_rate = falloff_rate(ii)
spp(ii)%shear_factor = shear_factor(ii)
spp(ii)%sigmaR = sigmaR(ii)
spp(ii)%sigmaZ = sigmaZ(ii)
spp(ii)%theta_gauss = theta_gauss(ii)
spp(ii)%psi_max = psi_max(ii)
spp(ii)%Spong_w = Spong_w(ii)
spp(ii)%Spong_b = Spong_b(ii)
spp(ii)%Spong_dlam = Spong_dlam(ii)
spp(ii)%dth = dth(ii)
spp(ii)%dgam = dgam(ii)
spp(ii)%dR = dR(ii)
spp(ii)%dZ = dZ(ii)
! * * These values can change in initial_energy_pitch_dist * * !
spp(ii)%Eo = Eo(ii)*C_E
spp(ii)%Eo_lims = Eo_lims((ii-1_idef)*2_idef + 1_idef:2_idef*ii)*C_E
spp(ii)%etao = etao(ii)
spp(ii)%etao_lims = etao_lims((ii-1_idef)*2_idef + 1_idef:2_idef*ii)
! * * These values can change in initial_energy_pitch_dist * * !
if (spp(ii)%spatial_distribution.eq.'TRACER') &
spp(ii)%Xtrace = Xtrace((ii-1_idef)*3_idef + 1_idef:3_idef*ii)
ALLOCATE( spp(ii)%vars%X(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%V(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%Rgc(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%Y(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%E(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%B(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%PSI_P(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%ne(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%Te(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%Zeff(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%g(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%eta(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%mu(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%Prad(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%Pin(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%flag(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%wt(spp(ii)%ppp) )
! write(6,'("0 size of PSI_P: ",I16)') size(spp(ii)%vars%PSI_P)
spp(ii)%vars%X = 0.0_rp
spp(ii)%vars%V = 0.0_rp
spp(ii)%vars%Rgc = 0.0_rp
spp(ii)%vars%Y = 0.0_rp
spp(ii)%vars%E = 0.0_rp
spp(ii)%vars%B = 0.0_rp
spp(ii)%vars%PSI_P = 0.0_rp
spp(ii)%vars%ne = 0.0_rp
spp(ii)%vars%Te = 0.0_rp
spp(ii)%vars%Zeff = 0.0_rp
spp(ii)%vars%g = 0.0_rp
spp(ii)%vars%eta = 0.0_rp
spp(ii)%vars%mu = 0.0_rp
spp(ii)%vars%Prad = 0.0_rp
spp(ii)%vars%Pin = 0.0_rp
spp(ii)%vars%flag = 1_is
spp(ii)%vars%wt = 0.0_rp
if (params%orbit_model(1:2).eq.'GC') then
ALLOCATE( spp(ii)%vars%Y0(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%V0(spp(ii)%ppp) )
ALLOCATE( spp(ii)%vars%k1(spp(ii)%ppp,4) )
ALLOCATE( spp(ii)%vars%k2(spp(ii)%ppp,4) )
ALLOCATE( spp(ii)%vars%k3(spp(ii)%ppp,4) )
ALLOCATE( spp(ii)%vars%k4(spp(ii)%ppp,4) )
ALLOCATE( spp(ii)%vars%k5(spp(ii)%ppp,4) )
ALLOCATE( spp(ii)%vars%k6(spp(ii)%ppp,4) )
if (params%orbit_model(3:5)=='pre') then
ALLOCATE( spp(ii)%vars%gradB(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%curlb(spp(ii)%ppp,3) )
spp(ii)%vars%gradB = 0.0_rp
spp(ii)%vars%curlb = 0.0_rp
else if (params%orbit_model(3:6)=='grad') then
ALLOCATE( spp(ii)%vars%BR(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%BPHI(spp(ii)%ppp,3) )
ALLOCATE( spp(ii)%vars%BZ(spp(ii)%ppp,3) )
spp(ii)%vars%BR = 0.0_rp
spp(ii)%vars%BPHI = 0.0_rp
spp(ii)%vars%BZ = 0.0_rp
end if
ALLOCATE( spp(ii)%vars%RHS(spp(ii)%ppp,5) )
spp(ii)%vars%Y0 = 0.0_rp
spp(ii)%vars%V0 = 0.0_rp
spp(ii)%vars%k1 = 0.0_rp
spp(ii)%vars%k2 = 0.0_rp
spp(ii)%vars%k3 = 0.0_rp
spp(ii)%vars%k4 = 0.0_rp
spp(ii)%vars%k5 = 0.0_rp
spp(ii)%vars%k6 = 0.0_rp
spp(ii)%vars%RHS = 0.0_rp
end if
end do
P%R0_RE=spp(1)%Ro
P%Z0_RE=spp(1)%Zo
P%n_REr0=max(sqrt(spp(1)%psi_max*2*spp(1)%sigmaR**2), &
sqrt(spp(1)%psi_max*2*spp(1)%sigmaZ**2))
call initial_energy_pitch_dist(params,spp)
DEALLOCATE(ppp)
DEALLOCATE(q)
DEALLOCATE(m)
DEALLOCATE(Eo)
DEALLOCATE(etao)
DEALLOCATE(Eo_lims)
DEALLOCATE(etao_lims)
DEALLOCATE(runaway)
DEALLOCATE(spatial_distribution)
DEALLOCATE(energy_distribution)
DEALLOCATE(pitch_distribution)
DEALLOCATE(Ro)
DEALLOCATE(PHIo)
DEALLOCATE(Zo)
DEALLOCATE(r_inner)
DEALLOCATE(r_outter)
DEALLOCATE(falloff_rate)
DEALLOCATE(shear_factor)
DEALLOCATE(sigmaR)
DEALLOCATE(sigmaZ)
DEALLOCATE(theta_gauss)
DEALLOCATE(psi_max)
DEALLOCATE(Spong_b)
DEALLOCATE(Spong_w)
DEALLOCATE(Spong_dlam)
DEALLOCATE(dth)
DEALLOCATE(dgam)
DEALLOCATE(dR)
DEALLOCATE(dZ)
DEALLOCATE(Xtrace)
end subroutine initialize_particles