subroutine initialize_collision_params(params,spp,P,F,init)
TYPE(KORC_PARAMS), INTENT(INOUT) :: params
TYPE(SPECIES), DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: spp
TYPE(PROFILES), INTENT(INOUT) :: P
TYPE(FIELDS), INTENT(IN) :: F
LOGICAL, INTENT(IN) :: init
INTEGER :: ii
REAL(rp) :: p_crit,gam_crit,maxEinterp,minEinterp
if (params%collisions.or.((TRIM(params%field_model).eq.'M3D_C1'.or. &
TRIM(params%field_model).eq.'NIMROD').and. &
params%radiation)) then
if (params%mpi_params%rank .EQ. 0) then
write(output_unit_write,'(/,"* * * * * * * INITIALIZING COLLISIONS * * * * * * *")')
end if
if (init) then
SELECT CASE (TRIM(params%collisions_model))
CASE (MODEL1)
call load_params_ss(params)
SELECT CASE(TRIM(params%bound_electron_model))
CASE ('NO_BOUND')
call load_params_ms(params)
cparams_ms%Ec=cparams_ss%Ec
if (.not.(P%ne_profile(1:6).eq.'RE-EVO')) then
cparams_ms%Ec_min=cparams_ms%Ec
else
cparams_ms%Ec_min=cparams_ms%Ec* &
cparams_ms%Gammac_min/cparams_ss%Gammac
endif
CASE('HESSLOW')
call load_params_ms(params)
cparams_ms%Ec=cparams_ss%Ec* &
(1._rp+sum((cparams_ms%Zo-cparams_ms%Zj)* &
cparams_ms%nz)/cparams_ss%ne)
if (.not.(P%ne_profile(1:6).eq.'RE-EVO')) then
cparams_ms%Ec_min=cparams_ms%Ec
else
cparams_ms%Ec_min=cparams_ms%Ec* &
cparams_ms%Gammac_min/cparams_ss%Gammac
end if
CASE('ROSENBLUTH')
call load_params_ms(params)
cparams_ms%Ec=cparams_ss%Ec* &
(1._rp+sum((cparams_ms%Zo-cparams_ms%Zj)* &
cparams_ms%nz)/cparams_ss%ne)
if (.not.(P%ne_profile(1:6).eq.'RE-EVO')) then
cparams_ms%Ec_min=cparams_ms%Ec
else
cparams_ms%Ec_min=cparams_ms%Ec* &
cparams_ms%Gammac_min/cparams_ss%Gammac
end if
CASE DEFAULT
write(output_unit_write,'("Default case")')
END SELECT
do ii=1_idef,params%num_species
ALLOCATE( spp(ii)%vars%nimp(spp(ii)%ppp, &
cparams_ms%num_impurity_species) )
spp(ii)%vars%nimp = 0.0_rp
end do
#ifdef FIO
if (TRIM(params%field_model) .eq. 'M3D_C1') then
call initialize_m3d_c1_imp(params,F,P, &
cparams_ms%num_impurity_species,.true.)
endif
#endif
CASE (MODEL2)
call load_params_ms(params)
CASE DEFAULT
write(output_unit_write,'("Default case")')
END SELECT
end if
if (params%LargeCollisions) then
if (params%mpi_params%rank .EQ. 0) then
write(output_unit_write,'(/,"* * * * * * * LARGE ANGLE COLLISIONS * * * * * * *")')
end if
if (TRIM(params%field_model) .eq. 'ANALYTICAL') then
!write(6,*) 'Eo',F%Eo
!write(6,*) 'Ec',cparams_ss%Ec
!write(6,*) 'Ec_min',cparams_ms%Ec_min
cparams_ss%avalanche=.TRUE.
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
if (abs(F%Eo).lt.cparams_ss%Ec) then
cparams_ss%avalanche=.FALSE.
end if
else
if (abs(F%Eo).lt.cparams_ms%Ec_min) then
cparams_ss%avalanche=.FALSE.
end if
end if
if (cparams_ss%avalanche) then
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
p_crit=1/sqrt(abs(F%Eo)/cparams_ss%Ec-1._rp)
else
p_crit=1/sqrt(abs(F%Eo)/cparams_ms%Ec_min-1._rp)
end if
end if
else if ((TRIM(params%field_model) .eq. 'EXTERNAL-PSI')) then
if (F%ReInterp_2x1t) then
maxEinterp=maxval(F%E_3D%PHI(:,F%ind_2x1t,:)* &
F%FLAG3D(:,F%ind_2x1t,:))
minEinterp=minval(F%E_3D%PHI(:,F%ind_2x1t,:)* &
F%FLAG3D(:,F%ind_2x1t,:))
else if (F%E_profile.eq.'MST_FSA') then
maxEinterp=F%E_dyn
end if
cparams_ss%avalanche=.TRUE.
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
if ((abs(maxEinterp).lt.cparams_ss%Ec).and. &
(abs(minEinterp).lt.cparams_ss%Ec)) &
cparams_ss%avalanche=.FALSE.
else
if ((abs(maxEinterp).lt.cparams_ms%Ec_min).and. &
(abs(minEinterp).lt.cparams_ms%Ec_min)) &
cparams_ss%avalanche=.FALSE.
end if
!write(6,*) 'maxEinterp',maxEinterp,'minEinterp',minEinterp, &
! 'E_c',cparams_ms%Ec_min,cparams_ss%avalanche
if (cparams_ss%avalanche) then
if (abs(maxEinterp).gt.abs(minEinterp)) then
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
p_crit=1/sqrt(abs(maxEinterp)/cparams_ss%Ec-1._rp)
else
p_crit=1/sqrt(abs(maxEinterp)/cparams_ms%Ec_min-1._rp)
end if
else
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
p_crit=1/sqrt(abs(minEinterp)/cparams_ss%Ec-1._rp)
else
p_crit=1/sqrt(abs(minEinterp)/cparams_ms%Ec_min-1._rp)
end if
end if
end if
else
write(6,*) 'Need to set p_crit!'
call korc_abort(25)
end if
if (cparams_ss%avalanche) then
cparams_ss%p_crit=p_crit
gam_crit=sqrt(1+p_crit*p_crit)
cparams_ss%gam_crit=gam_crit
cparams_ss%gam_therm=(gam_crit+1._rp)/2._rp
cparams_ss%p_therm=sqrt(cparams_ss%gam_therm*cparams_ss%gam_therm-1)
if(cparams_ss%min_secRE.eq.'THERM') then
cparams_ss%p_min=cparams_ss%p_therm
cparams_ss%gam_min=cparams_ss%gam_therm
else
cparams_ss%p_min=p_crit
cparams_ss%gam_min=gam_crit
end if
!write(6,*) p_crit,gam_crit,cparams_ss%p_therm,cparams_ss%gam_therm,cparams_ss%p_min,cparams_ss%gam_min
if (params%mpi_params%rank .EQ. 0) then
write(output_unit_write,*) 'Minimum energy of secondary RE is ',&
cparams_ss%min_secRE
write(output_unit_write,*) 'p_crit/(me*c) is: ',p_crit
write(output_unit_write,*) 'gam_min is: ',cparams_ss%gam_min
if(.not.init) then
if (TRIM(params%field_model) .eq. 'ANALYTICAL') then
write(output_unit_write,*) 'Maximum E_PHI : ',F%Eo*params%cpp%Eo,'V/m'
else if ((TRIM(params%field_model) .eq. 'EXTERNAL-PSI') &
.AND.(F%ReInterp_2x1t)) then
if (abs(maxEinterp).gt.abs(minEinterp)) then
write(output_unit_write,*) 'Maximum E_PHI : ',maxEinterp*params%cpp%Eo,'V/m'
else
write(output_unit_write,*) 'Maximum E_PHI : ',minEinterp*params%cpp%Eo,'V/m'
end if
endif
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
write(output_unit_write,*) 'E_CH is: ',cparams_ss%Ec*params%cpp%Eo,'V/m'
else
write(output_unit_write,*) 'E_CH is: ',cparams_ms%Ec_min*params%cpp%Eo,'V/m'
end if
write(output_unit_write,*) 'tau_c,rel is: ',cparams_ss%Tau*params%cpp%time,'s'
else
if (TRIM(params%field_model) .eq. 'ANALYTICAL') then
write(output_unit_write,*) 'Maximum E_PHI : ',F%Eo,'V/m'
else if ((TRIM(params%field_model) .eq. 'EXTERNAL-PSI') &
.AND.(F%ReInterp_2x1t)) then
if (abs(maxEinterp).gt.abs(minEinterp)) then
write(output_unit_write,*) 'Maximum E_PHI : ',maxEinterp,'V/m'
else
write(output_unit_write,*) 'Maximum E_PHI : ',minEinterp,'V/m'
end if
end if
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
write(output_unit_write,*) 'E_CH is: ',cparams_ss%Ec,'V/m'
else
write(output_unit_write,*) 'E_CH is: ',cparams_ms%Ec_min,'V/m'
end if
write(output_unit_write,*) 'tau_c,rel is: ',cparams_ss%Tau,'s'
end if
write(output_unit_write,'("* * * * * * * * * * * * * * * * * * * * * * * * * *",/)')
end if
else
if (params%mpi_params%rank .EQ. 0) then
if(.not.init) then
if (TRIM(params%field_model) .eq. 'ANALYTICAL') then
write(output_unit_write,*) 'Maximum E_PHI : ',F%Eo*params%cpp%Eo,'V/m'
else if ((TRIM(params%field_model) .eq. 'EXTERNAL-PSI') &
.AND.(F%ReInterp_2x1t)) then
if (abs(maxEinterp).gt.abs(minEinterp)) then
write(output_unit_write,*) 'Maximum E_PHI : ',maxEinterp*params%cpp%Eo,'V/m'
else
write(output_unit_write,*) 'Maximum E_PHI : ',minEinterp*params%cpp%Eo,'V/m'
end if
end if
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
write(output_unit_write,*) 'E_CH is: ',cparams_ss%Ec*params%cpp%Eo,'V/m'
else
write(output_unit_write,*) 'E_CH is: ',cparams_ms%Ec_min*params%cpp%Eo,'V/m'
end if
else
if (TRIM(params%field_model) .eq. 'ANALYTICAL') then
write(output_unit_write,*) 'Maximum E_PHI : ',F%Eo,'V/m'
else if ((TRIM(params%field_model) .eq. 'EXTERNAL-PSI') &
.AND.(F%ReInterp_2x1t)) then
if (abs(maxEinterp).gt.abs(minEinterp)) then
write(output_unit_write,*) 'Maximum E_PHI : ',maxEinterp,'V/m'
else
write(output_unit_write,*) 'Maximum E_PHI : ',minEinterp,'V/m'
end if
end if
if (TRIM(params%collisions_model).eq.'NO_BOUND') then
write(output_unit_write,*) 'E_CH is: ',cparams_ss%Ec,'V/m'
else
write(output_unit_write,*) 'E_CH is: ',cparams_ms%Ec_min,'V/m'
end if
end if
write(output_unit_write,*) 'No secondary REs will be calculated in this interval'
write(output_unit_write,*) 'p_min from initial or last time interval used'
write(output_unit_write,*) 'to calculate collision time scales'
write(output_unit_write,'("* * * * * * * * * * * * * * * * * * * * * * * * * *",/)')
end if
end if
end if
! for passing minimum gamma to Hollmann sampling routines
params%gam_min=sqrt(1+cparams_ss%p_min*cparams_ss%p_min* &
cparams_ss%pmin_scale*cparams_ss%pmin_scale)
!write(6,*) params%gam_min
if (params%mpi_params%rank .EQ. 0) then
write(output_unit_write,'("* * * * * * * * * * * * * * * * * * * * * * * * * *",/)')
end if
end if
end subroutine initialize_collision_params