You can not select more than 25 topics Topics must start with a chinese character,a letter or number, can include dashes ('-') and can be up to 35 characters long.

c_cblat1c.c 46 kB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289
  1. #include <math.h>
  2. #include <stdlib.h>
  3. #include <string.h>
  4. #include <stdio.h>
  5. #include <complex.h>
  6. #ifdef complex
  7. #undef complex
  8. #endif
  9. #ifdef I
  10. #undef I
  11. #endif
  12. #if defined(_WIN64)
  13. typedef long long BLASLONG;
  14. typedef unsigned long long BLASULONG;
  15. #else
  16. typedef long BLASLONG;
  17. typedef unsigned long BLASULONG;
  18. #endif
  19. #ifdef LAPACK_ILP64
  20. typedef BLASLONG blasint;
  21. #if defined(_WIN64)
  22. #define blasabs(x) llabs(x)
  23. #else
  24. #define blasabs(x) labs(x)
  25. #endif
  26. #else
  27. typedef int blasint;
  28. #define blasabs(x) abs(x)
  29. #endif
  30. typedef blasint integer;
  31. typedef unsigned int uinteger;
  32. typedef char *address;
  33. typedef short int shortint;
  34. typedef float real;
  35. typedef double doublereal;
  36. typedef struct { real r, i; } complex;
  37. typedef struct { doublereal r, i; } doublecomplex;
  38. #ifdef _MSC_VER
  39. static inline _Fcomplex Cf(complex *z) {_Fcomplex zz={z->r , z->i}; return zz;}
  40. static inline _Dcomplex Cd(doublecomplex *z) {_Dcomplex zz={z->r , z->i};return zz;}
  41. static inline _Fcomplex * _pCf(complex *z) {return (_Fcomplex*)z;}
  42. static inline _Dcomplex * _pCd(doublecomplex *z) {return (_Dcomplex*)z;}
  43. #else
  44. static inline _Complex float Cf(complex *z) {return z->r + z->i*_Complex_I;}
  45. static inline _Complex double Cd(doublecomplex *z) {return z->r + z->i*_Complex_I;}
  46. static inline _Complex float * _pCf(complex *z) {return (_Complex float*)z;}
  47. static inline _Complex double * _pCd(doublecomplex *z) {return (_Complex double*)z;}
  48. #endif
  49. #define pCf(z) (*_pCf(z))
  50. #define pCd(z) (*_pCd(z))
  51. typedef int logical;
  52. typedef short int shortlogical;
  53. typedef char logical1;
  54. typedef char integer1;
  55. #define TRUE_ (1)
  56. #define FALSE_ (0)
  57. /* Extern is for use with -E */
  58. #ifndef Extern
  59. #define Extern extern
  60. #endif
  61. /* I/O stuff */
  62. typedef int flag;
  63. typedef int ftnlen;
  64. typedef int ftnint;
  65. /*external read, write*/
  66. typedef struct
  67. { flag cierr;
  68. ftnint ciunit;
  69. flag ciend;
  70. char *cifmt;
  71. ftnint cirec;
  72. } cilist;
  73. /*internal read, write*/
  74. typedef struct
  75. { flag icierr;
  76. char *iciunit;
  77. flag iciend;
  78. char *icifmt;
  79. ftnint icirlen;
  80. ftnint icirnum;
  81. } icilist;
  82. /*open*/
  83. typedef struct
  84. { flag oerr;
  85. ftnint ounit;
  86. char *ofnm;
  87. ftnlen ofnmlen;
  88. char *osta;
  89. char *oacc;
  90. char *ofm;
  91. ftnint orl;
  92. char *oblnk;
  93. } olist;
  94. /*close*/
  95. typedef struct
  96. { flag cerr;
  97. ftnint cunit;
  98. char *csta;
  99. } cllist;
  100. /*rewind, backspace, endfile*/
  101. typedef struct
  102. { flag aerr;
  103. ftnint aunit;
  104. } alist;
  105. /* inquire */
  106. typedef struct
  107. { flag inerr;
  108. ftnint inunit;
  109. char *infile;
  110. ftnlen infilen;
  111. ftnint *inex; /*parameters in standard's order*/
  112. ftnint *inopen;
  113. ftnint *innum;
  114. ftnint *innamed;
  115. char *inname;
  116. ftnlen innamlen;
  117. char *inacc;
  118. ftnlen inacclen;
  119. char *inseq;
  120. ftnlen inseqlen;
  121. char *indir;
  122. ftnlen indirlen;
  123. char *infmt;
  124. ftnlen infmtlen;
  125. char *inform;
  126. ftnint informlen;
  127. char *inunf;
  128. ftnlen inunflen;
  129. ftnint *inrecl;
  130. ftnint *innrec;
  131. char *inblank;
  132. ftnlen inblanklen;
  133. } inlist;
  134. #define VOID void
  135. union Multitype { /* for multiple entry points */
  136. integer1 g;
  137. shortint h;
  138. integer i;
  139. /* longint j; */
  140. real r;
  141. doublereal d;
  142. complex c;
  143. doublecomplex z;
  144. };
  145. typedef union Multitype Multitype;
  146. struct Vardesc { /* for Namelist */
  147. char *name;
  148. char *addr;
  149. ftnlen *dims;
  150. int type;
  151. };
  152. typedef struct Vardesc Vardesc;
  153. struct Namelist {
  154. char *name;
  155. Vardesc **vars;
  156. int nvars;
  157. };
  158. typedef struct Namelist Namelist;
  159. #define abs(x) ((x) >= 0 ? (x) : -(x))
  160. #define dabs(x) (fabs(x))
  161. #define f2cmin(a,b) ((a) <= (b) ? (a) : (b))
  162. #define f2cmax(a,b) ((a) >= (b) ? (a) : (b))
  163. #define dmin(a,b) (f2cmin(a,b))
  164. #define dmax(a,b) (f2cmax(a,b))
  165. #define bit_test(a,b) ((a) >> (b) & 1)
  166. #define bit_clear(a,b) ((a) & ~((uinteger)1 << (b)))
  167. #define bit_set(a,b) ((a) | ((uinteger)1 << (b)))
  168. #define abort_() { sig_die("Fortran abort routine called", 1); }
  169. #define c_abs(z) (cabsf(Cf(z)))
  170. #define c_cos(R,Z) { pCf(R)=ccos(Cf(Z)); }
  171. #ifdef _MSC_VER
  172. #define c_div(c, a, b) {Cf(c)._Val[0] = (Cf(a)._Val[0]/Cf(b)._Val[0]); Cf(c)._Val[1]=(Cf(a)._Val[1]/Cf(b)._Val[1]);}
  173. #define z_div(c, a, b) {Cd(c)._Val[0] = (Cd(a)._Val[0]/Cd(b)._Val[0]); Cd(c)._Val[1]=(Cd(a)._Val[1]/Cd(b)._Val[1]);}
  174. #else
  175. #define c_div(c, a, b) {pCf(c) = Cf(a)/Cf(b);}
  176. #define z_div(c, a, b) {pCd(c) = Cd(a)/Cd(b);}
  177. #endif
  178. #define c_exp(R, Z) {pCf(R) = cexpf(Cf(Z));}
  179. #define c_log(R, Z) {pCf(R) = clogf(Cf(Z));}
  180. #define c_sin(R, Z) {pCf(R) = csinf(Cf(Z));}
  181. //#define c_sqrt(R, Z) {*(R) = csqrtf(Cf(Z));}
  182. #define c_sqrt(R, Z) {pCf(R) = csqrtf(Cf(Z));}
  183. #define d_abs(x) (fabs(*(x)))
  184. #define d_acos(x) (acos(*(x)))
  185. #define d_asin(x) (asin(*(x)))
  186. #define d_atan(x) (atan(*(x)))
  187. #define d_atn2(x, y) (atan2(*(x),*(y)))
  188. #define d_cnjg(R, Z) { pCd(R) = conj(Cd(Z)); }
  189. #define r_cnjg(R, Z) { pCf(R) = conjf(Cf(Z)); }
  190. #define d_cos(x) (cos(*(x)))
  191. #define d_cosh(x) (cosh(*(x)))
  192. #define d_dim(__a, __b) ( *(__a) > *(__b) ? *(__a) - *(__b) : 0.0 )
  193. #define d_exp(x) (exp(*(x)))
  194. #define d_imag(z) (cimag(Cd(z)))
  195. #define r_imag(z) (cimagf(Cf(z)))
  196. #define d_int(__x) (*(__x)>0 ? floor(*(__x)) : -floor(- *(__x)))
  197. #define r_int(__x) (*(__x)>0 ? floor(*(__x)) : -floor(- *(__x)))
  198. #define d_lg10(x) ( 0.43429448190325182765 * log(*(x)) )
  199. #define r_lg10(x) ( 0.43429448190325182765 * log(*(x)) )
  200. #define d_log(x) (log(*(x)))
  201. #define d_mod(x, y) (fmod(*(x), *(y)))
  202. #define u_nint(__x) ((__x)>=0 ? floor((__x) + .5) : -floor(.5 - (__x)))
  203. #define d_nint(x) u_nint(*(x))
  204. #define u_sign(__a,__b) ((__b) >= 0 ? ((__a) >= 0 ? (__a) : -(__a)) : -((__a) >= 0 ? (__a) : -(__a)))
  205. #define d_sign(a,b) u_sign(*(a),*(b))
  206. #define r_sign(a,b) u_sign(*(a),*(b))
  207. #define d_sin(x) (sin(*(x)))
  208. #define d_sinh(x) (sinh(*(x)))
  209. #define d_sqrt(x) (sqrt(*(x)))
  210. #define d_tan(x) (tan(*(x)))
  211. #define d_tanh(x) (tanh(*(x)))
  212. #define i_abs(x) abs(*(x))
  213. #define i_dnnt(x) ((integer)u_nint(*(x)))
  214. #define i_len(s, n) (n)
  215. #define i_nint(x) ((integer)u_nint(*(x)))
  216. #define i_sign(a,b) ((integer)u_sign((integer)*(a),(integer)*(b)))
  217. #define pow_dd(ap, bp) ( pow(*(ap), *(bp)))
  218. #define pow_si(B,E) spow_ui(*(B),*(E))
  219. #define pow_ri(B,E) spow_ui(*(B),*(E))
  220. #define pow_di(B,E) dpow_ui(*(B),*(E))
  221. #define pow_zi(p, a, b) {pCd(p) = zpow_ui(Cd(a), *(b));}
  222. #define pow_ci(p, a, b) {pCf(p) = cpow_ui(Cf(a), *(b));}
  223. #define pow_zz(R,A,B) {pCd(R) = cpow(Cd(A),*(B));}
  224. #define s_cat(lpp, rpp, rnp, np, llp) { ftnlen i, nc, ll; char *f__rp, *lp; ll = (llp); lp = (lpp); for(i=0; i < (int)*(np); ++i) { nc = ll; if((rnp)[i] < nc) nc = (rnp)[i]; ll -= nc; f__rp = (rpp)[i]; while(--nc >= 0) *lp++ = *(f__rp)++; } while(--ll >= 0) *lp++ = ' '; }
  225. #define s_cmp(a,b,c,d) ((integer)strncmp((a),(b),f2cmin((c),(d))))
  226. #define s_copy(A,B,C,D) { int __i,__m; for (__i=0, __m=f2cmin((C),(D)); __i<__m && (B)[__i] != 0; ++__i) (A)[__i] = (B)[__i]; }
  227. #define sig_die(s, kill) { exit(1); }
  228. #define s_stop(s, n) {exit(0);}
  229. #define z_abs(z) (cabs(Cd(z)))
  230. #define z_exp(R, Z) {pCd(R) = cexp(Cd(Z));}
  231. #define z_sqrt(R, Z) {pCd(R) = csqrt(Cd(Z));}
  232. #define myexit_() break;
  233. #define mycycle_() continue;
  234. #define myceiling_(w) {ceil(w)}
  235. #define myhuge_(w) {HUGE_VAL}
  236. //#define mymaxloc_(w,s,e,n) {if (sizeof(*(w)) == sizeof(double)) dmaxloc_((w),*(s),*(e),n); else dmaxloc_((w),*(s),*(e),n);}
  237. #define mymaxloc_(w,s,e,n) dmaxloc_(w,*(s),*(e),n)
  238. /* procedure parameter types for -A and -C++ */
  239. #define F2C_proc_par_types 1
  240. #ifdef __cplusplus
  241. typedef logical (*L_fp)(...);
  242. #else
  243. typedef logical (*L_fp)();
  244. #endif
  245. #if 0
  246. static float spow_ui(float x, integer n) {
  247. float pow=1.0; unsigned long int u;
  248. if(n != 0) {
  249. if(n < 0) n = -n, x = 1/x;
  250. for(u = n; ; ) {
  251. if(u & 01) pow *= x;
  252. if(u >>= 1) x *= x;
  253. else break;
  254. }
  255. }
  256. return pow;
  257. }
  258. static double dpow_ui(double x, integer n) {
  259. double pow=1.0; unsigned long int u;
  260. if(n != 0) {
  261. if(n < 0) n = -n, x = 1/x;
  262. for(u = n; ; ) {
  263. if(u & 01) pow *= x;
  264. if(u >>= 1) x *= x;
  265. else break;
  266. }
  267. }
  268. return pow;
  269. }
  270. #ifdef _MSC_VER
  271. static _Fcomplex cpow_ui(complex x, integer n) {
  272. complex pow={1.0,0.0}; unsigned long int u;
  273. if(n != 0) {
  274. if(n < 0) n = -n, x.r = 1/x.r, x.i=1/x.i;
  275. for(u = n; ; ) {
  276. if(u & 01) pow.r *= x.r, pow.i *= x.i;
  277. if(u >>= 1) x.r *= x.r, x.i *= x.i;
  278. else break;
  279. }
  280. }
  281. _Fcomplex p={pow.r, pow.i};
  282. return p;
  283. }
  284. #else
  285. static _Complex float cpow_ui(_Complex float x, integer n) {
  286. _Complex float pow=1.0; unsigned long int u;
  287. if(n != 0) {
  288. if(n < 0) n = -n, x = 1/x;
  289. for(u = n; ; ) {
  290. if(u & 01) pow *= x;
  291. if(u >>= 1) x *= x;
  292. else break;
  293. }
  294. }
  295. return pow;
  296. }
  297. #endif
  298. #ifdef _MSC_VER
  299. static _Dcomplex zpow_ui(_Dcomplex x, integer n) {
  300. _Dcomplex pow={1.0,0.0}; unsigned long int u;
  301. if(n != 0) {
  302. if(n < 0) n = -n, x._Val[0] = 1/x._Val[0], x._Val[1] =1/x._Val[1];
  303. for(u = n; ; ) {
  304. if(u & 01) pow._Val[0] *= x._Val[0], pow._Val[1] *= x._Val[1];
  305. if(u >>= 1) x._Val[0] *= x._Val[0], x._Val[1] *= x._Val[1];
  306. else break;
  307. }
  308. }
  309. _Dcomplex p = {pow._Val[0], pow._Val[1]};
  310. return p;
  311. }
  312. #else
  313. static _Complex double zpow_ui(_Complex double x, integer n) {
  314. _Complex double pow=1.0; unsigned long int u;
  315. if(n != 0) {
  316. if(n < 0) n = -n, x = 1/x;
  317. for(u = n; ; ) {
  318. if(u & 01) pow *= x;
  319. if(u >>= 1) x *= x;
  320. else break;
  321. }
  322. }
  323. return pow;
  324. }
  325. #endif
  326. static integer pow_ii(integer x, integer n) {
  327. integer pow; unsigned long int u;
  328. if (n <= 0) {
  329. if (n == 0 || x == 1) pow = 1;
  330. else if (x != -1) pow = x == 0 ? 1/x : 0;
  331. else n = -n;
  332. }
  333. if ((n > 0) || !(n == 0 || x == 1 || x != -1)) {
  334. u = n;
  335. for(pow = 1; ; ) {
  336. if(u & 01) pow *= x;
  337. if(u >>= 1) x *= x;
  338. else break;
  339. }
  340. }
  341. return pow;
  342. }
  343. static integer dmaxloc_(double *w, integer s, integer e, integer *n)
  344. {
  345. double m; integer i, mi;
  346. for(m=w[s-1], mi=s, i=s+1; i<=e; i++)
  347. if (w[i-1]>m) mi=i ,m=w[i-1];
  348. return mi-s+1;
  349. }
  350. static integer smaxloc_(float *w, integer s, integer e, integer *n)
  351. {
  352. float m; integer i, mi;
  353. for(m=w[s-1], mi=s, i=s+1; i<=e; i++)
  354. if (w[i-1]>m) mi=i ,m=w[i-1];
  355. return mi-s+1;
  356. }
  357. #endif
  358. static inline void cdotc_(complex *z, integer *n_, complex *x, integer *incx_, complex *y, integer *incy_) {
  359. integer n = *n_, incx = *incx_, incy = *incy_, i;
  360. #ifdef _MSC_VER
  361. _Fcomplex zdotc = {0.0, 0.0};
  362. if (incx == 1 && incy == 1) {
  363. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  364. zdotc._Val[0] += conjf(Cf(&x[i]))._Val[0] * Cf(&y[i])._Val[0];
  365. zdotc._Val[1] += conjf(Cf(&x[i]))._Val[1] * Cf(&y[i])._Val[1];
  366. }
  367. } else {
  368. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  369. zdotc._Val[0] += conjf(Cf(&x[i*incx]))._Val[0] * Cf(&y[i*incy])._Val[0];
  370. zdotc._Val[1] += conjf(Cf(&x[i*incx]))._Val[1] * Cf(&y[i*incy])._Val[1];
  371. }
  372. }
  373. pCf(z) = zdotc;
  374. }
  375. #else
  376. _Complex float zdotc = 0.0;
  377. if (incx == 1 && incy == 1) {
  378. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  379. zdotc += conjf(Cf(&x[i])) * Cf(&y[i]);
  380. }
  381. } else {
  382. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  383. zdotc += conjf(Cf(&x[i*incx])) * Cf(&y[i*incy]);
  384. }
  385. }
  386. pCf(z) = zdotc;
  387. }
  388. #endif
  389. static inline void zdotc_(doublecomplex *z, integer *n_, doublecomplex *x, integer *incx_, doublecomplex *y, integer *incy_) {
  390. integer n = *n_, incx = *incx_, incy = *incy_, i;
  391. #ifdef _MSC_VER
  392. _Dcomplex zdotc = {0.0, 0.0};
  393. if (incx == 1 && incy == 1) {
  394. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  395. zdotc._Val[0] += conj(Cd(&x[i]))._Val[0] * Cd(&y[i])._Val[0];
  396. zdotc._Val[1] += conj(Cd(&x[i]))._Val[1] * Cd(&y[i])._Val[1];
  397. }
  398. } else {
  399. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  400. zdotc._Val[0] += conj(Cd(&x[i*incx]))._Val[0] * Cd(&y[i*incy])._Val[0];
  401. zdotc._Val[1] += conj(Cd(&x[i*incx]))._Val[1] * Cd(&y[i*incy])._Val[1];
  402. }
  403. }
  404. pCd(z) = zdotc;
  405. }
  406. #else
  407. _Complex double zdotc = 0.0;
  408. if (incx == 1 && incy == 1) {
  409. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  410. zdotc += conj(Cd(&x[i])) * Cd(&y[i]);
  411. }
  412. } else {
  413. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  414. zdotc += conj(Cd(&x[i*incx])) * Cd(&y[i*incy]);
  415. }
  416. }
  417. pCd(z) = zdotc;
  418. }
  419. #endif
  420. static inline void cdotu_(complex *z, integer *n_, complex *x, integer *incx_, complex *y, integer *incy_) {
  421. integer n = *n_, incx = *incx_, incy = *incy_, i;
  422. #ifdef _MSC_VER
  423. _Fcomplex zdotc = {0.0, 0.0};
  424. if (incx == 1 && incy == 1) {
  425. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  426. zdotc._Val[0] += Cf(&x[i])._Val[0] * Cf(&y[i])._Val[0];
  427. zdotc._Val[1] += Cf(&x[i])._Val[1] * Cf(&y[i])._Val[1];
  428. }
  429. } else {
  430. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  431. zdotc._Val[0] += Cf(&x[i*incx])._Val[0] * Cf(&y[i*incy])._Val[0];
  432. zdotc._Val[1] += Cf(&x[i*incx])._Val[1] * Cf(&y[i*incy])._Val[1];
  433. }
  434. }
  435. pCf(z) = zdotc;
  436. }
  437. #else
  438. _Complex float zdotc = 0.0;
  439. if (incx == 1 && incy == 1) {
  440. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  441. zdotc += Cf(&x[i]) * Cf(&y[i]);
  442. }
  443. } else {
  444. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  445. zdotc += Cf(&x[i*incx]) * Cf(&y[i*incy]);
  446. }
  447. }
  448. pCf(z) = zdotc;
  449. }
  450. #endif
  451. static inline void zdotu_(doublecomplex *z, integer *n_, doublecomplex *x, integer *incx_, doublecomplex *y, integer *incy_) {
  452. integer n = *n_, incx = *incx_, incy = *incy_, i;
  453. #ifdef _MSC_VER
  454. _Dcomplex zdotc = {0.0, 0.0};
  455. if (incx == 1 && incy == 1) {
  456. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  457. zdotc._Val[0] += Cd(&x[i])._Val[0] * Cd(&y[i])._Val[0];
  458. zdotc._Val[1] += Cd(&x[i])._Val[1] * Cd(&y[i])._Val[1];
  459. }
  460. } else {
  461. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  462. zdotc._Val[0] += Cd(&x[i*incx])._Val[0] * Cd(&y[i*incy])._Val[0];
  463. zdotc._Val[1] += Cd(&x[i*incx])._Val[1] * Cd(&y[i*incy])._Val[1];
  464. }
  465. }
  466. pCd(z) = zdotc;
  467. }
  468. #else
  469. _Complex double zdotc = 0.0;
  470. if (incx == 1 && incy == 1) {
  471. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  472. zdotc += Cd(&x[i]) * Cd(&y[i]);
  473. }
  474. } else {
  475. for (i=0;i<n;i++) { /* zdotc = zdotc + dconjg(x(i))* y(i) */
  476. zdotc += Cd(&x[i*incx]) * Cd(&y[i*incy]);
  477. }
  478. }
  479. pCd(z) = zdotc;
  480. }
  481. #endif
  482. /* Common Block Declarations */
  483. struct {
  484. integer icase, n, incx, incy, mode;
  485. logical pass;
  486. } combla_;
  487. #define combla_1 combla_
  488. /* Table of constant values */
  489. static integer c__1 = 1;
  490. static integer c__5 = 5;
  491. static real c_b43 = (float)1.;
  492. /* Main program */ int main()
  493. {
  494. /* Initialized data */
  495. static real sfac = (float)9.765625e-4;
  496. /* Local variables */
  497. extern /* Subroutine */ int check1_(), check2_();
  498. static integer ic;
  499. extern /* Subroutine */ int header_();
  500. /* Test program for the COMPLEX Level 1 CBLAS. */
  501. /* Based upon the original CBLAS test routine together with: */
  502. /* F06GAF Example Program Text */
  503. /* .. Parameters .. */
  504. /* .. Scalars in Common .. */
  505. /* .. Local Scalars .. */
  506. /* .. External Subroutines .. */
  507. /* .. Common blocks .. */
  508. /* .. Data statements .. */
  509. /* .. Executable Statements .. */
  510. printf("Complex CBLAS Test Program Results\n");
  511. for (ic = 1; ic <= 10; ++ic) {
  512. combla_1.icase = ic;
  513. header_();
  514. /* Initialize PASS, INCX, INCY, and MODE for a new case. */
  515. /* The value 9999 for INCX, INCY or MODE will appear in the */
  516. /* detailed output, if any, for cases that do not involve */
  517. /* these parameters. */
  518. combla_1.pass = TRUE_;
  519. combla_1.incx = 9999;
  520. combla_1.incy = 9999;
  521. combla_1.mode = 9999;
  522. if (combla_1.icase <= 5) {
  523. check2_(&sfac);
  524. } else if (combla_1.icase >= 6) {
  525. check1_(&sfac);
  526. }
  527. /* -- Print */
  528. if (combla_1.pass) {
  529. printf(" ----- PASS -----\n");
  530. }
  531. /* L20: */
  532. }
  533. exit(0);
  534. } /* MAIN__ */
  535. /* Subroutine */ int header_()
  536. {
  537. /* Initialized data */
  538. static char l[15][13] = {"CBLAS_CDOTC " , "CBLAS_CDOTU " , "CBLAS_CAXPY " ,
  539. "CBLAS_CCOPY " , "CBLAS_CSWAP " , "CBLAS_SCNRM2" , "CBLAS_SCASUM" , "CBLAS_CSCAL " ,
  540. "CBLAS_CSSCAL" , "CBLAS_ICAMAX" };
  541. /* Format strings */
  542. /* Builtin functions */
  543. integer s_wsfe(), do_fio(), e_wsfe();
  544. /* .. Parameters .. */
  545. /* .. Scalars in Common .. */
  546. /* .. Local Arrays .. */
  547. /* .. Common blocks .. */
  548. /* .. Data statements .. */
  549. /* .. Executable Statements .. */
  550. printf("Test of subprogram number %3d %15s\n", combla_1.icase, l[combla_1.icase - 1]);
  551. return 0;
  552. } /* header_ */
  553. /* Subroutine */ int check1_(sfac)
  554. real *sfac;
  555. {
  556. /* Initialized data */
  557. static real strue2[5] = { (float)0.,(float).5,(float).6,(float).7,(float)
  558. .7 };
  559. static real strue4[5] = { (float)0.,(float).7,(float)1.,(float)1.3,(float)
  560. 1.7 };
  561. static complex ctrue5[80] /* was [8][5][2] */ = { {(float).1,(float).1},
  562. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  563. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  564. {(float)1.,(float)2.},{(float)-.16,(float)-.37},{(float)3.,(float)
  565. 4.},{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)
  566. 4.},{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)
  567. 4.},{(float)-.17,(float)-.19},{(float).13,(float)-.39},{(float)5.,
  568. (float)6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float)5.,
  569. (float)6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float)
  570. .11,(float)-.03},{(float)-.17,(float).46},{(float)-.17,(float)
  571. -.19},{(float)7.,(float)8.},{(float)7.,(float)8.},{(float)7.,(
  572. float)8.},{(float)7.,(float)8.},{(float)7.,(float)8.},{(float).19,
  573. (float)-.17},{(float).32,(float).09},{(float).23,(float)-.24},{(
  574. float).18,(float).01},{(float)2.,(float)3.},{(float)2.,(float)3.},
  575. {(float)2.,(float)3.},{(float)2.,(float)3.},{(float).1,(float).1},
  576. {(float)4.,(float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},
  577. {(float)4.,(float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},
  578. {(float)4.,(float)5.},{(float)-.16,(float)-.37},{(float)6.,(float)
  579. 7.},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)6.,(float)
  580. 7.},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)6.,(float)
  581. 7.},{(float)-.17,(float)-.19},{(float)8.,(float)9.},{(float).13,(
  582. float)-.39},{(float)2.,(float)5.},{(float)2.,(float)5.},{(float)
  583. 2.,(float)5.},{(float)2.,(float)5.},{(float)2.,(float)5.},{(float)
  584. .11,(float)-.03},{(float)3.,(float)6.},{(float)-.17,(float).46},{(
  585. float)4.,(float)7.},{(float)-.17,(float)-.19},{(float)7.,(float)
  586. 2.},{(float)7.,(float)2.},{(float)7.,(float)2.},{(float).19,(
  587. float)-.17},{(float)5.,(float)8.},{(float).32,(float).09},{(float)
  588. 6.,(float)9.},{(float).23,(float)-.24},{(float)8.,(float)3.},{(
  589. float).18,(float).01},{(float)9.,(float)4.} };
  590. static complex ctrue6[80] /* was [8][5][2] */ = { {(float).1,(float).1},
  591. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  592. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  593. {(float)1.,(float)2.},{(float).09,(float)-.12},{(float)3.,(float)
  594. 4.},{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)
  595. 4.},{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)
  596. 4.},{(float).03,(float)-.09},{(float).15,(float)-.03},{(float)5.,(
  597. float)6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float)5.,(
  598. float)6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float).03,
  599. (float).03},{(float)-.18,(float).03},{(float).03,(float)-.09},{(
  600. float)7.,(float)8.},{(float)7.,(float)8.},{(float)7.,(float)8.},{(
  601. float)7.,(float)8.},{(float)7.,(float)8.},{(float).09,(float).03},
  602. {(float).03,(float).12},{(float).12,(float).03},{(float).03,(
  603. float).06},{(float)2.,(float)3.},{(float)2.,(float)3.},{(float)2.,
  604. (float)3.},{(float)2.,(float)3.},{(float).1,(float).1},{(float)4.,
  605. (float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},{(float)4.,
  606. (float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},{(float)4.,
  607. (float)5.},{(float).09,(float)-.12},{(float)6.,(float)7.},{(float)
  608. 6.,(float)7.},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)
  609. 6.,(float)7.},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)
  610. .03,(float)-.09},{(float)8.,(float)9.},{(float).15,(float)-.03},{(
  611. float)2.,(float)5.},{(float)2.,(float)5.},{(float)2.,(float)5.},{(
  612. float)2.,(float)5.},{(float)2.,(float)5.},{(float).03,(float).03},
  613. {(float)3.,(float)6.},{(float)-.18,(float).03},{(float)4.,(float)
  614. 7.},{(float).03,(float)-.09},{(float)7.,(float)2.},{(float)7.,(
  615. float)2.},{(float)7.,(float)2.},{(float).09,(float).03},{(float)
  616. 5.,(float)8.},{(float).03,(float).12},{(float)6.,(float)9.},{(
  617. float).12,(float).03},{(float)8.,(float)3.},{(float).03,(float)
  618. .06},{(float)9.,(float)4.} };
  619. static integer itrue3[5] = { 0,1,2,2,2 };
  620. static real sa = (float).3;
  621. static complex ca = {(float).4,(float)-.7};
  622. static complex cv[80] /* was [8][5][2] */ = { {(float).1,(float).1},
  623. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  624. {(float)1.,(float)2.},{(float)1.,(float)2.},{(float)1.,(float)2.},
  625. {(float)1.,(float)2.},{(float).3,(float)-.4},{(float)3.,(float)4.}
  626. ,{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)4.}
  627. ,{(float)3.,(float)4.},{(float)3.,(float)4.},{(float)3.,(float)4.}
  628. ,{(float).1,(float)-.3},{(float).5,(float)-.1},{(float)5.,(float)
  629. 6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float)5.,(float)
  630. 6.},{(float)5.,(float)6.},{(float)5.,(float)6.},{(float).1,(float)
  631. .1},{(float)-.6,(float).1},{(float).1,(float)-.3},{(float)7.,(
  632. float)8.},{(float)7.,(float)8.},{(float)7.,(float)8.},{(float)7.,(
  633. float)8.},{(float)7.,(float)8.},{(float).3,(float).1},{(float).1,(
  634. float).4},{(float).4,(float).1},{(float).1,(float).2},{(float)2.,(
  635. float)3.},{(float)2.,(float)3.},{(float)2.,(float)3.},{(float)2.,(
  636. float)3.},{(float).1,(float).1},{(float)4.,(float)5.},{(float)4.,(
  637. float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},{(float)4.,(
  638. float)5.},{(float)4.,(float)5.},{(float)4.,(float)5.},{(float).3,(
  639. float)-.4},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)6.,
  640. (float)7.},{(float)6.,(float)7.},{(float)6.,(float)7.},{(float)6.,
  641. (float)7.},{(float)6.,(float)7.},{(float).1,(float)-.3},{(float)
  642. 8.,(float)9.},{(float).5,(float)-.1},{(float)2.,(float)5.},{(
  643. float)2.,(float)5.},{(float)2.,(float)5.},{(float)2.,(float)5.},{(
  644. float)2.,(float)5.},{(float).1,(float).1},{(float)3.,(float)6.},{(
  645. float)-.6,(float).1},{(float)4.,(float)7.},{(float).1,(float)-.3},
  646. {(float)7.,(float)2.},{(float)7.,(float)2.},{(float)7.,(float)2.},
  647. {(float).3,(float).1},{(float)5.,(float)8.},{(float).1,(float).4},
  648. {(float)6.,(float)9.},{(float).4,(float).1},{(float)8.,(float)3.},
  649. {(float).1,(float).2},{(float)9.,(float)4.} };
  650. /* System generated locals */
  651. integer i__1, i__2, i__3;
  652. real r__1;
  653. complex q__1;
  654. /* Local variables */
  655. static integer i__;
  656. extern /* Subroutine */ int cscal_(), ctest_();
  657. static complex mwpcs[5], mwpct[5];
  658. extern /* Subroutine */ int itest1_(), stest1_();
  659. static complex cx[8];
  660. extern real scnrm2test_();
  661. static integer np1;
  662. extern integer icamaxtest_();
  663. extern /* Subroutine */ int csscaltest_();
  664. extern real scasumtest_();
  665. static integer len;
  666. /* .. Parameters .. */
  667. /* .. Scalar Arguments .. */
  668. /* .. Scalars in Common .. */
  669. /* .. Local Scalars .. */
  670. /* .. Local Arrays .. */
  671. /* .. External Functions .. */
  672. /* .. External Subroutines .. */
  673. /* .. Intrinsic Functions .. */
  674. /* .. Common blocks .. */
  675. /* .. Data statements .. */
  676. /* .. Executable Statements .. */
  677. for (combla_1.incx = 1; combla_1.incx <= 2; ++combla_1.incx) {
  678. for (np1 = 1; np1 <= 5; ++np1) {
  679. combla_1.n = np1 - 1;
  680. len = f2cmax(combla_1.n,1) << 1;
  681. /* .. Set vector arguments .. */
  682. i__1 = len;
  683. for (i__ = 1; i__ <= i__1; ++i__) {
  684. i__2 = i__ - 1;
  685. i__3 = i__ + (np1 + combla_1.incx * 5 << 3) - 49;
  686. cx[i__2].r = cv[i__3].r, cx[i__2].i = cv[i__3].i;
  687. /* L20: */
  688. }
  689. if (combla_1.icase == 6) {
  690. /* .. SCNRM2TEST .. */
  691. r__1 = scnrm2test_(&combla_1.n, cx, &combla_1.incx);
  692. stest1_(&r__1, &strue2[np1 - 1], &strue2[np1 - 1], sfac);
  693. } else if (combla_1.icase == 7) {
  694. /* .. SCASUMTEST .. */
  695. r__1 = scasumtest_(&combla_1.n, cx, &combla_1.incx);
  696. stest1_(&r__1, &strue4[np1 - 1], &strue4[np1 - 1], sfac);
  697. } else if (combla_1.icase == 8) {
  698. /* .. CSCAL .. */
  699. cscal_(&combla_1.n, &ca, cx, &combla_1.incx);
  700. ctest_(&len, cx, &ctrue5[(np1 + combla_1.incx * 5 << 3) - 48],
  701. &ctrue5[(np1 + combla_1.incx * 5 << 3) - 48], sfac);
  702. } else if (combla_1.icase == 9) {
  703. /* .. CSSCALTEST .. */
  704. csscaltest_(&combla_1.n, &sa, cx, &combla_1.incx);
  705. ctest_(&len, cx, &ctrue6[(np1 + combla_1.incx * 5 << 3) - 48],
  706. &ctrue6[(np1 + combla_1.incx * 5 << 3) - 48], sfac);
  707. } else if (combla_1.icase == 10) {
  708. /* .. ICAMAXTEST .. */
  709. i__1 = icamaxtest_(&combla_1.n, cx, &combla_1.incx);
  710. itest1_(&i__1, &itrue3[np1 - 1]);
  711. } else {
  712. fprintf(stderr,"Shouldn't be here in CHECK1\n");
  713. exit(0);
  714. }
  715. /* L40: */
  716. }
  717. /* L60: */
  718. }
  719. combla_1.incx = 1;
  720. if (combla_1.icase == 8) {
  721. /* CSCAL */
  722. /* Add a test for alpha equal to zero. */
  723. ca.r = (float)0., ca.i = (float)0.;
  724. for (i__ = 1; i__ <= 5; ++i__) {
  725. i__1 = i__ - 1;
  726. mwpct[i__1].r = (float)0., mwpct[i__1].i = (float)0.;
  727. i__1 = i__ - 1;
  728. mwpcs[i__1].r = (float)1., mwpcs[i__1].i = (float)1.;
  729. /* L80: */
  730. }
  731. cscal_(&c__5, &ca, cx, &combla_1.incx);
  732. ctest_(&c__5, cx, mwpct, mwpcs, sfac);
  733. } else if (combla_1.icase == 9) {
  734. /* CSSCALTEST */
  735. /* Add a test for alpha equal to zero. */
  736. sa = (float)0.;
  737. for (i__ = 1; i__ <= 5; ++i__) {
  738. i__1 = i__ - 1;
  739. mwpct[i__1].r = (float)0., mwpct[i__1].i = (float)0.;
  740. i__1 = i__ - 1;
  741. mwpcs[i__1].r = (float)1., mwpcs[i__1].i = (float)1.;
  742. /* L100: */
  743. }
  744. csscaltest_(&c__5, &sa, cx, &combla_1.incx);
  745. ctest_(&c__5, cx, mwpct, mwpcs, sfac);
  746. /* Add a test for alpha equal to one. */
  747. sa = (float)1.;
  748. for (i__ = 1; i__ <= 5; ++i__) {
  749. i__1 = i__ - 1;
  750. i__2 = i__ - 1;
  751. mwpct[i__1].r = cx[i__2].r, mwpct[i__1].i = cx[i__2].i;
  752. i__1 = i__ - 1;
  753. i__2 = i__ - 1;
  754. mwpcs[i__1].r = cx[i__2].r, mwpcs[i__1].i = cx[i__2].i;
  755. /* L120: */
  756. }
  757. csscaltest_(&c__5, &sa, cx, &combla_1.incx);
  758. ctest_(&c__5, cx, mwpct, mwpcs, sfac);
  759. /* Add a test for alpha equal to minus one. */
  760. sa = (float)-1.;
  761. for (i__ = 1; i__ <= 5; ++i__) {
  762. i__1 = i__ - 1;
  763. i__2 = i__ - 1;
  764. q__1.r = -cx[i__2].r, q__1.i = -cx[i__2].i;
  765. mwpct[i__1].r = q__1.r, mwpct[i__1].i = q__1.i;
  766. i__1 = i__ - 1;
  767. i__2 = i__ - 1;
  768. q__1.r = -cx[i__2].r, q__1.i = -cx[i__2].i;
  769. mwpcs[i__1].r = q__1.r, mwpcs[i__1].i = q__1.i;
  770. /* L140: */
  771. }
  772. csscaltest_(&c__5, &sa, cx, &combla_1.incx);
  773. ctest_(&c__5, cx, mwpct, mwpcs, sfac);
  774. }
  775. return 0;
  776. } /* check1_ */
  777. /* Subroutine */ int check2_(sfac)
  778. real *sfac;
  779. {
  780. /* Initialized data */
  781. static complex ca = {(float).4,(float)-.7};
  782. static integer incxs[4] = { 1,2,-2,-1 };
  783. static integer incys[4] = { 1,-2,1,-2 };
  784. static integer lens[8] /* was [4][2] */ = { 1,1,2,4,1,1,3,7 };
  785. static integer ns[4] = { 0,1,2,4 };
  786. static complex cx1[7] = { {(float).7,(float)-.8},{(float)-.4,(float)-.7},{
  787. (float)-.1,(float)-.9},{(float).2,(float)-.8},{(float)-.9,(float)
  788. -.4},{(float).1,(float).4},{(float)-.6,(float).6} };
  789. static complex cy1[7] = { {(float).6,(float)-.6},{(float)-.9,(float).5},{(
  790. float).7,(float)-.6},{(float).1,(float)-.5},{(float)-.1,(float)
  791. -.2},{(float)-.5,(float)-.3},{(float).8,(float)-.7} };
  792. static complex ct8[112] /* was [7][4][4] */ = { {(float).6,(float)-.6}
  793. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  794. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  795. ,{(float).32,(float)-1.41},{(float)0.,(float)0.},{(float)0.,(
  796. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  797. float)0.},{(float)0.,(float)0.},{(float).32,(float)-1.41},{(float)
  798. -1.55,(float).5},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  799. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  800. float).32,(float)-1.41},{(float)-1.55,(float).5},{(float).03,(
  801. float)-.89},{(float)-.38,(float)-.96},{(float)0.,(float)0.},{(
  802. float)0.,(float)0.},{(float)0.,(float)0.},{(float).6,(float)-.6},{
  803. (float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{
  804. (float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{
  805. (float).32,(float)-1.41},{(float)0.,(float)0.},{(float)0.,(float)
  806. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)
  807. 0.},{(float)0.,(float)0.},{(float)-.07,(float)-.89},{(float)-.9,(
  808. float).5},{(float).42,(float)-1.41},{(float)0.,(float)0.},{(float)
  809. 0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)
  810. .78,(float).06},{(float)-.9,(float).5},{(float).06,(float)-.13},{(
  811. float).1,(float)-.5},{(float)-.77,(float)-.49},{(float)-.5,(float)
  812. -.3},{(float).52,(float)-1.51},{(float).6,(float)-.6},{(float)0.,(
  813. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  814. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float).32,
  815. (float)-1.41},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)
  816. 0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)
  817. 0.,(float)0.},{(float)-.07,(float)-.89},{(float)-1.18,(float)-.31}
  818. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  819. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float).78,(float)
  820. .06},{(float)-1.54,(float).97},{(float).03,(float)-.89},{(float)
  821. -.18,(float)-1.31},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  822. float)0.,(float)0.},{(float).6,(float)-.6},{(float)0.,(float)0.},{
  823. (float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{
  824. (float)0.,(float)0.},{(float)0.,(float)0.},{(float).32,(float)
  825. -1.41},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  826. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  827. float)0.},{(float).32,(float)-1.41},{(float)-.9,(float).5},{(
  828. float).05,(float)-.6},{(float)0.,(float)0.},{(float)0.,(float)0.},
  829. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float).32,(float)
  830. -1.41},{(float)-.9,(float).5},{(float).05,(float)-.6},{(float).1,(
  831. float)-.5},{(float)-.77,(float)-.49},{(float)-.5,(float)-.3},{(
  832. float).32,(float)-1.16} };
  833. static complex ct7[16] /* was [4][4] */ = { {(float)0.,(float)0.},{(
  834. float)-.06,(float)-.9},{(float).65,(float)-.47},{(float)-.34,(
  835. float)-1.22},{(float)0.,(float)0.},{(float)-.06,(float)-.9},{(
  836. float)-.59,(float)-1.46},{(float)-1.04,(float)-.04},{(float)0.,(
  837. float)0.},{(float)-.06,(float)-.9},{(float)-.83,(float).59},{(
  838. float).07,(float)-.37},{(float)0.,(float)0.},{(float)-.06,(float)
  839. -.9},{(float)-.76,(float)-1.15},{(float)-1.33,(float)-1.82} };
  840. static complex ct6[16] /* was [4][4] */ = { {(float)0.,(float)0.},{(
  841. float).9,(float).06},{(float).91,(float)-.77},{(float)1.8,(float)
  842. -.1},{(float)0.,(float)0.},{(float).9,(float).06},{(float)1.45,(
  843. float).74},{(float).2,(float).9},{(float)0.,(float)0.},{(float).9,
  844. (float).06},{(float)-.55,(float).23},{(float).83,(float)-.39},{(
  845. float)0.,(float)0.},{(float).9,(float).06},{(float)1.04,(float)
  846. .79},{(float)1.95,(float)1.22} };
  847. static complex ct10x[112] /* was [7][4][4] */ = { {(float).7,(float)-.8}
  848. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  849. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  850. ,{(float).6,(float)-.6},{(float)0.,(float)0.},{(float)0.,(float)
  851. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)
  852. 0.},{(float)0.,(float)0.},{(float).6,(float)-.6},{(float)-.9,(
  853. float).5},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  854. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float).6,(
  855. float)-.6},{(float)-.9,(float).5},{(float).7,(float)-.6},{(float)
  856. .1,(float)-.5},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  857. float)0.,(float)0.},{(float).7,(float)-.8},{(float)0.,(float)0.},{
  858. (float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{
  859. (float)0.,(float)0.},{(float)0.,(float)0.},{(float).6,(float)-.6},
  860. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},
  861. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},
  862. {(float).7,(float)-.6},{(float)-.4,(float)-.7},{(float).6,(float)
  863. -.6},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  864. float)0.},{(float)0.,(float)0.},{(float).8,(float)-.7},{(float)
  865. -.4,(float)-.7},{(float)-.1,(float)-.2},{(float).2,(float)-.8},{(
  866. float).7,(float)-.6},{(float).1,(float).4},{(float).6,(float)-.6},
  867. {(float).7,(float)-.8},{(float)0.,(float)0.},{(float)0.,(float)0.}
  868. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  869. ,{(float)0.,(float)0.},{(float).6,(float)-.6},{(float)0.,(float)
  870. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)
  871. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)-.9,(
  872. float).5},{(float)-.4,(float)-.7},{(float).6,(float)-.6},{(float)
  873. 0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)
  874. 0.,(float)0.},{(float).1,(float)-.5},{(float)-.4,(float)-.7},{(
  875. float).7,(float)-.6},{(float).2,(float)-.8},{(float)-.9,(float).5}
  876. ,{(float).1,(float).4},{(float).6,(float)-.6},{(float).7,(float)
  877. -.8},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  878. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  879. float)0.},{(float).6,(float)-.6},{(float)0.,(float)0.},{(float)0.,
  880. (float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,
  881. (float)0.},{(float)0.,(float)0.},{(float).6,(float)-.6},{(float)
  882. .7,(float)-.6},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  883. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  884. float).6,(float)-.6},{(float).7,(float)-.6},{(float)-.1,(float)
  885. -.2},{(float).8,(float)-.7},{(float)0.,(float)0.},{(float)0.,(
  886. float)0.},{(float)0.,(float)0.} };
  887. static complex ct10y[112] /* was [7][4][4] */ = { {(float).6,(float)-.6}
  888. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  889. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  890. ,{(float).7,(float)-.8},{(float)0.,(float)0.},{(float)0.,(float)
  891. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)
  892. 0.},{(float)0.,(float)0.},{(float).7,(float)-.8},{(float)-.4,(
  893. float)-.7},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,
  894. (float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float).7,
  895. (float)-.8},{(float)-.4,(float)-.7},{(float)-.1,(float)-.9},{(
  896. float).2,(float)-.8},{(float)0.,(float)0.},{(float)0.,(float)0.},{
  897. (float)0.,(float)0.},{(float).6,(float)-.6},{(float)0.,(float)0.},
  898. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},
  899. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float).7,(float)-.8}
  900. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  901. ,{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.}
  902. ,{(float)-.1,(float)-.9},{(float)-.9,(float).5},{(float).7,(float)
  903. -.8},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  904. float)0.},{(float)0.,(float)0.},{(float)-.6,(float).6},{(float)
  905. -.9,(float).5},{(float)-.9,(float)-.4},{(float).1,(float)-.5},{(
  906. float)-.1,(float)-.9},{(float)-.5,(float)-.3},{(float).7,(float)
  907. -.8},{(float).6,(float)-.6},{(float)0.,(float)0.},{(float)0.,(
  908. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  909. float)0.},{(float)0.,(float)0.},{(float).7,(float)-.8},{(float)0.,
  910. (float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,
  911. (float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)
  912. -.1,(float)-.9},{(float).7,(float)-.8},{(float)0.,(float)0.},{(
  913. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  914. float)0.,(float)0.},{(float)-.6,(float).6},{(float)-.9,(float)-.4}
  915. ,{(float)-.1,(float)-.9},{(float).7,(float)-.8},{(float)0.,(float)
  916. 0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float).6,(float)
  917. -.6},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  918. float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(
  919. float)0.},{(float).7,(float)-.8},{(float)0.,(float)0.},{(float)0.,
  920. (float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,
  921. (float)0.},{(float)0.,(float)0.},{(float).7,(float)-.8},{(float)
  922. -.9,(float).5},{(float)-.4,(float)-.7},{(float)0.,(float)0.},{(
  923. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  924. float).7,(float)-.8},{(float)-.9,(float).5},{(float)-.4,(float)
  925. -.7},{(float).1,(float)-.5},{(float)-.1,(float)-.9},{(float)-.5,(
  926. float)-.3},{(float).2,(float)-.8} };
  927. static complex csize1[4] = { {(float)0.,(float)0.},{(float).9,(float).9},{
  928. (float)1.63,(float)1.73},{(float)2.9,(float)2.78} };
  929. static complex csize3[14] = { {(float)0.,(float)0.},{(float)0.,(float)0.},
  930. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},
  931. {(float)0.,(float)0.},{(float)0.,(float)0.},{(float)1.17,(float)
  932. 1.17},{(float)1.17,(float)1.17},{(float)1.17,(float)1.17},{(float)
  933. 1.17,(float)1.17},{(float)1.17,(float)1.17},{(float)1.17,(float)
  934. 1.17},{(float)1.17,(float)1.17} };
  935. static complex csize2[14] /* was [7][2] */ = { {(float)0.,(float)0.},{(
  936. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  937. float)0.,(float)0.},{(float)0.,(float)0.},{(float)0.,(float)0.},{(
  938. float)1.54,(float)1.54},{(float)1.54,(float)1.54},{(float)1.54,(
  939. float)1.54},{(float)1.54,(float)1.54},{(float)1.54,(float)1.54},{(
  940. float)1.54,(float)1.54},{(float)1.54,(float)1.54} };
  941. /* System generated locals */
  942. integer i__1, i__2;
  943. /* Local variables */
  944. static complex cdot[1];
  945. static integer lenx, leny, i__;
  946. static complex ctemp;
  947. extern /* Subroutine */ int ctest_();
  948. static integer ksize;
  949. extern /* Subroutine */ int cdotctest_(), ccopytest_(), cdotutest_(),
  950. cswaptest_(), caxpytest_();
  951. static integer ki, kn;
  952. static complex cx[7], cy[7];
  953. static integer mx, my;
  954. /* .. Parameters .. */
  955. /* .. Scalar Arguments .. */
  956. /* .. Scalars in Common .. */
  957. /* .. Local Scalars .. */
  958. /* .. Local Arrays .. */
  959. /* .. External Functions .. */
  960. /* .. External Subroutines .. */
  961. /* .. Intrinsic Functions .. */
  962. /* .. Common blocks .. */
  963. /* .. Data statements .. */
  964. /* .. Executable Statements .. */
  965. for (ki = 1; ki <= 4; ++ki) {
  966. combla_1.incx = incxs[ki - 1];
  967. combla_1.incy = incys[ki - 1];
  968. mx = abs(combla_1.incx);
  969. my = abs(combla_1.incy);
  970. for (kn = 1; kn <= 4; ++kn) {
  971. combla_1.n = ns[kn - 1];
  972. ksize = f2cmin(2,kn);
  973. lenx = lens[kn + (mx << 2) - 5];
  974. leny = lens[kn + (my << 2) - 5];
  975. /* .. initialize all argument arrays .. */
  976. for (i__ = 1; i__ <= 7; ++i__) {
  977. i__1 = i__ - 1;
  978. i__2 = i__ - 1;
  979. cx[i__1].r = cx1[i__2].r, cx[i__1].i = cx1[i__2].i;
  980. i__1 = i__ - 1;
  981. i__2 = i__ - 1;
  982. cy[i__1].r = cy1[i__2].r, cy[i__1].i = cy1[i__2].i;
  983. /* L20: */
  984. }
  985. if (combla_1.icase == 1) {
  986. /* .. CDOTCTEST .. */
  987. cdotctest_(&combla_1.n, cx, &combla_1.incx, cy, &
  988. combla_1.incy, &ctemp);
  989. cdot[0].r = ctemp.r, cdot[0].i = ctemp.i;
  990. ctest_(&c__1, cdot, &ct6[kn + (ki << 2) - 5], &csize1[kn - 1],
  991. sfac);
  992. } else if (combla_1.icase == 2) {
  993. /* .. CDOTUTEST .. */
  994. cdotutest_(&combla_1.n, cx, &combla_1.incx, cy, &
  995. combla_1.incy, &ctemp);
  996. cdot[0].r = ctemp.r, cdot[0].i = ctemp.i;
  997. ctest_(&c__1, cdot, &ct7[kn + (ki << 2) - 5], &csize1[kn - 1],
  998. sfac);
  999. } else if (combla_1.icase == 3) {
  1000. /* .. CAXPYTEST .. */
  1001. caxpytest_(&combla_1.n, &ca, cx, &combla_1.incx, cy, &
  1002. combla_1.incy);
  1003. ctest_(&leny, cy, &ct8[(kn + (ki << 2)) * 7 - 35], &csize2[
  1004. ksize * 7 - 7], sfac);
  1005. } else if (combla_1.icase == 4) {
  1006. /* .. CCOPYTEST .. */
  1007. ccopytest_(&combla_1.n, cx, &combla_1.incx, cy, &
  1008. combla_1.incy);
  1009. ctest_(&leny, cy, &ct10y[(kn + (ki << 2)) * 7 - 35], csize3, &
  1010. c_b43);
  1011. } else if (combla_1.icase == 5) {
  1012. /* .. CSWAPTEST .. */
  1013. cswaptest_(&combla_1.n, cx, &combla_1.incx, cy, &
  1014. combla_1.incy);
  1015. ctest_(&lenx, cx, &ct10x[(kn + (ki << 2)) * 7 - 35], csize3, &
  1016. c_b43);
  1017. ctest_(&leny, cy, &ct10y[(kn + (ki << 2)) * 7 - 35], csize3, &
  1018. c_b43);
  1019. } else {
  1020. fprintf(stderr,"Shouldn't be here in CHECK2\n");
  1021. exit(0);
  1022. }
  1023. /* L40: */
  1024. }
  1025. /* L60: */
  1026. }
  1027. return 0;
  1028. } /* check2_ */
  1029. /* Subroutine */ int stest_(len, scomp, strue, ssize, sfac)
  1030. integer *len;
  1031. real *scomp, *strue, *ssize, *sfac;
  1032. {
  1033. /* System generated locals */
  1034. integer i__1;
  1035. real r__1, r__2, r__3, r__4, r__5;
  1036. /* Local variables */
  1037. static integer i__;
  1038. extern doublereal sdiff_();
  1039. static real sd;
  1040. /* ********************************* STEST ************************** */
  1041. /* THIS SUBR COMPARES ARRAYS SCOMP() AND STRUE() OF LENGTH LEN TO */
  1042. /* SEE IF THE TERM BY TERM DIFFERENCES, MULTIPLIED BY SFAC, ARE */
  1043. /* NEGLIGIBLE. */
  1044. /* C. L. LAWSON, JPL, 1974 DEC 10 */
  1045. /* .. Parameters .. */
  1046. /* .. Scalar Arguments .. */
  1047. /* .. Array Arguments .. */
  1048. /* .. Scalars in Common .. */
  1049. /* .. Local Scalars .. */
  1050. /* .. External Functions .. */
  1051. /* .. Intrinsic Functions .. */
  1052. /* .. Common blocks .. */
  1053. /* .. Executable Statements .. */
  1054. /* Parameter adjustments */
  1055. --ssize;
  1056. --strue;
  1057. --scomp;
  1058. /* Function Body */
  1059. i__1 = *len;
  1060. for (i__ = 1; i__ <= i__1; ++i__) {
  1061. sd = scomp[i__] - strue[i__];
  1062. r__4 = (r__1 = ssize[i__], dabs(r__1)) + (r__2 = *sfac * sd, dabs(
  1063. r__2));
  1064. r__5 = (r__3 = ssize[i__], dabs(r__3));
  1065. if (sdiff_(&r__4, &r__5) == (float)0.) {
  1066. goto L40;
  1067. }
  1068. /* HERE SCOMP(I) IS NOT CLOSE TO STRUE(I). */
  1069. if (! combla_1.pass) {
  1070. goto L20;
  1071. }
  1072. /* PRINT FAIL MESSAGE AND HEADER. */
  1073. combla_1.pass = FALSE_;
  1074. printf(" FAIL\n");
  1075. printf("CASE N INCX INCY MODE I COMP(I) TRUE(I) DIFFERENCE SIZE(I)\n");
  1076. L20:
  1077. printf("%4d %3d %5d %5d %5d %3d %36.8e %36.8e %12.4e %12.4e\n",combla_1.icase, combla_1.n, combla_1.incx, combla_1.incy,
  1078. combla_1.mode, i__, scomp[i__], strue[i__], sd, ssize[i__]);
  1079. L40:
  1080. ;
  1081. }
  1082. return 0;
  1083. } /* stest_ */
  1084. /* Subroutine */ int stest1_(scomp1, strue1, ssize, sfac)
  1085. real *scomp1, *strue1, *ssize, *sfac;
  1086. {
  1087. static real scomp[1], strue[1];
  1088. extern /* Subroutine */ int stest_();
  1089. /* ************************* STEST1 ***************************** */
  1090. /* THIS IS AN INTERFACE SUBROUTINE TO ACCOMMODATE THE FORTRAN */
  1091. /* REQUIREMENT THAT WHEN A DUMMY ARGUMENT IS AN ARRAY, THE */
  1092. /* ACTUAL ARGUMENT MUST ALSO BE AN ARRAY OR AN ARRAY ELEMENT. */
  1093. /* C.L. LAWSON, JPL, 1978 DEC 6 */
  1094. /* .. Scalar Arguments .. */
  1095. /* .. Array Arguments .. */
  1096. /* .. Local Arrays .. */
  1097. /* .. External Subroutines .. */
  1098. /* .. Executable Statements .. */
  1099. /* Parameter adjustments */
  1100. --ssize;
  1101. /* Function Body */
  1102. scomp[0] = *scomp1;
  1103. strue[0] = *strue1;
  1104. stest_(&c__1, scomp, strue, &ssize[1], sfac);
  1105. return 0;
  1106. } /* stest1_ */
  1107. doublereal sdiff_(sa, sb)
  1108. real *sa, *sb;
  1109. {
  1110. /* System generated locals */
  1111. real ret_val;
  1112. /* ********************************* SDIFF ************************** */
  1113. /* COMPUTES DIFFERENCE OF TWO NUMBERS. C. L. LAWSON, JPL 1974 FEB 15 */
  1114. /* .. Scalar Arguments .. */
  1115. /* .. Executable Statements .. */
  1116. ret_val = *sa - *sb;
  1117. return ret_val;
  1118. } /* sdiff_ */
  1119. /* Subroutine */ int ctest_(len, ccomp, ctrue, csize, sfac)
  1120. integer *len;
  1121. complex *ccomp, *ctrue, *csize;
  1122. real *sfac;
  1123. {
  1124. /* System generated locals */
  1125. integer i__1, i__2;
  1126. /* Builtin functions */
  1127. // double r_imag();
  1128. /* Local variables */
  1129. static integer i__;
  1130. static real scomp[20], ssize[20], strue[20];
  1131. extern /* Subroutine */ int stest_();
  1132. /* **************************** CTEST ***************************** */
  1133. /* C.L. LAWSON, JPL, 1978 DEC 6 */
  1134. /* .. Scalar Arguments .. */
  1135. /* .. Array Arguments .. */
  1136. /* .. Local Scalars .. */
  1137. /* .. Local Arrays .. */
  1138. /* .. External Subroutines .. */
  1139. /* .. Intrinsic Functions .. */
  1140. /* .. Executable Statements .. */
  1141. /* Parameter adjustments */
  1142. --csize;
  1143. --ctrue;
  1144. --ccomp;
  1145. /* Function Body */
  1146. i__1 = *len;
  1147. for (i__ = 1; i__ <= i__1; ++i__) {
  1148. i__2 = i__;
  1149. scomp[(i__ << 1) - 2] = ccomp[i__2].r;
  1150. scomp[(i__ << 1) - 1] = r_imag(&ccomp[i__]);
  1151. i__2 = i__;
  1152. strue[(i__ << 1) - 2] = ctrue[i__2].r;
  1153. strue[(i__ << 1) - 1] = r_imag(&ctrue[i__]);
  1154. i__2 = i__;
  1155. ssize[(i__ << 1) - 2] = csize[i__2].r;
  1156. ssize[(i__ << 1) - 1] = r_imag(&csize[i__]);
  1157. /* L20: */
  1158. }
  1159. i__1 = *len << 1;
  1160. stest_(&i__1, scomp, strue, ssize, sfac);
  1161. return 0;
  1162. } /* ctest_ */
  1163. /* Subroutine */ int itest1_(icomp, itrue)
  1164. integer *icomp, *itrue;
  1165. {
  1166. /* Local variables */
  1167. static integer id;
  1168. /* ********************************* ITEST1 ************************* */
  1169. /* THIS SUBROUTINE COMPARES THE VARIABLES ICOMP AND ITRUE FOR */
  1170. /* EQUALITY. */
  1171. /* C. L. LAWSON, JPL, 1974 DEC 10 */
  1172. /* .. Parameters .. */
  1173. /* .. Scalar Arguments .. */
  1174. /* .. Scalars in Common .. */
  1175. /* .. Local Scalars .. */
  1176. /* .. Common blocks .. */
  1177. /* .. Executable Statements .. */
  1178. if (*icomp == *itrue) {
  1179. goto L40;
  1180. }
  1181. /* HERE ICOMP IS NOT EQUAL TO ITRUE. */
  1182. if (! combla_1.pass) {
  1183. goto L20;
  1184. }
  1185. /* PRINT FAIL MESSAGE AND HEADER. */
  1186. combla_1.pass = FALSE_;
  1187. printf(" FAIL\n");
  1188. printf(" CASE N INCX INCY MODE COMP TRUE DIFFERENCE\n");
  1189. L20:
  1190. id = *icomp - *itrue;
  1191. printf("%4d %3d %5d %5d %5d %36d %36d %12d\n",combla_1.icase, combla_1.n, combla_1.incx, combla_1.incy,
  1192. combla_1.mode, *icomp, *itrue, id);
  1193. L40:
  1194. return 0;
  1195. } /* itest1_ */