Skip to content

Commit 29680dd

Browse files
a-ma72a-ma72
committed
Python rainflow counting package pyrfc
1 parent ed2ecf4 commit 29680dd

30 files changed

+21055
-125
lines changed

‎.gitignore

Lines changed: 7 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,11 @@
11
/build
2-
32
rainflow\.sublime-project
4-
53
rainflow\.sublime-workspace
64
/sf
5+
/python/__pycache__
6+
/python/build
7+
/python/venv
8+
/python/dist
9+
/python/pyrfc*
10+
/python/MANIFEST
11+
/src/config.h

‎CMakeLists.txt

Lines changed: 9 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@ endif()
2222
# Project name and version
2323
project( rainflow )
2424
set( RFC_VERSION_MAJOR "0" )
25-
set( RFC_VERSION_MINOR "4" )
25+
set( RFC_VERSION_MINOR "5" )
2626

2727

2828
# Build type and compiler selection
@@ -53,6 +53,7 @@ if( CMAKE_CURRENT_SOURCE_DIR STREQUAL CMAKE_SOURCE_DIR )
5353
option( RFC_USE_DELEGATES "Use delegates (functors)" ON )
5454
option( RFC_GLOBAL_EXTREMA "Always calculate global extrema" ON )
5555
option( RFC_HCM_SUPPORT "Support HCM (Clormann/Seeger) algorithm" ON )
56+
option( RFC_ASTM_SUPPORT "Support ASTM E 1049-85 algorithm" ON )
5657
option( RFC_TP_SUPPORT "Support turning points" ON )
5758
option( RFC_DH_SUPPORT "Support \"spread damage\" over turning points and damage history" ON )
5859
option( RFC_AT_SUPPORT "Support amplitude transformation regarding mean load influence on fatigue strength" ON )
@@ -64,11 +65,14 @@ else()
6465
message( STATUS "Build ${PROJECT_NAME} as subsequent project")
6566
endif()
6667

68+
# Save options in configuration file
69+
add_definitions( -DRFC_HAVE_CONFIG_H )
70+
configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/src/config.h.in ${CMAKE_CURRENT_SOURCE_DIR}/src/config.h )
71+
6772

6873
if( CMAKE_BUILD_TYPE MATCHES DEBUG )
6974
endif()
7075

71-
7276
# Compiler dependencies
7377
if( MSVC )
7478
# Turn off misguided "secure CRT" warnings in MSVC.
@@ -85,9 +89,6 @@ if( NOT LIBM_LIBRARY )
8589
endif()
8690

8791

88-
# Configuration file
89-
configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/config.h.in ${CMAKE_CURRENT_SOURCE_DIR}/config.h )
90-
9192
# MATLAB
9293
if( RFC_EXPORT_MEX )
9394
set( MATLAB_FIND_DEBUG 1 )
@@ -124,20 +125,20 @@ if( RFC_EXPORT_MEX )
124125
include_directories( ${Matlab_INCLUDE_DIRS} )
125126

126127
# MEX function (MATLAB)
127-
matlab_add_mex( NAME ${PROJECT_NAME} SRC rainflow.c OUTPUT_NAME rfc )
128+
matlab_add_mex( NAME ${PROJECT_NAME} SRC src/rainflow.c OUTPUT_NAME rfc )
128129
target_compile_definitions( ${PROJECT_NAME} PRIVATE MATLAB_MEX_FILE _SCL_SECURE_NO_WARNINGS )
129130
target_link_libraries( ${PROJECT_NAME} ${Matlab_LIBRARIES} ${LIBM_LIBRARY} )
130131
# install to /bin by default
131132
install( TARGETS ${PROJECT_NAME} RUNTIME DESTINATION bin LIBRARY DESTINATION bin )
132133
endif()
133134

134135
# Static rainflow library
135-
add_library( rfc STATIC rainflow.c )
136+
add_library( rfc STATIC src/rainflow.c )
136137
target_link_libraries( rfc ${LIBM_LIBRARY} )
137138

138139
# Test application, start project for MSVC
139140
if( RFC_TEST )
140-
add_executable( rfc_test rainflow.c test/rfc_test.c test/rfc_wrapper_simple.cpp test/rfc_wrapper_advanced.cpp )
141+
add_executable( rfc_test src/rainflow.c test/rfc_test.c test/rfc_wrapper_simple.cpp test/rfc_wrapper_advanced.cpp )
141142
target_compile_definitions( rfc_test PRIVATE _SCL_SECURE_NO_WARNINGS )
142143
target_link_libraries( rfc_test ${LIBM_LIBRARY} )
143144
target_sources( rfc_test PUBLIC greatest/greatest.h )

‎README.md

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
## Rainflow Counting Algorithm (4-point-method), C99 compliant
1+
## Rainflow Counting Algorithm (4-point-methods), C99 compliant
22

33
### "Rainflow Counting" consists of four main steps:
44

@@ -17,7 +17,7 @@
1717

1818
These steps are fully documented in standards such as
1919
ASTM E1049 "Standard Practices for Cycle Counting in Fatigue Analysis" [1]
20-
This implementation uses the 4-point algorithm mentioned in [3,4] and the 3-point HCM method proposed in [2]
20+
This implementation uses the 4-point algorithm mentioned in [3,4] and the 3-point HCM method proposed in [2] as well as the ASTM E 1049 (2011) standard in [1].
2121
To take the residue into account, you may implement a custom method or use some
2222
predefined functions.
2323

‎matlab/validate.m

Lines changed: 87 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -23,21 +23,21 @@
2323
[class_width, ...
2424
class_offset] = class_param( x, class_count );
2525
hysteresis = class_width;
26-
enforce_margin = 0;
27-
use_hcm = 0;
26+
enforce_margin = 0; % First and last data point may be excluded in tp
27+
use_hcm = 0; % Use 4 point method, not HCM
28+
use_astm = 0; % Use 4 point method, not ASTM
2829
residual_method = 0;
2930
spread_damage = 0;
3031

3132
[~,re,rm] = rfc( 'rfc', x, class_count, class_width, class_offset, hysteresis, ...
32-
residual_method, enforce_margin, use_hcm, spread_damage );
33+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
3334

3435
assert( sum( sum( rm ) ) == 0 );
3536

3637
assert( isempty(re) );
3738

3839
save( name, 'rm', 're' );
3940

40-
4141
%% One single cycle (up)
4242
name = 'one_cycle_up';
4343
class_count = 4;
@@ -47,13 +47,14 @@
4747
[class_width, ...
4848
class_offset] = class_param( x, class_count );
4949
hysteresis = class_width * 0.99;
50-
enforce_margin = 0;
51-
use_hcm = 0;
50+
enforce_margin = 0; % First and last data point may be excluded in tp
51+
use_hcm = 0; % Use 4 point method, not HCM
52+
use_astm = 0; % Use 4 point method, not ASTM
5253
residual_method = 0;
5354
spread_damage = 0;
5455

5556
[~,re,rm] = rfc( 'rfc', x, class_count, class_width, class_offset, hysteresis, ...
56-
residual_method, enforce_margin, use_hcm, spread_damage );
57+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
5758

5859
assert( sum( sum( rm ) ) == 1 );
5960
assert( rm( 3,2 ) == 1 );
@@ -62,24 +63,23 @@
6263

6364
save( name, 'rm', 're' );
6465

65-
6666
%% One single cycle (down)
6767
name = 'one_cycle_down';
6868
class_count = 4;
6969
x = export_series( name, [4,2,3,1], class_count );
7070
x_max = 4;
7171
x_min = 1;
72-
class_count = 4;
7372
[class_width, ...
7473
class_offset] = class_param( x, class_count );
7574
hysteresis = class_width * 0.99;
76-
enforce_margin = 0;
77-
use_hcm = 0;
75+
enforce_margin = 0; % First and last data point may be excluded in tp
76+
use_hcm = 0; % Use 4 point method, not HCM
77+
use_astm = 0; % Use 4 point method, not ASTM
7878
residual_method = 0;
7979
spread_damage = 0;
8080

8181
[~,re,rm] = rfc( 'rfc', x, class_count, class_width, class_offset, hysteresis, ...
82-
residual_method, enforce_margin, use_hcm, spread_damage );
82+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
8383

8484
assert( sum( sum( rm ) ) == 1 );
8585
assert( rm( 2,3 ) == 1 );
@@ -88,7 +88,6 @@
8888

8989
save( name, 'rm', 're' );
9090

91-
9291
%% Small example, taken from url:
9392
% [https://community.plm.automation.siemens.com/t5/Testing-Knowledge-Base/Rainflow-Counting/ta-p/383093]
9493
name = 'small_example';
@@ -99,13 +98,14 @@
9998
[class_width, ...
10099
class_offset] = class_param( x, class_count );
101100
hysteresis = class_width;
102-
enforce_margin = 0;
103-
use_hcm = 0;
101+
enforce_margin = 0; % First and last data point may be excluded in tp
102+
use_hcm = 0; % Use 4 point method, not HCM
103+
use_astm = 0; % Use 4 point method, not ASTM
104104
residual_method = 0;
105105
spread_damage = 0;
106106

107107
[~,re,rm] = rfc( 'rfc', x, class_count, class_width, class_offset, hysteresis, ...
108-
residual_method, enforce_margin, use_hcm, spread_damage );
108+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
109109

110110
assert( sum( sum( rm ) ) == 7 );
111111
assert( rm( 5,3 ) == 2 );
@@ -118,7 +118,6 @@
118118

119119
save( name, 'rm', 're' );
120120

121-
122121
%% Long data series
123122
rng(5,'twister') % Init random seed
124123
name = 'long_series';
@@ -130,16 +129,17 @@
130129
x_int = int16(round(xx));
131130
x = double(x_int);
132131
hysteresis = class_width;
133-
enforce_margin = 1;
134-
use_hcm = 0;
132+
enforce_margin = 1; % Enforce first and last data point included in tp
133+
use_hcm = 0; % Use 4 point method, not HCM
134+
use_astm = 0; % Use 4 point method, not ASTM
135135
residual_method = 0; % 0=RFC_RES_NONE, 7=RFC_RES_REPEATED
136136
spread_damage = 1; % 0=RFC_SD_HALF_23, 1=RFC_SD_RAMP_AMPLITUDE_23
137137

138138
assert( sum(abs(xx-x)) < 1e-4 );
139139

140140
[pd,re,rm,rp,lc,tp,dh] = ...
141141
rfc( 'rfc', double(x_int), class_count, class_width, class_offset, hysteresis, ...
142-
residual_method, enforce_margin, use_hcm, spread_damage );
142+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
143143

144144
% With residuum: pd == 9.8934e-06 (repeated)
145145
% Without residuum: pd == 1.1486e-07
@@ -159,7 +159,7 @@
159159

160160
[pd,re,rm,rp,lc,tp,dh] = ...
161161
rfc( 'rfc', x, class_count, class_width, class_offset, hysteresis, ...
162-
residual_method, enforce_margin, use_hcm, spread_damage );
162+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
163163

164164
assert( abs( sum( dh ) / pd - 1 ) < 1e-10 );
165165
hold( ax(2), 'all' );
@@ -198,13 +198,79 @@
198198
assert( all( test < 1e-1 ))
199199
end
200200

201+
%% Compare with ASTM E 1049-85 (MATLAB)
202+
if 0
203+
enforce_margin = 1; % Enforce first and last data point included in tp
204+
use_hcm = 0; % Use 4 point method, not HCM
205+
use_astm = 0; % Use 4 point method, not ASTM
206+
residual_method = 4; % 4=ASTM related
207+
spread_damage = 0; % 0=RFC_SD_HALF_23, 1=RFC_SD_RAMP_AMPLITUDE_23
208+
class_count = 1000;
209+
class_width = 5;
210+
class_offset = -2025;
211+
hysteresis = class_width;
212+
213+
if 1
214+
load long_series_csv.mat
215+
else
216+
x_int = [2,5,3,6,2,4,1,6,1,4,1,5,3,6,3,6,1,5,2];
217+
x_max = max(x_int)+0.5;
218+
x_min = min(x_int)-0.5;
219+
class_width = 1;
220+
class_offset = x_min;
221+
hysteresis = 0;
222+
class_count = (x_max - x_min) / class_width;
223+
end
224+
225+
residual_method = 7; % 7=repeated
226+
[pd7,re7,rm7,rp7,lc7,tp7] = ...
227+
rfc( 'rfc', double(x_int), class_count, class_width, class_offset, hysteresis, ...
228+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
229+
230+
residual_method = 4; % 4=ASTM related
231+
[pd4,re4,rm4,rp4,lc4,tp4] = ...
232+
rfc( 'rfc', double(x_int), class_count, class_width, class_offset, hysteresis, ...
233+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
234+
235+
use_astm = 1; % Use ASTM algorithm
236+
residual_method = 4; % 4=ASTM related
237+
[pd4a,re4a,rm4a,rp4a,lc4a,tp4a] = ...
238+
rfc( 'rfc', double(x_int), class_count, class_width, class_offset, hysteresis, ...
239+
residual_method, enforce_margin, use_hcm, use_astm, spread_damage );
240+
241+
% MATLAB - Rainflow counts for fatigue analysis (according to ASTM E 1049)
242+
c = rainflow( tp4a(:,2), 'ext' );
243+
edges = (0:class_count) .* class_width;
244+
[~,bin] = histc( c(:,2), edges );
245+
N = accumarray( bin, c(:,1) );
246+
N(class_count+1) = 0;
247+
N = N(1:class_count);
248+
Range = (edges(1:end-1)+edges(2:end)) / 2;
249+
pd_astm = sum( N(:)' ./ (1e7*(Range/2/1e3).^-5) );
250+
251+
figure
252+
plot( cumsum(rp7, 'reverse'), edges(1:end-1), 'g-', 'Disp', '4-point method, res=repeated' ), hold all
253+
plot( cumsum(rp4, 'reverse'), edges(1:end-1), 'k--', 'Disp', '4-point method, res=half cycles' ), hold all
254+
plot( cumsum(rp4a, 'reverse'), edges(1:end-1), 'b--', 'Disp', '3-point method, ASTM E 1049-85, res=half cycles' ), hold all
255+
plot( cumsum(N, 'reverse'), edges(1:end-1), 'r-', 'Disp', '3-point method, ASTM E 1049-85 (MATLAB)' ), hold all
256+
set( gca, 'XScale', 'log' );
257+
set( gca, 'YScale', 'log' );
258+
xlim( [0.9 1e3] );
259+
xlabel( 'Counts' );
260+
ylabel( 'Range (normalized, prepared)' );
261+
grid
262+
legend
263+
fprintf( 'Damage ratio 4pt,repeated vs. ASTM: %g%%\n', pd7/pd_astm );
264+
end
265+
201266
%% Long series, turning points only
202267
y = rfc( 'turningpoints', x, class_width*2, enforce_margin );
203268
figure
204269
plot( x, 'k-', 'DisplayName', 'time series' );
205270
hold all
206271
plot( y(2,:), y(1,:), 'r--', 'DisplayName', 'turning points' );
207272
legend( 'show' );
273+
208274
end
209275

210276
function rounded_data = export_series( filename, data, class_count, format )

‎python/CMakeLists.txt

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,18 @@
1+
project(pyrfc)
2+
3+
# Find Python 3
4+
find_package(Python3 COMPONENTS Interpreter Development)
5+
6+
# We also need to use the NumPy C API
7+
exec_program(${_Python3_EXECUTABLE}
8+
ARGS "-c \"import numpy; print(numpy.get_include())\""
9+
OUTPUT_VARIABLE NUMPY_INCLUDE_DIR
10+
RETURN_VALUE NUMPY_NOT_FOUND
11+
)
12+
if(NUMPY_NOT_FOUND)
13+
message(FATAL_ERROR "NumPy headers not found")
14+
endif()
15+
16+
add_library(pyrfc SHARED pyrfc.cpp)
17+
target_include_directories(pyrfc PRIVATE ${_Python3_INCLUDE_DIR})
18+
target_link_libraries(pyrfc PRIVATE ${_Python3_LIBRARY_RELEASE})

‎python/MANIFEST.in

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,3 @@
1+
include src/rainflow.c
2+
include src/rainflow.h
3+
include src/rainflow.hpp

‎python/README.md

Lines changed: 53 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,53 @@
1+
## **pyrfc** a rainflow counting algorithm Python package
2+
3+
### "Rainflow Counting" consists of four main steps:
4+
5+
1. Hysteresis Filtering
6+
2. Peak-Valley Filtering
7+
3. Discretization
8+
4. Four Point Counting Method:
9+
10+
* D
11+
/ \ Closed, if min(B,C) >= min(A,D) && max(B,C) <= max(A,D)
12+
B *<--/ Slope B-C is counted and removed from residue
13+
/ \ /
14+
/ * C
15+
\ /
16+
* A
17+
18+
These steps are fully documented in standards such as
19+
ASTM E1049 "Standard Practices for Cycle Counting in Fatigue Analysis" [1]
20+
This implementation uses the 4-point algorithm mentioned in [3,4] and the 3-point HCM method proposed in [2] as well as the ASTM E 1049 (2011) standard in [1].
21+
To take the residue into account, you may implement a custom method or use some
22+
predefined functions.
23+
24+
### Install
25+
pip install _{packagename}_.tar.gz
26+
where _{packagename}_ is the current package release, for example:
27+
28+
pip install pyrfc-0.1.0.tar.gz
29+
30+
### Test
31+
_pyrfc_ packages include some unit tests, which can be run:
32+
33+
python run_tests.py
34+
35+
### Examples
36+
For a quick introduction you can run and inspect a small example:
37+
python run_examples.py
38+
39+
![](jupyter_screenshot.png)
40+
41+
---
42+
### References:
43+
[1] "Standard Practices for Cycle Counting in Fatigue Analysis."
44+
ASTM Standard E 1049, 1985 (2011).
45+
West Conshohocken, PA: ASTM International, 2011.
46+
[2] "Rainflow - HCM / Ein Hysteresisschleifen-Zaehlalgorithmus auf werkstoffmechanischer Grundlage"
47+
U.H. Clormann, T. Seeger
48+
1985 TU Darmstadt, Fachgebiet Werkstoffmechanik
49+
[3] "Zaehlverfahren zur Bildung von Kollektiven und Matrizen aus Zeitfunktionen"
50+
FVA-Richtlinie, 2010.
51+
[https://fva-net.de/fileadmin/content/Richtlinien/FVA-Richtlinie_Zaehlverfahren_2010.pdf]
52+
[4] Siemens Product Lifecycle Management Software Inc., 2018.
53+
[https://community.plm.automation.siemens.com/t5/Testing-Knowledge-Base/Rainflow-Counting/ta-p/383093]

‎python/jupyter_screenshot.png

108 KB
Loading

‎python/requirements.txt

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,4 @@
1+
matplotlib==3.4.2
2+
numpy==1.21.1
3+
pandas==1.3.1
4+
seaborn==0.11.1

0 commit comments

Comments
 (0)