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""" 

9import numpy as np 

10from scipy import signal 

11import matplotlib.pyplot as plt 

12 

13from pySDC.projects.parallelSDC_reloaded.utils import solutionExact 

14 

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

16 

17 

18muVals = [0.1, 2, 10] 

19muPeriods = [] 

20 

21tEnd = 20 

22nSteps = 200 

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

24 

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

26for mu in muVals: 

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

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

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

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

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

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

33 

34 x = uExact[:, 0] 

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

36 period = tVals[idx] 

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

38 muPeriods.append(period) 

39 

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

41for mu, tEnd in zip(muVals, muPeriods): 

42 nSteps = 200 

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

44 

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

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

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

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

49 print(' -- done') 

50 

51# Figure settings 

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

53 plt.figure(figName) 

54 plt.legend() 

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

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

57 plt.tight_layout()