colrad Subroutine

  1. fidasim.f90
  2. libfida
  3. colrad

public subroutine colrad(plasma, ab, vn, dt, states, dens, photons)

Evolve density of states in time dt via collisional radiative model

Arguments

TypeIntentOptionalAttributesName
type(LocalProfiles), intent(in) :: plasma

Plasma parameters
real(kind=Float64), intent(in) :: ab

Ion/Neutral mass [amu]
real(kind=Float64), intent(in), dimension(:):: vn

Neutral velocitiy [cm/s]
real(kind=Float64), intent(in) :: dt

Time interval [s]
real(kind=Float64), intent(inout), dimension(:):: states

Density of states
real(kind=Float64), intent(out), dimension(nlevs):: dens

Density of neutrals
real(kind=Float64), intent(out) :: photons

Emitted photons

Calls

proc~~colrad~~CallsGraph proc~colrad colrad proc~get_rate_matrix get_rate_matrix proc~colrad->proc~get_rate_matrix proc~linsolve linsolve proc~colrad->proc~linsolve proc~eigen eigen proc~colrad->proc~eigen interface~interpol_coeff interpol_coeff proc~get_rate_matrix->interface~interpol_coeff dgetrs dgetrs proc~linsolve->dgetrs proc~matinv matinv proc~linsolve->proc~matinv dgetrf dgetrf proc~linsolve->dgetrf proc~elmhes elmhes proc~eigen->proc~elmhes proc~balance balance proc~eigen->proc~balance proc~hqr2 hqr2 proc~eigen->proc~hqr2 proc~balback balback proc~eigen->proc~balback proc~elmtrans elmtrans proc~eigen->proc~elmtrans proc~rswap RSWAP proc~elmhes->proc~rswap proc~lubksb lubksb proc~matinv->proc~lubksb proc~ludcmp ludcmp proc~matinv->proc~ludcmp proc~balance->proc~rswap proc~hqrvec hqrvec proc~hqr2->proc~hqrvec proc~balback->proc~rswap proc~comdiv Comdiv proc~hqrvec->proc~comdiv proc~outerprod outerprod proc~ludcmp->proc~outerprod proc~swap swap proc~ludcmp->proc~swap

Called by

proc~~colrad~~CalledByGraph proc~colrad colrad proc~ndmc ndmc proc~ndmc->proc~colrad proc~fida_weights_los fida_weights_los proc~fida_weights_los->proc~colrad proc~dcx dcx proc~dcx->proc~colrad proc~fida_f fida_f proc~fida_f->proc~colrad proc~pfida_f pfida_f proc~pfida_f->proc~colrad proc~read_equilibrium read_equilibrium proc~read_equilibrium->proc~colrad proc~attenuate attenuate proc~attenuate->proc~colrad proc~halo halo proc~halo->proc~colrad proc~pfida_mc pfida_mc proc~pfida_mc->proc~colrad proc~fida_weights_mc fida_weights_mc proc~fida_weights_mc->proc~colrad proc~fida_mc fida_mc proc~fida_mc->proc~colrad program~fidasim fidasim program~fidasim->proc~ndmc program~fidasim->proc~fida_weights_los program~fidasim->proc~dcx program~fidasim->proc~fida_f program~fidasim->proc~pfida_f program~fidasim->proc~read_equilibrium program~fidasim->proc~halo program~fidasim->proc~pfida_mc program~fidasim->proc~fida_weights_mc program~fidasim->proc~fida_mc proc~pnpa_mc pnpa_mc program~fidasim->proc~pnpa_mc proc~pnpa_f pnpa_f program~fidasim->proc~pnpa_f proc~npa_mc npa_mc program~fidasim->proc~npa_mc proc~npa_f npa_f program~fidasim->proc~npa_f proc~npa_weights npa_weights program~fidasim->proc~npa_weights proc~pnpa_mc->proc~attenuate proc~pnpa_f->proc~attenuate proc~npa_mc->proc~attenuate proc~npa_f->proc~attenuate proc~npa_weights->proc~attenuate

Contents

Source Code

colrad

Source Code

subroutine colrad(plasma,ab,vn,dt,states,dens,photons)
    !+ Evolve density of states in time `dt` via collisional radiative model
    type(LocalProfiles), intent(in)              :: plasma
        !+ Plasma parameters
    real(Float64), intent(in)                    :: ab
        !+ Ion/Neutral mass [amu]
    real(Float64), dimension(:), intent(in)      :: vn
        !+ Neutral velocitiy [cm/s]
    real(Float64), intent(in)                    :: dt
        !+ Time interval [s]
    real(Float64), dimension(:), intent(inout)   :: states
        !+ Density of states
    real(Float64), dimension(nlevs), intent(out) :: dens
        !+ Density of neutrals
    real(Float64), intent(out)                   :: photons
        !+ Emitted photons

    real(Float64), dimension(nlevs,nlevs) :: matrix  !! Matrix
    real(Float64) :: vnet_square    !! net velocity of neutrals squared
    real(Float64) :: eb             !! Energy of the fast neutral

    real(Float64), dimension(nlevs,nlevs) :: eigvec
    real(Float64), dimension(nlevs) :: eigval, coef
    real(Float64), dimension(nlevs) :: exp_eigval_dt
    real(Float64) :: iflux !!Initial total flux
    integer :: n

    photons=0.d0
    dens=0.d0

    iflux=sum(states)

    if(.not.plasma%in_plasma) then
        dens = states*dt
        return
    endif

    vnet_square=dot_product(vn-plasma%vrot,vn-plasma%vrot)  ![cm/s]
    eb = v2_to_E_per_amu*ab*vnet_square ![kev]
    call get_rate_matrix(plasma, ab, eb, matrix)

    call eigen(nlevs,matrix, eigvec, eigval)
    call linsolve(eigvec,states,coef) !coeffs determined from states at t=0
    exp_eigval_dt = exp(eigval*dt)   ! to improve speed (used twice)
    do n=1,nlevs
        if(eigval(n).eq.0.0) eigval(n)=eigval(n)+1 !protect against dividing by zero
    enddo

    states = matmul(eigvec, coef * exp_eigval_dt)  ![neutrals/cm^3/s]!
    dens   = matmul(eigvec,coef*(exp_eigval_dt-1.d0)/eigval)

    where (states.lt.0)
        states = 0.d0
    endwhere

    where (dens.lt.0)
        dens = 0.d0
    endwhere

    photons=dens(initial_state)*tables%einstein(final_state,initial_state) !! - [Ph/(s*cm^3)] - !!

end subroutine colrad