2323/** 
2424* Computes the arithmetic mean of a single-precision floating-point strided array using Welford's algorithm. 
2525* 
26- * @param N        number of indexed elements 
27- * @param X        input array 
28- * @param strideX  stride length 
29- * @return         output value 
30- */ 
31- float  API_SUFFIX (stdlib_strided_smeanwd )( const  CBLAS_INT  N , const  float  * X , const  CBLAS_INT  strideX  ) {
32- 	const  CBLAS_INT  ox  =  stdlib_strided_stride2offset ( N , strideX  );
33- 	return  API_SUFFIX (stdlib_strided_smeanwd_ndarray )( N , X , strideX , ox  );
34- }
35- 
36- /** 
37- * Computes the arithmetic mean of a single-precision floating-point strided array using Welford's algorithm and alternative indexing semantics. 
38- * 
3926* ## Method 
4027* 
4128* -   This implementation uses Welford's algorithm for efficient computation, which can be derived as follows 
@@ -57,6 +44,19 @@ float API_SUFFIX(stdlib_strided_smeanwd)( const CBLAS_INT N, const float *X, con
5744* @param N        number of indexed elements 
5845* @param X        input array 
5946* @param strideX  stride length 
47+ * @return         output value 
48+ */ 
49+ float  API_SUFFIX (stdlib_strided_smeanwd )( const  CBLAS_INT  N , const  float  * X , const  CBLAS_INT  strideX  ) {
50+ 	const  CBLAS_INT  ox  =  stdlib_strided_stride2offset ( N , strideX  );
51+ 	return  API_SUFFIX (stdlib_strided_smeanwd_ndarray )( N , X , strideX , ox  );
52+ }
53+ 
54+ /** 
55+ * Computes the arithmetic mean of a single-precision floating-point strided array using Welford's algorithm and alternative indexing semantics. 
56+ * 
57+ * @param N        number of indexed elements 
58+ * @param X        input array 
59+ * @param strideX  stride length 
6060* @param offsetX  starting index for X 
6161* @return         output value 
6262*/ 
@@ -67,18 +67,18 @@ float API_SUFFIX(stdlib_strided_smeanwd_ndarray)( const CBLAS_INT N, const float
6767	float  mu ;
6868
6969	if  ( N  <= 0  ) {
70- 		return  0.0   / 0.0  ; // NaN 
70+ 		return  0.0f   / 0.0f  ; // NaN 
7171	}
7272	if  ( N  ==  1  ||  strideX  ==  0  ) {
7373		return  X [ offsetX  ];
7474	}
7575	ix  =  offsetX ;
76- 	mu  =  0.0  ;
77- 	n  =  0 ;
76+ 	mu  =  0.0f  ;
77+ 	n  =  0.0  ;
7878	for  ( i  =  0 ; i  <  N ; i ++  ) {
79- 		n  +=  1 ;
80- 		mu  +=  ( X [ix ]- mu  ) / ( float ) n ;
81- 		ix  +=  strideX ;
79+ 		n  +=  1.0  ;
80+ 		mu  +=  (float )(( double )(  X [ix ]- mu  ) / n ) ;
81+ 		ix  +=  stride ;
8282	}
8383	return  mu ;
8484}
0 commit comments