Skip to content

Commit b934988

Browse files
pancettapancetta
committed
starting to merge boussinesq into FT framework
1 parent f817ac7 commit b934988

File tree

5 files changed

+149
-13
lines changed

5 files changed

+149
-13
lines changed

‎examples/boussinesq_2d_imex/ProblemClass.py‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -5,9 +5,9 @@
55

66
from pySDC.Problem import ptype
77
from pySDC.datatype_classes.mesh import mesh, rhs_imex_mesh
8-
from build2DFDMatrix import get2DMesh
9-
from buildBoussinesq2DMatrix import getBoussinesq2DMatrix, getBoussinesqBCHorizontal, getBoussinesqBCVertical, getBoussinesq2DUpwindMatrix
10-
from unflatten import unflatten
8+
from examples.boussinesq_2d_imex.build2DFDMatrix import get2DMesh
9+
from examples.boussinesq_2d_imex.buildBoussinesq2DMatrix import getBoussinesq2DMatrix, getBoussinesqBCHorizontal, getBoussinesqBCVertical, getBoussinesq2DUpwindMatrix
10+
from examples.boussinesq_2d_imex.unflatten import unflatten
1111

1212
class logging(object):
1313

‎examples/boussinesq_2d_imex/build2DFDMatrix.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@
22
sys.path.append('../')
33
import numpy as np
44
import scipy.sparse as sp
5-
from buildFDMatrix import getMatrix, getUpwindMatrix, getBCLeft, getBCRight
5+
from examples.boussinesq_2d_imex.buildFDMatrix import getMatrix, getUpwindMatrix, getBCLeft, getBCRight
66

77
#
88
#

‎examples/boussinesq_2d_imex/buildBoussinesq2DMatrix.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
import numpy as np
22
import scipy.sparse as sp
3-
from build2DFDMatrix import get2DMatrix, getBCHorizontal, getBCVertical, get2DUpwindMatrix
3+
from examples.boussinesq_2d_imex.build2DFDMatrix import get2DMatrix, getBCHorizontal, getBCVertical, get2DUpwindMatrix
44

55
def getBoussinesq2DUpwindMatrix(N, dx, u_adv, order):
66

‎examples/boussinesq_2d_imex/playground.py‎

Lines changed: 6 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -3,7 +3,7 @@
33

44
import numpy as np
55

6-
from ProblemClass import boussinesq_2d_imex
6+
from examples.boussinesq_2d_imex.ProblemClass import boussinesq_2d_imex
77
from examples.boussinesq_2d_imex.TransferClass import mesh_to_mesh_2d
88
from examples.boussinesq_2d_imex.HookClass import plot_solution
99

@@ -28,22 +28,20 @@
2828

2929
# This comes as read-in for the level class
3030
lparams = {}
31-
lparams['restol'] = 5E-6
31+
lparams['restol'] = 1E-06
3232

3333
swparams = {}
3434
swparams['collocation_class'] = collclass.CollGaussLobatto
3535
swparams['num_nodes'] = 3
3636
swparams['do_LU'] = True
3737

3838
sparams = {}
39-
sparams['maxiter'] = 16
39+
sparams['maxiter'] = 50
4040

4141
# setup parameters "in time"
4242
t0 = 0
43-
Tend = 384
44-
Nsteps = 128
45-
#Tend = 30
46-
#Nsteps = 5
43+
Tend = 960
44+
Nsteps = 320
4745
dt = Tend/float(Nsteps)
4846

4947
# This comes as read-in for the problem class
@@ -62,7 +60,7 @@
6260

6361
# This comes as read-in for the transfer operations
6462
tparams = {}
65-
tparams['finter'] = False
63+
tparams['finter'] = True
6664

6765
# Fill description dictionary for easy hierarchy creation
6866
description = {}

‎examples/fault_tolerance/boussinesq_example.py‎

Lines changed: 138 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,138 @@
1+
2+
from pySDC import CollocationClasses as collclass
3+
4+
from examples.boussinesq_2d_imex.ProblemClass import boussinesq_2d_imex
5+
from examples.boussinesq_2d_imex.TransferClass import mesh_to_mesh_2d
6+
from examples.boussinesq_2d_imex.HookClass import plot_solution
7+
8+
from pySDC.datatype_classes.mesh import mesh, rhs_imex_mesh
9+
from pySDC.sweeper_classes.imex_1st_order import imex_1st_order
10+
import pySDC.PFASST_stepwise as mp
11+
12+
from pySDC import Log
13+
from pySDC.Stats import grep_stats, sort_stats
14+
15+
import numpy as np
16+
17+
import pySDC.Plugins.fault_tolerance as ft
18+
19+
20+
if __name__ == "__main__":
21+
22+
# set global logger (remove this if you do not want the output at all)
23+
logger = Log.setup_custom_logger('root')
24+
25+
num_procs = 16
26+
27+
# This comes as read-in for the level class
28+
lparams = {}
29+
lparams['restol'] = 1E-06
30+
31+
swparams = {}
32+
swparams['collocation_class'] = collclass.CollGaussLobatto
33+
swparams['num_nodes'] = 3
34+
swparams['do_LU'] = True
35+
36+
sparams = {}
37+
sparams['maxiter'] = 50
38+
39+
# setup parameters "in time"
40+
t0 = 0
41+
Tend = 960
42+
Nsteps = 320
43+
dt = Tend/float(Nsteps)
44+
45+
# This comes as read-in for the problem class
46+
pparams = {}
47+
pparams['nvars'] = [(4,450,30), (4,450,30)]
48+
pparams['u_adv'] = 0.02
49+
pparams['c_s'] = 0.3
50+
pparams['Nfreq'] = 0.01
51+
pparams['x_bounds'] = [(-150.0, 150.0)]
52+
pparams['z_bounds'] = [( 0.0, 10.0)]
53+
pparams['order'] = [4, 2] # [fine_level, coarse_level]
54+
pparams['order_upw'] = [5, 1]
55+
pparams['gmres_maxiter'] = [50, 50]
56+
pparams['gmres_restart'] = [10, 10]
57+
pparams['gmres_tol'] = [1e-10, 1e-10]
58+
59+
# This comes as read-in for the transfer operations
60+
tparams = {}
61+
tparams['finter'] = True
62+
63+
# Fill description dictionary for easy hierarchy creation
64+
description = {}
65+
description['problem_class'] = boussinesq_2d_imex
66+
description['problem_params'] = pparams
67+
description['dtype_u'] = mesh
68+
description['dtype_f'] = rhs_imex_mesh
69+
description['sweeper_params'] = swparams
70+
description['sweeper_class'] = imex_1st_order
71+
description['level_params'] = lparams
72+
description['hook_class'] = plot_solution
73+
description['transfer_class'] = mesh_to_mesh_2d
74+
description['transfer_params'] = tparams
75+
76+
ft.hard_random = 0.03
77+
78+
strategies = ['NOFAULT']
79+
# strategies = ['NOFAULT','SPREAD','INTERP','INTERP_PREDICT','SPREAD_PREDICT']
80+
81+
for strategy in strategies:
82+
83+
print('------------------------------------------ working on strategy ',strategy)
84+
ft.strategy = strategy
85+
86+
# quickly generate block of steps
87+
MS = mp.generate_steps(num_procs,sparams,description)
88+
89+
# get initial values on finest level
90+
P = MS[0].levels[0].prob
91+
uinit = P.u_exact(t0)
92+
93+
cfl_advection = pparams['u_adv']*dt/P.h[0]
94+
cfl_acoustic_hor = pparams['c_s']*dt/P.h[0]
95+
cfl_acoustic_ver = pparams['c_s']*dt/P.h[1]
96+
print ("CFL number of advection: %4.2f" % cfl_advection)
97+
print ("CFL number of acoustics (horizontal): %4.2f" % cfl_acoustic_hor)
98+
print ("CFL number of acoustics (vertical): %4.2f" % cfl_acoustic_ver)
99+
100+
# call main function to get things done...
101+
uend,stats = mp.run_pfasst(MS,u0=uinit,t0=t0,dt=dt,Tend=Tend)
102+
103+
P.report_log()
104+
105+
# stats magic: get all residuals
106+
extract_stats = grep_stats(stats,level=-1,type='residual')
107+
108+
# find boundaries
109+
maxsteps = 0
110+
maxiter = 0
111+
minres = 0
112+
maxres = -99
113+
for k,v in extract_stats.items():
114+
maxsteps = max(maxsteps,getattr(k,'step'))
115+
maxiter = max(maxiter,getattr(k,'iter'))
116+
minres = min(minres,np.log10(v))
117+
maxres = max(maxres,np.log10(v))
118+
119+
# fill residual array
120+
residual = np.zeros((maxiter,maxsteps+1))
121+
residual[:] = 0
122+
for k,v in extract_stats.items():
123+
step = getattr(k,'step')
124+
iter = getattr(k,'iter')
125+
if iter is not -1:
126+
residual[iter-1,step] = v
127+
128+
# stats magic: get niter (probably redundant with maxiter)
129+
extract_stats = grep_stats(stats,iter=-1,type='niter')
130+
iter_count = np.zeros(maxsteps+1)
131+
for k,v in extract_stats.items():
132+
step = getattr(k,'step')
133+
iter_count[step] = v
134+
print(iter_count)
135+
136+
# np.savez('SDC_GRAYSCOTT_stats_hf_'+ft.strategy+'_new',residual=residual,iter_count=iter_count,hard_stats=ft.hard_stats)
137+
np.savez('PFASST_BOUSSINESQ_stats_hf_'+ft.strategy+'_P'+str(num_procs),residual=residual,iter_count=iter_count,hard_stats=ft.hard_stats)
138+

0 commit comments

Comments
 (0)