TL: script update

Author: tlunet

pySDC/playgrounds/dedalus/scripts/stab_rbc.py

@@ -7,14 +7,22 @@
 """
 import os
 import matplotlib.pyplot as plt
+import numpy as np
 
 from pySDC.playgrounds.dedalus.timestepper import SDCIMEX
 from pySDC.playgrounds.dedalus.problems.rbc import RBCProblem3D, OutputFiles
 
+tEnd = 10
 
-dtBase = 2.50e-03
-nSteps = 10
-nCycles = 2
+nStepsMax = 50
+nStepsMin = 8
+nVals = 19
+
+dtVals = 1/np.arange(nStepsMax, nStepsMin-1, -1)
+intervals = np.linspace(dtVals.min(), dtVals.max(), num=nVals)
+dtVals = np.unique([dtVals[max(np.argwhere(dtVals <= dt))] for dt in intervals])
+
+nStepsVals = [int(n) for n in 1/dtVals]
 
 Rayleigh = 1.5e5
 timeScheme = "SDC"
@@ -23,22 +31,21 @@
 
 SDCIMEX.setParameters(
     nNodes=4, nodeType="LEGENDRE", quadType="RADAU-RIGHT",
-    nSweeps=4, initSweep="COPY", explSweep="PIC", implSweep="MIN-SR-FLEX",
+    nSweeps=4, initSweep="COPY", explSweep="PIC", implSweep="MIN-SR-S",
     )
 
 os.makedirs(stabDir, exist_ok=True)
 
-dtFactors = [f*10**e for e in range(nCycles) for f in [1, 2, 5]]
-fmtSuffix = f":0{nCycles}d"
+fmtSuffix = f":0{len(str(nStepsMax))}d"
 
 plt.figure("spectrum")
-for fac in dtFactors:
-    dtRun = dtBase*fac
-    runDir = f"{stabDir}/dt_f" + ("{"+fmtSuffix+"}").format(fac)
+for nSteps in nStepsVals:
+    dtRun = 1/nSteps
+    runDir = f"{stabDir}/dt_N" + ("{"+fmtSuffix+"}").format(nSteps)
 
     prob = RBCProblem3D.runSimulation(
-        runDir, dtRun*nSteps, dtRun, timeScheme=timeScheme,
-        dtWrite=dtRun*nSteps, initField=initSol,
+        runDir, tEnd, dtRun, timeScheme=timeScheme,
+        dtWrite=tEnd, initField=initSol,
         aspectRatio=4, meshRatio=1, resFactor=1,
         Rayleigh=Rayleigh, Prandtl=0.7,
         )
@@ -49,6 +56,8 @@
     output = OutputFiles(runDir)
     spectrum = output.getSpectrum(which="all", start=1)
 
-#     plt.loglog(spectrum["kappa"], spectrum["p"], label="f"+("{"+fmtSuffix+"}").format(fac))
+    if np.any(np.isnan(spectrum["u"])):
+        break
+    plt.loglog(spectrum["kappa"], spectrum["u"], label="N"+("{"+fmtSuffix+"}").format(nSteps))
 
-# plt.legend()
+plt.legend()