Working on raiju config and output adjustments

This commit is contained in:
Anthony
2025-04-02 09:46:31 -06:00
parent 35b66c5012
commit 9b3655f1f2
5 changed files with 61 additions and 80 deletions

View File

@@ -69,7 +69,7 @@ module raijudefs
real(rp), parameter :: def_pdmb = 0.75
real(rp), parameter :: def_cfl = 0.3
real(rp), parameter :: cflMax = 0.3
logical, parameter :: def_doUseVelLRs = .false.
logical, parameter :: def_doUseVelLRs = .true.
! Domain limits
real(rp), parameter :: def_maxTail_buffer = 15.0 ! [Rp]

View File

@@ -219,12 +219,12 @@ module raijutypes
!! For debug
logical :: doClockConsoleOut
!! If we are driving, output clock info
logical :: doFatOutput
logical :: doOutput_potGrads
!! Output extra 3D arrays
logical :: doDebugOutput
logical :: doOutput_debug
!! Dump lots of otherwise unnecessary stuff
logical :: doLossOutput
!! Dump 3D loss variables like Tau
logical :: doOutput_3DLoss
!! Dump extra 3D loss variables
logical :: doOwnCorot
!! If true, we calculate corortation potential ourselves and include it in our velocity calculation
!! If false, we assume its already included in whatever ionospheric potential we're using

View File

@@ -284,26 +284,7 @@ module raijuIO
call AddOutVar(IOVars,"activeShells",outTmp2D,uStr="[Ni, Nk]")
deallocate(outTmp2D)
! Moments
!allocate(outTmp3D(is:ie,js:je,0:Grid%nSpc))
!outTmp3D = 0.0
!do s=0,Grid%nSpc
! outTmp3D(:,:,s) = State%Press(s)%data(is:ie,js:je)
!enddo
!call AddOutVar(IOVars,"Pressure",outTmp3D(:,:,:),uStr="nPa")
! Add density moment as #/cc
!outTmp3D = 0.0
!do s=1, Grid%nSpc
! ! Convert amu/cc to #/cc
! outTmp3D(:,:,s) = State%Den(s)%data(is:ie,js:je)/Grid%spc(s)%amu
! ! Don't include electrons to total number density
! if(Grid%spc(s)%spcType .ne. RAIJUELE) then
! outTmp3D(:,:,0) = outTmp3D(:,:,0) + outTmp3D(:,:,s)
! endif
!enddo
!call AddOutVar(IOVars,"Density",outTmp3D(:,:,:),uStr="#/cc")
!deallocate(outTmp3D)
! Moments
call AddOutSGV(IOVars, "Pressure", State%Press, outBndsO=outBnds2D, uStr="nPa" , dStr="Pressure from RAIJU flavors", doWriteMaskO=.false.)
call AddOutSGV(IOVars, "Density" , State%Den , outBndsO=outBnds2D, uStr="#/cc", dStr="Density from RAIJU flavors", doWriteMaskO=.false.)
@@ -325,44 +306,48 @@ module raijuIO
call AddOutVar(IOVars,"FTEntropy",outTmp2D,uStr="nPa*(Rp/nT)^(5/3)")
deallocate(outTmp2D)
if (Model%doLossOutput) then
!--- Losses ---!
! Calculate accumulated precipitation for each species
allocate(outPrecipN (is:ie,js:je,Grid%nSpc))
allocate(outPrecipE (is:ie,js:je,Grid%nSpc))
allocate(outPrecipAvgE(is:ie,js:je,Grid%nSpc))
allocate(outCCHeatFlux(is:ie,js:je,Grid%nSpc))
do s=1,Grid%nSpc
ks = Grid%spc(s)%kStart
ke = Grid%spc(s)%kEnd
outPrecipN(:,:,s) = sum(State%precipNFlux(is:ie,js:je,kS:kE), dim=3)
outPrecipE(:,:,s) = sum(State%precipEFlux(is:ie,js:je,kS:kE), dim=3)
outCCHeatFlux(:,:,s) = sum(State%CCHeatFlux (is:ie,js:je,kS:kE), dim=3)
where (outPrecipN(:,:,s) > TINY)
outPrecipAvgE(:,:,s) = outPrecipE(:,:,s)/outPrecipN(:,:,s) * erg2kev ! Avg E [keV]
elsewhere
outPrecipAvgE(:,:,s) = 0.0
end where
enddo
call AddOutVar(IOVars, "precipNFlux", outPrecipN , uStr="#/cm^2/s")
call AddOutVar(IOVars, "precipEFlux", outPrecipE , uStr="erg/cm^2/s")
call AddOutVar(IOVars, "precipAvgE" , outPrecipAvgE, uStr="keV")
call AddOutVar(IOVars, "CCHeatFlux" , outCCHeatFlux, uStr="eV/cm^2/s")
deallocate(outPrecipN)
deallocate(outPrecipE)
deallocate(outPrecipAvgE)
deallocate(outCCHeatFlux)
! (Ni, Nj, Nk)
call AddOutVar(IOVars, "lossRate", State%lossRates(is:ie,js:je,:), uStr="1/s")
call AddOutVar(IOVars, "lossRatePrecip", State%lossRatesPrecip(is:ie,js:je,:), uStr="1/s")
if (Model%doOutput_3DLoss) then
call AddOutVar(IOVars, "dEta_dt" , State%dEta_dt(is:ie,js:je,:), uStr="eta_units/s")
call AddOutVar(IOVars, "lossRate", State%lossRates(is:ie,js:je,:), uStr="1/s")
call AddOutVar(IOVars, "lossRatePrecip", State%lossRatesPrecip(is:ie,js:je,:), uStr="1/s")
! Calculate accumulated precipitation for each species
allocate(outPrecipN (is:ie,js:je,Grid%nSpc))
allocate(outPrecipE (is:ie,js:je,Grid%nSpc))
allocate(outPrecipAvgE(is:ie,js:je,Grid%nSpc))
allocate(outCCHeatFlux(is:ie,js:je,Grid%nSpc))
do s=1,Grid%nSpc
ks = Grid%spc(s)%kStart
ke = Grid%spc(s)%kEnd
outPrecipN(:,:,s) = sum(State%precipNFlux(is:ie,js:je,kS:kE), dim=3)
outPrecipE(:,:,s) = sum(State%precipEFlux(is:ie,js:je,kS:kE), dim=3)
outCCHeatFlux(:,:,s) = sum(State%CCHeatFlux (is:ie,js:je,kS:kE), dim=3)
where (outPrecipN(:,:,s) > TINY)
outPrecipAvgE(:,:,s) = outPrecipE(:,:,s)/outPrecipN(:,:,s) * erg2kev ! Avg E [keV]
elsewhere
outPrecipAvgE(:,:,s) = 0.0
end where
enddo
call AddOutVar(IOVars, "precipNFlux", outPrecipN , uStr="#/cm^2/s")
call AddOutVar(IOVars, "precipEFlux", outPrecipE , uStr="erg/cm^2/s")
call AddOutVar(IOVars, "precipAvgE" , outPrecipAvgE, uStr="keV")
call AddOutVar(IOVars, "CCHeatFlux" , outCCHeatFlux, uStr="eV/cm^2/s")
deallocate(outPrecipN)
deallocate(outPrecipE)
deallocate(outPrecipAvgE)
deallocate(outCCHeatFlux)
call AddOutVar(IOVars, "precipNFlux_Nk" , State%precipNFlux(is:ie,js:je,:), uStr="#/cm^2/s")
call AddOutVar(IOVars, "precipEFlux_Nk" , State%precipEFlux(is:ie,js:je,:), uStr="erg/cm^2/s")
call AddOutVar(IOVars, "CCHeatFlux_Nk" , State%CCHeatFlux (is:ie,js:je,:), uStr="eV/cm^2/s")
endif
if (Model%doFatOutput) then
if (Model%doOutput_potGrads) then
! (Ni, Nj, 2)
call AddOutVar(IOVars, "gradPotE" , State%gradPotE (is:ie+1,js:je+1,:), uStr="V/m")
call AddOutVar(IOVars, "gradPotCorot" , State%gradPotCorot(is:ie+1,js:je+1,:), uStr="V/m")
@@ -370,10 +355,6 @@ module raijuIO
call AddOutVar(IOVars, "gradPotE_cc" , State%gradPotE_cc (is:ie,js:je,:) , uStr="V/m")
call AddOutVar(IOVars, "gradPotCorot_cc", State%gradPotCorot_cc(is:ie,js:je,:) , uStr="V/m")
call AddOutVar(IOVars, "gradVM_cc" , State%gradVM_cc (is:ie,js:je,:) , uStr="V/m/lambda")
! (Ni, Nj, Nk)
call AddOutVar(IOVars, "precipNFlux_Nk" , State%precipNFlux(is:ie,js:je,:), uStr="#/cm^2/s")
call AddOutVar(IOVars, "precipEFlux_Nk" , State%precipEFlux(is:ie,js:je,:), uStr="erg/cm^2/s")
call AddOutVar(IOVars, "CCHeatFlux_Nk" , State%CCHeatFlux (is:ie,js:je,:), uStr="eV/cm^2/s")
! Calc pEffective based on current state
! Make full ghost size since that's what the subroutine expects
@@ -383,7 +364,7 @@ module raijuIO
deallocate(outTmp3D)
endif
if (Model%doDebugOutput) then
if (Model%doOutput_debug) then
call AddOutVar(IOVars, "eta_half" , State%eta_half (is:ie ,js:je ,:) , uStr="#/cm^3 * Rx/T")
call AddOutVar(IOVars, "iVel" , State%iVel (is:ie+1,js:je+1,:,:), uStr="m/s")
call AddOutVar(IOVars, "iVelL" , State%iVelL (is:ie+1,js:je+1,:,:), uStr="m/s")

View File

@@ -457,7 +457,7 @@ module raijuRecon
! ReconFaces(Model, Grid, isG, Qcc, QfaceL, QfaceR, QreconLO, QreconRO)
if (Model%doDebugOutput) then
if (Model%doOutput_debug) then
call ReconFaces(Model, Grid, isG, Qcc, QfaceL, QfaceR, State%etaFaceReconL(:,:,k,:), State%etaFaceReconR(:,:,k,:))
!call ReconFaces(Model, Grid, isG, Qcc, QfaceL, QfaceR, State%etaFaceReconL, State%etaFaceReconR)
else
@@ -504,7 +504,7 @@ module raijuRecon
! Now that valid faces are decided, we can apply our invalid/buffer boundary conditions
!call calcBoundaryFluxes(Grid%shGrid, State%active, Qflux)
if (Model%doDebugOutput) then
if (Model%doOutput_debug) then
State%etaFacePDML(:,:,k,:) = QfaceL
State%etaFacePDMR(:,:,k,:) = QfaceR
State%etaFlux (:,:,k,:) = Qflux

View File

@@ -114,7 +114,7 @@ module raijustarter
stop
endif
! Restart time
!--- Restarts ---!
if (Model%isSA) then
tmpStr = "restart/doRes"
else
@@ -138,12 +138,8 @@ module raijustarter
call iXML%Set_Val(tmpResId, trim(tmpStr), Model%RunID)
call genResInFname(Model, Model%ResF, runIdO=tmpResId) ! Determine filename to read from
endif
call iXML%Set_Val(Model%isLoud, "debug/isLoud",.false.)
call iXML%Set_Val(Model%writeGhosts, "debug/writeGhosts",.false.)
call iXML%Set_Val(Model%doDebugOutput, "debug/debugOutput",.false.)
! Plasmasphere settings
!--- Plasmasphere ---!
call iXML%Set_Val(Model%doPlasmasphere, "plasmasphere/doPsphere",.false.)
call iXML%Set_Val(Model%doPsphEvol, 'plasmasphere/doEvol',.true.)
! Determine number of species. First set default, then read from xml to overwrite if present
@@ -161,7 +157,7 @@ module raijustarter
call iXML%Set_Val(Model%nSpc, "prob/nSpc",Model%nSpc)
! Domain constraints
!--- Domain ---!
call iXML%Set_Val(Model%n_bndLim , "domain/bndLim" , 3)
call iXML%Set_Val(Model%maxTail_buffer, "domain/tail_buffer", def_maxTail_buffer)
call iXML%Set_Val(Model%maxSun_buffer , "domain/sun_buffer" , def_maxSun_buffer)
@@ -173,7 +169,7 @@ module raijustarter
Model%maxTail_active = abs(Model%maxTail_active)
Model%maxSun_active = abs(Model%maxSun_active)
! Solver params
!---Solver ---!
call iXML%Set_Val(Model%doUseVelLRs,'sim/useVelLRs',def_doUseVelLRs)
call iXML%Set_Val(Model%maxItersPerSec,'sim/maxIter',def_maxItersPerSec)
call iXML%Set_Val(Model%maxOrder,'sim/maxOrder',7)
@@ -199,9 +195,10 @@ module raijustarter
call iXML%Set_Val(Model%tiote, "prob/tiote",4.0)
! Active shell settings
call iXML%Set_Val(Model%doActiveShell, "activeShell/doAS",.true.)
call iXML%Set_Val(Model%worthyFrac, "prob/worthyFrac",fracWorthyDef)
call iXML%Set_Val(Model%doActiveShell, "activeShell/doAS",.false.)
call iXML%Set_Val(Model%worthyFrac, "activeShell/worthyFrac",fracWorthyDef)
! Corotation
if (Model%isSA) then
call iXML%Set_Val(Model%doOwnCorot, "prob/doCorot",.true.)
else
@@ -230,16 +227,19 @@ module raijustarter
end select
! Losses
!--- Losses ---!
call iXML%Set_Val(Model%doLosses, "losses/doLosses",.true.)
call iXML%Set_Val(Model%doSS , "losses/doSS" ,.true. )
call iXML%Set_Val(Model%doCC , "losses/doCC" ,.true. )
call iXML%Set_Val(Model%doCX , "losses/doCX" ,.true.)
call iXML%Set_Val(Model%doFLC , "losses/doFLC",.false.)
call iXML%Set_Val(Model%doLossOutput, "losses/doOutput",.false.) ! Several (Ni,Nj,Nk) arrays
call iXML%Set_Val(Model%doFatOutput, "output/doFat",.false.)
! TODO: Add flags to output certain data, like coupling information
!--- Output ---!
call iXML%Set_Val(Model%isLoud , "output/loudConsole",.false.)
call iXML%Set_Val(Model%writeGhosts , "output/writeGhosts",.false.)
call iXML%Set_Val(Model%doOutput_potGrads, "output/doFat" ,.false.)
call iXML%Set_Val(Model%doOutput_3DLoss , "output/lossExtras" ,.false.) ! Several (Ni,Nj,Nk) arrays
call iXML%Set_Val(Model%doOutput_debug , "output/doDebug" ,.false.)
! Model Hax?
call iXML%Set_Val(Model%doHack_rcmEtaPush, "hax/rcmEtaPush",.false.)
@@ -452,7 +452,7 @@ module raijustarter
enddo
! Only bother allocating persistent versions of debug stuff if we ned them
if (Model%doDebugOutput) then
if (Model%doOutput_debug) then
allocate( State%etaFaceReconL(sh%isg:sh%ieg+1, sh%jsg:sh%jeg+1, Grid%Nk, 2) )
allocate( State%etaFaceReconR(sh%isg:sh%ieg+1, sh%jsg:sh%jeg+1, Grid%Nk, 2) )
allocate( State%etaFacePDML (sh%isg:sh%ieg+1, sh%jsg:sh%jeg+1, Grid%Nk, 2) )