DSPLIB User Guide
DSPLIB_cholesky_common.cpp
Go to the documentation of this file.
1 /*******************************************************************************
2  * *
3  * module name :DSPLIB *
4  * *
5  * module descripton :Digital Signal Processing Library module for C7x+MMA *
6  * *
7  * Copyright (C) 2017-2018 Texas Instruments Incorporated - https://www.ti.com/ *
8  * ALL RIGHTS RESERVED *
9  * *
10  ******************************************************************************/
11 
23 /*******************************************************************************
24  *
25  * INCLUDES
26  *
27  ******************************************************************************/
28 #include "DSPLIB_cholesky_common.h"
29 
30 /*******************************************************************************
31  *
32  * INITIALIZATION
33  *
34  ******************************************************************************/
35 
36 template <typename dataType>
37 DSPLIB_STATUS DSPLIB_cholesky_inplace_isPosDefinite_init(int32_t order, int32_t colStride, uint8_t *pBlock)
38 {
39  DSPLIB_DEBUGPRINTFN(0, "%s\n", "Entering function");
40 
41  DSPLIB_STATUS status = DSPLIB_SUCCESS;
42 
43  typedef typename c7x::make_full_vector<dataType>::type vec;
44  int32_t eleCount = c7x::element_count_of<vec>::value;
45 
46  __SE_ELETYPE SE_ELETYPE = c7x::se_eletype<vec>::value;
47  __SE_VECLEN SE_VECLEN = c7x::se_veclen<vec>::value;
48  // __SA_VECLEN SA_VECLEN = c7x::sa_veclen<vec>::value;
49 
50  /* FOR POSITIVE DEFINITE TEST */
51  int32_t nVec = DSPLIB_ceilingDiv(order, eleCount);
52  int32_t nVecPosCheck2 = nVec / 2;
53  int32_t nVecPosCheck1 = nVec - nVecPosCheck2;
54  __SE_TEMPLATE_v1 sePosCheck1Params = __gen_SE_TEMPLATE_v1();
55  sePosCheck1Params.ICNT0 = eleCount;
56  sePosCheck1Params.ICNT1 = order;
57  sePosCheck1Params.DIM1 = colStride;
58  sePosCheck1Params.ICNT2 = nVecPosCheck1;
59  sePosCheck1Params.DIM2 = eleCount * 2;
60  sePosCheck1Params.DECDIM1 = __SE_DECDIM_DIM2;
61  sePosCheck1Params.DECDIM1_WIDTH = order;
62  sePosCheck1Params.DIMFMT = __SE_DIMFMT_3D;
63  sePosCheck1Params.VECLEN = SE_VECLEN;
64  sePosCheck1Params.ELETYPE = SE_ELETYPE;
65  *(__SE_TEMPLATE_v1 *) ((uint8_t *) pBlock + SE_SE0_PARAM_OFFSET) = sePosCheck1Params;
66 
67  __SE_TEMPLATE_v1 sePosCheck2Params = __gen_SE_TEMPLATE_v1();
68  sePosCheck2Params.ICNT0 = eleCount;
69  sePosCheck2Params.ICNT1 = order;
70  sePosCheck2Params.DIM1 = colStride;
71  sePosCheck2Params.ICNT2 = nVecPosCheck2;
72  sePosCheck2Params.DIM2 = eleCount * 2;
73  sePosCheck2Params.DECDIM1 = __SE_DECDIM_DIM2;
74  sePosCheck2Params.DECDIM1_WIDTH = order - eleCount;
75  sePosCheck2Params.DIMFMT = __SE_DIMFMT_3D;
76  sePosCheck2Params.VECLEN = SE_VECLEN;
77  sePosCheck2Params.ELETYPE = SE_ELETYPE;
78 
79  *(__SE_TEMPLATE_v1 *) ((uint8_t *) pBlock + SE_SE1_PARAM_OFFSET) = sePosCheck2Params;
80 
81  DSPLIB_DEBUGPRINTFN(0, "%s\n", "Exiting function");
82  return status;
83 }
84 
85 template DSPLIB_STATUS
86 DSPLIB_cholesky_inplace_isPosDefinite_init<float>(int32_t order, int32_t colStride, uint8_t *pBlock);
87 template DSPLIB_STATUS
88 DSPLIB_cholesky_inplace_isPosDefinite_init<double>(int32_t order, int32_t colStride, uint8_t *pBlock);
89 
90 /*******************************************************************************
91  *
92  * IMPLEMENTATION
93  *
94  ******************************************************************************/
95 
96 template <typename dataType>
97 dataType DSPLIB_cholesky_inplace_isPosDefinite(dataType *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock)
98 {
99  DSPLIB_DEBUGPRINTFN(0, "%s\n", "Entering function");
100 
101  __SE_TEMPLATE_v1 sePosCheck1Params, sePosCheck2Params;
102  sePosCheck1Params = *(__SE_TEMPLATE_v1 *) ((uint8_t *) pBlock + SE_SE0_PARAM_OFFSET);
103  sePosCheck2Params = *(__SE_TEMPLATE_v1 *) ((uint8_t *) pBlock + SE_SE1_PARAM_OFFSET);
104 
105  typedef typename c7x::make_full_vector<dataType>::type vec;
106 
107  /* test A for positive definite matrix: */
108  /* z_transpose*A*z>0 where z=1,2,...order */
109  int32_t nVec = DSPLIB_ceilingDiv(order, eleCount);
110 
111  vec vZTrans1 = (vec) 1.0;
112  vec vZTrans2 = (vec) 2.0;
113  vec vZTrans3 = (vec) 3.0;
114  vec vZTrans4 = (vec) 4.0;
115 
116  vec vOne = (vec) 1.0;
117  vec vTwo = (vec) 2.0;
118  vec vThree = (vec) 3.0;
119  vec vFour = (vec) 4.0;
120  vec vZero = (vec) 0.0;
121  vec vZ1 = vZero;
122  vec sum1 = vZero;
123  vec sum2 = vZero;
124  vec sum3 = vZero;
125  vec sum4 = vZero;
126  vec sum5 = vZero;
127  vec sum6 = vZero;
128  vec sum7 = vZero;
129  vec sum8 = vZero;
130 
131  vec finalSum1 = vZero;
132  vec finalSum2 = vZero;
133 
134  for (int32_t i = 0; i < eleCount; i++) {
135  vZ1.s[i] = (dataType) i;
136  }
137  vZ1 += vOne;
138  vec vZ2 = vZ1 + eleCount;
139 
140  __SE0_OPEN(A, sePosCheck1Params);
141  if (nVec > 1) {
142  __SE1_OPEN(A + eleCount, sePosCheck2Params);
143  }
144 
145  int32_t horizontal = 0;
146  int32_t vertical = 0;
147  for (horizontal = 0; horizontal < nVec - 1; horizontal += 2) {
148  for (vertical = 0; vertical < order - 3; vertical += 4) {
149 
150  vec v1 = c7x::strm_eng<0, vec>::get_adv();
151  vec v2 = c7x::strm_eng<1, vec>::get_adv();
152  vec v3 = c7x::strm_eng<0, vec>::get_adv();
153  vec v4 = c7x::strm_eng<1, vec>::get_adv();
154  vec v5 = c7x::strm_eng<0, vec>::get_adv();
155  vec v6 = c7x::strm_eng<1, vec>::get_adv();
156  vec v7 = c7x::strm_eng<0, vec>::get_adv();
157  vec v8 = c7x::strm_eng<1, vec>::get_adv();
158 
159  sum1 += v1 * vZTrans1;
160  sum2 += v2 * vZTrans1;
161  sum3 += v3 * vZTrans2;
162  sum4 += v4 * vZTrans2;
163  sum5 += v5 * vZTrans3;
164  sum6 += v6 * vZTrans3;
165  sum7 += v7 * vZTrans4;
166  sum8 += v8 * vZTrans4;
167 
168  vZTrans1 += vFour;
169  vZTrans2 += vFour;
170  vZTrans3 += vFour;
171  vZTrans4 += vFour;
172  }
173 
174  for (; vertical < order; vertical++) {
175  vec v1 = c7x::strm_eng<0, vec>::get_adv();
176  vec v2 = c7x::strm_eng<1, vec>::get_adv();
177 
178  sum1 += v1 * vZTrans1;
179  sum2 += v2 * vZTrans1;
180 
181  vZTrans1 += vOne;
182  }
183  sum1 += sum3;
184  sum1 += sum5;
185  sum1 += sum7;
186 
187  sum2 += sum4;
188  sum2 += sum6;
189  sum2 += sum8;
190 
191  finalSum1 += sum1 * vZ1;
192  finalSum2 += sum2 * vZ2;
193 
194  vZ1 += 2 * eleCount;
195  vZ2 += 2 * eleCount;
196 
197  vZTrans1 = vOne;
198  vZTrans2 = vTwo;
199  vZTrans3 = vThree;
200  vZTrans4 = vFour;
201 
202  sum1 = vZero;
203  sum2 = vZero;
204  sum3 = vZero;
205  sum4 = vZero;
206  sum5 = vZero;
207  sum6 = vZero;
208  sum7 = vZero;
209  sum8 = vZero;
210  }
211 
212 
213  if (horizontal != nVec) {
214  for (vertical = 0; vertical < order - 3; vertical += 4) {
215  vec v1 = c7x::strm_eng<0, vec>::get_adv();
216  vec v3 = c7x::strm_eng<0, vec>::get_adv();
217  vec v5 = c7x::strm_eng<0, vec>::get_adv();
218  vec v7 = c7x::strm_eng<0, vec>::get_adv();
219 
220  sum1 += v1 * vZTrans1;
221  sum3 += v3 * vZTrans2;
222  sum5 += v5 * vZTrans3;
223  sum7 += v7 * vZTrans4;
224 
225  vZTrans1 += vFour;
226  vZTrans2 += vFour;
227  vZTrans3 += vFour;
228  vZTrans4 += vFour;
229  }
230  for (; vertical < order; vertical++) {
231  vec v1 = c7x::strm_eng<0, vec>::get_adv();
232  sum1 += v1 * vZTrans1;
233  vZTrans1 += vOne;
234  }
235 
236  sum1 += sum3;
237  sum1 += sum5;
238  sum1 += sum7;
239 
240  finalSum1 += sum1 * vZ1;
241 
242  vZ1 += 2 * eleCount;
243  }
244 
245  vec finalSum = finalSum1 + finalSum2;
246  dataType sum = 0.0;
247  c7x_horizontal_add(finalSum, &sum);
248  __SE0_CLOSE();
249  __SE1_CLOSE();
250  DSPLIB_DEBUGPRINTFN(0, "Exiting function with return SUM: %7.5e\n", sum);
251  return sum;
252 }
253 template float
254 DSPLIB_cholesky_inplace_isPosDefinite<float>(float *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock);
255 
256 template double
257 DSPLIB_cholesky_inplace_isPosDefinite<double>(double *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock);
258 
259 /* ======================================================================== */
260 /* End of file: DSPLIB_cholesky_common.cpp */
261 /* ======================================================================== */
SE_SE0_PARAM_OFFSET
#define SE_SE0_PARAM_OFFSET
Definition: DSPLIB_add_ci.cpp:58
SE_SE1_PARAM_OFFSET
#define SE_SE1_PARAM_OFFSET
Definition: DSPLIB_bexp_ci.cpp:55
DSPLIB_cholesky_inplace_isPosDefinite_init< double >
template DSPLIB_STATUS DSPLIB_cholesky_inplace_isPosDefinite_init< double >(int32_t order, int32_t colStride, uint8_t *pBlock)
DSPLIB_cholesky_inplace_isPosDefinite< double >
template double DSPLIB_cholesky_inplace_isPosDefinite< double >(double *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock)
DSPLIB_cholesky_inplace_isPosDefinite_init
DSPLIB_STATUS DSPLIB_cholesky_inplace_isPosDefinite_init(int32_t order, int32_t colStride, uint8_t *pBlock)
Definition: DSPLIB_cholesky_common.cpp:37
DSPLIB_cholesky_inplace_isPosDefinite< float >
template float DSPLIB_cholesky_inplace_isPosDefinite< float >(float *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock)
DSPLIB_cholesky_inplace_isPosDefinite_init< float >
template DSPLIB_STATUS DSPLIB_cholesky_inplace_isPosDefinite_init< float >(int32_t order, int32_t colStride, uint8_t *pBlock)
DSPLIB_cholesky_inplace_isPosDefinite
dataType DSPLIB_cholesky_inplace_isPosDefinite(dataType *A, const int32_t order, const int32_t eleCount, uint8_t *pBlock)
Definition: DSPLIB_cholesky_common.cpp:97
DSPLIB_cholesky_common.h
DSPLIB_DEBUGPRINTFN
#define DSPLIB_DEBUGPRINTFN(N, fmt,...)
Definition: DSPLIB_types.h:83
DSPLIB_STATUS_NAME
DSPLIB_STATUS_NAME
The enumeration of all status codes.
Definition: DSPLIB_types.h:151
DSPLIB_SUCCESS
@ DSPLIB_SUCCESS
Definition: DSPLIB_types.h:152