Merge branch 'feature.thermionicEmission' into 'master.dev'

[feature.thermionicEmission] Thermionic emission modelling

See merge request piclas/piclas!658

docs/documentation/references.bib

@@ -777,3 +777,31 @@ @article{Zeng2020
   year={2020},
   publisher={AIP Publishing LLC}
 }
+
+@article{Coulombe1997,
+   author = {Sylvain Coulombe and Jean-Luc Meunier},
+   doi = {10.1088/0022-3727/30/5/009},
+   issn = {0022-3727},
+   issue = {5},
+   journal = {Journal of Physics D: Applied Physics},
+   month = {3},
+   pages = {776-780},
+   title = {Thermo-field emission: a comparative study},
+   volume = {30},
+   url = {https://iopscience.iop.org/article/10.1088/0022-3727/30/5/009},
+   year = {1997},
+}
+
+@article{Wu2022,
+   author = {Wu, Ernest Y. and Li, Baozhen},
+   doi = {10.1063/5.0087909},
+   issn = {0021-8979},
+   issue = {2},
+   journal = {Journal of Applied Physics},
+   month = {7},
+   pages = {025105},
+   title = {The Schottky emission effect: A critical examination of a century-old model},
+   volume = {132},
+   url = {https://aip.scitation.org/doi/10.1063/5.0087909},
+   year = {2022},
+}

docs/documentation/userguide/features-and-models/BC-particle-solver.md

@@ -1,3 +1,4 @@
+(sec:particle-boundary-conditions)=
 # Boundary Conditions - Particle Solver
 
 Within the parameter file it is possible to define different particle boundary conditions. The number of boundaries is defined by

docs/documentation/userguide/features-and-models/particle-initialization-and-emission.md

@@ -309,6 +309,26 @@ The velocity magnitude can be zero (per default) or a defined value (through `Ve
 be `open` or `reflective`. An example can be found in the regression test `regressioncheck/CHE_DSMC/SurfFlux_Tria_CurrentMassflow`
 For subsonic boundary conditions, where the velocity at the boundary is unknown, refer to Section {ref}`sec:particle-emission-adaptive`.
 
+### Thermionic Emission (including Schottky effect)
+
+The Richardson-Dushman equation including the Schottky effect is implemented and can be enabled to model thermionic emission
+
+$$j = A^* T_{\mathrm{w}} ^2 \exp\left(-\frac{W^*}{k_{\mathrm{B}} T_{\mathrm{w}}}\right), $$
+
+where the work function $W^*$ is defined by
+
+$$W^* = W - \Delta W \qquad \Delta W = \sqrt{\frac{q_{\mathrm{e}}^3 |\mathbf{E}|}{4\pi\epsilon_0}}.$$
+
+The magnitude of the electric field strength $|\mathbf{E}|$ is calculated with the average value of the interpolation points at the boundary. The material-specific properties such as the work function $W$ [eV] and the (modified) Richardson constant $A^*$ [A/cm²/K²] have to be provided as input. In addition to the surface flux parameters, a wall temperature $T_{\mathrm{w}}$ for the respective boundary has to be defined (as shown in Section {ref}`sec:particle-boundary-conditions`)
+
+    Part-Boundary1-WallTemp = 2000.
+    Part-Species1-Surfaceflux1-ThermionicEmission                    = TRUE
+    Part-Species1-Surfaceflux1-ThermionicEmission-SchottkyEffect     = TRUE
+    Part-Species1-Surfaceflux1-ThermionicEmission-WorkFunction       = 3
+    Part-Species1-Surfaceflux1-ThermionicEmission-RichardsonConstant = 120
+
+The provided temperature for the surface flux of the species determines the energy of emitted particles. While the thermionic emission can be enabled for PIC as well as DSMC simulations, the addition of the Schottky effect requires a field solver. An overview of the limitations of this modelling regarding the applied field strength, wall temperature and/or material is given by Ref. {cite}`Coulombe1997` and Ref. {cite}`Wu2022`. An example can be found in the regression test `regressioncheck/CHE_poisson/SurfFlux_ThermionicEmission_Schottky`.
+
 ### Circular Inflow
 
 The emission of particles from a surface flux can be limited to the area within a circle or a ring. The respective boundary has to
@@ -342,7 +362,7 @@ Multiple circular inflows can be defined on a single boundary through multiple s
 ### Adaptive/Subsonic Boundaries
 
 Different adaptive boundaries can be defined as a part of a surface flux to model subsonic in- and outflows, where the emission is
-adapted based on the prevalent conditions at the boundary. The modelling is based on the publications by {cite}`Farbar2014` and {cite}`Lei2017`.
+adapted based on the prevalent conditions at the boundary. The modelling is based on the publications by Ref. {cite}`Farbar2014` and Ref. {cite}`Lei2017`.
 
     Part-Species1-Surfaceflux1-Adaptive=TRUE
     Part-Species1-Surfaceflux1-Adaptive-Type=1

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/DSMC.ini

@@ -0,0 +1,8 @@
+! =============================================================================== !
+! Species1, e
+! =============================================================================== !
+Part-Species1-SpeciesName       = electron
+Part-Species1-InteractionID     = 4
+Part-Species1-Tref              = 273
+Part-Species1-dref              = 1E-15
+Part-Species1-omega             = 0.24

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/PartAnalyze_ref.csv

@@ -0,0 +1,12 @@
+001-TIME,002-Current-Spec-001-SF-001
+0.0000000000000000E+000,0.0000000000000000E+000
+0.1000000000000000E-009,0.1467593701480000E+001
+0.2000000000000000E-009,0.1467273266174000E+001
+0.3000000000000000E-009,0.1467593701480000E+001
+0.4000000000000000E-009,0.1467593701480000E+001
+0.5000000000000000E-009,0.1467593701480000E+001
+0.6000000000000000E-009,0.1467273266174000E+001
+0.7000000000000000E-009,0.1467593701480000E+001
+0.7999999999999999E-009,0.1467593701480000E+001
+0.8999999999999999E-009,0.1467273266174000E+001
+0.1000000000000000E-008,0.1467593701479998E+001

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/analyze.ini

@@ -0,0 +1,5 @@
+! compare the last line of PartAnalyze.csv with a reference file
+compare_data_file_name      = PartAnalyze.csv
+compare_data_file_reference = PartAnalyze_ref.csv
+compare_data_file_tolerance = 1e-3
+compare_data_file_tolerance_type = relative

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/command_line.ini

@@ -0,0 +1,2 @@
+MPI=1,2,4
+cmd_suffix=DSMC.ini

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/excludeBuild.ini

@@ -0,0 +1,2 @@
+! Skip N=1 hard coded
+PICLAS_POLYNOMIAL_DEGREE=1

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/hopr.ini

@@ -0,0 +1,39 @@
+ProjectName  = channel
+Debugvisu    = T
+DebugVisuLevel=1
+NVisu        =1
+Mode         =1
+
+DEFVAR = (REAL):   minus_x = 0.0
+DEFVAR = (REAL):   plus_x = 10.0
+
+DEFVAR = (REAL):   minus_y = 0.0
+DEFVAR = (REAL):   plus_y = 1.0
+
+DEFVAR = (REAL):   minus_z = 0.0
+DEFVAR = (REAL):   plus_z = 1.0
+
+Corner   =(/minus_x,minus_y,minus_z ,, plus_x,minus_y,minus_z ,, plus_x,plus_y,minus_z ,, minus_x,plus_y,minus_z ,, minus_x,minus_y,plus_z ,, plus_x,minus_y,plus_z ,, plus_x,plus_y,plus_z ,, minus_x,plus_y,plus_z /)
+nElems   =(/1,2,2/)
+elemtype =108
+
+BCIndex  =(/6 ,4 ,1 ,3 ,2 ,5/)
+!        =(/z-,y-,x+,y+,x-,z+/)
+nZones = 1
+nUserDefinedBoundaries=6
+BoundaryName=BC_Xplus
+BoundaryType=(/3,0,0,0/)
+BoundaryName=BC_Xminus
+BoundaryType=(/3,0,0,0/)
+BoundaryName=BC_Yplus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Yminus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Zplus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Zminus
+BoundaryType=(/4,0,0,0/)
+
+postscalemesh=true
+meshscale=1e-2
+jacobiantolerance=1e-20

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/parameter.ini

@@ -0,0 +1,86 @@
+! =============================================================================== !
+! EQUATION (linearscalaradvection)
+! =============================================================================== !
+IniExactFunc  = 0
+CFLscale   = 0.2  ! Scaling of theoretical CFL number
+! =============================================================================== !
+! DISCRETIZATION
+! =============================================================================== !
+N             = 1  ! Polynomial degree
+NAnalyze      = 1  ! Number of analyze points
+! =============================================================================== !
+! MESH
+! =============================================================================== !
+MeshFile      = channel_mesh.h5
+useCurveds    = F
+! if boundaries have to be changed (else they are used from Mesh directly):
+TrackingMethod = triatracking
+! =============================================================================== !
+! OUTPUT / VISUALIZATION
+! =============================================================================== !
+ProjectName      = SurfFlux_Tria_EmissionCurrent
+IterDisplayStep  = 10
+Part-AnalyzeStep = 1
+CalcSurfFluxInfo = T
+! =============================================================================== !
+! CALCULATION
+! =============================================================================== !
+tend            = 1.0E-9
+Analyze_dt      = 1.0E-9
+ManualTimeStep  = 1.0000E-10
+! =============================================================================== !
+! Load Balance
+! =============================================================================== !
+DoLoadBalance           = T
+DoInitialAutoRestart    = T
+Load-DeviationThreshold = 1e-9
+! =============================================================================== !
+! PARTICLES
+! =============================================================================== !
+Part-maxParticleNumber=500000
+Part-nSpecies=1
+Part-nBounds=6
+Part-Boundary1-SourceName = BC_Xplus
+Part-Boundary1-Condition  = open
+Part-Boundary2-SourceName = BC_Xminus
+Part-Boundary2-Condition  = reflective
+Part-Boundary2-WallTemp   = 2700.
+Part-Boundary3-SourceName = BC_Yplus
+Part-Boundary3-Condition  = reflective
+Part-Boundary4-SourceName = BC_Yminus
+Part-Boundary4-Condition  = reflective
+Part-Boundary5-SourceName = BC_Zplus
+Part-Boundary5-Condition  = symmetric
+Part-Boundary6-SourceName = BC_Zminus
+Part-Boundary6-Condition  = symmetric
+Part-FIBGMdeltas=(/1e-2,1e-2,1e-2/)
+! =============================================================================== !
+! Species1 - electron
+! =============================================================================== !
+Part-Species1-MassIC              = 9.11E-31
+Part-Species1-ChargeIC            = -1.60217653E-19
+Part-Species1-MacroParticleFactor = 2E5
+
+Part-Species1-nSurfaceFluxBCs=1
+Part-Species1-Surfaceflux1-BC=2
+
+Part-Species1-Surfaceflux1-VeloIC = 0.
+Part-Species1-Surfaceflux1-VeloVecIC = (/1,0,0/)
+Part-Species1-Surfaceflux1-velocityDistribution = maxwell_lpn
+Part-Species1-Surfaceflux1-MWTemperatureIC = 300
+Part-Species1-Surfaceflux1-ThermionicEmission = TRUE
+Part-Species1-Surfaceflux1-ThermionicEmission-WorkFunction = 4.54
+Part-Species1-Surfaceflux1-ThermionicEmission-RichardsonConstant = 60
+! =============================================================================== !
+! DSMC
+! =============================================================================== !
+Particles-HaloEpsVelo=2.0E+06
+UseDSMC=true
+Particles-DSMC-CollisMode=1 !(1:elast coll, 2: elast + rela, 3:chem)
+Part-NumberOfRandomSeeds=2
+Particles-RandomSeed1=5
+Particles-RandomSeed2=6
+Particles-DSMC-CalcSurfaceVal      = F
+Particles-DSMC-UseOctree           = F
+Particles-DSMC-UseNearestNeighbour = F
+Particles-DSMC-CalcQualityFactors  = F

regressioncheck/CHE_DSMC/SurfFlux_Tria_ThermionicEmission/readme.md

@@ -0,0 +1,5 @@
+# DSMC - Thermionic Emission
+* Particle emission using surface flux and thermionic emission
+* Comparing the calculated current to the expected value of the Richardson Dushman equation for Tungsten
+  * Input: W = 4.54 eV, A = 60 A/(cm^2 K^2), T_w = 2700 K
+  * Output: j = 1.47 A / cm^2 -> I = 1.4675 eV (A = 1 cm^2)

regressioncheck/CHE_poisson/SurfFlux_ThermionicEmission_Schottky/PartAnalyze_ref.csv

@@ -0,0 +1,12 @@
+001-TIME,002-Current-Spec-001-SF-001
+0.0000000000000000E+000,0.0000000000000000E+000
+0.1000000000000000E-011,0.1467593701480000E+001
+0.2000000000000000E-011,0.1467433483827000E+001
+0.3000000000000000E-011,0.1467433483827000E+001
+0.4000000000000000E-011,0.1467433483827000E+001
+0.5000000000000000E-011,0.1467593701480000E+001
+0.5999999999999999E-011,0.1467593701480000E+001
+0.6999999999999999E-011,0.1467433483827000E+001
+0.8000000000000000E-011,0.1467593701480000E+001
+0.9000000000000000E-011,0.1467433483827000E+001
+0.9999999999999999E-011,0.1467593701480000E+001

regressioncheck/CHE_poisson/SurfFlux_ThermionicEmission_Schottky/analyze.ini

@@ -0,0 +1,5 @@
+! compare the last line of PartAnalyze.csv with a reference file
+compare_data_file_name      = PartAnalyze.csv
+compare_data_file_reference = PartAnalyze_ref.csv
+compare_data_file_tolerance = 5e-3
+compare_data_file_tolerance_type = relative

regressioncheck/CHE_poisson/SurfFlux_ThermionicEmission_Schottky/command_line.ini

@@ -0,0 +1 @@
+MPI=1,2,4

regressioncheck/CHE_poisson/SurfFlux_ThermionicEmission_Schottky/hopr.ini

@@ -0,0 +1,39 @@
+ProjectName  = channel
+Debugvisu    = T
+DebugVisuLevel=1
+NVisu        =1
+Mode         =1
+
+DEFVAR = (REAL):   minus_x = 0.0
+DEFVAR = (REAL):   plus_x = 10.0
+
+DEFVAR = (REAL):   minus_y = 0.0
+DEFVAR = (REAL):   plus_y = 1.0
+
+DEFVAR = (REAL):   minus_z = 0.0
+DEFVAR = (REAL):   plus_z = 1.0
+
+Corner   =(/minus_x,minus_y,minus_z ,, plus_x,minus_y,minus_z ,, plus_x,plus_y,minus_z ,, minus_x,plus_y,minus_z ,, minus_x,minus_y,plus_z ,, plus_x,minus_y,plus_z ,, plus_x,plus_y,plus_z ,, minus_x,plus_y,plus_z /)
+nElems   =(/10,2,2/)
+elemtype =108
+
+BCIndex  =(/6 ,4 ,1 ,3 ,2 ,5/)
+!        =(/z-,y-,x+,y+,x-,z+/)
+nZones = 1
+nUserDefinedBoundaries=6
+BoundaryName=BC_Xplus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Xminus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Yplus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Yminus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Zplus
+BoundaryType=(/4,0,0,0/)
+BoundaryName=BC_Zminus
+BoundaryType=(/4,0,0,0/)
+
+postscalemesh=true
+meshscale=1e-2
+jacobiantolerance=1e-20