Coverage for /usr/share/miniconda3/envs/dolfin/lib/python3.8/site-packages/block/testing.py: 64%

1import atexit 

2import inspect 

3import os 

4import pathlib 

5import pickle 

6from urllib.parse import quote_plus 

7 

8import numpy as np 

9from dolfin import MPI 

10 

11 

12class RegressionFailure(Exception): 

13 pass 

14 

15 

16def _regr_root(): 

17 return pathlib.Path(__file__).parent.absolute() / 'data/regression' 

18 

19 

20_errs = 0 

21 

22 

23def _log_or_raise(msg): 

24 if os.environ.get('BLOCK_REGRESSION_ABORT'): 

25 raise RegressionFailure(msg) 

26 else: 

27 print('!', msg) 

28 global _errs 

29 _errs += 1 

30 if _errs == 1: 

31 @atexit.register 

32 def print_msg(): 

33 print(f'! {_errs} expected test(s) failed') 

34 print( 

35 f' Remove file(s) in {_regr_root()}/ and re-run to store new values as reference') 

36 

37 

38def check_expected(name, vec, show=False, rtol=1e-6, itol=0.1, prefix=None, expected=None): 

39 itol += 1 

40 if prefix is None: 40 ↛ 43line 40 didn't jump to line 43, because the condition on line 40 was never false

41 prefix = pathlib.Path(inspect.stack()[1].filename).stem 

42 

43 def _decimate(v): 

44 # Pickle not supported for block_vec/dolfin_vec, convert to numpy 

45 if hasattr(v, 'get_local'): 45 ↛ 47line 45 didn't jump to line 47, because the condition on line 45 was never false

46 v = v.get_local() 

47 if not isinstance(v, np.ndarray): 

48 v = np.concatenate(v) 

49 # To save disk&repo space, we decimate the vector and calculate mean+rms 

50 # within each chunk. This is quite arbitrary, but is intended to be robust 

51 # decimation which still catches most regression scenarios in practice. 

52 chunks = np.arange(0, len(v), 10) 

53 divisor = np.add.reduceat(np.ones_like(v), chunks) 

54 v_mean = np.add.reduceat(v, chunks) / divisor 

55 v_rms = np.sqrt(np.add.reduceat(v**2, chunks) / divisor) 

56 return (v_mean + v_rms) / 2 

57 

58 def _l2(v): 

59 if np.isscalar(v): 

60 return v 

61 elif hasattr(v, 'norm'): 

62 return v.norm('l2') 

63 else: 

64 return np.sqrt(np.sum(v**2)) 

65 

66 fname = _regr_root() / f'{quote_plus(prefix)}.{quote_plus(name)}.pickle' 

67 if fname.exists(): 67 ↛ 71line 67 didn't jump to line 71, because the condition on line 67 was never false

68 with open(fname, 'rb') as f: 

69 data = pickle.load(f) 

70 else: 

71 data = {} 

72 

73 cur_norm = _l2(vec) 

74 cur_iter = getattr(vec, '_regr_test_niter', None) 

75 ref_iter = data.get('iter') 

76 if expected is None: 

77 cur_vec = _decimate(vec) 

78 ref_vec = data.get('vec') 

79 ref_norm = data.get('norm') 

80 else: 

81 cur_vec = vec 

82 ref_vec = expected 

83 ref_norm = _l2(expected) 

84 is_serial = (MPI.size(MPI.comm_world) == 1) 

85 

86 if is_serial and ref_vec is not None: 

87 try: 

88 err_norm = _l2(cur_vec - ref_vec) 

89 except Exception as e: 

90 _log_or_raise(f'Failed to compute error norm: {e}') 

91 else: 

92 rdiff = err_norm / max(ref_norm, len(cur_vec)) 

93 if show: 

94 print(f'Norm of {name}: {cur_norm:.4f}, error: {err_norm:.4g}') 

95 if rdiff > rtol: 95 ↛ 96line 95 didn't jump to line 96, because the condition on line 95 was never true

96 _log_or_raise( 

97 f'Error in {name}: {err_norm:.4g} ({rdiff:.3g} > rtol {rtol:.3g}) ({prefix})') 

98 else: 

99 if show: 

100 print(f'Norm of {name}: {cur_norm:.4f}') 

101 if ref_norm is not None: 101 ↛ 106line 101 didn't jump to line 106, because the condition on line 101 was never false

102 rdiff = abs(cur_norm - ref_norm) / max(ref_norm, len(cur_vec)) 

103 if rdiff > rtol: 103 ↛ 104line 103 didn't jump to line 104, because the condition on line 103 was never true

104 _log_or_raise( 

105 f'Norm of {name} {cur_norm:.6g} != {ref_norm:.6g} ({rdiff:.3g} > rtol {rtol:.3g})') 

106 if is_serial and ref_iter is not None: 

107 if cur_iter is None or not ref_iter/itol <= cur_iter <= ref_iter*itol: 107 ↛ 108line 107 didn't jump to line 108, because the condition on line 107 was never true

108 _log_or_raise( 

109 f'Solver for {name} used {cur_iter} != ({ref_iter}/{itol}--{ref_iter}*{itol}) iterations ({prefix})') 

110 

111 if not fname.exists(): 111 ↛ 112line 111 didn't jump to line 112, because the condition on line 111 was never true

112 if is_serial: 

113 data = {} 

114 if cur_iter is not None: 

115 data['iter'] = cur_iter 

116 if expected is None: 

117 data['norm'] = cur_norm 

118 data['vec'] = cur_vec 

119 if data: 

120 with open(fname, 'wb') as f: 

121 pickle.dump(data, f) 

122 else: 

123 if MPI.rank(MPI.comm_world) == 0: 

124 print('Not saving regression results in parallel') 

125 

126 return vec