Coverage for pySDC/projects/parallelSDC/preconditioner_playground_MPI.py: 97%

1import os 

2import pickle 

3from collections import namedtuple 

4 

5import numpy as np 

6from mpi4py import MPI 

7 

8import pySDC.helpers.plot_helper as plt_helper 

9from pySDC.helpers.stats_helper import get_sorted 

10 

11from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI 

12from pySDC.implementations.problem_classes.AdvectionEquation_ND_FD import advectionNd 

13from pySDC.implementations.problem_classes.GeneralizedFisher_1D_FD_implicit import generalized_fisher 

14from pySDC.implementations.problem_classes.HeatEquation_ND_FD import heatNd_unforced 

15from pySDC.implementations.problem_classes.Van_der_Pol_implicit import vanderpol 

16from pySDC.implementations.sweeper_classes.generic_implicit_MPI import generic_implicit_MPI 

17 

18# from pySDC.projects.parallelSDC.BaseTransfer_MPI import base_transfer_mpi 

19 

20ID = namedtuple('ID', ['setup', 'qd_type', 'param']) 

21 

22 

23def main(comm=None): 

24 # initialize level parameters (part I) 

25 level_params = dict() 

26 level_params['restol'] = 1e-08 

27 

28 # initialize sweeper parameters (part I) 

29 sweeper_params = dict() 

30 sweeper_params['quad_type'] = 'RADAU-RIGHT' 

31 sweeper_params['num_nodes'] = comm.Get_size() 

32 sweeper_params['comm'] = comm 

33 

34 # initialize step parameters 

35 step_params = dict() 

36 step_params['maxiter'] = 100 

37 

38 # initialize controller parameters 

39 controller_params = dict() 

40 controller_params['logger_level'] = 30 

41 

42 # set up list of Q-delta types and setups 

43 qd_list = ['IEpar', 'Qpar', 'MIN', 'MIN3', 'MIN_GT'] 

44 setup_list = [ 

45 ('heat', 63, [10.0**i for i in range(-3, 3)]), 

46 ('advection', 64, [10.0**i for i in range(-3, 3)]), 

47 ('vanderpol', 2, [0.1 * 2**i for i in range(0, 10)]), 

48 ('fisher', 63, [2**i for i in range(-2, 3)]), 

49 ] 

50 # setup_list = [('fisher', 63, [2 * i for i in range(1, 6)])] 

51 

52 # pre-fill results with lists of setups 

53 results = dict() 

54 for setup, nvars, param_list in setup_list: 

55 results[setup] = (nvars, param_list) 

56 

57 # loop over all Q-delta matrix types 

58 for qd_type in qd_list: 

59 # assign implicit Q-delta matrix 

60 sweeper_params['QI'] = qd_type 

61 

62 # loop over all setups 

63 for setup, nvars, param_list in setup_list: 

64 # initialize problem parameters (part I) 

65 problem_params = dict() 

66 if setup != 'vanderpol': 

67 problem_params['nvars'] = nvars # number of degrees of freedom for each level 

68 

69 # loop over all parameters 

70 for param in param_list: 

71 # fill description for the controller 

72 description = dict() 

73 description['sweeper_class'] = generic_implicit_MPI # pass sweeper 

74 description['sweeper_params'] = sweeper_params # pass sweeper parameters 

75 description['step_params'] = step_params # pass step parameters 

76 # description['base_transfer_class'] = base_transfer_mpi 

77 

78 print('working on: %s - %s - %s' % (qd_type, setup, param)) 

79 

80 # decide which setup to take 

81 if setup == 'heat': 

82 problem_params['nu'] = param 

83 problem_params['freq'] = 2 

84 problem_params['bc'] = 'dirichlet-zero' # boundary conditions 

85 

86 level_params['dt'] = 0.1 

87 

88 description['problem_class'] = heatNd_unforced 

89 description['problem_params'] = problem_params 

90 description['level_params'] = level_params # pass level parameters 

91 

92 elif setup == 'advection': 

93 problem_params['c'] = param 

94 problem_params['order'] = 2 

95 problem_params['freq'] = 2 

96 problem_params['stencil_type'] = 'center' # boundary conditions 

97 problem_params['bc'] = 'periodic' # boundary conditions 

98 

99 level_params['dt'] = 0.1 

100 

101 description['problem_class'] = advectionNd 

102 description['problem_params'] = problem_params 

103 description['level_params'] = level_params # pass level parameters 

104 

105 elif setup == 'vanderpol': 

106 problem_params['newton_tol'] = 1e-09 

107 problem_params['newton_maxiter'] = 20 

108 problem_params['mu'] = param 

109 problem_params['u0'] = np.array([2.0, 0]) 

110 

111 level_params['dt'] = 0.1 

112 

113 description['problem_class'] = vanderpol 

114 description['problem_params'] = problem_params 

115 description['level_params'] = level_params 

116 

117 elif setup == 'fisher': 

118 problem_params['nu'] = 1 

119 problem_params['lambda0'] = param 

120 problem_params['newton_maxiter'] = 20 

121 problem_params['newton_tol'] = 1e-10 

122 problem_params['interval'] = (-5, 5) 

123 

124 level_params['dt'] = 0.01 

125 

126 description['problem_class'] = generalized_fisher 

127 description['problem_params'] = problem_params 

128 description['level_params'] = level_params 

129 

130 else: 

131 print('Setup not implemented..', setup) 

132 exit() 

133 

134 # instantiate controller 

135 controller = controller_nonMPI( 

136 num_procs=1, controller_params=controller_params, description=description 

137 ) 

138 

139 # get initial values on finest level 

140 P = controller.MS[0].levels[0].prob 

141 uinit = P.u_exact(0) 

142 

143 # call main function to get things done... 

144 uend, stats = controller.run(u0=uinit, t0=0, Tend=level_params['dt']) 

145 

146 # filter statistics by type (number of iterations) 

147 iter_counts = get_sorted(stats, type='niter', sortby='time') 

148 

149 # just one time-step, grep number of iteration and store 

150 niter = iter_counts[0][1] 

151 id = ID(setup=setup, qd_type=qd_type, param=param) 

152 results[id] = niter 

153 

154 assert len(results) == (6 + 6 + 10 + 5) * 5 + 4, 'ERROR: did not get all results, got %s' % len(results) 

155 

156 if comm.Get_rank() == 0: 

157 # write out for later visualization 

158 file = open('data/parallelSDC_iterations_precond_MPI.pkl', 'wb') 

159 pickle.dump(results, file) 

160 

161 assert os.path.isfile('data/parallelSDC_iterations_precond_MPI.pkl'), 'ERROR: pickle did not create file' 

162 

163 

164def plot_iterations(): 

165 """ 

166 Helper routine to plot iteration counts 

167 """ 

168 

169 file = open('data/parallelSDC_iterations_precond_MPI.pkl', 'rb') 

170 results = pickle.load(file) 

171 

172 # find the lists/header required for plotting 

173 qd_type_list = [] 

174 setup_list = [] 

175 for key in results.keys(): 

176 if isinstance(key, ID): 

177 if key.qd_type not in qd_type_list: 

178 qd_type_list.append(key.qd_type) 

179 elif isinstance(key, str): 

180 setup_list.append(key) 

181 print('Found these type of preconditioners:', qd_type_list) 

182 print('Found these setups:', setup_list) 

183 

184 assert len(qd_type_list) == 5, 'ERROR did not find four preconditioners, got %s' % qd_type_list 

185 assert len(setup_list) == 4, 'ERROR: did not find four setup, got %s' % setup_list 

186 

187 qd_type_list = ['IEpar', 'Qpar', 'MIN', 'MIN3', 'MIN_GT'] 

188 marker_list = ['s', 'o', '^', 'v', 'x'] 

189 color_list = ['r', 'g', 'b', 'c', 'm'] 

190 

191 plt_helper.setup_mpl() 

192 print('post setup') 

193 # loop over setups and Q-delta types: one figure per setup, all Qds in one plot 

194 for setup in setup_list: 

195 print('setup') 

196 plt_helper.newfig(textwidth=238.96, scale=0.89) 

197 

198 for qd_type, marker, color in zip(qd_type_list, marker_list, color_list): 

199 niter = np.zeros(len(results[setup][1])) 

200 for key in results.keys(): 

201 if isinstance(key, ID): 

202 if key.setup == setup and key.qd_type == qd_type: 

203 xvalue = results[setup][1].index(key.param) 

204 niter[xvalue] = results[key] 

205 ls = '-' 

206 lw = 1 

207 plt_helper.plt.semilogx( 

208 results[setup][1], 

209 niter, 

210 label=qd_type, 

211 lw=lw, 

212 linestyle=ls, 

213 color=color, 

214 marker=marker, 

215 markeredgecolor='k', 

216 ) 

217 

218 if setup == 'heat': 

219 xlabel = r'$\nu$' 

220 elif setup == 'advection': 

221 xlabel = r'$c$' 

222 elif setup == 'fisher': 

223 xlabel = r'$\lambda_0$' 

224 elif setup == 'vanderpol': 

225 xlabel = r'$\mu$' 

226 else: 

227 print('Setup not implemented..', setup) 

228 exit() 

229 

230 plt_helper.plt.ylim([0, 60]) 

231 plt_helper.plt.legend(loc=2, ncol=1) 

232 plt_helper.plt.ylabel('number of iterations') 

233 plt_helper.plt.xlabel(xlabel) 

234 plt_helper.plt.grid() 

235 

236 # save plot as PDF and PGF 

237 fname = 'data/parallelSDC_preconditioner_MPI_' + setup 

238 plt_helper.savefig(fname) 

239 

240 assert os.path.isfile(fname + '.pdf'), 'ERROR: plotting did not create PDF file' 

241 # assert os.path.isfile(fname + '.pgf'), 'ERROR: plotting did not create PGF file' 

242 assert os.path.isfile(fname + '.png'), 'ERROR: plotting did not create PNG file' 

243 

244 

245if __name__ == "__main__": 

246 comm = MPI.COMM_WORLD 

247 main(comm=comm) 

248 if comm.Get_rank() == 0: 

249 plot_iterations()