initialize_particles Subroutine

  1. korc_initialize.f90
  2. korc_initialize
  3. initialize_particles

public subroutine initialize_particles(params, F, P, spp)

Arguments

Type IntentOptional AttributesName
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.

Number of computational particles per MPI process used to simulate each electron species.

Charge of each species.

Mass of each species.

Initial energy of each electron species in case of using an initial mono-energetic distribution.

Initial pitch-angle of each electron species in case of using an initial mono-pitch-angle distribution.

Minimum and maximum energy limits of a given initial non-mono-energetic distribution.

Minimum and maximum pitch-angle limits of a given initial non-mono-pitch-angle distribution.

Flag to decide whether a given electron is a runaway (runaway=T) or not (runaway=F).

String describing the type of initial spatial distribution for each electron species.

String describing the type of initial energy distribution for each electron species.

String describing the type of initial pitch-angle distribution for each electron species.

Radial position of the center of the electrons' initial spatial distribution.

Azimuthal position of the electrons' initial spatial distribution, in case of using a disk at a certain poloidal section.

Height of the center of the electrons' initial spatial distribution.

Minimum minor radius of the electrons' initial spatial distribution.

Maximum minor radius of the electrons' initial spatial distribution.

Exponential falloff or standard deviation of a non-uniform radial distribution of electrons.

Shear factor used to generate an initial spatial distribution with an elliptic poloidal cross section. See Carbajal and del-Castillo-Negrete, Nuclear Fusion, submitted (2018).

Variance of the first dimension of a 2D Gaussian, spatial distribution function

Variance of the second dimension of a 2D Gaussian, spatial distribution function

Angle of counter-clockwise rotation (in degrees) of 2D Gaussian distribution relative to R,Z

Maximum value of the argument of the 2D gaussian exponential, used for an indicator function that limits the region of MH sampling

Calls

proc~~initialize_particles~~CallsGraph proc~initialize_particles initialize_particles sigmaz sigmaz proc~initialize_particles->sigmaz spong_dlam spong_dlam proc~initialize_particles->spong_dlam zo zo proc~initialize_particles->zo spong_b spong_b proc~initialize_particles->spong_b xtrace xtrace proc~initialize_particles->xtrace ro ro proc~initialize_particles->ro proc~initial_energy_pitch_dist initial_energy_pitch_dist proc~initialize_particles->proc~initial_energy_pitch_dist spong_w spong_w proc~initialize_particles->spong_w dr dr proc~initialize_particles->dr sigmar sigmar proc~initialize_particles->sigmar dz dz proc~initialize_particles->dz eo_lims eo_lims proc~initialize_particles->eo_lims phio phio proc~initialize_particles->phio eno eno proc~initialize_particles->eno mpi_barrier mpi_barrier proc~initial_energy_pitch_dist->mpi_barrier

Called by

proc~~initialize_particles~~CalledByGraph proc~initialize_particles initialize_particles program~main main program~main->proc~initialize_particles

Contents

Source Code

initialize_particles

Source Code

  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)%pinit = pinit(ii)
       spp(ii)%pRE = pinit(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 = Eno(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%ni(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%flagCon(spp(ii)%ppp) )       
       ALLOCATE( spp(ii)%vars%flagCol(spp(ii)%ppp) )
       ALLOCATE( spp(ii)%vars%initLCFS(spp(ii)%ppp) )
       ALLOCATE( spp(ii)%vars%flagRE(spp(ii)%ppp) )
       ALLOCATE( spp(ii)%vars%wt(spp(ii)%ppp) )
#ifdef FIO
       ALLOCATE( spp(ii)%vars%hint(spp(ii)%ppp))
#endif
       
       !     write(output_unit_write,'("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%ni = 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%flagCon(1:spp(ii)%pinit) = 1_is
       spp(ii)%vars%flagCol(1:spp(ii)%pinit) = 1_is
       spp(ii)%vars%flagRE(1:spp(ii)%pinit) = 1_is
       spp(ii)%vars%initLCFS(1:spp(ii)%pinit) = 1_is
       if (spp(ii)%pinit.lt.spp(ii)%ppp) then
          spp(ii)%vars%flagCon(spp(ii)%pinit+1:spp(ii)%ppp) = 0_is
          spp(ii)%vars%flagCol(spp(ii)%pinit+1:spp(ii)%ppp) = 0_is          
          spp(ii)%vars%flagRE(spp(ii)%pinit+1:spp(ii)%ppp) = 0_is
          spp(ii)%vars%initLCFS(spp(ii)%pinit+1:spp(ii)%ppp) = 0_is
       endif
       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'.or. &
               TRIM(params%field_model)=='M3D_C1'.or. &
               TRIM(params%field_model)=='NIMROD') 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(pinit)
    DEALLOCATE(q)
    DEALLOCATE(m)
    DEALLOCATE(Eno)
    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