-- -- Patricia Fasel -- Los Alamos National Laboratory -- 1990 August -- module TransPort (transPort) where import GamtebType import Consts import Utils import Distance import RoulSplit import PhotoElec import Compton import Pair (=:) a b = (a,b) -- transport a particle transPort :: Particle -> Probability -> ([Result], [Stat]) transPort p prob = if (dColl < dSurf) then (let pos' = transPos pos dir dColl p' = Part pos' dir w e eIndx cell seed doCompton = (r1 < (pComp / (pTot-pPhot))) (res, stat) = collision p' prob doCompton in (res, [nc=:1]++stat)) else (let pos' = transPos pos dir dSurf p' = Part pos' dir w e eIndx cell seed (res, stat) = noCollision p' prob surf in (res, [nnc=:1]++stat)) where (Part pos dir w e eIndx cell seed) = p (pComp,pPair,pPhot,pTot) = prob (r, r1) = genRand seed (dSurf, surf) = distSurf pos dir dColl = -((log r)/ pTot) -- no collision inside cylinder noCollision :: Particle -> Probability -> Int -> ([Result], [Stat]) noCollision p@(Part pos dir w e eIndx cell seed) prob surf = case surf of 1 -> ([(scatter, eIndx)=:w], [ns=:1]) 2 -> ([(escape, eIndx)=:w], [ne=:1]) 4 -> ([(transit, eIndx)=:w], [nt=:1]) 3 -> -- cross internal surface -- particle will split, die in russian roulette, or continue -- cell = [1..] causes roulet or split to alternate if (cell == 1) then (let (p1, p2) = split p (r1, s1) = transPort p1 prob (r2, s2) = transPort p2 prob in (r1++r2, [nsp=:1]++s1++s2)) else (let (p', stat', roulKill) = roulet p in if (roulKill) then ([], stat') else (let (res, stat) = transPort p' prob in (res, stat++stat'))) -- collision is in cylinder, do collision physics collision :: Particle -> Probability -> Bool -> ([Result], [Stat]) collision p prob doCompton = if (wgtKill) then ([], [nwk=:1]) else if (doCompton) then -- compton scattering (let (p'', prob', comptonCut) = compton p' in if (comptonCut) then ([], [nek=:1]) else transPort p'' prob') else -- pair production (let (p'', prob', pairCut) = pair p' in if (pairCut) then ([], [nek=:1]) else transPort p'' prob') where (Part pos dir w e eIndx cell seed) = p (pComp,pPair,pPhot,pTot) = prob (p', absorb, wgtKill) = photoElec p prob -- translate a particle position transPos :: Point -> Point -> Value -> Point transPos (x,y,z) (u,v,w) dist = (x',y',z') where x' = x + u*dist y' = y + v*dist z' = z + w*dist