subroutine PreLoopSdot(TimePreMax) use particles use scattering use simulation implicit none integer step, count_bcCheck, count_Poisson, i, j, bin, err1 ! add Tempcheck 11/2012 real ( kind = 8 ):: PreLoopPCalTime, TimePre,TimePreEarly,TimePreMax,PreLoop_Time_period, factor real ( kind = 8 ):: PreLoop_T_eNum,PreLoop_T_vel,PreLoop_T_eNum_eb,PreLoop_T_vel_eb,PreLoop_BoundFactor real ( kind = 8 ), dimension(:), allocatable :: Temp_bin_ParDrftVel integer, dimension(:), allocatable :: Temp_bin_ParNumber allocate ( Temp_bin_ParNumber(num_Bins) ) allocate ( Temp_bin_ParDrftVel(num_Bins) ) PreLoop_A_bin_eNum(1:num_Bins) = 0.0d0 PreLoop_A_bin_vel(1:num_Bins) = 0.0d0 Temp_bin_ParDrftVel(1:num_Bins) = 0.0d0 Temp_bin_ParNumber(1:num_Bins) = 0 TimePre = 0. ! in simulation module TimePreEarly = TimePreMax*0.5 step = 0 ! parameter setting eeHPABPre1 count_bcCheck = 0 count_Poisson = 0 PreLoop_Time_period = 0.0 do while(TimePre.le.TimePreMax) ! to find step = step + 1 TimePre = time_step*real(step,kind = 8) StepParOutL = 0 ! check Particles out StepParOutR = 0 ! do i = 1, num_maxCarriers if (Par_use(i).eq.1) then if ((Par_use(i).le.MaxE_eV).and.(Par_use(i).gt.0)) then err1 = 1 ! for debugging endif endif enddo call free_flight_scatter() ! Normal Using the minimum temperature count_bcCheck = count_bcCheck + 1 ! for boundary treat... AccumInjParL = AccumInjParL + num_InjLB ! elseif (Algo_bcCheck.eq.3) then AccumInjPar = AccumInjPar + InjParPerStep AccumInjParR = AccumInjParR + num_InjRB if (count_bcCheck.ge.num_CheckFreq) then count_bcCheck = 0 call BoundaryTreat() endif count_Poisson = count_Poisson + 1 ! Poisson Field Update -PoissonSetting, Algo_AssignCharge, Step_FreqFieldUpdate, time_Poisson, count_Poisson if ((PoissonSetting.eq.1).and.(time.ge.time_Poisson)) then if (count_Poisson.ge.Step_FreqFieldUpdate) then count_Poisson = 0 call ChargeAssign() ! call PotFromPoisson() ! endif endif if (TimePre.gt.TimePreEarly) then ! 0720 Start . PreLoop_Time_period = PreLoop_Time_period + 1.0 bin_ParNumber(:) = 0 ! call CalculateProperties() ! ! bin_ParDrftVel(:) = 0.0d0 num_Carriers = 0 do i = 1,num_maxCarriers if (Par_use(i).eq.1) then ! used particle num_Carriers = num_Carriers + 1 factor = 1./(1.+nonParabolPara2(Par_Al(i),Par_valley(i))*Par_ene(i)) do j = 1, DIM Par_vel(i,j) = hbar*Par_wvk(i,j)*factor/eff_mass(Par_Al(i),Par_valley(i)) enddo bin = num_Bins do j = 1,num_Bins ! ! Carrier density profile bin profile if ((Par_pos(i,1).ge.(loc_bin(j)-0.5*BinSize)).and.(Par_pos(i,1).lt.(loc_bin(j)+0.5*BinSize))) then bin = j endif enddo bin_ParNumber(bin) = bin_ParNumber(bin) + 1 bin_ParDrftVel(bin) = bin_ParDrftVel(bin) + Par_vel(i,1) endif enddo ! num_Carriers do j = 1,num_Bins ! if (bin_ParNumber(j).gt.0) then bin_ParDrftVel(j) = bin_ParDrftVel(j)/real(bin_ParNumber(j),kind = 8) else bin_ParDrftVel(j) = 0.0d0 endif enddo if (PropertyNode.eq.1) then ! if properties for node center bin_ParNumber(1) = 2*bin_ParNumber(1) bin_ParNumber(num_Bins) = 2*bin_ParNumber(num_Bins) endif do i = 1,num_Bins Temp_bin_ParNumber(i) = PreLoop_A_bin_eNum(i)*(PreLoop_Time_period-1.0) + real(bin_ParNumber(i),kind = 8) Temp_bin_ParDrftVel(i) = PreLoop_A_bin_vel(i)*(PreLoop_Time_period-1.0) + bin_ParDrftVel(i) PreLoop_A_bin_eNum(i) = Temp_bin_ParNumber(i)/PreLoop_Time_period PreLoop_A_bin_vel(i) = Temp_bin_ParDrftVel(i)/PreLoop_Time_period enddo else BinScatCount(:,:,:) = 0 endif enddo ! End of the time loop PreLoop_T_eNum = 0.0 PreLoop_T_vel = 0.0 PreLoop_T_eNum_eb = 0.0 PreLoop_T_vel_eb = 0.0 do i = 1, num_Bins PreLoop_BoundFactor = 1.0 if ((i.eq.1).or.(i.eq.num_Bins)) then PreLoop_BoundFactor = 0.5 endif PreLoop_T_eNum = PreLoop_T_eNum + PreLoop_A_bin_eNum(i)*PreLoop_BoundFactor PreLoop_T_vel = PreLoop_T_vel + PreLoop_A_bin_eNum(i)*PreLoop_A_bin_vel(i)*PreLoop_BoundFactor if ((i.gt.num_boundBin).and.(i.le.(num_Bins-num_boundBin))) then PreLoop_T_eNum_eb = PreLoop_T_eNum_eb + PreLoop_A_bin_eNum(i) PreLoop_T_vel_eb = PreLoop_T_vel_eb + PreLoop_A_bin_eNum(i)*PreLoop_A_bin_vel(i) endif enddo PreLoop_eNum = PreLoop_T_eNum PreLoop_vel = PreLoop_T_vel/PreLoop_eNum PreLoop_eNum_eb = PreLoop_T_eNum_eb PreLoop_vel_eb = PreLoop_T_vel_eb/PreLoop_eNum_eb PreLoop_je = n_d*PreLoop_eNum*PreLoop_vel*e_c*0.0001/real(num_Carriers,kind = 8) PreLoop_jeA = PreLoop_je*box_len(1) if (PoissonSetting.eq.1) then PreLoop_je_eb = n_d*PreLoop_eNum_eb*PreLoop_vel_eb*e_c*0.0001/real((num_Carriers-(num_boundBin+num_boundBin-1)*num_InitBin),kind = 8) PreLoop_jeA_eb = PreLoop_je_eb*(box_len(1)-BinSize*real((num_boundBin+num_boundBin-1),kind = 8)) else PreLoop_je_eb = n_d*PreLoop_eNum_eb*PreLoop_vel_eb*e_c*0.0001/real((num_Carriers-(num_boundBin+num_boundBin)*num_InitBin),kind = 8) PreLoop_jeA_eb = PreLoop_je_eb*(box_len(1)-BinSize*real((num_boundBin+num_boundBin),kind = 8)) endif PreLoopPCalTime = time_step*PreLoop_Time_period call PreLoopPhononCal(PreLoopPCalTime) ! final open (unit = 61, file = 'PreLoopOut.txt', status = 'replace', action = 'write' ) write(61,*)'Preloop for Sdot ','PreLoopOut.txt' write(61,'(a,a,a,a)'), 'Setting: ', 'TimePreMax ', 'averageDuration ', 'TimePreEarly ' write(61,'(2x,10e20.9)'),TimePreMax,PreLoop_Time_period,TimePreEarly write(61,'(a,a,a,a)'), 'Total: ', 'PreLoop_eNum, PreLoop_vel ', 'PreLoop_je, PreLoop_jeA ' , ' PreLoop_TP_EneWm2 PreLoop_TPLO_EneWm2 PreLoop_TPTO_EneWm2' write(61,'(2x,4f20.8,10e20.9)'), PreLoop_eNum, PreLoop_vel, PreLoop_je, PreLoop_jeA, PreLoop_TP_EneWm2,PreLoop_TPLO_EneWm2,PreLoop_TPTO_EneWm2 write(61,'(a,a,a,a)'), 'ExceptBound: ', 'PreLoop_eNum_eb, PreLoop_vel_eb ', 'PreLoop_je_eb, PreLoop_jeA_eb ' , 'PreLoop_TP_EneWm2_eb PreLoop_TPLO_EneWm2_eb PreLoop_TPTO_EneWm2_eb' write(61,'(2x,4f20.8,10e20.9)'), PreLoop_eNum_eb, PreLoop_vel_eb, PreLoop_je_eb, PreLoop_jeA_eb,PreLoop_TP_EneWm2_eb,PreLoop_TPLO_EneWm2_eb,PreLoop_TPTO_EneWm2_eb write(61,'(a)'), '----------------------------------------------------------------------------------------------------------------------------------------' write(61,'(a,a)'), 'Pos ', ' PreLoop_A_bin_eNum PreLoop_A_bin_vel PreLoop_Bin_EpOWm3 PreLoop_Bin_EpLOWm3 PreLoop_Bin_EpTOWm3 PreLoop_Bin_EpLOWm2 PreLoop_Bin_EpTOWm2' write(61,'(a)'), '----------------------------------------------------------------------------------------------------------------------------------------' do i = 1, num_Bins write(61,'(1x,e20.9,2f20.8,10e20.9)'),loc_bin(i),PreLoop_A_bin_eNum(i),PreLoop_A_bin_vel(i),PreLoop_Bin_EpOWm3(i),PreLoop_Bin_EpLOWm3(i),PreLoop_Bin_EpTOWm3(i),PreLoop_Bin_EpLOWm2(i),PreLoop_Bin_EpTOWm2(i) enddo close(61) return end subroutine PreLoopPhononCal(PreLoop_Time_period) ! LO phonon use simulation use scattering use particles ! real ( kind = 8 ), dimension(num_Valley,DIM) :: TempTotValleyVel integer:: j !num_ScatMech,num_Valley real ( kind = 8 )::PreLoop_Time_period PreLoop_Bin_EpLO(:) = 0.0d0 PreLoop_Bin_EpTO(:) = 0.0d0 PreLoop_Bin_EpO(:) = 0.0d0 PreLoop_Bin_EpLOWm3(:) = 0.0d0 PreLoop_Bin_EpTOWm3(:) = 0.0d0 PreLoop_Bin_EpOWm3(:) = 0.0d0 PreLoop_Bin_EpLOWm2(:) = 0.0d0 PreLoop_Bin_EpTOWm2(:) = 0.0d0 PreLoop_Bin_EpOWm2(:) = 0.0d0 PreLoop_TP_EneWm2 = 0.0d0 PreLoop_TP_EneWm2_eb = 0.0d0 PreLoop_TPLO_EneWm2 = 0.0d0 PreLoop_TPLO_EneWm2_eb = 0.0d0 PreLoop_TPTO_EneWm2 = 0.0d0 PreLoop_TPTO_EneWm2_eb = 0.0d0 do j = 1,num_Bins PreLoop_Bin_EpLO(j)= PreLoop_Bin_EpLO(j)+real((BinScatCount(j,2,1)-BinScatCount(j,3,1)),kind = 8)*polar_Ep_G+real((BinScatCount(j,2,2)-BinScatCount(j,3,2)),kind = 8)*polar_Ep_L+real((BinScatCount(j,2,3)-BinScatCount(j,3,3)),kind = 8)*polar_Ep_X PreLoop_Bin_EpTO(j)= PreLoop_Bin_EpTO(j)+real((BinScatCount(j,4,1)-BinScatCount(j,5,1)),kind = 8)*IntVal_Ep_G_L+real((BinScatCount(j,6,1)-BinScatCount(j,7,1)),kind = 8)*IntVal_Ep_G_X PreLoop_Bin_EpTO(j)= PreLoop_Bin_EpTO(j)+real((BinScatCount(j,4,2)-BinScatCount(j,5,2)),kind = 8)*IntVal_Ep_L_G+real((BinScatCount(j,6,2)-BinScatCount(j,7,2)),kind = 8)*IntVal_Ep_L_L+real((BinScatCount(j,8,2)-BinScatCount(j,9,2)),kind = 8)*IntVal_Ep_L_X PreLoop_Bin_EpTO(j)= PreLoop_Bin_EpTO(j)+real((BinScatCount(j,4,3)-BinScatCount(j,5,3)),kind = 8)*IntVal_Ep_X_G+real((BinScatCount(j,6,3)-BinScatCount(j,7,3)),kind = 8)*IntVal_Ep_X_L+real((BinScatCount(j,8,3)-BinScatCount(j,9,3)),kind = 8)*IntVal_Ep_X_X PreLoop_Bin_EpO(j) = PreLoop_Bin_EpLO(j) + PreLoop_Bin_EpTO(j) PreLoop_Bin_EpLOWm3(j) = (PreLoop_Bin_EpLO(j)*e_c)*n_d/real(num_InitNode,kind = 8)/PreLoop_Time_period ! J/s- BinSize PreLoop_Bin_EpTOWm3(j) = (PreLoop_Bin_EpTO(j)*e_c)*n_d/real(num_InitNode,kind = 8)/PreLoop_Time_period ! J/s- BinSize PreLoop_Bin_EpOWm3(j) = PreLoop_Bin_EpLOWm3(j) + PreLoop_Bin_EpTOWm3(j) if ((j.eq.1).or.(j.eq.num_Bins)) then PreLoop_Bin_EpLOWm2(j) = PreLoop_Bin_EpLOWm3(j)*BinSize*0.5 PreLoop_Bin_EpTOWm2(j) = PreLoop_Bin_EpTOWm3(j)*BinSize*0.5 else PreLoop_Bin_EpLOWm2(j) = PreLoop_Bin_EpLOWm3(j)*BinSize PreLoop_Bin_EpTOWm2(j) = PreLoop_Bin_EpTOWm3(j)*BinSize endif PreLoop_Bin_EpOWm2(j) = PreLoop_Bin_EpLOWm2(j) + PreLoop_Bin_EpTOWm2(j) PreLoop_TP_EneWm2 = PreLoop_TP_EneWm2 + PreLoop_Bin_EpOWm2(j) PreLoop_TPLO_EneWm2 = PreLoop_TPLO_EneWm2 + PreLoop_Bin_EpLOWm2(j) PreLoop_TPTO_EneWm2 = PreLoop_TPTO_EneWm2 + PreLoop_Bin_EpTOWm2(j) ! PreLoop_TPTO_Ene = PreLoop_TPTO_Ene + PreLoop_Bin_EpTO(j) if ((j.gt.num_boundBin).and.(j.le.(num_Bins-num_boundBin))) then PreLoop_TP_EneWm2_eb = PreLoop_TP_EneWm2_eb + PreLoop_Bin_EpOWm2(j) PreLoop_TPLO_EneWm2_eb = PreLoop_TPLO_EneWm2_eb + PreLoop_Bin_EpLOWm2(j) PreLoop_TPTO_EneWm2_eb = PreLoop_TPTO_EneWm2_eb + PreLoop_Bin_EpTOWm2(j) endif enddo return end subroutine InitGeneralParas() use simulation use scattering use particles ! real ( kind = 8 ), dimension(num_Valley,DIM) :: TempTotValleyVel integer:: i, j AccumParInNumBoundL = 0 AccumParOutNumBoundL = 0 AccumParOutNumBoundR = 0 AccumParNumBoundPlus = 0 AccumEneBoundL = 0.0d0 AccumEneBoundR = 0.0d0 AccumEneBound = 0.0d0 AccumInjParL = 0.0d0 AccumInjParR = 0.0d0 EneBoundL = 0.0d0 EneBoundR = 0.0d0 EneBound = 0.0d0 time = 0.0 time_fave_period = 0.0 bin_ParNumber(:) = 0 initbin_ParAveEne(:) = 0.0d0 bin_ParNumDen(:) = 0.0d0 bin_ParAveEne(:) = 0.0d0 bin_ParAveEneV(:) = 0.0d0 bin_ParDrftVel(:) = 0.0d0 bin_ParVal(:,:) = 0.0d0 binhist_ParNumDen(:,:) = 0.0d0 binhist_ParNumber(:,:) = 0 ParInNumBoundR = 0 ParInNumBoundL = 0 ParOutNumBoundL = 0 ParOutNumBoundR = 0 ParNumBoundPlus = 0 ! initialize Par_wvk(:,:) = 0.0d0 Par_pos(:,:) = 0.0d0 Par_vel(:,:) = 0.0d0 Par_valley(:) = 0 Par_use(:) = 0 Par_ene(:) = 0.0D+00 Par_tau(:) = 0.0D+00 TA_Node_Potential(:) = 0.0d0 TA_Node_PoiPot(:) = 0.0d0 TA_Node_Charge(:) = 0.0d0 TA_Node_EField(:) = 0.0d0 TA_bin_ParNumDen(:) = 0.0d0 TA_bin_ParAveEne(:) = 0.0d0 TA_bin_ParVal(:,:) = 0.0d0 TA_bin_ParNumber(:) = 0.0d0 TA_bin_ParDrftVel(:) = 0.0d0 TA_bin_PhononNum(:) = 0 TA_bin_PhononEne(:) = 0.0d0 TA_bin_PhononEneWm(:) = 0.0d0 TA_binhist_ParNumDen(:,:) = 0.0d0 TA_binhist_ParNumber(:,:) = 0.0d0 TABin_TemperAP(:) = 0.0d0 TABin_TemperOP(:) = 0.0d0 TABinQsupA(:) = 0.0d0 TABinQsupO(:) = 0.0d0 TABinQsupTotal(:) = 0.0d0 TAepInterEneWm2(:) = 0.0d0 TAppInterEneWm2(:) = 0.0d0 fave_bin_TpA(:) = 0.0d0 fave_bin_TpO(:) = 0.0d0 fave_bin_qsupA(:) = 0.0d0 fave_bin_qsupO(:) = 0.0d0 fave_bin_qsupT(:) = 0.0d0 fave_bin_epWm2(:) = 0.0d0 fave_bin_ppWm2(:) = 0.0d0 Tfave_InOut_B = 0.0 Tfave_InOut_LB = 0.0 Tfave_InOut_RB = 0.0 Tfave_G_ratio = 0.0 Tfave_bin_eNum(:) = 0.0 Tfave_bin_vel(:) = 0.0 Tfave_bin_PhonEWm(:) = 0.0 Tfave_bin_PoiPot(:) = 0.0 Tfave_bin_AllPot(:) = 0.0 Tfave_bin_Eneev(:) = 0.0 Tfave_bin_TpA(:) = 0.0d0 Tfave_bin_TpO(:) = 0.0d0 Tfave_bin_qsupA(:) = 0.0d0 Tfave_bin_qsupO(:) = 0.0d0 Tfave_bin_qsupT(:) = 0.0d0 Tfave_bin_epWm2(:) = 0.0d0 Tfave_bin_ppWm2(:) = 0.0d0 SysCurrent = 0.0d0 AllScatCount(:,:) = 0 ! 0709 Phonon BinScatCount(:,:,:) = 0 TempVar_BinScatCount(:,:,:) = 0 TempAddEneY = 0 do i = 1, num_Nodes loc_Node(i) = box_coor(1,1) + real(i,kind = 8)*NodeSize - NodeSize ! bin center position SecLow, SecMid, TempAddEneY pot_Node(i) = 0.0d0 ! additional Charge charge_Node(i) = 0.0d0 ! initial .. doping concentration field_Node(i) = 0.0d0 num_Node(i) = real(num_InitNode,kind = 8) if (num_BarTros.ge.1) then do j = 1, num_BarTros if (Par_pos(i,1).ge.loc_Bar(j)) then TempAddEneY = TempAddEneY + Bar_Height(j) endif enddo endif if (num_FieldSec.ge.1) then if (Par_pos(i,1).lt.loc_Start_Field(1)) then TempAddEneY = TempAddEneY + (Par_pos(i,1)-box_coor(1,1))*EFieldBg(1) else TempAddEneY = TempAddEneY + (loc_Start_Field(1)-box_coor(1,1))*EFieldBg(1) do j = 1, num_FieldSec if (Par_pos(i,1).ge.loc_Start_Field(j)) then if (Par_pos(i,1).ge.loc_End_Field(j)) then TempAddEneY = TempAddEneY + (loc_End_Field(j)-loc_Start_Field(j))*Field_Sec(j) else TempAddEneY = TempAddEneY + (Par_pos(i,1)-loc_Start_Field(j))*Field_Sec(j) endif endif enddo if (Par_pos(i,1).gt.loc_End_Field(num_FieldSec)) then TempAddEneY = TempAddEneY + (Par_pos(i,1)-loc_End_Field(num_FieldSec))*EFieldBg(1) endif endif else TempAddEneY = TempAddEneY + (Par_pos(i,1)-box_coor(1,1))*EFieldBg(1) endif InitPot_Node(i) = TempAddEneY enddo ! loc_Node, pot_Node, charge_Node num_Carriers = Initialnum_Carriers call initial_distribution() return end