Coverage for pySDC / projects / parallelSDC_reloaded / vanderpol_setup.py: 100%

1#!/usr/bin/env python3 

2# -*- coding: utf-8 -*- 

3""" 

4Created on Sun Nov 12 21:35:05 2023 

5 

6Script to numerically determine periods of Van der Pol oscillation for 

7different mu parameters. 

8""" 

9 

10import numpy as np 

11from scipy import signal 

12import matplotlib.pyplot as plt 

13 

14from pySDC.projects.parallelSDC_reloaded.utils import solutionExact 

15 

16script = __file__.split('/')[-1].split('.')[0] 

17 

18 

19muVals = [0.1, 2, 10] 

20muPeriods = [] 

21 

22tEnd = 20 

23nSteps = 200 

24tVals = np.linspace(0, tEnd, nSteps + 1) 

25 

26# Compute and plot unscaled solution to determined period for each mu 

27for mu in muVals: 

28 print(f"Computing exact solution up to t={tEnd} for mu={mu} ...") 

29 uExact = solutionExact(tEnd, nSteps, "VANDERPOL", mu=mu) 

30 plt.figure(f"{script}_traj") 

31 plt.plot(tVals, uExact[:, 0], '-', label=f"$\\mu={mu}$") 

32 plt.figure(f"{script}_accel") 

33 plt.plot(tVals, uExact[:, 1], '-', label=f"$\\mu={mu}$") 

34 

35 x = uExact[:, 0] 

36 idx = signal.find_peaks(x)[0][0] 

37 period = tVals[idx] 

38 print(f" -- done, found period={period:.1f}") 

39 muPeriods.append(period) 

40 

41# Compute and plot solution for each mu on one period, scale time with period 

42for mu, tEnd in zip(muVals, muPeriods, strict=True): 

43 nSteps = 200 

44 tVals = np.linspace(0, tEnd, nSteps + 1) 

45 

46 print(f"Computing exact solution up to t={tEnd:.1f} for mu={mu} ...") 

47 uExact = solutionExact(tEnd, nSteps, "VANDERPOL", mu=mu) 

48 plt.figure(f"{script}_traj_scaled") 

49 plt.plot(tVals / tEnd, uExact[:, 0], '-', label=f"$\\mu={mu}$") 

50 print(' -- done') 

51 

52# Figure settings 

53for figName in [f"{script}_traj", f"{script}_accel", f"{script}_traj_scaled"]: 

54 plt.figure(figName) 

55 plt.legend() 

56 plt.xlabel("time (scaled)" if "scaled" in figName else "time") 

57 plt.ylabel("trajectory" if "traj" in figName else "acceleration") 

58 plt.tight_layout()