From e1009a21df141e89f758eb2f57656886bf69ea1e Mon Sep 17 00:00:00 2001 From: Jorge Gonzalez Date: Tue, 18 Jul 2023 15:35:16 +0200 Subject: [PATCH] Correction due to angle Moving forward to a complete model for SEE. Now the yield has a correction for the incident angle on the surface (still to check if this is correct) Velocity distribution is still a dummy one. --- src/modules/common/moduleTable.f90 | 3 + src/modules/mesh/moduleMeshBoundary.f90 | 126 ++++++++++++++---------- src/modules/moduleBoundary.f90 | 12 +++ 3 files changed, 89 insertions(+), 52 deletions(-) diff --git a/src/modules/common/moduleTable.f90 b/src/modules/common/moduleTable.f90 index e4da335..ac5cbc7 100644 --- a/src/modules/common/moduleTable.f90 +++ b/src/modules/common/moduleTable.f90 @@ -30,6 +30,8 @@ MODULE moduleTable READ(id, '(A)', iostat = stat) dummy !If EOF or error, exit file IF (stat /= 0) EXIT + !If empty line, exit file + IF (dummy == "") EXIT !Skip comment IF (INDEX(dummy,'#') /= 0) CYCLE !Add row @@ -55,6 +57,7 @@ MODULE moduleTable !TODO: Make this a function IF (INDEX(dummy,'#') /= 0) CYCLE IF (stat /= 0) EXIT + IF (dummy == "") EXIT !Add data !TODO: substitute with extracting information from dummy BACKSPACE(id) diff --git a/src/modules/mesh/moduleMeshBoundary.f90 b/src/modules/mesh/moduleMeshBoundary.f90 index 90420c6..d130a8a 100644 --- a/src/modules/mesh/moduleMeshBoundary.f90 +++ b/src/modules/mesh/moduleMeshBoundary.f90 @@ -212,19 +212,24 @@ MODULE moduleMeshBoundary END SUBROUTINE symmetryAxis + !Secondary emission of electrons by impacting particle. SUBROUTINE secondaryEmission(edge, part) USE moduleSpecies USE moduleRandom + USE moduleConstParam IMPLICIT NONE CLASS(meshEdge), INTENT(inout):: edge CLASS(particle), INTENT(inout):: part REAL(8):: vRel, eRel + REAL(8), DIMENSION(1:3):: rElectron, XiElectron!Position of new electrons (impacting particle) + REAL(8), DIMENSION(1:3):: rIncident !Vector from imapcting particle position to particle position + REAL(8):: theta !incident angle REAL(8):: yield + REAL(8):: energy !incident energy corrected by threshold and INTEGER:: nNewElectrons REAL(8), ALLOCATABLE:: weight(:) INTEGER:: p - REAL(8), DIMENSION(1:3):: rElectron, XiElectron!Position of new electrons INTEGER:: cell TYPE(particle), POINTER:: secondaryElectron @@ -235,64 +240,81 @@ MODULE moduleMeshBoundary !Convert to relative energy eRel = part%species%m*vRel**2*5.D-1 - !Get number of secondary electrons macro-particles - yield = part%weight*bound%yield%get(eRel) - nNewElectrons = FLOOR(yield / bound%electron%weight) - IF (REAL(nNewElectrons) * bound%electron%weight < yield) THEN - nNewElectrons = nNewElectrons + 1 - ALLOCATE(weight(1:nNewElectrons)) - weight(1:nNewElectrons-1) = bound%electron%weight - weight(nNewElectrons) = yield - SUM(weight(1:nNewElectrons-1)) + !If energy is abound threshold calculate secondary electrons + IF (eRel >= bound%energyThreshold) THEN - ELSE - ALLOCATE(weight(1:nNewElectrons)) - weight(1:nNewElectrons) = bound%electron%weight - - END IF + !position of impacting ion + rElectron = edge%intersection(part%r) + XiElectron = mesh%cells(part%cell)%obj%phy2log(rElectron) - !position of impacting ion - rElectron = edge%intersection(part%r) - XiElectron = mesh%cells(part%cell)%obj%phy2log(rElectron) + !Calculate incident angle + rIncident = part%r - rElectron + theta = ACOS(DOT_PRODUCT(rIncident, edge%normal) / (NORM2(rIncident) * NORM2(edge%normal))) - !New cell of origin - IF (ASSOCIATED(edge%e1)) THEN - cell = edge%e1%n + !Calculate incident energy + energy = (eRel - bound%energyThreshold) * PI2 / (PI2 + theta**2) + bound%energyThreshold - ELSEIF (ASSOCIATED(edge%e2)) THEN - cell = edge%e2%n + !Get number of secondary electrons particles + yield = part%weight*bound%yield%get(eRel) * (1.D0 * theta**2 / PI2) !Check equation! + + !Convert the number to macro-particles + nNewElectrons = FLOOR(yield / bound%electron%weight) + + !If the weight of new macro-particles is below the yield, correct adding a particle + IF (REAL(nNewElectrons) * bound%electron%weight < yield) THEN + nNewElectrons = nNewElectrons + 1 + ALLOCATE(weight(1:nNewElectrons)) + weight(1:nNewElectrons-1) = bound%electron%weight + weight(nNewElectrons) = yield - SUM(weight(1:nNewElectrons-1)) + + ELSE + ALLOCATE(weight(1:nNewElectrons)) + weight(1:nNewElectrons) = bound%electron%weight + + END IF + + !New cell of origin + IF (ASSOCIATED(edge%e1)) THEN + cell = edge%e1%n + + ELSEIF (ASSOCIATED(edge%e2)) THEN + cell = edge%e2%n + + END IF + + !Create the new electron macro-particles + DO p = 1, nNewElectrons + !Create new macro-particle + ALLOCATE(secondaryElectron) + + !Assign species to electron + secondaryElectron%species => bound%electron + + !Assign position to particle + secondaryElectron%r = rElectron + secondaryElectron%Xi = XiElectron + + !Assign cell to electron + secondaryElectron%cell = cell + + !Assign weight + secondaryElectron%weight = weight(p) + + !Assume particle is inside the numerical domain + secondaryElectron%n_in = .TRUE. + + !Assign velocity + secondaryElectron%v = 2.D0 * edge%normal + 1.D-1 * (/ randomMaxwellian(), randomMaxwellian(), randomMaxwellian() /) + + !Add particle to list + CALL partSurfaces%setLock() + CALL partSurfaces%add(secondaryElectron) + CALL partSurfaces%unsetLock() + + END DO END IF - DO p = 1, nNewElectrons - !Create new macro-particle - ALLOCATE(secondaryElectron) - - !Assign species to electron - secondaryElectron%species => bound%electron - - !Assign position to particle - secondaryElectron%r = rElectron - secondaryElectron%Xi = XiElectron - - !Assign cell to electron - secondaryElectron%cell = cell - - !Assign weight - secondaryElectron%weight = weight(p) - - !Assume particle is inside the numerical domain - secondaryElectron%n_in = .TRUE. - - !Assign velocity - secondaryElectron%v = 2.D0 * edge%normal + 1.D-1 * (/ randomMaxwellian(), randomMaxwellian(), randomMaxwellian() /) - - !Add particle to list - CALL partSurfaces%setLock() - CALL partSurfaces%add(secondaryElectron) - CALL partSurfaces%unsetLock() - - END DO - !Absorb impacting particle CALL absorption(edge, part) diff --git a/src/modules/moduleBoundary.f90 b/src/modules/moduleBoundary.f90 index bfdbd70..b3c1d12 100644 --- a/src/modules/moduleBoundary.f90 +++ b/src/modules/moduleBoundary.f90 @@ -51,6 +51,7 @@ MODULE moduleBoundary !Yield as a function of ion energy TYPE(table1D):: yield CLASS(speciesGeneric), POINTER:: electron !Electron species for secondary emission + REAL(8):: energyThreshold CONTAINS END TYPE boundarySEE @@ -155,6 +156,7 @@ MODULE moduleBoundary CLASS(boundaryGeneric), ALLOCATABLE, INTENT(out):: boundary CHARACTER(:), ALLOCATABLE, INTENT(in):: tableFile INTEGER:: speciesID + INTEGER:: i ALLOCATE(boundarySEE:: boundary) SELECT TYPE(boundary) @@ -162,6 +164,16 @@ MODULE moduleBoundary CALL boundary%yield%init(tableFile) CALL boundary%yield%convert(eV2J/(m_ref*v_ref**2), 1.D0) boundary%electron => species(speciesID)%obj + !Search for the threshold energy in the table + DO i = 1, SIZE(boundary%yield%f) + IF (boundary%yield%f(i) > 0.D0) THEN + boundary%energyThreshold = boundary%yield%x(i) + + EXIT + + END IF + + END DO END SELECT